Цитата:
Сообщение от марианна
решить уравнение методом Ранге-Кутта 4-ого порядка
|
По-моему так, хотя результат получился несколько странный - почти чистая прямая. Но вроде ошибок не вижу.
Пояснения: исходное уравнение второго порядка приведено к системе уравнений первого порядка путём подстановки (dy/dt)=z. Тогда получаем систему:
dz/dt = (t³/5)*z - (t/5)*Sin(t)*y = f(t,y,z)
dy/dt = z = g(t,y,z);
Далее пишем программу и решаем, как это описано тут:
Метод Рунге — Кутты — Википедия
Код:
Const
t0=0.0;
y0=0.0;
z0=-3.41;
h=0.02;
t_max=0.44;
Var
T,Y,Z,Y1,Z1:real;
i,N:integer;
Function f(t,y,z:real):real;
begin
f:=t*t*t/5*z-t/5*sin(t)*y;
end;
Function g(t,y,z:real):real;
begin
g:=z;
end;
Procedure RK(t,y,z:real; var Uy:real; var Uz:real);
var k1,k2,k3,k4,q1,q2,q3,q4:real;
begin
k1:=h*f(t,y,z);
q1:=h*g(t,y,z);
k2:=h*f(t+h/2,y+q1/2,z+k1/2);
q2:=h*g(t+h/2,y+q1/2,z+k1/2);
k3:=h*f(t+h/2,y+q2/2,z+k2/2);
q3:=h*g(t+h/2,y+q2/2,z+k2/2);
k4:=h*f(t+h,y+q3,z+k3);
q4:=h*g(t+h,y+q3,z+q3);
Uz:=z+(k1+2*k2+2*k3+k4)/6;
Uy:=y+(q1+2*q2+2*q3+q4)/6;
end;
Begin
N:=Round((t_max-t0)/h);
T:=t0;
Y:=y0;
Z:=z0;
Writeln('t= ',T:4:2,' y= ',Y:8:5);
For i:=1 to N do
begin
RK(T,Y,Z,Y1,Z1);
T:=T+h;
Y:=Y1;
Z:=Z1;
Writeln('t= ',T:4:2,' y= ',Y:8:5);
end;
Readln
End.
P.S. Сейчас для проверки взял функцию y=-3.41t и подставил в исходное уравнение. Тогда слева - ноль, а оба слагаемых действительно совпали до пятого знака (при t=0.2 0.005456 и 0.005419, соответственно). Так что, похоже, всё правильно.