Improved scaling of hybrid Monte Carlo by using shadow Hamiltonian

Symplectic integrators integrate a modified Hamiltonian system that is close to the Hamiltonian of interest [78]. A cheap approximation to the modified Hamiltonian has been introduced recently. It is called the shadow Hamiltonian method [82] and it can achieve arbitrary accuracy in calculating the modified Hamiltonian using $k$ available values of forces and energies, with an accuracy $p=2k$, that is $\delta H=O(\delta t^{p})$.

The hybrid Monte Carlo partition function is  [56]

$$Z=\int[d\phi][d\pi]e^{-H(\phi,\pi)}.$$

Using any symplectic integrator to integrate the Hamiltonian's equation gives $<e^{-\delta H}>=1$. Expanding in powers of $\delta H$ gives

$$<\delta H>=\frac{1}{2}<\delta H^{2}>+O(\delta H^{3}).$$

Using the shadow Hamiltonian for the MD trajectory and the acceptance rule, the computational effort moving from one configuration to an approximately independent configuration is proportional to $N^{1+\frac{1}{2p}}$ (if we do not count the additional computational cost on calculating the shadow Hamiltonian). The asymptotic behavior can be improved towards linear in a systematic way. Using the shadow Hamiltonian, one would expect a speedup in the order of $N^{5/4}/N^{\frac{2p+1}{2p}}=N^{\frac{p-2}{4p}}$. As $k$ grows this method will gain a speedup of $N^{1/4}$. Because the shadow Hamiltonian method uses $k$ values of forces and energies which are already available during the integration of equation of motion, the additional computational cost is low and $c$ is only slightly larger than unity. There is an additional memory requirement to form the shadow Hamiltonian.

For $p=8$, we would expect a speedup over plain HMC of $N^{3/16}$. For $N=10^5$ the speedup would be $\approx9$. More importantly, the speedup will grow with $N$, which makes this method a scalable algorithm. For $N=10^{6}$ the speedup with an eigth order accurate shadow Hamiltonian will be $\approx13.$ Note that the use of the shadow Hamiltonian introduces a bias that has to be removed in the acceptance rule of Monte Carlo [69].