MATLAB System of ODEs Runge Kutta Fourth order

I have to solve the ODE y'' + 4y' + 3y = 0 y(0) = 2 y'(0) = -4
When tranformed into 2 first order ODE's it is dy1 = y2 and dy2 = -3y1 -4y2

I have to write a Matlab code which solves this using the fourth order Runge Kutta scheme on the domain [0,1]

This is what I have so far:

Code:

---

function [t,y] = System(f,n,h)
M = [0,1;-3,-4];
y(1, :) = [2, -4];
f =@(t,y)(M*y')';
t(1) = 0;

 
for j = 1 : n
k1 = h * f(t(j), y(j,:) );
k2 = h * f( t(j) + h/2, y(j,:) + k1/2 );
k3 = h * f( t(j) + h/2, y(j,:) + k2/2 );
k4 = h * f( t(j) + h, y(j,:) + k3 );
y(j+1,:) = y(j,:) + (k1 + 2*k2 + 2*k3 + k4) / 6;
end

---

I have tried so many different ways of writing this code and this one gives the least errors but does anyone have any idea where i can modify it so I don't get any errors and it solves properly?

Quote:

---

Originally Posted by CookieC

I have to solve the ODE y'' + 4y' + 3y = 0 y(0) = 2 y'(0) = -4
When tranformed into 2 first order ODE's it is dy1 = y2 and dy2 = -3y1 -4y2

I have to write a Matlab code which solves this using the fourth order Runge Kutta scheme on the domain [0,1]

This is what I have so far:

Code:

---

function [t,y] = System(f,n,h)
M = [0,1;-3,-4];
y(1, :) = [2, -4];
f =@(t,y)(M*y')';
t(1) = 0;

 
for j = 1 : n
k1 = h * f(t(j), y(j,:) );
k2 = h * f( t(j) + h/2, y(j,:) + k1/2 );
k3 = h * f( t(j) + h/2, y(j,:) + k2/2 );
k4 = h * f( t(j) + h, y(j,:) + k3 );
y(j+1,:) = y(j,:) + (k1 + 2*k2 + 2*k3 + k4) / 6;
end

---

I have tried so many different ways of writing this code and this one gives the least errors but does anyone have any idea where i can modify it so I don't get any errors and it solves properly?

---

I haven't checked the answer but:

Code:

---

function [t,y] = System(n,h)
  M = [0,1;-3,-4];
  y(1, :) = [2, -4];
  f =@(t,y)(M*y')';
  t(1) = 0;

 
  for j = 1 : n
    k1 = h * f(t(j), y(j,:) );
    k2 = h * f( t(j) + h/2, y(j,:) + k1/2 );
    k3 = h * f( t(j) + h/2, y(j,:) + k2/2 );
    k4 = h * f( t(j) + h, y(j,:) + k3 );
    y(j+1,:) = y(j,:) + (k1 + 2*k2 + 2*k3 + k4) / 6;
    t(j+1)=t(j)+h;
  end

---

runs

CB