clear all
close all
alpha = 1;
beta = 0.002;
gamma = 1;
delta = 0.001;

f=@(t,u) [alpha*u(1)-beta*u(1)*u(2); 
          -gamma*u(2)+delta*u(1)*u(2)];

u0=[500;50];   % initial population 500 rabbits, 50 foxes
Tspan=linspace(0,20,520); % time interval

[t,u]=ode45(f,Tspan,u0);

for i=1:length(Tspan)
  
  plot(t(1:i),u(1:i,1),'b',t(1:i),u(1:i,2),'r','LineWidth',2)
  axis([ 0 20 0 5000])
  xlabel('t')
  ylabel('number of animals', 'fontsize',16)
  legend('rabbits','foxes')
  
  pause(0.0003)
   
end

figure(2)

for i=1:length(Tspan)
  plot(u(1:i,1),u(1:i,2),'r','LineWidth',2)
  hold on
  plot(u(i,1),u(i,2),'bo','LineWidth',2)
  title(['Time t=', num2str(t(i))],'fontsize', 16)
  xlabel('rabbits','fontsize',16)
  ylabel('foxes','fontsize',16)
  
  axis([ 0 5000 0 5000])
  
  pause(0.0001)
    
end
