Symplectic integrator

In mathematics, a symplectic integrator (SI) is a numerical integration scheme for Hamiltonian systems. Symplectic integrators form the subclass of geometric integrators which, by definition, are canonical transformations. They are widely used in nonlinear dynamics, molecular dynamics, discrete element methods, accelerator physics, plasma physics, quantum physics, and celestial mechanics.
Introduction.
Symplectic integrators are designed for the numerical solution of Hamilton's equations, which read
where formula_2 denotes the position coordinates, formula_3 the momentum coordinates, and formula_4 is the Hamiltonian. 
The set of position and momentum coordinates formula_5 are called canonical coordinates.
The time evolution of Hamilton's equations is a symplectomorphism, meaning that it conserves the symplectic 2-form formula_6. A numerical scheme is a symplectic integrator if it also conserves this 2-form.
Symplectic integrators possess, as a conserved quantity, a Hamiltonian which is slightly perturbed from the original one. By virtue of these advantages, the SI scheme has been widely applied to the calculations of long-term evolution of chaotic Hamiltonian systems ranging from the Kepler problem to the classical and semi-classical simulations in molecular dynamics.
Most of the usual numerical methods, like the primitive Euler scheme and the classical Runge–Kutta scheme, are not symplectic integrators.
Methods for constructing symplectic algorithms.
Splitting methods for separable Hamiltonians.
A widely used class of symplectic integrators is formed by the splitting methods.
Assume that the Hamiltonian is separable, meaning that it can be written in the form
This happens frequently in Hamiltonian mechanics, with "T" being the kinetic energy and "V" the potential energy.
For the notational simplicity, let us introduce the symbol formula_8 to denote the canonical coordinates 
including both the position and momentum coordinates. Then, the set of the Hamilton's equations given in the introduction can be expressed in a single expression as
where formula_10 is a Poisson bracket. Furthermore, by introducing an operator formula_11, which returns a Poisson bracket of the operand with the Hamiltonian, the expression of the Hamilton's equation can be further simplified to
The formal solution of this set of equations is given as a matrix exponential:
Note the positivity of formula_14 in the matrix exponential.
When the Hamiltonian has the form of eq. (1), the solution (3) is equivalent to
The SI scheme approximates the time-evolution operator formula_16 in the formal solution (4) by a product of operators as
where formula_18 and formula_19 are real numbers, formula_20 is an integer, which is called the order of the integrator, and where formula_21. Note that each of the operators formula_22 and formula_23 provides a symplectic map, so their product appearing in the right-hand side of (5) also constitutes a symplectic map.
Since formula_24 for all formula_25, we can conclude that
By using a Taylor series, formula_27 can be expressed as
where formula_29 is an arbitrary real number. Combining (6) and (7), and by using the same reasoning for formula_30 as we have used for formula_31, we get
In concrete terms, formula_22 gives the mapping
and formula_23 gives
Note that both of these maps are practically computable.
Examples.
The simplified form of the equations (in executed order) are:
After converting into Lagrangian coordinates:
Where formula_41 is the force vector at formula_42, formula_43 is the acceleration vector at formula_42, and formula_45 is the scalar quantity of mass.
Several symplectic integrators are given below. An illustrative way to use them is to consider 
a particle with position formula_2 and velocity formula_3.
To apply a timestep with values formula_48 to the particle, carry out the following steps:
Iteratively:
A first-order example.
The symplectic Euler method is the first-order integrator with formula_54 and coefficients
Note that the algorithm above does not work if time-reversibility is needed.
The algorithm has to be implemented in two parts, one for positive time steps, one for negative time steps.
A second-order example.
The Verlet method is the second-order integrator with formula_56 and coefficients
Since formula_58, the algorithm above is symmetric in time. There are 3 steps to the algorithm, and step 1 and 3 are exactly the same, so the positive time version can be used for negative time.
A third-order example.
A third-order symplectic integrator (with formula_59) was discovered by Ronald Ruth in 1983.
One of the many solutions is given by
A fourth-order example.
A fourth-order integrator (with formula_61) was also discovered by Ruth in 1983 and distributed privately to the particle-accelerator community at that time. This was described in a lively review article by Forest.
This fourth-order integrator was published in 1990 by Forest and Ruth and also 
independently discovered by two other groups around that same time.
To determine these coefficients, the Baker–Campbell–Hausdorff formula can be used. Yoshida, in particular, gives an elegant derivation of coefficients for higher-order integrators. Later on, Blanes and Moan further developed partitioned Runge–Kutta methods for the integration of systems with separable Hamiltonians with very small error constants.
Splitting methods for general nonseparable Hamiltonians.
General nonseparable Hamiltonians can also be explicitly and symplectically integrated.
To do so, Tao introduced a restraint that binds two copies of phase space together to enable an explicit splitting of such systems.
The idea is, instead of formula_63, one simulates formula_64, whose solution agrees with that of formula_63 in the sense that formula_66.
The new Hamiltonian is advantageous for explicit symplectic integration, because it can be split into the sum of three sub-Hamiltonians, formula_67, formula_68, and formula_69. Exact solutions of all three sub-Hamiltonians can be explicitly obtained: both formula_70 solutions correspond to shifts of mismatched position and momentum, and formula_71 corresponds to a linear transformation. To symplectically simulate the system, one simply composes these solution maps.
Applications.
In plasma physics.
In recent decades symplectic integrator in plasma physics has become an active research topic, because straightforward applications of the standard symplectic methods do not suit the need of large-scale plasma simulations enabled by the peta- to exa-scale computing hardware. Special symplectic algorithms need to be customarily designed, tapping into the special structures of physics problem under investigation. One such example is the charged particle dynamics in an electromagnetic field. With the canonical symplectic structure, the Hamiltonian of the dynamics is formula_72 whose formula_73-dependence and formula_74-dependence are not separable, and standard explicit symplectic methods do not apply. For large-scale simulations on massively parallel clusters, however, explicit methods are preferred. 
To overcome this difficulty, we can explore the specific way that the formula_73-dependence and formula_74-dependence are entangled in this Hamiltonian, and try to design a symplectic algorithm just for this or this type of problem. First, we note that the formula_73-dependence is quadratic, therefore the first order symplectic Euler method implicit in formula_73 is actually explicit. This is what is used in the canonical symplectic particle-in-cell (PIC) algorithm. To build high order explicit methods, we further note that the formula_73-dependence and formula_74-dependence in this formula_81 are product-separable, 2nd and 3rd order explicit symplectic algorithms can be constructed using generating functions, and arbitrarily high-order explicit symplectic integrators for time-dependent eletromagnetic fields can also be constructed using Runge-Kutta techniques.
A more elegant and versatile alternative is to look at the following non-canonical symplectic structure of the problem, formula_82 Here formula_83 is a non-constant non-canonical symplectic form. General symplectic integrator for non-constant non-canonical symplectic structure, explicit or implicit, is not known to exist. However, for this specific problem, a family of high-order explicit non-canonical symplectic integrators can be constructed using the He splitting method. Splitting formula_84 into 4 parts, formula_85 we find serendipitously that for each subsystem, e.g., formula_86 and formula_87 the solution map can be written down explicitly and calculated exactly. Then explicit high-order non-canonical symplectic algorithms can be constructed using different compositions. Let formula_88 and formula_89 denote the exact solution maps for the 4 subsystems. A 1st-order symplectic scheme is formula_90 A symmetric 2nd-order symplectic scheme is, formula_91 which is a customarily modified Strang splitting. A formula_92-th order scheme can be constructed from a formula_93-th order scheme using the method of triple jump, formula_94 The He splitting method is one of key techniques used in the structure-preserving geometric particle-in-cell (PIC) algorithms 