MPC Toolbox in Matlab -...

37
MPC Toolbox in Matlab Seminararbeit von Florian H¨ aberlein FAKULT ¨ A TF ¨ UR MATHEMATIK UND PHYSIK MATHEMATISCHES INSTITUT Datum: 4. Juli 2006 Betreuung: Prof. Dr. L. Gr¨ une Dipl.-Wirtschaftsmath. J. Pannek

Transcript of MPC Toolbox in Matlab -...

Page 1: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

MPC Toolbox in Matlab

Seminararbeit

von

Florian Haberlein

FAKULTAT FUR MATHEMATIK UND PHYSIK

MATHEMATISCHES INSTITUT

Datum: 4. Juli 2006 Betreuung:Prof. Dr. L. GruneDipl.-Wirtschaftsmath. J. Pannek

Page 2: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss
Page 3: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

Inhaltsverzeichnis

Tabellenverzeichnis III

Abbildungsverzeichnis V

1 MPC Toolbox in Matlab 1

1.1 Einfuhrung in die MPC Toolbox . . . . . . . . . . . . . . . . . . . 11.1.1 Die Rolle der MPC Toolbox in der Matlab-Familie . . . . 11.1.2 Grundlegende Fahigkeiten der MPC Toolbox . . . . . . . . 2

1.2 Beispiel: Leerlauffullungsregelung . . . . . . . . . . . . . . . . . . 41.2.1 Problemstellung . . . . . . . . . . . . . . . . . . . . . . . . 41.2.2 Problemmodellierung . . . . . . . . . . . . . . . . . . . . . 41.2.3 Implementierung des Systems . . . . . . . . . . . . . . . . 51.2.4 Implementierung des Reglers . . . . . . . . . . . . . . . . . 61.2.5 Simulation des closed-loop-Verhaltens . . . . . . . . . . . . 91.2.6 Vergleich verschiedener Regler . . . . . . . . . . . . . . . . 111.2.7 Import des Reglers in Simulink . . . . . . . . . . . . . . . 12

1.3 Der Blick hinter die Kulissen der MPC Toolbox . . . . . . . . . . 141.3.1 Algorithmen fur unbeschrankte Probleme oder Probleme

mit saturierten Beschrankungen . . . . . . . . . . . . . . . 141.3.2 Algorithmen fur beschrankte Probleme . . . . . . . . . . . 15

A Plots und Quelltext 19

A.1 Plots . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19A.2 Quelltext zu Kapitel 1.2 . . . . . . . . . . . . . . . . . . . . . . . 26

Literaturverzeichnis 29

I

Page 4: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

II INHALTSVERZEICHNIS

Page 5: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

Tabellenverzeichnis

1.1 Ubersicht uber verschiedene Storungsarten . . . . . . . . . . . . . 3

III

Page 6: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

IV TABELLENVERZEICHNIS

Page 7: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

Abbildungsverzeichnis

1.1 Signalflussplan des cmpc-Reglers . . . . . . . . . . . . . . . . . . 71.2 Simulation fur die Fahrstufe

”D“ mit Modellfehler . . . . . . . . . 10

1.3 Simulation fur die Fahrstufe”N“ mit Modellfehler . . . . . . . . . 11

1.4 Vergleich von Simulationen fur die Fahrstufe”D“ mit beschranktem

Pradiktionshorizont . . . . . . . . . . . . . . . . . . . . . . . . . . 121.5 Signalflussplan des linearen Kontrollsystems in Simulink . . . . . 131.6 Signalflussplan der nichtlinearen Regelstrecke mit zusatzlicher Sto-

rung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 141.7 Prinzip der saturierten Beschrankungen . . . . . . . . . . . . . . . 151.8 Dantzig-Wolfe-Methode der aktiven Mengen . . . . . . . . . . . . 171.9 Beispiel zur Dantzig-Wolfe-Methode der aktiven Mengen . . . . . 18

A.1 Simulation fur die Fahrstufe”D“ ohne Modellfehler . . . . . . . . 19

A.2 Simulation fur die Fahrstufe”N“ ohne Modellfehler . . . . . . . . 20

A.3 Simulation fur die Fahrstufe”D“ mit Modellfehler . . . . . . . . . 20

A.4 Simulation fur die Fahrstufe”N“ mit Modellfehler . . . . . . . . . 21

A.5 Simulation fur die Fahrstufe”D“ ohne Modellfehler mit erhohten

Gewichtsfaktoren . . . . . . . . . . . . . . . . . . . . . . . . . . . 21A.6 Simulation fur die Fahrstufe

”N“ mit Modellfehler mit erhohten

Gewichtsfaktoren . . . . . . . . . . . . . . . . . . . . . . . . . . . 22A.7 Vergleich von Simulationen fur die Fahrstufe

”D“ mit unbeschrank-

tem Pradiktionshorizont . . . . . . . . . . . . . . . . . . . . . . . 23A.8 Vergleich von Simulationen fur die Fahrstufe

”D“ mit beschranktem

Pradiktionshorizont . . . . . . . . . . . . . . . . . . . . . . . . . . 24A.9 Simulation des linearen cmpc-Reglers an der nichtlinearen Regel-

strecke, Ausgangsvariablen . . . . . . . . . . . . . . . . . . . . . . 25A.10 Simulation des linearen cmpc-Reglers an der nichtlinearen Regel-

strecke, Kontrollvariablen . . . . . . . . . . . . . . . . . . . . . . 25

V

Page 8: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

VI ABBILDUNGSVERZEICHNIS

Page 9: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

Kapitel 1

MPC Toolbox in Matlab

In diesem Kapitel wird die MPC (Model Predictive Control) Toolbox in Matlabvorgestellt. Im ersten Abschnitt wird eine allgemeine Einfuhrung gegeben. Imzweiten Abschnitt erfolgt eine ausfuhrliche Behandlung eines Beispiels. Im drit-ten Abschnitt wird auf die verwendeten Algorithmen der in der MPC Toolboxvorhandenen Funktionen eingegangen.

1.1 Einfuhrung in die MPC Toolbox

In diesem Abschnitt wird die MPC Toolbox allgemein vorgestellt. Zuerst wird dieEinordnung der Toolbox in die Matlab-Produktfamilie erlautert, im Anschlusswerden die speziellen Fahigkeiten der MPC Toolbox aufgelistet.

1.1.1 Die Rolle der MPC Toolbox in der Matlab-Familie

Das Programmpaket Matlab bietet neben den im Kernprogramm vorhandenennumerischen Standardalgorithmen eine Reihe von Erweiterungspaketen — sog.Toolboxen — mit numerischen Algorithmen zu speziellen Themenbereichen derMathematik, Physik und den Ingenieurwissenschaften. Der in allen drei Wis-senschaften behandelte weitlaufige Bereich der Kontrolltheorie wird dabei unteranderem durch Simulink, die Control System Toolbox und die MPC Toolboxabgedeckt. Diese Toolboxen beschaftigen sich mit der Modellbildung, dem Reg-lerentwurf und den theoretischen Tests der entwickelten Regelstrecken. WeitereToolboxen ermoglichen dann den Entwurf von Prototypen des Reglers (z. B.mit Real-Time Workshop) und bieten Schnittstellen zwischen den realen Kon-trollsystemen und den entwickelten Reglern (z. B. die OPC Toolbox, OpenessProductivity Collaboration).Fur die Arbeit mit der MPC Toolbox in Matlab spielt die Control System Tool-box eine fundamentale Rolle. Mit ihr ist die Erstellung von Kontrollsystemen inden verschiedenen Modellarten moglich. Desweiteren bietet sie Algorithmen zur

1

Page 10: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

2 KAPITEL 1. MPC TOOLBOX IN MATLAB

Transformation zwischen den verschiedenen Modellen, zu Stabilitatstests der je-weiligen Modelle sowie zu grundlegenden Problemen der Kontrolltheorie. Da dieControl System Toolbox die existentielle Grundlage der MPC Toolbox bildet, istLetztere nur in Zusammenhang mit Ersterer zu betreiben.Im bereits erwahnten Programm Simulink finden sich eine Reihe von Moglich-keiten zur Modellierung von linearen und nichtlinearen Kontrollsystemen. Die indiesem Programm erstellten Modelle lassen sich in linearisierter Form ebenfallsin der MPC Toolbox zum Reglerentwurf heranziehen. Umgekehrt konnten die inder MPC Toolbox fur lineare Systeme entwickelten modellpradiktiven Regler furdie nichtlinearen Systeme in Simulink verwendet werden. Desweiteren bietet dieSystem Identification Toolbox die Moglichkeit, durch Tests an realen Systemenapproximierte Modelle am PC zu bestimmen, die wiederum sowohl in Simulinkals auch in der MPC Toolbox verwendet werden konnen.Zusammenfassend kann man sagen, dass die MPC Toolbox einen machtigen Bau-stein in der Matlab-Familie bildet. Die optimale Zusammenarbeit zwischen denverschiedenen Toolboxen und dem Kernprogramm sowie die Objektorientierungder entwickelten modellpradiktiven Regler begrunden diese Schlussfolgerung.

1.1.2 Grundlegende Fahigkeiten der MPC Toolbox

Die MPC Toolbox bietet Algorithmen zum Entwurf, zur Analyse und zur Si-mulation von modellpradiktiven Reglern. Als Benutzeroberflachen bietet sie zumEinen den sogenannten Control and Estimation Tools Manager, mit dem mangraphisch Regelstrecken und Regler erstellen, analysieren und testen kann. DieBenutzung ist intuitiv und ubersichtlich, jedoch ist diese Benutzerschnittstelle imLeistungsumfang gegenuber der herkommlichen kommandozeilenorientierten Be-nutzung eingeschrankt. Eine bessere Modularisierung und Wartbarkeit der Reg-ler und Kontrollsysteme und die Wiederverwendbarkeit von Teilsystemen undReglern, sowie die Benutzung von Scriptfiles bietet nur die klassische komman-dozeilenorientierte Benutzeroberflache. Es stehen alle im Control and EstimationTools Manager vorhandenen Funktionen auch durch kommandozeilenorientier-te Funktionen zur Verfugung, weswegen diese Arbeit sich nicht naher mit dergraphischen Oberflache beschaftigt. Eine detaillierte Anleitung zur graphischenOberflache findet sich in der Bedienungsanleitung zur Version 2 der MPC Tool-box [1, Kapitel 3].Das Anlegen der Modelle von Kontrollsystemen lasst sich neben der oben be-schriebenen Moglichkeit des Imports aus anderen Matlab-Toolboxen naturlichauch direkt durch Eingabe der Modelle mittels Kommandozeilenfunktionen be-werkstelligen. Dabei kann die MPC Toolbox MIMO-Systeme (Multi Input, Mul-ti Output, Kontroll- und Ausgangsvariablen konnen mehrdimensional sein) mitBeschrankungen in Kontroll- und Ausgangsvariablen verarbeiten. Grundsatzlichkonnen jedoch nur lineare, zeitinvariante Kontrollsysteme (sog. LTI-Systeme, li-

Page 11: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.1. EINFUHRUNG IN DIE MPC TOOLBOX 3

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5−0.5

0

0.5

1

1.5

2

t

z(t)

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5−0.2

0

0.2

0.4

0.6

0.8

1

t

z(t)

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 50

1

2

3

4

5

6

t

z(t)

z(t) = a z(t) =

(

0 t < t0

1 t ≥ t0z(t) = at + b

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5−0.8

−0.6

−0.4

−0.2

0

0.2

0.4

0.6

0.8

t

z(t)

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5−2.5

−2

−1.5

−1

−0.5

0

0.5

1

1.5

2

t

z(t)

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5−0.5

0

0.5

1

1.5

t

z(t)

z(t) = a sin(bt + t0)

z(t) = an,

t ∈ [nts, (n + 1)ts[,(an) ∼ N (0, 1)

z(t) =

(

a t ∈ [nts, (n + 1)ts

b t ∈ [(n + 1)ts , (n + 2)ts[

n := 2i

Tabelle 1.1: Ubersicht uber verschiedene Storungsarten. Obere Zeile v. l. n. r.:konstant, Sprung, linear; Untere Zeile v. l. n. r.: sinusformig, gaußverteilt, gepulst.a, b ∈ R, ts, t0 ∈ R

+0 , n, i ∈ N

+0

near, time invariant) der Form

x(k + 1) = Ax(k) + Bu(k) + Md(k) + Nw(k)

y(k) = y(k) + z(k)

verarbeitet werden. Dabei bezeichnet x(k) den k-ten Ausgang, u(k) die k-te Kon-trollanweisung, d(k) bzw. w(k) die gemessene bzw. ungemessene Storung; y(k)bezeichnet den Ausgang des Systems, z(k) das Messrauschen und y(k) den ge-messenen Ausgang. Die Matrizen A, B, M und N sind hierbei konstant, alsozeitinvariant. Desweiteren mussen die Systeme in diskretisierter Form vorliegen,die Abtastzeit kann jedoch selbst festgelegt werden. Neben dem Einfluss vonbekannten Storungen kann man auch unbekannte Storungen auf das System ein-wirken lassen. Dabei kann die Storung konstant, als Sprung, gepulst, als lineareFunktion, sinusformig oder gaußverteilt gewahlt werden. Eine Ubersicht uber dieverschiedenen Storungsarten gibt Tabelle 1.1. Bei der Reglererstellung wird dasKontrollziel als Setpoint-Regelung implementiert, eine Referenztrajektorie wirdals Folge von Setpoints zu den Abtastzeitpunkten realisiert. Als Reglerparame-ter lassen sich innerhalb der Problemstellung Pradiktions- und Kontrollhorizont

Page 12: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

4 KAPITEL 1. MPC TOOLBOX IN MATLAB

sowie Gewichtsfaktoren fur das Regelziel bei Ein- und Ausgangsgroßen festlegen.Die erstellten Regler lassen sich im Anschluss exportieren, so dass sie z. B. inSimulink fur die zugehorige nichtlineare Regelstrecke getestet werden konnen.

1.2 Beispiel: Leerlauffullungsregelung

In diesem Abschnitt wird die sog. Leerlauffullungsregelung (LFR) bei Motoren inAutos mit Automatikgetrieben mittels MPC behandelt. Der algorithmische Teildieses Abschnitts orientiert sich dabei an der Bedienungsanleitung zur Version 1der MPC Toolbox [8, Kapitel 2, Seite 2-22 bis 2-30], eine ausfuhrliche theoretischeBehandlung des Themas findet sich in der Arbeit von Hrovat und Bodenheimer [6,Seite 1778 - 1783]. Weiterfuhrende Literatur zum Einsatz von modellgestutztenReglern in Verbrennungsmotorsystemen findet sich im Buch von Guzzella undOnder [4].

1.2.1 Problemstellung

Das Ziel der Leerlauffullungsregelung ist es, einen Motor, der im Leerlauf be-trieben wird, durch Einfluss des Zundzeitpunkts und der Menge an zugefuhrterLuft wahrend des Zundzeitpunkts zu einer konstanten Mindestdrehzahl zu regeln.Das Kontrollsystem unterliegt dabei unbekannten Storungen: Ein Lenkeinschlagim Stand fuhrt beispielsweise dazu, dass der Oldruck der Servolenkung erhohtwerden muss, was wiederum zu einer plotzlichen Last auf den Motor fuhrt, eben-so konnen plotzliche Lastanderungen am Motor durch die Klimaanlage oder dieLichtmaschine hervorgerufen werden. Aufgrund unterschiedlicher Systemkompo-nenten und -eigenschaften (verschiedene Motorvarianten, verschiedene zusatzlichean den Motor gebundene Gerate) soll ein Regler fur die Massenproduktion erstelltwerden, der robust und stabil die auftretenden Storungen ausregelt. Der modell-pradiktive Ansatz bietet sich deswegen hier besonders an. Unabhangig von denveranderlichen Systemparametern mussen jedoch gewisse Schranken in den Ein-und Ausgangsvariablen eingehalten werden: Um einen sicheren Motorbetrieb zugewahrleisten, darf der Zundzeitpunkt vom voreingestellten Zundzeitpunkt umnicht mehr als 20 ◦ abweichen, desweiteren sollte er nach einer Veranderung zumvoreingestellten Zundzeitpunkt zuruckgeregelt werden um einen optimalen Kraft-stoffverbrauch zu gewahrleisten.

1.2.2 Problemmodellierung

Aufgrund der betrachteten Problemstellung liegt es nahe, den Zundzeitpunkt so-wohl als veranderbare Eingangsvariable als auch als zu regelnde Ausgangsvariablemit einer Gleichgewichtsregelung zu behandeln. Das allgemeine lineare Kontroll-

Page 13: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.2. BEISPIEL: LEERLAUFFULLUNGSREGELUNG 5

system wird dabei wie folgt im Frequenzbereich modelliert:

(

y1

y2

)

=

(

G1 G2

0 1

)(

u1

u2

)

+

(

Gd

0

)

w (1.1)

Dabei beschreibt y1 die Drehzahl des Motors, y2 und u2 beschreiben beide denZundzeitpunkt, u1 beschreibt die Menge der zugefuhrten Luft und w stellt dieunbekannte Storung in Form eines Lastdrehmoments dar.Fur das Modell betrachten wir zwei Betriebszustande, zum Einen den Betriebin der Fahrstufe

”D“ (dies entspricht bei Schaltgetrieben dem Betrieb mit einge-

legtem ersten Gang und einem”Schleifenlassen“ der Kupplung, so dass auf den

Motor ein konstantes Drehmoment von 30 Nm wirkt, welches das Auto langsamvorantreibt), zum Anderen den Betrieb in der Fahrstufe

”N“ (dies entspricht bei

Schaltgetrieben dem Leerlauf ohne eingelegten Gang, das Lastdrehmoment ist indiesem Fall 0 Nm). In beiden Fallen soll die Motordrehzahl auf konstant 800 rpmgehalten werden.Modelliert man fur das System die Ubertragungsfunktionen so, dass die Gleichge-wichtszustande im Nullpunkt liegen, ergibt sich als Beschrankung fur den Zund-zeitpunkt |u2| ≤ 0,7. Fur den Betrieb in der Fahrstufe

”D“ ergeben sich die

Ubertragungsfunktionen wie folgt:

G1 =9,62e−0,16s

s2 + 2,4s + 5,05

G2 =15,9(s + 3)e−0,04s

s2 + 2,4s + 5,05

Gd =−19,1(s + 3)

s2 + 2,4s + 5,05

(1.2)

Fur den Betrieb in der Fahrstufe”N“ ergeben sich folgende Ubertragungsfunktio-

nen:

G1 =20,5e−0,16s

s2 + 2,2s + 12,8

G2 =47,6(s + 3,5)e−0,04s

s2 + 2,2s + 12,8

Gd =−19,1(s + 3,5)

s2 + 2,2s + 12,8

(1.3)

1.2.3 Implementierung des Systems

Nun werden die Ubertragungsfunktionen fur (1.2) und (1.3) angelegt:

G1D = [0 0 9.62;1 2.4 5.05;0 0.16 0];

G2D = [0 15.9 15.9*3;1 2.4 5.05;0 0.04 0];

Page 14: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

6 KAPITEL 1. MPC TOOLBOX IN MATLAB

G1N = [0 0 20.5; 1 2.2 12.8; 0 0.16 0];

G2N = [0 47.6 47.6*3.5; 1 2.2 12.8; 0 0.04 0];

G0 = [0 0;0 0;0 0];

GI = [0 1;0 1;0 0];

Diese Ubertragungsfunktionen sind fur kontinuierliche Frequenzen formuliert (manbeachte, dass das System im Frequenz- und nicht im Zeitbereich formuliert ist),die MPC Toolbox kann jedoch nur diskrete Systeme verarbeiten. Deswegen wer-den nun die beiden Kontrollsysteme gemaß (1.1) mit einer diskreten Approxi-mation der Sprungantworten weiterverwendet. Dabei benutzen wir eine Endzeitder Approximation von tf = 4 und eine Abtastrate von ts = 0,1, dies bedeutet,dass wir tf/ts = 40 Koeffizienten der Sprungfunktion approximieren. Zur Fest-legung der Anzahl an stabilen Ausgabevariablen geben wir als drittes Argumentder Funktion tfd2step 2 an. Fur Systeme mit integrierenden Ausgabevariablen istdie dritte Ubergabevariable ein Vektor, mit der Dimension des Ausgabevektors,eine 1 in der i-ten Komponente signalisiert dabei einen stabilen i-ten Ausgang,eine 0 einen integrierenden i-ten Ausgang. Die restlichen vier Ubergabevariablenstellen die Transferfunktionen des ungestorten Systems dar, wobei die Elementeder Matrix spaltenweise von oben nach unten angegeben werden.

tf = 4;

ts = 0.1;

systemD = tfd2step(tf,ts,2,G1D,G0,G2D,GI);

systemN = tfd2step(tf,ts,2,G1N,G0,G2N,GI);

Um die Kontrollsysteme zu vervollstandigen, implementiert man nun noch dieStorungen fur beide Szenarien, dabei geht man genauso vor wie oben:

GdD = [0 -19.1 -19.1*3; 1 2.4 5.05; 0 0 0];

dD = tfd2step(tf,ts,2,GdD,G0);

GdN = [0 -19.1 -19.1*3.5; 1 2.2 12.8; 0 0 0];

dN = tfd2step(tf,ts,2,GdN,G0);

Die Vorarbeit ist nun geleistet und man kann zum Reglerentwurf und zu closed-loop-Simulationen ubergehen.

1.2.4 Implementierung des Reglers

Wir benutzen hierfur den cmpc-Regler (constraint MPC), der sich besonders furSimulationen an closed-loop Systemen mit relativ starken Beschrankungen anEin- und/oder Ausgangsvariablen eignet. Auf die implementierten Algorithmendieses und anderer MPC-Regler wird in Kapitel 1.3 eingegangen. Die Zusammen-hange dieses Reglers verdeutlicht Abbildung 1.1.Zunachst legen wir die Reglerparameter fest:

Page 15: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.2. BEISPIEL: LEERLAUFFULLUNGSREGELUNG 7

Abbildung 1.1: Signalflussplan des cmpc-Reglers

• Gewichtsfaktoren fur Ausgange: In der Grundeinstellung sind alle zu re-gelnden Ausgange gleich bewertet, dies bedeutet, Abweichungen von denSetpoints werden fur alle Ausgange gleich

”bestraft“. Im vorliegenden Mo-

dell ist die Drehzahlregelung viel wichtiger, als die Zundzeitpunktregelung,letztere spielt eine untergeordnete Rolle. Wir gewichten also y1 5-fach undy2 1-fach.

ywt = [5 1];

• Gewichtsfaktoren fur Regelung: In der Grundeinstellung sind alle Regelein-griffe unbestraft in der Optimierungsformel. In unserem Fall bestrafen wirzunachst beide 0,5-fach.

uwt = [.5 .5];

• Kontrollhorizont: Die Standardeinstellung fur den Kontrollhorizont ist M =1, dies bedeutet, es wird jeweils nur der erste Zeitschritt des aktuellen Ho-rizonts mit der gerade berechneten Kontrollfolge geregelt. Wir definierenunseren Kontrollhorizont mit N = 10, das heißt, es wird bei einer Abta-strate von ts = 0,1 immer genau eine Sekunde mit der aktuell berechnetenSteuerfolge geregelt. Es besteht die Moglichkeit, sogenannte

”Blockierfak-

toren“ zu definieren, diese geben die Anzahl der Zeitschritte an, uber diedie Kontrollfolge unverandert bleibt, also ∆u = 0. Dafur definiert man Mals Zeilenvektor, wobei die Summe uber alle

”Blockierfaktoren“ maximal so

groß sein darf, wie der Pradiktionshorizont.

Page 16: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

8 KAPITEL 1. MPC TOOLBOX IN MATLAB

m = 10;

• Pradiktionshorizont: Wir definieren den Pradiktionshorizont mit p = ∞.

p = inf;

• Simulationszeit: Die Endzeit der Simulation definieren wir mit tfsim = 20.

tfsim = 20;

• Referenztrajektorie: Die Standardeinstellung fur die Referenztrajektorie isteine Nullregelung fur alle Komponenten. Fur eine Regelung mit Setpointsverschieden von 0 definiert man die Referenztrajektorie als Spaltenvektor,wobei in der i-ten Komponente der Setpoint fur die i-te Ausgangsvariablesteht. Mochte man beispielsweise einer Referenztrajektorie folgen, so defi-niert man die Referenztrajektorie als Matrix. Dabei bezeichnet das Elementin der i-ten Zeile und der j-ten Spalte, den Setpoint fur die i-te Komponentezur Zeit t = jts. Hat die Matrix j Spalten mit jts < tfsim, so gilt ab der Zeitt = jts die letzte Spalte als konstante Setpointdefinition fur den Zeitbe-reich [jts, tfsim]. In unserem Fall wunschen wir eine konstante Nullregelung(beachte, das System wurde dementsprechend normiert).

r = [0 0];

• Beschrankungen der Kontrolle: Um den Kontrolleingriff zu beschranken de-finiert man einen Spaltenvektor, bestehend aus drei Matrizen nach dem sel-ben Prinzip wie bei der Referenztrajektorie. Die erste Matrix beschreibt dieunteren Schranken fur u, die zweite Matrix beschreibt die oberen Schrankenan u und die dritte Matrix beschreibt die Beschrankungen an ∆u. Die Stan-dardeinstellungen sind umin = −∞, umax = ∞ sowie (∆u)max = 106 fur alleZeitpunkte und alle Komponenten. Fur unser System ist u1 (Luftmenge) un-beschrankt sowohl in der oberen und unteren Grenze, als auch in der maxi-malen Anderung. Der Zundzeitpunkt ist beschrankt durch −0,7 ≤ u2 ≤ 0,7.Um die Motorfunktionalitat nicht zu beeintrachtigen, beschranken wir diemaximale Anderung des Zundzeitpunkts auf (∆u2)max ≤ 1,4.

ulim = [[-inf -0.7] [inf 0.7] [1e6 1.4]];

• Beschrankungen des Ausgangs: Die Definition der Beschrankungen des Aus-gangs entspricht der der Kontrolle, lediglich der Beschrankung an die An-derung der Ausgangsvariablen (∆y)max fallt weg. In unserem Fall liegenkeine Beschrankungen des Ausgangs vor, die Anderung des Zundzeitpunktswurde bereits durch u2 beschrankt.

Page 17: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.2. BEISPIEL: LEERLAUFFULLUNGSREGELUNG 9

ylim = [];

• Filterung und Verzogerung fur Storungen: Fur die Definition eines Entstor-filters und fur die Definition von Verzogerungen fur die Storungen definiertman je eine Matrix mit Zeitkonstanten. Die Standardeinstellungen sind kei-ne Filterung und keine Verzogerung (also ein unverzogerter Einheitssprungin der Storung). Wir definieren keine Filterung und keine Verzogerung furunser Modell.

tfilter = [];

• Storungen: Fur die Storungen ergibt sich das selbe Format wie fur die Re-ferenztrajektorie. Fur reine Storungen an den Ausgabevariablen (kein Sto-rungssystem und kein Storungsmodell vorhanden), stehen in der Matrix dieadditiven Storungen. Ist ein Storungssystem oder sowohl ein Storungsmo-dell als auch ein Storungssystem vorhanden, so werden in der Matrix dieSprungfunktionen definiert. In unserem Fall haben wir eine ungemesseneStorung, die durch eine Ubertragungsfunktion in das System eingeht (alsonur ein Storungssystem ohne Storungsmodell). Das zugehorige Storungssys-tem haben wir bereits angelegt. Desweiteren benotigen wir nun noch einenEinheitssprung auf das Storungssystem.

dstep = 1;

Nun sind alle System-, Modell- und Reglerparameter definiert und man kann zuden ersten Simulationen ubergehen.

1.2.5 Simulation des closed-loop-Verhaltens

Als erstes testen wir die Regler unter der Voraussetzung, dass wir keinen Modell-fehler gemacht haben, das modellierte System entspricht also dem realen System.Fur den Betrieb in der Fahrstufe

”D“ ergibt sich damit:

[y1,u1] = cmpc(systemD,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

figure;

plotall(y1,u1,ts);

Nun die gleiche Simulation fur den Betrieb in der Fahrstufe”N“:

[y2,u2] = cmpc(systemN,systemN,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dN,[],dstep);

figure;

plotall(y2,u2,ts);

Page 18: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

10 KAPITEL 1. MPC TOOLBOX IN MATLAB

Man beachte, dass wir nun zwei verschiedene Regler erstellt haben, die, jeder fursich, das jeweilige System stabilisieren (vergleiche die Plots in Abbildung A.1 undA.2 im Anhang A.1). Was man wirklich will, ist ein Regler, der das System inbeiden Fahrstufen regelt. Wir testen deswegen nun die Regler indem wir fur dasSystem

”D“ das Modell

”N“ und umgekehrt verwenden.

[y3,u3] = cmpc(systemD,systemN,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

figure;

plotall(y3,u3,ts);

[y4,u4] = cmpc(systemN,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dN,[],dstep);

figure;

plotall(y4,u4,ts);

0 1 2 3 4 5 6 7 8−15

−10

−5

0

5

10

15

Time

Outputs

y1

y2

Abbildung 1.2: Simulationfur die Fahrstufe

”D“ mit

Modellfehler

Man sieht, dass nun keiner der Regler das Systemstabilisieren kann, wenn ein Modellfehler vorliegt(vergleiche die Plots in Abbildung A.3 und A.4 imAnhang A.1 sowie in Abbildung 1.2 rechts). Manmuss also die Reglerparameter verandern um einebessere Robustheit zu erhalten. Ein Versuch kannsein, den Regeleingriff in der Optimierungsfunkti-on starker zu bestrafen um ein zu starkes Schwin-gen zu verhindern. Wir erhohen deswegen den Ge-wichtsfaktor fur die Bestrafung der Veranderungder Luftzufuhr von 0,5 auf 10 und den Gewichtsfak-tor fur die Bestrafung der Veranderung des Zund-zeitpunkts von 0,5 auf 20.

uwt = [10 20];

[y5,u5] = cmpc(systemD,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

figure;

plotall(y5,u5,ts);

[y6,u6] = cmpc(systemN,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dN,[],dstep);

figure;

plotall(y6,u6,ts);

Page 19: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.2. BEISPIEL: LEERLAUFFULLUNGSREGELUNG 11

0 1 2 3 4 5 6 7 8−3.5

−3

−2.5

−2

−1.5

−1

−0.5

0

0.5

1

Time

Outputs

y1

y2

Abbildung 1.3: Simulation fur dieFahrstufe

”N“ mit Modellfehler

Man sieht, dass fur das System ohne Mo-dellfehler der Regler weiterhin stabilisierendwirkt (vergleiche den Plot in Abbildung A.5im Anhang A.1), naturlich dauert es nunlanger, bis das System stabilisiert ist, da jaein starker Regeleingriff, der zu einer schnel-leren Nullregelung fuhrt, bestraft wird. Zumgleichen Ergebnis kommt man fur das Sys-tem in der Fahrstufe

”N“ wenn kein Modell-

fehler vorliegt. Fur das System mit Modell-fehler zeigt sich nun auch ein stabilisieren-des Regelverhalten (vergleiche den Plot inAbbildung A.6 im Anhang A.1 sowie in Ab-

bildung 1.3 links), das gleiche Ergebnis erhalt man fur das System”D“ mit feh-

lerbehaftetem Modell”N“.

Die Ergebnisse sollten jedoch nicht dazu verleiten, den erstellten Regler als globalasymptotisch stabil fur das System mit Modellfehler zu betrachten. Eine globaleasymptotische Stabilitat kann namlich nur fur die Regler ohne Modellfehler ga-rantiert werden, jedoch sollen uns in diesem Beispiel die erhaltenen Ergebnissegenugen.

1.2.6 Vergleich verschiedener Regler

Wir vergleichen nun den cmpc-Regler, der das quadratische Kostenfunktionalin jedem Schritt on-line minimiert, mit einem MPC-Regelansatz, der zuerst dieGain-Matrix berechnet und dann das System simuliert.

Kmpc = mpccon(systemD,ywt,uwt,m,p);

[y7,u7] = mpcsim(systemD,systemD,Kmpc,tfsim,r,ulim,tfilter,...

dD,[],dstep);

figure;

plot([0:ts:tfsim],y5(:,1),’-b’,[0:ts:tfsim],y5(:,2),’-k’,...

[0:ts:tfsim],y7(:,1),’-.b’,[0:ts:tfsim],y7(:,2),’-.k’);

legend(’y_1 mit cmpc’,’y_2 mit cmpc’,...

’y1 mit mpcsim’,’y_2 mit mpcsim’,0);

xlabel(’Time’);

title(’Outputs’);

Man sieht, dass sich die Ergebnisse nicht stark unterscheiden (vergleiche den Plotin Abbildung A.7 im Anhang A.1), obwohl wir zwei unterschiedliche Wege inder Optimierung des Kostenfunktionals gegangen sind. Der Vorteil der on-line-Optimierung lasst sich jedoch verdeutlichen, wenn wir den Pradiktionshorizontauf einen beschrankten Wert setzen, z. B. p = 80 und beide Simulationen nochein-mal durchfuhren.

Page 20: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

12 KAPITEL 1. MPC TOOLBOX IN MATLAB

p = 80;

[y8,u8] = cmpc(systemD,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

Kmpc = mpccon(systemD,ywt,uwt,m,p);

[y9,u9] = mpcsim(systemD,systemD,Kmpc,tfsim,r,ulim,tfilter,...

dD,[],dstep);

figure;

plot([0:ts:tfsim],y8(:,1),’-b’,[0:ts:tfsim],y8(:,2),’-k’,...

[0:ts:tfsim],y9(:,1),’-.b’,[0:ts:tfsim],y9(:,2),’-.k’);

legend(’y_1 mit cmpc’,’y_2 mit cmpc’,...

’y1 mit mpcsim’,’y_2 mit mpcsim’,0);

xlabel(’Time’);

title(’Outputs’);

Die Ergebnisse zeigen, dass in diesem Fall die on-line-Optimierung des Kosten-funktionals zu einer schnelleren Nullregelung fuhrt (vergleiche den Plot in Abbil-dung A.8 im Anhang A.1 sowie in Abbildung 1.4).

0 1 2 3 4 5 6 7 8−7

−6

−5

−4

−3

−2

−1

0

1

Time

Outputs

y1 mit cmpc

y2 mit cmpc

y1 mit mpcsimy

2 mit mpcsim

Abbildung 1.4: Vergleich von Simulationen fur die Fahrstufe”D“ mit beschrank-

tem Pradiktionshorizont

1.2.7 Import des Reglers in Simulink

Die gute Zusammenarbeit zwischen der MPC Toolbox und anderen Toolboxenin Matlab sowie Simulink erlaubt es, einen modellpradiktiven Regler fur lineareSysteme in Simulink an einer nichtlinearen Regelstrecke zu testen. Wir benutzen

Page 21: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.2. BEISPIEL: LEERLAUFFULLUNGSREGELUNG 13

hier den zuletzt betrachteten Regler fur beschrankte Systeme (cmpc, Pradikti-onshorizont p = 80, erhohte Gewichtsfaktoren). Simulink bietet zum Import eineS-Block-Funktion, die diesen Regler anspricht, die notigen Modelldaten und Reg-lereinstellungen mussen bereits im Arbeitsbereich geladen sein. Abbildung 1.5zeigt die lineare Regelstrecke und den modellpradiktiven Regler in Simulink.

Abbildung 1.5: Signalflussplan des linearen Kontrollsystems in Simulink.

Nun wird der Strecke eine leichte Nichtlinearitat sowie als eine weitere Storung einMessrauschen hinzugefugt. Die roten Blocke in Abbildung 1.6 zeigen die Nichtli-nearitaten und die zusatzliche Storung in der Regelstrecke. Die Abbildungen A.9und A.10 im Anhang A.1 zeigen, dass der Regler auch fur die nun nichtlineareStrecke gute Ergebnisse liefert.

Page 22: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

14 KAPITEL 1. MPC TOOLBOX IN MATLAB

Abbildung 1.6: Signalflussplan der nichtlinearen Regelstrecke mit zusatzlicherStorung.

1.3 Der Blick hinter die Kulissen der MPC Tool-

box

In diesem Abschnitt wird kurz auf die in der MPC Toolbox implementiertenAlgorithmen eingegangen.

1.3.1 Algorithmen fur unbeschrankte Probleme oder Pro-

bleme mit saturierten Beschrankungen

Fur Modelle (modelliert im Frequenzbereich oder im Zeitbereich) ohne Beschran-kungen oder mit saturierten Beschrankungen bietet die MPC Toolbox einen GPC-Ansatz (Generalized Predictive Control). Dabei wird zuerst iterativ uber die al-gebraische Riccati-Gleichung die sog. Gain-Matrix berechnet, im Anschluss wirdeine closed-loop-Simulation durchgefuhrt. Ist das Modell unbeschrankt, so lie-fert dieses Vorgehen gute Ergebnisse. Liegen nur schwache Beschrankungen vor,so kann man das System mit

”saturierten Beschrankungen“ betreiben. Dies be-

deutet, dass die Gain-Matrix, die ja fur unbeschrankte Probleme gilt, fur dasbeschrankte System eingesetzt wird. Sollte nun eine Kontrollanweisung oder eineAusgangsvariable eine Beschrankung verletzen, so wird die Kontrollanweisung soverandert, dass Beschrankungen nicht verletzt werden. Der

”zu große“ Anteil der

Page 23: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.3. DER BLICK HINTER DIE KULISSEN DER MPC TOOLBOX 15

Steuerung wird sozusagen”abgeschnitten“, so dass Beschrankungen erfullt blei-

ben, Abbildung 1.7 verdeutlicht dieses Vorgehen. Verstandlicherweise funktioniert

Abbildung 1.7: Prinzip der saturierten Beschrankungen

dieses Vorgehen nur fur schwach beschrankte Systeme und stellt in der Realitatsomit nicht das Verfahren der Wahl. Fur die Regelung uber den GPC wird aufKapitel 2 verwiesen.

1.3.2 Algorithmen fur beschrankte Probleme

Fur Modelle (modelliert im Frequenzbereich oder im Zeitbereich) mit Beschran-kungen wird das Kostenfunktional

min∆u(k)...∆u(k+m−l)

(

p∑

l=1

|Q (y(k + l|k) − r(k + l))|2 +

m∑

l=1

|R∆u(k + l − 1)|2)

in jedem Schritt online unter Beachtung der Beschrankungen gelost. Dabei be-zeichnen Q und R die Gewichtsmatrizen, y(k + l|k) den l-ten vorausgesagtenAusgangswert zum k-ten Abtastschritt und u(k + m − l) die k + m − l-te Kon-trollanweisung.Zur Losung dieses Problems verwendet die MPC Toolbox die Dantzig1-Wolfe2-Methode der aktiven Mengen (ausfuhrlich erlautert im Buch von Fletcher [3,Kapitel 10.3 und 10.6]). Dabei wird das quadratische Programm

minx

(

1

2xT Hx + fT x

)

=: minx

(q(x))

Ax ≤ b, x ≥ xmin

iterativ gelost. Voraussetzung ist hierfur, dass die Hessematrix H positiv definitist, dies garantiert, dass das Problem eindeutig losbar ist.

1George Bernard Dantzig, 1914 - 2005, US-amerikanischer Mathematiker2Philip Wolfe, US-amerikanischer Mathematiker

Page 24: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

16 KAPITEL 1. MPC TOOLBOX IN MATLAB

Die Dantzig-Wolfe-Methode der aktiven Mengen ist eine Erweiterung des be-kannten Simplex-Algorithmus, bei dem entlang des Beschrankungsrandes gesuchtwird. Da hier eine quadratische Funktion vorliegt, muss das Minimum aber nichtauf dem Rand des zu untersuchenden Gebiets liegen. Man unterscheidet deswegenzwei Falle:Zuerst sucht man das globale Minimum. Verletzt dies keine Beschrankungen, soist dies auch das Minimum unter Beschrankungen.Verletzt dieses globale Minimum aber Beschrankungen, so liegt das Minimum un-ter Beschrankungen auf dem Rand der zulassigen Menge. Man beachtet nun einegewisse Auswahlmenge (die aktive Menge) an Beschrankungen, die uber Gleich-heit anstatt Ungleichheit erfullt werden mussen. Im ersten Schritt sind dies dieBeschrankungen, welche vom globalen Minimum verletzt werden. Man sucht nunein Minimum auf dieser aktiven Menge an Restriktionen und uberpruft, ob nichtaktive Beschrankungen verletzt werden. Ist dies der Fall, werden diese mit indie aktive Menge der Beschrankungen aufgenommen (wieder mittels Gleichheitanstatt Ungleichheit) und es wird erneut gesucht. Werden aber keine inaktivenBeschrankungen verletzt, so ist ein potentieller Kandidat fur das Minimum unterBeschrankungen gefunden. Verletzt der Kandidat die Karush3-Kuhn4-Tucker5-Bedingung (s. u.) fur eine aktive Beschrankung, so wird diese inaktiv und eskann weitergesucht werden. Erfullt aber der potentielle Kandidat die Karush-Kuhn-Tucker-Bedingung, so ist er auch Losung des Problems.Die Karush-Kuhn-Tucker-Bedingung lasst sich folgendermaßen formulieren: Istx∗ ein lokales Minimum einer skalaren Funktion q(x) := 1

2xT Hx + fT x unter den

Nebenbedingungen g(x) := b − Ax ≤ 0 und gilt fur die Gradienten der aktivenGleichungen gi(x

∗) = 0 fur i ∈ I (aktive Menge), dass sie linear unabhangig sind,so existieren eindeutige Lagrange-Multiplikatoren λi ≥ 0, so dass

∇q(x∗) =∑

i∈I

λi∇gi(x∗).

Abbildung 1.8 zeigt die Vorgehensweise des Algorithmus.Zur genaueren Beschreibung des Algorithmus wird noch folgendes vermerkt: Imersten Schritt muss nicht notwendigerweise das globale Minimum fur das un-beschrankte System gesucht werden, ein

”Warmstart“ innerhalb der zulassigen

Menge ist moglich. Sollte in einem Suchschritt keine Beschrankung aktiv sein, sobestimmt sich die Suchrichtung als negativer Gradient der Zielfunktion, also als−Hx.Die Implementierung geschieht dabei analog zum Simplex-Algorithmus uber einTableau, in dem ein transformiertes und erweitertes System gespeichert wird. Injedem Iterationsschritt verandern Zeilenoperationen das Tableau so, dass am En-

3William Karush, 1917 - 1997, US-amerikanischer Mathematiker4Harold W. Kuhn, geb. 1925, US-amerikanischer Mathematiker5Albert William Tucker, 1905 - 1995, US-amerikanischer Mathematiker

Page 25: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

1.3. DER BLICK HINTER DIE KULISSEN DER MPC TOOLBOX 17

Abbildung 1.8: Dantzig-Wolfe-Methode der aktiven Mengen, aus [7].

de eine spezielle Form erhalten wird und die Losung abgelesen werden kann. Einegenauere Besprechung der Implementierung findet sich im Buch von Fletcher [3,Kapitel 10.6].Es lasst sich zeigen, dass der Algorithmus in einer endlichen Schrittzahl ein Er-gebnis liefert. Der Aufwand fur die Losung betragt O(n + p), wobei n = dim(x)und p = dim(b), nur selten werden dabei mehr als 3(n + p) Iterationen benotigt(eine ausfuhrliche Begrundung dieses Resultats findet sich im Buch von Dantzig[2, Kapitel 7.3]).Zur Verdeutlichung des Algorithmus folgt ein einfaches Beispiel im R

2: Abbildung1.9 zeigt die Hohenlinien der quadratischen Zielfunktion, der Punkt x stellt dabeidas globale Minimum ohne Beschrankungen dar, die weiße Flache beschreibt diezulassige Menge (beschrankt durch g1 ≥ 0, g2 ≥ 0 und g3 ≤ 0). Man startet miteinem Punkt in der zulassigen Menge, hier ist es der Punkt x(0). Dieser erfulltdie Beschrankung g1, weswegen diese in der aktiven Menge der Beschrankungenist. Man sucht nun das Minimum unter der Voraussetzung, dass x2 = 0, dies istx(1). Nun wird die zweite Beschrankung g2 aktiv: Uberpruft man nun die Karush-Kuhn-Tucker-Bedingung, so stellt man fest, dass

∇f

‖∇f‖= −

∇g1

‖∇g1‖

Page 26: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

18 KAPITEL 1. MPC TOOLBOX IN MATLAB

gilt, da der Gradient von f in negative x2-Richtung, der Gradient der Neben-bedingung aber in positive x2-Richtung zeigt. Die Beschrankung g1 wird alsoinaktiv.Man muss nun das Minimum unter der Bedingung g2 = 0 suchen. Dies findet sichin x(2). Betrachtet man wieder den Gradienten der Zielfunktion (dieser zeigt innegative x1-Richtung), so ergibt sich wieder ein negativer Lagrange-Multiplikator(der Gradient der Bedingung g2 zeigt in positive x1-Richtung). Die Beschrankungg2 wird also auch inaktiv.Beachte, nun sind keine Bedingungen aktiv und es wird in Richtung des aktu-ellen negativen Gradienten der Zielfunktion gesucht. Man gelangt schließlich zux(3), was nun die Beschrankung g3 aktiv werden lasst. Eine Suche unter der Be-schrankung g3 = 0 liefert dann das Minimum x(4). Hier zeigt der Gradient derZielfunktion in die selbe Richtung wie der Gradient der einzigen aktiven Neben-bedingung g3. Wir haben also die Karush-Kuhn-Tucker-Bedingung erfullt undsomit das Minimum unter den Beschrankungen gefunden.

Abbildung 1.9: Beispiel zur Dantzig-Wolfe-Methode der aktiven Mengen.

Page 27: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

Anhang A

Plots und Quelltext

A.1 Plots

0 2 4 6 8 10 12 14 16 18 20−5

−4

−3

−2

−1

0

1Outputs

Time

0 2 4 6 8 10 12 14 16 18 20−2

0

2

4

6

8Manipulated Variables

Time

y1

y2

u1

u2

Abbildung A.1: Simulation fur die Fahrstufe”D“ ohne Modellfehler. Man sieht,

dass der Regler das System asymptotisch stabil zur 0 regelt.

19

Page 28: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

20 ANHANG A. PLOTS UND QUELLTEXT

0 2 4 6 8 10 12 14 16 18 20−2.5

−2

−1.5

−1

−0.5

0

0.5

1Outputs

Time

0 2 4 6 8 10 12 14 16 18 20−1

0

1

2

3

4

5Manipulated Variables

Time

y2

y1

u1

u2

Abbildung A.2: Simulation fur die Fahrstufe”N“ ohne Modellfehler. Ebenso beim

zweiten Modell liefert der Regler ein asymptotisch stabiles Verhalten.

0 2 4 6 8 10 12 14 16 18 20−20

−10

0

10

20Outputs

Time

0 2 4 6 8 10 12 14 16 18 20−10

−5

0

5

10

15

20Manipulated Variables

Time

y1

y2

u1

u2

Abbildung A.3: Simulation fur die Fahrstufe”D“, wobei als Systemmodell die

Fahrstufe”N“ gewahlt wurde. Man sieht, dass der Regler das System nicht sta-

bilisieren kann.

Page 29: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

A.1. PLOTS 21

0 2 4 6 8 10 12 14 16 18 20−3

−2

−1

0

1

2Outputs

Time

0 2 4 6 8 10 12 14 16 18 20−1

0

1

2

3

4Manipulated Variables

Time

u1

u2

Abbildung A.4: Simulation fur die Fahrstufe”N“, wobei als Systemmodell die

Fahrstufe”D“ gewahlt wurde. Auch hier kann der Regler das System nicht stabi-

lisieren.

0 2 4 6 8 10 12 14 16 18 20−6

−5

−4

−3

−2

−1

0

1Outputs

Time

0 2 4 6 8 10 12 14 16 18 20−2

0

2

4

6

8Manipulated Variables

Time

y1

y2

u1

u2

Abbildung A.5: Simulation fur die Fahrstufe”D“ ohne Modellfehler mit erhohten

Gewichtsfaktoren fur die Bestrafung des Regeleingriffs in der Kostenfunktion.

Page 30: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

22 ANHANG A. PLOTS UND QUELLTEXT

0 2 4 6 8 10 12 14 16 18 20−4

−3

−2

−1

0

1Outputs

Time

0 2 4 6 8 10 12 14 16 18 200

0.5

1

1.5

2

2.5

3

3.5Manipulated Variables

Time

u1

u2

u1

u2

Abbildung A.6: Simulation fur die Fahrstufe”N“, wobei als Systemmodell die

Fahrstufe”D“ gewahlt wurde mit erhohten Gewichtsfaktoren fur die Bestrafung

des Regeleingriffs in der Kostenfunktion. Nun zeigt sich, dass der Regler auch furdas System mit Modellfehler ein stabilisierendes Ergebnis liefert.

Page 31: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

A.1. PLOTS 23

0 2 4 6 8 10 12 14 16 18 20−6

−5

−4

−3

−2

−1

0

1

Time

Outputs

y1 mit cmpc

y2 mit cmpc

y1 mit mpcsimy

2 mit mpcsim

Abbildung A.7: Vergleich von Simulationen fur die Fahrstufe”D“ ohne Modell-

fehler und Pradiktionshorizont p = ∞. Durchgezogene Linien: cmpc mit on-line-Optimierung des Kostenfunktionals. Strich-Punkt-Linien: mpcsim mit im Vorausberechneter Gain-Matrix. Es zeigt sich kein wesentlicher Unterschied zwischenden verschiedenen Verfahren.

Page 32: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

24 ANHANG A. PLOTS UND QUELLTEXT

0 2 4 6 8 10 12 14 16 18 20−7

−6

−5

−4

−3

−2

−1

0

1

Time

Outputs

y1 mit cmpc

y2 mit cmpc

y1 mit mpcsimy

2 mit mpcsim

Abbildung A.8: Vergleich von Simulationen fur die Fahrstufe”D“ ohne Modell-

fehler und Pradiktionshorizont p = 80. Durchgezogene Linien: cmpc mit on-line-Optimierung des Kostenfunktionals. Strich-Punkt-Linien: mpcsim mit im Vorausberechneter Gain-Matrix. Die on-line-Optimierung ist hier klar im Vorteil.

Page 33: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

A.1. PLOTS 25

0 2 4 6 8 10 12 14 16 18 20−8

−6

−4

−2

0

2

y1

Time

0 2 4 6 8 10 12 14 16 18 20−0.2

0

0.2

0.4

0.6

0.8

y2

Time

Abbildung A.9: Simulation des linearen cmpc-Reglers an der nichtlinearen Regel-strecke, Ausgangsvariablen.

0 2 4 6 8 10 12 14 16 18 20−1

0

1

2

3

4

5

6

u1

Time

0 2 4 6 8 10 12 14 16 18 20−0.2

0

0.2

0.4

0.6

0.8

u2

Time

Abbildung A.10: Simulation des linearen cmpc-Reglers an der nichtlinearen Re-gelstrecke, Kontrollvariablen.

Page 34: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

26 ANHANG A. PLOTS UND QUELLTEXT

A.2 Quelltext zu Kapitel 1.2

function bsp();

% Definition der Systemuebertragungsfunktionen

G1D = [0 0 9.62;1 2.4 5.05;0 0.16 0];

G2D = [0 15.9 15.9*3;1 2.4 5.05;0 0.04 0];

G1N = [0 0 20.5; 1 2.2 12.8; 0 0.16 0];

G2N = [0 47.6 47.6*3.5; 1 2.2 12.8; 0 0.04 0];

G0 = [0 0;0 0;0 0];

GI = [0 1;0 1;0 0];

% Approximation der Systemsprungfunktionen

tf = 4; % Endzeit der Approximation

ts = 0.1; % Abtastzeit

systemD = tfd2step(tf,ts,2,G1D,G0,G2D,GI);

systemN = tfd2step(tf,ts,2,G1N,G0,G2N,GI);

% Definition der Stoerungsuebertragungsfunktionen

% und Approximation der Stoerungssprungfunktionen

GdD = [0 -19.1 -19.1*3; 1 2.4 5.05; 0 0 0];

dD = tfd2step(tf,ts,2,GdD,G0);

GdN = [0 -19.1 -19.1*3.5; 1 2.2 12.8; 0 0 0];

dN = tfd2step(tf,ts,2,GdN,G0);

% Definition der Reglerparameter

ywt = [5 1]; % Gewichte fuer die Ausgaenge

uwt = [.5 .5]; % Gewichte fuer die Kontrolleingriffe

m = 10; % Kontrollhorizont

p = inf; % Praediktionshorizont

tfsim = 20; % Dauer der Simulation

r = [0 0]; % Setpunkte

ulim = [[-inf -0.7] [inf 0.7] [1e6 1.4]];

% Beschraenkungen an die Kontrolle

ylim = []; % Beschraenkungen an die Ausgaenge

tfilter = []; % Stoerungsfilter, -verzoegerung

dstep = 1; % Sprungfunktion fuer Stoerung

%================================================================

% 1. Simulation ohne Modellfehler

%===============================================================

% Fahrstufe D

Page 35: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

A.2. QUELLTEXT ZU KAPITEL 1.2 27

[y1,u1] = cmpc(systemD,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

figure; % neues Plotfenster

plotall(y1,u1,ts); % plottet Ausgaenge und Kontrollen

% Fahrstufe N

[y2,u2] = cmpc(systemN,systemN,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dN,[],dstep);

figure;

plotall(y2,u2,ts);

%===============================================================

% 2. Simulation mit Modellfehler

%===============================================================

% Fahrstufe D

[y3,u3] = cmpc(systemD,systemN,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

figure;

plotall(y3,u3,ts);

%Fahrstufe N

[y4,u4] = cmpc(systemN,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dN,[],dstep);

figure;

plotall(y4,u4,ts);

%===============================================================

% 3. Simulation mit erhoehtem Gewichtsfaktor fuer die Kontrolle

%===============================================================

uwt = [10 20]; % Gewichtsfaktor fuer die Kontrolle

% Fahrstufe D ohne Modellfehler

[y5,u5] = cmpc(systemD,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

figure;

plotall(y5,u5,ts);

% Fahrstufe N mit Modellfehler

[y6,u6] = cmpc(systemN,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dN,[],dstep);

figure;

Page 36: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

28 ANHANG A. PLOTS UND QUELLTEXT

plotall(y6,u6,ts);

%===============================================================

% 4. Simulation - Vergleich zwischen cmpc und mpcsim

%===============================================================

Kmpc = mpccon(systemD,ywt,uwt,m,p); % Erstellung der Gain-Matrix

[y7,u7] = mpcsim(systemD,systemD,Kmpc,tfsim,r,ulim,tfilter,...

dD,[],dstep);

figure;

plot([0:ts:tfsim],y5(:,1),’-b’,[0:ts:tfsim],y5(:,2),’-k’,...

[0:ts:tfsim],y7(:,1),’-.b’,[0:ts:tfsim],y7(:,2),’-.k’);

legend(’y_1 mit cmpc’,’y_2 mit cmpc’,...

’y1 mit mpcsim’,’y_2 mit mpcsim’,0);

xlabel(’Time’);

title(’Outputs’);

%===============================================================

% 5. Simulation - Vergleich zwischen cmpc und mpcsim

% mit beschraenktem Praediktionshorizont

%===============================================================

p = 80; % neuer Praediktionshorizont

[y8,u8] = cmpc(systemD,systemD,ywt,uwt,m,p,tfsim,r,ulim,ylim,...

tfilter,dD,[],dstep);

Kmpc = mpccon(systemD,ywt,uwt,m,p);

[y9,u9] = mpcsim(systemD,systemD,Kmpc,tfsim,r,ulim,tfilter,...

dD,[],dstep);

figure;

plot([0:ts:tfsim],y8(:,1),’-b’,[0:ts:tfsim],y8(:,2),’-k’,...

[0:ts:tfsim],y9(:,1),’-.b’,[0:ts:tfsim],y9(:,2),’-.k’);

legend(’y_1 mit cmpc’,’y_2 mit cmpc’,...

’y1 mit mpcsim’,’y_2 mit mpcsim’,0);

xlabel(’Time’);

title(’Outputs’);

Page 37: MPC Toolbox in Matlab - uni-bayreuth.denum.math.uni-bayreuth.de/de/teaching/archive/ss_2006/01062/ausar… · Einordnung der Toolbox in die Matlab-Produktfamilie erl¨autert, im Anschluss

Literaturverzeichnis

[1] A. Bemporad, M. Morari, N. L. Ricker, Getting Started with the Model Pre-

dictive Control Toolbox ; Version 2, The MathWorks, 2004 - 2006.

[2] G. B. Dantzig, Linear Programming and Extensions; Princeton UniversityPress, Princeton, 1963.

[3] R. Fletcher, Practical Methods of Optimization; John Wiley & Sons, Chiche-ster, UK, 1986.

[4] L. Guzzella, C. H. Onder, Introduction to Modeling and Control of Internal

Combustion Engine Systems ; Springer, 2004.

[5] D. Homberg, Nichtlineare Optimierung, Script zur Vorlesung; Tech-nische Universitat Berlin, 2005/2006, http://www.math.tu-berlin.de/

Vorlesungen/WS05/NonLinOpt/skript/nonlin-opt.pdf.

[6] D. Hrovat and B. Bodenheimer, Robust Automotive Idle Speed Control De-

sign Based on µ-Synthesis ; Proceedings of American Control Conference,1993.

[7] J. Matthes, L. Groll, R. Mikut, J. Jakel, Optimale Fuhrung von

Endoskopen mit redundanter Kinematik ; Automatisierungstechnik 48,R. Oldenbourg Verlag, 2000, http://www.sfb588.uni-karlsruhe.de/

publikationen/2001_06.pdf.

[8] M. Morari, N. L. Ricker, Model Predictive Control Toolbox User’s Guide;Version 1, The MathWorks, 1995 - 1998.

[9] MPC Toolbox Website: http://www.mathworks.com/products/mpc/.

29