hi I am trying to use this code to solve RK4 equation and I keep receiving this error . what should I do?

Question

ela mti 2020-12-24

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/701587-hi-i-am-trying-to-use-this-code-to-solve-rk4-equation-and-i-keep-receiving-this-error-what-should

评论： James Tursa 2020-12-25

    function [tout,Aout] = ivp_RK4(t0,A0,h,n_out,i_out)
    tout = zeros(n_out+1,1); tout(1) = 0;
    Aout = zeros(n_out+1,1); Aout(1) = 1;
    t = t0; A = A0;
    for j=2:n_out+1;
        i=1:i_out;
        k1 = A(i)-exp(-2*t);
        k2 =(A(i)+h/2)*k1-exp(-2*(t+h/2)) ;
        k3 =(A(i)+h/2)*k2-exp(-2*(t+h/2));
        k4 = (A(i)+h/2)*k3-exp(-2*(t+h/2));
        A = A + (h/6)*(k1 + 2*k2 + 2*k3 + k4);
        t = t + h;
    end
    tout(j) = t;
    Aout(j) = A;
    plot(t,A)
    end

Index exceeds the number of array elements (1).

Error in ivp_RK4 (line 7)

k1 = A(i+1)-exp(-2*t);

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

请先登录，再进行评论。

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

Answer 1

James Tursa 2020-12-24

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/701587-hi-i-am-trying-to-use-this-code-to-solve-rk4-equation-and-i-keep-receiving-this-error-what-should#answer_583387

在 MATLAB Online 中打开

Your implementation has several errors:

Looks like you are creating i as a vector and then using that for indexing in A(i)
You are not using the k's properly
You are not calculating k4 correctly (it should be a full step not a half step)
You are using A as both a vector and a scalar when updating
You are not saving the intermediate results in Aout and tout inside the loop

To make things simpler and easier to read and debug, I would suggest making a function handle for the derivative. Then construct your code using this function handle, indexing both Aout and tout inside the loop. E.g.,

f = @(t,A) A-exp(-2*t);
tout = zeros(n_out+1,1); tout(1) = t0;
Aout = zeros(n_out+1,1); Aout(1) = A0;
for i=1:n_out
    k1 = f(tout(i)    , Aout(i)       );
    k2 = f(tout(i)+h/2, Aout(i)+k1*h/2);
    k3 = f(tout(i)+h/2, Aout(i)+k2*h/2);
    k4 = f(tout(i)+h  , Aout(i)+k3*h  );
    Aout(i+1) = Aout(i) + (h/6)*(k1 + 2*k2 + 2*k3 + k4);
    tout(i+1) = tout(i) + h;
end