Calculating Thermodynamics Properties of Quantum Systems by a non-Markovian Monte Carlo Procedure

We present a history-dependent Monte Carlo scheme for the efficient calculation of the free-energy of quantum systems, inspired by the Wang-Landau sampling and metadynamics method. When embedded in a path integral formulation, it is of general applicability to a large variety of Hamiltonians. In the two-dimensional quantum Ising model, chosen here for illustration, the accuracy of free energy, critical temperature, and specific heat is demonstrated as a function of simulation time, and successfully compared with the best available approaches, particularly the Wang-Landau method over two different Monte Carlo procedures.

💡 Research Summary

The paper introduces a novel history‑dependent Monte Carlo algorithm for calculating the free‑energy of quantum many‑body systems. It merges two powerful ideas that have shaped modern sampling techniques: the Wang‑Landau (WL) flat‑histogram scheme and the collective‑variable biasing of metadynamics (MD). By embedding this hybrid approach within a path‑integral Monte Carlo (PIMC) framework, the authors obtain a method that is both generally applicable to a wide class of quantum Hamiltonians and capable of delivering the full free‑energy surface as a function of temperature, external fields, or any chosen collective variable.

The algorithm proceeds as follows. First, the quantum partition function is expressed as a discretized imaginary‑time path integral using the Trotter decomposition. The configuration space therefore consists of a set of spin (or boson) world‑lines defined on a lattice of spatial sites and imaginary‑time slices. A collective variable ξ – for example the total magnetisation, the energy, or any observable of interest – is defined on each configuration. A bias potential V(ξ) is introduced, initially zero, and is updated during the simulation in the spirit of metadynamics: each time the system visits a particular ξ‑bin, a small increment ΔV is added to V(ξ) for that bin. Simultaneously, a WL‑type density‑of‑states factor g(ξ) is maintained on a logarithmic scale. When the current ξ is visited, the algorithm adds ln f to g(ξ) and increments a histogram H(ξ). Once H(ξ) satisfies a predefined flatness criterion (e.g., 80 % of the average count), the modification factor f is reduced (commonly f←√f) and the histogram is reset. The bias potential V(ξ) is also smoothed, typically by convolution with a Gaussian kernel, to avoid abrupt jumps that could destabilise the Markov chain.

The combined update rule ensures that the sampling probability is proportional to exp