All Publications

In this work, we report potential energy surfaces (PESs) of the sodium dimer calculated by variational (VMC) and lattice-regularized diffusion Monte Carlo (LRDMC). The VMC calculation is accurate for determining the equilibrium distance and the qualitative shape of the experimental PES. Remarkably, after the application of the LRDMC projection to this single determinant ansatz, namely, the Jastrow Antisymmetrized Geminal Power (JAGP), chemical accuracy (∼1 kcal/mol) is reached in the binding energy, and the obtained equilibrium internuclear distance and harmonic vibrational frequency are in very good agreement with the experimental ones. This outcome is crucially dependent on the quality of the optimization used to determine the best possible trial function within the chosen ansatz. The strategy adopted in this work is to minimize the variational energy by initializing the trial function with the density functional theory (DFT) single determinant ansatz expanded exactly in the same atomic basis used for the corresponding VMC and LRDMC calculations. This atomic basis is reshaped ad-hoc for QMC calculations. Indeed, we multiply the standard Gaussian-type atomic orbitals by a one-body Jastrow factor, satisfying, in this way, the electron-ion cusp conditions. In order to achieve these important advantages, we have defined a very efficient DFT algorithm in the mentioned basis, by estimating the corresponding matrix elements on a mesh, and by using a much finer mesh grid in the vicinity of nuclei.

Within the ground-state auxiliary-field quantum Monte Carlo technique, we introduce discrete Hubbard-Stratonovich transformations (HSTs) that are also suitable for spatially inhomogeneous trial functions. The discrete auxiliary fields introduced here are coupled to local spin or charge operators fluctuating around their Hartree-Fock values. The formalism can be considered a generalization of the discrete HSTs by J. E. Hirsch [Phys. Rev. B 28, 4059 (1983)PRBMDO0163-182910.1103/PhysRevB.28.4059] or a compactification of the shifted-contour auxiliary-field Monte Carlo formalism by N. Rom et al. [Chem. Phys. Lett. 270, 382 (1997)CHPLBC0009-261410.1016/S0009-2614(97)00370-9]. An improvement of the acceptance ratio is found for a real auxiliary field, while an improvement of the average sign is found for a purely imaginary auxiliary field. Efficiencies of the different HSTs are tested in the single-band Hubbard model at and away from half filling by studying the staggered magnetization and energy expectation values, respectively.

An unbiased zeroerature auxiliary-field quantum Monte Carlo method is employed to analyze the nature of the semimetallic phase of the two-dimensional Hubbard model on the honeycomb lattice at half filling. It is shown that the quasiparticle weight Z of the massless Dirac fermions at the Fermi level, which characterizes the coherence of zero-energy single-particle excitations, can be evaluated in terms of the long-distance equal-time single-particle Green's function. If this quantity remains finite in the thermodynamic limit, the low-energy single-particle excitations of the correlated semimetallic phase are described by a Fermi-liquid-type single-particle Green's function. Based on the unprecedentedly large-scale numerical simulations on finite-size clusters containing more than 10 000 sites, we show that the quasiparticle weight remains finite in the semimetallic phase below a critical interaction strength. This is also supported by the long-distance algebraic behavior (∼r-2, where r is distance) of the equal-time single-particle Green's function that is expected for the Fermi liquid. Our result thus provides a numerical confirmation of Fermi-liquid theory in two-dimensional correlated metals.

We report a quantum Monte Carlo study, on a very simple but nevertheless very instructive model system of four hydrogen atoms, recently proposed in Gasperich et al. [J. Chem. Phys. 147, 074106 (2017)]. We find that the Jastrow correlated Antisymmetrized Geminal Power (JAGP) is able to recover most of the correlation energy even when the geometry is symmetric and the hydrogens lie on the edges of a perfect square. Under such conditions, the diradical character of the molecule ground state prevents a single determinant Ansatz to achieve an acceptable accuracy, whereas the JAGP performs very well for all geometries. Remarkably, this is obtained with a similar computational effort. Moreover, we find that the Jastrow factor is fundamental in promoting the correct resonances among several configurations in the JAGP, which cannot show up in the pure Antisymmetrized Geminal Power (AGP). We also show the extremely fast convergence of this approach in the extension of the basis set. Remarkably, only the simultaneous optimization of the Jastrow and the AGP part of our variational Ansatz is able to recover an almost perfect nodal surface, yielding therefore state of the art energies, almost converged in the complete basis set limit, when the so called diffusion Monte Carlo is applied.

By employing unbiased numerical methods, we show that pulse irradiation can induce unconventional superconductivity even in the Mott insulator of the Hubbard model. The superconductivity found here in the photoexcited state is due to the η-pairing mechanism, characterized by staggered pair-density-wave oscillations in the off-diagonal long-range correlation, and is absent in the ground-state phase diagram; i.e., it is induced neither by a change of the effective interaction of the Hubbard model nor by simple photocarrier doping. Because of the selection rule, we show that the nonlinear optical response is essential to increase the number of η pairs and thus enhance the superconducting correlation in the photoexcited state. Our finding demonstrates that nonequilibrium many-body dynamics is an alternative pathway to access a new exotic quantum state that is absent in the ground-state phase diagram, and also provides an alternative mechanism for enhancing superconductivity.

We show that, by an appropriate choice of auxiliary fields and exact integration of the phonon degrees of freedom, it is possible to define a "sign-free" path integral for the so-called Hubbard-Holstein model at half filling. We use a statistical method, based on an accelerated and efficient Langevin dynamics, for evaluating all relevant correlation functions of the model. Preliminary calculations at U/t=4 and U/t=1, for ω0/t=1, indicate a region around U≃g2ω0 without either antiferromagnetic or charge-density-wave orders that is much wider compared to previous approximate calculations. The elimination of the sign problem in a model without explicit particle-hole symmetry may open different perspectives for strongly correlated models, even away from the purely attractive or particle-hole symmetric cases.

Angle-resolved photoemission spectroscopy allows one to visualize in momentum space the probability weight maps of electrons subtracted from molecules deposited on a substrate. The interpretation of these maps usually relies on the plane wave approximation through the Fourier transform of single particle orbitals obtained from density functional theory. Here we propose a first-principle many-body approach based on quantum Monte Carlo (QMC) to directly calculate the quasi-particle wave functions (also known as Dyson orbitals) of molecules in momentum space. The comparison between these correlated QMC images and their single particle counterpart highlights features that arise from many-body effects. We test the QMC approach on the linear C2H2, CO2, and N2 molecules, for which only small amplitude remodulations are visible. Then, we consider the case of the pentacene molecule, focusing on the relationship between the momentum space features and the real space quasi-particle orbital. Eventually, we verify the correlation effects present in the metal CuCl42- planar complex.

By using variational Monte Carlo and auxiliary-field quantum Monte Carlo methods, we perform an accurate finite-size scaling of the s-wave superconducting order parameter and the pairing correlations for the negative-U Hubbard model at zero temperature in the square lattice. We show that the twist-averaged boundary conditions (TABCs) are extremely important to control finite-size effects and to achieve smooth and accurate extrapolations to the thermodynamic limit. We also show that TABCs are much more efficient in the grand-canonical ensemble rather than in the standard canonical ensemble with fixed number of electrons. The superconducting order parameter as a function of the doping is presented for several values of |U|/t and is found to be significantly smaller than the mean-field BCS estimate already for moderate couplings. This reduction is understood by a variational ansatz able to describe the low-energy behavior of the superconducting phase by means of a suitably chosen Jastrow factor including long-range density-density correlations.

The phase diagram of isotropically expanded graphene cannot be correctly predicted by ignoring either electron correlations, or mobile carbons, or the effect of applied stress, as was done so far. We calculate the ground state enthalpy (not just energy) of strained graphene by an accurate off-lattice quantum Monte Carlo correlated ansatz of great variational flexibility. Following undistorted semimetallic graphene at low strain, multideterminant Heitler-London correlations stabilize between ≃8.5% and ≃15% strain an insulating Kekulé-like dimerized (DIM) state. Closer to a crystallized resonating-valence bond than to a Peierls state, the DIM state prevails over the competing antiferromagnetic insulating state favored by density-functional calculations which we conduct in parallel. The DIM stressed graphene insulator, whose gap is predicted to grow in excess of 1 eV before failure near 15% strain, is topological in nature, implying under certain conditions 1D metallic interface states lying in the bulk energy gap.

We investigate a semimetal-superconductor phase transition of two-dimensional Dirac electrons at zero temperature by large-scale and essentially unbiased quantum Monte Carlo simulations for the half-filled attractive Hubbard model on the triangular lattice, in the presence of alternating magnetic π flux, that is introduced to construct two Dirac points in the one-particle bands at the Fermi level. This phase transition is expected to describe quantum criticality of the chiral XY class in the framework of the Gross-Neveu model, where, in the ordered phase, the U(1) symmetry is spontaneously broken and a mass gap opens in the excitation spectrum. We compute the order parameter of the s-wave superconductivity and estimate the quasiparticle weight from the long-distance behavior of the single-particle Green's function. These calculations allow us to obtain the critical exponents of this transition in a reliable and accurate way. Our estimate for the critical exponents is in good agreement with those obtained for a transition to a Kekulé valence bond solid, where an emergent U(1) symmetry is proposed [Z.-X. Li, Nat. Commun. 8, 314 (2017)2041-172310.1038/s41467-017-00167-6].

The dynamical spin structure factor is computed within a variational framework to study the one-dimensional J1-J2 Heisenberg model. Starting from Gutzwiller-projected fermionic wave functions, the low-energy spectrum is constructed from two-spinon excitations. The direct comparison with Lanczos calculations on small clusters demonstrates the excellent description of both gapless and gapped (dimerized) phases, including incommensurate structures for J2/J1>0.5. Calculations on large clusters show how the intensity evolves when increasing the frustrating ratio and give an unprecedented accurate characterization of the dynamical properties of (nonintegrable) frustrated spin models.

Understanding planetary interiors is directly linked to our ability of simulating exotic quantum mechanical systems such as hydrogen (H) and hydrogen-helium (H-He) mixtures at high pressures and temperatures. Equation of state (EOS) tables based on density functional theory are commonly used by planetary scientists, although this method allows only for a qualitative description of the phase diagram. Here we report quantum Monte Carlo (QMC) molecular dynamics simulations of pure H and H-He mixture. We calculate the first QMC EOS at 6000 K for a H-He mixture of a protosolar composition, and show the crucial influence of He on the H metallization pressure. Our results can be used to calibrate other EOS calculations and are very timely given the accurate determination of Jupiter's gravitational field from the NASA Juno mission and the effort to determine its structure.

Fifty years ago Walter Kohn speculated that a zero-gap semiconductor might be unstable against the spontaneous generation of excitons-electron-hole pairs bound together by Coulomb attraction. The reconstructed ground state would then open a gap breaking the symmetry of the underlying lattice, a genuine consequence of electronic correlations. Here we show that this excitonic insulator is realized in zero-gap carbon nanotubes by performing first-principles calculations through many-body perturbation theory as well as quantum Monte Carlo. The excitonic order modulates the charge between the two carbon sublattices opening an experimentally observable gap, which scales as the inverse of the tube radius and weakly depends on the axial magnetic field. Our findings call into question the Luttinger liquid paradigm for nanotubes and provide tests to experimentally discriminate between excitonic and Mott insulators.

Over the past several decades, computational approaches to studying strongly-interacting systems have become increasingly varied and sophisticated. This book provides a comprehensive introduction to state-of-the-art quantum Monte Carlo techniques relevant for applications in correlated systems. Providing a clear overview of variational wave functions, and featuring a detailed presentation of stochastic samplings including Markov chains and Langevin dynamics, which are developed into a discussion of Monte Carlo methods. The variational technique is described, from foundations to a detailed description of its algorithms. Further topics discussed include optimisation techniques, real-time dynamics and projection methods, including Green's function, reptation and auxiliary-field Monte Carlo, from basic definitions to advanced algorithms for efficient codes, and the book concludes with recent developments on the continuum space. Quantum Monte Carlo Approaches for Correlated Systems provides an extensive reference for students and researchers working in condensed matter theory or those interested in advanced numerical methods for electronic simulation.

By using variational wave functions and quantum Monte Carlo techniques, we investigate the interplay between electron-electron and electron-phonon interactions in the two-dimensional Hubbard-Holstein model. Here, the ground-state phase diagram is triggered by several energy scales, i.e., the electron hopping t, the on-site electron-electron interaction U, the phonon energy ω0, and the electron-phonon coupling g. At half filling, the ground state is an antiferromagnetic insulator for U2g2/ω0, while it is a charge-density-wave (or bipolaronic) insulator for U2g2/ω0. In addition to these phases, we find a superconducting phase that intrudes between them. For ω0/t=1, superconductivity emerges when both U/t and 2g2/tω0 are small; then, by increasing the value of the phonon energy ω0, it extends along the transition line between antiferromagnetic and charge-density-wave insulators. Away from half filling, phase separation occurs when doping the charge-density-wave insulator, while a uniform (superconducting) ground state is found when doping the superconducting phase. In the analysis of finite-size effects, it is extremely important to average over twisted boundary conditions, especially in the weak-coupling limit and in the doped case.

We introduce a novel approach for a fully quantum description of coupled electron-ion systems from first principles. It combines the variational quantum Monte Carlo solution of the electronic part with the path integral formalism for the quantum nuclear dynamics. On the one hand, the path integral molecular dynamics includes nuclear quantum effects by adding a set of fictitious classical particles (beads) aimed at reproducing nuclear quantum fluctuations via a harmonic kinetic term. On the other hand, variational quantum Monte Carlo can provide Born-Oppenheimer potential energy surfaces with a precision comparable to the most-Advanced post-Hartree-Fock approaches, and with a favorable scaling with the system size. In order to cope with the intrinsic noise due to the stochastic nature of quantum Monte Carlo methods, we generalize the path integral molecular dynamics using a Langevin thermostat correlated according to the covariance matrix of quantum Monte Carlo nuclear forces. The variational parameters of the quantum Monte Carlo wave function are evolved during the nuclear dynamics, such that the Born-Oppenheimer potential energy surface is unbiased. Statistical errors on the wave function parameters are reduced by resorting to bead grouping average, which we show to be accurate and well-controlled. Our general algorithm relies on a Trotter breakup between the dynamics driven by ionic forces and the one set by the harmonic interbead couplings. The latter is exactly integrated, even in the presence of the Langevin thermostat, thanks to the mapping onto an Ornstein-Uhlenbeck process. This framework turns out to be also very efficient in the case of noiseless (deterministic) ionic forces. The new implementation is validated on the Zundel ion (H5O2+) by direct comparison with standard path integral Langevin dynamics calculations made with a coupled cluster potential energy surface. Nuclear quantum effects are confirmed to be dominant over thermal effects well beyond room temperature, giving the excess proton an increased mobility by quantum tunneling.

We propose an ab initio molecular dynamics method, capable of dramatically reducing the autocorrelation time required for the simulation of classical and quantum particles at finite temperatures. The method is based on an efficient implementation of a first order Langevin dynamics modified by means of a suitable, position dependent acceleration matrix S. Here, we apply this technique to both Lennard-Jones models, to demonstrate the accuracy and speeding-up of the sampling, and within a quantum Monte Carlo based wave function approach, for determining the phase diagram of high-pressure hydrogen with simulations much longer than the autocorrelation time. With the proposed method, we are able to equilibrate in a few hundred steps even close to the liquid-liquid phase transition (LLT). Within our approach, we find that the LLT transition is consistent with recent density functionals predicting a much larger transition pressure when the long range dispersive forces are taken into account.

We present numerical results for the equation of state of an infinite chain of hydrogen atoms. A variety of modern many-body methods are employed, with exhaustive cross-checks and validation. Approaches for reaching the continuous space limit and the thermodynamic limit are investigated, proposed, and tested. The detailed comparisons provide a benchmark for assessing the current state of the art in many-body computation, and for the development of new methods. The ground-state energy per atom in the linear chain is accurately determined versus bond length, with a confidence bound given on all uncertainties.

We present a systematic investigation of the special twist method introduced by Rajagopal et al. [Phys. Rev. B 51, 10591 (1995)PRBMDO0163-182910.1103/PhysRevB.51.10591] for reducing finite-size effects in correlated calculations of periodic extended systems with Coulomb interactions and Fermi statistics. We propose a procedure for finding special twist values which, at variance with previous applications of this method, reproduce the energy of the mean-field infinite-size limit solution within an adjustable (arbitrarily small) numerical error. This choice of the special twist is shown to be the most accurate single-twist solution for curing one-body finite-size effects in correlated calculations. For these reasons we dubbed our procedure "exact special twist" (EST). EST only needs a fully converged independent-particles or mean-field calculation within the primitive cell and a simple fit to find the special twist along a specific direction in the Brillouin zone. We first assess the performances of EST in a simple correlated model such as the three-dimensional electron gas. Afterwards, we test its efficiency within ab initio quantum Monte Carlo simulations of metallic elements of increasing complexity. We show that EST displays an overall good performance in reducing finite-size errors comparable to the widely used twist average technique but at a much lower computational cost since it involves the evaluation of just one wave function. We also demonstrate that the EST method shows similar performances in the calculation of correlation functions, such as the ionic forces for structural relaxation and the pair radial distribution function in liquid hydrogen. Our conclusions point to the usefulness of EST for correlated supercell calculations; our method will be particularly relevant when the physical problem under consideration requires large periodic cells.

Clay minerals are ubiquitous in nature, and the manner in which they interact with their surroundings has important industrial and environmental implications. Consequently, a molecular-level understanding of the adsorption of molecules on clay surfaces is crucial. In this regard computer simulations play an important role, yet the accuracy of widely used empirical force fields (FF) and density functional theory (DFT) exchange-correlation functionals is often unclear in adsorption systems dominated by weak interactions. Herein we present results from quantum Monte Carlo (QMC) for water and methanol adsorption on the prototypical clay kaolinite. To the best of our knowledge, this is the first time QMC has been used to investigate adsorption at a complex, natural surface such as a clay. As well as being valuable in their own right, the QMC benchmarks obtained provide reference data against which the performance of cheaper DFT methods can be tested. Indeed using various DFT exchange-correlation functionals yields a very broad range of adsorption energies, and it is unclear a priori which evaluation is better. QMC reveals that in the systems considered here it is essential to account for van der Waals (vdW) dispersion forces since this alters both the absolute and relative adsorption energies of water and methanol. We show, via FF simulations, that incorrect relative energies can lead to significant changes in the interfacial densities of water and methanol solutions at the kaolinite interface. Despite the clear improvements offered by the vdW-corrected and the vdW-inclusive functionals, absolute adsorption energies are often overestimated, suggesting that the treatment of vdW forces in DFT is not yet a solved problem.