I'm using a Runge-Kutta of order 4 to resolve a system of differential equations. My codes are:
1
function [t,y]=RK(fun,tspan,y0,n,A,b,c)
h=(tspan(2)-tspan(1))/n;
y(:,1)=y0;
t(1)=tspan(1);
s=length(b);
for j=1:n
k(:,1)=feval(fun,t(j),y(:,j));
for i=2:s
k(:,i)=feval(fun,t(j)+c(i)*h,y(:,j)+h*k(:,1:i-1)*A(i,1:i-1)');
end
y(:,j+1)=y(:,j)+h*k*b;
t(j+1)=t(j)+h;
end
end
then
2
function dy=pred_prey(t,y)
a=2/3;
d=4/3;
r=(2/3)*(1/(exp(-0.5*t)+1)-0.5)*(1+(sin(t)/t));
mu=13/20-(3/5)*exp(-(3/t));
dy(1)=y(1)*r-a*y(1)*y(2);
dy(2)=-mu*y(2)+d*y(1)*y(2);
dy=dy';
end
and
3
function []=driver_pred_prey()
close all;
tspan=[0 5000];
n=1000;
fun='pred_prey';
y0=[2 5]';
A=[0 0 0 0; 1/2 0 0 0 ; 0 1/2 0 0; 0 0 1 0];
b=[1/6 1/3 1/3 1/6]';
c=[0 1/2 1/2 1]';
[t,y]=RK(fun, tspan, y0, n, A, b, c);
figure (1)
plot(t, y(1,:), 'k-'), hold
plot(t, y(2,:), 'k-.')
legend('resources','population')
figure (2)
plot(y(1,:),y(2,:))
title('orbits')
end
When I run driver_pred_prey everything is okay, except that no lines are showed in the plots.
