Combining Harmonic Sampling with the Worm Algorithm to Improve the Efficiency of Path Integral Monte Carlo

TL;DR

The paper tackles the low acceptance and slow sampling in Path Integral Monte Carlo (PIMC) for solids and dense systems by splitting the potential into harmonic and anharmonic parts. It introduces Harmonic PIMC (H-PIMC), which samples harmonic paths exactly and accepts based on the anharmonic residue, and its generalization Mixed PIMC (M-PIMC), which confines harmonic sampling to a local minimum vicinity and uses standard PIMC elsewhere. Benchmark results show H-PIMC dramatically improves acceptance and autocorrelation times for weakly to moderately anharmonic systems, accelerates energy convergence with fewer beads, and, for strong anharmonicity, M-PIMC optimizes autocorrelation by tuning the harmonic domain; M-PIMC-PBC extends these gains to periodic systems. The methods are further combined with the worm algorithm to achieve comparable efficiency gains for indistinguishable particles, suggesting substantial speedups for simulations of quantum solids, confined liquids, and complex periodic systems.

Abstract

We propose an improved Path Integral Monte Carlo (PIMC) algorithm called Harmonic PIMC (H-PIMC) and its generalization, Mixed PIMC (M-PIMC). PIMC is a powerful tool for studying quantum condensed phases. However, it often suffers from a low acceptance ratio for solids and dense confined liquids. We develop two sampling schemes especially suited for such problems by dividing the potential into its harmonic and anharmonic contributions. In H-PIMC, we generate the imaginary time paths for the harmonic part of the potential exactly and accept or reject it based on the anharmonic part. In M-PIMC, we restrict the harmonic sampling to the vicinity of local minimum and use standard PIMC otherwise, to optimize efficiency. We benchmark H-PIMC on systems with increasing anharmonicity, improving the acceptance ratio and lowering the auto-correlation time. For weakly to moderately anharmonic systems, at $β\hbar ω=16$, H-PIMC improves the acceptance ratio by a factor of 6-16 and reduces the autocorrelation time by a factor of 7-30. We also find that the method requires a smaller number of imaginary time slices for convergence, which leads to another two- to four-fold acceleration. For strongly anharmonic systems, M-PIMC converges with a similar number of imaginary time slices as standard PIMC, but allows the optimization of the auto-correlation time. We extend M-PIMC to periodic systems and apply it to a sinusoidal potential. Finally, we combine H- and M-PIMC with the worm algorithm, allowing us to obtain similar efficiency gains for systems of indistinguishable particles.

Figures (8)

• Figure 1: A schematic of the H-PIMC (upper panel) and M-PIMC (lower panel) algorithms. The full potential and its harmonic approximation are shown as solid and dashed lines, respectively (left). The shaded region indicates the harmonic domain. In H-PIMC, all the bead locations are proposed using a harmonic update. In M-PIMC, for a given system configuration, some beads fall within the harmonic domain (shaded region). For these, trial positions are sampled using the harmonic update, as in H-PIMC. For the remaining beads, trial positions are sampled using the free particle density matrix, as in PIMC.
• Figure 2: Comparison between PIMC and H-PIMC for harmonic potential. The energy convergence with the number of beads is shown for both PIMC and H-PIMC at $\beta \hbar \omega = 8$. The shaded region indicates a 0.5$\%$ deviation from the exact energy. The acceptance ratio and the energy autocorrelation time data are obtained using the number of beads required for energy convergence at each temperature. For example, at $\beta \hbar \omega =8$, the number of beads used in PIMC and H-PIMC are 48 and 8, respectively. Details of the parameters are given in Appendix \ref{['app:Params']}.
• Figure 3: A comparison between PIMC and H-PIMC results is shown for a weakly anharmonic system (upper panel), a moderately anharmonic system (middle panel), and a strongly anharmonic system (lower panel). The anharmonic potential and its harmonic approximation are shown in solid and dashed lines, respectively (first column). The average energy convergence is shown in the second column. The shaded region indicates a $0.5 \%$ deviation from the exact energy. The acceptance ratio (third column) and the energy autocorrelation time data (fourth column) are obtained using the number of beads required for energy convergence at each temperature. For the strongly anharmonic case (lower panel), the M-PIMC results for the optimized harmonic domain are shown as green squares.
• Figure 4: A comparison between PIMC and H-PIMC results for N = 2 indistinguishable particles is shown for a weakly anharmonic system (upper panel), a moderately anharmonic system (middle panel), and a strongly anharmonic system (lower panel). The anharmonic potential and its harmonic approximation are shown in solid and dashed lines, respectively (first column). The average energy convergence is shown in the second column. The shaded region indicates a $0.5 \%$ deviation from the exact energy. The acceptance ratio (third column) and the energy autocorrelation time data (fourth column) are obtained using the number of beads required for energy convergence at each temperature.
• Figure 5: M-PIMC results for strong anharmonicity at $\beta \hbar \omega = 16$. The M-PIMC calculations use 96 beads. The speedup in acceptance ratio is defined as $\frac{\text{M-PIMC value}}{\text{PIMC value}}$ and the energy autocorrelation time ratio is defined as $\frac{\text{M-PIMC value}}{\text{PIMC value}}$. The shaded area indicates the optimal harmonic domain.