Anfang Präsentation 2. Februar, 2005 Behandlung von Unstetigkeiten Wir wollen uns heute dem Problem...
26
2. Februar, 2005 Anfang Präsentation Behandlung von Unstetigkeiten • Wir wollen uns heute dem Problem der Behandlung von Unstetigkeiten in der Modellbeschreibung zuwenden. • Modelle aus den Ingenieurbereichen weisen häufig Unste-tigkeiten auf, die Schaltvorgänge, Begrenzer, trockene Reibung, Stösse oder ähnliche Vorgänge repräsentieren mögen. • Die Modellierungsumgebung muss sich dieser Vorgänge speziell annehmen, da sie in starkem Masse die Numerik der Differentialgleichungslöser beeinflussen.
Transcript of Anfang Präsentation 2. Februar, 2005 Behandlung von Unstetigkeiten Wir wollen uns heute dem Problem...
2. Februar, 2005Anfang Präsentation
Behandlung von Unstetigkeiten• Wir wollen uns heute dem Problem der Behandlung von
Unstetigkeiten in der Modellbeschreibung zuwenden.• Modelle aus den Ingenieurbereichen weisen häufig Unste-
tigkeiten auf, die Schaltvorgänge, Begrenzer, trockene Reibung, Stösse oder ähnliche Vorgänge repräsentieren mögen.
• Die Modellierungsumgebung muss sich dieser Vorgänge speziell annehmen, da sie in starkem Masse die Numerik der Differentialgleichungslöser beeinflussen.
2. Februar, 2005Anfang Präsentation
Übersicht• Die numerische Lösung von Differentialgleichungen• Unstetigkeiten in Zustandsgleichungen• Integration über Unstetigkeiten• Zustandsereignisse• Ereignisbearbeitung• Mehrwertige Kennlinien• Der elektrische Schalter• Die ideale Diode• Die Reibungskennlinie
2. Februar, 2005Anfang Präsentation
Numerische Differentialgleichungslöser• Alle heute im Einsatz stehenden numerischen
Differentialgleichungslöser arbeiten mit Polynom-extrapolation.
• Der Wert einer Zustandsvariablen x zum Zeitpunkt t+h, wobei h die momentane Schrittlänge ist, wird ermittelt, indem durch bereits bekannte Stützwerte von x und dx/dt der jetzigen Zeit t sowie der Vergangenheit ein Polynom nter Ordnung hindurchgelegt wird. Der Wert dieses Polynoms zum Zeitpunkt t+h repräsentiert die angenäherte Lösung der Differentialgleichung.
• Bei impliziten Verfahren wird die Zustandsableitung zum Zeitpunkt t+h ebenfalls als Stützwert beigezogen.
2. Februar, 2005Anfang Präsentation
Beispiele
Explizites Eulerverfahren erster Ordnung:
x(t+h) x(t) + h · x(t)
·
Implizites Eulerverfahren erster Ordnung:
x(t+h) x(t) + h · x(t+h)
·
2. Februar, 2005Anfang Präsentation
Unstetigkeiten in der Zustandsbeschreibung• Polynome sind immer stetig und stetig differenzierbar.• Wenn also die Zustandsraumbeschreibung des Systems:
• eine Unstetigkeit aufweist, ist die Polynomextrapolation eine sehr schlechte Approximation der Wirklichkeit.
• Somit weisen Verfahren mit fester Schrittlänge einen grossen Integrationsfehler auf, während Verfahren mit variabler Schrittlänge den Schritt in der Umgebung der Unstetigkeit sehr stark reduzieren.
x(t) = f(x(t),t)·
2. Februar, 2005Anfang Präsentation
Integration über Unstetigkeiten• Ein Verfahren variabler Schrittlänge reduziert die
Schrittlänge bei jeder Unstetigkeit.• Nach der Unstetigkeit wird die Schrittlänge nur langsam
wieder vergrössert, da das Verfahren nicht unterscheiden kann zwischen einer Unstetigkeit und einem Punkt grosser lokaler Steifigkeit (mit grossem absolutem Wert der Ableitung). h
t
UnstetigkeitenDie Schrittlänge ist dauernd viel zu klein. Das Verfahren ist zu-mindest sehr ineffizient, falls nicht ungenau.
2. Februar, 2005Anfang Präsentation
Das Zustandsereignis• Diese Probleme können vermieden werden, wenn dem
System explizit mitgeteilt wird, wo es Unstetigkeiten in der Modellbeschreibung gibt.
Beispiel: Begrenzerfunktion
x
f(x)
xp
fm
xm
fp
m = tg()1
f = fm
2f = m·x
3
f = fp
f = if x < xm then fm else if x < xp then m*x else fp ;
2. Februar, 2005Anfang Präsentation
Die Ereignisbearbeitung I
x
f(x)
xp
fm
xm
fp
xxp
t
Schwellwert
h
Iteration
h
xxp
t
Modellumschaltung
h
h
t
EreignisSchrittlängenreduktion während Iteration
2. Februar, 2005Anfang Präsentation
Die Ereignisbearbeitung II
h
t
h
t
Schrittlänge ohne Ereignis-behandlung
Schrittlänge mit Ereignis-behandlung
2. Februar, 2005Anfang Präsentation
Repräsentation von Unstetigkeiten
• In Modelica wird die Unstetigkeit durch ein if-Statement repräsentiert.
• Bei der Übersetzung werden diese Statements in korrekte Ereignisbeschreibungen (Sätze von Modellen mit Um-schaltbedingungen) abgebildet.
• Der Benützer braucht sich um die Ereignisbeschreibungs-mechanismen selbst nicht zu kümmern. Diese bleiben hinter den if-Statements verborgen.
f = if x < xm then fm else if x < xp then m*x else fp ;
2. Februar, 2005Anfang Präsentation
Probleme• Der Benützer sollte beachten, dass während der Iteration
die unstetige Kennlinie temporär verlassen wird.
• ist gefährlich, da absp temporär negativ werden kann.
• löst dieses Problem.
q = | p | p = p1 – p2 ;absp = if p > 0 then p else –p ;q = sqrt(absp) ;
p = p1 – p2 ;absp = noEvent( if p > 0 then p else –p ) ;q = sqrt(absp) ;
2. Februar, 2005Anfang Präsentation
Die „noEvent“-Klausel
• Die noEvent-Klausel sorgt dafür, dass if-Statements oder Boole’sche Ausdrücke, die normalerweise in Simulations-code übersetzt würden, bei welchem die Diskontinuitäten durch Ereignisse abgefangen werden, so stehengelassen werden, wie sie im Modellierungscode formuliert wurden.
• Dadurch wird die korrekte Simulation dieser Diskontinuitäten der Schrittsteuerung des numerischen Integrationsverfahrens überlassen.
p = p1 – p2 ;absp = noEvent( if p > 0 then p else –p ) ;q = sqrt(absp) ;
2. Februar, 2005Anfang Präsentation
Mehrwertige Kennlinien I• Die bisher eingeführten Sprachelemente reichen nicht aus,
um damit mehrwertige Kennlinien, wie zum Beispiel Hysteresefunktionen zu beschreiben.
• Wenn x grösser wird als xp, muss f von fm nach fp umgeschaltet werden.
• Wenn x kleiner wird als xm, muss f von fp nach fm umgeschaltet werden.
x
f(x)
xpxm
fp
fm
2. Februar, 2005Anfang Präsentation
Mehrwertige Kennlinien II
x
f(x)
xpxm
fp
fm
when initial() then reinit(f , fp);end when;when x > xp or x < xm then f = if x > 0 then fp else fm;end when;
}Diese Statements werden nur ausgeführt, wenn entweder x grösser wird als xp oder aber wenn x kleiner wird als xm.
Ausgeführt am Anfang der Simulation.
wird grösser
ist grösser
2. Februar, 2005Anfang Präsentation
Der elektrische Schalter I
i
u
Wenn der Schalter offen ist, ist der Strom i=0. Wenn der Schalter zu ist, ist die Spannung u=0.
0 = if offen then i else u ;
Das if-Statement in Modelica ist akausal. Es wird sortiert, genau wie alle anderen Statements.
2. Februar, 2005Anfang Präsentation
Der elektrische Schalter II
Mögliche Implementation: Schalter offen: s = 1 Schalter geschlossen: s = 0
0 = s · i + ( 1 – s ) · u
Schalter offen:
SF
Schalter geschlossen:SE
f = 0
e = 0
Sws ef
Die Kausalität des Switch Elements ist eine Funktion des Werts des Steuersignals s.
2. Februar, 2005Anfang Präsentation
Die ideale Diode I
u
i
Schalter geschlossen
Schalter offen
iu
Wenn u < 0 ist, ist der Schalter offen. Kein Strom fliesst.
Wenn u > 0 ist, ist der Schalter geschlossen. Strom kann fliessen. Die ideale Diode verhält sich wie ein Kurzschluss.
offen = u < 0 ; 0 = if offen then i else u ;
D ef
2. Februar, 2005Anfang Präsentation
Die ideale Diode II• Da Strom, wenn er durch die Diode fliesst, nicht einfach
unterbrochen werden kann, muss das Diodenmodell noch leicht modifiziert werden.
• Die Variable offen muss als Boolean deklariert werden. Ihr wird der Wert des logischen Ausdrucks rechts vom Gleichheitszeichen zugewiesen.
offen = u <= 0 and not i > 0 ; 0 = if offen then i else u ;
2. Februar, 2005Anfang Präsentation
Die Reibungskennlinie I• Kompliziertere Vorgänge, wie zum Beispiel die Reibungs-
kennlinie müssen sorgfältig analysiert werden.• Der Vorgang wird hier an Hand der Reibung erklärt.
fB
v
R0Rm
-Rm
-R0
Gleitreibung
Haftreibung
Wenn v ist die Rei-bungskraft eine Funktion der Geschwindigkeit.
Wenn v , wird die Reibungskraft so gerech-net, dass die Geschwin-digkeit 0 bleibt.
2. Februar, 2005Anfang Präsentation
Die Reibungskennlinie II• Wir unterscheiden fünf Fälle:
Haftung: Die Reibungskraft kompensiert die Summe der angrei-fenden Kräfte, es sei denn |f | > R0 .
Vorwärtsbewegung: Die Reibungskraft berechnet sich zu:fB = Rv · v + Rm.
Rückwärtsbewegung: Die Reibungskraft berechnet sich zu:fB = Rv · v Rm.
Einsetzen der Vorwärtsbewegung:
Die Reibungskraft berechnet sich zu:fB = Rm.
Einsetzen der Rückwärtsbewegung:
Die Reibungskraft berechnet sich zu:fB = Rm.
v = 0 a = 0
v = 0 a > 0
v = 0 a < 0
v > 0
v < 0
2. Februar, 2005Anfang Präsentation
Das Zustandsübergangsdiagramm• Der genaue Vorgang kann in einem Zustandsübergangs-
diagramm dargestellt werden.
Anfang
Rückwärts-bewegung
(v < 0)
Rückwärts-beschleunigung
(a < 0)Haftung (a = 0)
Vorwärts-beschleunigung
(a > 0)
Vorwärts-bewegung
(v > 0)
v < 0 v > 0v = 0
f < R0 f > R0
a 0 and not v < 0 a 0 and not v > 0
v < v >
v 0
v 0
2. Februar, 2005Anfang Präsentation
Das Reibungsmodell Imodel Reibung; parameter Real R0, Rm, Rv; parameter Boolean ic=false; Real fB, fc; Boolean Haft(final start = ic); Boolean Vor(final start = ic), Rück(final start = ic); Boolean StartVor(final start = ic), StartRück(final start = ic);
fB = if Vor then Rv*v + Rm else if Rück then Rv*v - Rm else if StartVor then Rm else if StartRück then -Rm else fc; 0 = if Haft or initial() then a else fc;
2. Februar, 2005Anfang Präsentation
Das Reibungsmodell II
when Haft and not initial() then reinit(v,0); end when;
Vor = initial() and v > 0 or pre(StartVor) and v > 0 or pre(Vor) and not v <= 0; Rück = initial() and v < 0 or pre(StartRück) and v < 0 or pre(Rück) and not v >= 0;
2. Februar, 2005Anfang Präsentation
Das Reibungsmodell III
StartVor = pre(Haft) and fc > R0 or pre(StartVor) and not (v > 0 or a <= 0 and not v > 0); StartRück = pre(Haft) and fc < -R0 or pre(StartRück) and not (v < 0 or a >= 0 and not v < 0); Haft = not (Vor or Rück or StartVor or StartRück); end Reibung;
2. Februar, 2005Anfang Präsentation
Referenzen I• Cellier, F.E. (1979), Combined Continuous/Discrete
System Simulation by Use of Digital Computers: Techniques and Tools, Swiss Federal Institute of Technology, ETH Zürich, Switzerland.
• Elmqvist, H., F.E. Cellier, and M. Otter (1993), “Object-oriented Modeling of Hybrid Systems,” Proc. ESS'93, SCS European Simulation Symposium, Delft, The Netherlands, pp.xxxi-xli.
• Cellier, F.E., M. Otter, and H. Elmqvist (1995), “Bond Graph Modeling of Variable Structure Systems,” Proc. ICBGM'95, 2nd SCS Intl. Conf. on Bond Graph Modeling and Simulation, Las Vegas, NV, pp.49-55.
2. Februar, 2005Anfang Präsentation
Referenzen II• Elmqvist, H., F.E. Cellier, and M. Otter (1994), “Object-
oriented Modeling of Power-electronic Circuits Using Dymola,” Proc. CISS'94, First Joint Conference of International Simulation Societies, Zurich, Switzerland, pp. 156-161.
• Glaser, J.S., F.E. Cellier, and A.F. Witulski (1995), “Object-oriented Switching Power Converter Modeling Using Dymola with Event-handling,” Proc. OOS'95, SCS Object-Oriented Simulation Conference, Las Vegas, NV, pp. 141-146.
