Diferencia entre revisiones de «Modelo Lokta-Volterra Prey-Predator. Grupo 6»

Revisión del 21:10 3 mar 2013

1 1 INTERPRETACION DEL MODELO

El modelo Lokta-Volterra,asumiendo las hipótesis dadas en el ejercicio, expresa: Sin zorros, los conejos se reproducen siguiendo el modelo de Maltus, por lo que su tasa de crecimiento es proporcional a su tamaño.

Maltus: R'(t) = aR(t)

Sin conejos que comer, la población de zorros disminuye proporcionalmente a su tamaño.
F'(t)=-bF(t)

La tasa a la cual los conejos son comidos por zorros es proporcional a la tasa de interacción entre zorros y conejos.
R'(t)=-cF(t)R(t)

La tasa por la cual los zorros nacen es proporcional a la tasa por la cual los conejos son comidos.
F'(t)= dF(t)R(t)

Por lo tanto, este modelo nos da el ritmo de crecimiento de las poblaciones de los depredadores y las presas (conejos y zorros) dado el número de miembros de cada una. El sistema que resulta es:

1.1 Ecuaciones:

[math] \frac{\mathrm{d} R}{\mathrm{d} t}=aR-cRF [/math]

[math] \frac{\mathrm{d} F}{\mathrm{d} t}=-bF+dRF [/math]

En estas ecuaciones los parámetros tienen los siguientes significados:

• a: el factor proporcional relacionado con el ritmo de reproducción de los conejos cuando no hay zorros
• b: es el factor con el que los zorros decrecen a un ritmo proporcional a su tamaño sin conejos que comer.
• c: es la proporción a la cual las interacciones entre zorros y conejos hacen decrecer a la población de los conejos.
• d: es el factor proporcional que maneja el ritmo al cual los zorros nacen. Este ritmo es proporcional al ritmo al que los conejos son comidos por los zorros.

2 2 LOKTA-VOLTERRA SEGUN EULER

Para una población de conejos y zorros ( R(t) y F(t)), particularizamos los parámetros arriba indicados (a=0.4, b=0.37, c=0.3, d=0.05) y consideramos que en un tiempo tЄ[0,100] tenemos una población de 3000 conejos y de 1000 zorros.

2.1 Codigo Matlab

Resolvemos el sistema según el siguiente código MATLAB:

% Z es la variable y su primer elemento sera R y el segundo F%
a=0.4,b=0.37,c=0.3,d=0.05;
t=0;
i=1;
z=[3;1];
f1(i)=Z(1),f2(i)=Z(2);
x(i)=t;
h=0.01;
while t<100
  i=i+1;
  D=Z;
  t=t+h;
  F=[a*Z(1)-c*Z(1)*Z(2),d*Z(1)*D(2)-b*Z(2)];
  p=h*F
  X=t;
  Z=D+h*F;
  f1(i)=Z(1),f2(i)=Z(2);
  x(i)=t;
end
figure(1)
hold on
plot(x;f1,'b',x,f2,'r')
figure(2)
hold on
plot(fi,f2,'g')
hold off

2.2 Graficas

Distribución de conejos y zorros en el tiempo. Siendo la línea de arriba (azul) los conejos y la de abajo (rojo) los zorros.

De la gráfica se deduce que el aumento de la población de presas conlleva a un aumento de la población de depredadores, del mismo modo que un aumento de la población de depredadores reduce la de presas a mínimos.Con esta gráfica concluimos que las densidades poblacionales de ambas especies están estrechamente relacionadas.

Distribución de conejos-zorros en el tiempo en 2D. En el eje de abscisas se marcan los conejos y en el eje de ordenadas los zorros.

Esta gráfica representa la interacción de depredadores y presas en el tiempo, reafirmando la relación previamente explicada. Por lo tanto, concluimos que el ecosistema es estable, oscilando alrededor del punto de equilibrio:

[math] R=\frac{c}{d}=\frac{0.3}{0.05}=6 \qquad F=\frac{a}{b}=\frac{0.4}{0.37}=1.081 [/math]

sería R=c/d y F=a/b. En nuestro problema por tanto: R=0.3/0.05= 6 yF=0.4/0.37=1.081.

3 3 EXTINCIÓN PARA DISTINTAS POBLACIONES

3.1 Codigo Matlab

sprintf('conejos 1500')
a=0.4,b=0.37,c=0.3,d=0.05;
t=0;
i=1;
Z=[1.5;1];
p1(i)=Z(1),p2(i)=Z(2);
x(i)=t;
ext=0;
while t<100
  i=i+1;
  D=Z;
  t=t+h;
  F=[a*Z(1)-c*Z(1)*Z(2);d*Z(1)*D(2)-b*Z(2)];
  P=h*F;
  X=t;
  Z=D+h*F;
  p1(i)=Z(1),p2(i)=Z(2);
  x(i)=t;
  if Z(2)<(50/1000)
  ext=1;
  end
end
if ext>0
  sprintf('extincion zorros')
end
sprintf('maximo y minimo de los conejos')
max(p1)
min(p1)
sprintf('maximo y minimo de los zorros')
max(p2)
min(p2)
sprintf('conejos 1000')
a=0.4,b=0.37,c=0.3,d=0.05;
t=0;
i=1;
Z=[1;1];
p1(i)=Z(1),p2(i)=Z(2);
x(i)=t;
ext=0;
while t<100
  i=i+1;
  D=Z;
  t=t+h;
  F=[a*Z(1)-c*Z(1)*Z(2),d*Z(1)*D(2)-b*Z(2)];
  P=h*F;
  X=t;
  Z=D+h*F;
  p1(i)=Z(1),p2(i)=Z(2);
  x(i)=t;
  if Z(2)<=(50/1000)
  ext=1;
  end
end
if ext>0
  sprintf('extincion zorros')
end
sprintf('maximo y minimo de los conejos')
max(p1)
min(p1)
sprintf('maximo y minimo de los zorros')
max(p2)
min(p2)
sprintf('conejos 250')
a=0.4,b=0.37,c=0.3,d=0.05;
t=0;
i=1;
Z=[0.25;1];
p1(i)=Z(1),p2(i)=Z(2);
x(i)=t;
ext=0;
while t<100
  i=i+1;
  D=Z;
  t=t+h;
  F=[a*Z(1)-c*Z(1)*Z(2),d*Z(1)*D(2)-b*Z(2)];
  P=h*F;
  X=t;
  Z=D+h*F;
  p1(i)=Z(1),p2(i)=Z(2);
  x(i)=t;
  if Z(2)<=(50/1000)
  ext=1;
  end
end
if ext>0
  sprintf('extincion zorros')
end
sprintf('maximo y minimo de los conejos')
max(p1)
min(p1)
sprintf('maximo y minimo de los zorros')
max(p2)
min(p2)

Programando en MATLAB para estas tres poblaciones se obtienen los siguientes resultados:

De lo que se observa que las proporciones quedan:

Concluimos, por tanto, que si se disminuye la población de conejos la proporción conejos/zorros se mantiene constante, es decir, si disminuye la población de la presa también lo hará proporcionalmente la del depredador (como se observa en el pico de la gráfica remarcado más arriba).

Reduciendo la población de conejos a 250 obtendríamos un máximo de 43 zorros, inferior a los 50 que marcan el límite de extinción, lo que nos indica que la extinción sería algo inevitable.

4 4 LOKTA-VOLTERRA POR RUNGE-KUTTA

Resolver la ecuación por el método de Runge-Kutta de orden 4. Compara los resultados con los de Euler para los diferentes tiempos de las h fijadas.¿Cuál es la dificultad cuando se usa el método del trapecio?.

4.1 Código Matlab

a=0.4,b=0.37,c=0.3,d=0.05;
t=0;
i=1;
z=[3;1];
f3(i)=Z(1),f4(i)=Z(2);
x(i)=t;
while t<100 
  i=i+1;
  D=Z;
  t=t+h;
  K1=[a*D(1)-c*D(1)*D(2);d*D(1)*D(2)-b*D(2)];
  K2=[a*(D(1)+h/2*K1(1))-c*(D(1)+h/2*K1(1))*(D(2)+h/2*K1(2));d*(D(1)+h/2*K1(1))*(D(2)+h/2*K1(2))-b*(D(2)+h/2*K1(2))];
  K3=[a*(D(1)+h/2*K2(1))-c*(D(1)+h/2*K2(1))*(D(2)+h/2*K2(2));d*(D(1)+h/2*K2(1))*(D(2)+h/2*K2(2))-b*(D(2)+h/2*K2(2))];
  K4=[a*(D(1)+h*K3(1))-c*(D(1)+h*K3(1))*(D(2)+h*K3(2));d*(D(1)+h*K3(1))*(D(2)+h*K3(2))-b*(D(2)+h*K3(2))];
  z=D+h/6*(K1+2.*K2+2.*K3+K4);
  f3(i)=z(1);f4(i)=z(2);
  x(i)=t;
end
figure(3)
hold on
plot(f1,f2,'b',f3,f4,'r')

4.2 Gráfica