Startseite   Übersicht


Numerische Integrationsverfahren (I)

Einleitung

Bei Simulationen am Computer kommt man nicht drum herum, die Formeln in ein Stück'chen Code zu gießen. In der (allgemein gesagt) Wissenschaft werden die Zusammenhänge meist mit Differentialgleichungen beschrieben und diese müssen mit Integralen berechnet werden. Diese Integrale werden dann aber nicht analytisch gelöst (also nicht die Integrale formelmäßig bestimmen), sondern die Zahlenwerte werden bestimmt. Dies ist ein Themengebiet der Numerik. Das ist also - einfach gesagt - ein Rechnen mit Zahlen und nicht mit Formeln (was aber zum Teil sehr kompliziert sein kann).

Ich möchte mal einige, sehr einfache numerische Integrationsverfahren vorstellen, die zur Simulation mit einem Computer eingesetzt werden können. Dabei beschränke ich mich erst mal auf das Wesentliche, um nicht zu weit abzuschweifen. Die Differentialgleichungen werden nur angeschnitten und darauf aufbauend die Integrale grafisch erläutert. Die Rechenvorschriften für die numerische Lösung sind erst mal nur eine einfache Auswahl.

Einfach...damit diese von dem ambitionierten Programmierer auch leicht umgesetzt werden können. Einfach...damit der Sinn der Algorithmen noch nachvollzogen werden kann. Einfach...damit der ein oder andere auch weiter wie bis hier hin liest. Einfach...nur so.

Differentialgleichungen

...jetzt kommt ein wenig was aus der Mathematik...

In der technischen Welt werden Vorgänge sehr häufig mit Differentialgleichungen beschrieben. Gewöhnliche Differentialgleichungen sind Gleichungen in der Form

(es heißt dny(x), aber der Formeleditor wollte nicht...)

oder kurz geschrieben

Differentialgleichungen sind also einfach Gleichungen, in denen die (mehrfache) Ableitungen einer Funktion (und die Funktion selbst) vorkommen. y(n)' ist die hierbei n-te Ableitung von y. Die Notation mit dem Strich ist nur eine Kurzschreibweise für die sonst ausführliche Notation mit dy nach dx.

Allgemein (implizit) formuliert gilt

Das heißt, das F[...] von x, y und den Ableitungen von y abhängt (erste, zweite, ... n-te Ableitung). Die Ableitung einer Funktion wird i.d.R. in der Schule im Thema Differentialrechnung (und damit eng Verbunden die Integralrechnung) behandelt.

Differentialgleichungen erster Ordnung werden

oder kurz

beschrieben. Hierbei ist dann y' (die ersten Ableitung von y) nur von x und y abhängig.

Die Umkehrung der Differentialrechnung ist die Integralrechnung. Betrachtet man hierzu nochmal eine Differentialgleichung erster Ordnung

und integriert beide Seiten

so erhält man ja links vom Gleichheitszeichen gerade y(x), also

Ganz allgemein gilt ja, das f(x) über Integration von f '(x) bestimmt wird.

Das heißt ja nichts anderes, als das man Differentialgleichungen über Integrale löst (die Ableitung einer Funktion ist ja beschrieben und man will die Funktion selbst haben). Also für die Simulation am Computer heißt das, das man sich mit (numerischen) Integrationsverfahren herumplagen muss. Man interessiert sich ja nicht für die n-te Ableitung einer Funktion, sondern für die Funktion selbst. Die Ableitung ist ja nur (das notwendige Übel für) die einfache Beschreibung eines technischen Vorgangs/Zusammenhangs.

Hat man beispielsweise die siebte Ableitung einer Funktion (woher auch immer), so muss man ja nur sieben mal integrieren, um auf die Funktion selbst zu kommen. Wenn man also schon einmal integrieren kann, so sollte es leicht fallen, dies auch sieben mal (das selbe) hintereinander zu machen.

Beispiele für Differentialgleichungen sind aus der Elektrotechnik Spule und Kondensator.



mit i(t) dem zeitlichen Verlauf des Stroms, u(t) dem zeitlichen Verlauf der Spannung, di(t)/dt der zeitlichen Änderung des Stroms und du(t)/dt der zeitlichen Änderung der Spannung. L ist die Induktivität der Spule und C die Kapazität des Kondensators.

Aus der Mechanik sind zum Beispiel Kraft = Masse * Beschleunigung bekannt:

Die Beschleunigung a (Aufgrund der Krafteinwirkung F auf die Masse m) ist die 1. Ableitung der Geschwindigkeit, oder anders herum, die Geschindigkeit erhält man über Integration der Beschleunigung.

Und daraus ermittelt man den Weg über das Integral der Geschwindigkeit.

Für eine gleichbleibende Geschwindigkeit (z.B. beim Auto oder beim Raumfahrzeug) ergeben diese Rechenvorschriften ja gerade die bekannte Formeln x=v*t (die Beschleunigung ist ja 0).

Damit wir auch auf die Astronomie kommen: Die Bewegungsgleichungen (nicht nur die o.a. Grundgleichung F=m*a) sind genau das, was man zur Positionsbestimmung. Und über diese Bewegungsgleichungen können dann z.B. mit dem Zwei- oder Dreikörperproblem die Umlaufbahnen von Planeten etc. (idealisiert) bestimmt werden.

Grafische Interpretation der Integration

Betrachtet man eine Funktion y=f(x) und man will das Integral berechnen, so kann man sich das ganze mal grafisch vor Augen halten. Das Integral der Funktion f(x) ist ja nichts anderes als die Fläche unter dem Kurvenverlauf der Funktion f(x).

Hier im Bild ist das Integral der Kurve zwischen A und B die gefüllte rosa Fläche. Was hält uns denn noch von der Berechnung des Flächeninhalts ab? Hmm, eigentlich nichts, aber der Verlauf der Kurve ist so ein kleines Problem. Versuchen wir es mal mit einer Näherung der Fläche. Im unteren Bereich ist alles rechtwinklig. Da kann man also die Fläche im ersten Schritt über Quadrate rechnen. Schaut man sich das erst mal in der Grafik an:

Durch das Rechteck unten (rosa) wird die Fläche unter der Kurve (ja einigermaßen genau) angenähert. Es fehlt ja (nur) die blaue Fläche, um auf den exakten Wert zu kommen.

Genauer gehts mit Rechtecken, wenn das rosa Quadrat nicht so breit ist. Um dann die Fläche von A nach B zu bestimmen, kann man ja auch mehrere nebeneinanderliegende Quadrate in die Gesamtfläche einzechnen. Dann wird ja die blaue Fläche immer kleiner und die rosa Fläche nähert sich damit immer mehr der exakten Fläche an. Hier kann man das für drei Quadrate dann sehen:

Das Quadrat kann man auch anders ansetzen. Bei den bisherigen Quadraten ist die Fläche ja etwas zu klein, da das Quadrat unterhalb der Kurve angesetzt wurde. Selbstverständlich kann man dies auch oberhalb ansetzen. Im Ergebnis heißt das, das die angenäherte Fläche dann um den Betrag der blauen Fläche zu gross ist.

Ein weiterer Ansatz ist nicht ein Quadrat, sondern ein rechtwinkliges Trapez. Damit sieht die angenäherte Fläche folgendermaßen aus:

Was direkt auffällt ist, das die blaue Fläche sehr klein ist (viel kleiner als bei den Quadraten). Das heißt also, das der Fehler der angenäherten Fläche mit einem Trapez bei gleicher Breite AB sehr viel kleiner ist als bei den Quadraten. Daraus kann man direkt sagen, das bei der Flächenannäherung mit Trapezen die Gesamtfläche mit viel weniger Teilflächen angenähert werden kann, ohne das der Fehler größer als bei der Annäherung über viele kleine Quadrate ist. Anders ausgedrückt heißt das, das...

Wenn man nun die Breite der Teilflächen immer kleiner macht, so wird auch der Fehler immer kleiner. Daraus folgt, das man die echte Fläche entsprechend genau bestimmen kann. Schlußendlich heißt das, das wir das Integral der Funktion f(x) recht genau berechnen können.

Integration mit Einschrittverfahren

Um es einfach zu halten, werden einfach nur mal einige Einschrittverfahren vorgestellt. Einschrittverfahren bestimmen allgemein gesagt nur aus den Daten des aktullen Schritt die Daten des nächsten Schritt (kommen wir gleich noch dazu). Die oben gezeigten grafischen Interpretationen der Integration sind allesamt Einschrittverfahren.

Zuerst einmal muss man sich im Klaren darüber sein, das wir im Computer nicht mit kontinuierlichen Werten rechnen können, es sind immer diskrete Werte. Diskrete Werte heißt, das man die Funktion f(x) nicht an beliebigen x-Stellen bestimmt, sondern immer an bestimmten, nah beieinanderliegenden Stützstellen (x1, x2, ..., xn). Dazwischen haben wir keine Werte. Wählt man die Stützstellen aber dicht genug beieinander, so wird man kaum einen Unterschied zur kontinuierlichen Funktion feststellen. Die Dichte der Stützstellen ist ein anderes Thema (ich sag nur Schrittweitensteuerung).

Euler (vorwärts)

Ziel der Berechnung ist ja die Bestimmung des Integrals einer Funktion. Das Integral kann man - wie oben gesehen - grafisch interpretieren. Darüber hinaus benutzen wir diskrete Stützstellen.

Beim Euler-Verfahren nimmt man die diskreten Stützstellen äquidistant an, also eine konstante Schrittweite. Als Schrittweite definiert man z.B. Δt (siehe in der Grafik oben die Strecke AB) und rechnet wie mit den oben gezeigten Quadraten den Wert für den nächsten Stützpunkt aus dem bisherigen Wert aus.

mit

oder anders ausgedrückt

Da hier zu dem aktuellen Integralwert in ein Quadrat hinzuaddiert (Grundseitenlänge Δt, Höhe f(tn, in)), entspricht das der grafischen Interpretation mit dem Quadrat unter der Kurve.

Dieses Verfahren nennt man explizit, da der Wert von in+1 direkt aus den Werten zum Zeitpunkt n bestimmt werden. Das in+1 steht alleine auf der linken Seite, die rechte Seite der Gleichung hängt nur von Werten zum Zeitpunkt n ab.

Euler (rückwärts)

Bei dem zweiten Verfahren wird ähnlich vorgegangen, jedoch ergeben sich kleine Unterschiede. Die Äquidistanz der Stützstellen bleibt, jedoch wird hierbei das Quadrat rückwärts gebildet, d.h. aus dem Funktionswert des nächsten Stützpunktes wird die Höhe bestimmt.

mit demselben



Man sieht, das das grafische Gegenstück oben das Quadrat über der Kurve ist. Aber man sieht auch, das der Integralwert in+1 nicht alleine in der Formel steht. Links vom Gleichheitszeichen steht dieser zwar alleine, aber in der Funktion f fließt dieser mit ein. Das heißt, das der Wert nicht explizit, sondern implizit gegeben ist, weswegen dieser Ansatz auch implizites Euler-Verfahren genannt wird. Die Bestimmung des Integralwertes zum Zeitpunkt n+1 ist dann i.d.R. erst mit einer Zwischenrechnug möglich.

Den genauen Unterschied dieser beiden Methoden lass ich mal bei Seite, dieser kann in den Quellen nachgelesen werden.

Trapez

Auch das Trapez kann nach obiger Grafik recht einfach geschrieben werden:

Auch hier sieht man, das der Integralwert in+1 nicht explizit angegeben ist, sondern links und rechts in der Gleichung auftaucht.

Runge-Kutta

Das Runge-Kutta Verfahren (4. Ordnung) ist ein recht bekannter Algorithmus zur numerischen Bestimmung des Integrals.

...das beschreib ich dann irgend wann mal, damit ich jetzt nicht zu weit ausholen muss...

Zum Abschluss ein Beispiel

Ich nehm mal ein Beispiel aus der Elektrotechnik. Es ist eigentlich egal, da die Berechnungsgrundlage ja überall angewendet werden kann (hehe, Zwei- und Dreikörperproblem in der Astronomie). Betrachtet werden soll hier mal der Spannungsverlauf an einem Kondensator (und Widerstand, RC-Glied) bei aufgeprägtem Stromverlauf. Erinnern wir uns nochmal:


und für den Widerstand gilt:


...der rest folgt noch...


Weiterführende Links

Wikipedia - Differentialgleichungen
Wikipedia - Explizites Eulerverfahren
Wikipedia - Implizites Eulerverfahren
Wikipedia - Trapezintegration
Wikipedia - Runge-Kutta Verfahren
Wikipedia - Integralrechnung

Stand: 27. Februar 2010