Symplectic matrices

ExtendedSymplecticMetric <: AbstractMetric

The extension of the RealSymplecticMetric at a point $p \in \mathrm{Sp}(2n)$ as an inner product over the embedding space $ℝ^{2n×2n}$, i.e.

\[ ⟨x, y⟩_p = ⟨p^{-1}x, p^{-1}⟩_{\mathrm{Fr}} = \operatorname{tr}(x^{\mathrm{T}}(pp^{\mathrm{T}})^{-1}y), \text{ for all } x, y \in ℝ^{2n×2n}.\]

RealSymplecticMetric <: RiemannianMetric

The canonical Riemannian metric on the symplectic manifold, defined pointwise for $p \in \mathrm{Sp}(2n)$ by [Fio11]]

\[\begin{align*} & g_p \colon T_p\mathrm{Sp}(2n)×T_p\mathrm{Sp}(2n) → ℝ, \\ & g_p(Z_1, Z_2) = \operatorname{tr}((p^{-1}Z_1)^{\mathrm{T}} (p^{-1}Z_2)). \end{align*}\]

This metric is also the default metric for the SymplecticMatrices manifold.

SymplecticElement{T}

A lightweight structure to represent the action of the matrix representation of the canonical symplectic form,

\[J_{2n}(λ) = λ\begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix} ∈ ℝ^{2n×2n},\]

where we write $J_{2n} = J_{2n}(1)$ for short. The canonical symplectic form is represented by

\[\omega_{2n}(x, y) = x^{\mathrm{T}}J_{2n}y, \quad x, y ∈ ℝ^{2n}.\]

The entire matrix is however not instantiated in memory, instead a scalar $λ$ of type T is stored, which is used to keep track of scaling and transpose operations applied to each SymplecticElement. This type acts similar to I from LinearAlgeba.

Constructor

SymplecticElement(λ=1)

Generate the sumplectic matrix with scaling $1$.

SymplecticMatricesMatrices{T, 𝔽} <: AbstractEmbeddedManifold{𝔽, DefaultIsometricEmbeddingType}

The symplectic manifold consists of all $2n×2n$ matrices which preserve the canonical symplectic form over $𝔽^{2n×2n}×𝔽^{2n×2n}$,

\[ \omega\colon 𝔽^{2n×2n}×𝔽^{2n×2n} → 𝔽, \quad \omega(x, y) = p^{\mathrm{T}} J_{2n} q, \ x, y \in 𝔽^{2n×2n},\]

where $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

The symplectic manifold consists of

\[\mathrm{Sp}(2n, ℝ) = \bigl\{ p ∈ ℝ^{2n×2n} \, \big| \, p^{\mathrm{T}}J_{2n}p = J_{2n} \bigr\},\]

The tangent space at a point $p$ is given by [BZ21]

\[\begin{align*} T_p\mathrm{Sp}(2n) &= \{X \in ℝ^{2n×2n} \ |\ p^{T}J_{2n}X + X^{T}J_{2n}p = 0 \}, \\ &= \{X = pJ_{2n}S \ \mid\ S ∈ R^{2n×2n}, S^{\mathrm{T}} = S \}. \end{align*}\]

Constructor

SymplecticMatrices(2n, field=ℝ; parameter::Symbol=:type)

Generate the (real-valued) symplectic manifold of $2n×2n$ symplectic matrices. The constructor for the SymplecticMatrices manifold accepts the even column/row embedding dimension $2n$ for the real symplectic manifold, $ℝ^{2n×2n}$.

exp(M::SymplecticMatrices, p, X)
exp!(M::SymplecticMatrices, q, p, X)

The Exponential mapping on the Symplectic manifold with the RealSymplecticMetric Riemannian metric.

For the point $p \in \mathrm{Sp}(2n)$ the exponential mapping along the tangent vector $X \in T_p\mathrm{Sp}(2n)$ is computed as [WSF18]

\[ \operatorname{exp}_p(X) = p \operatorname{Exp}((p^{-1}X)^{\mathrm{T}}) \operatorname{Exp}(p^{-1}X - (p^{-1}X)^{\mathrm{T}}),\]

where $\operatorname{Exp}(⋅)$ denotes the matrix exponential.

inv(::SymplecticMatrices, A)
inv!(::SymplecticMatrices, A)

Compute the symplectic inverse $A^+$ of matrix $A ∈ ℝ^{2n×2n}$. See symplectic_inverse for details.

rand(::SymplecticStiefel; vector_at=nothing, σ::Real=1.0)

Generate a random point on $\mathrm{Sp}(2n)$ or a random tangent vector $X \in T_p\mathrm{Sp}(2n)$ if vector_at is set to a point $p \in \mathrm{Sp}(2n)$.

A random point on $\mathrm{Sp}(2n)$ is constructed by generating a random Hamiltonian matrix $Ω \in \mathfrak{sp}(2n,F)$ with norm σ, and then transforming it to a symplectic matrix by applying the Cayley transform

\[ \operatorname{cay}: \mathfrak{sp}(2n,F) → \mathrm{Sp}(2n), \ \Omega \mapsto (I - \Omega)^{-1}(I + \Omega).\]

To generate a random tangent vector in $T_p\mathrm{Sp}(2n)$, this code employs the second tangent vector space parametrization of SymplecticMatrices. It first generates a random symmetric matrix $S$ by S = randn(2n, 2n) and then symmetrizes it as S = S + S'. Then $S$ is normalized to have Frobenius norm of σ and X = pJS is returned, where J is the SymplecticElement.

gradient(M::SymplecticMatrices, f, p, backend::RiemannianProjectionBackend;
         extended_metric=true)
gradient!(M::SymplecticMatrices, f, p, backend::RiemannianProjectionBackend;
         extended_metric=true)

Compute the manifold gradient $\text{grad}f(p)$ of a scalar function $f \colon \mathrm{Sp}(2n) → ℝ$ at $p \in \mathrm{Sp}(2n)$.

The element $\text{grad}f(p)$ is found as the Riesz representer of the differential $\text{D}f(p) \colon T_p\mathrm{Sp}(2n) → ℝ$ with respect to the Riemannian metric inner product at $p$ [Fio11]]. That is, $\text{grad}f(p) \in T_p\mathrm{Sp}(2n)$ solves the relation

\[ g_p(\text{grad}f(p), X) = \text{D}f(p) \quad\forall\; X \in T_p\mathrm{Sp}(2n).\]

The default behaviour is to first change the representation of the Euclidean gradient from the Euclidean metric to the RealSymplecticMetric at $p$, and then we projecting the result onto the correct tangent tangent space $T_p\mathrm{Sp}(2n, ℝ)$ w.r.t the Riemannian metric $g_p$ extended to the entire embedding space.

Arguments:

• extended_metric = true: If true, compute the gradient $\text{grad}f(p)$ by first changing the representer of the Euclidean gradient of a smooth extension of $f$, $∇f(p)$, with respect to the RealSymplecticMetric at $p$ extended to the entire embedding space, before projecting onto the correct tangent vector space with respect to the same extended metric $g_p$. If false, compute the gradient by first projecting $∇f(p)$ onto the tangent vector space, before changing the representer in the tangent vector space to comply with the RealSymplecticMetric.

riemannian_gradient(M::SymplecticMatrices, p, Y)

Given a gradient $Y = \operatorname{grad} \tilde f(p)$ in the embedding $ℝ^{2n×2n}$ or at least around the SymplecticMatrices M where p (the embedding of) a point on M, we restrict $\tilde f$ to the manifold and denote that by $f$. Then the Riemannian gradient $X = \operatorname{grad} f(p)$ is given by

\[ X = Yp^{\mathrm{T}}p + J_{2n}pY^{\mathrm{T}}J_{2n}p,\]

where $J_{2n}$ denotes the SymplecticElement.

inv!(M::SymplecticMatrices, A)

Compute the symplectic_inverse of a suqare matrix A inplace of A

project_normal!(::MetricManifold{𝔽,<:Euclidean,ExtendedSymplecticMetric}, Y, p, X)

Project onto the normal of the tangent space $(T_p\mathrm{Sp}(2n))^{\perp_g}$ at a point $p ∈ \mathrm{Sp}(2n)$, relative to the riemannian metric $g$ RealSymplecticMetric.

That is,

\[(T_p\mathrm{Sp}(2n))^{\perp_g} = \{Y ∈ ℝ^{2n×2n} : g_p(Y, X) = 0 \test{ for all } X \in T_p\mathrm{Sp}(2n)\}.\]

The closed form projection operator onto the normal space is given by [GSAS21]

\[\operatorname{P}^{(T_p\mathrm{Sp}(2n))\perp}_{g_p}(X) = pJ_{2n}\operatorname{skew}(p^{\mathrm{T}}J_{2n}^{\mathrm{T}}X),\]

where $\operatorname{skew}(A) = \frac{1}{2}(A - A^{\mathrm{T}})$ and $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

This function is not exported.

symplectic_inverse(A)

Given a matrix

\[ A ∈ ℝ^{2n×2k},\quad A = \begin{bmatrix} A_{1,1} & A_{1,2} \\ A_{2,1} & A_{2, 2} \end{bmatrix}\]

the symplectic inverse is defined as:

\[A^{+} := J_{2k}^{\mathrm{T}} A^{\mathrm{T}} J_{2n},\]

where $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

The symplectic inverse of A can be expressed explicitly as:

\[A^{+} = \begin{bmatrix} A_{2, 2}^{\mathrm{T}} & -A_{1, 2}^{\mathrm{T}} \\[1.2mm] -A_{2, 1}^{\mathrm{T}} & A_{1, 1}^{\mathrm{T}} \end{bmatrix}.\]

symplectic_inverse_times(::SymplecticMatrices, p, q)
symplectic_inverse_times!(::SymplecticMatrices, A, p, q)

Directly compute the symplectic inverse of $p \in \mathrm{Sp}(2n)$, multiplied with $q \in \mathrm{Sp}(2n)$. That is, this function efficiently computes $p^+q = (J_{2n}p^{\mathrm{T}}J_{2n})q ∈ ℝ^{2n×2n}$, where $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

change_representer(MetMan::MetricManifold{<:Any, <:Euclidean, ExtendedSymplecticMetric},
                   EucMet::EuclideanMetric, p, X)
change_representer!(MetMan::MetricManifold{<:Any, <:Euclidean, ExtendedSymplecticMetric},
                    Y, EucMet::EuclideanMetric, p, X)

Change the representation of a matrix $ξ ∈ ℝ^{2n×2n}$ into the inner product space $(ℝ^{2n×2n}, g_p)$ where the inner product is given by $g_p(ξ, η) = \langle p^{-1}ξ, p^{-1}η \rangle = \operatorname{tr}(ξ^{\mathrm{T}}(pp^{\mathrm{T}})^{-1}η)$, as the extension of the RealSymplecticMetric onto the entire embedding space.

By changing the representation we mean to apply a mapping

\[ c_p : ℝ^{2n×2n} → ℝ^{2n×2n},\]

defined by requiring that it satisfy the metric compatibility condition

\[ g_p(c_p(ξ), η) = ⟨p^{-1}c_p(ξ), p^{-1}η⟩ = ⟨ξ, η⟩^{\text{Euc}} \;∀\; η ∈ T_p\mathrm{Sp}(2n, ℝ).\]

In this case, we compute the mapping

\[ c_p(ξ) = pp^{\mathrm{T}} ξ.\]

change_representer(::SymplecticMatrices, ::EuclideanMetric, p, X)
change_representer!(::SymplecticMatrices, Y, ::EuclideanMetric, p, X)

Compute the representation of a tangent vector $ξ ∈ T_p\mathrm{Sp}(2n, ℝ)$ s.t.

\[ g_p(c_p(ξ), η) = ⟨ξ, η⟩^{\text{Euc}} \text{for all } η ∈ T_p\mathrm{Sp}(2n, ℝ).\]

with the conversion function

\[ c_p : T_p\mathrm{Sp}(2n, ℝ) → T_p\mathrm{Sp}(2n, ℝ), \quad c_p(ξ) = \frac{1}{2} pp^{\mathrm{T}} ξ + \frac{1}{2} pJ_{2n} ξ^{\mathrm{T}} pJ_{2n},\]

where $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

Each of the terms $c_p^1(ξ) = p p^{\mathrm{T}} ξ$ and $c_p^2(ξ) = pJ_{2n} ξ^{\mathrm{T}} pJ_{2n}$ from the above definition of $c_p(η)$ are themselves metric compatible in the sense that

\[ c_p^i : T_p\mathrm{Sp}(2n, ℝ) → ℝ^{2n×2n}\quad g_p^i(c_p(ξ), η) = ⟨ξ, η⟩^{\text{Euc}} \;∀\; η ∈ T_p\mathrm{Sp}(2n, ℝ),\]

for $i \in {1, 2}$. However the range of each function alone is not confined to $T_p\mathrm{Sp}(2n, ℝ)$, but the convex combination

\[ c_p(ξ) = \frac{1}{2}c_p^1(ξ) + \frac{1}{2}c_p^2(ξ)\]

does have the correct range $T_p\mathrm{Sp}(2n, ℝ)$.

check_point(M::SymplecticMatrices, p; kwargs...)

Check whether p is a valid point on the SymplecticMatrices M=$\mathrm{Sp}(2n)$, i.e. that it has the right AbstractNumbers type and $p^{+}p$ is (approximately) the identity, where $A^+$ denotes the symplectic_inverse.

The tolerance can be set with kwargs....

check_vector(M::SymplecticMatrices, p, X; kwargs...)

Checks whether X is a valid tangent vector at p on the SymplecticMatrices M=$\mathrm{Sp}(2n)$, which requires that

\[p^{T}J_{2n}X + X^{T}J_{2n}p = 0\]

holds (approximately), where $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

The tolerance can be set with kwargs...

distance(M::SymplecticMatrices, p, q)

Compute an approximate geodesic distance between two Symplectic matrices $p, q \in \mathrm{Sp}(2n)$, as done in [WSF18].

\[ \operatorname{dist}(p, q) ≈ \lVert\operatorname{Log}(p^+q)\rVert_{\mathrm{Fr}},\]

where the $\operatorname{Log}(⋅)$ operator is the matrix logarithm.

This approximation is justified by first recalling the Baker-Campbell-Hausdorf formula,

\[\operatorname{Log}(\operatorname{Exp}(A)\operatorname{Exp}(B)) = A + B + \frac{1}{2}[A, B] + \frac{1}{12}[A, [A, B]] + \frac{1}{12}[B, [B, A]] + \ldots \;.\]

Then we write the expression for the exponential map from $p$ to $q$ as

\[ q = \operatorname{exp}_p(X) = p \operatorname{Exp}((p^{+}X)^{\mathrm{T}}) \operatorname{Exp}([p^{+}X - (p^{+}X)^{\mathrm{T}}]), X \in T_p\mathrm{Sp},\]

and with the geodesic distance between $p$ and $q$ given by

\[\operatorname{dist}(p, q) = \lVert X \rVert_p = \lVert p^+ X \rVert_{\mathrm{Fr}}\]

we see that

\[ \begin{align*} \lVert\operatorname{Log}(p^+q)\rVert_{\mathrm{Fr}} &=\Bigl\lVert \operatorname{Log}\bigl( \operatorname{Exp}((p^{+}X)^{\mathrm{T}}) \operatorname{Exp}(p^{+}X - (p^{+}X)^{\mathrm{T}}) \bigr) \Bigr\rVert_{\mathrm{Fr}} \\ &=\lVert p^{+}X + \frac{1}{2}[(p^{+}X)^{\mathrm{T}}, p^{+}X - (p^{+}X)^{\mathrm{T}}] + \ldots\lVert_{\mathrm{Fr}} \\ &≈\lVert p^{+}X\rVert_{\mathrm{Fr}} = \operatorname{dist}(p, q). \end{align*}\]

inner(::SymplecticMatrices{<:Any,ℝ}, p, X, Y)

Compute the canonical Riemannian inner product RealSymplecticMetric

\[ g_p(X, Y) = \operatorname{tr}((p^{-1}X)^{\mathrm{T}} (p^{-1}Y))\]

between the two tangent vectors $X, Y \in T_p\mathrm{Sp}(2n)$.

inverse_retract(M::SymplecticMatrices, p, q, ::CayleyInverseRetraction)

Compute the Cayley Inverse Retraction $X = \mathcal{L}_p^{\mathrm{Sp}}(q)$ such that the Cayley Retraction from $p$ along $X$ lands at $q$, i.e. $\mathcal{R}_p(X) = q$ [BZ21].

For $p, q ∈ \mathrm{Sp}(2n, ℝ)$ then, we can define the inverse cayley retraction as long as the following matrices exist.

\[ U = (I + p^+ q)^{-1}, \quad V = (I + q^+ p)^{-1},\]

where $(⋅)^+$ denotes the symplectic_inverse.

Then inverse cayley retraction at $p$ applied to $q$ is

\[\mathcal{L}_p^{\mathrm{Sp}}(q) = 2p\bigl(V - U\bigr) + 2\bigl((p + q)U - p\bigr) ∈ T_p\mathrm{Sp}(2n).\]

is_flat(::SymplecticMatrices)

Return false. SymplecticMatrices is not a flat manifold.

manifold_dimension(::SymplecticMatrices)

Returns the dimension of the symplectic manifold embedded in $ℝ^{2n×2n}$, i.e.

\[ \operatorname{dim}(\mathrm{Sp}(2n)) = (2n + 1)n.\]

project!(::MetricManifold{𝔽,<:Euclidean,ExtendedSymplecticMetric}, Y, p, X) where {𝔽}

Compute the projection of $X ∈ R^{2n×2n}$ onto $T_p\mathrm{Sp}(2n, ℝ)$ with respect to the RealSymplecticMetric $g$.

The closed form projection mapping is given by [GSAS21]

\[ \operatorname{P}^{T_p\mathrm{Sp}(2n)}_{g_p}(X) = pJ_{2n}\operatorname{sym}(p^{\mathrm{T}}J_{2n}^{\mathrm{T}}X),\]

where $\operatorname{sym}(A) = \frac{1}{2}(A + A^{\mathrm{T}})$ and and $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

project(::SymplecticMatrices, p, A)
project!(::SymplecticMatrices, Y, p, A)

Given a point $p \in \mathrm{Sp}(2n)$, project an element $A \in ℝ^{2n×2n}$ onto the tangent space $T_p\mathrm{Sp}(2n)$ relative to the euclidean metric of the embedding $ℝ^{2n×2n}$.

That is, we find the element $X \in T_p\operatorname{Sp}(2n)$ which solves the constrained optimization problem

\[ \operatorname{min}_{X \in ℝ^{2n×2n}} \frac{1}{2}\lVert X - A\rVert^2, \quad \text{such that}\; h(X) := X^{\mathrm{T}} J_{2n} p + p^{\mathrm{T}} J_{2n} X = 0,\]

where $h: ℝ^{2n×2n} → \operatorname{skew}(2n)$ denotes the restriction of $X$ onto the tangent space $T_p\operatorname{SpSt}(2n, 2k)$ and $J_{2n} = \begin{bmatrix} 0_n & I_n \\ -I_n & 0_n \end{bmatrix}$ denotes the SymplecticElement.

retract(::SymplecticMatrices, p, X, ::CayleyRetraction)
retract!(::SymplecticMatrices, q, p, X, ::CayleyRetraction)

Compute the Cayley retraction on $p ∈ \mathrm{Sp}(2n, ℝ)$ in the direction of tangent vector $X ∈ T_p\mathrm{Sp}(2n, ℝ)$, as defined in by Birtea et al in proposition 2 [BCC20].

Using the symplectic_inverse $A^+$ of a matrix $A \in ℝ^{2n×2n}$ the retraction $\mathcal{R}: T\mathrm{Sp}(2n) → \mathrm{Sp}(2n)$ is defined pointwise as

\[\begin{align*} \mathcal{R}_p(X) &= p \operatorname{cay}\left(\frac{1}{2}p^{+}X\right), \\ &= p \operatorname{exp}_{1/1}(p^{+}X), \\ &= p (2I - p^{+}X)^{-1}(2I + p^{+}X). \end{align*}\]

Here $\operatorname{exp}_{1/1}(z) = (2 - z)^{-1}(2 + z)$ denotes the Padé (1, 1) approximation to $\operatorname{exp}(z)$.