Línea de fase, ecuaciones diferenciales autónomas

Sigue leyendo

Anuncios

Métodos numéricos en ecuaciones diferenciales y MATLAB

Esta parte del blog está dedicada a los métodos numéricos aplicados a obtener aproximaciones de soluciones de ecuaciones diferenciales de primer orden. Para aclarar de modo más técnico este problema consideremos la ecuación diferencial con condición inicial

\frac{dy}{dx}=f(x,y), \hskip2cm y(a)=y_0

El problema que buscamos resolver consiste en determinar una aproximación para el valor y(b) para b\neq a.

 

Usted lector se preguntará seguramente, ¿porqué es necesario resolver este problema si ya conozco métodos que me permiten resolver muchas ecuaciones diferenciales?

Usted, amable lector de este blog, podrá imaginarse sin mucha dificultad ecuaciones diferenciales de primer orden, que al menos a primera vista  parezcan no tenener solución o al menos que hallar la solución no sea una tarea simple . Por ejemplo, pensemos en la ecuación diferencial.

 

\frac{dy}{dx}=(2x+y)x^2\sin^3(e^{x+y}).

Su solución pareciera no es fácil de determinar, en caso de que exista. Hablando de existencia de soluciones,   podemos responder afirmativamente a los problemas de existencia y unicidad de soluciones de la ecuaión diferencial para todo punto de en plano, usando el teorema de existencia y unicidad de ecuaciones diferenciales y que el lector encontrará en la sección correspondiente de este blog). Como es muy factible que no podamos resolver esta ecuación diferencia con condición inicial dada y(a)=y_0, es por ello que sería de gran utilidad hallar aproximaciones para valores y(b) para b\neq a.

 

En el documento siguiente, en formato pdf “SOLUCIÓN NUMÉRICA EDO’s”, podrá hallar información sobre los métodos más elementales que se usan para aproximar soluciones para una ecuación diferencial de primer orden y para sistemas de ecuaciones diferenciales de primer orden. Desafortunadamente no contiene explicaciones explícitas del desarrollo de los distintos métodos, pero espero le sea de ayuda si ha tenido un primer contacto con este problema. Conforme el tiempo me lo vaya permitiendo, iré desarrollando en esta sección las ideas contenidas en el referido documento.

Contiene información sobre

* Método de Euler

* Método de Euler mejorado, también conocido como Predictor-Corrector.

* Método del punto medio

* Método de Runge-Kutta de cuarto orden

* Ideas de las estructuras computaciones que sirven de modelos para desarrollar programas.

* Programas en Matlab

* Ejemplos de los distintos métodos desarrollados

* Obtención del método de Runge-Kutta de segundo orden y determinación de algunos casos particulares de la familia de métodos de Runge-Kutta de segundo orden obtenidos.

* Compraciones numéricas de los métodos usados.

* Desarrollos sobre como extender los métodos de Euler y punto medio  para sistemas de ecuaciones diferenciales de primer orden.

 

Link para Documento

SOLUCION NUMERICA EDO’s

El siguiente link tiene un conjunto de problemas para resolver numéricamente usando los métodos que se han revisado

PROBLEMAS NUMÉRICOS

 

 

Estimado lector, permítame sugerirle la entrada

ECUACIONES DIFERENCIALES Y MATLAB

que se encuentra en este blog y que espero también le sea de ayuda.

APROXIMACIONES NUMÉRICAS DE SOLUCIONES DE ECUACIONES DIFERENCIALES MEDIANTE DESARROLLOS DE TAYLOR

 

Un tema que no se trata es el de la relación entre soluciones numéricas de una ecuación diferencial y el desarrollo de Taylor para una función.

Antes de ver dicha relación retomaremos aspectos que se usarán sobre el desarrollo de Taylor de una función y aproximaciones de Taylor de orden1,2 y 3 para una función dada.

APROXIMACIONES DE TAYLOR DE ORDENES 1,2 Y 3.

Recordemos que si en un punto x_0 conocemos información de todos los órdenes para una función y; es decir,  conocemos y(x_0), y^{(1)}(x_0), \dots,y^{(n)}(x_0),\dots entonces podemos obtener los valores de  y para puntos x^* suficientemente cercanos a x_0. Más precisamente,

y(x^*)=y(x_0)+y^{(1)}(x_0)\cdot (x^*-x_0)+\frac{1}{2!}y^{(2)}(x_0)\cdot (x^*-x_0)^2+\frac{1}{3!}y^{(3)}(x_0)\cdot (x^*-x_0)^3+\dots+\frac{1}{n!}y^{(n)}(x_0)\cdot (x^*-x_0)^n+\dots

Al tomar h=x^*-x_0, la expresión anterior toma la forma

y(x^*)=y(x_0)+y^{(1)}(x_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2+\frac{1}{3!}y^{(3)}(x_0)\cdot h^3+\dots\linebreak+\frac{1}{n!}y^{(n)}(x_0)\cdot h^n+\dots.

Otra forma de escribir la expresión anterior, cuya obtención queda fuera de los alcances de esta sección (el lector interesado puede consultar pràcticamente cualquier texto de Análisis Matemático) está dada por

y(x^*)=y(x_0)+y^{(1)}(x_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2+\frac{1}{3!}y^{(3)}(x_0)\cdot h^3+\dots\linebreak+\frac{1}{n!}y^{(n)}(x_0)\cdot h^n+ \frac{1}{n!}\int_{x_0}^{x^*}(x^*-t)^n y^{(n+1)}(t)dt.

donde el último sumando \frac{1}{n!}\int_{x_0}^{x^*}(x^*-t)^n y^{(n+1)}(t)dt se conoce como el residuo.Este residuo no es más que la diferencia entre y(x^*) y la expresión

y(x_0)+y^{(1)}(x_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2+\frac{1}{3!}y^{(3)}(x_0)\cdot h^3+\dots\linebreak+\frac{1}{n!}y^{(n)}(x_0)\cdot h^n.

Esta última expresión es conocida como la aproximación de orden n del valor y(x^*). Esto lo que escribimos como

y(x^*)\cong y(x_0)+y^{(1)}(x_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2+\frac{1}{3!}y^{(3)}(x_0)\cdot h^3\linebreak+\dots+\frac{1}{n!}y^{(n)}(x_0)\cdot h^n.

De esta forma, tenemos las aproximaciones de los órdenes siguientes

y(x^*)\cong y(x_0)+y^{(1)}(x_0)\cdot h                                (aproximación de orden 1)

y(x^*)\cong y(x_0)+y^{(1)}(x_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2          (aproximación de orden 2)

y(x^*)\cong y(x_0)+y^{(1)}(x_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2+\frac{1}{3!}y^{(3)}(x_0)\cdot h^3          (aproximación de orden 3).

APROXIMACIONES NUMERICAS DE SOLUCIONES UANDO DESARROLLO DE TAYLOR DE ORDEN 1

A continuación mostraremos como dar la descripción del desarrollo de Taylor de orden 1

y(x^*)\cong y(x_0)+y^{(1)}(x_0)\cdot h          (1)

para obtener una aproximación del valor y(x^*)  (con x^* suficientemente cercano a x_0) que toma la solución y(x) de la ecuación diferencial

y'=f(x,y)     (2)

que satisface la condiciónn inicial y(x_0)=y_0.

Observemos que de (1), para obtener la aproximación al valor  y(x^*) se conoce el valor y(x_0) es cual es igual a y_0, y solamente resta por conocer el valor y'(x_0).

Al ser y(x) solución de (2), satisface la igualdad

y'(x_0)=f(x_0,y(x_0))

pero y(x_0)=y_0. Así,

y'(x_0)=f(x_0,y(x_0))=f(x_0,y_0). Es decir, y'(x_0)=f(x_0,y_0), y al sustituir esta igualdad en (1), obtenemos

y(x^*)\cong y(x_0)+y^{(1)}(x_0)\cdot h=y(x_0)+f(x_0,y_0)\cdot h,

con lo que se tiene la aproximación

y(x^*)\cong y(x_0)+f(x_0,y_0)\cdot h,,

la cual corresponde al método de Euler.

 

Es decir, el método de Euler es la aproximación para la solución de la ecuación diferencial que se obtiene usando el desarrollo de Taylor de orden 1.

 

APROXIMACIONES NUMERICAS DE SOLUCIONES UANDO DESARROLLO DE TAYLOR DE ORDEN $$

Del inciso anterior tenemos que al usar la igualdad  y'(x_0)=f(x_0,y_0) en el desarrollo de Taylor de orden 2

y(x^*)\cong y(x_0)+y^{(1)}(x_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2,

obtenemos

y(x^*)\cong y(x_0)+f(x_0,y_0)\cdot h+\frac{1}{2!}y^{(2)}(x_0)\cdot h^2,          (3)

y ahora solamente resta obtener una expresión para y^{(2)}(x_0), la cual se logra calcular a partir de la misma ecuación diferencial

y'(x)=f(x,y(x))

y el uso de la ragla de la cadena

Así,

y''(x_0)=\frac{df(x_0,y(x_0))}{dx}=\frac{\partial f \left(x_0,y(x_0)\right)}{\partial x}+\frac{\partial f\left(x_0,y(x_0)\right)}{\partial y} y'(x_0)=\frac{\partial f \left(x_0,y_0\right)}{\partial x}+\frac{\partial f\left(x_0,y_0\right)}{\partial y} f(x_0,y_0).

Es decir

y''(x_0)=\frac{\partial f \left(x_0,y_0\right)}{\partial x}+\frac{\partial f\left(x_0,y_0\right)}{\partial y} f(x_0,y_0).     (4).

 

Finalmente, sustituyendo (4) en (3), obtenemos

y(x^*)\cong y(x_0)+h\cdot f(x_0,y_0) +\frac{h^2}{2!}\left[\frac{\partial f \left(x_0,y_0\right)}{\partial x}+\frac{\partial f\left(x_0,y_0\right)}{\partial y} f(x_0,y_0)\right].       (5)

Esta expresión corresponde a la aproximación de la solución de la ecuación diferencial usando el desarrollo de Taylor de orden 2.

 

 

El algoritmo asociado a (5) para obtener la familia de aproximaciones \left\{y_i\right\}_{i=0,\dots, n} a los valores \left\{y(x_i)\right\}_{i=1,\dots, n} que la solución y(x) toma en los puntos a=x_0,x_1,\dots,x_n=b corresponde a

y_{i+1}= y(x_i)+h\cdot f(x_i,y_i) +\frac{h^2}{2!}\left[\frac{\partial f \left(x_i,y_i\right)}{\partial x}+\frac{\partial f\left(x_i,y_i\right)}{\partial y} f(x_i,y_i)\right].       (6)

Modelado y solucion de ecuaciones diferenciales de primer orden

LEY DE ENFRIAMIENTO DE NEWTON: La tasa de variación de la temperatura T (con respecto del tiempo) de un objeto es proporcional a la diferencia de la temperatura  del medio ambiente T_{amb} donde se encuentra inmerso y de la temperatura T(t) del objeto.

En un laboratorio se desea conocer la tasa de variación de temperatura de un cierto material cuando se introduce en un medio ambiente de temperatura T_{amb}=20^\circ \;C, ya que una vez determinada la tasa de variación, se puede conocer el material de que este está hecho. Así, se idea un experimento. Se calienta el objeto a una temperatura de 70^\circ \;, después de 20 minutos, su temperatura ha bajado a 50^\circ \, C.

Modele mediante una ecuación diferencial con condición inicial este problema y resuelve la ecuación diferencial con la condición inicial que ha establecido. De acuerdo con sus resultados, ¿diría usted que el objeto es un buen conductor de calor ? ¿sería preferible usar este material como aislante?

MODELOS DE CRECIMIENTO DE POBLACIÓN

MODELO MALTHUSIANO

En este modelo de crecimiento poblacional se supone que la tasa de variación instantánea de una población es proporcional, en todo momento, a la población existente. Este tipo de modelo se usan para estudiar, por ejemplo, la forma como cambia la población de un cultivo de bacterias cuando estas se repoducen por división celular. Se supone asimismo que las bacterias tienen condiciones adecuadas para su reproducción. La temperatura es apropiada, tienen espacio vital y alimento suficiente.

Si P(t) denota la cantidad de bacterias que hay en el cultivo en el instante t, se tiene que la variación instantánea de esta, dada por \frac{dP}{dt}, es proporcional a la población. Por lo que

\frac{dP}{dt}=kP(t).

Si la población inicial de bacterias es P_0,  usando el método de variables separables, el lector puede mostrar que al tomar P(0)=P_0 como condición inicial, entonces la solución particular está dada por la función

P(t)=P_0e^{kt}.

Es de notar que en este modelo no han sido tomados en cuenta muchos factores que pueden afectar el crecimiento de una población. Entre ellos podemos mencionar una tasa de mortandad por causas naturales, enfermedades, falta de alimento, sobrepoblación, entre otros.

Una primera mejora en nuestro modelo, olvidando que en particular hemos considerado una población de bacterias. Es tomar en cuenta que además de una tasa de generación de nuevos individuos por unidad de tiempo, también hay una tasa de fallecimientos (por causas naturales) por unidad de tiempo.

Si la tasa de nacimientos por unidad de tiempo es proporcional a la población, con constante de proporcionalidad \alpha y la tasa de muertes, en todo instante, también es proporcional a la cantidad de individuos, con constante de proporcionalidad \beta, entonces el modelo poblacional est\’a dao por

\frac{dP}{dt}=\alpha P -\beta P=(\alpha-\beta)P.

Un primer ejercicio consiste en mostrar que si \alpha>\beta entonces la población crecerá exponencialmente, mientras que en el caso \alpha<\beta, tendremos una extinción exponencial de la población.

MODELO LOGÍSTICO

Un modelo más realista que el malthusiano para el estudio de la evolución de una población es el conocido como modelo logístico, dado por la ecuación diferencial \frac{dP}{dt}=kP\left(1-\frac{P}{K}\right). Un modelo más realista que el malthusiano para el estudio de la evolución de una población es el conocido como modelo logístico, dado por la ecuación diferencial   \frac{dP}{dt}=kP\left(1-\frac{P}{K}\right)., donde la constante K se cono como la capacidad de carga del medio ambiente. Una forma de pensar esta ecuación es partir de un modelo malthusiano \frac{dP}{dt}=\alpha P y considerar que no existe espacio vital o alimento suficiente que permita continuar a la población con un crecimiento malthusiano o exponencial. Tomando en consideración estas hipótesis, se dará un enfrentamiento (por el alimento o por el espacio vital)   de un individuo de la población con otros, de estos enfrentamientos, se producirá la muerte de una parte proporcional de la población, \beta P.  Al tomar en cuenta a todos los individuos de la población y todos los posibles encuentros, se producirá la muerte de  \beta P^2 elementos de la población. Por lo que, el modelo que toma encuenta estas consideraciones sobre la restricción de alimento o de espacio vital está dado por \frac{dP}{dt}=\alpha P-\beta P^2.

  • Suponiendo que una población tiene originalmente un crecimiento malthusiano con constante de crecimiento dado por \alpha=0.2 y constante de decrecimiento \beta=0.00081 debido a encuentros de los individuos al luchar por alimento o espacio vital, obtenga una aproximación de la capacidad de carga del medio ambiente. Suponga que la cantidad de individuos en un instante t=0 es igual a 300 individuos. Establezca un modelo de población con estos datos en términos de una ecuación diferencial de primer orden con condición inicial. Resuelva la ecuación diferencial usando el comando dsolve y grafique la solución. obtenida. Para la condición inicial dada, ¿la población aumentará?, ¿disminuirá?.

Para lograr modelo más realistas de dinámica de poblaciones sería necesario considerar otros elementos que pueden influir en ella, como la aparición de enfermedades o epidemias, fenómenos de emigración o inmigración. Sólo por mencionar algunos.

Ecuaciones diferenciales de primer orden y ode45

Algunas veces las ecuaciones diferenciales que se obienen al modelar algún fenómeno no son fáciles de resolver o simplemente su presentación tiene una forma compleja, de modo que hallar soluciones explícitas o implícitas de manera analítica se convierte en un problema difícil de resolver o de plano la ecuación diferencial no se puede resolver mediante alguno de los métodos conocidos.

En una situación de este estilo solo queda el camino de los métodos numéricos. El lector se dirá :“¡pero, yo no sé nada de métodos numéricos!”. Aún en una situación así, Matlab tiene ya establecidas rutinas numéricas que podemos usar y nos permiten aproximar soluciones para un problema con condición inicial.

Por ejemplo, consideremos el problema un problema de dinámica de poblaciones que obedece a un modelo logístico

\frac{dP}{dt}=0.0025P(100-P), \hskip1cm P(0)=0.01

aún cuando podemos resover esta ecuación diferencial por el método de variables separables, usemos la instrucción ode45 que nos permitirá obtener una aproximación numérica de la solución a este problema y graficar la solución en un intervalo que consideremos nos puede dar información sobre la solución. Para lograr esto, podemos introducir la siguiente sucesión de instrucciones en la ventana de comandos

>> f = @(x,P) 0.0025*P*(100-P);

>>Intervalo=[0 70];

>>CondInic=0.01;

>>[x,y]=ode45(f,Intervalo,CondInic);

>>plot(x,y)

O en su caso, crear un programa en un archivo tipo script que contenga la misma sucesión de instrucciones

f = @(t,P) 0.0025*P*(100-P);

% función f(t,P)=0.0025*P*(100-P)
Intervalo=[0 70];

% intervalo sobre el cual se realizará la aproximación numérica
CondInic=0.01;

% Condición inicial de la ec’n diferencial P(0)=0.01
[t,P]=ode45(f,Intervalo,CondInic);

% gneración de vectores t y P
plot(t,P, ‘r’,’LineWidth’,2)

% se grafican los puntos (t,P) para obrener la grafica de la                                                %aproximación de la solucion del problema con condición inicial
grid on

% se grafica el enmallado
axis([0 70 -5 105])

% la grafica aparecerá en la ventana dada pr el rectangulo

% [0, 70]X[-5, 105]
print djpeg100 PoblacionLogistica.jpg 

% se guarda la gráfica en el  archivo  PoblacionLogistica.jpg

Recordemos que las expresiones del tipo

% Condición inicial de la ec’n diferencial P(0)=0.01

son comentarios sobre cada una las instrucciones y no forman parte del programa.

El programa anterior genera la gráfica siguiente

PoblacionLogistica

MÁS SOBRE DINÁMICA DE POBLACIONES

Para este modelo supongamos que en cierta región geográfica  hay una población de 60 millones habitantes para el año 2000 y que la tasa de crecimiento anual es de 2.5% anual. A partir de ese año se da un nuevo fenómeno de emigración e inmigración . cuya tasa de variación al año t debida a dicho fenómeno es igual a I(t).

Para modelar la tasa de variación de la población consideraremos el hecho de que la tasa de variación de la población al tiempo t  es igual a

la tasa de crecimiento +la tasa de variación debida a los fenómenos de emigración e inmigración

De esta forma, siendo \frac{dP}{dt} la tasa de variación de la población anual, tenemos que esta debe ser igual a 0.025P(t)+I(t), pues la población crece a una tasa anual de 2.5 \% anual e I(t) es la tasa de variación de la población debida a los fenómenos de emigración e inmnigración es igual a I(t).

Así, el modelo de la población con condición inicial está dado por la ecuación diferencial

\frac{dP}{dt}=0.025\; P(t)+\; I(t),\;\qquad P(2000)=6\cdot 10^7.

 Por ejemplo, si la tasa de variación de la población debida a los procesos de emigración e inmigración está dada por I(t)=0.02t, entonces tenemos el problema con valor inicial dado por

\frac{dP}{dt}=0.025\; P(t)+\; 0.02t,\;\qquad P(2000)=6\cdot 10^7.

Observemos que, por la forma que tiene I(t)=0.02t la inmigración es más mayor que la emigración de individuos de la población y el modelo es una ecuación lineal de primer orden no homogénea cuya solución puede obtenerse mediante el método del factor integrante..

EPIDEMIOLOGÍA:

Consideremos que en una población fija P_0, se expande una epidemia. Deseamos establecer hipótesis bajo las cuales sea posible dar un primer modelo para la forma en que crece la enfermedad entre la población.

Dividamos la población total en dos subconjuntos disjuntos, uno de ellos corresponde a una población que es susceptible de infectarse de la enfermedad. A esta población susceptible la denotaremos por S(t) que nos da el número de individuos susceptibles de contraer la enfermedad en el instante t . Por otra parte, el segundo sunconjunto de la población total corresponde a los individuos de la población que han contraido la enfermedad, es decir a la cantidad de individuos que han sido infectados. A esta población la denotaremos por I(t).

De esta forma, se tiene que la población total P_0 es la suma de la cantidad de individuos susceptibles S(t) y la cantidad de individuos infectados I(t); es decir,

P_0=S(t)+I(t). \hskip1cm (1)

De esta última ecuación tenemos que la cantidad de individuos que aín no han sdo inefcados es igual a la diferencia de la población total y la cantidad de sujetos infectados. En conscuencia,

S(t)=P_0-I(t) \hskip1cm (2).

Experimentalmente, se ha determinado que en todo instante,  la tasa a la cual crece la cantidad de individuos infectados es proporcional, tanto a la cantidad de personas que han contraido la enfermedad, I(t), como al número de personas susceptibles de contraerla S(t)=P_0-I(t) (por la ecuación (2)). Por lo tanto, existe una constante \alpha tal que

\frac{dI}{dt}=\alpha I(t)S(t)=\alpha I(t)(P_0-I(t)).

En resumen, el crecimiento de la cantidad de individuos infectados I(t) debe satisfacer la ecuación diferencial de primer orden

\frac{dI}{dt}=\alpha I(t)(P_0-I(t)).

Observemos que esta es una ecuación diferencial con variables separadas. Si denotamos por $I_0$ a la cantidad inicial de individuos infectados, obtenemos el problema con condición inicial

\frac{dI}{dt}=\alpha I(t)(P_0-I(t)), \hskip1cm I(0)=I_0.

Observe la gran semejanza con la ecuación logística de crecimiento de poblaciones.

Ejercicio. Resuelva este problema con condición inicial y encuentre una expresión explícita para la solución.

Ejercicio. Encuentre \lim_{t\to\infty} I(t). ¿Este límite es independiente de la cantidad inicial de individuos infectados? De una interpretación de su respuesta.

MEZCLAS:

Mezclador

Consideremos un reactor de mezclado continuo que contiene un volumen inicial V_0 de una mezcla obtenida con un soluto y solvente (pueden ser sal y agua pura, por ejemplo). Supongamos que el reactor tiene conectados dos tubos, A y B,  el primero permite entrada otra mezcla (salina) obtenida con los mismos tipos de soluto y solvente, mientras que el segundo permite la salida de la mezcla. También se considera que al momento de entrar una mezcla (salina) por el tubo A, este se combina homogéneamente  con la existente en el reactor. Esta mezcla homogénea es la que sale por el tubo B del reactor.

Sea S(t) la cantidad de soluto (medido en kg) en el tanque al tiempo t. La cantidad de soluto en el tanque variará dependiendo de la cantidad inicial S_0 de soluto contenida en el volumen V_0

Supongamos que la solución que entra por el tubo A tiene una concentración constante de C_1 Kg/litro y además entran L_1 litros/seg de solución al tanque por este tubo. No olvidemos que al entrar esta solución al tanque se mezcla de modo homogéneo por agitación con la mezcla que había previamente en el tanque, de modo que esta nueva mezcla homogénea sale del tanque de mezclado continuo a una razón de L_2 litros por segundo.

Si la cantidad de soluto contenida en el tanque al tiempo t  es de S(t) Kg y V(t) es el volumen de la solución, medida en litros, entonces para cualquier instante la concentración de soluto contenida en el tanque está dada por

\frac{S(t)}{V(t)} Kg/litro

Al entrar solución al tanque y salir solución del mismo, el volumen de solución V(t) contenido en el tanque en todo momento  será igual, mayor o menor al volumen inicial de solución V_0 contenido en el tanque, dependiendo de si la tasa de entrada L_1 de solución al tanque es igual a la tasa de salida L_2 de solución (L_1=L_2), si es mayor a esta (L_1>L_2) o menor a ella (L_1<L_2).

En todo caso, el volumen V(t) de solución  contenido en el tanque, en todo instante está dado por

V(t)=V_0+(L_1-L_2)t

Por otra parte, la variación de la cantidad de soluto en el tanque está dada por la diferencia de la cantidad de soluto que entra por unidad de tiempo y la cantidad de soluto que sale, también por unidad de tiempo.

Pero

      cantidad de soluto que entra por unidad de tiempo=C_1 Kg/litro\cdot L_1 litro/min=C_1L_1Kg/min

y

cantidad de soluto que sale por unidad de tiempo=C_2 Kg/litro\cdot L_2 litro/min=C_2L_2Kg/min

Por lo que

\frac{dS(t)}{dt}=C_1L_1-C_2L_2.

pero la concentración C_2(t) de soluto que sale por unidad de tiempo está dada por el cociente de la cantidad de soluto que hay en el tanque al tiempo t y el volumen de solución contenido en el tanque en el mismo instante t; es decir,

C_2(t)=\frac{S(t)}{V(t)}=\frac{S(t)}{V_0+(L_1-L2)t}.

Por lo tanto

\frac{dS(t)}{dt}=C_1L_1-\frac{S(t)}{V_0+(L_1-L2)t}L_2=C_1L_1-\frac{L_2}{V_0+(L_1-L2)t}S(t)

o equivalentemente,

\frac{dS(t)}{dt}+\frac{L_2}{V_0+(L_1-L2)t}S(t)=C_1L_1.

De esta forma, el modelo matemático está dado por el problema con condición inicial

\frac{dS(t)}{dt}+\frac{L_2}{V_0+(L_1-L2)t}S(t)=C_1L_1, \hskip1cm S(0)=S_0.

En las siguiente liga, estimad@s lector@s de este blog, encontrarán una lista de problemas de aplicaciones de ecuaciones de primer orden

ProbsAplicaciones1erOrden

Existencia y unicidad de soluciones de Ecuaciones Diferenciales Ordinarias de primer orden

Una herramienta importante de las ecuaciones diferenciales ordinarias es el llamado teorema de existencia y unicidad de soluciones para ecuaciones diferenciales ordinarias.

Supongamos que la  población P(t) de una especie obedece a un modelo dado por un problema con condición inicial dado por

\frac{dP}{dt}=f(P,t), \hskip1cm P(t_0)=P_0.

Si el problema tuviera dos soluciones P_1(t) y P_1(t), tendríamos dos comportamientos diferentes para la población, es decir para tiempos diferentes de t=0, tendríamos dos posible valores valores de la población. ¿ cuál de ellos sería la posible?.  Una situación de esta forma no sería deseable. Sería mejor poder hallar condiciones sobre la función f(P,t) que nos garantizara que para la condición inicial,

  • El problema con condición inicial \frac{dP}{dt}=f(P,t), \hskip1cm P(t_0)=P_0 tenga al menos una solución \phi(t)
  • que la solución \phi(t) sea única.

De esta forma debemos responder a dos interrogantes. La respuesta a la primera cuestión sobre la existencia de  al menos una solución para el problema con condición inicial la da el siguiente resultado

TEOREMA DE EXISTENCIA: Dado el problema con condición inicial \frac{dy}{dt}=f(y,t), \hskip1cm y(t_0)=y_0, si la función f es continua en una región \Omega del plano y (t_0,y_0)\in \Omega, entonces existen \epsilon>0 y una funciòn \phi:(t_0-\epsilon,t_0+\epsilon)\rightarrow R que resuelve el problema con condición inicial. Es decir

\frac{d\phi(t)}{dt}=f(t,\phi(t)) para cada t\in (t_0-\epsilon,t_0+\epsilon)
y
\phi(t_0)=y_0.

Ejemplo: Consideremos la ecuación diferencial con condición inicial

\frac{dy}{dt}=\frac{y}{t}\hskip1cm y(2)=3.

En este caso f(t,y)=\frac{y}{t} y (t_0,y_0)=(2,3). Además la función f es continua en la región \Omega dada por el semiplano superior y   el punto (t_0,y_0)=(2,3)\in \Omega. De esta forma se satisface la condición del teorema de existencia y podemos garantizar que el problema con condición inicial tiene al menos una solución. Observemos que el teorema no nos garantiza que la solución, que hemos asegurado existe, sea única.

El siguiente resultado nos da condiciones que debe satisfacer la función f(t,y) en el punto (t_0,y_0) para que no solo podamos asegurar la existencia de al menos una solución al problema con condición inicial, sino que también podamos garantizar que la solución existente sea única.

TEOREMA DE UNICIDAD: Para el problema con condición inicial \frac{dy}{dt}=f(y,t), \hskip1cm y(t_0)=y_0, si las funciones f y \frac{\partial f}{\partial y} son continuas en una región \Omega del plano y (t_0,y_0)\in \Omega, entonces existen

\epsilon>0

y una única función

\phi:(t_0-\epsilon,t_0+\epsilon)\rightarrow R

que resuelve el problema con condición inicial. Es decir, tal que

\frac{d\phi(t)}{dt}=f(t,\phi(t))

para cada t\in (t_0-\epsilon,t_0+\epsilon)

y

\phi(t_0)=y_0.

Ejemplo: Consideremos el problema con condición inicial

\frac{dy}{dt}=\sqrt{y-t}, \hskip1cm y(1)=5.

Es de notar que f(t,y)=\sqrt{t+y} y que el dominio donde esta función es continua es la región del plano

\Omega=\left\{ (t,y) \; | \; y-t\> 0 \right\}=\left\{ (t,y) \; |\; y> t \right\}

el cual es el conjunto de puntos de puntos del plano que se encuentran por encima de la recta identidad y=t.

Además la función \frac{\partial f}{\partial y}=\frac{1}{2\sqrt{y-t}} también es continua en la misma región \Omega.

De esta forma ambas funciones f y \frac{\partial f}{\partial y} son continuas en la región \Omega. Por lo tanto, para cada condición inicial (t_0,y_0)\in \Omega, el problema con condición inicial

\frac{dy}{dt}=\sqrt{y-t}, \hskip1cm y(t_0)=y_0.

tiene solución y esta es única. En particular, existe solución al problema con condición inicial

\frac{dy}{dt}=\sqrt{y-t}, \hskip1cm y(1)=5,

y la solución es única.

EJERCICIOS: CONSIDERE CADA UNA DE LAS SIGUIENTES ECUACIONES DIFERENCIALES Y DETERMINE LOS SUBCONJUNTOS DEL PLANO (t,y) donde podamos asegurar

  1. Para cada condición inicial en tal subconjunto la ecuación diferencial tiene al menos una solución que satisface tal condición inicial
  2. Dado un punto (condición inicial) en tal subconjunto, la solución de la ecuación diferencial que satisfae la condición inicial es única.
  • y'=\frac{1}{2t+y}.
  • y'=t^{1/3}.
  • y'=y^{1/3}.

Consideremos la ecuación diferencial y'=y^{1/5}. En este caso f(t,y)=y^{1/5} y esta función es continua en todo el plano R^2. Por el teorema de existencia de soluciones de ecuaciones diferenciles,  para cualquier punto (condición inicial) dado por (t_0,y_0) en el plano existe al menos una solución de la ecuación diferencial que pasa por el punto (t_0,y_0). Observemos que la continuidad de la función f no es suficiente para garantizar que la solución existente sea única. Solamente podemos asegurar la existenca de al menos una solución que satisface la condición inicial.

Para asegurar la unicidad de la solución necesitamos verificar la condición que la función \frac{\partial f}{\partial y} sea también continua en el punto (t_0,y_0).

Pero

\frac{\partial f}{\partial y}= \frac{\partial y^{1/5}}{\partial y}=\frac{1}{5} \frac{1}{y^{4/5}}

es continua en el plano R^2 salvo cuando y=0. Es decir,  en el plano, excepto el eje x.

En consecuencia de lo anterior, si la condición inicial (t_0,y_0) no se encuentra sobre el eje x,por ejemplo (t_0,y_0)=(2,3), entonces al ser  las funciones f y \frac{\partial f}{\partial y}=\frac{1}{5} \frac{1}{y^{4/5}} continuas en el punto (t_0,y_0) podemos garantizar no solo que al menos existe una solución que pasa por ese punto, sino también que dicha solución debe ser única.

ProbsExistenciaYunicidad