Métodos numéricos para ecuaciones diferenciales ordinarias

Los métodos numéricos para ecuaciones diferenciales ordinarias son procedimientos utilizados para encontrar aproximaciones numericas a las soluciones de ecuaciones diferenciales ordinarias (EDO). Su uso también se conoce como integración numérica, aunque este término a veces se toma para significar el cálculo de una integración.

Muchas ecuaciones diferenciales no pueden resolverse usando funciones típicas ("análisis"). Sin embargo, a efectos prácticos, como en ingeniería, una aproximación numérica a la solución suele ser suficiente. Los algoritmos estudiados aquí pueden usarse para calcular tal aproximación. Un método alternativo es utilizar técnicas de cálculo infinitesimal para obtener una expansión en serie de la solución.

Las ecuaciones diferenciales ordinarias se presentan en muchas disciplinas científicas, por ejemplo, en física, química, biología y economía. Además, algunos métodos en ecuaciones diferenciales parciales numéricas convierten una ecuación diferencial parcial en una ecuación diferencial ordinaria, que luego debe resolverse. (solve by: Wellington Castillo)

Una ecuación diferencial de primer orden es un problema de valor inicial (PVI) de la forma,^[1]​

donde f es una función que asigna [t₀,∞) × R^d a R^d, con la condición inicial y₀ ∈ R^d es un vector dado. Primer orden significa que solo la primera derivada de y aparece en la ecuación, y las derivadas más altas están ausentes.

Sin pérdida de generalidad en los sistemas de orden superior, en este artículo se restringe la explicación a las ecuaciones diferenciales de "primer orden", porque un EDO de orden superior se puede convertir en un sistema más grande de ecuaciones de primer orden mediante la introducción de variables adicionales. Por ejemplo, la ecuación de segundo orden y'' = -y puede reescribirse como dos ecuaciones de primer orden: y' = z y z' = -y.

En esta sección, se describen métodos numéricos para los PVI, teniendo en cuenta que los problemas de condición de frontera (PCF) requieren un conjunto diferente de herramientas. En un PCF, se definen valores o componentes de la solución y en más de un punto. Debido a esto, se deben usar diferentes métodos para resolverlo. Por ejemplo, el método de disparo (y sus variantes) o métodos globales como las diferencias finitas, el método de Galerkin o el método de colocación son apropiados para esa clase de problemas.

El teorema de Picard-Lindelöf establece que existe una solución única, siempre que f sea lipschitzianamente continua.

Los métodos numéricos para resolver PVI de primer orden a menudo se dividen en una de estas dos grandes categorías: método lineal multipaso o método de Runge-Kutta. Se puede lograr una separación adicional dividiendo los métodos en aquellos que son explícitos y aquellos que son implícitos. Por ejemplo, los métodos lineales multipaso implícitos incluyen el método de Adams-Moulton y la fórmula de diferenciación hacia atrás (FDA), mientras que el método de Runge-Kutta^[2]​ incluye Runge-Kutta diagonalmente implícito (RKDI), Runge-Kutta diagonalmente implícito simple (RKDIS) y Gauss-Radau (basado en la cuadratura gaussiana). Los ejemplos explícitos de la familia lineal multipaso incluyen el método lineal multipaso, y cualquier método de Runge-Kutta con una diagonal inferior es explícito. Una regla general suelta dicta que las ecuaciones diferenciales rígidas requieren el uso de esquemas implícitos, mientras que los problemas no rígidos se pueden resolver de manera más eficiente con esquemas explícitos.

Los llamados métodos lineales generales (MLG) son una generalización de las dos grandes clases de métodos anteriores.

Desde cualquier punto de una curva, se puede encontrar una aproximación de otro punto cercano en la curva moviéndose una corta distancia sobre una línea tangente a la curva.

Comenzando con la ecuación diferencial (1), se reemplaza la derivada y' por la aproximación respecto a una diferencia finita

que cuando se reorganiza produce la siguiente fórmula

y usando (1) da:

Esta fórmula generalmente se aplica de la manera que se explica a continuación.

Se elige el tamaño de paso h y se construye la secuencia t₀, t₁ = t₀ + h, t₂ = t₀ + 2h, ... Denotando por y_n una estimación numérica de la solución exacta y(t_n).

De acuerdo con (3), se calculan estas estimaciones mediante el siguiente esquema recursivo:${displaystyle y_{n+1}=y_{n}+hf(t_{n},y_{n}).qquad qquad (4)}$ Este es el método de Euler (en contraste con el método de Euler hacia atrás, que se describe a continuación). El método lleva el nombre de Leonhard Euler que lo describió en 1768.

Es un ejemplo de un método explícito. Esto significa que el nuevo valor y_n+1 se define en términos de datos que ya se conocen, como y_n.

Si, en lugar de (2), se usa la aproximación

se obtiene el "método de Euler hacia atrás":

Esto implica que se trata de un método implicito, lo que significa que previamente se debe resolver otra ecuación con el fin de encontrar y_n+1. Para ello, a menudo se usa el método del punto fijo o alguna modificación del método de Newton-Raphson.

Sin embargo, normalmente cuesta más tiempo resolver esta ecuación que los cálculos de los métodos explícitos. Este costo debe tenerse en cuenta cuando se selecciona el método a utilizar. La ventaja de los métodos implícitos como (6) es que generalmente son más estables para resolver una ecuación rígida, lo que significa que se puede usar un tamaño de paso h más grande.

Los integradores exponenciales describen una gran clase de métodos que han experimentado un gran desarrollo.^[3]​ Si origen se remonta al menos a la década de 1960.

En lugar de (1), se asume que la ecuación diferencial es cualquiera de la forma

o se ha linealizado localmente sobre una forma original para producir un término lineal ${displaystyle -Ay}$ y un término no lineal ${displaystyle {mathcal {N}}(y)}$.

Los integradores exponenciales se construyen multiplicando (7) por ${displaystyle e^{At}}$ e integrando exactamente el resultado sobre un intervalo de tiempo ${displaystyle [t_{n},t_{n+1}=t_{n}+h]}$:

Esta ecuación integral es exacta, pero no define la integral.

El integrador exponencial de primer orden se puede determinar manteniendo ${displaystyle {mathcal {N}}(y(t_{n}+ au ))}$ constante durante todo el intervalo:

El método de Euler a menudo no es lo suficientemente exacto. En términos más precisos, solo tiene orden uno (el concepto de orden se explica a continuación). Esto hizo que los matemáticos buscaran métodos de orden superior.

Una posibilidad es usar no solo el valor previamente calculado y_n para determinar y_n+1, sino hacer que la solución dependa de más valores calculados previamente. Esto produce el llamado "método de varios pasos". Quizás el más simple es el método del salto de rana, que es de segundo orden y (más o menos) se basa en dos valores previos cada vez.

Casi todos los métodos prácticos de varios pasos pertenecen a la familia del método lineal multipaso, que tienen la forma

Otra posibilidad es usar más puntos en el intervalo [t_n, t_n+1]. Esto lleva a la familia del método de Runge-Kutta, llamada así por Carl Runge y Martin Wilhelm Kutta. Uno de sus métodos de cuarto orden es especialmente popular.

Una buena implementación de uno de estos métodos para resolver una EDO a menudo implica algo más que elegir una fórmula y un paso de intervalo adecuados.

A menudo es ineficiente usar el mismo tamaño de paso todo el tiempo, por lo que se han desarrollado "métodos de tamaño de paso variable". Por lo general, el tamaño del paso se elige de modo que el error (local) por paso esté por debajo de cierto nivel de tolerancia. Esto significa que los métodos también deben calcular un "indicador de error", una estimación del error local.

Una extensión de esta idea es elegir dinámicamente entre diferentes métodos de diferentes órdenes (esto se llama un "método de orden variable"). Los métodos basados en la extrapolación de Richardson, como el algoritmo de Bulirsch-Stoer, a menudo se utilizan para construir varios métodos de diferentes órdenes.

Otras características deseables incluyen:

Algunos otros métodos no se detallan el presente artículo. Entre estos métodos alternativos, se encuentran:

Para aplicaciones que requieren computación paralela en supercomputadoras, el grado de complejidad requerido por un método numérico se vuelve relevante. En vista de los desafíos de los sistemas informáticos que trabajan a exaescala, se están estudiando métodos numéricos para problemas de valor inicial que puedan proporcionar concurrencia en la dirección temporal.^[4]​ Parareal es un ejemplo relativamente conocido de tal método de integración paralelo en el tiempo, pero las primeras ideas al respecto se remontan a la década de 1960.^[5]​

El análisis numérico no es solo el diseño de métodos numéricos, sino que también implica el estudio de cómo funcionan. Tres conceptos son centrales en este estudio:

Se dice que un método numérico es convergente si la solución numérica se aproxima a la solución exacta cuando el tamaño de paso h tiende a 0. Más precisamente, se requiere que para cada EDO (1) con una función lipschitziana f y cada t^* > 0,

Todos los métodos mencionados anteriormente son convergentes.

Supóngase que el método numérico es

El "error local (de truncamiento)" del método es el error cometido en un paso del método. Es decir, es la diferencia entre el resultado dado por el método, suponiendo que no se cometió ningún error en los pasos anteriores, y la solución exacta:

Se dice que el método es consistente si

El método tiene orden ${displaystyle p}$ si

Por lo tanto, un método es consistente si tiene un orden mayor que 0. El método de Euler (hacia adelante) (4) y el método de Euler hacia atrás (6) introducido anteriormente tienen orden 1, por lo que son consistentes. La mayoría de los métodos utilizados en la práctica alcanzan un orden superior. La consistencia es una condición necesaria para la convergencia ^{[cita requerida]}, pero no suficiente; Para que un método sea convergente, debe ser consistente y cero-estable.

Un concepto relacionado es el de "error global (de truncamiento)", el error acumulado en todos los pasos necesarios para alcanzar un tiempo fijo t. Explícitamente, el error global en el tiempo t es y_N − y(t) donde N = (t−t₀)/h. El error global de un método de un solo paso p es O (h^p); en particular, tal método es convergente. Esta afirmación no es necesariamente cierta para los métodos de varios pasos.

Para algunas ecuaciones diferenciales, la aplicación de métodos estándar, como el método de Euler, el método de Runge-Kutta explícito o el método lineal multipaso (por ejemplo, el método de Adams-Bashforth), exhiben inestabilidad en las soluciones, aunque otros métodos pueden producir soluciones estables. Este "comportamiento difícil" en la ecuación (que puede no ser necesariamente complejo en sí mismo) se describe como "rigidez", y a menudo es causado por la presencia de diferentes escalas de tiempo en el problema subyacente. Por ejemplo, una colisión en un sistema mecánico como en un oscilador de impacto generalmente ocurre en una escala de tiempo mucho más pequeña que la necesaria para el caracterizar el movimiento de los objetos. Esta discrepancia genera "virajes bruscos" en las curvas de los parámetros de estado.

Los problemas rígidos son ubicuos en cinética química, teoría del control, mecánica de sólidos deformables, pronóstico del tiempo, biología, plasma y electrónica. Una forma de superar la rigidez es extender la noción de ecuación diferencial a la de inclusión diferencial, que permite y modela la falta de suavidad.^[6]​^[7]​

A continuación se muestra una cronología de algunos desarrollos importantes en este campo:

Los problemas de valor límite (PVL) generalmente se resuelven numéricamente solucionando un problema de matriz aproximadamente equivalente, obtenido al discretizar el PVL original. El procedimiento más utilizado para resolver numéricamente estos problemas en una dimensión se llama método de las diferencias finitas. Este método aprovecha las combinaciones lineales de valores de puntos para construir coeficientes de diferencias finitas que describen derivadas de la función. Por ejemplo, la aproximación en diferencia finita de segundo orden a la primera derivada viene dada por:

y la diferencia finita de segundo orden para la segunda derivada viene dado por:

En ambas fórmulas, ${displaystyle h=x_{i}-x_{i-1}}$ es la distancia entre los valores vecinos x en el dominio discretizado. Luego se construye un sistema lineal que luego se puede resolver con el estándar matrix methods. Por ejemplo, supongamos que la ecuación a resolver es:

El siguiente paso sería discretizar el problema y usar aproximaciones lineales derivadas como

y resolver el sistema resultante de ecuaciones lineales. Esto llevaría a ecuaciones como:

A primera vista, este sistema de ecuaciones parece tener dificultades asociadas con el hecho de que la ecuación no implica términos que no se multipliquen por variables, pero en realidad esto es falso. En i=1 y nminus&;1 existe un término que involucra los valores límite ${displaystyle u(0)=u_{0}}$ y ${displaystyle u(1)=u_{n}}$ y dado que estos dos valores son conocidos, simplemente se puede sustituirlos en esta ecuación y, como resultado, obtener un sistema lineal de ecuaciones no homogéneo que tiene soluciones no triviales.

---