Path integral molecular dynamics for fermions: Alleviating the sign problem with the Bogoliubov inequality

Ab Initio Path Integral Monte Carlo Simulations of the Uniform Electron Gas on Large Length Scales

The accurate description of non-ideal quantum many-body systems is of prime importance for a host of applications within physics, quantum chemistry, materials science, and related disciplines. At finite temperatures, the gold standard is given by ab initio path integral Monte Carlo (PIMC) simulations, which do not require any empirical input but exhibit an exponential increase in the required computation time for Fermionic systems with an increase in system size N. Very recently, computing Fermionic properties without this bottleneck based on PIMC simulations of fictitious identical particles has been suggested. In our work, we use this technique to perform very large (N ≤ 1000) PIMC simulations of the warm dense electron gas and demonstrate that it is capable of providing a highly accurate description of the investigated properties, i.e., the static structure factor, the static density response function, and the local field correction, over the entire range of length scales.

Fermionic physics from ab initio path integral Monte Carlo simulations of fictitious identical particles

The ab initio 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 sign problem. Very recently, Xiong and Xiong [J. Chem. Phys. 157, 094112 (2022)] have suggested to partially circumvent the sign problem 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 sign problem. 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).

Quadratic scaling bosonic path integral molecular dynamics

Bosonic exchange symmetry leads to fascinating quantum phenomena, from exciton condensation in quantum materials to the superfluidity of liquid 4He. Unfortunately, path integral molecular dynamics (PIMD) simulations of bosons are computationally prohibitive beyond ∼100 particles, due to a cubic scaling with the system size. We present an algorithm that reduces the complexity from cubic to quadratic, allowing the first simulations of thousands of bosons using PIMD. Our method is orders of magnitude faster, with a speedup that scales linearly with the number of particles and the number of imaginary time slices (beads). Simulations that would have otherwise taken decades can now be done in days. In practice, the new algorithm eliminates most of the added computational cost of including bosonic exchange effects, making them almost as accessible as PIMD simulations of distinguishable particles.

Electronic density response of warm dense hydrogen on the nanoscale

The properties of hydrogen at warm dense matter (WDM) conditions are of high importance for the understanding of astrophysical objects and technological applications such as inertial confinement fusion. In this work, we present extensive ab initio path integral Monte Carlo results for the electronic properties in the Coulomb potential of a fixed ionic configuration. This gives us unique insights into the complex interplay between the electronic localization around the protons with their density response to an external harmonic perturbation. We find qualitative agreement between our simulation data and a heuristic model based on the assumption of a local uniform electron gas model, but important trends are not captured by this simplification. In addition to being interesting in their own right, we are convinced that our results will be of high value for future projects, such as the rigorous benchmarking of approximate theories for the simulation of WDM, most notably density functional theory.

Thermodynamics of fermions at any temperature based on parametrized partition function

In this work, we study the recently developed parametrized partition function formulation and show how we can infer the thermodynamic properties of fermions based on numerical simulation of bosons and distinguishable particles at various temperatures. In particular, we show that in the three-dimensional space defined by energy, temperature, and the parameter characterizing parametrized partition function, we can map the energies of bosons and distinguishable particles to fermionic energies through constant-energy contours. We apply this idea to both noninteracting and interacting Fermi systems and show it is possible to infer the fermionic energies at all temperatures, thus providing a practical and efficient approach to obtain thermodynamic properties of Fermi systems with numerical simulation. As an example, we present energies and heat capacities for 10 noninteracting fermions and 10 interacting fermions and show good agreement with the analytical result for the noninteracting case.

On the thermodynamic properties of fictitious identical particles and the application to fermion sign problem

By generalizing the recently developed path integral molecular dynamics for identical bosons and fermions, we consider the finite-temperature <p>thermodynamic properties of fictitious identical particles with a real parameter ξ interpolating continuously between bosons (ξ=1) and fermions (ξ=-1). Through general analysis and numerical experiments we find that the average energy may have good analytical property as a function of this real parameter ξ, which provides the chance to calculate the thermodynamical properties of identical fermions by an extrapolation with a simple polynomial function after accurately calculating the thermodynamic properties of the fictitious particles for ξ{greater than or equal to}0. Using several examples, it is shown that our method can efficiently give accurate energy values for finite-temperature fermionic systems. Our work provides a chance to circumvent the fermion sign problem for some quantum systems.

Path-integral molecular dynamics for anyons, bosons, and fermions

In this article we develop a general method to numerically calculate physical properties for a system of anyons with path integral molecular dynamics. We provide a unified method to calculate the thermodynamics of identical bosons, fermions, and anyons. Our method is tested and applied to systems of anyons, bosons, and fermions in a two-dimensional harmonic trap. We also consider a method to calculate the energy for fermions as an application of the path integral molecular dynamics to simulate the anyon model.

Path integral molecular dynamics for thermodynamics and Green's function of ultracold spinor bosons

Most recently, the path integral molecular dynamics has been successfully used to consider the thermodynamics of single-component identical bosons and fermions. In this work, the path integral molecular dynamics is developed to simulate the thermodynamics, Green's function and momentum distribution of two-component bosons in three dimensions. As an example of our general method, we consider the thermodynamics of up to sixteen bosons in a three-dimensional harmonic trap. For noninteracting spinor bosons, our simulation shows a bump in the heat capacity. As the repulsive interaction strength increases, however, we find the gradual disappearance of the bump in the heat capacity. We believe this simulation result can be tested by ultracold spinor bosons with optical lattices and magnetic-field Feshbach resonance to tune the inter-particle interaction. We also calculate Green's function and momentum distribution of spinor bosons. Our work facilitates the exact numerical simulation of spinor bosons, whose property is one of the major problems in ultracold Bose gases.

Numerical calculation of Green's function and momentum distribution for spin-polarized fermions by path integral molecular dynamics

Most recently, path integral molecular dynamics (PIMD) has been successfully applied to perform simulations of identical bosons and fermions by B. Hirshberg et al.. In this work, we demonstrate that PIMD can be developed to calculate Green's function and extract momentum distribution for spin-polarized fermions. In particular, we show that the momentum distribution calculated by PIMD has potential application to numerous quantum systems, e.g. ultracold fermionic atoms in optical lattices.

Path integral molecular dynamics simulations for Green's function in a system of identical bosons

Path integral molecular dynamics (PIMD) has been successfully applied to perform simulations of large bosonic systems in a recent study [Hirshberg et al., Proc. Natl. Acad. Sci. U. S. A. 116, 21445 (2019)]. In this work, we extend PIMD techniques to study Green's function for bosonic systems. We demonstrate that the development of the original PIMD method enables us to calculate Green's function and extract momentum distribution from our simulations. We also apply our method to systems of identical interacting bosons to study Berezinskii-Kosterlitz-Thouless transition around its critical temperature.

Prediction of a Supersolid Phase in High-Pressure Deuterium

Supersolid is a mysterious and puzzling state of matter whose possible existence has stirred a vigorous debate among physicists for over 60 years. Its elusive nature stems from the coexistence of two seemingly contradicting properties, long-range order and superfluidity. We report computational evidence of a supersolid phase of deuterium under high pressure (p>800  GPa) and low temperature (T<1.0  K). In our simulations, that are based on bosonic path integral molecular dynamics, we observe a highly concerted exchange of atoms while the system preserves its crystalline order. The exchange processes are favored by the soft core interactions between deuterium atoms that form a densely packed metallic solid. At the zero temperature limit, Bose-Einstein condensation is observed as the permutation probability of N deuterium atoms approaches 1/N with a finite superfluid fraction. Our study provides concrete evidence for the existence of a supersolid phase in high-pressure deuterium and could provide insights on the future investigation of supersolid phases in real materials.

The Sign Problem in Density Matrix Quantum Monte Carlo

Density matrix quantum Monte Carlo (DMQMC) is a recently developed method for stochastically sampling the N-particle thermal density matrix to obtain exact-on-average energies for model and ab initio systems. We report a systematic numerical study of the sign problem in DMQMC based on simulations of atomic and molecular systems. In DMQMC, the density matrix is written in an outer product basis of Slater determinants. In principle, this means that DMQMC needs to sample a space that scales in the system size, N, as O[(exp(N))²]. In practice, removing the sign problem requires a total walker population that exceeds a system-dependent critical walker population (N_c), imposing limitations on both storage and compute time. We establish that N_c for DMQMC is the square of N_c for FCIQMC. By contrast, the minimum N_c in the interaction picture modification of DMQMC (IP-DMQMC) is only linearly related to the N_c for FCIQMC. We find that this difference originates from the difference in propagation of IP-DMQMC versus canonical DMQMC: the former is asymmetric, whereas the latter is symmetric. When an asymmetric mode of propagation is used in DMQMC, there is a much greater stochastic error and is thus prohibitively expensive for DMQMC without the interaction picture adaptation. Finally, we find that the equivalence between IP-DMQMC and FCIQMC seems to extend to the initiator approximation, which is often required to study larger systems with large basis sets. This suggests that IP-DMQMC offers a way to ameliorate the cost of moving between a Slater determinant space and an outer product basis.

Isotopic separation of helium through graphyne membranes: a ring polymer molecular dynamics study

Microscopic-level understanding of the separation mechanism for two-dimensional (2D) membranes is an active area of research due to potential implications of this class of membranes for various technological processes. Helium (He) purification from the natural resources is of particular interest due to the shortfall in its production. In this work, we applied the ring polymer molecular dynamics (RPMD) method to graphdiyne (Gr2) and graphtriyne (Gr3) 2D membranes having variable pore sizes for the separation of He isotopes, and compare for the first time with rigorous quantum calculations. We found that the transmission rate through Gr3 is many orders of magnitude greater than Gr2. The selectivity of either isotope at low temperatures is a consequence of a delicate balance between the zero-point energy effect and tunneling of ⁴He and ³He. In particular, a remarkable tunneling effect is reported on the Gr2 membrane at 10 K, leading to a much larger permeation of the lighter species as compared to the heavier isotope. RPMD provides an efficient approach for studying the separation of He isotopes, taking into account quantum effects of light nuclei motions at low temperatures, which classical methods fail to capture.

Finite-size effects in the reconstruction of dynamic properties from ab initio path integral Monte Carlo simulations

We systematically investigate finite-size effects in the dynamic structure factor S(q,ω) of the uniform electron gas obtained via the analytic continuation of ab initio path integral Monte Carlo data for the imaginary-time density-density correlation function F(q,τ). Using the recent scheme by Dornheim et al. [Phys. Rev. Lett. 121, 255001 (2018)PRLTAO0031-900710.1103/PhysRevLett.121.255001], we find that the reconstructed spectra are not afflicted with any finite-size effects for as few as N=14 electrons both at warm dense matter (WDM) conditions and at the margins of the strongly correlated electron liquid regime. Our results further corroborate the high quality of our current description of the dynamic density response of correlated electrons, which is of high importance for many applications in WDM theory and beyond.