Simplemente quiero saber si el método numérico Dormand-Prince o el método numérico Cash-Karp es más preciso.

Como acabo de terminar de optimizar muchos de ellos en un software, DifferentialEquations.jl , decidí presentar una comparación de los métodos principales de Order 4/5. El método Fehlberg se omitió porque comúnmente se sabe que es menos eficiente que el método DP5.

Historias de fondo

Dormand-Prince 4/5

El método Dormand-Prince fue desarrollado para ser preciso como un par 4/5 con uso de extrapolación local (es decir, paso con el par de orden 5. Esto se debe a que fue diseñado para estar cerca del coeficiente de error de truncamiento del principio óptimo (es decir, mínimo) (bajo la restricción de tener también el número mínimo de pasos para lograr el orden 5). Tiene una interpolación de orden 4 que es gratuita, pero necesita pasos adicionales para una interpolación de orden 5.

Cash-Karp 4/5

El método Cash-Karp se desarrolló para satisfacer diferentes restricciones, es decir, para tratar mejor los problemas no uniformes. Eligieron que , el porcentaje del paso de tiempo en el ésimo paso (es decir, es el tiempo en el que se calcula el ésimo paso) para que sea lo más uniforme posible, y aún así alcanzar el orden 5. Entonces también se derivó para incorporar métodos de primer, segundo, tercer y cuarto orden con esta uniformidad de$c_i$$i$$t+c_i \Delta t$$i$$c_i$. Están espaciados de tal manera que puede averiguar dónde comienza una parte rígida por la cual la diferencia es grande. Además, tenga en cuenta que cuanto más rígida es la ecuación, peor es un método de orden superior (porque necesita límites en derivados más altos). Entonces, desarrollan una estrategia que utiliza los 5 métodos integrados para "dejar de fumar temprano": es decir, si detecta rigidez, deténgase en la etapa $i<6$para disminuir el número de llamadas a funciones y ahorrar tiempo. Entonces, al final, este "par" se desarrolló con muchas otras restricciones en mente, por lo que no hay razón para esperar que sea "más preciso", al menos como un par 4/5. Si agrega toda esta otra maquinaria, entonces, en problemas (semi) rígidos, será más preciso (pero en ese caso es posible que desee utilizar un método diferente como el método W-Rosenbrock). Esta es una razón por la que este par no se ha convertido en estándar sobre el par DP5, pero aún puede ser útil (¿tal vez sería bueno para un método híbrido que cambia a un solucionador rígido cuando se encuentra rigidez?).

Bogacki y Shampine 4/5

Para completar la respuesta, analicemos el par Bogacki y Shampine que se mencionó en el comentario. El método BS5 elimina la restricción de "usar la menor cantidad de llamadas a funciones" (usa 8 en lugar de 6) para hacer 2 cosas:

• Obtenga coeficientes de error de truncamiento de principio realmente bajos.

• Produzca una interpolación de orden 5 con coeficientes de error más bajos.

Estos coeficientes son tan bajos que para muchos problemas con tolerancias que los usuarios probablemente usan, mide como si fuera el sexto orden. Su artículo muestra que para llamadas de función baratas, esto puede ser más eficiente que DP5, ya que aproximadamente la misma cantidad en DP5 fue sobre RKF5 (el método Fehlberg).

Puede poner dos y dos juntos y ver: espere un segundo, Shampine es la misma persona que desarrolló la suite MATLAB ODE, esto fue después de que se publicó el documento del par BS5, ¿por qué MATLAB no ode45usa el par BS5? Una razón es que se realizó principalmente antes de que el par BS5 se retransmitiera. La otra razón es porque la ode45función se desarrolló para minimizar el tiempo. Si bien el par BS5 es más eficiente (es decir, obtiene una precisión menor), el propósito deode45es tener un error lo suficientemente bueno como para hacer una trama lo suficientemente buena. Esto significa que, para lidiar con los pasos grandes, también produce dos soluciones interpoladas adicionales entre cada paso. Para el método DP5, hay una interpolación de orden 4 "libre", por lo que esto es mucho más rápido que usar BS5. Dado que también es "lo suficientemente preciso" con tolerancias moderadas, este método se establece como estándar porque brinda una mejor experiencia de usuario estándar que BS5 al hacer computación interactiva (por lo que esta elección fue específica del contexto).

Tsitorous 4/5

Aquí hay uno menos que la gente conoce. Se deriva en este documento . Se deriva utilizando menos suposiciones que el método DP5 e intenta obtener un par con coeficientes de error de truncamiento de principio más bajos. En sus pruebas, afirma que logra esto. También tiene una interpolación de orden libre 4 como el método DP5.

Pruebas numéricas

Escribí el paquete numérico DifferentialEquations.jl para ser un conjunto bastante completo de solucionadores para Julia. A lo largo del camino de guerra, implementé más de 100 métodos de Runge-Kutta y optimicé mucho a mano. Tres de los integradores optimizados a mano son los métodos DP5, BS5 y Tsit5 (no hice CK5 porque, como se señaló en la historia de fondo, su caso principal es para problemas que son un poco rígidos. Creo que la mejor manera de manejarlos es usar DP5 / BS5 y cambiar a solucionadores rígidos según sea necesario de una manera como LSODE, pero esa es una historia para un momento diferente) (una forma de ver que están cerca de lo óptimo es que estos métodos son más rápidos que las dopri5implementaciones de Hairer , por lo que son al menos implementaciones decentes). Pruebas entre muchos métodos de Runge-Kutta en ecuaciones no rígidasse puede encontrar en la carpeta de puntos de referencia . A medida que avanzo, agrego más, pero como pueden ver en los diagramas de precisión de trabajo del ODE lineal y del problema de tres cuerpos, mido los métodos DP5 y Tsit5 para tener una eficiencia casi idéntica, superando el método BS5 en el ODE lineal, mientras que DP5 y BS5 son casi idénticos en el problema de tres cuerpos con Tsit5 detrás. A partir de esta información, al menos por ahora, me he decidido por el método DP5 como predeterminado, que coincide con las recomendaciones anteriores. Eso puede cambiar con futuras pruebas (¡o podría agregar puntos de referencia! Siéntase libre de contribuir o marque el repositorio para brindar más apoyo a este esfuerzo).

Conclusión

En conclusión, los pares de la Orden 5 son así:

• El par Dormand-Prince 4/5 es un buen par, ya que está bien optimizado en términos de coeficiente de error de truncamiento principal y tiene una interpolación de orden 4 barata, lo que lo hace rápido para producir parcelas decentes.

• El par Cash-Karp tiene más restricciones para manejar mejor las ecuaciones rígidas. Sin embargo, para obtener el beneficio completo, querrá usar el algoritmo completo con los 5 métodos integrados.

• El método Bogacki & Shampine Order 5 puede ser el más eficiente en términos de producción de llamadas de error por función (tiene un estimador de doble error, por lo que en problemas más difíciles probablemente sea mejor), pero eso le permite tomar pasos de tiempo más largos. Sin embargo, si solo desea producir una trama uniforme, debe contrarrestar este método: use una tolerancia más baja (por lo que tomará más tiempo que DP5 pero con menos error) o use pasos más interpolados. Al final, esto significaba que podría no ser mejor para aplicaciones interactivas, aunque podría ser mejor para algunas aplicaciones informáticas científicas.

• El Tsitorous 4/5. Fue desarrollado bastante recientemente (2011) para vencer al DP5 en una comparación directa. Mis pruebas no me dan una razón para creer que es mucho mejor que DP5 que ahora debería considerarse como el nuevo método estándar, pero las pruebas futuras pueden comenzar a ponerse a su favor.

Editar

Mejoré la implementación de Tsit5. Ahora funciona mejor que DP5 en la mayoría de las pruebas, tanto las implementaciones DifferentialEquations.jl como Hairer dopri (aunque uno podría sorprenderse de que las implementaciones DifferentialEquations.jl sean realmente más rápidas, lo que, por supuesto, ayuda a la implementación Tsit5). Ahora lo recomiendo como el método predeterminado de orden 4/5.

Si está interesado en comparar dos integradores, resuelva

$(q,p) = (1,0)$$dt = 0.1$

$t$

• PJ Channell y C Scovel. Integración simpléctica de sistemas hamiltonianos . Journal of Computational Physics 73, 468 (1973).