Introduction

Isospectral flows

• $Y_0\in \mathbb{M}(N,\mathbb{C})$, $\mathcal{L}:\mathbb{M}(N,\mathbb{C})\rightarrow\mathbb{M}(N,\mathbb{C})$ linear
• $[A,B] = AB - BA$
• The spectrum of $Y$ is time-invariant
• $\mathcal{L}$ self-adjoint w.r.t. Frobenious inner product $\Rightarrow$ flow Hamiltonian w.r.t. $H=\frac{1}{2}Tr(\mathcal{L}YY^*)$

Geometric integration of isospectral flows

• Second order in time ($h$ is the time-step)
• Spectrum of $Y$ exactly preserved
• Hamiltonian conserved up to $\mathcal{O}(h^2)$

Cubic equation

• $X$ unknown
• $Y\in \mathbb{M}(N,\mathbb{R}) \ or \ \mathbb{M}(N,\mathbb{C})$ given

Thm. Given $Y\in\mathbb{M}(N,\mathbb{C})$, equation (*) has a unique solution for sufficiently small $h>0$ in some neighbourhood of $Y$.
Furthermore, when equation (*) takes place in $\mathfrak{su}(N)$, there exists such a neighbourhood for any $h\leq\frac{1}{3\|\mathcal{L}\|_{op}\|Y\|}$.

Numerical schemes

Explicit fixed point iteration

• Convergence is guaranteed by the Theorem above, for  $h$ sufficiently small and $X_0:=Y$
• The complexity is $\mathcal{O}(N^3)$
• Only two matrices multiplications are required when working in $\mathfrak{su}(N)$

Linear iteration

• Convergence can be derived by the Theorem above, for  $h$ sufficiently small and $X_0:=Y$
• The complexity is $\mathcal{O}(N^3)$
• To solve the second equation we take the LU-factorization of $I\pm hP_k$

Inexact Newton iteration

$G(X):=X - h[\mathcal{L} X,X]- h^2(\mathcal{L} X)X(\mathcal{L} X) - Y$

$DG(X)=I - h(\mathcal{B}_1 + h\mathcal{B}_2)$

$\mathcal{B}_2=(\mathcal{L}\cdot)X(\mathcal{L}X)+(\mathcal{L}X)\cdot(\mathcal{L}X)+(\mathcal{L}X)X(\mathcal{L}\cdot)$

$\mathcal{B}_1=[\mathcal{L}\cdot,X]+[\mathcal{L}X,\cdot]$

Inexact Newton iteration

Need an approximation for $DG(X)^{-1}$

• $DG(X)^{-1}\approx I + h\mathcal{B}_1$ [Exact 1st order]
• $DG(X)^{-1}\approx I + h\mathcal{B}_1 + h^2\mathcal{B}_2$
• $DG(X)^{-1}\approx I + h\mathcal{B}_1 + h^2\mathcal{B}_1^2$
• $DG(X)^{-1}\approx I + h\mathcal{B}_1 + h^2(\mathcal{B}_1^2+\mathcal{B}_2)$ [Exact 2nd order]

Rmk. the second one is the best compromise between exactness and computational cost

• The total complexity is $\mathcal{O}(N^3)$

Numerical experiments

Euler Equations

Spatial semidiscretization of the Euler equations on a compact orientable surface $S$, $\omega_t\in C^\infty(S)$

$\dot{\omega}=\nabla\omega\cdot\nabla^\perp\Delta^{-1}\omega$

Give the Hamiltonian isospectral flow on $\mathfrak{su}(N)$, for $N\geq 1$

$\dot{W}=[W,\Delta_N^{-1}W]$

Rmk. $\Delta_N$ is tridiagonal, the evaluation of $\Delta_N^{-1}W$ (via LU-factorization done once) requires $\mathcal{O}(N^2)$ operations

Euler Equations

Time-step $h=0.5$, normalized initial data

Rmk. Explicit fixed point performs the best in terms of CPU time

Heisenberg spin chain

Spatial semidiscretization of the Landau–Lifshitz–Gilbert Hamiltonian PDE, $\sigma_t:S^1\rightarrow S^2$

$\dot{\sigma}=\sigma\times\partial_{xx}\sigma$

Give the Hamiltonian isospectral flow on $\mathfrak{su}(2)^{\oplus_N}$, for $N\geq 1$

$\dot{S}_i=[S_i,S_{i-1}+S_{i+1}]$

for $i=1,...,N$ with $S_1\equiv S_N$

Time-step $h=0.5$, normalized initial data

Rmk. Explicit fixed point performs the best in terms of stabilty and CPU time. 

Heisenberg spin chain

An unpractical scheme

Quadratic iteration

• The first equation is CARE (Continuous Algebraic Riccati Equation) in $P_k$
  
  
   
• Non-uniqueness issues

Conclusions

Conclusions

• Iterative schemes for the solution of cubic matrix equations arising in the numerical solution of certain conservative PDEs
• Geometrical integrators: preservation of important physical features
• Both fixed point iterations and inexact Newton work well
• Fixed point iterations more efficient for large time-step and large dimensions
• Implicit fixed point iteration is more stable

Thanks for your attention