% LEP P1-3: Predator-Prey Relation
%

Vspan = [0 500]; 
y0 = [500 200]; 
k1=0.02;
k2=0.00004;
k3=0.0004;
k4=0.04;
kd=0.005;
[v,y]=ode45(@(v,y)ODEfun(v,y,k1,k2,k3,k4,kd),Vspan,y0);
tiledlayout(2,2)
set (gcf,'Position',[0,0.1,600,400])
nexttile
plot(v,y(:,1),'LineWidth',2)
grid on
xlabel('time (days)')
ylabel('No of Rabbits')
title('No. of Rabbits vs. time')
xlim([0 500])
ylim([0 1000])
nexttile
plot(v,y(:,2),'LineWidth',2)
grid on
xlabel('time (days)')
ylabel('No of Foxes')
title('No. of Foxes vs. time')
xlim([0 500])
ylim([0 8000])
nexttile
plot(y(:,1),y(:,2),'LineWidth',2)
grid on
xlabel('No of Rabbits')
ylabel('No of Foxes')
title('No. of Rabbits vs. Foxes')
xlim([0 550])
ylim([0 10000])

%

function dYfuncvecdV = ODEfun(V,Yfuncvec,k1,k2,k3,k4,kd)
x = Yfuncvec(1); 
y = Yfuncvec(2);
% Differential equations
f1 = k1*x - k2*x*y-kd*x;
f2 = k3*x*y - k4*y;
dYfuncvecdV =[f1; f2]; 
end
