Regelungstechnik 2 für EIT Teil 10/11: Zustandsraum-Darstellung für Mehrgrössen-Systeme (MIMO)
Version 1.0 Prof. Dr. David Zogg Institut für Automation IA Hochschule für Technik Fachhochschule Nordwestschweiz Windisch, Mai 2012
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
1/35
Dokumentenkontrolle
Änderungen
Version
Datum
Autoren
Bemerkung
1.0
02.05.2012
David Zogg
Erstellung
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
2/35
1.
Zweck
4
2.
Referenzen
4
3.
Symbol- und Abkürzungsverzeichnis
5
3.1.
Allgemein
5
3.2.
Gasturbinen-Modell
6
4.
Einführung
7
4.1.
Lernziele
7
4.2.
Von der klassischen Regeltechnik zur modellbasierten Regelung
7
4.3.
Zustandsraum-Darstellung (Repetition)
8
4.4.
Praxisbeispiel Gasturbine
9
4.5.
Teilmodelle der Gasturbine
11
4.6.
Serieschaltung der Teilmodelle in MATLAB
14
4.7.
Analyse der Regelstrecke in MATLAB
15
4.8.
Zustandsregler mit LQR-Entwurf auslegen (Repetition)
18
4.9.
Beobachter mit LQE-Entwurf auslegen
19
4.10.
Synthese zum kombinierten Beobachter/Zustandsregler (LQG-Regulator)
21
4.11.
Der Schritt zum Folgeregler mit Soll/Ist-Vergleich
23
4.12. Die kombinierte Auslegungsmethode für Beobachter und Zustandsregler (LQG/LTR)
24
4.13.
Anwendung des LQG/LTR-Entwurfs auf die Gasturbine
27
4.14.
Simulation des entworfenen Reglers unter realen Bedingungen
31
4.15.
Zusammenfassung
35
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
3/35
1.
Zweck
Das vorliegende Skript dient als Grundlage für das Modul „Regelungstechnik 2“ (rt2) im 6. Semester des Studiengangs Elektro- und Informationstechnik (EIT). Im Modul „Regelungstechnik 1“ (rt1) des 5. Semesters wurde die Auslegung von zeitkontinuierlichen PID- und Zustandsreglern behandelt. Im vorliegenden 10./11. Teil wird die Zustandsraum-Darstellung auf massive Mehrgrössen-Systeme (MIMO = Multiple Input Multiple Output) erweitert. Neben der bekannten LQR-Auslegungsmethode für Zustandsregler wird die verwandte LQEAuslegungsmethode für Beobachter eingeführt. Die Theorie wird anhand eines realen Beispiels einer modernen Gasturbine für ein Kombikraftwerk erklärt.
2.
Referenzen
[ 1 ] J.P. Keller: Zustandsraummethoden, Skript Systemtechnik RT2, Windisch 2009 [ 2 ] J. Zellweger, Zustandsregelkreise, Kapitel 5, Skript, Windisch [ 3 ] H. Gassmann: Theorie der Regelungstechnink, 2. Auflage, Verlag Harri Deutsch, Rapperswil/Frankfurt 2003 [ 4 ] H.P. Geering: Regelungstechnik, 3. Auflage, Springer Verlag, Zürich/Berlin 1994 [ 5 ] K. Müller: Entwurf robuster Regelungen, B.G. Teubner Stuttgart, 1996
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
4/35
3.
Symbol- und Abkürzungsverzeichnis
3.1. Allgemein A
Zustandsraumdarstellung: Systemmatrix, Dynamikmatrix
B
Zustandsraumdarstellung: Eingangsmatrix für Stellgrösse u
Bw, Bw Zustandsraumdarstellung: Eingangsmatrix für Prozessrauschen w C
Zustandsraumdarstellung: Ausgangsmatrix
D
Zustandsraumdarstellung: Durchgangsmatrix (Direct Feedthrough)
G
Regelstrecke (Übertragungsfunktion)
Gu, Gu Regelstrecke mit Stellgrössen u am Eingang Gd, Gd Regelstrecke mit Störgrössen d am Eingang H
Beobachter-Rückführmatrix (Beobachter-Verstärkungsmatrix)
I, In
Identitäts-Matrix (Dimension n x n)
K
Zustands-Rückführmatrix
L
Offener Loop (Kreisverstärkung G0)
Lf, Lf
Offener Loop des LQ-Estimators (Kalman-Filter, Beobachter)
Lr, Lr
Offener Loop des LQ-Regulators (Zustandsregler)
Qx, Qx
Gewichtungsmatrix für Zustandsgrössen (LQ-Regulator)
Qy, Qy Gewichtungsmatrix für Ausgangsgrössen (LQ-Regulator) Qw, Qw Gewichtungsmatrix für Prozessrauschen bzw. Schätzfehler (LQ-Estimator) Ru, Ru Gewichtungsmatrix für Stellgrössen (LQ-Regulator) Rv, Rv
Gewichtungsmatrix für Messrauschen (LQ-Estimator)
d
Störgrössen
e
Regelfehler
i
Imaginär-Anteil
s
Laplace-Operator (zeitkontinuierlich)
u
Eingangsgrössen (Stellgrössen, Steuergrösse)
w
Führungsgrössen (Sollwert)
x
Zustandsgrössen
y
Ausgangsgrössen (Messgrössen, Istwert)
ρ
Auslegungsparameter für LQ-Regulator (Methode LQG/LTR)
µ
Auslegungsparameter für LQ-Estimator (Methode LQG/LTR)
ω
Kreisfrequenz (rad/s)
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
5/35
LQ
Linear-Quadratisch
LQR
Linear-Quadratischer Regulator
LQE
Linear-Quadratischer Estimator
LQG
Linear-Quadratisch mit Gauss’scher Verteilung (der Rauschsignale)
LTR
Loop Transfer Recovery
3.2. Gasturbinen-Modell Aij
Wärmeübergangsfläche zwischen Bilanzelement i und j
cv,i
Spezifische isochore Wärmekapazität des Mediums im Bilanzelement i
cp
Spezifische isobare Wärmekapazität des Mediums (Stoffstrom, Massenstrom)
hBr
Spezifische Enthalpie des Brennstoffes
mi
Masse des Mediums im Bilanzelement i
kij
Wärmedurchgangszahl zwischen Bilanzelement i und j
mi
Masse des Mediums im Bilanzelement i
m*i
Massenstrom des Mediums durch das Bilanzelement i
m*Br,i
Brennstoff-Massenstrom im Bilanzelement i (Zufuhr in Brennkammer)
Ti
Temperatur des Mediums im Bilanzelement i (= Temperatur am Austritt)
Tin,i
Temperatur des Mediums am Eintritt des Bilanzelementes i
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
6/35
4.
Einführung
4.1. Lernziele Teil 10 (Zustandsraum-Darstellung) Lernziel Repetition: Schritt von der klassischen Regelung zur modellbasierten Regelung Zustandsraum-Darstellung auf massive Mehrgrössensysteme mit vielen Ein-/Ausgängen anwenden können
Taxonomiestufe (Bloom) Verständnis Anwendung
Teil 11 (Auslegungsmethoden) Lernziel Repetition: Zustandsregler nach LQR-Methode auslegen können Neu: Optimalen Beobachter nach LQE-Methode auslegen können Mehrgrössen-Folgeregelung mit Beobachter und Zustandsregler verstehen Kombinierte Auslegungsmethode LQR/LTR anwenden können
Taxonomiestufe (Bloom) Anwendung Anwendung Verständnis Anwendung
4.2. Von der klassischen Regeltechnik zur modellbasierten Regelung In der klassischen Regeltechnik wird der Prozess (Regelstrecke) zwar analysiert, für die Reglerauslegung ist aber kein explizites Prozessmodell notwendig. In einem Grossteil der praktischen Anwendungen ist deshalb der PID-Regler vorherrschend (Abbildung 1). Dieser Regler kann auch dann eingestellt werden, wenn die Regelstrecke nicht genau bekannt ist.
Abbildung 1: Regelkreis mit klassischem PID-Regler
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
7/35
Die klassische Regeltechnik stösst aber bei folgenden Situationen an ihre Grenzen: • • • • •
Massive Mehrgrössensysteme mit vielen Ein- und Ausgangsgrössen Starke Kopplungen zwischen den Signalen (bzw. Zustandsgrössen) Nicht alle Grössen messbar (Zustands-Beobachter notwendig) Hochdynamische („schnelle“) Regelstrecken, welche hohe Reglerperformance benötigen (z.B. aktive Schwingungsdämpfung in mechanischen Systemen) Sicherheitskritische oder sehr langsame Regelstrecken, welche nicht von Hand oder durch Einstellregeln (Sprungantwort, Ziegler-Nichols) eingestellt werden können
In solchen Situationen ist es ratsam, zunächst ein Prozessmodell aufzustellen. Dieses Modell kann anschliessend für die Auslegung des Reglers verwendet werden (Abbildung 2). Dazu sind folgende Möglichkeiten vorhanden: • •
Verwendung des Prozessmodells zur mathematischen Auslegung des Reglers ( modellbasierter Zustandsregler) Direkte „Einpflanzung“ des Prozessmodells in den Regler ( Verwendung eines Beobachters mit Zustandsregler)
Abbildung 2: Regelkreis mit modellbasiertem Regler
4.3. Zustandsraum-Darstellung (Repetition) Die Zustandsraum-Darstellung ist durch folgendes Differentialgleichungssystem gegeben: (1)
x& = Ax + Bu y = Cx + Du
Die Stellgrössen sind im Eingangsvektor u zusammengefasst, die (nicht messbaren) Zustandsgrössen im Zustandsvektor x und die (messbaren) Ausgangsgrössen im Ausgangsvektor y. Die Systemmatrix A ist für die Dynamik der Regelstrecke zuständig, die Eingangsmatrix B beschreibt den Einfluss der Stellgrössen, die Ausgangsmatrix C selektiert die messbaren Ausgangsgrössen und die Durchgangsmatrix D definiert die (statischen) direkten Zusammenhänge zwischen Eingang u und Ausgang y (oft ist D=0). Das zugehörige Singalflussbild (Blockdiagramm) zeigt Abbildung 3
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
8/35
Abbildung 3: Allgemeines Signalflussbild für Zustandsraum-Darstellung
Folgend ist die Zustandsraum-Darstellung für ein System 3. Ordnung mit 3 Eingangsgrössen u1..u3, drei Zustandsgrössen x1..x3 und einer gemessenen Ausgangsgrösse y dargestellt (Beispiel):
x&1 a11 a12 x& = a a 2 21 22 x&3 a31 a32 y = [c1 c2 (2)
a13 x1 b11 b12 b13 u1 a23 ⋅ x2 + b21 b22 b23 ⋅ u2 a33 x3 b31 b32 b33 u3 x1 u1 c3 ]⋅ x2 + [d1 d 2 d 3 ]⋅ u2 x3 u3
Falls auch der Einfluss der Störgrössen d modelliert werden soll, wird folgende erweiterte Variante gewählt (Beispiel mit 2 Störgrössen d1 und d2 am Eingang):
(3)
Bu A 448 6447 644 47 444 8 647 4Bd 48 4 x&1 a11 a12 a13 x1 bu11 bu12 bu13 u1 bd 11 bd 12 x& = a d1 2 21 a22 a23 ⋅ x2 + bu 21 bu 22 bu 23 ⋅ u 2 + bd 21 bd 22 ⋅ d x&3 a31 a32 a33 x3 bu 31 bu 32 bu 33 u3 bd 31 bd 32 2 Du Dd C 48 x1 647 6447 448 u1 647 48 d y = [c1 c2 c3 ] ⋅ x2 + [d u1 d u 2 d u 3 ] ⋅ u 2 + [d d 1 d d 2 ] ⋅ 1 d 2 x3 u3
Dabei wird die Eingangsmatrix in zwei Teile Bu (Einfluss der Stellgrössen u) und Bd (Einfluss der Störgrössen d) aufgeteilt. Ebenso wird die Durchgangsmatrix in Du und Dd aufgeteilt.
4.4. Praxisbeispiel Gasturbine Gasturbinen dienen in modernen Kombikraftwerken der Stromerzeugung (Abbildung 4). Durch eine Kopplung mit einer Dampfturbine (Nutzung der Abwärme der Gasturbine) werden heute Kombi-Wirkungsgrade von bis zu 60% erreicht. Speziell dazu wurden von ABB (heute Alstom) Gasturbinen mit sequentieller Verbrennung entwickelt. Zwei hintereinandergeschaltete Brennkammern ermöglichen eine vollständige Verbrennung mit kleiner NOx-Emission bei gleichzeitig hoher Austrittstemperatur, was wiederum für den nachgeschalteten Dampfkreislauf von Vorteil ist.
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
9/35
Abbildung 4: Gasturbine GT24/26 von ABB/Alstom mit 240 MW Leistung
Diese Gasturbinen hatten in den Anfängen allerdings eine „Kinderkrankheit“, es entstanden nämlich unerklärliche Pulsationen in den Brennkammern. Das verursachte hohe Kosten durch Ausfälle. Natürlich wurde fieberhaft nach Lösungen gesucht. Ein Ansatz betraf die Regeltechnik. Die konventionelle Regelung der Abgastemperatur war über einen einschleifigen PID-Regelkreis gelöst (Abbildung 5). Damit wurden die Werte der Temperatursensoren gemittelt, welche im Austrittgehäuse über dem Umfang verteilt waren. Die gemittelte Temperatur wurde anschliessend über einen PID-Regler auf die Brennstoffzufuhr geführt, welche wiederum alle Brenner zusammen ansteuerte. Damit war es nicht möglich, Brenner einzeln anzusteuern bzw. Störungen in den Temperaturverteilungen auszuregeln.
Abbildung 5: Klassische Regelung über einschleifigen PID-Regler
In einem Zwischenschritt wurde der einschleifige Regelkreis durch mehrere parallel angeordnete PID-Regelkreise ersetzt, welcher jeder für sich eine Abgastemperatur auf die entsprechende Brennstoffzufuhr zurückführt (Abbildung 6). Dabei muss berücksichtigt werden, dass in der Turbine eine Ablenkung der Strömung durch die Leitund Laufschaufeln entsteht. Die Sensoren müssen also den richtigen Aktuatoren zugeordnet werden. Zudem kann man sich einfach vorstellen, dass die Brennstoffzufuhr an einem Brenner auch die Temperatur der seitlich benachbarten Kammern beeinflusst,
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
10/35
da eine Vermischung stattfindet. Diese Kopplungen können mit den (getrennt) parallel laufenden PID-Reglern nicht berücksichtigt werden.
Abbildung 6: Klassische Regelung über mehrere parallele PID-Regler (hellgelbe Pfeile: Kopplungen)
Als nächster Schritt wurde deshalb der modellbasierte Ansatz gewählt. Damit werden die Kopplungen automatisch berücksichtigt, weil sie bereits im Prozessmodell vorhanden sind. Auch die Ablenkung in der Turbine wird automatisch berücksichtigt.
Abbildung 7: Modellbasierte Regelung (Prozessmodell mit Kopplungen)
4.5. Teilmodelle der Gasturbine Nun wollen wir das Prozessmodell für die Regelstrecke der Gasturbine aufstellen. Dazu teilen wir das System in mehrere Bilanzelemente auf (Abbildung 8). Für die Teilmodelle Brennkammer, Turbine und Austrittsgehäuse sind je n Bilanzelemente über den Umfang der Maschine verteilt (n ≥ Anzahl Temperatursensoren am Austritt). Das sind insgesamt 3 x n Bilanzelemente (für n = 6 also 3 x 6 = 18 Bilanzelemente).
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
11/35
Abbildung 8: Bilanzelemente, verteilt über den Umfang
Für jedes Bilanzelement wird die thermodynamische Gleichung der Wärmebilanz aufgestellt. Folgende Differentialgleichung gilt für ein Bilanzelement der Brennkammer:
(4)
Die Symbole zu ( 4 ) sind im Abschnitt 3.2 definiert. Die Gleichung ( 4 ) ist wie folgt zu lesen: „Die Veränderung der inneren Energie des Bilanzelementes (Dynamik) ist gleich der zugeführten Wärme von den benachbarten Elementen j und k (Wärmeübertragung) plus der konvektiv zugeführten Energie (Konvektion) plus der zugeführten Brennstoffenergie (Brennstoff)“. Gleichung ( 4 ) ist in allgemeiner Form dargestellt und kann nun für alle Bilanzelemente einzeln angewandt werden. Für das Teilmodell der Brennkammer wird die Gleichung in die Form ( 5 ) gebracht:
(5)
Wird das System nun für n gekoppelte Bilanzelemente aufgestellt, resultiert folgende Zustandsraum-Darstellung für die Brennkammer (Beispiel n=6, Übung!):
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
12/35
(6)
Die Brennstoffmassenströme m*Br,i gehen als bekannte Stellgrössen u ein, während dem die unbekannten Temperaturen Tin,i am Eintritt der Brennkammer als Störgrössen ud eingehen. Dabei sind die Koeffizienten in ( 6 ) wie folgt definiert:
(7)
In einem weiteren Schritt wird das Teilmodell der Turbine betrachtet. Hier wird nur die Ablenkung durch die Leit- und Laufschaufeln berücksichtigt (kein Wärmeaustausch), was durch die „Shift-Matrix“ ( 8 ) definiert werden kann. Damit wird die Temperatur des zweiten Elementes am Eintritt dem ersten Element am Austritt zugeordnet, die Temperatur des dritten Elementes dem zweiten, usw.
(8)
ud y K 4448 67 678 64447 8 T1,out 0 k 0 0 0 0 T1,in T 2,out 0 0 k 0 0 0 T2,in T3, out 0 0 0 k 0 0 T3,in = ⋅ T4,out 0 0 0 0 k 0 T4,in T5,out 0 0 0 0 0 k T5,in T6,out k 0 0 0 0 0 T6,in
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
13/35
Als drittes Teilmodell folgt nun noch das Austrittsgehäuse. Hier wird wieder die Gleichung ( 5 ) bzw. die Zustandsraum-Darstellung ( 6 ) als Basis genommen, wobei hier keine Brennstoffmassenströme als Stellgrössen vorhanden sind (u = 0, also Bu = 0):
(9)
Für das betrachtete Teilmodell sind die Temperaturen Ti,in in ( 9 ) als Störgrössen ud zusammengefasst, wobei diese gleich den Austrittstemperaturen des vorgeschalteten Teilmodelles sind (also der Turbine). Damit ist der nächste Schritt schon definiert: Die drei Teilmodelle müssen nun in Serie hintereinander geschaltet werden, was in folgendem Abschnitt mit Hilfe von MATLAB durchgeführt wird.
4.6. Serieschaltung der Teilmodelle in MATLAB Die Teilmodelle der Brennkammer, Turbine und des Austrittsgehäuses werden nun in mit Hilfe von MATLAB in Serie geschaltet. Wir verwenden dazu folgende MATLAB-Befehle: Funktion Teilmodell Brennkammer ( 6 ) mit Koeffizienten ( 7 )
Teilmodell Turbine ( 8 ) mit Koeffizient k Teilmodell Austrittsgehäuse ( 9 ) mit Koeffizienten ( 7 )
keine Brennstoffzufuhr: Serieschaltung der Teilmodelle zum Gesamtsystem „Gasturbine“ Auswahl der Eingangsgrössen: nur Stellgrössen u nur Störgrössen d
MATLAB-Befehl GBrK = ss(A, [Bu Bd], C, [Du Dd]) aik = 1; bid = 10; aii = 2*aik+bid; biu = 500; GTrb = K k = 0.8 GAus = ss(A, Bd, C, Dd) aik = 0.5; bid = 1; aii = 2*aik+bid; biu = 0; G = GAus*GTrb*GBrK
Gu = G(:, 1:6); Gd = G(:, 7:12);
File: RT2_Strecke_GT.m
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
14/35
4.7. Analyse der Regelstrecke in MATLAB Die Analyse der Regelstrecke erfolgt wiederum in MATLAB. Es ist zu beachten, dass wir es hier mit einem Mehrgrössensystem zu tun haben, welches 6 Eingänge und 6 Ausgänge hat (MIMO: Multiple Input Multiple Output). Ein Bode-Diagramm ist aber nur für ein System mit einem Eingang und einem Ausgang definiert (SISO: Single Input Single Output). Anstelle des Bode-Diagrammes verwenden wir deshalb hier sogenannte „Singularwertverläufe“. Das ist im Prinzip der Amplitudengang für ein Mehrgrössensystem (Phasengang entfällt). Die MATLAB-Funktion „sigma“ stellt die Singularwerte über der Frequenz dar. Definition der Singularwerte [ 5 ], S.167 (Zusatzinfo für Interessierte): Die singulären Werte σi geben an, um welchen Faktor sich ein Vektor x durch die Abbildung an der Matrix Ax in der Länge ändert. Im Gegensatz zu den Eigenwerten λi kann sich dabei die Richtung des Vektors ändern. Wie auch bei den Eigenwerten gilt: ( 10 )
det [ A] =
n
n
∏σ i =∏ λi i =1
i =1
Die singulären Werte charakterisieren das „Verstärkungsverhalten“ einer Matrix besser als die Eigenwerte, denn die grösste bzw. kleinste Längenänderung erfolgt nicht zwangsläufig auch in der Richtung des Eingangsvektors x. Dies kommt durch die Eigenschaft ( 11 )
σ max ≥ λmax
σ min ≤ λmin
zum Ausdruck. σmax ist der grösste Singularwert. Er beschreibt die maximal mögliche Längenänderung eines Vektors durch Abbildung an einer Matrix und kann als Verallgemeinerung des Betrags von Skalaren aufgefasst werden. Für die Analyse der Regelstrecke verwenden wir folgende MATLAB-Befehle: Funktion Pole und Nullstellen Impulsantwort auf Stellgrössen u Bode-Diagramm Singularwertverläufe
MATLAB-Befehl pzmap(Gu) impulse(Gu) bode(Gu) sigmaplot(Gu)
File: RT2_Strecke_GT.m Abbildung 9 zeigt
die Impulsantwort für das MIMO-System mit 6 Stellgrössen (horizontal) und 6 Ausgangsgrössen (vertikal). Es ist eindeutig die Struktur der Matrix A in ( 6 ) zusammen mit der „Shift-Matrix“ ( 8 ) ersichtlich. Die grösste Impulsantwort liegt in der Kammer mit dem Brennstoff-Impuls vor, während kleinere Impulsantworten in den benachbarten Kammern vorliegen. Die Reaktionen in den benachbarten Kammern erfolgen aufgrund der seitlichen Wärmeübertragung bzw. Kopplungen. Abbildung 10 zeigt das Bode-Diagramm für das MIMO-System.
Da das Bode-Diagramm nur für SISO-Systeme definiert ist (siehe Bemerkung oben), werden hier 6x6 BodeDiagramme für jede Eingangs-/Ausgangskombination dargestellt. Daraus können keine sinnvollen Aussagen für das Gesamtsystem gemacht werden. Abbildung 11 zeigt die Singularwertverläufe, welche nun anstelle des Bode-Diagrammes
treten. Die Singularwertverläufe verhalten sich ähnlich wie Amplitudengänge. Hier sind die zwei kleinsten und zwei grössten Singularwertverläufe dargestellt (kleinste und grösste Verstärkungen des Systems).
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
15/35
Abbildung 9: Impulsantwort auf die Stellgrössen. In(1)..In(6) = Stellgrössen u1..u6, Out(1)..Out(6) = Ausgangsgrössen y1..y6
Abbildung 10: Bode-Diagramm für das MIMO-System. In(1)..In(6) = Stellgrössen u1..u6, Out(1)..Out(6) = Ausgangsgrössen y1..y6
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
16/35
Abbildung 11: Singularwertverläufe des MIMO-Systems (Durchtrittsfrequenz bei ca. 20 rad/s).
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
17/35
4.8. Zustandsregler mit LQR-Entwurf auslegen (Repetition) Ist die Regelstrecke in der Zustandsraum-Darstellung gegeben, kann daraus der Zustandsregler berechnet werden. Es resultiert das geregelte System in Abbildung 12. Dieser Regelkreis bewirkt nichts anderes, als die Zustandsgrössen mit den Anfangsbedingungen x0 auf Null zu bringen (x 0). Es ist noch kein Soll-/Istvergleich wie in einem konventionellen Regelkreis vorhanden. Dieser wird später eingeführt.
Abbildung 12: Strecke im Zustandsraum (A,B,C,D) mit Zustandsrückführung (-K)
Der Zustandsregler tritt als Matrix K auf, welche die Zustandsgrössen x auf die Eingangsgrössen (Stellgrössen) u zurückführt. Es ist also ein reiner P-Regler, welcher allerdings die Kopplungen zwischen den Zustandsgrössen berücksichtigen kann (ausserdiagonale Elemente der K-Matrix). Die Matrix K wird aufgrund einer Optimierung mit folgendem Gütekriterium berechnet: T
( 12 )
J x ,u (u ) =
∫ [x
T
]
(t ) ⋅ Qx ⋅ x(t ) + u T (t ) ⋅ Ru ⋅ u (t ) ⋅ dt
0
Das Gütekriterium ( 12 ) wird minimiert, was in einer Minimierung der Zustandsgrösse x(t) über der Zeit sowie der Stellgrösse u(t) über der Zeit resultiert. Mit den Gewichtungsmatrizen Qx und Ru kann das Verhältnis von Zustands- zu Stellgrössenminimierung gewählt werden. Damit kann das Verhältnis der Regelgüte (Performance) zur verwendeten Energie (Stellgrösse) eingestellt werden (wie eine Art Schieberegler!). Anstatt der (nicht messbaren) Zustandsgrössen x(t) können auch die (messbaren) Ausgangsgrössen y(t) gewichtet werden: T
( 13 )
J y ,u (u ) =
∫ [y
T
]
(t ) ⋅ Q y ⋅ y (t ) + u T (t ) ⋅ Ru ⋅ u (t ) ⋅ dt
0
In ( 13 ) tritt anstelle von Qx die Gewichtungsmatrix Qy. Zwischen Qx und Qy besteht folgender Zusammenhang: ( 14 )
Qx = C T ⋅ Q y ⋅ C
Dies kann durch Einsetzen von ( 14 ) in ( 12 ) einfach nachvollzogen werden. Das Optimierungsproblem wird nun wie folgt definiert: Minimierung von ( 12 ) unter Berücksichtigung der folgenden Nebenbedingungen: ( 15 )
x& = ( A − BK ) ⋅ x
x(0) = x0
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
18/35
Die Nebenbedingungen ( 15 ) bedeuten nichts anderes, dass die Differentialgleichung des Regelkreises mit Zustandsrückführung sowie die Anfangsbedingungen x(0) eingehalten werden. Damit wird auch die Bedeutung von LQR klar, was nichts anderes als „Linear Quadratic Regulator“ bedeutet. Linear ist die Regelstrecke, während das Gütekriterium quadratisch ist. Ohne auf die Herleitung einzugehen, wird als Lösung des Optimierungsproblems die Matrix-Riccatigleichung gefunden: ( 16 )
−1
AT P + PA − PBRu B T P + Qx = 0
welche folgende Zustandsrückführung K ergibt: ( 17 )
−1
K = Ru B T P
Der entsprechende MATLAB-Befehl lautet: [K,P,E] = lqr(A,B,Qx,Ru) mit K = Zustandsrückführmatrix P = Lösung der Matrix-Riccatigleichung E = Eigenwerte von A-BK (Dynamik Regelkreis mit Zustandsrückführung)
4.9. Beobachter mit LQE-Entwurf auslegen Wir wollen nun den Schritt vom Zustandsregler zum Beobachter wagen. Der Beobachter kann mit einer ähnlichen Optimierungsmethode ausgelegt werden. Doch zunächst zur allgemeinen Struktur eines Beobachter, welche in Abbildung 13 dargestellt ist. Der Beobachter dient dazu, die Zustandsgrössen zu schätzen. Dies ist insbesondere dann wichtig, wenn nicht alle Zustandsgrössen gemessen werden können. Ohne Beobachter könnte in diesem Fall kein Zustandsregler ausgelegt werden. Im Beobachter ist im Prinzip eine „Kopie“ der Regelstrecke abgebildet, wobei die Zustandsraum-Matrizen (A,B,C,D) im Beobachter nie genau den realen Matrizen im Prozess entsprechen. Zudem ist der Anfangswert x0 im Allgemeinen nicht genau bekannt. Aus diesen Gründen müssen die Schätzwerte xˆ im Beobachter laufend korrigiert werden, was durch einen permanenten Vergleich der geschätzten Ausgangswerte yˆ mit den gemessenen Ausgangsgangswerten y geschieht. Der Fehler wird über die BeobachterRückführung H auf die Ableitung der Schätzwerte x&ˆ geführt. Auch hier ist also wieder ein reiner P-Regler im Spiel (mit entsprechenden Kopplungen in den ausserdiagonalen Elementen von H). Den Beobachter kann man auch als Ersatz für einen „Differentiator“ ansehen, welcher aus den bekannten Ausgangsgrössen y die Ableitungen x&ˆ berechnet. Nur macht er dies „besser“ als ein Differentiator, in dem er vorwärts integriert und damit wie ein Filter wirkt. Die Beobachter-Rückführung H definiert, wie stark das Filter auf sich verändernde Messgrössen y reagiert. Damit ist der Beobachter bei entsprechender Auslegung auch für verrauschte Messgrössen geeignet (im Gegensatz zum Differentiator).
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
19/35
Abbildung 13: Beobachter Struktur mit Beobachter-Rückführung L
Der Beobachter bzw. dessen Rückführmatrix H kann auf verschiedene Arten ausgelegt werden, u.A. über Polvorgabe. Wir wollen hier aber die optimale Auslegung analog zum Zustandsregler betrachten. Folgende Berechnungsmethode führt auch zum sogenannten „Kalman-Filter“. Dazu werden zwei zusätzliche (künstliche) Rauschsignale definiert:
•
•
Prozessrauschen w am Eingang des Prozesses. Diese Störgrösse modelliert Rauschvorgänge im Prozess selber und kann als zusätzliche Stör-Eingangsgrösse betrachtet werden. Der Einfluss auf den Prozess wird über die Matrix G definiert. Man kann die Grösse Bw⋅w auch als Abweichung zwischen x& und x&ˆ bzw. w als Mass für den Schätzfehler ( e = x − xˆ ) betrachten. Messrauschen v am Ausgang des Prozesses. Diese Störgrösse modelliert Rauschvorgänge in Sensoren und wirkt direkt auf die Ausgangsgrösse y . Damit wird also der Einfluss des Messrauschens modelliert.
Die Matrix H kann aufgrund einer Optimierung mit folgendem Gütekriterium berechnet werden: T
( 18 )
J w,v ( w) =
∫ [w
T
]
(t ) ⋅ Qw ⋅ w(t ) + v T (t ) ⋅ Rv ⋅ v (t ) ⋅ dt
0
Das Gütekriterium ( 18 ) wird minimiert, was in einer Minimierung des Einflusses des Prozessrauschens w(t) über der Zeit sowie des Einflusses des Messrauschens v(t) über der Zeit resultiert. Mit den Gewichtungsmatrizen Qw und Rv kann das Verhältnis von Schätzfehler-Minimierung zu Messrausch-Unterdrückung gewählt werden (also wieder eine Art Schieberegler!). Die Matrizen Qw und Rv können auch durch folgende Kovarianzmatrizen beschrieben werden: ( 19 )
E ( ww T ) = Qw
E (vvT ) = Rv
E ( w) = E (v) = E ( wv T ) = 0
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
20/35
In ( 19 ) wird der Erwartungswert des quadrierten Prozessrauschens w mit Qw bezeichnet, der Erwartungswert des quadrierten Messrauschens v mit Rv und die beiden Rauschsignale sind unkorreliert (Kreuzvarianzmatrix = 0). Das Optimierungsproblem wird nun wie folgt definiert: Minimierung von ( 18 ) unter Berücksichtigung der folgenden Nebenbedingungen: ( 20 )
x&ˆ = ( A − HC ) ⋅ xˆ + ( B − HD) ⋅ u + H ⋅ y
xˆ (0) = xˆ0
Die Nebenbedingungen ( 20 ) bedeuten hier, dass die Differentialgleichung des Beobachters mit Beobachterückführung sowie die Anfangsbedingungen eingehalten werden. Der Begriff LQE bedeutet „Linear Quadratic Estimator“. Linear ist die Regelstrecke (bzw. die Beobachtergleichung), während das Gütekriterium quadratisch ist. Der resultierende Beobachter wird engl. mit „Estimator“ („Schätzer“) bezeichnet. Ohne auf die Herleitung einzugehen, wird als Lösung des Optimierungsproblems die Matrix-Riccatigleichung gefunden, nun in der Beobachter-Form: ( 21 )
−1
AΣ + ΣAT − ΣC T Rv CΣ + BQw B T = 0
welche folgende Beobachterrückführung H ergibt: ( 22 )
H = ΣC T Rv
−1
Der entsprechende MATLAB-Befehl lautet: [H,S,E] = lqe(A,Bw,C,Qw,Rv) mit H = Beobachter-Verstärkungsmatrix S = Lösung der Matrix-Riccatigleichung E = Eigenwerte von A-LC (Beobachter-Dynamik)
4.10. Synthese zum kombinierten Beobachter/Zustandsregler (LQGRegulator) Wirklich interessant wird es nun, wenn der optimale Beobachter mit dem optimalen Zustandsregler zu einem vollständigen Regelkreis geschlossen wird, was Abbildung 14 zeigt. Die Zustandsgrössen x werden durch den Beobachter geschätzt und anschliessend über den Zustandsregler auf die Eingangsgrössen u zurückgeführt. Dies funktioniert auch dann, wenn nicht alle Zustandsgrössen gemessen werden. Wir haben es hier mit dem LQG-Regulator zu tun, wobei L für die lineare Regelstrecke, Q für das quadratische Gütekriterium und G für das Gauss’sche Rauschen (w und v) steht. Die Gleichungen des vollständigen Beobachters mit Zustandsregler lauten wie folgt: x&ˆ = ( A − HC − BK + HDK ) ⋅ xˆ + H ⋅ y ( 24 ) u = − K ⋅ xˆ ( 23 )
Gleichung ( 23 ) beschreibt die Schätzung der Zustandsgrössen xˆ aus den Messgrössen y und ( 24 ) beschreibt die Zustandsrückführung auf die Stellgrössen u .
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
21/35
Abbildung 14: Kombination von Beobachter und Zustandsrückführung
In MATLAB kann das System wie folgt ausgelegt und zusammengesetzt werden: K = lqr(A,B,Qx,Ru) H = lqe(A,Bw,C,Qw,Rv) Die kombinierte Zustandsraumdarstellung von Beobachter und Zustandsregler lautet wie folgt (Eingang y , Ausgang u , Zustandsgrössen xˆ ): Alqg Blqg Clqg Dlqg
= = = =
A-HC-BK+HDK H -K 0
Klqg = ss(Alqg, Blqg, Clqg, Dlqg) Als Alternative können in MATLAB auch folgende Befehle verwendet werden: Ge = ss(A, [B Bw], C, [D H]) K = lqr(A,B,Qx,Ru) Kest = kalman(Ge,Qw,Rv) Klqg = lqgreg(Kest, K) Die Zusammensetzung des LQG Regulators erfolgt gemäss folgender Grafik:
Abbildung 15: Synthese zum LQG-Regulator (Quelle: MATLAB Control Toolbox)
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
22/35
4.11. Der Schritt zum Folgeregler mit Soll/Ist-Vergleich In Abbildung 14 haben wir noch keinen Soll/Ist-Vergleich, welchen wir nun einführen wollen. Dies erfolgt mit einem einfachen Schritt, wie Abbildung 16 zeigt. Anstelle des Messwertes y wird nun die Differenz zwischen Sollwert w und Istwert y, also der Regelfehler e in den Beobachter geführt. Man beachte die Vorzeichen, welche hier umgekehrt definiert werden, um keinen Vorzeichenwechsel im Regelkreis herbeizuführen (negative Rückführung –K muss erhalten bleiben).
Abbildung 16: Einführung des Sollwertes w und Regelfehlers e (Soll/Ist-Vergleich)
Aus obiger Abbildung wird klar, dass der Beobachter nicht mehr die absoluten Zustandsgrössen x , sondern nur noch deren Abweichungen vom Betriebspunkt schätzt. Bei einem Regelfehler null sind auch die geschätzten Zustandsgrössen xˆ null. Dies ist grundsätzlich nichts Neues, da wir hier sowieso von einem linearen Modell der Regelstrecke ausgehen, welches in der Realität in einem bestimmten Betriebspunkt linearisiert wird. Durch Umzeichnen wird nun die bekannte, klassische Darstellung des Regelkreises mit Regler und Strecke in Abbildung 17 erreicht (die Durchgangsmatrix D wurde hier 0 gesetzt). Der Regler besteht nun aus den Teilen des Beobachters und Zustandsreglers. Damit ist auch sofort ersichtlich, wie die Strecke im Regler abgebildet wird (Pfeil). Wir sind also beim Ziel von Abbildung 2 angelangt.
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
23/35
Abbildung 17: Folgeregler in klassischer Darstellung
4.12. Die kombinierte Auslegungsmethode für Beobachter und Zustandsregler (LQG/LTR) Wir haben nun einen kombinierten Regler, bestehend aus Beobachter und Zustandsregler. Wir kennen auch die Auslegungsmethoden, mit welchen der Beobachter und der Zustandsregler einzeln ausgelegt werden:
•
LQR-Methode für die Auslegung des Zustandsreglers mit quadratischem GüteT
kriterium J x ,u (u ) =
∫ [x
T
]
(t ) ⋅ Qx ⋅ x(t ) + u T (t ) ⋅ Ru ⋅ u (t ) ⋅ dt
0
•
LQE-Methode für die Auslegung des Beobachter bzw. Kalman-Filters mit quadT
ratischem Gütekriterium J w,v ( w) =
∫ [w
T
]
(t ) ⋅ Qw ⋅ w(t ) + v T (t ) ⋅ Rv ⋅ v(t ) ⋅ dt
0
Damit haben wir insgesamt vier Gewichtungsmatrizen für das „Design“ des Reglers (Beispiel mit m=3 Eingängen und p=3 Ausgängen, Ausserdiagonale Elemente der Matrizen auf 0 gesetzt):
•
Gewichtungsmatrizen für den LQ-Regulator (Zustandsregler):
•
Gewichtungsmatrizen für den LQ-Estimator (Beobachter):
Damit haben wir für obiges Beispiel insgesamt 12 Auslegungsparameter! Die Frage ist, wie man diese am besten einstellt. Für den Zustandsregler alleine war dies noch relativ einfach durch Ausprobieren möglich. Man konnte z.B. mit folgenden Gewichtungen starten:
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
24/35
Schritt 1:
1 0 0 Ru = 0 1 0 0 0 1
1 0 0 Qx = 0 1 0 0 0 1
um dann sukzessive eine Gewichtungsmatrix um Zehnerpotenten zu erhöhen, z.B.: Schritt 2:
1 0 0 Ru = 0 1 0 0 0 1
10 0 0 Qx = 0 10 0 0 0 10
Man hat dann den Effekt auf die Sprungantwort betrachtet und bei erreichter Performance die Gewichtungsmatrizen „eingefroren“. Geübte konnten sich auch darin versuchen, die Diagonalelemente verschieden hoch zu gewichten, usw. Man hatte damit zwar viele Möglichkeiten, die Einstellung brauchte aber einiges an Übung und Erfahrung. Wenn der Beobachter nun noch dazukommt, wird das Einstellen von vier statt nur zwei Gewichtungsmatrizen noch wesentlich anspruchsvoller. Deshalb wollen wir dieses Vorgehen nun vereinfachen, indem wir auf ein kombiniertes Auslegungsverfahren zurückgreifen. Dieses kombinierte Auslegungsverfahren heisst LQG/LTR [ 4 ]: • •
LQG = Linear Quadratic Gaussian LTR = Loop Transfer Recovery
Den Namen des ersten Schrittes LQG kennen wir bereits, während wir den Namen des zweiten Schrittes bald begreifen werden. Wir wollen das Ganze praktisch an den Gewichtungsmatrizen erklären. Das Ziel des Verfahrens muss ja sein, die Wahl der Gewichtungsmatrizen massiv zu vereinfachen. Dies erreichen wir mit folgendem Ansatz (Beispiel wieder mit m=3, p=3): •
Gewichtungsmatrizen für den LQ-Regulator (Zustandsregler):
•
Gewichtungsmatrizen für den LQ-Estimator (Beobachter):
Wir sehen sofort, dass wir nur noch 2 Auslegungsparameter haben statt 12 wie vorher! Wir haben also nur noch einen Parameter ρ, welcher die Auslegung des Zustandsreglers beeinflusst, und einen Parameter µ, welcher die Auslegung des Beobachter beeinflusst. Wir erinnern uns auch an die Idee der „Schieberegler“, mit welchem wir die Performance des Zustandsreglers und des Beobachters festlegen konnten. Die Parameter ρ und µ kann man sich also je durch einen „Schieberegler“ vorstellen. Damit können wir alles einstellen! Der aufmerksame Leser fragt sich, wie die Wahl der Gewichtungsmatrizen Qx und Qw zustande kommt. Auch dies ist ganz einfach. Für Qx wird das Gütekriterium der Ausgangsgewichtung ( 13 ) mit ( 14 ) genommen. Die Ausgangsgewichtung Qy wird dabei auf die Identitätsmatrix Ip gesetzt, also auf eine Matrix der Dimension p x p mit den Diagonalelementen 1. Die Ausgangsgewichtung hat den Vorteil, dass wirklich nur beeinfluss-
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
25/35
bare Grössen gewichtet werden (diese werden ja auch gemessen und mit dem Sollwert verglichen). Bei der Zustandsgewichtung kann es vorkommen, dass einzelne Zustandsgrössen zu stark gewichtet werden, welche in der Realität gar nicht im gewünschten Masse beeinflusst werden können (diese werden ja auch nicht geregelt). Analog wird beim Beobachter eine spezielle Wahl der Gewichtungsmatrix Qw getroffen. Vereinfachend wird angenommen, dass die Eingangsmatrix des Prozessrauschens Bw gleich der Eingangsmatrix B sei (siehe Abschnitt 4.9), also: ( 25 )
T
Qw = Bw ⋅ I m ⋅ Bw = B T ⋅ I m ⋅ B
Auch hier kann also wieder mit der Identitätsmatrix Im als Gewichtung gearbeitet werden, also einer Diagonalmatrix der Dimension m x m mit den Elementen 1. Da die Eingangsmatrix B und die Ausgangsmatrix C bekannt sind, ist die Wahl von Qx und Qw trivial. Basierend auf den Parametern ρ und µ wird nun das Vorgehen der LQG/LTRMethode beschrieben: •
1. Schraube = Parameter ρ. Mit diesem Parameter wird das „Loop Shaping“ beschrieben. Der LQ-Regulator gibt also die Form des offenen Loops durch Lr vor, wobei ρ die Bandbreite bestimmt. Dabei werden die guten Robustheitseigenschaften des LQ-Regulators übernommen (garantierte Phasenreserve von 60°). In MATLAB werden dazu folgende Befehle ausgeführt: K = lqr(A,B,Qx,Ru) Zustands-Rückführung K Lr = K*inv(s*eye(n)-A)*B Loop Shape Lr des LQ-Regulators
•
2. Schraube = Parameter µ. Mit diesem Parameter wird der „Loop Transfer Recovery“-Schritt (LTR) vorgenommen, was bedeutet, dass der offene Loop L des Gesamtsystems an das Loop Shaping Lr angeglichen wird. Der Parameter µ bestimmt dabei, wie stark der Loop angeglichen wird. In MATLAB werden dazu folgende Befehle ausgeführt: H = lqe(A,Bw,C,Qw,Rv) Beobachter-Rückführung H Alqg = A-H*C-B*K+H*D*K Synthese LQG-Regulator Blqg = H Clqg = -K Dlqg = 0 Klqg = ss(Alqg,Blqg,Clqg,Dlqg) L = Gu*Klqg offener Loop Gesamtsystem
Wir haben also ein sehr einfaches Vorgehen in zwei Schritten. Obige Reihenfolge gilt allerdings nur für den Fall, dass die Anzahl p der Ausgangsgrössen grösser als die Anzahl m der Stellgrössen ist. Ein solcher Fall ist beispielsweise der Servoantrieb mit einer Stellgrösse (Strom bzw. Kraft/Moment) und mehreren voneinander abhängigen Ausgangsgrössen (Geschwindigkeit, Position). Hier können auch die Sollwerte für die einzelnen Ausgangsgrössen nicht beliebig vorgegeben werden, sondern die Abhängigkeiten müssen berücksichtigt werden (Position = Integrierte Geschwindigkeit). Ein anderes Beispiel ist der Mehrmassenschwinger (Laborversuch!), bei welchem nur eine Stellgrösse, aber mehrere Ausgangsgrössen (Positionen der einzelnen Massen) vorhanden sind. Im Laborversuch wird nur eine Position geregelt, ansonsten müssten auch hier die Abhängigkeiten berücksichtigt werden.
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
26/35
Im umgekehrten Fall, wenn also die Anzahl m der Stellgrössen grösser oder gleich der Anzahl m der Ausgangsgrössen ist, gilt die „duale Methode“. Diese wird für die Mehrheit der Anwendungen verwendet, da im Normalfall die Ausgangsgrössen praktisch unabhängig sind und deshalb mindestens gleich viel Stellgrössen vorhanden sein müssen, um alle Ausgangsgrössen beeinflussen zu können. Die „duale Methode“ wird folgend beschrieben: •
1. Schraube = Parameter µ. Mit diesem Parameter wird das „Loop Shaping“ beschrieben. Der LQ-Estimator gibt also die Form des offenen Loops durch Lf vor, wobei µ die Bandbreite bestimmt. Dabei werden die guten Robustheitseigenschaften des LQ-Estimators übernommen (garantierte Phasenreserve von 60°). In MATLAB werden dazu folgende Befehle ausgeführt: H = lqe(A,Bw,C,Qw,Rv) Beobachter-Rückführung H Lf = C*inv(s*eye(n)-A)*H Loop Shape Lf des LQ-Estimators
•
2. Schraube = Parameter ρ. Mit diesem Parameter wird der „Loop Transfer Recovery“-Schritt (LTR) vorgenommen, was bedeutet, dass der offene Loop L des Gesamtsystems an das Loop Shaping Lf angeglichen wird. Der Parameter ρ bestimmt dabei, wie stark der Loop angeglichen wird. In MATLAB werden dazu folgende Befehle ausgeführt: K = lqr(A,B,Qx,Ru) Zustands-Rückführung K Alqg = A-H*C-B*K+H*D*K Synthese LQG-Regulator Blqg = H Clqg = -K Dlqg = 0 Klqg = ss(Alqg,Blqg,Clqg,Dlqg) L = Gu*Klqg offener Loop Gesamtsystem
Wie bei jeder Methode muss darauf hingewiesen werden, dass der Regler beim Entwurf nicht „übertunt“ werden sollte. Speziell ist auf folgende Punkte zu achten: • •
„Loop Shaping“: Die Bandbreite von L sollte gegenüber der Regelstrecke G nicht zu stark erhöht werden. „Loop Transfer Recovery“: Der Parameter µ bzw. ρ sollte nicht zu gross gewählt werden, weil sonst ein „high gain“-Regler mit hohem „D-Anteil“ entstehen kann, welcher bei verrauschten Signalen nicht mehr realisierbar ist. Es darf also durchaus noch eine Abweichung von Lg und L übrig bleiben.
4.13. Anwendung des LQG/LTR-Entwurfs auf die Gasturbine Wir wollen nun die LQG/LTR-Methode aus Abschnitt 4.12 auf das Gasturbinen-Modell in Abschnitt 4.4 bis 4.7 anwenden. Dazu erstellen wir ein MATLAB-File, welches die einzelnen Schritte des Reglerentwurfes durchführt. Die MATLAB-Befehle sind im Wesentlichen in Abschnitt 4.12 aufgeführt. Das MATLAB-File soll folgende Schritte ausführen: • • •
Anwendung der LQG/LTR-Methode auf das Gasturbinen-Modell (mit 6 Stellgrössen, 12 Zustandsgrössen und 6 Ausgangsgrössen) Darstellung der Singularwertverläufe für das Loop Shaping (Befehl sigmaplot verwenden) Darstellung der Singularwertverläufe des offenen Loops für das Gesamtsystem und Vergleich mit dem Verlauf des Loop Shaping
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
27/35
• •
Darstellung der Singularwertverläufe der Regelstrecke und des resultierenden Reglers (als Kontrolle) Tuning der Parameter µ und ρ, bis ein optimaler, realisierbarer Regler herauskommt
File: RT2_LQGLTR_GT.m Folgend sind ein paar typische Outputs (Singularwertverläufe) von MATLAB gezeigt. Beim Loop Shaping in Abbildung 18 sollte darauf geachtet werden, dass die Durchtrittsfrequenz gegenüber der Regelstrecke (bei ca. 20 rad/s, siehe Abbildung 11) nicht zu stark erhöht wird.
Singular Values 50 Lf 40 30
Singular Values (dB)
20 10 0 -10 -20 -30 -40 -50 -2 10
10
-1
0
10
10
1
2
10
10
3
4
10
Frequency (rad/s)
Abbildung 18: Singularwertverläufe Loop Shaping Lf (Durchtrittsfrequenz bei ca. 50 rad/s)
Im LTR-Schritt muss für eine gute Auslegung darauf geachtet werden, dass der GesamtLoop L nicht zu nahe an den Shaping Loop Lf herangedrückt wird, was Abbildung 19 zeigt. Der resultierende Regler ist in Abbildung 20 dargestellt (rote Kurve Klqg). Hier ist ein PD-Verhalten mit beschränktem D-Anteil ersichtlich (Amplitudenabfall bei höheren Frequenzen, also PD-T2). Im Gegensatz dazu ist auch eine schlechte Auslegung mit „übertuntem“ LTR-Parameter dargestellt. In Abbildung 21 ist ersichtlich, wie der Gesamt-Loop L künstlich auf den Shaping Loop Lf gedrückt wird. Abbildung 22 zeigt den resultierenden Regler (rote Kurve Klqg), welcher bei hohen Frequenzen noch eine grosse Verstärkung zeigt („high gain“). Bei verrauschten Signalen kann dies zu Problemen führen, wie die Simulation in folgendem Abschnitt zeigen wird. Die Regelstrecke Gu kann durch den Regler Klqg also nicht auf eine beliebige Vorgabe „gedrückt“ werden! Natürlich ist dieses Verhalten auch abhängig von den Eigenschaften der Regelstrecke Gu, hier haben wir es immerhin mit einem gekoppelten System 12. Ordnung zu tun!
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
28/35
Singular Values 60 Lf 40
L
20
Singular Values (dB)
0 -20 -40 -60 -80 -100 -120 -140 -2 10
10
-1
0
10
10
1
2
10
10
3
4
10
Frequency (rad/s)
Abbildung 19: Loop Shaping Lf und Loop Transfer Recovery L (gute Auslegung)
Singular Values 40 Gu 20
L Klqg
0
Singular Values (dB)
-20 -40 -60 -80 -100 -120 -140 -2 10
10
-1
0
10
10
1
2
10
10
3
4
10
Frequency (rad/s)
Abbildung 20: Regelstrecke Gu, Loop Transfer Recovery L, Regler Klqg (gute Auslegung)
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
29/35
Singular Values 60 Lf 40
L
20
Singular Values (dB)
0 -20 -40 -60 -80 -100 -120 -2 10
10
-1
0
10
10
1
2
10
10
3
4
10
Frequency (rad/s)
Abbildung 21: Loop Shaping Lf und Loop Transfer Recovery L („übertunt“)
Singular Values 60 Gu L
40
Klqg 20
Singular Values (dB)
0 -20 -40 -60 -80 -100 -120 -2 10
10
-1
0
10
10
1
2
10
10
3
4
10
Frequency (rad/s)
Abbildung 22: Regelstrecke Gu, Loop Transfer Recovery L, Regler Klqg („übertunt“)
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
30/35
4.14. Simulation des entworfenen Reglers unter realen Bedingungen Wichtig ist, dass der entworfene Regler unter möglichst „realen“ Bedingungen ausgetestet wird, bevor er implementiert wird. Dies ist unabhängig von der Entwurfsmethode! Hier wird eine Simulation in MATLAB/Simulink durchgeführt, um reale Effekte wie Rauschen auf den Messsignalen oder Stellgrössenbeschränkungen zu berücksichtigen. Dazu wird folgendes Simulink-Modell aufgebaut (File: RT2_Regelkreis_GT.mdl):
Abbildung 23: Simulink-Modell zum Austesten des Reglers mit Rauschen am Ausgang, Stellgrössenbeschränkung und zusätzlichen Störgrössen am Eingang
Als Regelstrecke G wird in Abbildung 23 das Modell aus Abschnitt 4.6 (Serieschaltung der Teilmodelle) verwendet. Es hat 12 Eingänge, wobei die ersten 6 Eingänge die Stellgrössen sind und die weiteren 6 Eingänge die Störgrössen. Dies wird genutzt, um zusätzliche Störgrössen am Eingang auf das System zu geben. Zudem kann ein Rauschen mit einstellbarer Varianz am Ausgang des Systems überlagert werden, wozu der Block „Random Number“ verwendet wird. Die Stellgrössenbeschränkung am Eingang ist durch einen Block „Saturation“ realisiert. Schliesslich können verschiedene Sollwertvorgaben gefahren werden. Der Regler ist wie die Strecke mit einem LTI-Block definiert, wobei die Übertragungsmatrix –Klqg aus Abschnitt 4.13 verwendet wird (Achtung negatives Vorzeichen!). Die Ein- und Ausgangsgrössen sind skaliert, um physikalisch einigermassen sinnvolle Werte zu haben (Stellgrössen sind Brennstoffmassenströme in kg/s, Ausgangsgrössen sind Temperaturen in °C). Folgend sind ein paar typische Simulationen gezeigt. Zunächst wird eine Sollwertvorgabe gemacht, es wird also das Verhalten als Folgeregler untersucht. Abbildung 24 zeigt, wie der Istwert in Kammer 1 (y1) schön dem Sollwert folgt, die Temperaturen der anderen Kammern (y2..y6) bleiben auf dem ursprünglichen Wert, trotz Querkopplungen im System. Der Regler hat also optimal entkoppelt! Dies ist auch klar in den Stellgrössenverläufen auf Abbildung 25 ersichtlich. Es reagiert nicht nur eine Stellgrösse, sondern auch die Stellgrössen der benachbarten Kammern. Es ist auch ersichtlich, dass Stellgrösse u2 am stärksten reagiert statt Stellgrösse u1. Wir erinnern uns an die Ablenkung in der Turbine, durch die „Shift-Matrix“ in Abschnitt 4.5 beschrieben. Auch dieses Verhalten wird im Regler automatisch übernommen, da eine Kopie des Modells im Beobachter abgebildet ist. Der Regler steuert also automatisch die richtige Brennkammer an!
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
31/35
Abbildung 24: Verlauf der Ausgangsgrössen y1..y6 (Temperaturen) auf einen Sollwertsprung (600°C in Kammer 1 bei 0..2 sec)
Abbildung 25: Verlauf der Stellgrössen u1..u6 (Brennstoffmassenströme) beim Sollwertsprung
In Abbildung 24 ist aber auch eine kleine bleibende Sollwertabweichung ersichtlich, 600 °C wird nicht ganz erreicht. Warum dies? Wir betrachten nochmals den Regler in Abbildung 20 (rote Kurve) – dies ist ein PD-Regler, der I-Anteil fehlt also! Mit einem kleinen Trick könnte auch dieses Problem noch gelöst werden. Man müsste das System mit einem künstlichen I-Anteil (auf allen Kanälen) erweitern. Für die Gasturbine spielt dies aber keine Rolle und wir lassen die Erweiterung weg. Nun wollen wir austesten, was bei einem „übertunten“ Regler geschieht, gemäss Abbildung 22. Da wir Rauschen auf dem Messsignal haben, sollte dies nun deutlich verstärkt werden. Und tatsächlich, Abbildung 26 rechts zeigt eine stark verrauschte Stellgrösse. Diese Reglereinstellung ist also nicht brauchbar! Natürlich könnte man in der Realität ein zusätzliches Filter einbauen, um die Messdaten zu filtern, bevor sie in den Regler kommen. Aber ein solches Filter würde durch dessen Zeitkonstante die Performance des Systems ebenfalls verschlechtern. Somit können wir auch den „trägeren“ Regler verwenden. Nebenbei sehen wir auch, wie der Sollwert in Abbildung 26 links etwas besser erreicht wird. Dies ist auf die grössere Verstärkung des Reglers zurückzuführen.
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
32/35
Abbildung 26: Verlauf von Ausgangsgrössen y1..y6 (links) und Stellgrössen u1..u6 (rechts) bei „übertuntem“ Regler
Als letzte Simulation lassen wir nun die Störgrössen am Eingang „spielen“. Diese sollen zufällige Temperaturvariationen am Eingang der Brennkammer simulieren. Die Sollwerte werden konstant gehalten, der Regler hat also die Aufgabe der reinen Störunterdrückung. Der Regler erfüllt diese Aufgabe optimal, wie folgende Abbildungen zeigen (PlotSkalierungen vergrössert).
Abbildung 27: „Zufällig“ variierende Störgrössen d1..d6 (Temperaturen ändern jede Sekunde)
Abbildung 28: Resultierende Ausgangstemperaturen y1..y6 (Betriebspunkt 500 °C)
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
33/35
Abbildung 29: Stellgrössen u1..u6 (Betriebspunkt 200 kg/s)
Damit wurde gezeigt, dass die Auslegung nach LQG/LTR optimale Regler sowohl für die Folgeregelung wie auch die Störunterdrückung liefert. Wichtig ist – wie bei allen Einstellmethoden – dass der Regler nicht „übertunt“ wird und unter möglichst realen Bedingungen ausgetestet wird.
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
02.05.2012
34/35
4.15. Zusammenfassung Teil 10 (Zustandsraum-Darstellung) Gelerntes Klassische vs. modellbasierte Regelung Klassisch: Kopplungen im System werden nicht berücksichtigt, Black-Box-Ansatz Modellbasiert (Prozessmodell als Basis): White-Box-Ansatz Zustandsraum-Darstellung MIMO-Systeme (Multiple Input Multiple Output) Kopplungen und Dynamik bei Modellierung berücksichtigt Serieschaltung/Zusammensetzen von Teilmodellen in MATLAB möglich Analyse von MIMO-Systemen Singularwertverläufe treten an die Stelle des Bode-Diagrammes Analog zum Amplitudengang, aber auf Mehrgrössensysteme erweitert
Lernziel erreicht? ☺☺☺
☺☺☺
☺☺
Teil 11 (Auslegungsmethoden) Gelerntes Zustandsregler nach LQR-Methode auslegen können „Performance“ vs. „Energieverbrauch“ Optimalen Beobachter nach LQE-Methode auslegen können „Minimaler Schätzfehler“ vs. „Rauschunterdrückung“ Mehrgrössen-Folgeregelung mit Beobachter und Zustandsregler verstehen Wie klassischer Regelkreis Kopie der Regelstrecke im Regler! Alle Kopplungen usw. automatisch berücksichtigt! Kombinierte Auslegungsmethode LQR/LTR anwenden können Nur noch 2 Auslegungsparameter! Schrittweises Vorgehen mit Loop Shaping und Loop Transfer Recovery Automatisch robuste Regler! Achtung „Übertuning“! Immer mit Rauschanteilen und Stellgrössenbeschränkungen simulieren!
RT2, Teil 10/11, Zustandsraum-Darstellung für Mehrgrössen-Systeme – Prof. Dr. D. Zogg
Lernziel erreicht? ☺☺☺ ☺☺☺ ☺☺☺
☺☺
02.05.2012
35/35