Path Integral Monte Carlo Method Research Articles

On accurate calculations of equilibrium H/D fractionations at super-cold conditions (≤200 K) I: Full Partition Function Ratio (FPFR) vs. Path Integral Monte Carlo (PIMC)

Accurate estimations of Hydrogen/Deuterium (H/D) isotope effects have been a long-standing problem since the beginning of stable isotope geochemistry, primarily due to their strong nuclear quantum effects (NQEs). Currently, the D/H ratio plays a key role in deciphering the origin and evolution of water ice and organics in the solar system. Additionally, the Big-bang origin of Deuterium necessitates that the observed D/H variations be linked to isotope exchange reactions, such as HD + H2O ↔ H2 + HDO. Therefore, an accurate benchmark dataset describing the equilibrium H/D isotope fractionations at low temperatures (≤ 200 K) is essential for interpreting rapidly accumulating astrochemical observation data. In this study, we implemented and compared two methods developed for evaluating the NQEs of H/D isotopes: 1) the Full Partition Function Ratio (FPFR) method and 2) the Path Integral Monte Carlo (PIMC) method. Both were assessed for their accuracy and efficiency in estimating equilibrium isotope effects (EIEs) using 10 common interstellar small molecules (H2, HF, HCl, HCN, HNC, H2O, H2S, NH3, HCHO, and CH4) as examples. Most of them are calculated at high precision for the first time. For the FPFR method, several higher-order corrections to the Urey-Bigeleisen-Mayer model were considered, including anharmonicity, quantum mechanical rotation with nuclear-spin statistics, and corrections to the Born-Oppenheimer approximation. For the PIMC method, four different actions related to the thermal density matrix are used. These include the primitive and three fourth-order actions (Takahashi and Imada, Suzuki-Chin and Chin approximations), based on the direct scaled-coordinates transformation. Three main conclusions are obtained: 1) Under the framework of Born-Oppenheimer approximation, the PIMC method holds the same and even slightly better accuracy (more comprehensive descriptions of quantum anharmonicities) in estimating EIEs for H/D isotopes as the FPFR method. However, both methods suffer from the neglection of the DBOC when using the Born-Oppenheimer potential energy surfaces. For the PIMC method, additional corrections for particle distinguishability must be considered when nuclear spin statistics exhibit its effects at low temperatures (e.g., < 160 K for H2). 2) Among the four actions, three fourth-order actions show better convergence and efficiency (2.2 to 26.5 times speedups at 50 K) than the primitive action, making them superior candidates for PIMC simulations at super-cold (≤ 200 K) conditions. 3) A combination of accurate potential energy surface with the PIMC method serves as an alternative way for theoretical estimations of H/D isotope effects.

Bounds to electron spin qubit variability for scalable CMOS architectures

Spins of electrons in silicon MOS quantum dots combine exquisite quantum properties and scalable fabrication. In the age of quantum technology, however, the metrics that crowned Si/SiO2 as the microelectronics standard need to be reassessed with respect to their impact upon qubit performance. We chart spin qubit variability due to the unavoidable atomic-scale roughness of the Si/SiO2 interface, compiling experiments across 12 devices, and develop theoretical tools to analyse these results. Atomistic tight binding and path integral Monte Carlo methods are adapted to describe fluctuations in devices with millions of atoms by directly analysing their wavefunctions and electron paths instead of their energy spectra. We correlate the effect of roughness with the variability in qubit position, deformation, valley splitting, valley phase, spin-orbit coupling and exchange coupling. These variabilities are found to be bounded, and they lie within the tolerances for scalable architectures for quantum computing as long as robust control methods are incorporated.

Effect of zero-point motion on properties of quantum particles adsorbed on a substrate

We qualitatively investigate the effect of zero-point motion (ZPM) on the structure and properties of a film composed of quantum particles adsorbed on a graphite substrate. The amplitude of ZPM is controlled by a change of the particle mass while keeping the interactions fixed. In that sense it is assumed that the interactions can be controlled by future doping methods. The worm-algorithm path integral Monte Carlo (WAPIMC) method is applied to simulate this system in the grand-canonical ensemble, where particles can be exchanged with the external particle reservoir. Another method, namely the multiconfigurational time-dependent Hartree method for bosons is additionally applied to verify some of the WAPIMC results and to provide further information on the entropy and the condensate fraction. Several important findings are reported. It is found that ZPM plays an important role in defining order and disorder in the crystalline structure of the adsorbed film. The total energy of the film drops with a reduction in the amplitude of ZPM, that is, it becomes more negative which is an indication to stronger adsorption. For a few particle numbers, a significant condensate fraction is detected that however drops sharply at critical values of the ZPM amplitude. Most importantly, a connection is established between chaos, in coordinate as well as momentum space, and the Heisenberg uncertainty principle. The importance of the present study lies in the fact that adsorbed two-dimensional films serve as an excellent experimental testbed for demonstrating low-dimensional quantum phenomena in the ground state. The present examination contributes also to a further understanding of the properties of heavy quantum particles adsorbed on substrates.

Fermionic physics from abinitio path integral Monte Carlo simulations of fictitious identical particles.

The abinitio path integral Monte Carlo (PIMC) method is one of the most successful methods in statistical physics, quantum chemistry and related fields, but its application to quantum degenerate Fermi systems is severely hampered by an exponential computational bottleneck: the notorious fermion signproblem. Very recently, Xiong and Xiong [J. Chem. Phys. 157, 094112 (2022)] have suggested to partially circumvent the signproblem by carrying out simulations of fictitious systems guided by an interpolating continuous variable ξ ∈ [-1, 1], with the physical Fermi- and Bose-statistics corresponding to ξ = -1 and ξ = 1. It has been proposed that information about the fermionic limit might be obtained by calculations within the bosonic sector ξ > 0 combined with an extrapolation throughout the fermionic sector ξ < 0, essentially bypassing the signproblem. Here, we show how the inclusion of the artificial parameter ξ can be interpreted as an effective penalty on the formation of permutation cycles in the PIMC simulation. We demonstrate that the proposed extrapolation method breaks down for moderate to high quantum degeneracy. Instead, the method constitutes a valuable tool for the description of large Fermi-systems of weak quantum degeneracy. This is demonstrated for electrons in a 2D harmonic trap and for the uniform electron gas (UEG), where we find excellent agreement (∼0.5%) with exact configuration PIMC results in the high-density regime while attaining a speed-up exceeding 11 orders of magnitude. Finally, we extend the idea beyond the energy and analyze the radial density distribution (2D trap), as well as the static structure factor and imaginary-time density-density correlation function (UEG).

Correlation properties of a one-dimensional repulsive Bose gas at finite temperature

We present a comprehensive study shedding light on how thermal fluctuations affect correlations in a Bose gas with contact repulsive interactions in one spatial dimension. The pair correlation function, the static structure factor, and the one-body density matrix are calculated as a function of the interaction strength and temperature with the exact ab-initio Path Integral Monte Carlo method. We explore all possible gas regimes from weak to strong interactions and from low to high temperatures. We provide a detailed comparison with a number of theories, such as perturbative (Bogoliubov and decoherent classical), effective (Luttinger liquid) and exact (ground-state and thermal Bethe Ansatz) ones. Our Monte Carlo results exhibit an excellent agreement with the tractable limits and provide a fundamental benchmark for future observations which can be achieved in atomic gases, cavity quantum-electrodynamic and superconducting-circuit platforms.

Unusual fluid of one-dimensional disordered bosons at finite temperature

In previous important examinations by Aleiner et al. [Nat. Phys. 6, 900 (2010)] and Michal et al. [Proc. Natl. Acad. Sci. (U.S.A.) 113, E4455 (2016)], the fluid insulator transition has been explored in one-dimensional (1D) disordered bosons by analytical methods. However, the nature of this ``fluid'' is still not quite understood. Here we shed more light on this ``fluid'' by studying its properties under changes in the average particle number $\ensuremath{\langle}N\ensuremath{\rangle}$ and temperature $T$ in the absence of superfluidity. For this purpose, the worm algorithm path-integral Monte Carlo method is utilized in the grand-canonical ensemble. A measure for the width of path fluctuations is introduced that we call ${\ensuremath{\sigma}}_{J}$. It is found that by increasing $\ensuremath{\langle}N\ensuremath{\rangle}$, ${\ensuremath{\sigma}}_{J}$ rises as the bosons become less repulsive and gain more mobility. The compressibility $\ensuremath{\chi}(\ensuremath{\mu})$, $\ensuremath{\mu}$ being the chemical potential, reveals fluid behavior at the lower $\ensuremath{\gamma}$ and a proposed Bose-glass phase at the larger $\ensuremath{\gamma}$ (the Lieb-Liniger interaction parameter). The kinetic energy ${E}_{K}$ versus $\ensuremath{\gamma}$ at various $T$ signals a transition by its convergence towards a common value. ${E}_{K}$ versus $\ensuremath{\mu}$ reveals a peculiar inversion as one goes from small to larger values of the speckle strength $\ensuremath{\langle}V\ensuremath{\rangle}$ demonstrating a complex influence of the disorder. A tipover of ${E}_{K}$ versus $\ensuremath{\gamma}$ from an increasing to a decreasing trend at the lower $\ensuremath{\gamma}$ also signals a transition. Further, Luttinger's liquid theory cannot be applied to describe 1D bosons in strong disorder at finite $T$ because of the absence of superfluidity. Moreover, $\ensuremath{\langle}V\ensuremath{\rangle}$ largely controls the behavior of $\ensuremath{\chi}(\ensuremath{\mu})$, the pair-correlation function at the origin ${g}_{2}(0)$, and ${E}_{K}$. Surprisingly, the diffusion constant ${D}_{0}$ does not obey the classical isomorphism as it rises with $T$. We also emphasize the importance of this work via its applicability to quantum wires.

A path integral Monte Carlo (PIMC) method based on Feynman-Kac formula for electrical impedance tomography

A path integral Monte Carlo method (PIMC) based on a Feynman-Kac formula for the Laplace equation with mixed boundary conditions is proposed to solve the forward problem of the electrical impedance tomography (EIT). The forward problem is an important part of iterative algorithms of the inverse EIT problem, and the proposed PIMC provides a local solution to find the potentials and currents on individual electrodes. Improved techniques are proposed to compute with better accuracy both the local time of reflecting Brownian motions (RBMs) and the Feynman-Kac formula for mixed boundary problems of the Laplace equation. Accurate voltage-to-current maps on the electrodes of a model 3-D EIT problem with eight electrodes are obtained by solving a mixed boundary problem with the proposed PIMC method.

Swarm-like ground-state searching based on a path-integral quantum Monte Carlo method for automatically fluctuation-controlled simulated annealing

The path-integral quantum Monte Carlo (PIQMC) method is widely used as a classical simulation algorithm for quantum annealing. Replicas represent different points in imaginary time. We propose and demonstrate a swarm-like ground-state searching algorithm based on the PIQMC method for automatically fluctuation-controlled simulated annealing. Replicas help one another to search for the ground state of an Ising system as if forming a swarm and working cooperatively. Their interactions are local, simple and fluctuate to a certain degree. Such fluctuations are necessary to escape local minima, but the fluctuations do not need to be explicitly controlled. The size of fluctuations is automatically adjusted as annealing proceeds. We solve max-cut problems for algorithm evaluation, each of which corresponds to a graph of 100 vertices. We also solve the same problems using a conventional method based on the Metropolis algorithm for comparison.

Accounting for long–range interaction in the Kelbg pseudopotential

AbstractAn expression for the high–temperature density matrix taking into account long–range interaction effects is obtained. The angular averaged Ewald (AAE) potential and Kelbg's solution of the Bloch equation in the first order of perturbation theory are used for this purpose. The expression generalizes the high–temperature density matrix for the Coulomb potential derived earlier by G. Kelbg. An analytical formula for the diagonal Kelbg‐AAE pseudopotential is obtained and analysed. The results can be useful in path integral Monte Carlo methods.

Path-integral quantum Monte Carlo calculations of light nuclei

We describe a path-integral ground-state quantum Monte Carlo method for light nuclei in continuous space. We show how to efficiently update and sample the paths with spin-isospin dependent and spin-orbit interactions. We apply the method to the triton and alpha particle using both local chiral interactions with next-to-next-to-leading-order %(N$^2$LO) and the Argonne interactions. For operators, like the total energy, that commute with the Hamiltonian, our results agree with Green's function Monte Carlo and auxiliary field diffusion Monte Carlo calculations. For operators that do not commute with the Hamiltonian and for Euclidean response functions, the path-integral formulation allows straightforward calculation without forward walking or the increased variance typical of diffusion methods. We demonstrate this by calculating density distributions, root mean square radii, and Euclidean response functions for single-nucleon couplings.

Quantum computation of Restricted Boltzmann Machines by Monte Carlo Methods

In recent years, the diversification of problems that require computers to solve has attracted attention to the construction of meta-heuristics that can be applied to a wide range of problems, and to specialized computers that implement these meta-heuristics in their devices. The representative meta-heuristics are Simulated Annealing (SA) and its extension to quantum computation, Quantum Annealing (QA), and its path-integral Monte Carlo method for classical simulation Crosson and Harrow showed that for certain problems where QA outperformed SA, SQA achieved performance close to that of QA, and SQA sometimes outperformed SA by an exponential time factor. On the other hand, it remains unclear whether SQA can work efficiently on a wide range of other problems. In this study, we experimentally compared SA and SQA on instances of the restricted Boltzmann machine RBM, known as a fundamental building block in deep learning, and 3SAT, a fundamental combinatorial optimization problem. The results show that SQA gives slightly better solutions than SA as the problem size increases for RBM in terms of both accuracy and computation time in our setting, but the opposite trend is observed for 3SAT, indicating that there is no significant difference between the two methods. From the viewpoint of artificial intelligence research, it is necessary to further examine whether deep learning can be made more efficient by applying QA and SQA to RBM.

Hole-induced anomaly in the thermodynamic behavior of a one-dimensional Bose gas

We reveal an intriguing anomaly in the temperature dependence of the specific heat of a one-dimensional Bose gas.~The observed peak holds for arbitrary interaction and remembers a superfluid-to-normal phase transition in higher dimensions, but phase transitions are not allowed in one dimension.~The presence of the anomaly signals a region of unpopulated states which behaves as an energy gap and is located below the hole branch in the excitation spectrum.~The anomaly temperature is found to be of the same order of the energy of the maximum of the hole branch.~We rely on the Bethe Ansatz to obtain the specific heat exactly and provide interpretations of the analytically tractable limits.~The dynamic structure factor is computed with the Path Integral Monte Carlo method for the first time.~We notice that at temperatures similar to the anomaly threshold, the energy of the thermal fluctuations become comparable with the maximal hole energy, leading to a qualitative change in the structure of excitations.~This excitation pattern experiences the breakdown of the quasiparticle description for any value of the interaction strength at the anomaly, similarly to any superfluid phase transition at the critical temperature.~We provide indications for future observations and how the hole anomaly can be employed for in-situ thermometry, identifying different collisional regimes and understanding other anomalies in atomic, solid-state, electronic, spin-chain and ladder systems.

Momentum Distribution Functions and Pair Correlation Functions of Unpolarized Uniform Electron Gas in Warm Dense Matter Regime

In this paper we continued our research of the uniform electron gas in a warm dense matter regime, focusing on the momentum distribution functions and pair correlation functions. We use the single–momentum path integral Monte Carlo method, based on the Wigner formulation of quantum statistics to calculate both momentum- and coordinate-depending distributions and average values of quantum operators for many-fermion Coulomb systems. We discovered that the single-particle momentum distribution function deviates from the ideal Fermi distribution and forms the so-called “quantum tails” at high momenta, if non-ideality is strong enough in both degenerate and non-degenerate cases. This effect is always followed by the appearance of the short-range order on pair correlation functions and can be explained by the tunneling through the effective potential wells surrounding the electrons. Furthermore, we calculated the average kinetic and potential energies in the wide range of states, expanding our previous results significantly.

The density–temperature range of exchange–correlation exciton existence by the fermionic path integral Monte Carlo method

A recently developed fermionic path integral Monte Carlo approach [Filinov et al., Phys. Rev. E 102, 033203 (2020) and Filinov et al., J. Phys. A 55, 035001 (2021)] has been applied for the estimation of the density–temperature range of exchange–correlation exciton existence in a strongly coupled degenerate uniform electron gas. The approach allows us to reduce the “fermionic sign problem” taking into account the interference effects of the Coulomb and exchange interaction of electrons in the basic Monte Carlo cell and its periodic images. Our results for radial distribution functions demonstrate the formation and decay of a short-range quantum ordering of electrons associated with exchange–correlation excitons in the literature. Such excitons have never been observed earlier in standard path integral Monte Carlo simulations.

Density Functional Theory Perspective on the Nonlinear Response of Correlated Electrons across Temperature Regimes.

We explore a new formalism to study the nonlinear electronic density response based on Kohn–Sham density functional theory (KS-DFT) at partially and strongly quantum degenerate regimes. It is demonstrated that the KS-DFT calculations are able to accurately reproduce the available path integral Monte Carlo simulation results at temperatures relevant for warm dense matter research. The existing analytical results for the quadratic and cubic response functions are rigorously tested. It is demonstrated that the analytical results for the quadratic response function closely agree with the KS-DFT data. Furthermore, the performed analysis reveals that currently available analytical formulas for the cubic response function are not able to describe simulation results, neither qualitatively nor quantitatively, at small wavenumbers q < 2qF, with qF being the Fermi wavenumber. The results show that KS-DFT can be used to describe warm dense matter that is strongly perturbed by an external field with remarkable accuracy. Furthermore, it is demonstrated that KS-DFT constitutes a valuable tool to guide the development of the nonlinear response theory of correlated quantum electrons from ambient to extreme conditions. This opens up new avenues to study nonlinear effects in a gamut of different contexts at conditions that cannot be accessed with previously used path integral Monte Carlo methods.

Zonal Estimators for Quasiperiodic Bosonic Many-Body Phases.

In this work, we explore the relevant methodology for the investigation of interacting systems with contact interactions, and we introduce a class of zonal estimators for path-integral Monte Carlo methods, designed to provide physical information about limited regions of inhomogeneous systems. We demonstrate the usefulness of zonal estimators by their application to a system of trapped bosons in a quasiperiodic potential in two dimensions, focusing on finite temperature properties across a wide range of values of the potential. Finally, we comment on the generalization of such estimators to local fluctuations of the particle numbers and to magnetic ordering in multi-component systems, spin systems, and systems with nonlocal interactions.

Bound States of the Exchange—Correlation Excitons in the Uniform Electron Gas by the Monte Carlo Simulations

The modified path integral representation of Wigner functions and the new Monte Carlo approach has been suggested to account for the impact of the interparticle interaction on the Pauli exclusion principle of fermions. This approach also allows to calculate the momentum distribution functions and to reduce the “sign problem” that is inaccessible to the standard path integral Monte Carlo methods. The obtained pair electron–electron distribution functions for the “uniform electron gas” demonstrate the short-range quantum ordering of electrons associated with exchange-correlation excitons. The exchange-correlation exciton is caused by the interaction of electrons with positively charged exchange holes and the excluded volume effect. The developed approach allows one to study the density–temperature range of the exciton arising, existence, and decay. Using the potential of the mean force and semiclassical Bohr–Sommerfeld quantization condition, we have demonstrated the existence of bound states disturbing the Maxwellian distribution and estimated their average energy levels. The exchange-correlation excitons have not been observed earlier in the standard path integral Monte Carlo simulations.

Negative string tension of a higher-charge Schwinger model via digital quantum simulation

Abstract We study some properties of generalized global symmetry for the charge-q Schwinger model in the Hamiltonian formalism, which is the (1 + 1)D quantum electrodynamics with a charge-q Dirac fermion. This model has the $\mathbb {Z}_q\, 1$-form symmetry, which is a remnant of the electric $U(1)\, 1$-form symmetry in the pure Maxwell theory. It is known that, if we put the theory on closed space, then the Hilbert space is decomposed into q distinct sectors, called universes, and some states with higher energy density do not decay to the ground state due to the selection rule of the 1-form symmetry. Even with open boundaries, we can observe the stability of such states by seeing a negative string tension behavior, meaning that opposite charges repel each other. In order to see negative string tensions, the vacuum angle θ has to be large enough and the standard path-integral Monte Carlo method suffers from the sign problem. We develop a method based on the adiabatic state preparation to see this feature with digital quantum simulation and confirm it using a classical simulator of quantum devices. In particular, we measure the local energy density and see how it jumps between the inside and outside of the insertion of the probe charges. We explicitly see that the energy density inside is lower than that outside. This is a clear signature of the negative string tension.

Path-integral calculation of the third dielectric virial coefficient of noble gases.

We present a rigorous framework for fully quantum calculation of the third dielectric virial coefficient Cɛ(T) of noble gases, including exchange effects. The quantum effects are taken into account with the path-integral Monte Carlo method. Calculations employing state-of-the-art pair and three-body potentials and pair polarizabilities yield results generally consistent with the few scattered experimental data available for helium, neon, and argon, but rigorous calculations with well-described uncertainties will require the development of surfaces for the three-body nonadditive polarizability and the three-body dipole moment. The framework, developed here for the first time, will enable new approaches to primary temperature and pressure metrology based on first-principles calculations of gas properties.

Single-momentum path integral Monte Carlo simulations of uniform electron gas in warm dense matter regime

In this paper, the single-momentum path integral Monte Carlo method, previously developed for simple quantum systems and hydrogen plasma, is adapted to simulations of the uniform electron gas. The developed method is based on the combination of Wigner formalism and the path integral approach and is able to calculate various thermodynamic values and distribution functions without differentiation of the partition function. Since the exchange interaction between electrons is taken into account by the Gram determinants of the exchange matrix, the fermionic sign problem is reduced significantly, and in the case of coordinate-depending variables, is completely eliminated. The method was applied to study thermodynamic properties of the uniform electron gas in warm dense matter regime. Average kinetic, potential, and exchange-correlation energy were calculated in a wide range of states.