Integration von Differentialgleichungen

Beispiele für Differentialgleichungen:
- radioaktiver Zerfall
  - Anzahl N der Kerne mit Zerfallskonstante λ
- radioaktive Zerfallskette
  - Teilchensorte 1 zerfällt in Sorte 2 mit λ₁
  - Sorte 2 zerfällt weiter mit λ₂
- angeregte Schwingung mit viskoser Reibung
  - Masse m schwingt an Feder mit Federkonstante c
  - Reibung ist proportional zur Geschwindigkeit v
  - äußere zeitabhängige Erregerkraft F(t)
Lösungsverfahren von Euler:
- einfachstes (und ungenauestes) Verfahren
- Grundidee
  - Ableitung ersetzen durch (endliche) Differenzen
  - ergibt normale Gleichungen für f
- im Beispiel "radioaktiver Zerfall"
  - bei gegebenem Startwert N(0) schrittweise lösen
  - hier sogar direkt auflösbar
- in Matlab
  - ```
  >> N0 = 10000;        % Anfangsmenge
  >> lambda = 0.1;      % Zerfallskonstante
  >> h = 0.1;           % Schrittweite
  >> n = 0:1000;        % 1000 Schritte
  >> t = h*n;           % Zeitwerte
  >> N = N0*(1-lambda*h).^n;   % Werte fuer N(t)
  >> plot(t, N);
```
Lösungsverfahren in Matlab:
- viele sehr weitentwickelte Verfahren (ode solver)
- einige Solver für schwierige (steife) Differentialgleichungen
- Auswahl des Solvers u.a. abhängig von
  - Genauigkeitsforderung
  - "Glattheit" der Gleichungen
  - Schwierigkeit der Gleichung (steif oder nicht steif)
- Daumenregel
  - erste Wahl: ode45
  - falls der versagt oder zu lange braucht: ode15s
  - wenn der auch nicht klappt: Experten fragen
Vorgehensweise am Beispiel "radioaktiver Zerfall"
- Gleichungen schreiben als y' = f(t, y)
- rechte Seite als Matlab-Funktion definieren (radioaktiv.m)
- Solver aufrufen
  - >> [t,y] = ode45(@radioaktiv,[0 100],[10000]);
- Parameter
  - betrachtetes Zeitintervall [0 100]
  - Anfangswert [10000]
- Ergebniswerte y und zugehörige Zeiten t
- Ergebnis in gewohnter Form plotten
  - >> plot(t, y);
Beispiel "Zerfallskette":
- System von Differentialgleichungen
- wie vorher, aber alle Größen als Vektoren
- Gleichungen als y' = f(t, y)
- rechte Seiten als Matlab-Funktion radiokette.m
- Anfangsbedingungen z.B.
  - 10000 Teilchen von Sorte 1
  - keine Teilchen von Sorte 2
- Lösen
  - >> [t,y] = ode45(@radiokette,[0 100],[10000 0]);
- Spalten 1, 2 von y enthalten Zeitwerte von y₁, y₂, daher plotten mit
  - ```
  >> plot(t, y(:,1), t, y(:,2))
  >> legend('N_1(t)', 'N_2(t)', 0);
```
Beispiel "angeregte Schwingung mit viskoser Reibung"
- Problem: Gleichung zweiter Ordnung
- Trick: Geschwindigkeit als zweite Variable einführen
- damit Umschreiben in 2 Gleichungen 1. Ordnung
- rechte Seiten als Matlab-Funktion erzwungen.m
- Anfangsbedingungen z.B. Masse in Ruhe am Ursprung
- Lösen
  - >> [t, y] = ode45(@erzwungen, [0 30], [0 0]);
- Darstellung
  - zweite Spalte ist Geschwindigkeit
  - zusätzlich äußere Anregung plotten
  - ```
  >> F = cos(t);
  >> plot(t, y(:,1), t, y(:,2), t, F)
  >> legend('Auslenkung', 'Geschwindigkeit', 'Kraft', 'Location', 'NorthWest');
```
Übergabe von Parametern:
- Problem:
  - Parameter (m, b, c, A, omega) fest in Funktion erzwungen
  - sollen direkt angebbar sein
- Parameter als zusätzliche Argumente der Funktion (erzwungenp.m)
- Erzeugen einer Hilfsfunktion myfunc, die nur von t und y abhängt
  - >> myfunc = @(t, y) erzwungenp(t, y, 1 , 0.286, 4.68, 0.706, 1);
- Aufruf des Solvers mit der Hilfsfunktion
  - >> [t, y] = ode45(myfunc, [0 40], [1 0]);
Aufgaben:
- Aufgabe 12