Skip to content
Snippets Groups Projects
Commit 42fedde2 authored by Adam Friedrich Schrey's avatar Adam Friedrich Schrey
Browse files

Merge branch 'develop' into 'master'

merge develop into master

See merge request !1
parents 814e4267 68879533
No related branches found
No related tags found
1 merge request!1merge develop into master
Showing
with 3772 additions and 1340 deletions
%% Cell type:markdown id: tags:
# Introduction to JupyterHub
# Systemtheorie 2 und JupyterHub
## 1st Step
Login in the page jupyter.rwth-aachen.de/
Please use your RWTH account
## Wie greife ich auf JupyterHub zu?
Loggen Sie sich unter <a href="https://jupyter.rwth-aachen.de/hub/spawn?profile=sys2"> https://jupyter.rwth-aachen.de/hub/spawn?profile=sys2 </a> mit Ihrem RWTH-Account ein.
## 2.step
Select the system theory 2 profile
## Wie nutze ich JupyterHub?
Auf der Weboberfläche sollten auf der linken Seite verschiedene Ordner angezeigt werden. Unter dem Pfad "lecture-tutorials/notebooks" sind Dateien (Notebooks) zum Thema "Systemtheorie 2".
## 3. Step
Direct into the folder in which we provide you with existing notebooks for your lecture.
In diesen Notebooks können Beispiele in Form von Code-Blöcken ausgeführt werden. Meist geben diese einen statischen Graphen aus. Jedoch gibt es auch oft interaktive Elemente, bei denen Werte über einen Schieberegler oder Knopf verändert werden können.
## 4. Step
The notebooks for each lecture are mainly based on the octave kernel. A switch to scilab or python kernels could be needed but is indicated within the respective noteooks. You can switch the kernel on the top right corner of your browser window.
Die Code-Blöcke können im Notebook geändert werden und diese Änderung kann unter dem Unterpunkt "Git" auch wieder rückgängig gemacht werden.
"V01_zeitdiskrete_systeme.ipynb" bis "V10_kalmanfilter_demonstration.ipynb" sind Notebooks zu den Themen der Vorlesung. In "Sandkasten.ipynb" kann man beliebige Funktionen untersuchen, eigene Simulationen erstellen und Verschiedenes ausprobieren.
### Ausführen der Code-Blöcke in den Folien:
<img src="introduction_pictures/run.png"/>
### Beispiel eines interaktiven Elementes:
<img src="introduction_pictures/interact.png"/>
### Änderungen rückgängig machen:
<img src="introduction_pictures/git.png"/>
<img src="introduction_pictures/discard_changes.png"/>
## Python und Matlab
Ältere Notebooks, welche Matlab-Code verwenden, können unter dem Pfad "lecture-tutorials/notebooks/alt" gefunden werden.Aber auch andere Notebooks können Matlab-Code beinhalten.Wenn man zwischen verschiedenen Notebooks mit Python-Code und Matlab-Code wechselt, muss man den Kernel womöglichmanuell ändern (Octave Kernel für Matlab):
### Ändern des Kernels:
<img src="introduction_pictures/switch_kernel.png"/>
## Fehlerbehebung
Sollte es zu Fehlern kommen oder das System in eine Endlosschleife geraten, so kann man unter dem Menü-Punkt "Kernel" >> "Restart Kernel" die aktuellen Berechnungen und gespeicherten Variablen zurücksetzen.
Sollten die Notebooks oder zugehöriger Python-Code fehlerhaft sein, können Sie dies gerne melden.
......
File moved
introduction_pictures/discard_changes.png

72.8 KiB

introduction_pictures/git.png

73 KiB

introduction_pictures/interact.png

120 KiB

introduction_pictures/run.png

123 KiB

introduction_pictures/switch_kernel.png

40.4 KiB

This diff is collapsed.
# -*- coding: utf-8 -*-
import numpy
class LinProb(object):
def __init__(self, nodeNum):
self.G = numpy.zeros((nodeNum,nodeNum), dtype=numpy.complex)
self.b = numpy.zeros((nodeNum,1), dtype=numpy.complex)
def getSize(self):
return self.G.shape[0]
def resetbVector(self):
self.b = numpy.zeros((self.b.shape[0],1), dtype=numpy.complex)
def solve(self):
result = numpy.dot(numpy.linalg.inv(self.G), self.b)
return result
def show(self):
print(self.G)
print(self.b)
class Circuit(object):
def __init__(self, nnode):
self.n = nnode
self.NetList = list()
def addComponent(self,h):
print("adding component" + str(h))
self.nelements = len(self.NetList)+1
self.NetList.append(h)
def initialize(self,dt):
self.Pn = LinProb(self.n)
self.Pn.show()
self.nelements = len(self.NetList)
print("number of components " + str(self.nelements))
for i in range(0,self.nelements):
print("initializing component " + str(i+1))
self.NetList[i].initialize(self.Pn,dt)
return self.Pn
def initializeDp(self,dt,om):
self.Pn = LinProb(self.n)
self.nelements = len(self.NetList)
for i in range(0,self.nelements):
self.NetList[i].initialize(self.Pn,dt,om)
def solveLinTrans(self,tini,tfin,dt):
print("show problem")
self.Pn.show()
t = tini
k = 0
while t < tfin:
#print("solving for t = " + str(t))
self.Pn.resetbVector()
for i in range(0,self.nelements):
self.Pn = self.NetList[i].step(self.Pn,dt,t)
if t == 0:
self.vout = self.Pn.solve()
else:
self.vout = numpy.append(self.vout, self.Pn.solve(), 1)
for i in range(0,self.nelements):
self.NetList[i].postStep(self.vout[:,k],dt)
t = t+dt
k = k+1
return self.vout
\ No newline at end of file
import numpy
from Block import Block
inv = numpy.linalg.inv
dot = numpy.dot
#########################################################################
#### Abstract class for discrete time block in a block diagram ########
#########################################################################
class DTBlock(Block):
def __init__(self):
self.ninput = 0
self.noutput = 0
self.nstate = 0
self.nparam = 0
self.inpos = numpy.zeros((self.ninput),dtype=numpy.int)
self.outpos = numpy.zeros((self.noutput),dtype=numpy.int)
self.x = numpy.zeros((self.nstate))
self.xo = numpy.zeros((self.nstate))
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
self.Ts = 0.001
self.lastt = -10
def Step(self,t,dt):
if t-self.lastt >= self.Ts:
self.updatediscrete()
self.lastt = t
def updatediscrete(self):
self.y[0] = 0
def Reset(self):
self.lastt = -10;
#########################################################################
################################ FIR ####################################
#########################################################################
class FIR(DTBlock):
def __init__(self,inp,out,a,ts):
super(FIR, self).__init__()
self.ninput = 1
self.inpos = numpy.zeros((self.ninput),dtype=numpy.int)
self.inpos[0] = inp
self.noutput = 1
self.outpos = numpy.zeros((self.noutput),dtype=numpy.int)
self.outpos[0] = out
self.Ts = ts
self.pt = len(a)
self.a = numpy.zeros((self.pt))
self.a = a
self.buf = numpy.zeros((self.pt))
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
self.y[0] = 0
def updatediscrete(self):
for i in range (self.pt-1,0,-1):
self.buf[i] = self.buf[i-1]
self.buf[0] = self.u[0]
self.y[0] = 0
for i in range (0,self.pt):
self.y[0] = self.y[0]+self.a[i]*self.buf[i]
def Reset(self):
self.buf = zeros(self.pt,1)
self.lastt = -10
#########################################################################
#########################Discrete Integral###############################
#########################################################################
class DTIntegral(DTBlock):
def __init__(self,inp,out,ts,xo):
super(DTIntegral, self).__init__()
self.ninput = 1
self.inpos = numpy.zeros((self.ninput),dtype=numpy.int)
self.inpos[0] = inp
self.noutput = 1
self.outpos = numpy.zeros((self.noutput),dtype=numpy.int)
self.outpos[0] = out
self.Ts = ts
self.xo = xo
self.x = xo
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
self.y[0] = xo
def updatediscrete(self):
self.y[0] = self.x
self.x = self.x + self.u[0]*self.Ts
def Reset(self):
self.x = self.xo
#########################################################################
######################### ZOH #################################
#########################################################################
class ZOH(DTBlock):
def __init__(self,inp,out,ts):
super(ZOH, self).__init__()
self.ninput = 1
self.inpos = numpy.zeros((self.ninput),dtype=numpy.int)
self.inpos[0] = inp
self.noutput = 1
self.outpos = numpy.zeros((self.noutput),dtype=numpy.int)
self.outpos[0] = out
self.Ts = ts
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
self.y[0] = 0
def updatediscrete(self):
self.y[0] = self.u[0]
#########################################################################
####################### Delay #################################
#########################################################################
class DTDelay(DTBlock):
def __init__(self,inp,out,init,ts):
super(DTDelay, self).__init__()
self.ninput = 1
self.inpos = numpy.zeros((self.ninput),dtype=numpy.int)
self.inpos[0] = inp
self.noutput = 1
self.outpos = numpy.zeros((self.noutput),dtype=numpy.int)
self.outpos[0] = out
self.Ts = ts
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
self.y[0] = 0
self.oldv = init
self.init = init
def updatediscrete(self):
self.y[0] = self.oldv
self.oldv = self.u[0]
def ChangeParameters(self,value):
self.Ts = value
def Reset(self):
self.oldv = self.init
#########################################################################
################################ IIR ####################################
#########################################################################
class IIR(DTBlock):
def __init__(self,inp,out,a,b,ts):
super(IIR, self).__init__()
self.ninput = 1
self.inpos = numpy.zeros((self.ninput),dtype=numpy.int)
self.inpos[0] = inp
self.noutput = 1
self.outpos = numpy.zeros((self.noutput),dtype=numpy.int)
self.outpos[0] = out
self.Ts = ts
self.pn = len(a)
self.a = numpy.zeros((self.pn))
self.a = a
self.pd = len(b)
self.b = numpy.zeros((self.pd))
self.b = b
self.bufn = numpy.zeros((self.pn))
self.bufd = numpy.zeros((self.pd))
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
self.y[0] = 0
def updatediscrete(self):
for i in range (self.pn-1,0,-1):
self.bufn[i] = self.bufn[i-1]
self.bufn[0] = self.u[0]
self.y[0] = 0
for i in range (0,self.pn):
self.y[0] = self.y[0]+self.a[i]*self.bufn[i]
for i in range (self.pd-1,0,-1):
self.bufd[i] = self.bufd[i-1]
self.bufd[0] = self.y[0]
for i in range (0,self.pd):
self.y[0] = self.y[0]-self.b[i]*self.bufd[i]
def Reset(self):
self.bufn = zeros(self.pn,1)
self.bufd = zeros(self.pd,1)
self.lastt = -10
#########################################################################
############### Discrete Time Transfer Function ########################
#########################################################################
class DTTransferFunction(DTBlock):
def __init__(self,inp,out,num,den,Ts):
super(DTTransferFunction, self).__init__()
self.Ts = Ts
self.ninput = 1
self.noutput = 1
self.inpos = numpy.zeros((self.ninput),dtype=numpy.int)
self.outpos = numpy.zeros((self.noutput),dtype=numpy.int)
self.inpos[0] = inp
self.outpos[0] = out
self.nstate = len(den)-1
denord = len(den)
numord = len(num)
self.A = numpy.zeros((self.nstate,self.nstate))
self.B = numpy.zeros((self.nstate))
self.C = numpy.zeros((self.nstate))
self.D = 0
scale = den[0];
for i in range(0,denord):
den[i]= den[i]/scale
for i in range(0,numord):
num[i]= num[i]/scale
if numord==denord:
self.D = num[0]
for i in range(0,numord-1):
self.C[i]=num[numord-i-1]-num[0]*den[numord-i-1]
else:
for i in range(0,numord):
self.C[i]=num[numord-i-1]
for i in range(0,self.nstate):
self.A[self.nstate-1,i] = -den[denord-i-1]
for i in range(0,self.nstate-1):
self.A[i,i+1]=1
self.B[self.nstate-1] = 1
self.x = numpy.zeros((self.nstate))
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
for i in range(0,self.nstate):
self.x[i] = 0
def updatediscrete(self):
self.y[0] = self.D*self.u[0]
for i in range (0,self.nstate):
self.y[0] = self.y[0]+self.C[i]*self.x[i]
xn = numpy.zeros((self.nstate))
for i in range (0,self.nstate):
xn[i] = self.B[i]*self.u[0]
for j in range (0,self.nstate):
xn[i] = xn[i] + self.A[i,j]*self.x[j]
for i in range (0,self.nstate):
self.x[i] = xn[i]
def Reset(self):
for i in range(0,self.nstate):
self.x[i] = 0
#########################################################################
########################### DT State Space #############################
#########################################################################
class DTStateSpace(DTBlock):
def __init__(self,inp,out,A,B,C,D,xo,Ts):
super(DTStateSpace, self).__init__()
self.Ts = Ts
self.ninput = len(inp)
self.noutput = len(out)
self.nstate = len(A)
self.A = A
self.B = B
self.C = C
self.D = D
self.inpos = inp
self.outpos = out
self.x = numpy.zeros((self.nstate))
self.xo = numpy.zeros((self.nstate))
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
for i in range(0,self.nstate):
self.x[i] = xo[i]
self.xo[i] = xo[i]
def updatediscrete(self):
for i in range(0,self.noutput):
self.y[i] = 0
for j in range(0,self.nstate):
if self.C.ndim == 1:
self.y[i] = self.y[i] + self.C[j]*self.x[j]
elif self.C.ndim == 2:
#prevent [ValueError: setting an array element with a sequence.] for 2 dimensional self.C:
self.y[i] = self.y[i] + self.C[i][j]*self.x[j]
else:
print("WARNING: Can not handle self.C with dimension "+str(self.C.ndim))
for j in range(0,self.ninput):
if self.D.ndim == 1:
self.y[i] = self.y[i] + self.D[j]*self.u[j]
elif self.D.ndim == 2:
#prevent [ValueError: setting an array element with a sequence.] for 2 dimensional self.C:
self.y[i] = self.y[i] + self.D[i][j]*self.u[j]
else:
print("WARNING: Can not handle self.D with dimension "+str(self.C.ndim))
self.x = numpy.matmul(self.B,self.u.T) + numpy.matmul(self.A,self.x)
def Reset(self):
for i in range(0,self.nstate):
self.x[i] = 0
#########################################################################
########################### Kalman-Filter ##############################
#########################################################################
class KalmanFilter(DTBlock):
def __init__(self,inp,out,nin,nmeas,F,G,C,Q,R,xo,Pk,Ts):
self.F = F
self.G = G
self.C = C
self.Q = Q
self.R = R
self.x = xo
self.Pk = Pk
self.Ts = Ts
self.nin = nin
self.nmeas = nmeas
self.ninput = len(inp)
self.noutput = len(out)
self.inpos = inp;
self.outpos = out;
self.nstate = len(F);
self.u = numpy.zeros((self.ninput))
self.y = numpy.zeros((self.noutput))
self.lastt = -10
def updatediscrete(self):
uu = numpy.zeros((self.nstate))
m = numpy.zeros((self.nstate))
if (self.nin != 0):
for i in range(0,self.nin):
uu[i] = self.u[i];
for i in range(0,self.nmeas):
m[i] = self.u[self.nin+i]
# Zustandsvorhersage:
# x = F*x + G*u
self.x = dot(self.F,self.x)+dot(self.G,uu)
# Kovarianz der Zustandsvorhersage:
# Pk = F*P*F.T + Q
self.Pk = dot(self.F,dot(self.Pk,self.F.T)) + self.Q
# Kovarianz der Innovation:
# S = C*Pk*C.T + R
S = dot(self.C, dot(self.Pk, self.C.T))+ self.R
# Kalman Gain:
# Kk = Pk*self.C.T*inv(S)
Kk = dot(dot(self.Pk, self.C.T), inv(S))
# Messvorhersage:
# z = C*x
z = dot(self.C,self.x)
# Innovation:
v = m - z
# Kovarianzkorrektur:
# Pk = (I-(Kk*C))*Pk
self.Pk= dot((numpy.eye(self.nstate)-dot(Kk,self.C)),self.Pk)
# Zustandskorrektur:
self.x = self.x+dot(Kk,v)
# Widergabe des Zustandsvektors:
self.y = self.x
def Reset(self):
self.x = xo
self.Pk = Pk
# -*- coding: utf-8 -*-
import numpy
class RCComponent(object):
def __init__(self, p, n):
self.Pos = p
self.Neg = n
self.linear = 1
self.hyin = 0
self.hyout = 0
def applyMatrixStamp(self, P,dt):
return P
def initialize(self, P,dt):
return P
def step(self,P,dt,t):
return P
def postStep(self,vsol,dt):
return vsol
class Resistance(RCComponent):
def __init__(self,p,n,par):
super(Resistance, self).__init__(p,n)
self.r = par
if par > 0:
self.Gc = 1/par
else:
print("Error: Resistance value not correct")
def initialize(self,P,dt):
print("initializing resistance")
return Resistance.applyMatrixStamp(self,P,dt)
def applyMatrixStamp(self,P,dt):
if self.Pos > -1:
P.G[self.Pos,self.Pos] = P.G[self.Pos,self.Pos] + self.Gc
if self.Neg > -1:
P.G[self.Neg,self.Neg] = P.G[self.Neg,self.Neg] + self.Gc
if (self.Pos > -1) and (self.Neg > -1):
P.G[self.Pos,self.Neg] = P.G[self.Pos,self.Neg] - self.Gc
P.G[self.Neg,self.Pos] = P.G[self.Neg,self.Pos] - self.Gc
P.show()
return P
class Capacitor(RCComponent):
def __init__(self,p,n,par):
super(Capacitor, self).__init__(p,n)
self.C = par
self.Gl = 0
def initialize(self,P,dt):
print("initializing capacitor")
self.Bc = 0
self.Vc = 0
self.Ic = 0
return Capacitor.applyMatrixStamp(self,P,dt)
def applyMatrixStamp(self,P,dt):
self.Gc = 2*self.C/dt
if self.Pos > -1:
P.G[self.Pos,self.Pos] = P.G[self.Pos,self.Pos] + self.Gc
if self.Neg > -1:
P.G[self.Neg,self.Neg] = P.G[self.Neg,self.Neg] + self.Gc
if (self.Pos > -1) and (self.Neg > -1):
P.G[self.Pos,self.Neg] = P.G[self.Pos,self.Neg] - self.Gc
P.G[self.Neg,self.Pos] = P.G[self.Neg,self.Pos] - self.Gc
P.show()
return P
def step(self,P,dt,t):
self.Bc = self.Gc*self.Vc+self.Ic
if self.Pos > -1:
P.b[self.Pos] = P.b[self.Pos] + self.Bc
if self.Neg > -1:
P.b[self.Neg] = P.b[self.Neg] - self.Bc
return P
def postStep(self,vsol,dt):
if self.Pos > -1:
if self.Neg > -1:
self.Vc = vsol[self.Pos] - vsol[self.Neg]
else:
self.Vc = vsol[self.Pos]
else:
self.Vc = -vsol[self.Neg];
self.Ic = self.Gc*self.Vc-self.Bc;
class Inductor(RCComponent):
def __init__(self,p,n,par):
super(Inductor, self).__init__(p,n)
self.L = par
def initialize(self,P,dt):
print("initializing inductor")
self.Bl = 0
self.Vl = 0
self.Il = 0
return Inductor.applyMatrixStamp(self,P,dt)
def applyMatrixStamp(self,P,dt):
self.Gl = dt/(2*self.L)
if self.Pos > -1:
P.G[self.Pos,self.Pos] = P.G[self.Pos,self.Pos] + self.Gl
if self.Neg > -1:
P.G[self.Neg,self.Neg] = P.G[self.Neg,self.Neg] + self.Gl
if (self.Pos > -1) and (self.Neg > -1):
P.G[self.Pos,self.Neg] = P.G[self.Pos,self.Neg] - self.Gl
P.G[self.Neg,self.Pos] = P.G[self.Neg,self.Pos] - self.Gl
P.show()
return P
def step(self,P,dt,t):
self.Bl = self.Gl*self.Vl+self.Il
if self.Pos > -1:
P.b[self.Pos] = P.b[self.Pos] - self.Bl
if self.Neg > -1:
P.b[self.Neg] = P.b[self.Neg] + self.Bl
return P
def postStep(self,vsol,dt):
if self.Pos > -1:
if self.Neg > -1:
self.Vl = vsol[self.Pos] - vsol[self.Neg]
else:
self.Vl = vsol[self.Pos]
else:
self.Vl = -vsol[self.Neg];
self.Il = self.Gl*self.Vl+self.Bl;
class IdealACVoltageSource(RCComponent):
def __init__(self,p,n,V,omega,phase):
super(IdealACVoltageSource, self).__init__(p,n)
self.Vs = V
self.om = omega
self.ph = phase
def initialize(self,P,dt):
print("initializing ideal ac voltage source")
return IdealACVoltageSource.applyMatrixStamp(self,P,dt)
def applyMatrixStamp(self,P,dt):
self.ntot = P.getSize()
self.possource = self.ntot
print("position source at " + str(self.possource))
P.G = numpy.append(P.G, numpy.zeros((self.ntot,1), dtype=numpy.complex), 1)
P.G = numpy.append(P.G, numpy.zeros((1,self.ntot+1), dtype=numpy.complex), 0)
P.b = numpy.append(P.b, numpy.zeros((1,1), dtype=numpy.complex), 0)
if self.Pos > -1:
P.G[self.possource,self.Pos] = 1
P.G[self.Pos,self.possource] = 1
if self.Neg > -1:
P.G[self.possource,self.Neg] = -1
P.G[self.Neg,self.possource] = -1
P.show()
return P
def step(self,P,dt,t):
P.b[self.possource,0] = self.Vs * numpy.sin(self.om*t+self.ph)
return P
\ No newline at end of file
# -*- coding: utf-8 -*-
import numpy as np
class Schema(object):
def __init__(self, ti,tf,dt,numflows):
self.tini = ti
self.tfin = tf
self.nmodels = 0
self.nflows = numflows
self.dt = dt
self.flows = np.zeros(numflows)
self.BlockList = list()
def AddComponent(self,h):
self.nmodels = len(self.BlockList)+1
self.BlockList.append(h)
def AddListComponents(self,l):
p = len(l)
for i in range(0,p):
self.BlockList.append(l[i])
self.nmodels = self.nmodels+1
def Init(self):
self.t = self.tini
for i in range(0,self.nflows):
self.flows[i] = 0
def Reset(self):
for i in range(0,self.nmodels):
self.BlockList[i].Reset()
def Run(self,listout):
self.Init()
qmax = listout.size
npoints = (self.tfin-self.tini)/self.dt
npt = int(npoints)
out = np.zeros((qmax+1,npt-1))
p=0
while self.Step() and p<npt-1:
out[0,p] = self.GetTime()
for q in range(0,qmax):
out[q+1, p] = self.GetFlow(listout[q])
p=p+1
return out
def ChangeTimeStep(self,dt):
self.dt = dt
def ChangeParameters(self,pos,val):
self.ModelList[pos].ChangeParameters(val)
def Step(self):
for i in range(0,self.nmodels):
self.BlockList[i].GetInput(self)
self.BlockList[i].Step(self.t,self.dt)
self.BlockList[i].WriteOutput(self)
self.t = self.t+self.dt
if self.t < self.tfin:
flag = 1
else:
flag = 0
return flag
def GetFlow(self,n):
val = self.flows[n]
return val
def GetTime(self):
val = self.t
return val
\ No newline at end of file
%% Cell type:markdown id:sacred-child tags:
# <span style='color:OrangeRed'>Sandkasten</span>
%% Cell type:code id:adaptive-optimum tags:
``` python
from systheo2functions import *
%matplotlib inline
```
%% Cell type:code id:outdoor-confidence tags:
``` python
# Untersuche Eigenschaften einer beliebigen Regelstrecke:
Gs = 80/((1+2*s)**7) #Gebe hier eine beliebige Funktion ein; s**N steht für s^N
print("Gs:"+str(Gs))
plt_bode(Gs)
margin(Gs)
#plt_nyquist(G,x-Skalierung,y-Skalierung,x-Anfang,x-Ende):
plt_nyquist(Gs,0.5,0.4,-20,0) #Benutze plt_nyquist(Gs,1,1,0,0) um das ganze Diagramm zu sehen
```
%% Cell type:code id:fresh-reviewer tags:
``` python
# Erstellen einer einfachen Simulation mit einem Eingangs- und einem Ausgangs-Signal:
tini = 0 # Start time
tfinal = 1.5 # End time
dt = 0.001 # Time Step
nflows = 3 # Number of data flows in the schematic
Ts = 0.1 # Sampling time for discrete time
sc = Schema(tini,tfinal,dt,nflows) # Instance of the simulation schematic
c1 = SinusoidalSignalSource(1,0,1,2*pi,0)#SinusoidalSignalSource(out,startv,Am,om,phi)
c2 = TransferFunction(1,2,[1],[1,1]);#TransferFunction(inp,out,num,den)
sc.AddListComponents(np.array([c1,c2]));
#Run the schematic and plot:
out = sc.Run(np.array([1, 2]))
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
fig.set_dpi(120)
ax.plot(out[0,:],out[1,:],out[0,:],out[2,:])
ax.grid()
ax.legend(['Eingangssignal (SinusoidalSignalSource)','Ausgangssignal durch TransferFunction'])
plt.show()
```
%% Cell type:code id:hungry-newton tags:
``` python
```
%% Cell type:code id:elementary-lambda tags:
``` python
```
%% Cell type:code id:moving-knock tags:
``` python
```
%% Cell type:markdown id:sapphire-member tags:
# Liste nützlicher Funktionen für Matrizen und mathematische Funktionen:
print_det(M)
print_rank(M)
print_eig(M)
matmul_loop(mats)
plt_bode(G)
plt_nyquist(G,scale_x,scale_y,x0,x1)
margin(G)
# Liste nützlicher Klassen für Simulationen:
### ⚠️ Achtung: Die zur Verfügung stehenden Klassen wurden nicht ausgiebig genug getestet um fehlerfreie Simulationen zu garantieren, außerdem müssen "inp" und "out" bei manchen Klassen als list-Objekt und bei anderen als integer übergeben werden.
StepSource(out,startv,endv,ts)
SinusoidalSignalSource(out,startv,Am,om,phi)
SquareSignal(out,hi,lo,f,duty)
Constant(out,value)
Division(in1,in2,out)
Product(in1,in2,out)
Sum(in1,in2,out,sign1,sign2)
Saturation(inp,out,minout,maxout)
Gain(inp,out,k)
PI(inp,out,Kp,Ki,ini)
PID(inp,out,Kp,Ki,Kd,ini)
Integrator(inp,out,ini)
Not(inp,out)
Noise(inp,out,sigma)
StateSpace(inp,out,A,B,C,D,xo)
TransferFunction(inp,out,num,den)
### -------
FIR(inp,out,a,ts)
DTIntegral(inp,out,ts,xo)
ZOH(inp,out,ts)
DTDelay(inp,out,init,ts)
IIR(inp,out,a,b,ts)
DTTransferFunction(inp,out,num,den,Ts)
DTStateSpace(inp,out,A,B,C,D,xo,Ts)
KalmanFilter(inp,out,N,M,F,G,C,Q,R,xo,Pk,Ts)
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
%% Cell type:markdown id: tags:
# <span style='color:OrangeRed'>V7 - Steuerbarkeit und Beobachtbarkeit</span>
%% Cell type:code id: tags:
``` python
from systheo2functions import *
%matplotlib inline
```
%% Cell type:markdown id: tags:
## <span style='color:Gray'>Beispiel #1 </span>
%% Cell type:markdown id: tags:
<div style="font-family: 'times'; font-size: 13pt; text-align: justify">
Es wird das folgende Eingrößensystem betrachtet:
<br>$\dot{x}(t)$= $Ax(t)$ + $bu(t)$
<br>$y(t)$=$c^Tx(t)$
<br><br> $x$ = $\left[ \begin{array}{rrrr}
x_1 \\
x_2 \\
x_3 \\
\end{array}\right] $
<br><br> $A$ = $\left[ \begin{array}{rrrr}
-1 & -2 & -2 \\
0 & -1 & 1 \\
1 & 0 & -1 \\
\end{array}\right] $
%% Cell type:code id: tags:
``` octave
A = [-1 -1 -2; 0 -1 1; 1 0 -1];
``` python
A = np.array([[-1,-1,-2],[0,-1,1],[1,0,-1]]);
```
%% Cell type:markdown id: tags:
<div style="font-family: 'times'; font-size: 13pt; text-align: justify">
$b$ = $\left[ \begin{array}{rrrr}
2 \\
0 \\
1 \\
\end{array}\right] $
%% Cell type:code id: tags:
``` octave
b = [2; 0; 1];
``` python
b = np.array([[2],[0],[1]]);
```
%% Cell type:markdown id: tags:
<div style="font-family: 'times'; font-size: 13pt; text-align: justify">
$c^T$ = $\left[ \begin{array}{rrrr}
1 & 1 & 0 \\
\end{array}\right] $
%% Cell type:code id: tags:
``` octave
c = [1 1 0];
``` python
c = np.array([1,1,0]);
```
%% Cell type:markdown id: tags:
<div style="font-family: 'times'; font-size: 13pt; text-align: justify">
<b>Zielsetzung</b>: Überprüfung der Steuerbarkeit und Beobachtbarkeit.
<br>Zuerst betrachten wir die <span style='color:red'>Steuerbarkeit</span>:
<br>Ein System 3. Ordnung ist vollständig zustandssteuerbar, wenn gilt:
$Rang$[$B\,|\,AB\,|\,A^2B$]=$Rang$($S_C$) = $n$ = $3$
<br><b>Mit</b>:
<br><br> $A$ = $\left[ \begin{array}{rrrr}
-1 & -2 & -2 \\
0 & -1 & 1 \\
1 & 0 & -1 \\
\end{array}\right] $
<br><br> $B$ = $b$ = $\left[ \begin{array}{rrrr}
2 \\
0 \\
1 \\
\end{array}\right] $
<br><b>Wird</b>:
<br><br> $Ab$ = $\left[ \begin{array}{rrrr}
-4 \\
1 \\
1 \\
\end{array}\right] $
<br><br> $A^2b$ = $\left[ \begin{array}{rrrr}
0 \\
0 \\
-5 \\
\end{array}\right] $
<br>Damit erhalten wir die Steuerbarkeitsmatrix:
<br>$Rang$[$b\,|\,Ab\,|\,A^2b$] = $\left[ \begin{array}{rrrr}
2 & -4 & 0 \\
0 & 1 & 0 \\
1 & 1 & -5 \\
\end{array}\right] $
%% Cell type:code id: tags:
``` octave
S_C = [b A*b A*A*b]
``` python
column1 = b
column2 = np.matmul(A,b)
column3 = np.matmul(np.matmul(A,A),b)
S_C = np.c_[column1,column2,column3]
print("S_C:\n"+str(S_C))
```
%% Output
S_C =
2 -4 1
0 1 0
1 1 -5
S_C:
[[ 2 -4 1]
[ 0 1 0]
[ 1 1 -5]]
%% Cell type:code id: tags:
``` octave
det(S_C)
``` python
print_det(S_C)
```
%% Output
ans = -11
Determinante: -11.0
%% Cell type:code id: tags:
``` octave
rank(S_C)
``` python
print_rank(S_C)
```
%% Output
ans = 3
Rang: 3
%% Cell type:markdown id: tags:
<div style="font-family: 'times'; font-size: 13pt; text-align: justify">
<code>det(S_C)</code> $\ne$ 0 und <code>rank(S_C)</code> = $n$ = $3$
<br>Damit ist das gegebene System <span style='color:red'>vollständig steuerbar</span>!
%% Cell type:markdown id: tags:
<div style="font-family: 'times'; font-size: 13pt; text-align: justify">
<br>Jetzt betrachten wir die <span style='color:red'>Beobachtbarkeit</span>:
<br>Ein Eingrößensystem 3. Ordnung ist genau dann beobachtbar, wenn: $Rang$$\left[ \begin{array}{rrrr}
c^T \\
c^TA\\
c^TA^2\\
\end{array}\right] $ =$Rang$($S_O$) = $n$ = $3$
<br>Im vorliegenden Fall erhält man mit:
$c^T$ = [$1$ $1$ $0$], $c^TA$ = [$-1$ $-3$ $-1$], $c^TA^2$ = [$0$ $5$ $0$]
<br> Damit erhalten wir die Beobachtbarkeitsmatrix:
<br><br> $S_O$ = $\left[ \begin{array}{rrrr}
1 & 1 & 0 \\
-1 & -3 & -1 \\
0 & 5 & 0 \\
\end{array}\right] $
%% Cell type:code id: tags:
``` octave
S_O = [c; c*A; c*A*A]
``` python
S_O = np.c_[c, np.matmul(c,A), np.matmul(np.matmul(c,A),A)]
print("S_C:\n"+str(S_O))
```
%% Output
S_O =
1 1 0
-1 -2 -1
0 3 1
S_C:
[[ 1 -1 0]
[ 1 -2 3]
[ 0 -1 1]]
%% Cell type:code id: tags:
``` octave
det(S_O)
``` python
print_det(S_O)
```
%% Output
ans = 2
Determinante: 2.0
%% Cell type:code id: tags:
``` octave
rank(S_O)
``` python
print_rank(S_O)
```
%% Output
ans = 3
Rang: 3
%% Cell type:markdown id: tags:
<div style="font-family: 'times'; font-size: 13pt; text-align: justify">
<code>det(S_O)</code> $\ne$ 0 und <code>rank(S_O)</code> = $n$ = $3$
<br>Damit ist das gegebene System <span style='color:red'>vollständig beobachtbar</span>!
......
This diff is collapsed.
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment