Klas Modin PRO
Mathematician at Chalmers University of Technology and the University of Gothenburg
Lie-Poisson integrators for
quasi-geostrophic equations
Klas Modin
Collaborator: Milo Viviani
Quasi-geostrophic equation
v˙+∇vv=−∇p−2Ω~×v
\dot v + \nabla_v v = -\nabla p - 2 \tilde\Omega\times v
divv=0
\operatorname{div} v = 0
3D Navier-Stokes equation
Shallow water equation
Quasi-geostrophic equation
Quasi-geostrophic equation
v˙+∇vv=−∇p−2Ω~×v
\dot v + \nabla_v v = -\nabla p - 2 \tilde\Omega\times v
divv=0
\operatorname{div} v = 0
Jule Charney
3D Navier-Stokes equation
Shallow water equation
Quasi-geostrophic equation
Quasi-geostrophic equation
v˙+∇vv=−∇p−2Ω~×v
\dot v + \nabla_v v = -\nabla p - 2 \tilde\Omega\times v
divv=0
\operatorname{div} v = 0
3D Navier-Stokes equation
Shallow water equation
Quasi-geostrophic equation
horizontal scale ≫ vertical scale
Quasi-geostrophic equation
v˙+∇vv=−∇p−2Ω~×v
\dot v + \nabla_v v = -\nabla p - 2 \tilde\Omega\times v
divv=0
\operatorname{div} v = 0
3D Navier-Stokes equation
Shallow water equation
Quasi-geostrophic equation
horizontal scale ≫ vertical scale
Quasi-geostrophic equation
v˙+∇vv=−∇p−2Ω~×v
\dot v + \nabla_v v = -\nabla p - 2 \tilde\Omega\times v
divv=0
\operatorname{div} v = 0
3D Navier-Stokes equation
Shallow water equation
Quasi-geostrophic equation
Coriolis + pressure ≫ inertial forces
Vorticity formulation
v˙+∇vv=−∇p−2Ω~×v
\dot v + \nabla_v v = -\nabla p - 2 \tilde\Omega\times v
divv=0
\operatorname{div} v = 0
Apply curl to v
ω˙−{ψ,ω}=0
\dot\omega - \{\psi,\omega \} = 0
Δψ=ω−f
\Delta\psi = \omega - f
ω˙−Lvω=0
\dot\omega - L_v\omega = 0
ω:S2→Rψ:S2→R
\omega: S^2\to \mathbb{R} \qquad \psi: S^2\to \mathbb{R}
level-sets of ω
Geometry of QGE
Lie-Poisson system on Xμ(S2)∗≃C0∞(S2)
G=Diffμ(S2)
Te∗G≃g∗
Casimir functions:
Cf(ω)=∫S2f(ω)μ
\displaystyle\mathcal C_f(\omega) = \int_{S^2}f(\omega)\mu
Finite-dim (weak) co-adjoint orbits:
ω=k=1∑NΓkδqk
\displaystyle\omega = \sum_{k=1}^N \Gamma_k \delta_{q^k}
Statistical mechanics theories for 2D Euler
Idea by Onsager (1949):
- Approximate ω by PV for large N
- Find invariant measure and presume ergodicity
Miller (1990) and Robert & Sommeria (1991): (MRS)
- Minimize microcanonical entropy under energy and Casimir constraints
Is the MRS prediction correct?
No!
2D Euler equations are not ergodic
...but perhaps MRS is "generically" correct
Flow ergodic except at "KAM islands"
Poincaré section of finite dimensional Hamiltonian system
Ultimate problem
for geometric integration
To test MRS we need to:
- Run long simulations
- Preserve the Casimirs
(energy + enstrophy alone not enough) - Preserve the Lie-Poisson structure
(criterion in MRS)
On T2 such discretization exists (sine-bracket)
[Zeitlin 1991, McLachlan 1993]
based on quantization theory by Hoppe (1989)
Numerical simulations support MRS on T2
[Abramov & Majda 2003]
A torus is not a sphere
MRS generally assumed valid also on S2
However, non-structure preserving simulations by Dritschel, Qi, & Marston (2015) contradict MRS on S2
DQM simulation yield persistent unsteadiness
Our mission: construct trustworthy discretization on S2
2D Euler to isospectral flow via Berezin-Toeplitz quantization
C∞(M)∋f↦TfN∈gN
C^\infty(M)\ni f \mapsto T^N_f \in \mathfrak{g}_N
Exists if M compact quantizable Kähler manifold
Idea: approximate Poisson algebra with matrix algebras
{f,g}→[TfN,TgN]N→∞
\{f,g\} \to [T^N_f,T^N_g] \quad N\to\infty
ω˙={Δ−1ω,ω}
\displaystyle\dot \omega = \left\{\Delta^{-1}\omega,\omega \right\}
W˙=[ΔN−1W,W]
\displaystyle\dot W = [\Delta_N^{-1}W,W]
From 2D Euler
To isospectral
ω↦W
\omega \mapsto W
Lie-Poisson isospectral flows
W˙=[B(W),W]
\dot W = [B(W),W]
Let B:g→g
isospectral flow
Analytic function f yields first integral
Cf(W):=tr(f(W))
C_f(W) := \mathrm{tr}(f(W))
Casimir function
Hamiltonian case
B(W)=∇H(W)†
B(W) = \nabla H(W)^\dagger
Hamiltonian function
Note: Non-canonical Poisson structure (Lie-Poisson)
Explicit B-T quantization on S2
C∞(S2)∋ω↦W∈su(N)
C^\infty(\mathbb S^2)\ni \omega\mapsto W \in \mathfrak{su}(N)
[Hoppe, 1989]
- Express ω in spherical harmonics expansion ω=l=1∑∞m=−l∑lωlmYlm
- Truncate at lmax=N−1
- For fixed m, linear map between (ωlm)l=1N−m) and m:th diagonal of W
- Gives N linear maps
Complicated coefficients, expressed by Wigner 3-j symbols of very high order
~2 weeks to compute coefficients for N=1025
Discrete S2 Laplacian on su(N)
- "Magic" formula [Hoppe & Yau, 1998]
ΔNW=2N2−1([XN,[XN,W]]−21[X+N,[X−N,W]]−21[X−N,[X+N,W]])
banded matrices
W˙=[ΔN−1W,W]
\displaystyle\dot W = [\Delta_N^{-1}W,W]
Recall
What is ΔN and how compute ΔN−1W ?
(Naive approach requires O(N3) operations with large constant)
O(N2) operations
- ΔN admits sparse LU-factorization with O(N2 non-zeros
Spatial discretization obtained!
W˙=[ΔN−1W,W]
\displaystyle \dot W = [\Delta_N^{-1}W,W]
Note: corresponds to
N2 spherical harmonics
O(N2) operations
O(N3) operations
Isospectral flow ⇒ discrete Casimirs
CfN(W)=tr(f(W))
\displaystyle C^N_f(W) = \operatorname{tr}(f(W))
In progress: spatial convergence
Conjecture [M. & Viviani, 2020]
Fixed interval [0,T], constant C=C(T,ω0) s.t. ∥ω(t,⋅)−TN−1(WN(t))∥∞≤C/N2
Classical global existence and uniqueness in L∞ setting
⟹ L∞-norm conserved (it's a Casimir)
Toeplitz-Berezin quantization theory gives
∥TN−1([TN(ω1),TN(ω2)])−{ω1,ω2}∥∞=O(1/N2)
\lVert T_N^{-1}([T_N(\omega_1),T_N(\omega_2)]) - \{\omega_1,\omega_2 \} \rVert_{\infty} = O(1/N^2)
Time discretization
Aim: numerical integrator that is
-
isospectral, Wk→Wk+1 an isospectral map
necessary to preserve Casimirs
-
symplectic, Wk→Wk+1 a Lie-Poisson map su(N)∗→su(N)∗
necessary to (nearly) preserve energy and phase space structure
What about symplectic Runge-Kutta methods (SRK)?
- Not Lie-Poisson preserving!
- Not isospectral!
Isospectral Symplectic
Runge-Kutta methods
[M. & Viviani 2019]
W˙=[B(W),W]
\dot W = [B(W),W]
Given s-stage Butcher tableau (aij,bi) for SRK
XiYiKijW~iWk+1=−(Wn+j=1∑saijXj)hB(W~i)=hB(W~i)(Wk+j=1∑saijYj=hB(W~i)(j′=1∑s(aij′Xj′+ajj′Kij′))=Wk+j=1∑saij(Xj+Yj+Kij=Wk+i=1∑s[hB(W~i),W~i]
\begin{aligned}
X_i &= - (W_n + \sum_{j=1}^s a_{ij}X_j)h B(\tilde W_i) \\
Y_i &= hB(\tilde W_i)(W_k + \sum_{j=1}^sa_{ij}Y_j \\
K_{ij} &= h B(\tilde W_i)(\sum_{j'=1}^s(a_{ij'}X_{j'}+a_{jj'}K_{ij'})) \\
\tilde W_i &= W_k + \sum_{j=1}^s a_{ij}(X_j + Y_j + K_{ij} \\
W_{k+1} &= W_k + \sum_{i=1}^s [hB(\tilde W_i),\tilde W_i]
\end{aligned}
Theorem: method is isospectral and Lie-Poisson preserving on any reductive Lie algebra
Numerical results for QGE
Zonal jet and vortex structures on Jupiter
Copyright: NASA, Cassini Imaging Team
Simulation of unstable Rossby-Haurwitz wave
Thank you!
M. & Viviani
A Casimir preserving scheme for long-time simulation of spherical ideal hydrodynamics
J. Fluid Mech., 2020
M. & Viviani
Lie-Poisson methods for isospectral flows
Found. Comp. Math., 2020
- Slides at slides.com/kmodin
- Videos at bitbucket.org/kmodin/euler-sphere-quantization
Lie-Poisson integrators for quasi-geostrophic equations Klas Modin
Berezin-Toeplitz quantization and Lie-Poisson integrators for quasi-geostrophic equations
By Klas Modin
Berezin-Toeplitz quantization and Lie-Poisson integrators for quasi-geostrophic equations
Presentation given 2020-06 at the online-only FoCM Workshop "Geometric Integration and Computational Mechanics" June, 15-18, 2020.
- 1,396
