kramann.info
© Guido Kramann

Realisierung eines Zustandsreglers mit Beobachter

pol.zip - Scilab-Skripte zu den nachfolgenden Beispielen.

Nachtrag: Musterlösung zu vorangegangener Übung Aufgabe 1:

Das folgende Scilab-Skript stellt eine mögliche Lösung zu Aufgabe 1 des vorangegangenen Kapitels dar:

clear;
clc;
clearglobal;
D=0.0018;
J=0.00045;
k=0.015;
R=1.0;
L=0.003;
//r1 und r2 so wählen, das eine doppelte reelle Polstelle bei -10 entsteht:
doppelpolstelle = -10
r2 = (2*(-doppelpolstelle)-(D/J))*L-R
r1 = J*L*(1/k)*(-(D/J)*(R+r2)/L + 0.25*((R+r2)/L+D/J)^2)-k
//Stationäre Werte für i und omega bestimmen für u=5Volt:
u = 5;
A = [-D/J,k/J;-k/L,-R/L];
b = [0;u/L]
y_stat = inv(A)*(-b);
w_stat = y_stat(1,1)
i_stat = y_stat(2,1)
//Überprüfen der Eigenwerte im *-System:
Astern = [-D/J,k/J;-(k+r1)/L,-(R+r2)/L];
spec(Astern)
//Rückführung auf tatsächliche Systemmatrix und "Rest" und
//Simulation mit und ohne Regler:
//ungeregelt:
function f = rechteSeite(t,y)
    uu=5;
    AA = [-D/J,k/J;-k/L,-R/L];
    BB=[0;1/L]
    f=AA*y+BB*uu;
endfunction
//mit Zustandsregler:
function f = rechteSeite2(t,y)
    u=-r1*(y(1,1)-w_stat)-r2*(y(2,1)-i_stat)+k*w_stat+R*i_stat;
    B=[0;1/L]
    f=A*y+B*u;
endfunction
//oder alternativ:
function f = rechteSeite3(t,y)
    RR=[r1,r2];
    BB=[0;1/L];
    CC=[0;1];
    pp=(((k+r1)*w_stat+(R+r2)*i_stat)/L);
    f=(A-BB*RR)*y + CC*pp;
endfunction
t = 0:0.01:1;
y0 = [0,0]';
t0 = 0;
y  = ode(y0,t0,t,rechteSeite);
y2  = ode(y0,t0,t,rechteSeite3);
plot(t,y(1,:)','blu',t,y2(1,:)','gre');
a = gca();
a.x_label.text = 't/s';
a.y_label.text = 'y';
a.title.text = 'Vergleich Sprungantwort (u=5Volt) mit Zustandsregler';
hl=legend(['Sprungantwort';'Zustandsregler'],a=1);

Code 0-1: Bestimmung der Regelmatrix R über die Vorgabe der Eigenwerte des geregelten Systems.

Bild 0-1: Plot: Sprungantwort und Zustandsregler aus obigem Skript.

Nachtrag: Musterlösung zu vorangegangener Übung Aufgabe 2:

Vorgabe der Eigenwerte und Bestimmung der passensten Regelmatrix mit der Scilab-Funktion ppol:

D=0.0018;        
J=0.00045;        
k=0.015;
R=1.0;
L=0.003;
    A = [
            [-D/J,k/J];
            [-k/L,-R/L]
        ];
    B = [
            [0];
            [1/L]
        ]
//RR=ppol(A,B,[-10+10*%i,-10-10*%i])
RR=ppol(A,B,[-10+2*%i,-10-2*%i])
//RR=ppol(A,B,[-10,-10])
//Ergebnis überprüfen:
spec(A-B*RR)
//ungeregelt:
function f = rechteSeite(t,y)
    u=5;
    B=[0;1/L]
    f=A*y+B*u;
endfunction
//geregelt:
function f = rechteSeite3(t,y)
    BB=[0;1/L];
    CC=[0;1];
    pp=(((k+r1)*w_stat+(R+r2)*i_stat)/L);
    f=(A-BB*RR)*y + CC*pp;
endfunction
t = 0:0.01:1;
y0 = [0,0]';
t0 = 0;
y  = ode(y0,t0,t,rechteSeite);
y2  = ode(y0,t0,t,rechteSeite3);
plot(t,y(1,:)','blu',t,y2(1,:)','gre');
a = gca();
a.x_label.text = 't/s';
a.y_label.text = 'y';
a.title.text = 'Vergleich Sprungantwort (u=5Volt) mit Zustandsregler';
hl=legend(['Sprungantwort';'Zustandsregler'],a=1);

Code 0-2: Vorgabe der Eigenwerte und Bestimmung der passensten Regelmatrix mit der Scilab-Funktion ppol.

Bild 0-2: Plot: Sprungantwort und Zustandsregler aus obigem Skript über ppol bestimmt.

Nachtrag: Umsetzung eines Riccati-optimierten Reglers mit Scilab

clear;
clc;
clearglobal;
D=0.0018;
J=0.00045;
k=0.015;
R=1.0;
L=0.003;
//r1 und r2 so wählen, das eine doppelte reelle Polstelle bei -10 entsteht:
doppelpolstelle = -20
r2 = (2*(-doppelpolstelle)-(D/J))*L-R
r1 = J*L*(1/k)*(-(D/J)*(R+r2)/L + 0.25*((R+r2)/L+D/J)^2)-k
//Stationäre Werte für i und omega bestimmen für u=5Volt:
u = 5;
A = [-D/J,k/J;-k/L,-R/L];
b = [0;u/L]
y_stat = inv(A)*(-b);
w_stat = y_stat(1,1)
i_stat = y_stat(2,1)
//Systemmatrix:
AA = [[-D/J,k/J];[-k/L,-R/L]];
//Gewichtungsmatrix für Zustandsgrößen:
QQ = [[1,0];[0,1]];
//Ausgangsmatrix:
CC = [[1,0];[0,1]];
//Durchgangsmatrix:
//DD = [[1,0];[0,1]];
DD = [0;0];
//Eingangsmatrix / Rückführungsmatrix
BB = [0;1/L];
//[X,RCOND,FERR] = riccsl(AA,CC,DD,'cont','schr')
GG = syslin('c',AA,BB,CC,DD);
[XX,FF] = gcare(GG);
XX
RR = BB'*XX
r1_ricc=RR(1,1);
r2_ricc=RR(1,2);
//ungeregelt:
function f = rechteSeite(t,y)
    u=5;
    B=[0;1/L]
    f=A*y+B*u;
endfunction
//-20 als doppelte Polstellen:
function f = rechteSeite3(t,y)
    u = -r1*(y(1,1)-w_stat)-r2*(y(2,1)-i_stat)+k*w_stat+R*i_stat;
    f(1,1)=(-D/J)*y(1,1)+(k/J)*y(2,1);
    f(1,2)=(-k/L)*y(1,1)+(-R/L)*y(2,1)+(1/L)*u;
endfunction
//Riccati-Lösung:
function f = rechteSeite4(t,y)
    u = -r1_ricc*(y(1,1)-w_stat)-r2_ricc*(y(2,1)-i_stat)+k*w_stat+R*i_stat;
    f(1,1)=(-D/J)*y(1,1)+(k/J)*y(2,1);
    f(1,2)=(-k/L)*y(1,1)+(-R/L)*y(2,1)+(1/L)*u;
endfunction
t = 0:0.01:1;
y0 = [0,0]';
yb0 = [0,0,0,0]';
ybb0 = [0,0,20,20]';
t0   = 0;
y    = ode(y0,t0,t,rechteSeite);
y3   = ode(y0,t0,t,rechteSeite3);
y4   = ode(y0,t0,t,rechteSeite4);
plot(t,y(1,:)','blu',t,y3(1,:)','gre',t,y4(1,:)','red');
a = gca();
a.x_label.text = 't/s';
a.y_label.text = 'y';
a.title.text = 'Vergleich Sprungantwort (u=5Volt) mit Riccati-Zustandsregler';
hl=legend(['Sprungantwort';'reelle Pole';'Riccati optimiert'],a=1);

Code 0-3: Optimierung der Eigenwerte und Bestimmung der passensten Regelmatrix mit der Scilab-Funktion gcare (Riccati).

Bild 0-3: Plot: Sprungantwort und Zustandsregler aus obigem Skript über Riccati-Optimierung bestimmt.

Ergänzung: Standarddarstellung linearer Übertragungsglieder im Zeitbereich

Die folgende Darstellung ist die Standarddarstellung linearer Übertragungsglieder im Zeitbereich.
Die Bedeutung der einzelnen Matrizen ist die folgende:
A: Systemmatrix
B: Eingangsmatrix
C: Ausgangsmatrix
D: Durchgangsmatrix

Bild 0-4: Standarddarstellung linearer Übertragungsglieder im Zeitbereich als DGLS erster Ordnung.

Häufig ist bei den hier betrachteten Systemen die Ausgangsmatrix C die Einheitsmatrix und die Durchgangsmatrix die Nullmatrix.

Zustandsregler mit Beobachter

Beim Einsatz von Zustandseglern bei realen Regelstrecken tritt häufig das Problem auf, dass sich nicht alle der zurückzuführenden Zustandsgrößen auch an der Regelstrecke messen lassen.
Eine Lösung für dieses Problem stellt dar, parallel zur Regelung eine Simulation der geregelten Strecke mitlaufen zu lassen und von dort die fehlenden Zustandsgrößen abzugreifen.
Im (natürlich unrealistischen) Idealfall könnte man auch sämtliche Zustandsgrößen simulieren und diese auf die Simulation und das reale System multipliziert mit der Regelmatrix R negativ zurückführen:

Bild 0-5: Idealfall Zustandsregler mit Beobachter-Modell: der simulierte Zustandsvektor wird auf das reale System und das simulierte zurückgeführt.

Ein solches Konzept kann aus folgenden Gründen nicht funktionieren:

Das Modell bildet nicht exakt das reale System ab. Beide Zustandsverläufe driften nach einiger Zeit auseinander.
Beim realen System treten Störungen auf, die im Modell nicht mit abgebildet werden, weshalb der Regler falsch reagieren würde.
Die Anfangsbedingungen des realen Systems sind i.d.R. nicht genau genug bekannt. Die Simulation würde deshalb an einem vom realen System abweichenden Zustand starten.

Das Konzept wird aber praktisch anwendbar, wenn die Differenz zwischen realem und simuliertem Zustand y zusätzlich auf den Eingang der Simulation zurückgeführt wird.
Das entsprechende Schaubild sieht dann folgendermaßen aus:

Zustandsregler mit Beobachtermodell: Die Differenz zwischen realem Systemzustand und simuliertem Systemzustand wird auf den Eingang der simulierten Regelstrecke rückgekoppelt.

Bild 0-6: Zustandsregler mit Beobachtermodell: Die Differenz zwischen realem Systemzustand und simuliertem Systemzustand wird auf den Eingang der simulierten Regelstrecke rückgekoppelt.

Um die Wirksamkeit der Rückkopplung der Differenz zwischen realer und simulierter Zustandsgröße einstellen zu können, wird noch eine Gewichtsmatrix L eingeführt.
Sinnvollerweise werden nur diejenigen Zustände aus der Simulation verwendet, die sich nicht direkt messen lassen.

Umsetzung mit Scilab

Anhand des vorangehend behandelten Zustandsreglers für den Elektromotor soll das Konzept des Zustandsreglers mit Beobachter besser verdeutlicht werden:

1. Idealfall:

Das Modell entspricht exakt dem realen System, die Anfangsbedingungen sind vollständig bekannt, es treten keine Störungen auf.
In der folgenden Scilab-Simulation werden zwei Systeme simuliert.
Die Gleichungen 1 und 2 sollen das reale System repräsentieren.
Die Gleichungen 3 und 4 repräsentieren den Beobachter, der hier idealisiert exakt synchron mit dem realen System läuft:

function f = rechteSeite5(t,y)
    u = -r1*(y(3,1)-w_stat)-r2*(y(4,1)-i_stat)+k*w_stat+R*i_stat;
//"Reales" System:
    f(1,1)=(-D/J)*y(1,1)+(k/J)*y(2,1);
    f(1,2)=(-k/L)*y(1,1)+(-R/L)*y(2,1)+(1/L)*u;
//"simuliertes" System (hier identisch, AB bekannt!, das ist unrealistisch):
    f(1,3)=(-D/J)*y(3,1)+(k/J)*y(4,1);
    f(1,4)=(-k/L)*y(3,1)+(-R/L)*y(4,1)+(1/L)*u;
endfunction

Code 0-4: Simulationsmodell für den Idealfall eines Zustandsreglers mit Beobachter in Scilab-Formulierung.

2. Idealfall mit falschen Anfangsbedingungen:

Selbst bei diesem Idealfall funktioniert die Regelung nicht, wenn das zweite, die Simulation/ den Beobachter repräsentierende Modell mit Anfangsbedingungen integriert wird, die nicht denen des realen Systems entsprechen:

t = 0:0.01:1;
y0 = [0,0]';
ybb0 = [0,0,20,20]';
t0   = 0;
y5b  = ode(ybb0,t0,t,rechteSeite5); //falsche AB

Code 0-5: Simulation mit differerierenden Anfangsbedingungen zwischen Beobachter und realem System.

3. Rückkopplung der Differenz zwischen realem und simuliertem Zustand

Bei folgender Scilab-Funktion liegt ein Beobachter vor, bei dem die Differenz zwischen realem und simuliertem Zustand zusätzlich auf den Eingang gegeben wird.

//Beobachter 1: Rückführung der Zustandsdifferenz y-ysim:
function f = rechteSeite6(t,y)
    u = -r1*(y(3,1)-w_stat)-r2*(y(4,1)-i_stat)+k*w_stat+R*i_stat;
//"Reales" System:
    f(1,1)=(-D/J)*y(1,1)+(k/J)*y(2,1);
    f(1,2)=(-k/L)*y(1,1)+(-R/L)*y(2,1)+(1/L)*u;
//"simuliertes" System (hier identisch, AB verkehrt, aber Rückführung der Zustandsdifferenz):
    l1=5;
    l2=5;
    f(1,3)=(-D/J)*y(3,1)+(k/J)*y(4,1)            +l1*(y(1,1)-y(3,1));
    f(1,4)=(-k/L)*y(3,1)+(-R/L)*y(4,1)+(1/L)*u   +l2*(y(2,1)-y(4,1));
endfunction

Code 0-6: Beobachter mit Rückkopplung der Differenz zwischen realem und simuliertem Zustand.

4. Einer der Systemzustände ist nicht meßbar

Im folgenden Fall wird nun konsequent berücksichtigt, dass eine der Zustandsgrößen (y2==i, Stromstärke) nicht am realen System gemessen werden kann.
Sinnvollweise wird dann die meßbare Zustandsgröße y1 im Regler benutzt und nicht dessen simuliertes Pendant y3.
Außerdem kann die nicht meßbare Zustandsgröße auch nicht benutzt werden, um die Differenz zwischen realem und simuliertem Zustand zu bilden.
Demzufolge wird das Gewicht l2 der Matrix L zu Null gesetzt.

function f = rechteSeite7(t,y)
    u = -r1*(y(1,1)-w_stat)-r2*(y(4,1)-i_stat)+k*w_stat+R*i_stat;
//"Reales" System:
    f(1,1)=(-D/J)*y(1,1)+(k/J)*y(2,1);
    f(1,2)=(-k/L)*y(1,1)+(-R/L)*y(2,1)+(1/L)*u;
//"simuliertes" System (hier identisch, AB verkehrt, aber Rückführung der Zustandsdifferenz):
    l1=5;
    l2=0; //zu Null gesetzt, weil nicht meßbar.
    f(1,3)=(-D/J)*y(3,1)+(k/J)*y(4,1)            +l1*(y(1,1)-y(3,1));
    f(1,4)=(-k/L)*y(3,1)+(-R/L)*y(4,1)+(1/L)*u   +l2*(y(2,1)-y(4,1));
endfunction

Code 0-7: Realistischer Fall: y2 ist nicht meßbar.

Bild 0-7: Vergleich diverser Varianten von Beobachtern (vergl. Text).

Die blaue Kurve zeigt die Sprungantwort des Systems bei u=5.
Der grünen Kurve entspricht ein Zustandsregler, bei dem von Hand eine doppelte reelle Polstelle bei -20 vorgegeben wurde.
Bei alle weiteren Kurven wird das gleiche geregelte System (doppelte Polstelle bei -20) sowohl im Beobachter, als auch im "realen System" benutzt.
Rote Kurve: Man sieht, dass hier das Systemverhalten trotz identischer Modelle im realen und simulierten System von der grünen Kurve abweicht, da die Anfangsbedingungen der "Simulation" nicht mit denen des realen Systrems übereinstimmen.
Hellblaue Kurve: Durch Rückführung der Differenz zwischen realem und simuliertem Zustand kann das Problem der falschen AB teilweise behoben werden.
Violette Kurve: Hier ist y2 nicht meßbar, aber es wird das "reale" y1 benutzt.