Logística con Umbral (15-A)
| Trabajo realizado por estudiantes | |
|---|---|
| Título | Logística con Umbral (Grupo 15-A) |
| Asignatura | Ecuaciones Diferenciales |
| Curso | Curso 2014-15 |
| Autores | • Javier Colorado Martínez (1278) • Álvaro Llera Fernández (501) • Francisco Javier Alcaraz de Amuriza (567) • Daniel Ballesteros Gálvez (1376) • Alberto Rodríguez Fernández (357) • Francisco Javier Barral González (632) |
| Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura | |
Contenido
1 Introducción
Nuestro Problema de Valor Inicial (P.V.I.) viene definido por la siguiente ecuación: [math]\left\{\begin{matrix}y’=-ry(1-{y\over M_1})(1-{y\over M_2}), t\gt0\\\ltbr/\gt\\y=y_0\end{matrix}\right.[/math] Donde \(M_1=30\), \(M_2=100\), \(r=0.04\) e \(y_0=60\). La ecuación a estudiar está definida en el intervalo \(I=[0,100]\).
2 Resolución numérica
A continuación, se va a proceder a modelizar dicha ecuación mediante los métodos de Euler, Heun y Runge-Kutta de orden 4. El procedimiento por el que se ha optado ha consistido en introducir en MATLAB cada uno de los métodos con los diferentes pasos, a saber, \(h=1\), \(h=0.1\) y \(h=0.01\).
2.1 Euler
%MÉTODO DE EULER
hold on
%Datos conocidos e intervalo de estudio.
t0=0; tN=100; r=0.04; M1=30; M2=100; y0=60;
%-------------------------------------------------------
%CASO 1: paso de 1
h1=1;
N1=(tN-t0)/h1;
t1=t0:h1:tN;
y1=zeros(1,N1+1);
y1(1)=y0;
for n=1:N1
y1(n+1)=y1(n)+h1*(-r*y1(n)*(1-y1(n)/M1)*(1-y1(n)/M2));
end
%-------------------------------------------------------
%CASO 2: paso de 0.1
h2=0.1;
N2=(tN-t0)/h2;
t2=t0:h2:tN;
y2=zeros(1,N2+1);
y2(1)=y0;
for n=1:N2
y2(n+1)=y2(n)+h2*(-r*y2(n)*(1-y2(n)/M1)*(1-y2(n)/M2));
end
%------------------------------------------------------
%CASO 3: PASO DE 0.01
h3=0.01;
N3=(tN-t0)/h3;
t3=t0:h3:tN;
y3=zeros(1,N3+1);
y3(1)=y0;
for n=1:N3
y3(n+1)=y3(n)+h3*-r*y3(n)*(1-y3(n)/M1)*(1-y3(n)/M2);
end
%-----------------------------------------------------
%Representación en gráficas.
plot(t1,y1,'r') %Caso 1.
plot(t2,y2,'g') %Caso 2.
plot(t3,y3,'b') %Caso 3.
title('Método de Euler con todos los casos')
xlabel('Tiempo')
ylabel('Evolución de la y')
legend('Caso 1 Euler','Caso 2 Euler','Caso 3 Euler','Location','best')
hold off
Como podemos apreciar en las gráficas
2.2 Heun
%MÉTODO DE HEUN
hold on
%Datos conocidos e intervalo de estudio.
t0=0; tN=100; r=0.04; M1=30; M2=100; y0=60;
%--------------------------------------------------------------
%CASO 1: paso de 1
h1=1;
N1=(tN-t0)/h1;
t1=t0:h1:tN;
y1=zeros(1,N1+1);
y1(1)=y0;
for n=1:N1
k1=-r*y1(n)*(1-y1(n)/M1)*(1-y1(n)/M2);
k2=-r*(y1(n)+k1*h1)*(1-(y1(n)+k1*h1)/M1)*(1-(y1(n)+k1*h1)/M2);
y1(n+1)=y1(n)+h1/2*(k1+k2);
end
%---------------------------------------------------------------
%CASO 2: PASO DE 0.1
h2=0.1;
N2=(tN-t0)/h2;
t2=t0:h2:tN;
y2=zeros(1,N2+1);
y2(1)=y0;
for n=1:N2
k1=-r*y2(n)*(1-y2(n)/M1)*(1-y2(n)/M2);
k2=-r*(y2(n)+k1*h2)*(1-(y2(n)+k1*h2)/M1)*(1-(y2(n)+k1*h2)/M2);
y2(n+1)=y2(n)+h2/2*(k1+k2);
end
%---------------------------------------------------------------
%CASO 3: PASO DE 0.01
h3=0.01;
N3=(tN-t0)/h3;
t3=t0:h3:tN;
y3=zeros(1,N3+1);
y3(1)=y0;
for n=1:N3
k1=-r*y3(n)*(1-y3(n)/M1)*(1-y3(n)/M2);
k2=-r*(y3(n)+k1*h3)*(1-(y3(n)+k1*h3)/M1)*(1-(y3(n)+k1*h3)/M2);
y3(n+1)=y3(n)+h3/2*(k1+k2);
end
%---------------------------------------------------------------
%Representación en gráficas.
plot(t1,y1,'r') %Caso 1.
plot(t2,y2,'g') %Caso 2.
plot(t3,y3,'b') %Caso 3.
title('Método de Heun con todos los casos')
xlabel('Tiempo')
ylabel('Evolución de la y')
legend('Caso 1 Heun','Caso 2 Heun','Caso 3 Heun','Location','best')
hold off
Como podemos apreciar
2.3 Runge-Kutta de orden 4
%MÉTODO DE RUNGE-KUTTA DE ORDEN 4
hold on
%Datos conocidos e intervalo de estudio.
t0=0; tN=100; r=0.04; M1=30; M2=100; y0=60;
%--------------------------------------------------------------------------
%CASO 1: paso de 1
h1=1;
N1=(tN-t0)/h1;
t1=t0:h1:tN;
y1=zeros(1,N1+1);
y1(1)=y0;
for n=1:N1
k1=-r*y1(n)*(1-y1(n)/M1)*(1-y1(n)/M2);
k2=-r*(y1(n)+1/2*k1*h1)*(1-(y1(n)+1/2*k1*h1)/M1)*(1-(y1(n)+1/2*k1*h1)/M2);
k3=-r*(y1(n)+1/2*k2*h1)*(1-(y1(n)+1/2*k2*h1)/M1)*(1-(y1(n)+1/2*k2*h1)/M2);
k4=-r*(y1(n)+k3*h1)*(1-(y1(n)+k3*h1)/M1)*(1-(y1(n)+k3*h1)/M2);
y1(n+1)=y1(n)+h1/6*(k1+2*k2+2*k3+k4);
end
%--------------------------------------------------------------------------
%CASO 2: paso de 0.1
h2=0.1;
N2=(tN-t0)/h2;
t2=t0:h2:tN;
y2=zeros(1,N2+1);
y2(1)=y0;
for n=1:N2
k1=-r*y2(n)*(1-y2(n)/M1)*(1-y2(n)/M2);
k2=-r*(y2(n)+1/2*k1*h2)*(1-(y2(n)+1/2*k1*h2)/M1)*(1-(y2(n)+1/2*k1*h2)/M2);
k3=-r*(y2(n)+1/2*k2*h2)*(1-(y2(n)+1/2*k2*h2)/M1)*(1-(y2(n)+1/2*k2*h2)/M2);
k4=-r*(y2(n)+k3*h2)*(1-(y2(n)+k3*h2)/M1)*(1-(y2(n)+k3*h2)/M2);
y2(n+1)=y2(n)+h2/6*(k1+2*k2+2*k3+k4);
end
%--------------------------------------------------------------------------
%CASO 3: paso de 0.01
h3=0.01;
N3=(tN-t0)/h3;
t3=t0:h3:tN;
y3=zeros(1,N3+1);
y3(1)=y0;
for n=1:N3
k1=-r*y3(n)*(1-y3(n)/M1)*(1-y3(n)/M2);
k2=-r*(y3(n)+1/2*k1*h3)*(1-(y3(n)+1/2*k1*h3)/M1)*(1-(y3(n)+1/2*k1*h3)/M2);
k3=-r*(y3(n)+1/2*k2*h3)*(1-(y3(n)+1/2*k2*h3)/M1)*(1-(y3(n)+1/2*k2*h3)/M2);
k4=-r*(y3(n)+k3*h3)*(1-(y3(n)+k3*h3)/M1)*(1-(y3(n)+k3*h3)/M2);
y3(n+1)=y3(n)+h3/6*(k1+2*k2+2*k3+k4);
end
%--------------------------------------------------------------------------
%Representación en gráficas.
plot(t1,y1,'r') %Caso 1.
plot(t2,y2,'g') %Caso 2.
plot(t3,y3,'b') %Caso 3.
title('Método de Runge-Kutta con todos los casos')
xlabel('Tiempo')
ylabel('Evolución de la y')
legend('Caso 1 RK4','Caso 2 RK4','Caso 3 RK4','Location','best')
hold off
Como podemos en la gráfica
2.4 Interpretaciones
De todas las gráficas representadas anteriormente, se puede apreciar que aquellas que corresponden a un paso de \(h=1\) (el caso 1, considerado como un paso amplio), se alejan de la solución exacta, a diferencia de lo que sucede con pasos menores (casos 2 (\(h=0.1\)) y 3 (\(h=0.01\))). A su vez, la solución más cercana a la exacta es la obtenida mediante el método de Runge-Kutta de orden 4.
2.5 Dinámica de poblaciones
Desde un punto de vista de la dinámica de poblaciones, dado que se trata de un método logístico, se puede apreciar que la población no crece indefinidamente, sino que llega a un punto en el que el ecosistema no admite más individuos de la especie en cuestión.
3 Método de Heun para diferentes poblaciones iniciales
%MÉTODO DE HEUN CON DISTINTOS VALORES INICIALES
clf; clear;
% Datos Iniciales.
t0=0; tN=100; r=0.04; M1=30; M2=100;
y0=120; z0=20;
%Paso e intervalo de estudio.
h=0.1;
N=(tN-t0)/h;
t=t0:h:tN;
%--------------------------------------------------------------------------
%Heun con valor inicial 120.
y1=zeros(1,N+1);
y1(1)=y0;
for n=1:N
k1=-r*y1(n)*(1-y1(n)/M1)*(1-y1(n)/M2);
k2=-r*(y1(n)+k1*h)*(1-(y1(n)+k1*h)/M1)*(1-(y1(n)+k1*h)/M2);
y1(n+1)=y1(n)+h/2*(k1+k2);
end
%--------------------------------------------------------------------------
%Heun con valor incial 20.
z1=zeros(1,N+1);
z1(1)=z0;
for n=1:N
K1=-r*z1(n)*(1-z1(n)/M1)*(1-z1(n)/M2);
K2=-r*(z1(n)+K1*h)*(1-(z1(n)+K1*h)/M1)*(1-(z1(n)+k1*h)/M2);
z1(n+1)=z1(n)+h/2*(K1+K2);
end
%--------------------------------------------------------------------------
%Comparación de gráficas.
hold on
plot(t,y1,'r')
plot(t,z1,'g')
title('Comparativa del método de Heun')
xlabel('Tiempo')
ylabel('Evolución de la población')
legend('Heun 120','Heun 20','Location','best')
hold off3.1 Interpretación
4 Modelos de competencia
[math]\left\{\begin{matrix}x’=a_1x-b_1x^2-c_1xy, \\\ltbr/\gt\\ y’=a_2x-b_2x^2-c_2xy \end{matrix}\right.[/math]
5 Apartado 5 (comparación de especies)
Caso 1: En cuanto el animal 2 aumenta su población, el animal 1 comienza a aumentar su población a su costa, se trata de un parasitismo. Cuando el animal 2 está casi extinto, el animal 1 empieza también a desaparecer.
Caso 2: Aunque pudiera parecer un amensalismo ya que cuando una especie aumenta mientras que la otra desciende, se trata de un neutralismo ya que ninguna de las dos especies tiende a desaparecer.
Caso 3: Se trata de una simbiosis ya que en la parte final de la gráfica ambas especies tienden a crecer.
Caso 4: A primera vista se podría sospechar que se trata de un parasitismo. No obstante, si fuera así, al verse reducido el número de individuos de la especie 2, el animal 1 (supuesto parásito) debería tender a disminuir también su población al quedarse sin huésped. De ahí deducimos que no es un caso de parasitismo, porque al verse disminuida la especie 2, la especie 1 no tiende a disminuir su población, sino que ambas se estabilizan. Podría tratarse de un caso de competición, porque el crecimiento de ambas especies llega a un máximo, a partir del cual el animal 1 sigue en tendencia positiva de crecimiento, mientras que el animal 2 decrece, aunque no se ha llegado a un consenso claro en esta interpretación.
Caso 5: En este caso se puede apreciar un amensalismo, ya que la especie 2 se ve reducida a medida que aumenta el número de individuos de la especie 1. Como observación, en la gráfica se ve que cuantos más animales 1 hay, más rápido decrece la población del animal 2, y además, cuantos menos individuos hay de la especie 2, más despacio decrece ésta. Caso 6: Vemos que una de las especies ve incrementada su población en proporción al número de individuos de la segunda especie, mientras ésta no se ve afectada por la primera.