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

    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);

Answer 1

[James Tursa]

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

After the loop finishes, the answers are in the Aout and tout vectors. So you would

plot(tout,Aout)

I am confused what the i_out is supposed to be used for. Can you clarify?

[ela mti]

first of all thank you so much !

i_out is frequency to output x and y to [xout, yout]

[ela mti]

actually A should be reduced over the time I tried to fix the code with your suggestion code but it cant be right

[James Tursa]

What is the differential equation you are trying to solve?