Exchange Monte Carlo for continuous-space Path Integral Monte Carlo simulation

TL;DR

This work tackles the slow sampling of permutation sectors in finite-temperature bosonic PIMC by introducing Exchange Monte Carlo (EMC) updates that operate along the axis of interacting imaginary-time slices rather than temperature. Together with Stochastic Potential Switching (SPS) and a No-U-Turn Sampler (NUTS)–based local update strategy, the method accelerates exploration of winding-number sectors while maintaining correct thermodynamic statistics, reducing cost for long-range potentials. A nonuniform imaginary-time discretization is devised to flatten exchange rates and enable more frequent replica round-trips, with bond-refresh steps ensuring detailed balance in the expanded ensemble. Numerical results for $^4$He demonstrate reduced autocorrelations, accurate energy and observable estimates, and a finite-size scaling analysis yielding $T_c$ consistent with established worm algorithm results, highlighting the practical impact for scalable quantum simulations of bosonic systems.

Abstract

We present a novel Exchange Monte Carlo (EMC) method designed for application in continuous-space Path Integral Monte Carlo (PIMC) simulations at finite temperature. Traditional PIMC methods for bosonic systems suffer from long autocorrelation times, particularly when measuring observables affected by particle permutations, such as the winding number. To address this issue, we introduce an exchange update scheme that facilitates replica transitions between different interaction regimes, significantly accelerating Monte Carlo dynamics-especially for global observables sensitive to permutation effects. Furthermore, we incorporate Stochastic Potential Switching (SPS) to efficiently decompose interactions, substantially enhancing computational efficiency for long-range interatomic pair potentials such as the Lennard-Jones and Aziz potentials.

Figures (11)

• Figure 1: Stochastic potential switching treatment of particle interactions in two dimensions. The space is divided into bins of size $r_c$. For interactions of particle A, we only need to consider the nearest bins (orange region) and interactions that are switched to $\bar{U}$ (particle B and C).
• Figure 2: Autocorrelation time of observables vs. number of interacting imaginary slices $p'$. Autocorrelation time quantifies how long it takes for an observable in a time series to become statistically independent, and it is calculated using binning. Superfluid fraction, which is affected by the configuration's topology, can be sampled more efficiently during simulations for partially interacting systems.
• Figure 3: Exchange update between replica index $j$ and $j-1$, indicating the number of interacting slices in the configuration. If the proposal is accepted, replica 1(left replica) will turn on one more interacting slice, and replica 2(right replica) will turn off one more interacting slice. Curved lines represent an abbreviated path to the end of the worldlines.
• Figure 4: Trajectory of replica at system size $N = 128$ as it traverses between fully interacting system($p' = P = 94$) and partially interacting system ($p' = 76$). Measurements of superfluid fraction and longest cycle length are compared with the MCMC simulation, which has only one replica staying at $p' = 94$.
• Figure 5: Exchange rate at different $p'$ for uniform and nonuniform distribution of imaginary time slices. After modeling the slice distribution, the exchange rates are flattened, and the smallest exchange rate is increased.