Επίλυση διαφορικής εξίσωσης με τη μέθοδο Runge-Kutta Octave/Matlab
Η Runge-Kutta είναι μία προσεγγιστική μέθοδος επίλυσης διαφορικών εξισώσεων. Στην ουσία πρόκειται για μια μεγάλη οικογένεια μεθόδων με πολλές παραλλαγές.
Η μέθοδος είναι περισσότερο ακριβής από τη μέθοδο Euler (Δείτε Υπολογιστικές μέθοδοι για επίλυση διαφορικών εξισώσεων) και στην πράξη χρησιμοποιείται πολύ περισσότερο.
Εδώ θα εξετάσουμε ένα απλό παράδειγμα της κλασικής μεθόδου Runge-Kutta (RK4). Έστω η διαφορική εξίσωση:
$$ \frac{\mathrm{d}\,y}{\mathrm{d}\,x}= x\,y $$
με αρχικές συνθήκες x=0, y=1. Να υπολογιστεί πίνακας τιμών της συνάρτησης y=f(x) στο διάστημα [0,1].
Η μέθοδος στηρίζεται στον αναδρομικό τύπο:
$$ \begin{align} y_{n+1} &= y_n + \tfrac{1}{6} h\left(k_1 + 2k_2 + 2k_3 + k_4 \right)\\ x_{n+1} &= x_n + h \\ \end{align} $$
όπου:
$$ \begin{align} k_1 &= f(x_n, y_n) \\ k_2 &= f(x_n + \tfrac{1}{2}h , y_n + \tfrac{h}{2} k_1) \\ k_3 &= f(x_n + \tfrac{1}{2}h , y_n + \tfrac{h}{2} k_2) \\ k_4 &= f(x_n + h , y_n + hk_3). \end{align} $$
Ένα απλό πρόγραμμα που επιλύει τη διαφορική εξίσωση είναι:
clear;
x0 = input("initial value of x = ");
y0 = input("initial value of y = ");
x1 = input("upper limit of x = ");
N = input("number of steps = ");
function z = f1(x, y)
z = x*y;
end
x = zeros(N, 1);
y = zeros(N, 1);
h = (x1-x0) / N;
x(1) = x0;
y(1) = y0;
for i=1:N
k1 = h * f1(x(i),y(i));
k2 = h * f1(x(i) + h/2.0, y(i) + 0.5 * k1);
k3 = h * f1(x(i) + h/2.0, y(i) + 0.5 * k2);
k4 = h * f1(x(i) + h, y(i) + k3);
y(i+1) = y(i) + (k1 + 2.0*k2 +2.0*k3 + k4 ) / 6.0;
x(i+1) = x(i) + h;
end
[x y]
σχολιασμοί, εξωτερικοί σύνδεσμοι, βοήθεια, ψηφοφορίες, αρχεία, κτλ.
Εκπαιδευτικό υλικό από τον
Αθανάσιο Σταυρακούδη
σας παρέχετε κάτω από την άδεια
Creative Commons Attribution-NonCommercial-ShareAlike 4.0 License.
Σας παρακαλώ να ενημερωθείτε για κάποιους επιπλέον περιορισμούς
http://stavrakoudis.econ.uoi.gr/stavrakoudis/?iid=401.