clear all

a= 2;   % prirodni prirastaj plena
b=1.1;  % stopa umanjenja plena usled grabljivaca
c=1;    % natalitet grabljivaca
m=0.9;  % mortalitet grabljivaca
tend=50;
tp = 16; % u trenutku t_p desava se lov
%tp=15;
%tp=17;
alfa=0.5; % 50% ulovljenih jedinki jedne vrste (plena)
beta=0.5; % 50% ulovljenih jedinki druge vrste (grabljivaca)

P0=2;
G0=1;

% Stacionarno resenje:
%P0=m/c; % pocetna populacija plena (stac.)
%G0=a/b; % pocetna populacija grabljivaca (stac.)

f=@(t,x) [a*x(1)-b*x(1)*x(2); c*x(1)*x(2) - m*x(2)];

% Ako je logisticki rast zeceva
% K=5;
% f=@(t,x) [a*x(1)*(1-x(1)/K)-b*x(1)*x(2); c*x(1)*x(2) - m*x(2)];

[t,x]=ode45(f,[0,tend],[P0,G0]);

% Lov:
%[t1,x1]=ode45(f,[0,tp],[P0,G0]);
%[t2,x2]=ode45(f,[tp,tend],[x1(end,1)*(1-alfa),x1(end,2)*(1-beta)]);
% 
%t=[t1;t2];
%x=[x1;x2];

% Grafik promene broja jedinki obe vrste kroz vreme
figure(1)
plot(t,x(:,1),t,x(:,2))
xlabel('vreme')
ylabel('populacija')
legend('plen','grabljivice')

% Fazni dijagram
figure(2)
plot(x(:,1),x(:,2))
% plot(x(:,1),x(:,2),'o') % u slucaju stacionarnog resenja to ce biti samo jedna tacka, pa je crtamo kao kruzic
xlabel('plen')
ylabel('grabiljivice')
