¿Por qué los métodos Runge – Kutta de orden superior no se usan con más frecuencia?

17

Solo tenía curiosidad por saber por qué los métodos Runge – Kutta de alto orden (es decir, mayores de 4) casi nunca se discuten / emplean (al menos que yo sepa). Entiendo que requiere un mayor tiempo de cómputo por paso (por ejemplo, RK14 con paso incorporado de 12 ° orden ), pero ¿hay alguna otra desventaja de usar métodos Runge-Kutta de orden superior (por ejemplo, problemas de estabilidad)? Cuando se aplica a ecuaciones con soluciones altamente oscilantes en escalas de tiempo extremas, ¿no se preferirían tales métodos de orden superior?

ode runge-kutta Mathews24
fuente
2
Creo que esta es una pregunta muy subjetiva. El mayor inconveniente que ya ha notado es el costo de la computación. Generalmente tratamos de equilibrar entre precisión y tiempo computacional. En los PDE, cuando las personas hablan de un orden superior, generalmente piensan en el tercer o cuarto orden. Y el paso del tiempo también se mantiene en el mismo orden.
Vikram
3
En PDE, un esquema de precisión de alto orden para la dependencia temporal no tiene sentido si la precisión espacial es peor. De hecho, la precisión de la dependencia espacial es principalmente de 2 ° o 3 ° orden, especialmente cuando se trabaja en mallas no estructuradas. La gente necesita controlar el truncamiento de error global con el menor costo, por lo tanto, considera Runge-Kutta con un orden de precisión suficientemente alto en casos particulares.
tqviet
@tqviet Si se usan aproximaciones de diferencia hacia atrás o centrales de hasta 8 para las derivadas espaciales, RK8 sería adecuado, ¿no? En general, ¿hay problemas de precisión o estabilidad con el uso de aproximaciones de diferencias finitas de alto orden de las derivadas espaciales?
Mathews24
1
@ Mathews24: No mencioné la estabilidad, que depende en gran medida de la ecuación. Cuando se aplica un esquema altamente preciso a la dependencia espacial, adoptamos RK a la dependencia temporal con al menos el mismo orden de precisión, pero la condición de estabilidad puede requerir un valor menor de . Δt
tqviet

Respuestas:

17

Existen miles de documentos y cientos de códigos que utilizan métodos Runge-Kutta de quinto orden o más. Tenga en cuenta que el integrador explícito más utilizado en MATLAB es ODE45, que avanza la solución utilizando un método Runge-Kutta de quinto orden.

Ejemplos de métodos de Runge-Kutta de alto orden ampliamente utilizados

El documento de Dormand & Prince que da un método de quinto orden tiene más de 1700 citas según Google Scholar . La mayoría de ellos son documentos que utilizan su método para resolver algún problema. El documento del método Cash-Karp tiene más de 400 citas . Quizás el método de orden superior a 5 más utilizado es el método de octavo orden de Prince-Dormand, que tiene más de 400 citas en Google Scholar . Podría dar muchos otros ejemplos; y tenga en cuenta que muchas (si no la mayoría) de las personas que usan estos métodos nunca citan los documentos.

Tenga en cuenta también que los métodos de extrapolación de alto orden y corrección diferida son métodos Runge-Kutta .

Métodos de alto orden y error de redondeo

Si su precisión está limitada por errores de redondeo, entonces debe usar un método de orden superior . Esto se debe a que los métodos de orden superior requieren menos pasos (y menos evaluaciones de funciones, aunque haya más evaluaciones por paso), por lo que cometen menos errores de redondeo. Puede verificar esto fácilmente con experimentos simples; Es un buen problema de tarea para un primer curso de análisis numérico.

Los métodos de décimo orden son extremadamente útiles en la aritmética de doble precisión. Por el contrario, si todo lo que tuviéramos fuera el método de Euler, entonces el error de redondeo sería un problema importante y necesitaríamos números de punto flotante de muy alta precisión para muchos problemas donde los solucionadores de alto orden funcionan bien.

Los métodos de alto orden pueden ser igual de estables

UNsi

Métodos de alto orden en mecánica celeste

Usted pregunta

Cuando se aplica a ecuaciones con soluciones altamente oscilantes en escalas de tiempo extremas, ¿no se preferirían tales métodos de orden superior?

Tienes toda la razón! Un excelente ejemplo de esto es la mecánica celeste. No soy un experto en esa área. Pero este documento , por ejemplo, compara métodos para la mecánica celeste y ni siquiera considera un orden inferior a 5. Concluye que los métodos de orden 11 o 12 son a menudo los más eficientes (con el método Prince-Dormand de orden 8 también a menudo muy eficiente).

David Ketcheson
fuente
Ketchson: ¿podría proporcionar alguna evidencia o explicación sobre esta afirmación: "los métodos de extrapolación de alto orden y corrección diferida son métodos de Runge-Kutta"? Especialmente los "métodos de corrección diferida". Gracias.
tqviet
@David Ketcheson ¿Puede analizar cómo cambiaría su respuesta si utilizara técnicas informáticas validadas (verificadas), como el intervalo redondeado hacia afuera o la aritmética radial? ¿Qué tal si se usara un intervalo redondeado hacia afuera de precisión doble superior o una aritmética radial? ¿Qué sucederá con el ajuste y la dependencia a medida que aumente el orden de Runge-Kutta, y solo por diversión, digamos que el ODE es muy rígido?
Mark L. Stone el
@ MarkL.Stone Ese es un conjunto completamente diferente de preguntas. Si desea hacerles preguntas, publíquelas como preguntas separadas. Sin embargo, no soy un experto en esas cosas y no podré responder.
David Ketcheson
1
@tqviet Eche un vistazo a este documento para obtener una explicación.
David Ketcheson el
12

Mientras utilice la aritmética de coma flotante de precisión doble estándar, no se necesitan métodos de orden muy alto para obtener una solución con alta precisión en un número razonable de pasos. En la práctica, encuentro que la precisión de la solución normalmente está limitada a un error relativo de 1.0e-16 por la representación de coma flotante de doble precisión en lugar del número / longitud de los pasos que se toman con RKF45.

Si cambia a un esquema aritmético de coma flotante de precisión superior a la doble, vale la pena usar un método de décimo orden.

Brian Borchers
fuente
55
Creo que esta respuesta es engañosa. Los métodos de alto orden conducen a un error de redondeo mucho menor, mientras que los métodos de bajo orden sufren un error de redondeo dominante cuando la precisión requerida es grande o el intervalo de tiempo es largo; mira mi respuesta a continuación.
David Ketcheson
2
El punto es que en coma flotante de doble precisión ni siquiera puede representar una solución con una precisión relativa mejor que 1.0e-16. En muchas situaciones prácticas, el viejo RKF45 lo llevará a ese nivel de precisión durante el período que le interese sin requerir pequeños pasos. Puede que no sea una buena opción para sistemas rígidos o situaciones donde se requiere un integrador simpléctico, pero un método de Runge Kutta de orden superior tampoco es una gran solución en esas situaciones. Estoy de acuerdo en que durante períodos de tiempo muy largos, los métodos de Runge Kutta de orden superior pueden tener algún sentido.
Brian Borchers el
10

Solo para agregar a la excelente respuesta de Brian Borcher, muchas aplicaciones de la vida real admiten ODE o DAE muy rígidos. Intuitivamente, estos problemas experimentan cambios abruptos no suaves a lo largo del tiempo, por lo que se modelan mejor utilizando polinomios de bajo orden distribuidos finamente en pasos cortos, a diferencia de los polinomios de alto orden estirados en pasos largos. Además, la estabilidad a menudo requiere el uso de métodos implícitos. métodos , para los cuales la penalización computacional de los métodos de orden superior es mucho más pronunciada.

Más rigurosamente, los métodos de orden superior son menos estables que los métodos de orden inferior para problemas rígidos. Tenemos, por ejemplo, las barreras de Dahlquist para métodos lineales de varios pasos.

r≤2

Se pueden hacer declaraciones similares (pero mucho más complicadas) para la estabilidad de L en las fórmulas RK. En todos los casos, el aumento del orden a menudo no siempre conduce a soluciones más precisas. Lo siguiente es un extracto del artículo seminal de 1974 de Prothero y Robinson:

Al utilizar métodos de un paso estables en A para resolver grandes sistemas de ecuaciones diferenciales no lineales rígidas, hemos encontrado que
(a) algunos métodos estables en A proporcionan soluciones altamente inestables, y
(b) la precisión de las soluciones obtenidas cuando las ecuaciones son rígido a menudo parece no estar relacionado con el orden del método utilizado.

Para tratamientos aún más rigurosos de este tema, vea el texto clásico de Hairer & Wanner, "Resolviendo ecuaciones diferenciales ordinarias II: Rigidez y diferencia - Problemas algebraicos", 1991.

En la práctica, las ecuaciones rígidas casi siempre se resuelven utilizando la regla trapezoidal o la fórmula TR-BDF2 (funciones ode23t y ode23tb en MATLAB). Ambos son métodos implícitos de segundo orden. Por supuesto, donde la estabilidad no es un problema (es decir, en ecuaciones no rígidas) somos libres de elegir entre varias opciones; RK45 es la opción más común.

Richard Zhang
fuente
Muy interesante. ¿Hay alguna explicación (intuitiva) de por qué el orden tiene que ser menor o igual a 2 para que sea un método multipaso estable en A? Y solo para aclarar, cuando dice que se pueden hacer declaraciones similares para las fórmulas RK, ¿es de orden 2 una vez más?
Mathews24
Pero para los métodos Runge-Kutta, hay métodos A-estables de orden arbitrario.
David Ketcheson
@DavidKetcheson Sí, pero no son muy estables en A (es decir, L-estables). Tienen muchos problemas cuando se usan para resolver DAE, por ejemplo, simular circuitos de transistores simples. De hecho, TR es infame por causar zumbidos artificiales en SPICE, que es lo que motivó el desarrollo de TR-BDF2.
Richard Zhang
@DavidKetcheson Para referencia, consulte doi.org/10.1090/S0025-5718-1974-0331793-2 . La noción de estabilidad A no es lo suficientemente fuerte para DAE, y los métodos estables de alto orden a menudo producen resultados extraños cuando se usan para resolver DAE.
Richard Zhang el
Claro, pero la pregunta no se trata de DAE o de métodos de varios pasos.
David Ketcheson
9

La configuración de referencia

En el software Julia DifferentialEquations.jl implementamos muchos métodos de orden superior, incluidos los métodos Feagin. Puede verlo en nuestra lista de métodos , y luego hay toneladas de otros que puede usar como cuadros suministrados . Debido a que todos estos métodos se combinan, puede compararlos fácilmente. Puede ver los puntos de referencia que tengo en línea aquí , y ver que es muy simple comparar muchos algoritmos diferentes. Entonces, si desea tomarse unos minutos para ejecutar los puntos de referencia, hágalo. Aquí hay un resumen de lo que sale.

En primer lugar, es importante tener en cuenta que, si observa cada uno de los puntos de referencia, verá que nuestra DP5(Orden 5 de Dormand-Prince) y los DP8métodos son más rápidos que los códigos Hairer Fortran ( dopri5y dop853), por lo que estas implementaciones están muy bien optimizadas . Estos muestran que, como se señaló en otro hilo, el uso excesivo de los métodos de Dormand-Prince se debe a que los métodos ya están escritos, no porque sigan siendo los mejores. Entonces, la comparación real entre las implementaciones más optimizadas es entre los métodos Tsitorous, los métodos Verner y los métodos Feagin de DifferentialEquations.jl.

Los resultados

En general, los métodos de un orden superior a 7 tienen un costo computacional adicional que generalmente no se ve superado por el orden, dadas las tolerancias elegidas. Una razón para esto es que las opciones de coeficientes para los métodos de orden inferior están más optimizadas (tienen pequeños "coeficientes de error de truncamiento principales", que son más importantes cuando no es asimétricamente pequeño). Puede ver que en muchos problemas como aquí, los métodos Verner Efficient 6 y 7 funcionan extremadamente bien, pero métodos como el Verner Efficient 8 pueden tener una pendiente más baja. Esto se debe a que las "ganancias" de orden superior se combinan con tolerancias más bajas, por lo que siempre hay una tolerancia donde los métodos de orden superior serán más eficientes.

Sin embargo, la pregunta es entonces, ¿qué tan bajo? En una implementación bien optimizada, eso se vuelve bastante bajo por dos razones. La primera razón es porque los métodos de orden inferior implementan algo llamado FSAL (primero igual que el último). Esta propiedad significa que los métodos de orden inferior reutilizan una evaluación de función del paso anterior en el siguiente paso y, por lo tanto, tienen efectivamente una evaluación de función menos. Si esto se usa correctamente, entonces algo como un método de quinto orden (Tsitorous o Dormand-Prince) en realidad está tomando 5 evaluaciones de funciones en lugar de las 6 que sugieren los cuadros. Esto también es cierto para el método Verner 6.

La otra razón se debe a las interpolaciones. Una razón para usar un método de orden muy alto es tomar menos pasos y simplemente interpolar valores intermedios. Sin embargo, para obtener los valores intermedios, la función de interpolación puede necesitar más evaluaciones de funciones de las que se usaron para dar el paso.Si nos fijamos en los métodos de Verner, se necesitan 8 evaluaciones de funciones adicionales para el método de Orden 8 para obtener un interpolante de Orden 8. Muchas veces, los métodos de bajo orden proporcionan un interpolante "libre", por ejemplo, la mayoría de los métodos de quinto orden tienen una interpolación libre de cuarto orden (sin evaluaciones de funciones adicionales). Esto significa que si necesita valores intermedios (que necesitará para una buena gráfica si está utilizando un método de orden superior), hay algunos costos adicionales ocultos. Tenga en cuenta el hecho de que estos valores interpolados son realmente importantes para el manejo de eventos y la resolución de ecuaciones diferenciales de retardo y verá por qué el costo adicional de interpolación tiene en cuenta.

Entonces, ¿qué pasa con los métodos de Feagin?

Entonces verá que los métodos Feagin faltan sospechosamente en los puntos de referencia. Están bien, las pruebas de convergencia funcionan con números de precisión arbitrarios, etc., pero para que realmente funcionen bien, debes pedir algunas tolerancias bastante absurdamente bajas. Por ejemplo, encontré en puntos de referencia inéditos que el Feagin14rendimiento superior Vern9(el método eficiente Verner de noveno orden) en tolerancias como 1e-30. Para aplicaciones con dinámica caótica (como en los problemas de Astrofísica de 3 cuerpos o Pleides), es posible que desee esta cantidad de precisión debido a la dependencia sensible (los errores en los sistemas caóticos se componen rápidamente). Sin embargo, la mayoría de las personas probablemente están calculando con números de coma flotante de doble precisión, y no he encontrado un punto de referencia donde superen en este dominio de tolerancia.

Además, no hay interpolante para acompañar los métodos de Feagin. Entonces, lo que hago es simplemente ponerles una interpolación de Hermite de tercer orden para que exista uno (y funciona sorprendentemente bien). Sin embargo, si no hay una función de interpolación estándar, puede hacer el método recursivo de Hermite (use esta interpolación para obtener el punto medio, luego haga una interpolación de quinto orden, etc.) para obtener una interpolación de alto orden, pero esto es muy costoso y el resultado es la interpolación no necesariamente tiene un término de error de truncamiento de principio bajo (por lo que solo es bueno cuando dtes realmente pequeño, ¡exactamente lo opuesto al caso que queremos!). Entonces, si alguna vez necesita una interpolación realmente buena para que coincida con su precisión, necesita al menos volver a algo así Vern9.

Nota sobre extrapolación

Tenga en cuenta que los métodos de extrapolación son simplemente algoritmos para generar métodos arbitrarios de Runge-Kutta. Sin embargo, para su orden toman más pasos de los necesarios y tienen altos coeficientes de error de truncamiento de principio, por lo que no son tan eficientes como un método RK bien optimizado en un orden dado. Pero dado el análisis anterior, esto significa que existe un dominio de tolerancia extremadamente baja en el que estos métodos funcionarán mejor que los métodos RK "conocidos". Pero en cada punto de referencia que he ejecutado, parece que no he llegado tan bajo.

Nota sobre estabilidad

La elección realmente no tiene nada que ver con problemas de estabilidad. De hecho, si pasa por los cuadros DifferentialEquations.jl (puede hacerlo solo plot(tab)para las regiones de estabilidad) verá que la mayoría de los métodos tienen regiones de estabilidad sospechosamente similares. Esto es realmente una elección. Por lo general, al derivar los métodos, el autor generalmente hace lo siguiente:

  1. Encuentre los coeficientes de error de truncamiento del principio más bajo (es decir, los coeficientes para los siguientes términos de orden)
  2. Sujeto a las restricciones de orden
  3. Y acerque la región de estabilidad a la del método Dormand-Prince Order 5.

¿Por qué la última condición? Bueno, debido a que ese método tiende a ser siempre estable con la forma en que se realizan las elecciones de tamaño de paso adaptativo controlado por PI, por lo que es una buena barra para regiones de estabilidad "suficientemente buenas". Por lo tanto, no es casualidad que todas las regiones de estabilidad sean similares.

Conclusión

Hay compensaciones en cada elección de método. Los métodos de RK de orden superior simplemente no son tan eficientes con tolerancias más bajas, ya que es más difícil optimizar la elección de los coeficientes y porque el número de compuestos de evaluación de funciones (y crece aún más rápido cuando intervienen las interpolaciones). Sin embargo, si la tolerancia es lo suficientemente baja, ganan, pero las tolerancias que se requieren pueden estar muy por debajo de las aplicaciones "estándar" (es decir, solo se aplican a sistemas caóticos).

Chris Rackauckas
fuente