Local implicit midpoint step for a Hamiltonian vector field

Let $d\ge 1$, let $U\subseteq {\mathbb{R}}^{2d}$ be open, and let $H\in C^2(U)$. Write $z=(q,p)$ in canonical coordinates, set $J=\left(\begin{array}{cc}0& I\\ -I& 0\end{array}\right)$, and use the Euclidean norm and its induced operator norm. Define the finite-regularity coordinate field $F(q,p):=({\mathrm{\nabla }}_{p}H(q,p),-{\mathrm{\nabla }}_{q}H(q,p))=J\mathrm{\nabla }H(q,p)\mathrm{.}$

Given $z\in U$ and a signed step $h$, supply a candidate radius $\rho \ge 0$. If no radius is supplied or $\rho <0$, return branch-not-certified. If the closed ball $\bar{B}(z,\rho )$ is certified not to lie in $U$, return domain-failure; if its containment in $U$ cannot be certified, return branch-not-certified. For a ball contained in $U$, supply certified bounds $M,L\ge 0$ with $M\ge {\sup }_{W\in \bar{B}(z,\rho )}\mathrm{\parallel }F(W)\mathrm{\parallel },L\ge {\sup }_{W\in \bar{B}(z,\rho )}\mathrm{\parallel }DF(W)\mathrm{\parallel }\mathrm{.}$ If either bound or its certificate is unavailable, or if either admissibility inequality $\mathrm{\mid }h\mathrm{\mid }M\mathrm{/}2\le \rho$ or $\mathrm{\mid }h\mathrm{\mid }L\mathrm{/}2<1$ fails, return branch-not-certified. Otherwise the exact local implicit midpoint step selects the unique $Z\in \bar{B}(z,\rho )$ satisfying $Z=z+\frac{h}{2}F(Z)$ and returns $z^{+}=2Z-z$. The returned state is in ${\mathbb{R}}^{2d}$; a further step is admitted only when that state lies in $U$ and passes a new ball check.

For a finite iteration implementation of an admitted branch, input an output tolerance $\tau >0$ and an update budget $B\in {\mathbb{N}}_0$. Use $\theta =\mathrm{\mid }h\mathrm{\mid }L\mathrm{/}2$ and proceed as follows.

1. Set $Z_0=z$. If $B=0$, return inconclusive with current stage iterate $Z_0$, approximate output $z$, and output-error bound $2\rho$.
2. If $B\ge 1$, compute $Z_1=z+(h\mathrm{/}2)F(Z_0)$ and set $d_1=\mathrm{\parallel }{Z}_1-Z_0\mathrm{\parallel }$.
3. For $k=1,\dots ,B$, set $b_k={\theta }^k{d}_1\mathrm{/}(1-\theta )$ and the output radius $r_k=2b_k$. If $r_k\le \tau$, return the approximation $2Z_k-z$ with output-error certificate $r_k$. Otherwise, if $k<B$, compute $Z_{k+1}=z+(h\mathrm{/}2)F(Z_k)$ and continue. If the check fails at $k=B$, return inconclusive with current stage iterate $Z_B$, approximate output $2Z_B-z$, and output-error bound $2b_B$.

When $\theta =0$, use ${\theta }^k=0$ for $k\ge 1$; the procedure never evaluates $0^0$. If a finite implementation cannot certify an input or an exact operation needed for these instructions, return inconclusive without a tolerance claim. The certificates concern the displayed exact recurrence and exact output operations; inexact evaluations require additional certified error accounting. An approximate result is not asserted to be exactly symplectic. The exact map is the selected branch value $2Z-z$.