Diferencia entre revisiones de «Modelos de mezclas(Grupo 15)»

De MateWiki
Saltar a: navegación, buscar
 
(No se muestran 109 ediciones intermedias de 3 usuarios)
Línea 1: Línea 1:
−
 
+
{{Trabajo|Modelo de mezclas|[[:Categoría:Ecuaciones Diferenciales|Ecuaciones Diferenciales]]|[[:Categoría:Trabajos 2012-13|2012-13]]}}
−
Dos pantanos A y B con 100 Hm3 de agua cada uno están unidos por una presa que deja pasar agua de A a B.
+
Dos pantanos A y B con 100 <math>Hm^3</math> de agua cada uno están unidos por una presa que deja pasar agua de A a B.
−
El pantano A recive 3 Hm3/dia de agua limpia provniente de rios y el B 1.5 Hm3/dia.
+
El pantano A recibe 3 <math>Hm^3/dia</math> de agua limpia provniente de rios y el B 1.5 <math>Hm^3/dia</math>.
−
Para manatener estable la presa, de A a B se deja pasar una media de 3 Hm3/dia mientras que la presa al final de B desaloja 4.5 Hm3/dia. Se produce un vertido en el pantano A que deja 20 toneladas de un cierto contaminante.
+
Para manatener estable la presa, de A a B se deja pasar una media de 3 <math>Hm^3/dia</math> mientras que la presa al final de B desaloja 4.5 <math>Hm^3/dia</math>. Se produce un vertido en el pantano A que deja 20 toneladas de un cierto contaminante.
 
Suponemos las siguientes hipotesis:
 
Suponemos las siguientes hipotesis:
 
*El contaminante está disuelto de forma homogénea en el agua de los pantanos.
 
*El contaminante está disuelto de forma homogénea en el agua de los pantanos.
Línea 25: Línea 25:
  
 
:<math>\begin{pmatrix} X_A(0)  \\ X_B(0) \end{pmatrix}=\begin{pmatrix} 20 \\ 0 \end{pmatrix}</math>
 
:<math>\begin{pmatrix} X_A(0)  \\ X_B(0) \end{pmatrix}=\begin{pmatrix} 20 \\ 0 \end{pmatrix}</math>
−
La solución de este sistema es inmediata, por variables separadas se obtiene <math>X_A</math> y sustituyendo en la otra ecuación se obtiene <math>X_B</math> :
+
 
 +
*La solución general de este sistema es inmediata, por variables separadas se obtiene <math>X_A</math> y sustituyendo en la otra ecuación se obtiene <math>X_B</math> :
 
:<math> X_A(t)=20e^{-3t/100}</math>
 
:<math> X_A(t)=20e^{-3t/100}</math>
 
:<math> X_B(t)=40(e^{^-3t/100}-e^{-4.5t/100})</math>
 
:<math> X_B(t)=40(e^{^-3t/100}-e^{-4.5t/100})</math>
  
−
*Ahora,suponemos que hay un tercer pantano C unido a B por una segunda presa que suelta 6 Hm3/dia a B, reciviendo 1.5 Hm3/dia de agua limpia de rios.
+
*Ahora,suponemos que hay un tercer pantano C unido a B por una segunda presa que suelta 6 <math>Hm^3/dia</math>, reciviendo 1.5 <math>Hm^3/dia</math> de agua limpia de rios.
 +
El sistema pasaría a depender también de la velocidad en C, ahora como no conocemos el volumen de C consideramos VolA=VolB=VolC= 100 <math> Hm^3 </math>=cte, tampoco sabemos la velocidad de C, pero si se tiene que mantener el volumen constante tendrá que recivir 4.5 <math> Hm^3 </math> de B (que es lo que suelta) y él mismo soltará 6 <math> Hm^3 </math>, quedando el sistema de la siguiente forma:
  
 +
:<math>
 +
\left\{\begin{matrix} X'_A=\frac{-3X_A}{100}\\ X'_B=\frac{3X_A}{100}-\frac{4.5X_B}{100}\\ X'_C=\frac{4.5X_B}{100}-\frac{6X_C}{100}\\X_A(t_{0})=20, \quad X_B(t_{0})=0\end{matrix}\right.
 +
</math>
  
 
=='''SEGUNDO SISTEMA DE ECUACIONES'''==
 
=='''SEGUNDO SISTEMA DE ECUACIONES'''==
  
−
Supongamos que se activa un plan de limpieza que consiste en bombear 1 Hm3/dia de agua del pantano B al A ajustando las cantidades de agua que dejan pasar las presas para mantener estables los niveles de los pantanos.
+
Supongamos que se activa un plan de limpieza que consiste en bombear 1 <math>Hm^3/dia</math> de agua del pantano B al A ajustando las cantidades de agua que dejan pasar las presas para mantener estables los niveles de los pantanos.
 
El sistema de ecuaciones va a cambiar, ahora el pantano A recivirá agua contaminada de B y para que se mantenga estable, dejara pasar más agua a B reciviendo este más agua contaminada.
 
El sistema de ecuaciones va a cambiar, ahora el pantano A recivirá agua contaminada de B y para que se mantenga estable, dejara pasar más agua a B reciviendo este más agua contaminada.
 
Los cambios se pueden ver en las siguientes ecuaciones:
 
Los cambios se pueden ver en las siguientes ecuaciones:
 +
 +
:<math> X'_A(t)=-{4\over 100}X_A+{1\over 100}X_B </math>
 +
 +
:<math> X'_B(t)={4\over 100}X_A-{5.5\over 100}X_B </math>
 +
 +
En forma matricial tenemos el siguiente problema de Cauchy:
 +
:<math> \begin{pmatrix} X'_A \\ X'_B \end{pmatrix}={1\over 100}\begin{pmatrix} -4 & 1 \\ 4 & -5.5 \end{pmatrix}\begin{pmatrix} X_A \\ X_B \end{pmatrix} </math>
 +
 +
:<math> \begin{pmatrix} X_A(0) \\ X_B(0) \end{pmatrix}=\begin{pmatrix} 20 \\ 0 \end{pmatrix} </math>
 +
Para solucionarlo es útil pasarlo a una ecuación diferencial de segundo orden, siendo esta de coeficientes constantes y homogenea, resolviendola sacando soluciones con el polinomio caracteristico, considerando que ambas forman un espacio vectorial de dimension 2 y sustituyendo las condiciones iniciales obtenemos la solución general, una vez obtenida sustituimos en el sistema inicial y sacamos la concentración del otro pantano.
 +
 +
 +
:<math>
 +
\begin{matrix}X''_A+{19\over 200}X'_A+{9\over 5000}X_A=0 \\ X_A(0)=20 \\ X'_A(0)=-{4\over 100}X_A(0)+{1\over 100}X_B(0)=-{4\over 5}\end{matrix}</math>
 +
*Las soluciones son:
 +
:<math> X_A(t)=(10/\sqrt{73})e^{-19t/400}\left [ (\sqrt{73}+35)e^{\sqrt{73}t/400}+(\sqrt{73}-35)e^{-\sqrt{73}t/400}\right ]    </math>
 +
 +
:<math> X_B(t)=(320/\sqrt{73})e^{-19t/400}sinh(\sqrt{73}/400)</math>
 +
 +
=='''RESOLUCIÓN MEDIANTE EL MÉTODO DE EULER'''==
 +
 +
El método de Euler es un método númerico explicito de primer orden.Resolvemos los dos sistemas con MATLAB mediante el método de Euler y damos una aproximación a las soluciones generales obtenidas en los aparatados anteriores.
 +
 +
'''PRIMER SISTEMA'''
 +
{{matlab|codigo=
 +
 +
clear all
 +
format long
 +
% PRIMER SISTEMA
 +
t0=0; tN=6*100/4.5; % Extremos del intervalo temporal a estudiar
 +
N=100; h=(tN-t0)/N; % Número y amplitud de subintervalos
 +
t=t0:h:tN; % Valores de t
 +
x0=20;y0=0; % Asignacion de valores iniciales x->XA y->XB
 +
x(1)=x0;y(1)=y0; % Inicialización de los vectores
 +
%Algoritmo
 +
for n=1:N
 +
    x(n+1)=x(n)+h*(-3/100)*x(n); %Aplicación del método de Euler
 +
    y(n+1)=y(n)+h*(3/100)*x(n)-h*(4.5/100)*y(n);
 +
end
 +
% Dibujo de gráficas
 +
plot(t,x,'b',t,y,'r')
 +
}}
 +
 +
'''SEGUNDO SISTEMA'''
 +
{{matlab|codigo=
 +
 +
clear all
 +
format long
 +
% SEGUNDO SISTEMA
 +
t0=0; tN=6*100/4.5; % Extremos del intervalo temporal a estudiar
 +
N=100; h=(tN-t0)/N; % Número y amplitud de subintervalos
 +
t=t0:h:tN; % Valores de t
 +
x0=20;y0=0; % Asignacion de valores iniciales x->XA y->XB
 +
x(1)=x0;y(1)=y0; % Inicialización de los vectores
 +
%Algoritmo
 +
for n=1:N
 +
    x(n+1)=x(n)+h*(-4/100)*x(n)+h*(1/100)*y(n); %Aplicación del método de Euler
 +
    y(n+1)=y(n)+h*(4/100)*x(n)-h*(5.5/100)*y(n);
 +
end
 +
% Dibujo de gráficas
 +
plot(t,x,'b',t,y,'r')
 +
}}
 +
 +
 +
Gráficamente, las soluciones son:
 +
 +
[[Archivo:DIBUJO SISTEMA1 EULER.jpg|marco|centro]][[Archivo:DIBUJO SISTEMA2 EULER.jpg|marco|centro]]
 +
 +
Se observa que ambas soluciones son del tipo exponencial y se ve que para un tiempo grande la contaminación acaba desapareiendo. Tambien se observa que ambos pantanos llegan a tener la misma contaminación y que para el sistema 2, en el que se activa un plan de limpieza, la grafica del pantano B se desplaza hacia la izquierda lo que indica que tardaría menos en desaparecer toda la contaminación.
 +
 +
*La diferencia de tiempo que tarda en desaparecer la mitad del contaminante inicial (20 toneladas) en el pantano A cuando se activa el sistema de limpieza será la diferencia de valores medidos en ambas graficas, es decir, t2-t1=23-18=5 dias y para la tercera parte es t4-t3=13-10=3 dias.
 +
 +
=='''RESOLUCIÓN MEDIANTE EL MÉTODO DE RUNGE-KUTTA'''==
 +
 +
El método de Runge-Kutta es un método implicito de cuarto orden con lo que nos dára una mejor aproximación que el método de Euler al tener mayor orden, pero al ser implicito la programación será más difícl. Resolvemos los dos sistemas en MATLAB mediante el método de Runge-Kutta:
 +
 +
'''PRIMER SISTEMA'''
 +
 +
{{matlab|codigo=
 +
clear all
 +
format long
 +
t0=0; tN=6*100/4.5;
 +
N=100; h=(tN-t0)/N;
 +
t=t0:h:tN;
 +
x0=20;y0=0;
 +
x(1)=x0;y(1)=y0;
 +
for n=1:N
 +
% Aplicación del método de Runge-Kutta
 +
    p1=(-3/100)*x(n);
 +
    q1=(3/100)*x(n)-(4.5/100)*y(n);
 +
    p2=(-3/100)*(x(n)+h*p1/2);
 +
    q2=(3/100)*(x(n)+h*p1/2)-(4.5/100)*(y(n)+h*q1/2);
 +
    p3=(-3/100)*(x(n)+h*p2/2);
 +
    q3=(3/100)*(x(n)+h*p2/2)-(4.5/100)*(y(n)+h*q2/2);
 +
    p4=(-3/100)*(x(n)+h*p3);
 +
    q4=(3/100)*(x(n)+h*p3)-(4.5/100)*(y(n)+h*q3);
 +
    x(n+1)=x(n)+(h/6)*(p1+2*p2+2*p3+p4);
 +
    y(n+1)=y(n)+(h/6)*(q1+2*q2+2*q3+q4);
 +
end
 +
plot(t,x,'b',t,y,'r')
 +
}}
 +
 +
 +
'''SEGUNDO SISTEMA'''
 +
 +
{{matlab|codigo=
 +
clear all
 +
format long
 +
t0=0; tN=6*100/4.5;
 +
N=100; h=(tN-t0)/N;
 +
t=t0:h:tN;
 +
x0=20;y0=0;
 +
x(1)=x0;y(1)=y0;
 +
for n=1:N
 +
% Aplicación del método de Runge-Kutta
 +
    p1=(-4/100)*x(n)+(1/100)*y(n);
 +
    q1=(4/100)*x(n)-(5.5/100)*y(n);
 +
    p2=(-4/100)*(x(n)+h*p1/2)+(1/100)*(y(n)+h*q1/2);
 +
    q2=(4/100)*(x(n)+h*p1/2)-(5.5/100)*(y(n)+h*q1/2);
 +
    p3=(-4/100)*(x(n)+h*p2/2)+(1/100)*(y(n)+h*q2/2);
 +
    q3=(4/100)*(x(n)+h*p2/2)-(5.5/100)*(y(n)+h*q2/2);
 +
    p4=(-4/100)*(x(n)+h*p3)+(1/100)*(y(n)+h*q3);
 +
    q4=(4/100)*(x(n)+h*p3)-(5.5/100)*(y(n)+h*q3);
 +
    x(n+1)=x(n)+(h/6)*(p1+2*p2+2*p3+p4);
 +
    y(n+1)=y(n)+(h/6)*(q1+2*q2+2*q3+q4);
 +
end
 +
plot(t,x,'b',t,y,'r')
 +
}}
 +
Gráficamente, las soluciones son:
 +
 +
[[Archivo:S1RK.jpg|marco|centro]][[Archivo:S2RK.jpg|marco|centro]]
 +
*Para comparar el método de Runge-Kutta con el de Euler, se pueden comparar los gráficos de ambos métodos o bien como en nuestro caso, hemos cojido el 70 valor de la solucion x (por cojer uno) del Runge-Kutta y resulta que nos da el mismo valor que haciendolo con la solucion exacta(x=1,2658) y haciendo lo mismo con el de Euler nos da un valor diferente(x=1,1960) con lo que se verifica que el método de Runge-Kutta da una mejor aproximación de la solución.
 +
 +
=='''CONTAMINANTE INICIAL DESCONOCIDO'''==
 +
 +
 +
'''Si no conocemos la cantidad de contaminante inicial pero sabemos tras unos dias se redujo el contaminante a sólo una tonelada en A y dos en B,¿Cuánto contaminante se estima que se vertió inicialmente?'''
 +
En este caso desconocemos el contaminante inicial,pero para el tiempo tN=0 y otro tiempo inicial por ejmplo t0=100, conocemos los valores del contaminante <math>(X_A,X_B)=(1,2)</math>, con esto aplicando el método de Euler podemos estimar el contaminante inicial.
 +
 +
'''CÓDIGO EN MATLAB'''
 +
 +
{{matlab|codigo=
 +
t0=100; tN=0;
 +
x0=[1 2]'; %Valores iniciales para t0
 +
N=100; h=(tN-t0)/N; %Número y amplitud de subintervalos
 +
A=[-3/100,0;3/100,-4.5/100]; %Matriz de coeficientes del primer sistema
 +
x=x0;
 +
xA(1)=x(1); %Iniciamos el vector cantidad de contaminante en A
 +
xB(1)=x(2); %Iniciamos el vector cantidad de contaminante en B
 +
for n=1:N
 +
    x=x+h*(A*x); %Aplicamos el método de Euler
 +
    xA(n+1)=x(1);
 +
    xB(n+1)=x(2);
 +
end
 +
x=t0:h:tN;
 +
%Dibujo
 +
plot(x,xA,'b',x,xB,'r')
 +
}}
 +
 +
Gráficamente la solución es:
 +
[[Archivo:CID.jpg|marco|centro]]
 +
 +
*De la gráfica se puede estimar la cantidad de contaminate inicial y sera XA(0) y XB(0).
 +
 +
 +
 +
[[Categoría:Grado en Ingeniería Civil y Territorial]]
 +
[[Categoría:Ecuaciones Diferenciales]]
 +
[[Categoría:Trabajos 2012-13]]

Revisión actual del 00:28 3 jun 2013

Dos pantanos A y B con 100 [math]Hm^3[/math] de agua cada uno están unidos por una presa que deja pasar agua de A a B. El pantano A recibe 3 [math]Hm^3/dia[/math] de agua limpia provniente de rios y el B 1.5 [math]Hm^3/dia[/math]. Para manatener estable la presa, de A a B se deja pasar una media de 3 [math]Hm^3/dia[/math] mientras que la presa al final de B desaloja 4.5 [math]Hm^3/dia[/math]. Se produce un vertido en el pantano A que deja 20 toneladas de un cierto contaminante. Suponemos las siguientes hipotesis:

  • El contaminante está disuelto de forma homogénea en el agua de los pantanos.
  • Al entrar o salir agua en un pantano, está se mezcla con el agua del pantano de forma inmediata creando un mezcla homogénea.
  • La variación de contaminante en un lago es la diferencia ente el contaminate que entra y el que sale del lago.

1 PRIMER SISTEMA DE ECUACIONES

Llamamos [math]X_A(t)[/math] a la cantidad de contaminante del lago A y [math]X_B(t)[/math] a la cantidad de contaminante del lago B. Definimos la variación de contaminante en un lago como:

[math]X'_A(t)=velocidad(entrada)-velocidad(salida)[/math]

Conocemos la velocidad de entrada y salida de volumen de agua de cada pantano, está multiplicada por la concentración de contaminante en el instante t nos da la velocidad de entrada y salida de contaminante. Definimos como concentración de contaminante, siendo vol(t) el volumen del lago en un instante t:

[math]C(t)={X(t)\over vol(t)}[/math]

Aplicando las hipotesis y sabiendo que el contaminante se traspasa de lago a lago y que el agua de los rios es agua limpia, nos quedan las siguiente ecuaciones: \begin{array}{c} X'_A(t)=-{3\over 100}X_A \\ X'_B(t)={3\over 100}X_A-{4.5\over 100}X_B \end{array}

En forma matricial y considerando que [math]X_A(0)=20[/math] nos queda el siguiente problema de Cauchy:

[math]\begin{pmatrix} X'_A\\ X'_B\end{pmatrix}={1\over 100}\begin{pmatrix} -3 & 0 \\ 3 & -4.5 \end{pmatrix}\begin{pmatrix} X_A \\ X_B \end{pmatrix}[/math]
[math]\begin{pmatrix} X_A(0) \\ X_B(0) \end{pmatrix}=\begin{pmatrix} 20 \\ 0 \end{pmatrix}[/math]
  • La solución general de este sistema es inmediata, por variables separadas se obtiene [math]X_A[/math] y sustituyendo en la otra ecuación se obtiene [math]X_B[/math] :
[math] X_A(t)=20e^{-3t/100}[/math]
[math] X_B(t)=40(e^{^-3t/100}-e^{-4.5t/100})[/math]
  • Ahora,suponemos que hay un tercer pantano C unido a B por una segunda presa que suelta 6 [math]Hm^3/dia[/math], reciviendo 1.5 [math]Hm^3/dia[/math] de agua limpia de rios.

El sistema pasaría a depender también de la velocidad en C, ahora como no conocemos el volumen de C consideramos VolA=VolB=VolC= 100 [math] Hm^3 [/math]=cte, tampoco sabemos la velocidad de C, pero si se tiene que mantener el volumen constante tendrá que recivir 4.5 [math] Hm^3 [/math] de B (que es lo que suelta) y él mismo soltará 6 [math] Hm^3 [/math], quedando el sistema de la siguiente forma:

[math] \left\{\begin{matrix} X'_A=\frac{-3X_A}{100}\\ X'_B=\frac{3X_A}{100}-\frac{4.5X_B}{100}\\ X'_C=\frac{4.5X_B}{100}-\frac{6X_C}{100}\\X_A(t_{0})=20, \quad X_B(t_{0})=0\end{matrix}\right. [/math]

2 SEGUNDO SISTEMA DE ECUACIONES

Supongamos que se activa un plan de limpieza que consiste en bombear 1 [math]Hm^3/dia[/math] de agua del pantano B al A ajustando las cantidades de agua que dejan pasar las presas para mantener estables los niveles de los pantanos. El sistema de ecuaciones va a cambiar, ahora el pantano A recivirá agua contaminada de B y para que se mantenga estable, dejara pasar más agua a B reciviendo este más agua contaminada. Los cambios se pueden ver en las siguientes ecuaciones:

[math] X'_A(t)=-{4\over 100}X_A+{1\over 100}X_B [/math]
[math] X'_B(t)={4\over 100}X_A-{5.5\over 100}X_B [/math]

En forma matricial tenemos el siguiente problema de Cauchy:

[math] \begin{pmatrix} X'_A \\ X'_B \end{pmatrix}={1\over 100}\begin{pmatrix} -4 & 1 \\ 4 & -5.5 \end{pmatrix}\begin{pmatrix} X_A \\ X_B \end{pmatrix} [/math]
[math] \begin{pmatrix} X_A(0) \\ X_B(0) \end{pmatrix}=\begin{pmatrix} 20 \\ 0 \end{pmatrix} [/math]

Para solucionarlo es útil pasarlo a una ecuación diferencial de segundo orden, siendo esta de coeficientes constantes y homogenea, resolviendola sacando soluciones con el polinomio caracteristico, considerando que ambas forman un espacio vectorial de dimension 2 y sustituyendo las condiciones iniciales obtenemos la solución general, una vez obtenida sustituimos en el sistema inicial y sacamos la concentración del otro pantano.


[math] \begin{matrix}X''_A+{19\over 200}X'_A+{9\over 5000}X_A=0 \\ X_A(0)=20 \\ X'_A(0)=-{4\over 100}X_A(0)+{1\over 100}X_B(0)=-{4\over 5}\end{matrix}[/math]
  • Las soluciones son:
[math] X_A(t)=(10/\sqrt{73})e^{-19t/400}\left [ (\sqrt{73}+35)e^{\sqrt{73}t/400}+(\sqrt{73}-35)e^{-\sqrt{73}t/400}\right ] [/math]
[math] X_B(t)=(320/\sqrt{73})e^{-19t/400}sinh(\sqrt{73}/400)[/math]

3 RESOLUCIÓN MEDIANTE EL MÉTODO DE EULER

El método de Euler es un método númerico explicito de primer orden.Resolvemos los dos sistemas con MATLAB mediante el método de Euler y damos una aproximación a las soluciones generales obtenidas en los aparatados anteriores.

PRIMER SISTEMA 

clear all
format long
% PRIMER SISTEMA
t0=0; tN=6*100/4.5; % Extremos del intervalo temporal a estudiar
N=100; h=(tN-t0)/N; % Número y amplitud de subintervalos
t=t0:h:tN; % Valores de t
x0=20;y0=0; % Asignacion de valores iniciales x->XA y->XB
x(1)=x0;y(1)=y0; % Inicialización de los vectores
%Algoritmo
for n=1:N
    x(n+1)=x(n)+h*(-3/100)*x(n); %Aplicación del método de Euler 
    y(n+1)=y(n)+h*(3/100)*x(n)-h*(4.5/100)*y(n); 
end
% Dibujo de gráficas
plot(t,x,'b',t,y,'r')


SEGUNDO SISTEMA 

clear all
format long
% SEGUNDO SISTEMA
t0=0; tN=6*100/4.5; % Extremos del intervalo temporal a estudiar
N=100; h=(tN-t0)/N; % Número y amplitud de subintervalos
t=t0:h:tN; % Valores de t
x0=20;y0=0; % Asignacion de valores iniciales x->XA y->XB
x(1)=x0;y(1)=y0; % Inicialización de los vectores
%Algoritmo
for n=1:N
    x(n+1)=x(n)+h*(-4/100)*x(n)+h*(1/100)*y(n); %Aplicación del método de Euler 
    y(n+1)=y(n)+h*(4/100)*x(n)-h*(5.5/100)*y(n); 
end
% Dibujo de gráficas
plot(t,x,'b',t,y,'r')


Gráficamente, las soluciones son:

centro
centro

Se observa que ambas soluciones son del tipo exponencial y se ve que para un tiempo grande la contaminación acaba desapareiendo. Tambien se observa que ambos pantanos llegan a tener la misma contaminación y que para el sistema 2, en el que se activa un plan de limpieza, la grafica del pantano B se desplaza hacia la izquierda lo que indica que tardaría menos en desaparecer toda la contaminación.

  • La diferencia de tiempo que tarda en desaparecer la mitad del contaminante inicial (20 toneladas) en el pantano A cuando se activa el sistema de limpieza será la diferencia de valores medidos en ambas graficas, es decir, t2-t1=23-18=5 dias y para la tercera parte es t4-t3=13-10=3 dias.

4 RESOLUCIÓN MEDIANTE EL MÉTODO DE RUNGE-KUTTA

El método de Runge-Kutta es un método implicito de cuarto orden con lo que nos dára una mejor aproximación que el método de Euler al tener mayor orden, pero al ser implicito la programación será más difícl. Resolvemos los dos sistemas en MATLAB mediante el método de Runge-Kutta:

PRIMER SISTEMA 

clear all
format long
t0=0; tN=6*100/4.5; 
N=100; h=(tN-t0)/N; 
t=t0:h:tN; 
x0=20;y0=0;
x(1)=x0;y(1)=y0; 
for n=1:N
% Aplicación del método de Runge-Kutta
    p1=(-3/100)*x(n);
    q1=(3/100)*x(n)-(4.5/100)*y(n);
    p2=(-3/100)*(x(n)+h*p1/2);
    q2=(3/100)*(x(n)+h*p1/2)-(4.5/100)*(y(n)+h*q1/2);
    p3=(-3/100)*(x(n)+h*p2/2);
    q3=(3/100)*(x(n)+h*p2/2)-(4.5/100)*(y(n)+h*q2/2);
    p4=(-3/100)*(x(n)+h*p3);
    q4=(3/100)*(x(n)+h*p3)-(4.5/100)*(y(n)+h*q3);
    x(n+1)=x(n)+(h/6)*(p1+2*p2+2*p3+p4);
    y(n+1)=y(n)+(h/6)*(q1+2*q2+2*q3+q4);
end
plot(t,x,'b',t,y,'r')


SEGUNDO SISTEMA 

clear all
format long
t0=0; tN=6*100/4.5; 
N=100; h=(tN-t0)/N; 
t=t0:h:tN; 
x0=20;y0=0;
x(1)=x0;y(1)=y0; 
for n=1:N
% Aplicación del método de Runge-Kutta
    p1=(-4/100)*x(n)+(1/100)*y(n);
    q1=(4/100)*x(n)-(5.5/100)*y(n);
    p2=(-4/100)*(x(n)+h*p1/2)+(1/100)*(y(n)+h*q1/2);
    q2=(4/100)*(x(n)+h*p1/2)-(5.5/100)*(y(n)+h*q1/2);
    p3=(-4/100)*(x(n)+h*p2/2)+(1/100)*(y(n)+h*q2/2);
    q3=(4/100)*(x(n)+h*p2/2)-(5.5/100)*(y(n)+h*q2/2);
    p4=(-4/100)*(x(n)+h*p3)+(1/100)*(y(n)+h*q3);
    q4=(4/100)*(x(n)+h*p3)-(5.5/100)*(y(n)+h*q3);
    x(n+1)=x(n)+(h/6)*(p1+2*p2+2*p3+p4);
    y(n+1)=y(n)+(h/6)*(q1+2*q2+2*q3+q4);
end
plot(t,x,'b',t,y,'r')

Gráficamente, las soluciones son:

centro
centro
  • Para comparar el método de Runge-Kutta con el de Euler, se pueden comparar los gráficos de ambos métodos o bien como en nuestro caso, hemos cojido el 70 valor de la solucion x (por cojer uno) del Runge-Kutta y resulta que nos da el mismo valor que haciendolo con la solucion exacta(x=1,2658) y haciendo lo mismo con el de Euler nos da un valor diferente(x=1,1960) con lo que se verifica que el método de Runge-Kutta da una mejor aproximación de la solución.

5 CONTAMINANTE INICIAL DESCONOCIDO

Si no conocemos la cantidad de contaminante inicial pero sabemos tras unos dias se redujo el contaminante a sólo una tonelada en A y dos en B,¿Cuánto contaminante se estima que se vertió inicialmente? En este caso desconocemos el contaminante inicial,pero para el tiempo tN=0 y otro tiempo inicial por ejmplo t0=100, conocemos los valores del contaminante [math](X_A,X_B)=(1,2)[/math], con esto aplicando el método de Euler podemos estimar el contaminante inicial.

CÓDIGO EN MATLAB 

t0=100; tN=0; 
x0=[1 2]'; %Valores iniciales para t0
N=100; h=(tN-t0)/N; %Número y amplitud de subintervalos
A=[-3/100,0;3/100,-4.5/100]; %Matriz de coeficientes del primer sistema
x=x0; 
xA(1)=x(1); %Iniciamos el vector cantidad de contaminante en A
xB(1)=x(2); %Iniciamos el vector cantidad de contaminante en B
for n=1:N
    x=x+h*(A*x); %Aplicamos el método de Euler
    xA(n+1)=x(1);
    xB(n+1)=x(2); 
end
x=t0:h:tN;
%Dibujo
plot(x,xA,'b',x,xB,'r')


Gráficamente la solución es:

centro
  • De la gráfica se puede estimar la cantidad de contaminate inicial y sera XA(0) y XB(0).
Obtenido de «https://mat.caminos.upm.es/w/index.php?title=Modelos_de_mezclas(Grupo_15)&oldid=1954»