Dany jest układ dwóch zagadnień początkowych (Cauchy):
\(\displaystyle{ x^{'} = f(t,x(t),y(t)), \ \ x(a) = \alpha,}\)
\(\displaystyle{ y^{'} = g(t,x(t),y(t)), \ \ y(a) = \beta.}\)
Schemat Runge-Kutta rzędu 4 dla tego układu ma postać:
\(\displaystyle{ x_{n+1} = x_{n} +\frac{h}{6}(k_{1} + 2k_{2}+2k_{3}+ k_{4}),}\)
\(\displaystyle{ y_{n+1} = y_{n} +\frac{h}{6}(c_{1} + 2c_{2}+2c_{3}+ c_{4}),}\)
\(\displaystyle{ k_{1} = f(t_{n},x_{n},y_{n}), }\)
\(\displaystyle{ c_{1}= g(t_{n},x_{n},y_{n}), }\)
\(\displaystyle{ k_{2}= f\left(t_{n}+\frac{h}{2}, x_{n}+\frac{h}{2}k_{1}+\frac{h}{2}c_{1} \right),}\)
\(\displaystyle{ c_{2}= g \left (t_{n}+\frac{h}{2}, x_{n}+\frac{h}{2}k_{1} +\frac{h}{2}c_{1}\right)),}\)
\(\displaystyle{ k_{3} = f\left( t_{n} + \frac{h}{2}, x_{n}+\frac{h}{2}k_{2}, y_{n}+ \frac{h}{2}c_{2} \right),}\)
\(\displaystyle{ c_{3} = g\left( t_{n} + \frac{h}{2}, x_{n}+\frac{h}{2}k_{2}, y_{n}+ \frac{h}{2}c_{2} \right),}\)
\(\displaystyle{ k_{4}=f(t_{n}+h, x_{n}+ hk_{3}, y_{n}+hc_{3}),}\)
\(\displaystyle{ c_{4}=g(t_{n}+h, x_{n}+ hk_{3}, y_{n}+hc_{3}).}\)
Metodą możemy rozwiązywać numerycznie równania różniczkowe zwyczajne drugiego rzędu, sprowadzając je do układu równań rzędu pierwszego.
Przykład
Rozwiązać równanie
\(\displaystyle{ y^{''}(t) - 3y^{'}(t) +2y = 6e^{3t}, \ \ 0\leq t \leq 1}\)
z warunkami początkowymi: \(\displaystyle{ y(0)=1, \ \ y^{'}(0)=-1. }\)
Stosując podstawienia:
\(\displaystyle{ x^{'}(t) = 3x- 2y +6e^{3t}, \ \ x(0)=-1,}\)
\(\displaystyle{ y^{'}(t)= x, \ \ y(0) =1, \ \ 0\leq t \leq 1, }\)
sprowadzamy rozwiązanie równania rózniczkowego rzędu drugiego do rozwiązania układu równań pierwszego rzędu.
Program w MATLAB
Kod: Zaznacz cały
function ukl_rk4(f,g,a,b,x0,y0,n)
% Rozwiązanie układu równań: x'=f(t,x,y), x(a)=x0,
% y'=g(t,x,y), y(a)=y0 metodą Runge-Kutta rzędu 4.
fprintf('\n')
disp(' Rozwiązanie układu równań rózniczkowych metodą Runge-Kutta rzędu 4 ')
h=(b-a)/n;
x=x0;
y=y0;
disp('________________________________')
disp(' t x y ')
disp('________________________________')
fprintf('\n')
fprintf('%8.2f %10.3f %10.3f\n', a, x, y)
for i=1:n
% Runge-Kutta rzędu 4
t=a+(i-1)*h;
k1=feval(f,t,x,y);
c1=feval(g,t,x,y);
k2=feval(f,t+h/2,x+h*k1/2,y+h*c1/2);
c2=feval(g,t+h/2,x+h*k1/2,y+h*c1/2);
k3=feval(f,t+h/2,x+h*k2/2,y+h*c2/2);
c3=feval(g,t+h/2,x+h*k2/2,y+h*c2/2);
k4=feval(f,t+h,x+h*k3,y+h*c3);
c4=feval(g,t+h,x+h*k3,y+h*c3);
x=x+h*(k1+2*k2+2*k3+k4)/6;
y=y+h*(c1+2*c2+2*c3+c4)/6;
t=t+h;
fprintf('%8.2f %10.3f %10.3f\n',t,x,y)
end
Kod: Zaznacz cały
function f=f1(t,x,y)
f=3*x-2*y+6*exp(3*t);
function g=g1(t,x,y)
g=x;
>> ukl_rk4('f1','g1',0,1,-1,1,10)
Rozwiązanie układu równań rózniczkowych metodą Runge-Kutta rzędu 4
________________________________
t x y
________________________________
0.00 -1.000 1.000
0.10 -0.763 0.909
0.20 -0.142 0.860
0.30 1.082 0.901
0.40 3.223 1.107
0.50 6.735 1.591
0.60 12.257 2.520
0.70 20.694 4.139
0.80 33.311 6.798
0.90 51.878 10.998
1.00 78.850 17.452