One of the major obstacles in modern quantum chemistry is the need to employ exceedingly large one-electron basis sets in electron correlation methods to obtain reliably accurate results. The family of explicitly correlated wave function methods circumvents this issue by introducing terms that explicitly depend on the interelectronic distance, thereby enforcing the correct electron–electron cusp behavior, albeit at the cost of additional complicated many-electron integrals when evaluated without further approximation. The main focus of the present work is the efficient, low-scaling evaluation of explicitly correlated approaches, examined here through the explicitly correlated correction to second-order Møller–Plesset perturbation theory (MP2-F12) as a representative example. While these methods are highly effective in reducing the basis set incompleteness error, the cost associated with the corresponding correction often far exceeds that of the underlying correlation method. To overcome these limitations, several techniques and strategies are introduced to significantly reduce computational demand and resolve the main bottlenecks that arise. A particularly straightforward approach is provided by the scaled MP2-F12 ansatz, which introduces a single empirically determined factor that scales the direct-type geminal–geminal terms to closely reproduce the MP2-F12 correction, thereby omitting the most expensive exchange-type contributions. Furthermore, the computation of Fock matrix elements spanning multiple orbital spaces is drastically accelerated by extending techniques originally developed for Hartree–Fock and density functional theory. The direct contributions are evaluated using a modified resolution-of-the-identity Coulomb (RI-J) approach employing the J-engine ansatz for integral construction, while the exchange contributions are treated via three-dimensional numerical quadrature in the form of the seminumerical sn-LinK method. These approaches reduce the computational cost by more than three orders of magnitude compared to the approximation-free evaluation while maintaining high accuracy, with additional speedups achieved by exploiting the massive parallelism of graphics processing units. The central focus of this work is the efficient decomposition of complicated exchange-type multi-electron integrals, where numerical quadrature plays a key role in improving performance. The complementary auxiliary basis set resolution of the identity (CABS-RI) reduces the most expensive exchange-type integral to a three-electron integral, which is then efficiently decomposed via numerical quadrature. In this context, a highly efficient batch-wise, distance-dependent integral screening exploiting the short-range behavior of F12 operators is introduced, significantly improving performance and enabling near-linear-scaling evaluations while retaining high accuracy. These strategies are extended to all exchange-type contributions in F12 theory, leading to an alternative formulation of integrals better suited for efficient factorization. A novel combination of numerical quadrature with density fitting is presented for the efficient computation of products of four-center two-electron integrals, alongside algorithms covering all exchange-type contributions, including advantageous batching strategies. In general, the formal computational cost with system size M is reduced from O(M⁵) to O(M⁴), with further improvements appearing feasible through the use of localized molecular orbitals combined with advanced integral screening techniques.
Computational acceleration of orbital-invariant local correlation methods on graphics processing units (GPUs) has remained largely unexplored due to substantial algorithmic complexities. The runtime efficiency of GPU-implemented local correlation theories can be significantly constrained by the degree of parallelism in the orbital localization procedure, the iterative solution of the local wave function, and the adaptation of CUDA kernels to inherently local or sparse operations. Using the second-order Møller-Plesset perturbation (MP2) theory, we present a multi-GPU implementation for large-scale third-order many-body expansion orbital-specific virtual MP2 (MBE(3)-OSV-MP2) energy calculations. Accordingly, our algorithms and implementation address the GPU parallelization to maximize peak utilization and the parallelism of local MP2 computation in several aspects, including Jacobi-Pipek-Mezey localization, randomized OSV (rOSV) generation, direct MP2 integral regeneration, and CUDA kernel adaptation for local operations. The GPU-based MBE(3)-OSV-MP2 energy computation achieves an effective empirical N1.9 scaling up to (Gly)40/def2-TZVP and 84% parallel efficiency without localization up to 24 GPUs distributed on multiple nodes. The present implementation delivers 40-fold wall time speedup of the canonical RI-MP2 and 10-fold speedup of the central processing units (CPU)-based MBE(3)-OSV-MP2 for (H2O)128/cc-pVDZ and (H2O)190/cc-pVDZ, respectively. A large-scale computation of a 784-atom insulin peptide yields the full MBE(3)-OSV-MP2 energies in 24 min with cc-pVDZ (7571 basis functions) and in 6.4 h with cc-pVTZ (17,448 basis functions) on 8 NVIDIA A800 GPUs. Our work opens new possibilities for fast GPU-based local correlation calculations on real-world macromolecules.
Post-Hartree-Fock methods are wavefunction-based quantum chemistry approaches that explicitly account for the electron correlation beyond the mean-field approximation, allowing accurate determination of electron orbitals and energies in a given system. However, their computational cost scales as the fourth power or higher with respect to the number of atoms, severely limiting their practical applications. This paper proposes an efficient GPU implementation of the Møller-Plesset method, a post-Hartree-Fock approach based on perturbation theory. We focus on the computation of second- and thirdorder perturbation energies (MP2 and MP3) for closed-shell systems, which incorporate electron correlation on top of the Hartree-Fock energy. To accelerate the extensive contraction operations arising in the Møller-Plesset method, our implementation employs a formulation that maximizes the use of highly optimized matrix multiplication routines. Furthermore, we effectively utilize computational capabilities of GPU architecture through a hierarchical reduction using warp-level primitives and shared memory, as well as contiguous coalesced access to matrix elements. Performance evaluation on an NVIDIA A100 GPU shows up to 255.5 and 71.1 times speedups for the MP2 correlation energy computation, compared with the widely used quantum chemistry software package PySCF and its GPU plugin GPU4PySCF, respectively. In the MP3 calculations, our DGEMM-based implementation achieves a maximum speedup of 96.2 times over PySCF.
The computational acceleration of orbital-invariant local correlation methods on graphics processing units (GPUs) has remained largely unexplored due to substantial algorithmic complexities. The runtime efficiency of GPU-implemented local correlation theories can be significantly constrained by the parallelizable degree of the orbital localization procedure, the iterative solution of the local wave function, and the adaptation of CUDA kernels to inherently local or sparse operations. Using the second-order Møller-Plesset perturbation (MP2) theory, we present a multi-GPU implementation for large-scale third-order many-body expansion orbital-specific virtual MP2 (MBE(3)-OSV-MP2) energy calculations. Accordingly, our algorithms and implementation address the GPU parallelization ability for peak utilization and parallelism of local MP2 computation in several aspects, including Jacobi-Pipek-Mezey localization, randomized OSV generation, direct MP2 integral regeneration, as well as CUDA kernel adaptation to local operations. The GPU-based MBE(3)-OSV-MP2 energy computation achieves O(N1.9) scaling and 84\% parallel efficiency up to 24 GPUs distributed on multiple nodes. The present implementation delivers 40-fold wall-time speedup of the canonical RI-MP2 and 10-fold speedup of the CPU-based MBE(3)-OSV-MP2 for (H2O)128/cc-pVDZ and (H2O)190/cc-pVDZ, respectively. A large scale computation of 784-atom insulin peptide yields the full MBE(3)-OSV-MP2 energies in 24 minutes with cc-pVDZ (7571 basis functions) and 6.4 hours with cc-pVTZ (17448 basis functions) on 8 NVIDIA A800 GPUs. Our work opens up new possibilities for performing fast GPU-based local correlation calculations on real-life macromolecules.
A parallel algorithm is described for the coupled-cluster singles and doubles method augmented with a perturbative correction for triple excitations [CCSD(T)] using the resolution-of-the-identity (RI) approximation for two-electron repulsion integrals (ERIs). The algorithm bypasses the storage of four-center ERIs by adopting an integral-direct strategy. The CCSD amplitude equations are given in a compact quasi-linear form by factorizing them in terms of amplitude-dressed three-center intermediates. A hybrid MPI/OpenMP parallelization scheme is employed, which uses the OpenMP-based shared memory model for intranode parallelization and the MPI-based distributed memory model for internode parallelization. Parallel efficiency has been optimized for all terms in the CCSD amplitude equations. Two different algorithms have been implemented for the rate-limiting terms in the CCSD amplitude equations that entail and -scaling computational costs, where N O and N V denote the number of correlated occupied and virtual orbitals, respectively. One of the algorithms assembles the four-center ERIs requiring N V 4 and N O 2 N V 2 -scaling memory costs in a distributed manner on a number of MPI ranks, while the other algorithm completely bypasses the assembling of quartic memory-scaling ERIs and thus largely reduces the memory demand. It is demonstrated that the former memory-expensive algorithm is faster on a few hundred cores, while the latter memory-economic algorithm shows a better strong scaling in the limit of a few thousand cores. The program is shown to exhibit a near-linear scaling, in particular for the compute-intensive triples correction step, on up to 8000 cores. The performance of the program is demonstrated via calculations involving molecules with 24–51 atoms and up to 1624 atomic basis functions. As the first application, the complete basis set (CBS) limit for the interaction energy of the π-stacked uracil dimer from the S66 data set has been investigated. This work reports the first calculation of the interaction energy at the CCSD(T)/aug-cc-pVQZ level without local orbital approximation. The CBS limit for the CCSD correlation contribution to the interaction energy was found to be −8.01 kcal/mol, which agrees very well with the value −7.99 kcal/mol reported by Schmitz, Hättig, and Tew [ Phys. Chem. Chem. Phys. 2014, 16, 22167−22178]. The CBS limit for the total interaction energy was estimated to be −9.64 kcal/mol.