function preda_predatore(x0, y0, N, DT, A, B, C, D, E)
% programma preda_predatore

warning('off');

if nargin < 9, 
    disp(sprintf(' inserire x0, y0, N, DT, A, B, C, D, E'))
    return
end

P1=[0 0]'
P2=[B/D (A-E*B/D)/C]'
P3=[A/E 0]'
x=zeros(N/DT,1);
y=zeros(N/DT,1);
x(1)=x0;
y(1)=y0;

for i = 1:N/DT
    x(i+1,1)=x(i)+(A*x(i)-C*x(i)*y(i)-E*x(i)^2)*DT;
    y(i+1,1)=y(i)+(-B*y(i)+D*x(i)*y(i))*DT;
end


timef=size(x,1)*DT;
t=(0:DT:timef-DT)';

figure(1);clf;
 plot(x, y);
 hold on;
 plot(P1(1),P1(2), 'xr')
 hold on
 plot(P2(1),P2(2), '*r')
 hold on
 plot(P3(1),P3(2), '+r')
s_title=sprintf('Traiettoria popolazioni');
title(s_title);
ylabel('Predatori')
xlabel('Prede');
xlim([0 max(x)*1.2]);
ylim([0 max(y)*1.2]);


figure(2);clf;
  subplot(2,1,1)
  plot(t,x)
  ylabel('Popolazione prede')
  xlabel('tempo');
  xlim([t(1) t(end)]);
  ylim([0 max(x)*1.2]);
  subplot(2,1,2)
  plot(t,y)
  ylabel('Popolazione predatori')
  xlabel('tempo');
  xlim([t(1) t(end)]);
  ylim([0 max(y)*1.2]);
s_title=sprintf('Dinamica popolazione');
title(s_title);


return