## Runge-Kutta 4th order method

### Mariam Gasra (view profile)

on 5 May 2019
Latest activity Answered by David Wilson

on 6 May 2019

### David Wilson (view profile)

% It calculates ODE using Runge-Kutta 4th order method
% Author Ido Schwartz
clc; % Clears the screen
clear;
h=5; % step size
x = 0:h:100; % Calculates upto y(3)
Y = zeros(1,length(x));
y(1) = [-0.5;0.3;0.2];
% initial condition
F_xy = @(t,r) 3.*exp(-t)-0.4*r; % change the function as you desire
for i=1:(length(x)-1) % calculation loop
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; % main equation
end
display(Y(i+1));
if i run the programme i get answer =0;
how can i solve this problem if i have three initial condition -0.5 ,0.3 and 0.2
while x=0:5:100
and how i can plot the answer with respect to x?

### David Wilson (view profile)

on 6 May 2019

First up, you will need a much smaller step size to get an accurate solution using this explicit RK4 (with no error control). I suggest h = 0.05. Validate using say ode45 (which does have error control).
Then you will need to run your ode above three separate times, once starting from y(1) = 0.5, again with y(1) = 0.3, etc.
Then finally plot the result with plot(x,y,'o-').
h=0.5; % step size
x = 0:h:100; % Calculates upto y(3)
Y = zeros(1,length(x));
%y(1) = [-0.5;0.3;0.2];
y(1) = -0.5; % redo with other choices here.
% initial condition
F_xy = @(t,r) 3.*exp(-t)-0.4*r; % change the function as you desire
for i=1:(length(x)-1) % calculation loop
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; % main equation
end
% validate using a decent ODE integrator
tspan = [0,100]; y0 = -0.5;
[tx, yx] = ode45(F_xy, tspan, y0)
plot(x,y,'o-', tx, yx, '--')