Elegir un solver de ODE

Ecuaciones diferenciales ordinarias

Una ecuación diferencial ordinaria (ODE) contiene una o más derivadas de una variable dependiente y con respecto a una única variable independiente t, normalmente denominada tiempo. La notación que se utiliza aquí para representar derivadas de y con respecto a t es $$y\text{'}$$ para la primera derivada, $$y\text{'}\text{'}$$, etc. El orden de la ODE equivale a la derivada del orden más alto de y que aparece en la ecuación.

Por ejemplo, esta es una ODE de segundo orden:

$$y\text{'}\text{'}=9y$$

En un problema de valores iniciales, la ODE se resuelve empezando por el estado inicial. Con la condición inicial $$y_0$$ y un periodo de tiempo en el que se va a obtener la respuesta, $$\left(t_0,t_f\right)$$, la solución se consigue de forma iterativa. En cada paso, el solver aplica un algoritmo determinado a los resultados de los pasos anteriores. En el primer paso, la condición inicial proporciona la información necesaria para que proceda la integración. El resultado final es que el solver de ODE devuelve un vector de unidades de tiempo $$t=\left[t_0,t_1,t_2,\mathrm{...},t_f\right]$$ y la solución correspondiente en cada paso $$y=\left[y_0,y_1,y_2,\mathrm{...},y_f\right]$$.

Tipos de ODE

Los solvers de ODE de MATLAB^® resuelven estos tipos de ODE de primer orden:

• ODE explícitas con formato $$y\text{'}=f\left(t,y\right)$$.

• ODE implícitas de forma lineal con formato $$M\left(t,y\right)y\text{'}=f\left(t,y\right)$$, donde $$M\left(t,y\right)$$ es una matriz de masa no singular. La matriz de masa puede depender del estado o del tiempo, o puede ser una matriz constante. Las ODE implícitas de forma lineal conllevan combinaciones lineales de la primera derivada de y, que están codificadas en la matriz de masa.

  Las ODE implícitas de forma lineal siempre se pueden transformar a una forma explícita, $$y\text{'}=M^{-1}\left(t,y\right)f\left(t,y\right)$$. En cambio, al especificar la matriz de masa directamente en el solver de ODE se evita esta transformación, lo que no es práctico y puede ser computacionalmente costoso.

• Si faltan componentes de $$y\text{'}$$, entonces las ecuaciones se llaman ecuaciones algebraicas diferenciales o DAE, y el sistema de DAE contiene variables algebraicas. Las variables algebraicas son variables dependientes cuyas derivadas no aparecen en las ecuaciones. Un sistema de DAE se puede reescribir como un sistema equivalente de ODE de primer orden tomando las derivadas de las ecuaciones para eliminar las variables algebraicas. El número de derivadas necesario para reescribir una DAE como ODE se llama índice diferencial. Los solvers ode15s y ode23t pueden resolver DAE con índice -1.

• ODE totalmente implícitas con formato $$f\left(t,y,y\text{'}\right)=0$$. Las ODE totalmente implícitas no se pueden reescribir de forma explícita, y también pueden contener variables algebraicas. El solver ode15i está diseñado para problemas totalmente implícitos, como DAE con índice -1.

Sistemas de ODE

Puede especificar cualquier número de ecuaciones ODE emparejadas para resolver y, en principio, el número de ecuaciones solo está limitado por la memoria disponible del ordenador. Si el sistema de ecuaciones tiene n ecuaciones,

$$\left(\begin{array}{c}y{\text{'}}_1\\ y{\text{'}}_2\\ ⋮\\ y{\text{'}}_n\end{array}\right)=\left(\begin{array}{c}{f}_1\left(t,y_1,y_2,\mathrm{...},y_n\right)\\ f_2\left(t,y_1,y_2,\mathrm{...},y_n\right)\\ ⋮\\ f_n\left(t,y_1,y_2,\mathrm{...},y_n\right)\end{array}\right),$$

la función que codifica las ecuaciones devuelve un vector con n elementos, que corresponden a los valores de $$y{\text{'}}_1,y{\text{'}}_2,\dots ,y{\text{'}}_n$$. Por ejemplo, considere el sistema de dos ecuaciones

$$\{\begin{array}{l}y{\text{'}}_1=y_2\\ y{\text{'}}_2=y_1{y}_2-2.\end{array}$$

Una función que codifica estas ecuaciones es

function dy = myODE(t,y)
  dy(1) = y(2);
  dy(2) = y(1)*y(2)-2;
end

ODE de orden superior

Los solvers de ODE de MATLAB solo resuelven ecuaciones de primer orden. Debe reescribir las ODE de orden superior como un sistema equivalente de ecuaciones de primer orden usando las sustituciones genéricas

$$\begin{array}{l}{y}_1=y\\ y_2=y\text{'}\\ y_3=y\text{'}\text{'}\\ ⋮\\ y_n=y^{(n-1)}.\end{array}$$

El resultado de estas sustituciones es un sistema de n ecuaciones de primer orden

$$\{\begin{array}{l}y{\text{'}}_1=y_2\\ y{\text{'}}_2=y_3\\ ⋮\\ y{\text{'}}_n=f\left(t,y_1,y_2,\mathrm{...},y_n\right).\end{array}$$

Por ejemplo, considere la ODE de primer orden

$$y\text{'}\text{'}\text{'}-y\text{'}\text{'}y+1=0.$$

Con las sustituciones

$$\begin{array}{l}{y}_1=y\\ y_2=y\text{'}\\ y_3=y\text{'}\text{'}\end{array}$$

da como resultado el sistema equivalente de primer orden

$$\{\begin{array}{l}y{\text{'}}_1=y_2\\ y{\text{'}}_2=y_3\\ y{\text{'}}_3=y_1{y}_3-1.\end{array}$$

Entonces, el código para este sistema de ecuaciones es

function dydt = f(t,y)
  dydt(1) = y(2);
  dydt(2) = y(3);
  dydt(3) = y(1)*y(3)-1;
end

ODE complejas

Considere la ecuación ODE compleja

$$y\text{'}=f\left(t,y\right),$$

donde $$y=y_1+i{y}_2$$. Para resolverla, separe las partes real e imaginaria en distintos componentes de la solución y vuelva a combinar los resultados al final. De forma conceptual, sería

$$\begin{array}{l}yv=\left[\text{Real}\left(y\right)\text{Imag}\left(y\right)\right]\\ fv=\left[\text{Real}\left(f\left(t,y\right)\right)\text{Imag}\left(f\left(t,y\right)\right)\right].\end{array}$$

Por ejemplo, si la ODE es $$y\text{'}=yt+2i$$, puede representar la ecuación con el archivo de función:

function f = complexf(t,y)
  f = y.*t + 2*i;
end

Así, el código para separar las partes real e imaginaria es

function fv = imaginaryODE(t,yv)
  % Construct y from the real and imaginary components
  y = yv(1) + i*yv(2);            

  % Evaluate the function
  yp = complexf(t,y);             

  % Return real and imaginary in separate components
  fv = [real(yp); imag(yp)]; 
end     

Cuando ejecuta un solver para obtener la solución, la condición inicial y0 también se separa en las partes real e imaginaria para proporcionar una condición inicial a cada componente de la solución.

y0 = 1+i;
yv0 = [real(y0); imag(y0)];
tspan = [0 2];
[t,yv] = ode45(@imaginaryODE, tspan, yv0);

Cuando obtenga la solución, combine los componentes de la parte real e imaginaria para obtener el resultado final.

y = yv(:,1) + i*yv(:,2);

Selección básica de solver

ode45 funciona bien con la mayoría de problemas de ODE y generalmente debería ser nuestra primera opción de solver. En cambio, ode23, ode78, ode89 y ode113 pueden ser más eficientes que ode45 para problemas con requisitos de precisión más o menos estrictos.

Algunos problemas de ODE pueden ser rígidos o difíciles de evaluar. La rigidez es un término sin una definición precisa pero, en general, ocurre cuando hay una diferencia para escalar en algún lugar del problema. Por ejemplo, si la solución de una ODE tiene dos componentes que varían drásticamente en escalas de tiempo, la ecuación puede ser rígida. Puede identificar un problema como rígido si los solvers no rígidos (como ode45) no pueden resolver el problema o son extremadamente lentos. Si observa que un solver no rígido es muy lento, intente usar un solver rígido como ode15s. Con un solver rígido, puede mejorar la fiabilidad y la eficiencia proporcionando la matriz jacobiana o su patrón de escasez.

Puede usar los objetos ode para automatizar la selección de solver en función de las propiedades del problema. Si no tiene claro qué solver usar, esta tabla muestra directrices generales sobre cuándo usar cada solver.

Solver	Tipo de problema	Precisión	Cuándo usarlo
ode45	No rígido	Media	Casi siempre. El primer solver que debe intentar es ode45.
ode23		Baja	ode23 puede ser más eficiente que ode45 en problemas con tolerancias aproximadas o con rigidez moderada.
ode113		De baja a alta	ode113 puede ser más eficiente que ode45 en problemas con tolerancias de error estrictas o cuando sea caro evaluar la función de la ODE.
ode78		Alta	ode78 puede ser más eficiente que ode45 en problemas con soluciones sencillas y requisitos de precisión altos.
ode89		Alta	ode89 puede ser más eficiente que ode78 en problemas muy sencillos, cuando se integra en largos intervalos de tiempo o cuando las tolerancias son especialmente estrictas.
ode15s	Rígido	De baja a media	Pruebe ode15s cuando ode45 falle o sea ineficiente y sospeche que el problema es rígido. Use también ode15s para resolver ecuaciones algebraicas diferenciales (DAE).
ode23s		Baja	ode23s puede ser más eficiente que ode15s en problemas con tolerancias de error aproximadas. Resuelve algunos problemas rígidos para los que ode15s no es efectivo. ode23s calcula la jacobiana en cada paso, por lo que es bueno proporcionar la jacobiana a través de odeset para maximizar la eficiencia y la precisión. Si hay una matriz de masa, debe ser constante.
ode23t		Baja	Utilice ode23t si el problema solo es moderadamente rígido y necesita una solución sin amortiguamiento numérico. ode23t resuelve ecuaciones algebraicas diferenciales (DAE).
ode23tb		Baja	Al igual que ode23s, el solver ode23tb puede ser más eficiente que ode15s en problemas con tolerancias de error aproximadas.
ode15i	Totalmente implícito	Baja	Utilice ode15i para problemas totalmente implícitos f(t,y,y’) = 0 y para ecuaciones algebraicas diferenciales (DAE) con índice 1.

Para más detalles y recomendaciones sobre cuándo usar cada solver, consulte [5].

Resumen de ejemplos y archivos de ODE

Hay varios ejemplos de archivos disponibles que sirven como excelentes puntos de partida para la mayoría de problemas de ODE. Para ejecutar la app Ejemplos de ecuaciones diferenciales que sirve para que pueda explorar y ejecutar ejemplos, escriba

odeexamples

Para abrir un archivo de muestra particular y editarlo, escriba

edit exampleFileName.m

Para ejecutar un ejemplo, escriba

exampleFileName

Esta tabla contiene una lista de archivos de muestra de DAE y ODE, y de los solvers y las opciones que utilizan. Se incluyen enlaces para el subconjunto de ejemplos que también están publicados directamente en la documentación.

Archivo de ejemplo	Solver utilizado	Opciones especificadas	Descripción	Enlace a la documentación
amp1dae	ode23t	'Mass'	DAE rígida: circuito eléctrico con una matriz de masa singular y constante	Solve Stiff Transistor Differential Algebraic Equation
ballode	ode23	'Events' 'OutputFcn' 'OutputSel' 'Refine' 'InitialStep' 'MaxStep'	Ubicación de eventos simple: pelota botando	Ubicación de eventos de EDO
batonode	ode45	'Mass'	ODE con matriz de masa dependiente del estado y el tiempo: movimiento de un bastón	Solve Equations of Motion for Baton Thrown into Air
brussode	ode15s	'JPattern' 'Vectorized'	Problema grande y rígido: difusión en una reacción química (el Brusselator)	Solve Stiff ODEs
burgersode	ode15s	'Mass' 'MStateDependence' 'JPattern' 'MvPattern' 'RelTol' 'AbsTol'	ODE con matriz de masa muy dependiente del estado: ecuación de Burgers resuelta con una técnica de malla en movimiento	Solve ODE with Strongly State-Dependent Mass Matrix
fem1ode	ode15s	'Mass' 'MStateDependence' 'Jacobian'	Problema rígido con una matriz de masa dependiente del tiempo: método de elementos finitos	—
fem2ode	ode23s	'Mass'	Problema rígido con una matriz de masa constante: método de elementos finitos	—
hb1ode	ode15s	—	Problema de ODE rígido resuelto en un intervalo muy largo: reacción química de Robertson	—
hb1dae	ode15s	'Mass' 'RelTol' 'AbsTol' 'Vectorized'	DAE linealmente implícita rígida de una ley de conservación: reacción química de Robertson	Solve Robertson Problem as Semi-Explicit Differential Algebraic Equations (DAEs)
ihb1dae	ode15i	'RelTol' 'AbsTol' 'Jacobian'	DAE totalmente implícita rígida: reacción química de Robertson	Solve Robertson Problem as Implicit Differential Algebraic Equations (DAEs)
iburgersode	ode15i	'RelTol' 'AbsTol' 'Jacobian' 'JPattern'	Sistema de ODE implícita: ecuación de Burgers	—
kneeode	ode15s	'NonNegative'	El "problema de la rodilla" con restricciones de no negatividad	Nonnegative ODE Solution
orbitode	ode45	'RelTol' 'AbsTol' 'Events' 'OutputFcn'	Ubicación de eventos avanzada: problema de tres cuerpos restringido	Ubicación de eventos de EDO
rigidode	ode45	—	Problema no rígido: ecuaciones de Euler de un cuerpo rígido sin fuerzas externas	Resolver ecuaciones diferenciales ordinarias (EDO) no rígidas
vdpode	ode15s	'Jacobian'	Ecuación van der Pol parametrizable (rígido para grandes μ)	Solve Stiff ODEs

Referencias

[1] Shampine, L. F. and M. K. Gordon, Computer Solution of Ordinary Differential Equations: the Initial Value Problem, W. H. Freeman, San Francisco, 1975.

[2] Forsythe, G., M. Malcolm, and C. Moler, Computer Methods for Mathematical Computations, Prentice-Hall, New Jersey, 1977.

[3] Kahaner, D., C. Moler, and S. Nash, Numerical Methods and Software, Prentice-Hall, New Jersey, 1989.

[4] Shampine, L. F., Numerical Solution of Ordinary Differential Equations, Chapman & Hall, New York, 1994.

[5] Shampine, L. F. and M. W. Reichelt, “The MATLAB ODE Suite,” SIAM Journal on Scientific Computing, Vol. 18, 1997, pp. 1–22.

[6] Shampine, L. F., Gladwell, I. and S. Thompson, Solving ODEs with MATLAB, Cambridge University Press, Cambridge UK, 2003.

Consulte también

ode | odeset | odextend

Temas relacionados

• Resolver ecuaciones diferenciales ordinarias (EDO) no rígidas
• Solve Stiff ODEs
• Solve Differential Algebraic Equations (DAEs)

Sitios web externos

• Ecuaciones diferenciales ordinarias