%  duas_especies.m
% 
%  Modelo para populações que interagem (duas espécies)  
%
% O programa resolve um modelo de dinâmico de duas populações que interagem. 
% Os resultados são apresentados em gráficos com o plano fase e a evolução das 
% duas espécies ao longo do tempo.  



p=0;q=0;                           
%
%Definição do modelo
%
disp('Introduza as fórmulas que definem os tamanhos das duas espécies,') 
disp('usando a notação "p" para uma das espécies e "q" para a outra:')
disp('(Exemplo: next_p = p+.4*p*(1-p/10)-.05*p*q')
disp('          next_q = q-.2*q+.08*p*q)')
%
next_p=input ('next_p = ','s');
p=eval(next_p);                    
%
next_q=input ('next_q = ','s');
q=eval(next_q);                    
% 
disp(' ')
disp('Introduza os limites para "p" e "q" na construção dos graficos')
disp('(por defeito é considerado [pmin pmax qmin qmax]=[0 10 0 10])')
%
limits=input('[pmin pmax qmin qmax] = ');
if isempty(limits) limits=[0 10 0 10]; end;
%
disp(' ')
disp('Introduza o número de passos temporais: (por é considerado defeito n=100)')
n=input('n= ');
if isempty(n) n=100; end;
%
disp(' ')
disp('Prima qualquer tecla para iniciar.') 
disp('Posicione o cursor no plano fase para escolher os valores iniciais ')
disp('dos tamanhos das populações e clique no botão da esquerda do rato. ')
disp('Clique no botão da direita do rato para finalizar.')
%
pause
%
% início de uma nova figura no plano fase
pp=figure('Units','normalized','Position',[.01,.3,.485,.6],...   
   'Name','Plano fase');     
hold on;
axis(limits); grid;      
xlabel('Espécie  p');ylabel('Espécie q');
title(['next\_p=',next_p,',  next\_q=',next_q]);
%
% início de uma nova figura no gráfico temporal
tp=figure('Units','normalized','Position',[.505,.3,.485,.6],...   
   'Name','grafico temporal');                  
xlabel('tempo');ylabel('populacoes');         
title(['next\_p=',next_p,',  next\_q=',next_q]);
%
figure(pp)                 % vai para o plano fase
[p,q,button]=ginput(1);    % capta os valores iniciais das populações  
%
while button==1            % ciclo até que o botão da esquerda seja primido
   pops=[p;q];
   for i=1:n               %  ciclo até ao número máximo de iterações
      newp=eval(next_p);
      newq=eval(next_q);
      p=newp;
      q=newq;
      pops=[pops,[p;q]];
   end
   plot(pops(1,:),pops(2,:),pops(1,:),pops(2,:),'o');   % gráfico do plano fase
   %
   figure(tp)                                    % vai para o gráfico temporal
   plot([0:n], pops(1,:), [0:n], pops(2,:));     % gráfico da última simulação
   xlabel('tempo');ylabel('populações');         
   title(['next\_p=',next_p,' (azul),  next\_q=',next_q,' (verde)']);
   %
   figure(pp)                                    % regresso ao plano fase
   [p,q,button]=ginput(1);                       %  capta os  valores iniciais das populações seguintes 
end
%
hold off                                % fim do plano fase

