system of equation Runge-Kutta 4th order for system of equation using matlab

Question

I need to do matlab code to solve the system of equation by using Runge-Kutta method 4th order but in every try i got problem and can't solve the derivative is (d^2 y)/dx^(2) +dy/dx-2y=0 , h=0.1 Y(0)=1 , dy/dx (0)=-2

{clear all, close all, clc
%{
____________________TASK:______________________
Solve the system of differential equations below 
in the interval 0<x<1, with stepsize h = 0.1.
y= y1                 y(0)=0
y3= 2y1-y2            y2(0)=-2   

_______________________________________________

%}
h = 0.1;
x    = 0:h:1
N = length(x);
y   = zeros(N,1);
y3    = zeros(N,1);

g = @(x, y, y1, y2) y1;
f = @(x, y, y1, y2) 2*y1-y2;
y1(1) = 0;
y2(1) =-2;

for i = 1:(N-1)
   k_1 = x(i)+y(i)
   k_11=g(x(i),y,y(i))
    k_2 = (x(i)+h/2)+(y(i)+0.5*h*k_1)
    k_22=g((x(i)+0.5*h),y,(y(i)+0.5*h*k_11))
    k_3 = (x(i)+h/2)+(y(i)+0.5*h*k_2)
    k_33=g((X(i)+0.5*h),y,(y(i)+0.5*h*k_22));
    k_4 = (x(i)+h)+(y(i)+h*k_33)
    k_44=g((x(i)+h),y,(y(i)+k_33*h));


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


end
Answer_Matrix = [x' y3 ];}

Answer 1

You used functions, that's not really necessary, but it might be easier that way to see the formula more clearly. In your functions however, you list arguments that used present in the function. That's not needed, and creates unwanted overhead.

In your initial conditions you should use y and y3, since that are the ones you use in the loop. Also in the first condition you've made a typo.

In the loop you forget to call the function f, and to update the y vector.

Making these changes in your code results in the following:

h = 0.1;
x = 0:h:1;
N = length(x);
y = zeros(N,1);
y3 = zeros(N,1);

g = @(y2) y2;
f = @(y1, y2) 2*y1-y2;

y(1) = 1;
y3(1) = -2;

for i = 1:(N-1)
    k_1 = f(y(i), y3(i));
    k_11 = g(y3(i));

    k_2 = f(y(i)+0.5*h*k_1, y3(i) +0.5*h*k_11);
    k_22 = g((y3(i)+0.5*h*k_11));

    k_3 = f(y(i)+0.5*h*k_2, y3(i) +0.5*h*k_22);
    k_33 = g((y3(i)+0.5*h*k_22));

    k_4 = f(y(i)+h*k_3, y3(i) +h*k_33);
    k_44 = g((y3(i)+h*k_33));

    y3(i+1) = y3(i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*h  ;
    y(i+1) = y(i) + (1/6)*(k_11+2*k_22+2*k_33+k_44)*h  ;
end
Answer_Matrix = [x' y];
% solution of DE is exp(-2x) and is plotted as reference
plot(x,y,x,exp(-2*x))

---

As mentioned before, you can also solve this without the use of functions:

h = .1;
x = 0:h:1;
N = length(x);

% allocate memory
y = zeros(N,1);
z = zeros(N,1);

% starting values
y(1) = 1;
z(1) = -2;

for i=1:N-1
    ky1 = z(i);
    kz1 = -z(i) + 2*y(i);

    ky2 = z(i) + h/2*kz1;
    kz2 = -z(i) - h/2*kz1 + 2*y(i) + 2*h/2*ky1;

    ky3 = z(i) + h/2*kz2;
    kz3 = -z(i) - h/2*kz2 + 2*y(i) + 2*h/2*ky2;

    ky4 = z(i) + h*kz3;
    kz4 = -z(i) - h*kz3 + 2*y(i) + 2*h*ky3;

    y(i+1) = y(i) + h/6*(ky1 + 2*ky2 + 2*ky3 + ky4);
    z(i+1) = z(i) + h/6*(kz1 + 2*kz2 + 2*kz3 + kz4);
end

% exp(-2*x) is solution of DE and is plotted as reference 
plot(x,y,x,exp(-2*x))