Pliki do wykorzystania w poniższym ćwiczeniu można pobrać za pomocą poniższych linków:
Celem ćwiczenia jest zastosowanie metody Eulera oraz metody Rugego-Kutty 4 rzędu do numerycznego rozwiązania równań ruchu dynamiki Newtona. Jako przykład takiego zagadnienia posłuży nam wahadło matematyczne.
Ruch wahadła matematycznego najwygodniej jest opisać w układzie współrzędnych biegunowych związanych z jego punktem zaczepienia. Otrzymamy wtedy równanie różniczkowe wraz z warunkami początkowymi: \[ \left\{ \begin{array}{l} \dfrac{d^2\alpha}{dt^2} = -\dfrac{g}{l}\sin(\alpha) \\ \alpha(t_0) = \alpha_0 \\ \dfrac{d\alpha}{dt}(t_0) = \omega_0 \end{array} \right. \] gdzie:
Równanie to możemy sprowadzić do układu równań różniczkowych zwyczajnych pierwszego rzędu za pomocą podstawienia: \[ \dfrac{d\alpha}{dt} = \omega \] Układ równań ma teraz postać: \[ \tag{*} \left\{ \begin{array}{l} \dfrac{d\omega}{dt} = -\dfrac{g}{l}\sin(\alpha) \\ \dfrac{d\alpha}{dt} = \omega \\ \alpha(t_0) = \alpha_0 \\ \omega(t_0) = \omega_0 \end{array} \right. \]
Mamy układ dwóch równań różniczkowych pierwszego rzędu: \[ \left\{ \begin{array}{l} \dfrac{d\omega}{dt} = F_1(\alpha, \omega, t) \\ \dfrac{d\alpha}{dt} = F_2(\alpha, \omega, t) \\ \end{array} \right. \] Szukanymi funkcjami są \(\omega = \omega(t)\) oraz \(\alpha = \alpha(t)\). Układ ten można rozwiązać metodą Eulera. Jedna iteracja całkowania z krokiem \(h\) będzie miała postać: \[\begin{equation} \left\{ \begin{array}{l} \omega(t_{i+1}) = \omega(t_i) + h \cdot F_1(\alpha(t_i), \omega(t_i), t_i) \\ \alpha(t_{i+1}) = \alpha(t_i) + h \cdot F_2(\alpha(t_i), \omega(t_i), t_i) \\ \end{array} \right. \end{equation}\]
Dla wahadła opisanego układem równań (*):
void rhs_fun(double t, double *X, double *F);
która oblicza wartości prawych stron równań różniczkowych. Argumenty funkcji to:
t
- zmienna niezależna (czas),X
- tablica zmiennych zależnych (\(\alpha\) i \(\omega\)),F
- tablica do której zapisane zostaną obliczone prawe strony równań różniczkowych.void veuler(double t, double *X, double h, int n,
void (* fun)(double, double *,double *), double *X1);
która wykonuje jeden krok całkowania układu równań różniczkowych zwyczajnych pierwszego rzędu metodą Eulera. Argumenty funkcji to:
t
- zmienna niezależna,X
- tablica wartości zmiennych zależnych w kroku t
,h
- krok całkowania,n
- rozmiar tablicy,fun
- wskaźnik do funkcji obliczającej prawe strony równań,X1
- tablica do której zapisane zostaną wartości zmiennych zależnych w kroku t + h
.Napisz program, który używając metody Eulera wyznacza zależności kąta wychylenia wahadła \(\alpha\) oraz prędkości kątowej \(\omega\) od czasu dla \(t \in \left[0, \ldots 10\right]\).
Narysuj wykres trajektorii układu w przestrzeni fazowej (\(\alpha-\omega\)).
Powtórz obliczenia korzystając z metody Rungego-Kutty 4 rzędu, która jest zaimplementowana w bibliotece rk4.cpp
. Odpowiednia funkcja nazywa się vrk4
a jej nagłówek jest analogiczny do nagłówka funkcji veuler
.
Wyznacz zależność energii całkowitej wahadła od czasu \(E(t)\). Energia całkowita wahadła wyraża się wzorem: \[ E = \frac{ml^2}{2}\left(\frac{d\alpha}{dt}\right)^2 + mgl(1-\cos(\alpha)) \] Uwaga: Przy braku dyssypacji, energia mechaniczna powinna być stała.
Powtórz obliczenia dla różnych kroków czasowych.