Integrators

Explicit Euler Method.

\[q_{n+1} = q_{n} + h \, v (t_{n}, q_{n})\]

Explicit Euler integrator cache.

Implicit Euler Method.

Implicit Euler integrator cache.

Symplectic Euler-A Method for separable Hamiltonians.

The momentum is updated explicitly and the position is updated with the new momentum,

\[\begin{aligned} p_{n+1} &= p_{n} + h \, f (t_{n}, q_{n}) , \\ q_{n+1} &= q_{n} + h \, v (t_{n+1}, p_{n+1}) . \end{aligned}\]

This is the adjoint of SymplecticEulerB. See SymplecticEulerMethod for the separability assumption.

Symplectic Euler-B Method for separable Hamiltonians.

The position is updated explicitly and the momentum is updated with the new position,

\[\begin{aligned} q_{n+1} &= q_{n} + h \, v (t_{n}, p_{n}) , \\ p_{n+1} &= p_{n} + h \, f (t_{n+1}, q_{n+1}) . \end{aligned}\]

This is the adjoint of SymplecticEulerA. See SymplecticEulerMethod for the separability assumption.

Symplectic Euler integrator cache.

Abstract supertype of the symplectic Euler methods SymplecticEulerA and SymplecticEulerB.

Both methods are implemented for a separable Hamiltonian, that is

\[H (t, q, p) = T (t, p) + V (t, q) ,\]

so that the vector fields satisfy $v = v(t,p)$ and $f = f(t,q)$. Under this assumption the otherwise implicit partitioned scheme decouples into two explicit substeps and no nonlinear solver is required. Separability cannot be checked at runtime, so applying these methods to a non-separable Hamiltonian silently computes something that is neither symplectic nor consistent with the symplectic Euler method.

Implicit Midpoint Method.

For an ordinary differential equation $\dot{q} = v(t,q)$, that is an ODEProblem, the method reads

\[q_{n+1} = q_{n} + h \, v \bigg( t_{n} + \frac{h}{2} , \frac{q_{n} + q_{n+1}}{2} \bigg)\]

The nonlinear solver solves for the stage vector field $V = v(t_{n} + h/2, q_{n} + h V / 2)$, so that the update reads $q_{n+1} = q_{n} + h V$.

For a partitioned differential equation, that is a PODEProblem or HODEProblem,

\[\dot{q} = v (t, q, p) , \qquad \dot{p} = f (t, q, p) ,\]

the same quadrature is applied to both components. With the stage time $\tilde{t} = t_{n} + h/2$ and the midpoints $Q = q_{n} + h V / 2$ and $P = p_{n} + h F / 2$, the nonlinear solver solves

\[V = v (\tilde{t}, Q, P) , \qquad F = f (\tilde{t}, Q, P) ,\]

for the stage vector fields $V$ and $F$, and the updates read

\[q_{n+1} = q_{n} + h \, V , \qquad p_{n+1} = p_{n} + h \, F .\]

Both components are coupled through $v$ and $f$, so the solver solution vector holds both stage vector fields and the nonlinear system is twice the size of the one for an ordinary differential equation. Even for a separable Hamiltonian the method does not decouple into two explicit substeps the way the symplectic Euler methods do: it is the Gauss method with a single stage, applied to the partitioned system.

For an implicit differential equation, that is an IODEProblem or LODEProblem,

\[\begin{aligned} p &= \vartheta (t, q, v) , & \dot{p} &= f (t, q, v) , & \dot{q} &= v , \end{aligned}\]

the same quadrature is applied to the momentum map and to the force, which amounts to the Gauss method with a single stage. With the stage time $\tilde{t} = t_{n} + h/2$ and the midpoint $Q = q_{n} + h V / 2$, the nonlinear solver solves

\[\vartheta (\tilde{t}, Q, V) = p_{n} + \frac{h}{2} \, f (\tilde{t}, Q, V)\]

for the stage velocity $V$, and the updates read

\[\begin{aligned} q_{n+1} &= q_{n} + h \, V , & p_{n+1} &= p_{n} + h \, f (\tilde{t}, Q, V) . \end{aligned}\]

Whenever $\vartheta$ is regular, so that the implicit equation is equivalent to a partitioned ordinary differential equation, this is the same map as the implicit midpoint method applied to that equation. The solver solution vector holds $V$ in both cases, so the nonlinear system has the same size as for an ordinary differential equation.

Implicit midpoint integrator cache.

Fields

• x: nonlinear solver solution vector, holding the stage vector field $V$
• q: midpoint of the time step, $Q = q_{n} + h V / 2$
• v: stage vector field $V$

Implicit midpoint integrator cache for implicit differential equations.

Fields

• x: nonlinear solver solution vector, holding the stage velocity $V$
• q: midpoint of the time step, $Q = q_{n} + h V / 2$
• v: stage velocity $V$
• θ: momentum map $\vartheta$ evaluated at the stage
• f: force $f$ evaluated at the stage

Implicit midpoint integrator cache for partitioned differential equations.

Fields

• x: nonlinear solver solution vector, holding the stage vector fields $V$ and $F$
• q: midpoint of the time step, $Q = q_{n} + h V / 2$
• p: midpoint of the time step, $P = p_{n} + h F / 2$
• v: stage vector field $V$
• f: stage vector field $F$

Crank-Nicolson Method, also known as the trapezoidal rule.

For an ordinary differential equation $\dot{q} = v(t,q)$, that is an ODEProblem, the method reads

\[q_{n+1} = q_{n} + \frac{h}{2} \, \big[ v (t_{n}, q_{n}) + v (t_{n+1}, q_{n+1}) \big]\]

The nonlinear solver solves for the vector field at the new time step, $V = v(t_{n+1}, q_{n+1})$, while $\bar{V} = v(t_{n}, q_{n})$ is computed once per time step. The method is symmetric and second order, but not symplectic (it is conjugate to a symplectic method).

For a partitioned differential equation, that is a PODEProblem or HODEProblem,

\[\dot{q} = v (t, q, p) , \qquad \dot{p} = f (t, q, p) ,\]

the trapezoidal rule is applied to both components,

\[\begin{aligned} q_{n+1} &= q_{n} + \frac{h}{2} \, ( \bar{V} + V ) , & p_{n+1} &= p_{n} + \frac{h}{2} \, ( \bar{F} + F ) , \end{aligned}\]

where $\bar{V} = v(t_{n}, q_{n}, p_{n})$ and $\bar{F} = f(t_{n}, q_{n}, p_{n})$ are computed once per time step, as they are given by the solution at the beginning of the time step, while $V = v(t_{n+1}, q_{n+1}, p_{n+1})$ and $F = f(t_{n+1}, q_{n+1}, p_{n+1})$ are solved for. The solver solution vector holds both of them, so the nonlinear system is twice the size of the one for an ordinary differential equation.

For an implicit differential equation, that is an IODEProblem or LODEProblem,

\[\begin{aligned} p &= \vartheta (t, q, v) , & \dot{p} &= f (t, q, v) , & \dot{q} &= v , \end{aligned}\]

the trapezoidal rule is applied to the position and to the momentum alike, which amounts to the Lobatto IIIA method with two stages,

\[\begin{aligned} q_{n+1} &= q_{n} + \frac{h}{2} \, ( \bar{V} + V ) , & p_{n+1} &= p_{n} + \frac{h}{2} \, \big[ f (t_{n}, q_{n}, \bar{V}) + f (t_{n+1}, q_{n+1}, V) \big] . \end{aligned}\]

In contrast to the explicit case, the velocity at the beginning of the time step is not given by a function evaluation but implicitly by $\vartheta (t_{n}, q_{n}, \bar{V}) = p_{n}$. Both velocities are therefore solved for simultaneously, so that the nonlinear system reads

\[\begin{aligned} 0 &= \vartheta (t_{n}, q_{n}, \bar{V}) - p_{n} , \\ 0 &= \vartheta (t_{n+1}, q_{n+1}, V) - p_{n} - \frac{h}{2} \, \big[ f (t_{n}, q_{n}, \bar{V}) + f (t_{n+1}, q_{n+1}, V) \big] , \end{aligned}\]

and is twice the size of the one for an ordinary differential equation. Whenever $\vartheta$ is regular, so that the implicit equation is equivalent to a partitioned ordinary differential equation, this is the same map as the Crank-Nicolson method applied to that equation.

Crank-Nicolson integrator cache.

Fields

• x: nonlinear solver solution vector, holding the vector field $V = v(t_{n+1}, q_{n+1})$
• q: solution at the end of the time step
• v: vector field at the end of the time step
• v̄: vector field at the beginning of the time step, $v(t_{n}, q_{n})$, which is constant throughout the nonlinear solve

Crank-Nicolson integrator cache for implicit differential equations.

In contrast to CrankNicolsonCache, nothing is constant during the solve: the velocity at the beginning of the time step is part of the nonlinear solver solution vector, so every field is computed from it and has to be accessed through cache(int, ST).

Fields

• x: nonlinear solver solution vector, holding $\bar{V}$ and $V$
• q: solution at the end of the time step
• v, θ, f: velocity, momentum map and force at the end of the time step
• v̄, θ̄, f̄: velocity, momentum map and force at the beginning of the time step

Crank-Nicolson integrator cache for partitioned differential equations.

As in CrankNicolsonCache, the vector fields at the beginning of the time step are constant during the solve, so they are read from the cache at working precision.

Fields

• x: nonlinear solver solution vector, holding $V$ and $F$
• q, p: solution at the end of the time step
• v, f: vector fields at the end of the time step
• v̄, f̄: vector fields at the beginning of the time step, $v(t_{n}, q_{n}, p_{n})$ and $f(t_{n}, q_{n}, p_{n})$, which are constant throughout the nonlinear solve

Compute the stages of the Crank-Nicolson method for a partitioned differential equation from the nonlinear solver solution x.

Requires v̄ and f̄ in the cache at working precision to hold $v(t_{n}, q_{n}, p_{n})$ and $f(t_{n}, q_{n}, p_{n})$, which integrate_step! computes at the beginning of every time step. Calling this function before that would silently use stale or zero values.

Compute the stages of the Crank-Nicolson method from the nonlinear solver solution x.

Requires v̄ in the cache at working precision to hold $v(t_{n}, q_{n})$, which integrate_step! computes at the beginning of every time step. Calling this function before that would silently use a stale or zero v̄.