We fully develop the Monte Carlo many-body Green’s function (MC-GF) method with the following enhancements: (1) The truncation order of the perturbation expansion of the Dyson self-energy is raised from the second order (MC-GF2) to the third order (MC-GF3) with the aid of a computerized procedure to enumerate and transform all 84 third-order Goldstone diagrams into Monte Carlo integrable expressions and then into central processing unit (CPU)/graphical processing unit (GPU)-parallel computer codes. (2) An efficient algorithm is proposed that computes all off-diagonal and diagonal elements of the MC-GF2 and MC-GF3 self-energy matrices by common subexpression elimination. (3) The frequency-independent approximation is lifted by introducing a method that computes frequency derivatives of the MC-GF2 and MC-GF3 self-energies up to any arbitrarily high order at nearly no additional computational cost. (4) The imaginary-time integration in the Laplace-transformed expressions of the self-energy is carried out stochastically (instead of using a quadrature in the previous implementations), resulting in a 50- to 200-fold speedup. (5) The efficiency of the redundant-walker convergence acceleration scheme is analyzed numerically, and the guidelines are established to select an optimal number of walkers for maximal efficiency. When such an optimal number is used, the cost per sample is constant of molecular size on either many CPUs or many GPUs. (6) The computational cost to obtain a binding energy within a given statistical uncertainty is observed to increase as (tentatively) O( n 4 ) and O( n 5 ) of molecular size ( n is the number of electrons) with and without the diagonal approximation, respectively, as opposed to O( n 5 ) and O( n 6 ) of the corresponding deterministic algorithms. With this method applied to the electron binding energies of C 60, we show that the third-order corrections to the self-energies are much greater in electron binding energies than in ground-state energies. They display a sign of oscillatory convergence toward experimental results, not necessarily improving the agreement with increasing perturbation order, justifying MC-GF3 and motivating even higher-order methods.
We fully develop the Monte Carlo many-body Green’s function (MC-GF) method with the following enhancements: (1) The truncation order of the perturbation expansion of the Dyson self-energy is raised from the second order (MC-GF2) to the third order (MC-GF3) with the aid of a computerized procedure to enumerate and transform all 84 third-order Goldstone diagrams into Monte Carlo integrable expressions and then into central processing unit (CPU)/graphical processing unit (GPU)-parallel computer codes. (2) An efficient algorithm is proposed that computes all off-diagonal and diagonal elements of the MC-GF2 and MC-GF3 self-energy matrices by common subexpression elimination. (3) The frequency-independent approximation is lifted by introducing a method that computes frequency derivatives of the MC-GF2 and MC-GF3 self-energies up to any arbitrarily high order at nearly no additional computational cost. (4) The imaginary-time integration in the Laplace-transformed expressions of the self-energy is carried out stochastically (instead of using a quadrature in the previous implementations), resulting in a 50- to 200-fold speedup. (5) The efficiency of the redundant-walker convergence acceleration scheme is analyzed numerically, and the guidelines are established to select an optimal number of walkers for maximal efficiency. When such an optimal number is used, the cost per sample is constant of molecular size on either many CPUs or many GPUs. (6) The computational cost to obtain a binding energy within a given statistical uncertainty is observed to increase as (tentatively) O( n 4 ) and O( n 5 ) of molecular size ( n is the number of electrons) with and without the diagonal approximation, respectively, as opposed to O( n 5 ) and O( n 6 ) of the corresponding deterministic algorithms. With this method applied to the electron binding energies of C 60, we show that the third-order corrections to the self-energies are much greater in electron binding energies than in ground-state energies. They display a sign of oscillatory convergence toward experimental results, not necessarily improving the agreement with increasing perturbation order, justifying MC-GF3 and motivating even higher-order methods.
A new, alternative set of interpretation rules of Feynman-Goldstone diagrams for many-body perturbation theory is proposed, which translates diagrams into algebraic expressions suitable for direct Monte Carlo integrations. A vertex of a diagram is associated with a Coulomb interaction (rather than a two-electron integral) and an edge with the trace of a Green's function in real space and imaginary time. With these, 12 diagrams of third-order many-body perturbation (MP3) theory are converted into 20-dimensional integrals, which are then evaluated by a Monte Carlo method. It uses redundant walkers for convergence acceleration and a weight function for importance sampling in conjunction with the Metropolis algorithm. The resulting Monte Carlo MP3 method has low-rank polynomial size dependence of the operation cost, a negligible memory cost, and a naturally parallel computational kernel, while reproducing the correct correlation energies of small molecules within a few mEh after 10(6) Monte Carlo steps.
In the Monte Carlo second-order many-body perturbation (MC-MP2) method, the long sum-of-product matrix expression of the MP2 energy, whose literal evaluation may be poorly scalable, is recast into a single high-dimensional integral of functions of electron pair coordinates, which is evaluated by the scalable method of Monte Carlo integration. The sampling efficiency is further accelerated by the redundant-walker algorithm, which allows a maximal reuse of electron pairs. Here, a multitude of graphical processing units (GPUs) offers a uniquely ideal platform to expose multilevel parallelism: fine-grain data-parallelism for the redundant-walker algorithm in which millions of threads compute and share orbital amplitudes on each GPU; coarse-grain instruction-parallelism for near-independent Monte Carlo integrations on many GPUs with few and infrequent interprocessor communications. While the efficiency boost by the redundant-walker algorithm on central processing units (CPUs) grows linearly with the number of electron pairs and tends to saturate when the latter exceeds the number of orbitals, on a GPU it grows quadratically before it increases linearly and then eventually saturates at a much larger number of pairs. This is because the orbital constructions are nearly perfectly parallelized on a GPU and thus completed in a near-constant time regardless of the number of pairs. In consequence, an MC-MP2/cc-pVDZ calculation of a benzene dimer is 2700 times faster on 256 GPUs (using 2048 electron pairs) than on two CPUs, each with 8 cores (which can use only up to 256 pairs effectively). We also numerically determine that the cost to achieve a given relative statistical uncertainty in an MC-MP2 energy increases as O ( n 3 ) or better with system size n, which may be compared with the O ( n 5 ) scaling of the conventional implementation of deterministic MP2. We thus establish the scalability of MC-MP2 with both system and computer sizes.
A systematically converging series of ab initio, postdensity functional, size-consistent, electron-correlated approximations is desired for predictive computing of electronic band structures of insulating, semiconducting, and metallic solids. A series that meets all of these desiderata (except the applicability to metals) is ab initio many-body Green’s function theory based on Gaussian-type orbital (GTO) basis sets. Here, its leading-order approximation, the second-order Green’s function (GF2) method in the diagonal and frequency-independent approximations with the aug-cc-pVDZ basis set, is applied to the fundamental band gaps of three semiconductors (diamond, silicon, and silicon carbide in the zincblende structure) using cluster models. Corrections are made to the basis set incompleteness errors by the explicit correlation (F12) ansatz (GF2-F12) for the valence band edges. The crystals are modeled as surface-passivated clusters of increasing sizes, whose wave functions are expanded by up to 2709 GTO basis functions. Immense computational costs of these calculations are overcome by the highly scalable stochastic algorithm of the Monte Carlo GF2-F12 method, whose operation cost per state increases only as a cubic power of system size, which has a tiny memory footprint and easily achieves near-perfect parallel efficiency on thousands of CPUs or on hundreds of GPUs. The correlated, F12-corrected highest-occupied and lowest-unoccupied molecular orbital energy (HOMO–LUMO) gap is 5.78 ± 0.07 eV for C 87 H 76 as compared with the experimental value of the fundamental (indirect) band gap of bulk diamond at 5.48 eV. The correlated, F12-corrected HOMO–LUMO gaps for Si 75 H 76 and Si 32 C 43 H 76 are 2.56 ± 0.15 and 3.50 ± 0.12 eV, respectively, which are expected to decrease further with increasing cluster sizes. The experimental fundamental (indirect) band gaps of bulk silicon and silicon carbide are 1.17 and 2.42 eV, respectively.
A scalable, stochastic algorithm evaluating the fourth-order many-body perturbation (MP4) correction to energy is proposed. Three hundred Goldstone diagrams representing the MP4 correction are computer generated and then converted into algebraic formulas expressed in terms of Green’s functions in real space and imaginary time. They are evaluated by the direct (i.e., non-Markov, non-Metropolis) Monte Carlo (MC) integration accelerated by the redundant-walker and control-variate algorithms. The resulting MC-MP4 method is efficiently parallelized and is shown to display O(n5.3) size-dependence of cost, which is nearly two ranks lower than the O(n7) dependence of the deterministic MP4 algorithm. It evaluates the MP4/aug-cc-pVDZ energy for benzene, naphthalene, phenanthrene, and corannulene with the statistical uncertainty of 10 mEh (1.1% of the total basis-set correlation energy), 38 mEh (2.6%), 110 mEh (5.5%), and 280 mEh (9.0%), respectively, after about 109 MC steps.
We report an implementation of self-consistent Green's function many-body theory within a second-order approximation (GF2) for application with molecular systems. This is done by iterative solution of the Dyson equation expressed in matrix form in an atomic orbital basis, where the Green's function and self-energy are built on the imaginary frequency and imaginary time domain, respectively, and fast Fourier transform is used to efficiently transform these quantities as needed. We apply this method to several archetypical examples of strong correlation, such as a H32 finite lattice that displays a highly multireference electronic ground state even at equilibrium lattice spacing. In all cases, GF2 gives a physically meaningful description of the metal to insulator transition in these systems, without resorting to spin-symmetry breaking. Our results show that self-consistent Green's function many-body theory offers a viable route to describing strong correlations while remaining within a computationally tractable single-particle formalism.
We present an implementation of a fully self-consistent finite temperature second order Green's function perturbation theory (GF2) within the diagrammatic Monte Carlo framework. In contrast to the previous implementations of stochastic GF2 ({\it J. Chem. Phys.},{\bf 151}, 044144 (2019)), the current self-consistent stochastic GF2 does not introduce a systematic bias of the resulting electronic energies. Instead, the introduced implementation accounts for the stochastic errors appearing during the solution of the Dyson equation. We present an extensive discussion of the error handling necessary in a self-consistent procedure resulting in dressed Green's function lines. We test our method on a series of simple molecular examples.