Diamond Scheme

Mathematical Introduction

The multi-Hamiltonian PDE is $$K\mathbf{z}_t + L\mathbf{z}_x = \nabla S(\mathbf{z}), \label{eqn:hampde}$$ where $K$ and $L$ are constant $n \times n$ real skew-symmetric matrices, $\mathbf{z}: \mathbb{R}^2 \Rightarrow \mathbb{R}^n$, and $S: \mathbb{R}^n \rightarrow \mathbb{R}$. As an example, by introducing $\mathbf{z} = (u, v, w)$, $v = u_t$ and $w = u_x$ the one dimensional wave equation, $u_{tt} - u_{xx} = f'(u)$, can be written in this multi-Hamiltonian form with $$\label{eqn:1dwaveeqn} K = \begin{pmatrix} 0 & -1 & 0\\1 & 0 & 0\\0 & 0 & 0 \end{pmatrix}, \quad L = \begin{pmatrix} 0 & 0 & 1\\0 & 0 & 0\\-1 & 0 & 0 \end{pmatrix}, \quad \mathrm{and} \quad \nabla S(\mathbf{z}) = \begin{pmatrix} -f'(u)\\v\\-w\end{pmatrix}.$$ Other equations that can be written in this form include Schrödinger's equation, the Korteweg-de Vries equation, and Maxwell's electromagnetic wave equation. The solution to this PDE satisfies the multi-symplectic conservation law $$\omega_t + \kappa_x = 0, \label{eqn:symlaw}$$ where $\omega = \frac{1}{2} (d\mathbf{z} \wedge K d\mathbf{z})$ and $\kappa = \frac{1}{2} (d\mathbf{z} \wedge L d\mathbf{z})$. A numerical method that satisfies a discrete version of this conservation law is called a multi-symplectic integrator. The Euler-box and Preissman/Keller-box schemes are two such methods; both have order two. A higher order method is given by Reich. It only works when $n$ is odd, and for problems with periodic boundary conditions, it also requires an expensive non-linear solve at each time step. The method presented here, which is based on Reich's work, overcomes these problems.

Code Introduction

The Diamond Scheme numerically solves the multi-Hamiltonian PDE given above. Downloads are below. Currently there are uniprocessor and MPI python implementations, and a uniprocessor C version. Any problem with periodic boundary conditions is supported. There is a README in the latest version that describes how to run it.