Stochastic many-body perturbation theory for high-order calculations

High-order perturbative $\textit{ab initio}$ calculations are challenging due to the rapidly growing configuration space and the difficulty of assessing convergence. In this work, we introduce perturbation theory quantum Monte Carlo (PTQMC), a stochastic approach designed to compute high-order many-body perturbative corrections. By representing the perturbative wave function with random walkers in configuration space, PTQMC avoids the exponential scaling inherent to conventional constructions of high-rank excitation operators. Benchmark calculations for the Richardson pairing model demonstrate that PTQMC accurately reproduces exact many-body perturbation theory (MBPT) coefficients up to 16th order, even in strongly divergent regimes. We further show that combining PTQMC with series resummation techniques yields stable and precise energy estimates in cases where the straightforward perturbative series fails. Finally, we propose the effective number of configurations, $e^{S}$, as a global measure of perturbative wave-function complexity that can be directly extracted within PTQMC. We demonstrate that the saturation behavior of $e^{S}$ provides a more reliable indicator of the validity of perturbative expansions than energy convergence alone.

💡 Research Summary

The paper introduces a novel stochastic algorithm called perturbation theory quantum Monte Carlo (PTQMC) designed to overcome two long‑standing obstacles in many‑body perturbation theory (MBPT): the exponential growth of the configuration space with perturbative order and the difficulty of assessing convergence of high‑order series. The key idea is to represent the $n$‑th order MBPT wave‑function as a population of signed random walkers distributed over many‑body configurations. Walkers propagate from one configuration to another with a probability proportional to the magnitude of the Hamiltonian matrix element divided by the Møller–Plesset energy denominator, while the sign of the spawned walker follows the sign of the matrix element. This stochastic spawning implements the first term of the exact MBPT recursion relation, while the second term (involving lower‑order energy corrections) is added deterministically. After each propagation step, walkers occupying the same configuration are summed, allowing positive and negative contributions to annihilate, analogous to the annihilation step in full configuration‑interaction QMC (FCIQMC). The resulting signed population $w^{(n)}_I$ provides an unbiased estimator of the exact perturbative coefficient $c^{(n)}_I$, and the $n$‑th order energy correction is obtained as $E^{(n)}=\sum_I \langle\Phi_0|H_1|\Phi_I\rangle w^{(n-1)}_I$.

Because the algorithm samples directly the MBPT recursion, its computational cost scales as $O(n N_w)$, where $n$ is the perturbative order and $N_w$ the number of walkers, a dramatic reduction compared with the conventional $O(N^{2n})$ scaling of diagrammatic or operator‑based approaches. The authors benchmark PTQMC on the Richardson pairing Hamiltonian, a four‑level, half‑filled model that admits exact MBPT coefficients up to at least 16th order. Using a modest walker population ($N_w=10^4$), PTQMC reproduces the exact eighth‑order correlation energy at $g=-1.2$ within statistical error, and systematic convergence of the error with $N_w^{-1/2}$ is demonstrated. Extending to strongly divergent coupling regimes ($-1.2<g<-0.76$), PTQMC remains stable and yields unbiased estimates of the MBPT coefficients up to 16th order, even though the raw perturbative series exhibits large oscillations and fails to converge in the traditional sense.

Recognizing that a divergent raw series cannot be used directly for physical predictions, the authors combine the high‑order PTQMC data with Padé resummation. By constructing rational approximants $