%Parámetros del modelo de Gause, capacidades de carga y tasas de
%crecimiento
K1 = 1; 
K2 = 1; 
r1 = 1; 
r2 = 0.1;
a = 0.5; %coeficiente de competencia de 2 contra 1
b = 0.5; %coeficiente de competencia de 1 contra 2
T = 500; %numero de generaciones
N = [0.1 0.1]; %poblaciones iniciales
dN = @(t,N) [r1*N(1)*(K1 - N(1) - a*N(2))/K1 ; r2*N(2)*(K2 - N(2) - b*N(1))/K2];
dt = 1:T;
[t,N] = ode45(dN,dt,N);

figure(1) %definimos una nueva figura, que va a mostrar el curso temporal de las poblaciones

plot(dt,N(:,1))
title ('Poblaciones a lo largo del tiempo')

hold on %indica que vamos a superponer distintas gráficas
plot(dt,N(:,2))

line([1 T], [K1,K1], 'color', 'red') %marca las capacidades de carga para comparar con la población
line([1 T],[K2,K2], 'color', 'black')
hold off

xlabel('t')
legend('N1(t)','N2(t)','K1','K2')
axis([0 T 0 max(K1,K2)+0.5]) %determina el tamaño de los ejes, quedando 
                                %el eje x entre 0 y T (es decir, el tiempo 
                                %total del sistema) y el eje y entre 0 
                                %(mínimo tamaño poblacional) y un valor un 
                                %poco por arriba de K1 y K2 (los máximos 
                                %tamaños poblacionales posibles). 

figure(2)
plot(N(:,1),N(:,2))
title ('Diagrama de fases')
xlabel('N1')
ylabel('N2')

hold on
line([0 K1],[K1/a 0],'color','red') %esta linea y la siguiente marcan las nulclinas
line([0 K2/b],[K2 0],'color','black')
scatter(N(1,1),N(1,2),'filled') %marca el punto inicial, para que quede claro en qué sentido progresa el sistema
hold off

legend('N','dN1/dt=0','dN2/dt=0')
axis([0 max(K1,K2/b) 0 max(K2,K1/a)])