Symplectic Systems
Symplectic systems are ODEs whose vector field has a specific structure that is very restrictive. Before we introduce symplectic vector fields on manifolds we first have to define what a symplectic structure is:
A symplectic structure or symplectic 2-form $\Omega$ assigns to each $x\in\mathcal{M}$ a mapping
\[\Omega_x:T_x\mathcal{M}\times{}T_x\mathcal{M} \to \mathbb{R},\]
such that $\Omega_x$ is skew-symmetric and nondegenerate, and $\Omega$ is closed.
We forego the precise definition of closedness because it would require us to introduce differential forms [16, 27]. This property is also closely related to the Jacobi identity [28, Chapter 4.4]. After having defined a symplectic structure, we can introduce Hamiltonian vector fields[1]:
A Hamiltonian vector field at $x\in\mathcal{M}$ corresponding to the function $H:\mathcal{M}\to\mathbb{R}$ (called the Hamiltonian) is a vector field that has the following property:
\[\Omega(X_H, \dot{\gamma}(0)) = \frac{d}{dt}\bigg|_{t = 0}H(\gamma(t)), \Omega(X_H, \dot{\gamma}(0)) = \frac{d}{dt}\bigg|_{t = 0}H(\gamma(t)),\]
where $\gamma$ is a $C^\infty$ curve through $x$.
Of particular importance for us here are canonical Hamiltonian systems:
To obtain a canonical Hamiltonian system we take $\mathcal{M} = \mathbb{R}^{2d}$ and
\[\Omega_z = \mathbb{J}_{2d}^T = \begin{pmatrix} \mathbb{O} & -\mathbb{I} \\ \mathbb{I} & \mathbb{O} \end{pmatrix}\]
for all $z$. We call $\mathbb{J}_{2d}$ the Poisson tensor. In this case the vector field can be written as:
\[X_H(z) = \mathbb{J}_{2d}\nabla_z{}H,\]
where $\nabla{}H$ is the Euclidean gradient.
Typically we further split the $z$ coordinates on $\mathbb{R}^{2d}$ into $(q, p)$ coordinates, where $q$ and $p$ are the first $d$ components of $z$ and the second $d$ components of $z$ respectively, i.e.
\[z = \begin{pmatrix} q \\ p \end{pmatrix}.\]
We can then reformulate a Hamiltonian vector field as two separate vector fields:
\[\begin{aligned} \dot{q} & = \frac{\partial{}H}{\partial{}p} \quad\text{and} \\ \dot{p} & = - \frac{\partial{}H}{\partial{}q}. \end{aligned}\]
Solution of Symplectic Systems
The flow of a Hamiltonian ODE has very restrictive properties, the most important one of these is called symplecticity [1]. This property dramatically restricts the dynamically accessible states of the flow map. For a canonical Hamiltonian system symplecticity is defined as follows:
A map $\phi:\mathbb{R}^{2d}\to\mathbb{R}^{2d}$ is called symplectic on $U\subset\mathbb{R}^{2d}$ if
\[ (\nabla_z\phi)^T\mathbb{J}_{2d}\nabla_z\phi = \mathbb{J}_{2d}, (\nabla_z\phi)^T\mathbb{J}_{2d}\nabla_z\phi = \mathbb{J}_{2d},\]
for all $z\in{}U.$
A similar definition of symplecticity holds for the general case of symplectic manifolds $(\mathcal{M}, \Omega)$ [27]. The following holds:
The flow of a Hamiltonian system is symplectic
Proof
We proof this statement only for canonical Hamiltonian systems here. Consider the flow of a Hamiltonian ODE $\varphi^t:\mathbb{R}\to\mathbb{R}$. For this we have:
\[\begin{aligned} \frac{d}{dt}\left( (\nabla\varphi^t)^T\mathbb{J}\nabla\varphi^t \right) & = (\mathbb{J}\nabla^2H)^T\mathbb{J}\nabla\varphi^t + (\nabla\varphi^t)^T\mathbb{J}\mathbb{J}\nabla^2H \\ & = \nabla^2H\nabla\varphi^t - (\nabla\varphi^t)^T\nabla^2H, \end{aligned}\]
and this expression is zero at $t=0.$ This computation holds for all points on the domain on which $H$ is defined and $\varphi^t$ is symplectic.
The discipline of finding numerical approximations of flows $\varphi^t$ such that these numerical approximations also preserve certain properties of that flow (such as symplecticity) is referred to as structure-preserving numerical integration or geometric numerical integration [1]. The Julia library GeometricIntegrators [2] offers a wide array of such geometric numerical integrators for a broad class of systems (not just canonical Hamiltonian systems).
In addition to being symplectic, the flow of a Hamiltonian vector field is also energy-preserving:
The flow of a Hamiltonian vector field preserves the Hamiltonian $H,$ i.e.
\[ H(\varphi^t(z_0)) = H(z_0),\]
for all $t\in[0, T].$
Proof
We differentiate $H(\varphi^t(z_0))$ with respect to $t$:
\[ \frac{d}{dt}H(\varphi^t(z_0)) = (\nabla_{\varphi^t(z_0)}H)^T\frac{d}{dt}\varphi^t(z_0) = (\nabla_{\varphi^t(z_0)}H)^T\mathbb{J}\nabla_{\varphi^t(z_0)}H = \mathbb{O},\]
because $\mathbb{J}$ is skew-symmetric.
It is important to note that symplecticity is a very strong property[2] that may not be achievable in some practical applications. If preservation of symplecticity is not achievable, it may however still be advantageous to consider weaker properties such as volume preservation.
Library Functions
GeometricMachineLearning.PoissonTensor — TypePoissonTensor(n2)Returns a (canonical) Poisson tensor of size $2n\times2n$:
\[\mathbb{J}_{2n} = \begin{pmatrix} \mathbb{O} & \mathbb{I}_n \\ -\mathbb{I}_n & \mathbb{O} \\ \end{pmatrix}\]
Arguments
It can also be called with a backend and a type:
using GeometricMachineLearning
backend = CPU()
T = Float16
PoissonTensor(backend, 4, T)
# output
4×4 PoissonTensor{Float16, Matrix{Float16}}:
0.0 0.0 1.0 0.0
0.0 0.0 0.0 1.0
-1.0 0.0 0.0 0.0
0.0 -1.0 0.0 0.0GeometricMachineLearning.QPT — TypeQPTThe type for data in $(q, p)$ coordinates. It encompasses various array types.
Examples
using GeometricMachineLearning: QPT
# allocate two vectors
data1 = (q = rand(5), p = rand(5))
# allocate two matrices
data2 = (q = rand(5, 4), p = rand(5, 4))
# allocate two tensors
data3 = (q = rand(5, 4, 2), p = rand(5, 4, 2))
(typeof(data1) <: QPT, typeof(data2) <: QPT, typeof(data3) <: QPT)
# output
(true, true, true)We can also do:
using GeometricMachineLearning: QPT, PoissonTensor
𝕁 = PoissonTensor(4)
qp = (q = [1, 2], p = [3, 4])
𝕁 * qp
# output
(q = [3, 4], p = [-1, -2])GeometricMachineLearning.QPTOAT — TypeQPTOATA union of two types:
const QPTOAT = Union{QPT, AbstractArray}This could be data in $(q, p)\in\mathbb{R}^{2d}$ form or come from an arbitrary vector space.
References
- [27]
- V. I. Arnold. Mathematical methods of classical mechanics. Vol. 60 of Graduate Texts in Mathematics (Springer Verlag, Berlin, 1978).
- [16]
- S. I. Richard L. Bishop. Tensor Analysis on Manifolds (Dover Publications, Mineola, New York, 1980).
- [1]
- E. Hairer, C. Lubich and G. Wanner. Geometric Numerical integration: structure-preserving algorithms for ordinary differential equations (Springer, Heidelberg, 2006).
- 1Also compare this to the definition of the Riemannian gradient.
- 2Symplecticity imposes in general greater restrictions on the flow map than e.g. conservation of energy. The famous Ge-Marsden theorem [29] says that one cannot achieve preservation of energy and preservation of symplecticity at the same time, unless the numerical method is the exact integral of the Hamiltonian ODE. In practice we hence almost always choose a symplectic over an energy-preserving scheme.