Runge Kutta 4th order

Question

Melda Harlova 2019-5-8

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/460919-runge-kutta-4th-order

回答： Meysam Mahooti 2021-5-5

采纳的回答： Erivelton Gualter

Hello,

Here is the task that i have to solve:

y1' = y2

y2'=f(x,y1,y2) with y1(0)=0 and y2(0)=y20

where f(x,y1,y2) = -axy2-y1, a=0.03 and y20 = 0.25

and here is my matlab code:

--

function main
h = pi/10;
x = 0 : h : 2*pi;
n = length(x);
y(1,1)=0;  y(2,1) = 0,25;
for k = 1: (n-1)
                k1 = h * rung(x(k), y(1:2,k)); 
                k2 = h * rung(x(k) + h/2, y(1:2,k) + k1/2);
                k3 = h * rung(x(k) + h/2, y(1:2,k) + k2/2);
                k4 = h * rung(x(k) + h, y(1:2,k) + k3);
                y(1:2,k+1)=y(1:2,k) + (k1 + 2*k2 + 2*k3 + k4)/6;
end
figure(1), plot(x, y(2,:)) %everything from second row
figure(2), plot(x, y(1,:), x, ??, '*')
function f=rung(x,y)
f(1,1) = y(2); 
f(2,1) = -0,03*x*y(2) – y(1);

So, my question is - Is this correct or not? And what would be the best option for h and x.

And how can i found and use numerical solution? Because i know that instead of the question marks here figure(2), plot(x, y(1,:), x, ??, '*') should be the exact numerical solution.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Erivelton Gualter 2019-5-8

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/460919-runge-kutta-4th-order#answer_374137

Your code seems to have some issues. Try this one:

h = pi/10;
x = 0 : h : 2*pi;
n = length(x);
y = zeros(2,n);
y(:,1) = [0, 0.25];
for k = 1: (n-1)
                k1 = h * rung(x(k), y(:,k)); 
                k2 = h * rung(x(k)+h/2, y(:,k) + k1/2);
                k3 = h * rung(x(k)+h/2, y(:,k) + k2/2);
                k4 = h * rung(x(k)+h, y(:,k) + k3);
                y(:,k+1) = y(:,k) + (k1 + 2*k2 + 2*k3 + k4)/6;
end
subplot(211); plot(x, y(1,:))
subplot(212); plot(x, y(2,:))
function f = rung(x,y)
    f(1,1) = y(2); 
    f(2,1) = -0.03*x*y(2) - y(1);
end