¿Qué significa "simpléctico" en referencia a los integradores numéricos, y los utiliza la odeint de SciPy?

25

En este comentario escribí:

... integrador predeterminado de SciPy, que supongo que solo usa métodos simplécticos.

en el que me refiero a SciPy's odeint, que utiliza un "método no rígido (Adams)" o un "método rígido (BDF)". Según la fuente :

def odeint(func, y0, t, args=(), Dfun=None, col_deriv=0, full_output=0,
           ml=None, mu=None, rtol=None, atol=None, tcrit=None, h0=0.0,
           hmax=0.0, hmin=0.0, ixpr=0, mxstep=0, mxhnil=0, mxordn=12,
           mxords=5, printmessg=0):
    """
    Integrate a system of ordinary differential equations.

    Solve a system of ordinary differential equations using lsoda from the
    FORTRAN library odepack.

    Solves the initial value problem for stiff or non-stiff systems
    of first order ode-s::
        dy/dt = func(y, t0, ...)
    where y can be a vector.
    """

Aquí hay un ejemplo en el que propago la órbita de un satélite alrededor de la Tierra durante tres meses solo para demostrar que funciona como se esperaba.

Creo que los integradores no simplécticos tienen la propiedad indeseable de que tenderán a no conservar energía (u otras cantidades) y, por lo tanto, no son deseables en la mecánica orbital, por ejemplo. Pero no estoy exactamente seguro de qué es lo que hace que un integrador simpléctico sea simpléctico.

¿Es posible explicar cuál es la propiedad (que hace que un integrador simpléctico sea simpléctico) de una manera directa y (bastante) fácil de entender pero no inexacta? Pregunto desde el punto de vista de cómo funciona el integrador internamente , en lugar de cómo funciona en las pruebas.

¿Y es correcta mi sospecha de que odeintsolo usa integradores simplécticos?

UH oh
fuente
44
Como regla general, solo debe esperar que un integrador de caja negra sea simpléctico si requiere que separe las ecuaciones de posición y momento.
origimbo
@origimbo Gracias. Estos sí, y parece que odeintes una envoltura de Python para códigos fuente bastante antiguos, establecidos y bien referenciados (pregunta editada, referencias ODEPACK y LSODA) aunque ciertamente admito que lo usé en modo de caja negra. Mi ejemplo vinculado muestra que el vector de estado 6D consta de tres posiciones y tres velocidades.
uhoh
11
Los integradores ODE en ODEPACK y LSODA no son integradores simplécticos.
Brian Borchers
2
Aquí hay un ejemplo trabajado que compara dos solucionadores muy simples: Euler y Symplectic Euler: idontgetoutmuch.wordpress.com/2013/08/06/… .
idontgetoutmuch
2
El libro de Hairer, Nørsett y Wanner ofrece una buena explicación de los métodos simplécticos. Mire la Figura 16.1 en particular, y las figuras aquí .
JM

Respuestas:

47

Déjame comenzar con correcciones. No, odeintno tiene integradores simplécticos. No, la integración simpléctica no significa conservación de energía.

¿Qué significa symplectic y cuándo debe usarlo?

En primer lugar, ¿qué significa symplectic? Simpléctico significa que la solución existe en una variedad simpléctica. Una variedad simpléctica es un conjunto de soluciones que se define mediante una forma 2. Los detalles de las variedades simplécticas probablemente suenan como tonterías matemáticas, por lo que la esencia de esto es que hay una relación directa entre dos conjuntos de variables en tal variedad. La razón por la cual esto es importante para la física es porque las ecuaciones de Hamilton naturalmente tienen que las soluciones residen en una variedad simpléctica en el espacio de fase, siendo la división natural los componentes de posición y momento. Para la verdadera solución hamiltoniana, ese camino de fase espacial es energía constante.

O(Δtn)n

Lo que esto significa es que los integradores simplécticos tienden a capturar los patrones a largo plazo mejor que los integradores normales debido a esta falta de deriva y esto casi garantiza la periodicidad. Este cuaderno muestra bien esas propiedades en el problema de Kepler . La primera imagen muestra de lo que estoy hablando con la naturaleza periódica de la solución.

ingrese la descripción de la imagen aquí

Esto se resolvió utilizando el integrador simpléctico de sexto orden de Kahan y Li de DifferentialEquations.jl . Puede ver que la energía no está exactamente conservada, pero su variación depende de qué tan lejos esté el colector de solución perturbado del colector verdadero. Pero dado que la solución numérica en sí misma reside en una variedad simpléctica, tiende a ser casi exactamente periódica (con alguna deriva numérica lineal que puede ver), lo que lo hace muy bien para la integración a largo plazo. Si hace lo mismo con RK4, puede sufrir un desastre:

ingrese la descripción de la imagen aquí

Puede ver que el problema es que no hay una verdadera periodicidad en la solución numérica y, por lo tanto, las horas extraordinarias tienden a derivar.

Esto resalta la verdadera razón para elegir integradores simplécticos: los integradores simplécticos son buenos en integraciones a largo plazo en problemas que tienen la propiedad simpléctica (sistemas hamiltonianos) . Así que veamos algunas cosas. Tenga en cuenta que no siempre necesita integradores simplécticos incluso en un problema simpléctico. Para este caso, un método adaptativo de Runge-Kutta de quinto orden puede funcionar bien. Aquí está Tsit5:

ingrese la descripción de la imagen aquí

Noten dos cosas. Uno, obtiene una precisión lo suficientemente buena como para que no se pueda ver la deriva real en el diagrama del espacio de fases. Sin embargo, en el lado derecho puede ver que existe esta deriva de energía, por lo que si realiza una integración lo suficientemente larga, este método no funcionará tan bien como el método de solución con las propiedades periódicas. Pero eso plantea la pregunta, ¿cómo le va con respecto a la eficiencia frente a la simple integración con extrema precisión? Bueno, esto es un poco menos seguro. En DiffEqBenchmarks.jl puede encontrar algunos puntos de referencia que investigan esta pregunta. Por ejemplo, este cuadernoexamina el error de energía frente al tiempo de ejecución en un sistema de ecuaciones de Hamilton de un modelo Boson cuádruple y muestra que si desea una precisión realmente alta, incluso para tiempos de integración bastante largos, es más eficiente usar un RK de alto orden o Runge-Kutta Nystrom ( RKN) método. Esto tiene sentido porque para satisfacer la propiedad simpléctica, los integradores renuncian a cierta eficiencia y tienen que ser un paso de tiempo fijo (hay algunas investigaciones que avanzan en este último pero no está muy avanzado).

Además, observe en estos dos cuadernos que también puede tomar un método estándar y proyectarlo nuevamente en la solución múltiple en cada paso (o cada pocos pasos). Esto es lo que están haciendo los ejemplos que utilizan la devolución de llamada de DifferentialEquations.jl ManifoldProjection . Verá que las garantías de las leyes de conservación se mantienen pero con un costo adicional de resolver un sistema implícito en cada paso. También puede usar un solucionador de ODE totalmente implícito o matrices de masas singulares para agregar ecuaciones de conservación, pero el resultado final es que estos métodos son más costosos computacionalmente como una compensación.

Para resumir, la clase de problemas en los que desea alcanzar un integrador simpléctico son aquellos que tienen una solución en una variedad simpléctica (sistemas hamiltonianos) en los que no desea invertir los recursos computacionales para tener una tolerancia muy exacta (tolerancia <1e-12) solución y no necesita energía exacta / etc. conservación. Esto pone de relieve que se trata de propiedades de integración a largo plazo, por lo que no debería acudir a todas ellas como lo sugiere la literatura. Pero siguen siendo una herramienta muy importante en muchos campos, como la astrofísica, donde tienes integraciones a largo plazo que debes resolver lo suficientemente rápido sin tener una precisión absurda.

¿Dónde encuentro integradores simplécticos? ¿Qué tipo de integradores simplécticos existen?

Generalmente hay dos clases de integradores simplécticos. Existen los integradores simplécticos de Runge-Kutta (que son los que se muestran en los ejemplos anteriores) y hay métodos implícitos de Runge-Kutta que tienen la propiedad simpléctica. Como @origimbo menciona, los integradores simplécticos de Runge-Kutta requieren que les proporciones una estructura particionada para que puedan manejar la posición y las partes de momento por separado. Sin embargo, en contra del comentario, los métodos implícitos de Runge-Kutta son simplécticos sin requerir esto, sino que requieren resolver un sistema no lineal. Esto no es tan malo porque si el sistema no es rígido, este sistema no lineal puede resolverse con iteración funcional o aceleración de Anderson, pero los métodos RK simplécticos probablemente deberían preferirse por su eficiencia (

Dicho esto, odeint no tiene métodos de ninguna de estas familias , por lo que no es una buena opción si está buscando integradores simplécticos. En Fortran, el sitio de Hairer tiene un pequeño conjunto que puedes usar . Mathematica tiene algunos incorporados . Los solucionadores GSL ODE tienen integradores de puntos gaussianos RK implícitos que IIRC son simplécticos, pero esa es la única razón para usar los métodos GSL.

Pero el conjunto más completo de integradores simplécticos se puede encontrar en DifferentialEquations.jl en Julia (recuerde que esto se usó para los cuadernos anteriores). La lista de métodos simplécticos disponibles de Runge-Kutta se encuentra en esta página y notará que el método de punto medio implícito también es simpléctico (el método trapezoidal implícito de Runge-Kutta se considera "casi simpléctico" porque es reversible). No solo tiene el conjunto de métodos más grande, sino que también es de código abierto (puede ver el código y sus pruebas en un lenguaje de alto nivel) y tiene muchos puntos de referencia . Un buen cuaderno introductorio para usarlo para resolver problemas físicos es este cuaderno tutorial. Pero, por supuesto, se recomienda comenzar con el paquete a través del primer tutorial de ODE .

En general, puede encontrar un análisis detallado de las suites de ecuaciones diferenciales numéricas en esta publicación de blog . Es bastante detallado, pero dado que tiene que abarcar muchos temas, cada uno lo hace con menos detalles que este, así que no dude en solicitar que se amplíe de alguna manera.

Chris Rackauckas
fuente
10
¡Con esta respuesta, parece que he ganado el Jackpot Stack Exchange! Esta es la respuesta perfecta para mí, ya que entiendo algo de inmediato y las partes que no me dejan ansioso por leer más. Realmente aprecio el tiempo que ha tomado para obtener esta respuesta exhaustivamente, así como para incluir otros enlaces útiles e instructivos.
uhoh
Antes de entrar en detalles matemáticos, podríamos decir que simpléctico significa conservación del volumen , ¿no?
Miguel
2
FTR, la razón por la cual el Runge-Kutta de 5to orden adaptativo se desempeña mucho mejor que RK4 aquí no es que tenga un orden más alto sino que elija tamaños de pasos más adecuados. La razón por la que RK4 funciona tan mal es principalmente porque el tamaño del paso es inapropiadamente alto en el perigeo; el mismo solucionador con la mitad del tamaño del paso daría una solución mucho mejor. (Solo, perdería mucho tiempo resolviendo la órbita finamente alrededor del apogeo, donde esto no es necesario en absoluto.)
Leftaroundabout
1
Excelente exposición. Como pregunta adicional: el OP comienza con una referencia a Python: ¿hay tutoriales / paquetes de Python recomendados a lo largo de las líneas de los ejemplos de Julia vinculados?
Quetzalcoatl
1
El único paquete de Python que conozco para este tipo de integradores es diffeqpy , donde no está documentado en el archivo README, pero puede acceder a todos estos mismos métodos y volver a escribir esto en Python usando ese paquete.
Chris Rackauckas
14

pqH(p,q)

dqdt=+Hp
dpdt=Hq.
H
p(t),q(t)=ϕt(p(t0),q(t0))
dpdqpqson unidimensionales, se puede pensar que esto dice que el área dentro de las curvas cerradas en el espacio de fase se conserva. Esto asegura todo tipo de buenas propiedades de estabilidad, ya que las "bolas" de las trayectorias deben mantenerse "cercanas" entre sí.

En términos numéricos, un integrador simpléctico actúa de la misma manera, conservando también esta área / dos formas. A su vez, esto significa que hay un "hamiltoniano numérico" conservado (que puede no ser [leer 'no es'] el mismo que el exacto). Tenga en cuenta que la estabilidad no es lo mismo que la precisión, por lo que la mayoría de las ventajas de los métodos simplécticos se obtienen cuando se integra durante mucho tiempo (por ejemplo, su método puede colocar rápidamente un satélite en el lado equivocado de la Tierra, sin permitir que se descomponga en eso).

origimbo
fuente
¡Gracias por esto! Ahora usaré palabras por encima de mi calificación salarial. Las n-bolas de trayectorias están más en riesgo cuando están cerca de bifurcaciones, como las de las simulaciones de 3 cuerpos. cf. Doedel y col. 2007, int. J. Bifurcación y Caos, v 17, no. 8 (2007) 2625–2677 ¿Cómo lo hice? También ieec.cat/hosted/web-libpoint/papers/…
uhoh el
2
A menos que el lector conozca los detalles matemáticos, la mención de la estabilidad es engañosa, ya que la conservación del volumen no significa que las trayectorias individuales permanezcan cercanas.
Miguel
1
@Miguel Creo que esta es una de esas situaciones en las que el lector que no sigue los detalles matemáticos está condenado de cualquier manera, pero en términos de la troika habitual de numeración de precisión, estabilidad y eficiencia computacional, diría que enfatizar la estabilidad Los beneficios son útiles. Me complace aceptar sugerencias para una reescritura si puede pensar en algo mejor que no sea deliberadamente inexacto.
origimbo
22
1
@Miguel: Pero la masa de partículas se puede dividir en dos o más partes. Su volumen total solo tiene que mantenerse constante.
Wolfgang Bangerth