Using the Runge-Kutta integration method in a system

Question

h=0.005;                                            
x = 0:h:40;                                         
y = zeros(1,length(x)); 
y(1) = 0;                                         
F_xy = ;                    

for i=1:(length(x)-1)                              
    k_1 = F_xy(x(i),y(i));
    k_2 = F_xy(x(i)+0.5*h,y(i)+0.5*h*k_1);
    k_3 = F_xy((x(i)+0.5*h),(y(i)+0.5*h*k_2));
    k_4 = F_xy((x(i)+h),(y(i)+k_3*h));

    y(i+1) = y(i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*h;  
end

I have the following code, I think it's right. I know there's parts missing on the F_xy because this is my follow up question.

I have dx/dt = = −x(2 − y) with t_0 = 0, x(t_0) = 1

and dy/dt = y(1 − 2x) with t_0 = 0, y(t_0) = 2.

My question is that I don't know how to get these equations in to the code. All help appreciated

Lutz Lehmann Lutz Lehmann · Accepted Answer · 2018-03-18T14:16:57

You are using both t and x as independent variable in an inconsistent manner. Going from the actual differential equations, the independent variable is t, while the dependent variables of the 2-dimensional system are x and y. They can be combined into a state vector u=[x,y] Then one way to encode the system close to what you wrote is

h=0.005;                                            
t = 0:h:40;   
u0 = [1, 2]                                       
u = [ u0 ]                                         
function udot = F(t,u)
     x = u(1); y = u(2);
     udot = [ -x*(2 - y),  y*(1 - 2*x) ]
end                 

for i=1:(length(t)-1)                              
    k_1 = F(t(i)      , u(i,:)        );
    k_2 = F(t(i)+0.5*h, u(i,:)+0.5*h*k_1);
    k_3 = F(t(i)+0.5*h, u(i,:)+0.5*h*k_2);
    k_4 = F(t(i)+    h, u(i,:)+    h*k_3);

    u(i+1,:) = u(i,:) + (h/6)*(k_1+2*k_2+2*k_3+k_4);  
end

with a solution output

Using the Runge-Kutta integration method in a system

2 Answers