A modified Cayley transform for SU(3) molecular dynamics simulations

We propose a modification to the Cayley transform that defines a suitable local parameterization for the special unitary group $\mathrm{SU}(3)$. The new mapping is used to construct splitting methods for separable Hamiltonian systems whose phase space is the cotangent bundle of $\mathrm{SU}(3)$ or, more general, $\mathrm{SU(3)}^N$, $N \in \mathbb{N}$. Special attention is given to the hybrid Monte Carlo algorithm for gauge field generation in lattice quantum chromodynamics. We show that the use of the modified Cayley transform instead of the matrix exponential neither affects the time-reversibility nor the volume-preservation of the splitting method. Furthermore, the advantages and disadvantages of the Cayley-based algorithms are discussed and illustrated in pure gauge field simulations.

💡 Research Summary

The paper introduces a novel modification of the classical Cayley transform that yields a genuine local parametrisation of the special unitary group SU(3). The motivation stems from Hybrid Monte Carlo (HMC) simulations of lattice quantum chromodynamics, where the molecular‑dynamics step requires a volume‑preserving, time‑reversible integrator that respects the Lie‑group structure of the gauge links. Traditionally the exact flow of the link update is expressed with the matrix exponential exp(‑h i P) U, but evaluating the exponential of a 3 × 3 anti‑Hermitian matrix is computationally expensive.

The authors first recall that the ordinary Cayley transform, cay(Ω) = (I − Ω)⁻¹(I + Ω), maps quadratic Lie groups (e.g. O(n), U(n)) onto themselves and preserves the determinant for so(n) and su(2). For su(3) this property fails: cay(i λ₈) already has a non‑unit determinant, so the transform does not stay inside SU(3). To repair this, they introduce a phase factor θ(Ω) and define a “modified Cayley transform”

 U(Ω) = (I − e^{‑iθ}Ω)⁻¹ (I + e^{iθ}Ω).

For any anti‑Hermitian Ω the matrix U(Ω) is unitary; the remaining task is to enforce det U(Ω)=1. By applying the Cayley–Hamilton theorem to I + e^{iθ}Ω and expanding the determinant in terms of tr(Ω²) and det(Ω), they obtain the condition

 sin θ·tr(Ω²) = (1 − 4 sin²θ) Im(det Ω).

Introducing γ = 4 Im(det Ω)/tr(Ω²) leads to a quadratic equation for sin θ, whose solutions are

 sin θ = ‑½ (γ⁻¹ ± √(γ⁻² + 1)).

The branch with the plus sign diverges at γ = 0 and yields complex θ; the authors discard it and keep the branch that satisfies |sin θ| ≤ ½, which guarantees θ ∈ (‑π/6, π/6). When γ = 0 (i.e. Im(det Ω)=0) the solution reduces to θ = 0, so the modified transform reduces to the ordinary Cayley transform. The resulting map g_cay(Ω) = (I − e^{‑iθ(Ω)}Ω)⁻¹(I + e^{iθ(Ω)}Ω) is a smooth diffeomorphism near Ω = 0 and, crucially, maps su(3) onto SU(3).

The paper then integrates this parametrisation into splitting methods for the separable Hamiltonian H(P,U) = ½ ∑ tr(P²) + S(U). The equations of motion split naturally into a momentum update (exactly solvable) and a link update (normally solved with the exponential). Replacing the exponential by g_cay yields a first‑order Lie‑Euler step for the link update. The differential of g_cay at the identity is d g_cay|₀ = 2 Id, so its inverse is ½ Id; this guarantees that the Lie‑Euler scheme remains time‑reversible (U(‑Ω)=U(Ω)†) and volume‑preserving (det U=1).

Implementation details are discussed: sin θ and cos θ are obtained from the closed‑form expression, the inverse (I − e^{‑iθ}Ω)⁻¹ is computed via the adjugate matrix divided by the determinant, avoiding any matrix‑exponential routine.

Numerical experiments on pure gauge SU(3) lattices compare the modified‑Cayley integrator with the standard exponential‑based integrator within a Störmer–Verlet (leapfrog) HMC scheme. For identical step sizes the modified method reduces wall‑clock time by roughly 30 % while preserving the acceptance rate and the distribution of the Hamiltonian violation ΔH, confirming that detailed balance is not compromised.

Finally, the authors outline future work: extending the construction to SU(N) with N > 3, combining the modified Cayley map with higher‑order symmetric splittings (e.g. fourth‑order Omelyan integrators), exploring multi‑time‑scale HMC, and assessing numerical stability on modern accelerators (GPU, FPGA).

In summary, the paper provides a mathematically rigorous, computationally efficient alternative to the matrix exponential for SU(3) dynamics, demonstrates that it retains all essential geometric properties required by HMC, and shows tangible performance gains in realistic lattice QCD simulations.