Counterdiabatic Hamiltonian Monte Carlo

Generating samples from a continuous probability density is a central algorithmic problem across statistics, engineering, and the sciences. For high-dimensional settings, Hamiltonian Monte Carlo (HMC) is the default algorithm across mainstream software packages. However, despite the extensive line of work on HMC and its widespread empirical success, it remains unclear how many iterations of HMC are required as a function of the dimension $d$. On one hand, a variety of results show that Metropolized HMC converges in $O(d^{1/4})$ iterations from a warm start close to stationarity. On the other hand, Metropolized HMC is significantly slower without a warm start, e.g., requiring $Ω(d^{1/2})$ iterations even for simple target distributions such as isotropic Gaussians. Finding a warm start is therefore the computational bottleneck for HMC. We resolve this issue for the well-studied setting of sampling from a probability distribution satisfying strong log-concavity (or isoperimetry) and third-order derivative bounds. We prove that non-Metropolized HMC generates a warm start in $\tilde{O}(d^{1/4})$ iterations, after which we can exploit the warm start using Metropolized HMC. Our final complexity of $\tilde{O}(d^{1/4})$ is the fastest algorithm for high-accuracy sampling under these assumptions, improving over the prior best of $\tilde{O}(d^{1/2})$. This closes the long line of work on the dimensional complexity of MHMC for such settings, and also provides a simple warm-start prescription for practical implementations.

Learned dynamical proposals can generate configurations without providing a tractable endpoint density. Nonequilibrium path probabilities offer a way to correct such proposals, but correction alone does not determine their ability to connect separated regions. We introduce Neural Non-Equilibrium Hamiltonian Monte Carlo (NHMC), which combines conditional momentum distributions with reversible, volume preserving dynamics. The forward--reverse path ratio gives the work used for training, importance weighting, normalizer estimation, and Metropolis correction. We then construct a configuration space round-trip kernel whose reverse and forward paths share an intermediate configuration. Conditional on that configuration, its path record update is independence Metropolis--Hastings. We compare its stationary inter-region flow with the flow obtained under exact conditional matching and bound their difference by the conditional path mismatch. This separates errors in the conditional proposal from dependence between the endpoint region and the intermediate configuration. Many well calculations test normalizers and mode probabilities; a controlled lattice $φ^4$ calculation compares ideal and accepted inter-sector flows and measures sector relaxation. In two dimensional $U(1)$, $\mathrm{SU}(2)$, and $\mathrm{SU}(3)$ gauge theories, corrected path weight estimates reproduce finite volume Wilson loops; $U(1)$ topological susceptibilities are also compared with analytic predictions across several lattice volumes.

We develop a new low-accuracy sampler, called smoothed Picard Hamiltonian Monte Carlo, which combines Gaussian smoothing, Picard iteration, and higher-order discretization. For a log-concave target $π\propto \exp(-V)$ in dimension $d$ satisfying $0 \prec αI \preceq \nabla^2 V \preceq βI$, with condition number $κ:= β/α$, smoothed Picard HMC returns a sample with $\sqrt α\,W_2(\cdot,π) \le \varepsilon$ using $\widetilde O(κ^2 + κ^{7/6} d^{1/6}/\varepsilon^{1/3})$ gradient queries. We also prove stronger $W_q$ bounds, and then develop an algorithmic framework, the recursive warm start generator, to upgrade these $W_q$ bounds to stronger divergence guarantees. This produces a warm start for the proximal bouncy particle sampler, introduced in a companion work, leading to a high-accuracy log-concave sampler with complexity $\widetilde O((κ^{7/6} d^{1/6} + κ^{1/2} d^{1/4})\mathrm{polylog}(1/\varepsilon))$.