diff --git a/docs/methods/post_hartree_fock/mp2.md b/docs/methods/post_hartree_fock/mp2.md index 03496621fd..c08d9bc812 100644 --- a/docs/methods/post_hartree_fock/mp2.md +++ b/docs/methods/post_hartree_fock/mp2.md @@ -6,13 +6,20 @@ $$ E^{(2)} = - \sum_{ijab} \frac{(ia|jb)[2(ia|jb)-(ib|ja)]}{\epsilon_{a}+\epsilon_{b}-\epsilon_{i}-\epsilon_{j}} $$ -with the GPW method as described in [](#DelBen2012) and [](#DelBen2013). CP2K also implements double-hybrid functionals. +with the GPW method as described in [](#DelBen2012) and [](#DelBen2013). CP2K also implements +double-hybrid functionals. -CP2K provides three different implementations: A canonical implementation, a GPW-based implementation and a RI-based implementation. Different functionals require rescaling of different contributions (singulet, triplet). This is achieved using the keywords [SCALE_S](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_S) and [SCALE_T](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_T) for all implementations. +CP2K provides three different implementations: A canonical implementation, a GPW-based +implementation and a RI-based implementation. Different functionals require rescaling of different +contributions (singulet, triplet). This is achieved using the keywords +[SCALE_S](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_S) and +[SCALE_T](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_T) for all implementations. ## Canonical Implementation -The canonical implementation is available for GPW and GAPW reference calculations and thus allows for core-corrections. It is the most costly implementation. Gradients (Nuclear forces, stress tensors) are not analytically available. A simple WF_CORRELATION section looks like +The canonical implementation is available for GPW and GAPW reference calculations and thus allows +for core-corrections. It is the most costly implementation. Gradients (Nuclear forces, stress +tensors) are not analytically available. A simple WF_CORRELATION section looks like ```none &WF_CORRELATION @@ -30,7 +37,9 @@ The canonical implementation is available for GPW and GAPW reference calculation ## GPW-based Implementation -The MP2-based implementation is available only for GPW calculations.[](#DelBen2012) It is the cheaper than the canonical implementation. Gradients (Nuclear forces, stress tensors) are not analytically available. A simple WF_CORRELATION section looks like +The MP2-based implementation is available only for GPW calculations.[](#DelBen2012) It is the +cheaper than the canonical implementation. Gradients (Nuclear forces, stress tensors) are not +analytically available. A simple WF_CORRELATION section looks like ```none &WF_CORRELATION @@ -53,7 +62,10 @@ The MP2-based implementation is available only for GPW calculations.[](#DelBen20 ## RI-based Implementation -The RI-based implementation is the most affordable MP2-implementation.[](#DelBen2013) Analytical gradients are available for GPW-based integrals (not for GAPW reference calculations).[](#DelBen2015b) MME-integration does not support stress tensor calculations, Obara-Saika integration does not implement any gradients. A simple WF_CORRELATION section looks like +The RI-based implementation is the most affordable MP2-implementation.[](#DelBen2013) Analytical +gradients are available for GPW-based integrals (not for GAPW reference +calculations).[](#DelBen2015b) MME-integration does not support stress tensor calculations, +Obara-Saika integration does not implement any gradients. A simple WF_CORRELATION section looks like ```none &WF_CORRELATION @@ -71,13 +83,167 @@ The RI-based implementation is the most affordable MP2-implementation.[](#DelBen &RI_MP2 # Larger block sizes require more memory but reduce communication, -1 (default) let CP2K choose it BLOCK_SIZE 2 - # This keyword determines how many copies of the rank-three tensor are kept in the memory of all CPUs. Large number of groups increasethe memory demands but reduce communication (default: -1, automatic determination) + # This keyword determines how many copies of the rank-three tensor are kept in the memory of all CPUs. Large number of groups increase the memory demands but reduce communication (default: -1, automatic determination) NUMBER_INTEGRATION_GROUPS 2 &END &END ``` +## RI-MP2 Gradient Calculations + +Analytical gradients (nuclear gradients, stress tensors) are available for RI-MP2 and related +double-hybrid functionals.[](#DelBen2015b),[](#Stein2022) HF calculations may also be accelerated by +the ADMM approximation which is recommended if very diffuse basis functions are included. The +computational costs of the MP2 part is 2-3 times as high as for energy-only calculations. A typical +input file of a RI-MP2 calculation using ADMM-accelerated HF is + +```none +&GLOBAL + PRINT_LEVEL LOW + PROJECT example + RUN_TYPE FORCE_EVAL +&END GLOBAL + +&FORCE_EVAL + METHOD Quickstep + &DFT + # This file contains basis sets for common elements (H, C, N, O; for F, Al, Cl, Si, P, S only TZ basis sets) + BASIS_SET_FILE_NAME BASIS_RI_cc-TZ + BASIS_SET_FILE_NAME BASIS_ADMM + POTENTIAL_FILE_NAME POTENTIAL + # Plug in your preconverged wfn file here + WFN_RESTART_FILE_NAME ./example_HF.wfn + # Use if no RI basis set is available (possible values: SMALL, MEDIUM, LARGE, HUGE) + AUTO_BASIS RI_AUX LARGE + + &AUXILIARY_DENSITY_MATRIX_METHOD + # Purification is not implemented + ADMM_PURIFICATION_METHOD NONE + # Try different options (check with reference values or use the default) + EXCH_CORRECTION_FUNC PBEX + # other methods are not implemented + METHOD BASIS_PROJECTION + &END AUXILIARY_DENSITY_MATRIX_METHOD + &MGRID + # Adjust as usual + CUTOFF 600 + REL_CUTOFF 50 + &END MGRID + &POISSON + PERIODIC XYZ + POISSON_SOLVER WAVELET + &END POISSON + &QS + # Choose as tight as necessary + EPS_DEFAULT 1.0E-10 + METHOD GPW + &END QS + &SCF + # Choose as tight as necessary + EPS_SCF 1.0E-8 + MAX_SCF 100 + SCF_GUESS RESTART + &END SCF + &XC + &HF + FRACTION 1.0000000 + &INTERACTION_POTENTIAL + # Adjust the cutoff radius according to your cell + CUTOFF_RADIUS 1.5 + POTENTIAL_TYPE TRUNCATED + T_C_G_DATA t_c_g.dat + &END INTERACTION_POTENTIAL + &SCREENING + # Tune these parameters + EPS_SCHWARZ 1.0E-10 + EPS_SCHWARZ_FORCES 1.0E-5 + SCREEN_ON_INITIAL_P .FALSE. + &END SCREENING + &END HF + &WF_CORRELATION + MEMORY 500 + NUMBER_PROC 1 + &CANONICAL_GRADIENTS + # The default is usually good enough + EPS_CANONICAL 1E-6 + # This is the option which should always work + # Try to set it to .TRUE. (may segfault) + FREE_HFX_BUFFER .FALSE. + &CPHF + # Choose as tight as necessary + EPS_CONV 1.0E-6 + # Smaller values of EPS_CONV may require more iterations + MAX_ITER 10 + &END CPHF + &END CANONICAL_GRADIENTS + &INTEGRALS + &WFC_GPW + # Adjust these parameters to your needs of accuracy + CUTOFF 200 + EPS_FILTER 1.0E-12 + EPS_GRID 1.0E-8 + REL_CUTOFF 50 + &END WFC_GPW + &END INTEGRALS + &RI_MP2 + # Determine an automatic block size, if memory is low, set it to 1 + BLOCK_SIZE -1 + &END RI_MP2 + &END WF_CORRELATION + &XC_FUNCTIONAL NONE + &END XC_FUNCTIONAL + &END XC + &END DFT + &PRINT + &FORCES + &END FORCES + &END PRINT + &SUBSYS + # Adjust cell information and coordinations as necessary + &CELL + ABC [angstrom] 5.0 5.0 5.0 + PERIODIC XYZ + &END CELL + &COORD + O 0.000000 0.000000 -0.211000 + H 0.000000 -0.844000 0.495000 + H 0.000000 0.744000 0.495000 + &END COORD + &KIND H + # orbital and RI basis set should match + BASIS_SET cc-TZ + BASIS_SET RI_AUX RI_TZ + # Use the largest affordable one + BASIS_SET AUX_FIT cpFIT3 + POTENTIAL GTH-HF-q1 + &END KIND + &KIND O + BASIS_SET cc-TZ + BASIS_SET RI_AUX RI_TZ + BASIS_SET AUX_FIT cpFIT3 + POTENTIAL GTH-HF-q6 + &END KIND + &TOPOLOGY + &CENTER_COORDINATES + &END CENTER_COORDINATES + &END TOPOLOGY + &END SUBSYS +&END FORCE_EVAL +``` + +A few more words are required regarding systems with many degenerate occupied orbital pairs as they +are found in structures with many symmetry-equivalent atoms. In that case, non-diagonal elements of +the density matrix of almost degenerate occupied-orbital pairs have to be calculated explicitly for +numerical reason. This is significantly more expensive than for non-degenerate orbital pairs. This +behavior can be tuned using the +[EPS_CANONICAL](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS.EPS_CANONICAL) +keyword. + ## Performance Considerations -MP2-calculations are generally very expensive for large systems due to their quintical scaling with respect to the number of atoms. Thus, for larger systems, one should switch to the cheaper RPA or SOS-MP2 methods. -CP2K employs local matrix-matrix multiplications for contractions. These can be accelerated with the SpLA library if CP2K was linked against it. In that way, GPU acceleration is possible if no GPU-accelerated DGEMM implementation is available and SpLA was configured accordingly. +MP2-calculations are generally very expensive for large systems due to their quintical scaling with +respect to the number of atoms. Thus, for larger systems, one should switch to the cheaper RPA or +SOS-MP2 methods. CP2K employs local matrix-matrix multiplications for contractions. These can be +accelerated with the SpLA library if CP2K was linked against it. In that way, GPU acceleration is +possible if no GPU-accelerated DGEMM implementation is available and SpLA was configured +accordingly. diff --git a/docs/methods/post_hartree_fock/preliminaries.md b/docs/methods/post_hartree_fock/preliminaries.md index 42bf44add6..8d37c5df38 100644 --- a/docs/methods/post_hartree_fock/preliminaries.md +++ b/docs/methods/post_hartree_fock/preliminaries.md @@ -1,48 +1,99 @@ # Preliminaries -Here, we explain the relevant preparation steps and things to consider with all post-Hartree-Fock methods regarding choice of basis sets and pseudopotentials. Because of the different building blocks, we have to distinguish at least two kinds of basis sets: the primary basis set (PBS) to represent the orbital functions and the RI basis set to expand products of PBS functions into a set of auxiliary functions such that the four-center electron repulsion integrals (in Mulliken-notation) are approximated by +Here, we explain the relevant preparation steps and things to consider with all post-Hartree-Fock +methods regarding choice of basis sets and pseudopotentials. Because of the different building +blocks, we have to distinguish at least two kinds of basis sets: the primary basis set (PBS) to +represent the orbital functions and the RI basis set to expand products of PBS functions into a set +of auxiliary functions such that the four-center electron repulsion integrals (in Mulliken-notation) +are approximated by $$ (ia|jb)=\sum_{PQ}(ia|P) (P|Q)^{-1} (Q|jb) $$ -On top of the PBS and the RI basis set, the ADMM approximation of HF may introduce a third kind of basis set (see there for further information). + +On top of the PBS and the RI basis set, the ADMM approximation of HF may introduce a third kind of +basis set (see there for further information). ## Step 1: Choice of Primary Basis Sets -Because of the slow convergence of the calculated properties from post-Hartree-Fock calculations with respect to the size of the PBS, post-Hartree-Fock methods demand for larger PBSs than ordinary hybrid DFT or HF calculations and for very accurate results a basis set convergence scheme employing consistent PBS for accurate results. For this purpose, post-Hartree-Fock calculations are commonly performed with correlation-consistent (cc)-basis sets. Without an extrapolation scheme, we recommend basis sets of at least triple-zeta (cc-TZ) or augmented double-zeta (aug-cc-DZ) quality. With an extrapolation scheme, one needs a second basis set for each atom kind of different size. The usual extrapolation formula for correlation energies is given by +Because of the slow convergence of the calculated properties from post-Hartree-Fock calculations +with respect to the size of the PBS, post-Hartree-Fock methods demand for larger PBSs than ordinary +hybrid DFT or HF calculations and for very accurate results a basis set convergence scheme employing +consistent PBS for accurate results. For this purpose, post-Hartree-Fock calculations are commonly +performed with correlation-consistent (cc)-basis sets. Without an extrapolation scheme, we recommend +basis sets of at least triple-zeta (cc-TZ) or augmented double-zeta (aug-cc-DZ) quality. With an +extrapolation scheme, one needs a second basis set for each atom kind of different size. The usual +extrapolation formula for correlation energies is given by $$ E(X) = A + \frac{B}{X^3} $$ -with $X$ being 2, 3 or 4 in case of DZ, TZ or QZ zeta basis sets, respectively. If possible, one should refrain from non-augmented DZ basis sets as their results are usually not converged enough for extrapolation. Suitable basis sets can be found in the `BASIS_RI_cc-TZ` and `BASIS_ccGRB` basis set files in the data directory. +with $X$ being 2, 3 or 4 in case of DZ, TZ or QZ zeta basis sets, respectively. If possible, one +should refrain from non-augmented DZ basis sets as their results are usually not converged enough +for extrapolation. Suitable basis sets can be found in the `BASIS_RI_cc-TZ` and `BASIS_ccGRB` basis +set files in the data directory. ## Step 2: Choice of RI Basis Sets (Only for RI-based Methods) -The RI approximation asks for an additional auxiliary basis set. They are usually specified as [BASIS_SET RI_AUX](#CP2K_INPUT.FORCE_EVAL.SUBSYS.KIND.BASIS_SET) in the input file. CP2K provides several options to find suitable RI basis sets: -- For a few main group elements, the `BASIS_RI_cc-TZ` file contains RI basis sets for the given PBSs. The literature may help with missing elements. -- Optimize RI basis sets yourself. This is recommended if a larger number of calculations is required to reduce the computational costs significantly. Check the [OPT_RI_BASIS](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI.OPT_RI_BASIS) and the regtest suite (QS/regtest-ri-opt) fur further information. -- In all other cases, CP2K can generate RI basis sets automatically (see the [AUTO_BASIS RI_AUX](#CP2K_INPUT.FORCE_EVAL.DFT.AUTO_BASIS). These basis sets provide a decent accuracy but are roughly twice as large as optimized basis sets. +The RI approximation asks for an additional auxiliary basis set. They are usually specified as +[BASIS_SET RI_AUX](#CP2K_INPUT.FORCE_EVAL.SUBSYS.KIND.BASIS_SET) in the input file. CP2K provides +several options to find suitable RI basis sets: + +- For a few main group elements, the `BASIS_RI_cc-TZ` file contains RI basis sets for the given + PBSs. The literature may help with missing elements. +- Optimize RI basis sets yourself. This is recommended if a larger number of calculations is + required to reduce the computational costs significantly. Check the + [OPT_RI_BASIS](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI.OPT_RI_BASIS) and the regtest suite + (QS/regtest-ri-opt) fur further information. +- In all other cases, CP2K can generate RI basis sets automatically (see the + [AUTO_BASIS RI_AUX](#CP2K_INPUT.FORCE_EVAL.DFT.AUTO_BASIS). These basis sets provide a decent + accuracy but are roughly twice as large as optimized basis sets. ## Step 3: Pseudopotentials -Pseudopotentials should reflect the underlying hybrid or HF calculation used to determine orbitals and orbital energies (see the respective sections). +Pseudopotentials should reflect the underlying hybrid or HF calculation used to determine orbitals +and orbital energies (see the respective sections). ## Step 4: Preoptimize Reference Orbitals -It is recommended to restart the orbital calculation step from preoptimized HF or hybrid orbitals. This reduces computational costs because DFT and HF calculations do not perform well in case of the large number of CPU cores required for post-Hartree-Fock methods. +It is recommended to restart the orbital calculation step from preoptimized HF or hybrid orbitals. +This reduces computational costs because DFT and HF calculations do not perform well in case of the +large number of CPU cores required for post-Hartree-Fock methods. ## Step 5: Setup of Integral Calculation Step -The integration step is setup with the [INTEGRALS](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS) section. Because most densities occurring in the formalisms of the post-Hartree-Fock methods have a zero net-charge, such that the Coulomb operator does not need to be truncated it is required by HF. CP2K has its own integration routines for post-Hartree-Fock methods exploiting this circumstance. The available integration methods (compare the [ERI_METHOD](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS.ERI_METHOD) keyword) are Gaussian-Plane-Wave integration (GPW), Minimax-Ewald integration (MME) and Obara Saika integration (OS). +The integration step is setup with the +[INTEGRALS](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS) section. Because most densities +occurring in the formalisms of the post-Hartree-Fock methods have a zero net-charge, such that the +Coulomb operator does not need to be truncated it is required by HF. CP2K has its own integration +routines for post-Hartree-Fock methods exploiting this circumstance. The available integration +methods (compare the [ERI_METHOD](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS.ERI_METHOD) +keyword) are Gaussian-Plane-Wave integration (GPW), Minimax-Ewald integration (MME) and Obara Saika +integration (OS). -GPW is fully supported, MME does not implement stress tensors and does not support short-ranged operators (truncated Coulomb, erfc-Coulomb, ...), OS implements gradients only in case of low scaling methods and does not support automatically generated RI basis sets. The GPW and MME integration methods are further configured in the respective [WFC_GPW](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS.WFC_GPW) and [ERI_MME](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS.ERI_MME) sections. In case of the GPW integration method, the respective CUTOFF parameters should be tuned as known from ordinary G(A)PW integrations in CP2K although the primary cutoff parameter can be chosen much smaller (150-300 Ry) then usual. +GPW is fully supported, MME does not implement stress tensors and does not support short-ranged +operators (truncated Coulomb, erfc-Coulomb, ...), OS implements gradients only in case of low +scaling methods and does not support automatically generated RI basis sets. The GPW and MME +integration methods are further configured in the respective +[WFC_GPW](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS.WFC_GPW) and +[ERI_MME](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.INTEGRALS.ERI_MME) sections. In case of the +GPW integration method, the respective CUTOFF parameters should be tuned as known from ordinary +G(A)PW integrations in CP2K although the primary cutoff parameter can be chosen much smaller +(150-300 Ry) then usual. -The RI metric can be set by the user with the [RI](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI) section. In case of low scaling methods (with activated [LOW_SCALING](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.LOW_SCALING) section), the metric is by default the overlap metric, else it will be the potential operator will be used. If the default is not requested, the RI section should be set. +The RI metric can be set by the user with the [RI](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI) +section. In case of low scaling methods (with activated +[LOW_SCALING](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.LOW_SCALING) section), the metric is by +default the overlap metric, else it will be the potential operator will be used. If the default is +not requested, the RI section should be set. ## General Structure of Post-Hartree-Fock Methods -The general structure is very similar. Each post-Hartree-Fock method has its own section which has to be activated to start the respective calculation. Each section has subsections to further configure the respective kind of calculation. Most post-Hartree-Fock calculations require the following structure +The general structure is very similar. Each post-Hartree-Fock method has its own section which has +to be activated to start the respective calculation. Each section has subsections to further +configure the respective kind of calculation. Most post-Hartree-Fock calculations require the +following structure ```none &WF_CORRELATION @@ -80,7 +131,12 @@ The general structure is very similar. Each post-Hartree-Fock method has its own ## Gradient Calculations -Gradient calculations are available for RI-MP2, RI-RPA and RI-SOS-MP2 calculations. In case of low-scaling implementations, consult the [CPHF](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.LOW_SCALING.CPHF) section in the LOW_SCALING section, otherwise the [CANONICAL_GRADIENTS](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS) section for further information. The setup of the CANONICAL_GRADIENTS section islooks like +Gradient calculations are available for RI-MP2, RI-RPA and RI-SOS-MP2 calculations. In case of +low-scaling implementations, consult the +[CPHF](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.LOW_SCALING.CPHF) section in the LOW_SCALING +section, otherwise the +[CANONICAL_GRADIENTS](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS) section for +further information. The setup of the CANONICAL_GRADIENTS section islooks like ```none &CANONICAL_GRADIENTS @@ -98,4 +154,8 @@ Gradient calculations are available for RI-MP2, RI-RPA and RI-SOS-MP2 calculatio ## Basis Set Superposition Error (BSSE) -Energy calculations using post-Hartree-Fock methods suffer from severe basis set superposition error (BSSE). In CP2K, this is enhanced by the many different basis sets in use. Check the [BSSE](#CP2K_INPUT.FORCE_EVAL.BSSE) section for further information on the automatic setup of BSSE calculations and the [GHOST](#CP2K_INPUT.FORCE_EVAL.SUBSYS.KIND.GHOST) for the manual setup of BSSE calculations. +Energy calculations using post-Hartree-Fock methods suffer from severe basis set superposition error +(BSSE). In CP2K, this is enhanced by the many different basis sets in use. Check the +[BSSE](#CP2K_INPUT.FORCE_EVAL.BSSE) section for further information on the automatic setup of BSSE +calculations and the [GHOST](#CP2K_INPUT.FORCE_EVAL.SUBSYS.KIND.GHOST) for the manual setup of BSSE +calculations. diff --git a/docs/methods/post_hartree_fock/rpa.md b/docs/methods/post_hartree_fock/rpa.md index a8d9b0c025..d90cc2d507 100644 --- a/docs/methods/post_hartree_fock/rpa.md +++ b/docs/methods/post_hartree_fock/rpa.md @@ -14,20 +14,31 @@ E^{LT-RI-SOS-MP2} &= - \int_{0}^{\infty}d\tau\text{Tr}\left(\overline{Q}(\tau)^2 \overline{Q}_{RS}(\tau) &= \sum_{ia}B_{iaR}e^{-\left(\epsilon_a-\epsilon_i\right)\tau}B_{iaS} $$ -as implemented in [](#DelBen2013). The integration is performed numerically. Double-hybrid functionals are available in CP2K. +as implemented in [](#DelBen2013). The integration is performed numerically. Double-hybrid +functionals are available in CP2K. -CP2K provides two different implementations: a quartically-scaling and a low-scaling cubically-scaling implementation. Different functionals require rescaling which is achieved using the keywords [SCALE_RPA](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.SCALE_RPA) in case of RPA and [SCALE_S](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_S) in case of LT-RI-SOS-MP2. +CP2K provides two different implementations: a quartically-scaling and a low-scaling +cubically-scaling implementation. Different functionals require rescaling which is achieved using +the keywords [SCALE_RPA](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.SCALE_RPA) in case of +RPA and [SCALE_S](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_S) in case of LT-RI-SOS-MP2. Here, we will focus on the quartically-scaling implementation only. ## RI-dRPA -CP2K implements two quadrature schemes: Clenshaw-Curtis and Minimax. The first requires 30-40 quadrature points, the latter 6-8 quadrature points. Minimax quadrature rules have to be preoptimized such that not all possible numbers of quadrature points are available. +CP2K implements two quadrature schemes: Clenshaw-Curtis and Minimax. The first requires 30-40 +quadrature points, the latter 6-8 quadrature points. Minimax quadrature rules have to be +preoptimized such that not all possible numbers of quadrature points are available. -RPA correlation energies are usually combined with exact exchange energies. In CP2K, this is available by activating the [HF](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.HF) section which is setup like an ordinary HF section. +RPA correlation energies are usually combined with exact exchange energies. In CP2K, this is +available by activating the [HF](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.HF) section +which is setup like an ordinary HF section. -Several Beyond-RPA schemes are available: The Renormalized Screened Exchange (RSE) correction, the Approximate Exchange Kernel (AXK) and the Second-Order Screened Exchange corrections (SOSEX). These are enabled with the [RSE](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.RSE) keyword and the [EXCHANGE_CORRECTION](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.EXCHANGE_CORRECTION) section. -This results into the following possible RPA section +Several Beyond-RPA schemes are available: The Renormalized Screened Exchange (RSE) correction, the +Approximate Exchange Kernel (AXK) and the Second-Order Screened Exchange corrections (SOSEX). These +are enabled with the [RSE](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.RSE) keyword and the +[EXCHANGE_CORRECTION](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.EXCHANGE_CORRECTION) +section. This results into the following possible RPA section ```none &RI_RPA @@ -52,12 +63,24 @@ This results into the following possible RPA section ## RI-dRPA Gradient Calculations -Analytical gradients are only available for RPA calculations using a minimax grid, but not for the beyond-RPA methods (RSE, AXK, SOSEX). The general setup is similar to RI-MP2 gradients. In addition, there are two more relevant keywords: [DOT_PRODUCT_BLKSIZE](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS.DOT_PRODUCT_BLKSIZE) and [MAX_PARALLEL_COMM](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS.MAX_PARALLEL_COMM). The first splits the contraction along the auxiliary index to improve numerical stability. By default, this feature is turned off. The second keyword determines the number of parallel communication channels for non-blocking communication. Larger numbers allow more overlap but increase the memory requirements. Larger values than 3 are commonly not necessary. +Analytical gradients are only available for RPA calculations using a minimax grid, but not for the +beyond-RPA methods (RSE, AXK, SOSEX).[](#Stein2024) The general setup is similar to RI-MP2 +gradients. In addition, there are two more relevant keywords: +[DOT_PRODUCT_BLKSIZE](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS.DOT_PRODUCT_BLKSIZE) +and +[MAX_PARALLEL_COMM](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS.MAX_PARALLEL_COMM). +The first splits the contraction along the auxiliary index to improve numerical stability. By +default, this feature is turned off. The second keyword determines the number of parallel +communication channels for non-blocking communication. Larger numbers allow more overlap but +increase the memory requirements. Larger values than 3 are commonly not necessary. ## LT-RI-SOS-MP2 -CP2K implements two quadrature schemes only the Minimax scheme requiring usually 6-8 quadrature points. Minimax quadrature rules have to be preoptimized such that not all possible numbers of quadrature points are available. Scaling of the energy contribution is possible with the [SCALE_S](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_S) keyword which is also used by MP2 calculations. -The input is simpler than in the RPA-case +CP2K implements two quadrature schemes only the Minimax scheme requiring usually 6-8 quadrature +points.[](#DelBen2013) Minimax quadrature rules have to be preoptimized such that not all possible +numbers of quadrature points are available. Scaling of the energy contribution is possible with the +[SCALE_S](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.SCALE_S) keyword which is also used by MP2 +calculations. The input is simpler than in the RPA-case ```none &RI_SOS_MP2 @@ -70,10 +93,25 @@ The input is simpler than in the RPA-case ## LT-RI-SOS-MP2 Gradient Calculations -Analytical gradients are available and work similarly to RI-MP2 calculations and RI-RPA calculations. As in case of RI-MP2, it makes use of the [EPS_CANONICAL](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS.EPS_CANONICAL) keyword. Larger values improve the numerical accuracy but increase the computational costs significantly because CP2K assumes that the number of relevant pairs is negligible. +Analytical gradients are available and work similarly to RI-MP2 calculations and RI-RPA +calculations. As in case of RI-MP2, it makes use of the +[EPS_CANONICAL](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.CANONICAL_GRADIENTS.EPS_CANONICAL) +keyword. Larger values improve the numerical accuracy but increase the computational costs +significantly because CP2K assumes that the number of relevant pairs is negligible. ## Performance Considerations -Without employing the low-scaling implementation, RI-RPA and RI-SOS-MP2 calculations scale quartically with respect to the number of atoms. The implementation relies on parallel matrix-matrix multiplications. Their costs can be reduced with the COSMA library and accelerated on GPUs if COSMA was configured accordingly. +Without employing the low-scaling implementation, RI-RPA and RI-SOS-MP2 calculations scale +quartically with respect to the number of atoms. The implementation relies on parallel matrix-matrix +multiplications. Their costs can be reduced with the COSMA library and accelerated on GPUs if COSMA +was configured accordingly. -HF calculations may be accelerated using ADMM which is recommended if very diffuse basis functions are employed. The RSE correction employs the same HF section as for the calculation of the HF exchange energy. The costs of AXK and SOSEX depend on the implementation in use. The HF-based implementation is cheaper for larger systems. Its performance depends on the available memory per process group (see [GROUP_SIZE](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.GROUP_SIZE) and [MAX_MEMORY](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.HF.MEMORY.MAX_MEMORY)). The costs becomes neglegible for large systems. Without the HF implementation, CP2K relies on SpLA-accelerated local matrix-matrix multiplications. These become the bottleneck for larger systems but introduce less numerical noise than the HF-based implementation. +HF calculations may be accelerated using ADMM which is recommended if very diffuse basis functions +are employed. The RSE correction employs the same HF section as for the calculation of the HF +exchange energy. The costs of AXK and SOSEX depend on the implementation in use. The HF-based +implementation is cheaper for larger systems. Its performance depends on the available memory per +process group (see [GROUP_SIZE](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.GROUP_SIZE) and +[MAX_MEMORY](#CP2K_INPUT.FORCE_EVAL.DFT.XC.WF_CORRELATION.RI_RPA.HF.MEMORY.MAX_MEMORY)). The costs +becomes neglegible for large systems. Without the HF implementation, CP2K relies on SpLA-accelerated +local matrix-matrix multiplications. These become the bottleneck for larger systems but introduce +less numerical noise than the HF-based implementation. diff --git a/src/common/bibliography.F b/src/common/bibliography.F index d34ca18336..f3c99e3bfd 100644 --- a/src/common/bibliography.F +++ b/src/common/bibliography.F @@ -91,7 +91,8 @@ MODULE bibliography Rycroft2009, Thomas2015, Brehm2018, Brehm2020, Shigeta2001, Heinecke2016, & Brehm2021, Bussy2021a, Bussy2021b, Ditler2021, Ditler2022, Mattiat2019, & Mattiat2022, Belleflamme2023, Knizia2013, Musaelian2023, Eriksen2020, & - Bussy2023, Wang2018, Zeng2023, Graml2024, Solca2024, Freeman1977, Gruneis2009 + Bussy2023, Wang2018, Zeng2023, Graml2024, Solca2024, Freeman1977, Gruneis2009, & + Stein2022, Stein2024 CONTAINS @@ -5031,7 +5032,7 @@ CONTAINS CALL add_reference(key=Gruneis2009, ISI_record=s2a( & "TY JOUR", & - "AU Grüneis, Andreas", & + "AU Grueneis, Andreas", & " Marsman, Martijn", & " Harl, Judith", & " Schimka, Laurids", & @@ -5050,6 +5051,44 @@ CONTAINS "ER"), & DOI="10.1063/1.3250347") + CALL add_reference(key=Stein2022, ISI_record=s2a( & + "TY JOUR", & + "AU Stein, Frederick", & + "AU Hutter, Juerg", & + "TI Double-hybrid density functionals for the condensed phase: ", & + " Gradients, stress tensor, and auxiliary-density matrix method acceleration", & + "PY 2022", & + "Y1 2022/02/17", & + "JO The Journal of Chemical Physics", & + "JA J. Chem. Phys.", & + "VL 156", & + "IS 7", & + "SP 074107", & + "SN 0021-9606", & + "Y2 6/13/2024", & + "UR https://doi.org/10.1063/5.0082327", & + "ER"), & + DOI="10.1063/5.0082327") + + CALL add_reference(key=Stein2024, ISI_record=s2a( & + "TY JOUR", & + "AU Stein, Frederick", & + "AU Hutter, Juerg", & + "TI Massively parallel implementation of gradients within the random ", & + " phase approximation: Application to the polymorphs of benzene", & + "PY 2024", & + "Y1 2024/01/12", & + "JO The Journal of Chemical Physics", & + "JA J. Chem. Phys.", & + "VL 160", & + "IS 2", & + "SP 024120", & + "SN 0021-9606", & + "Y2 6/13/2024", & + "UR https://doi.org/10.1063/5.0180704", & + "ER"), & + DOI="10.1063/5.0180704") + END SUBROUTINE add_all_references END MODULE bibliography diff --git a/src/input_cp2k_mp2.F b/src/input_cp2k_mp2.F index d9beac99e6..16def9ee8f 100644 --- a/src/input_cp2k_mp2.F +++ b/src/input_cp2k_mp2.F @@ -14,7 +14,7 @@ MODULE input_cp2k_mp2 USE bibliography, ONLY: & Bates2013, DelBen2012, DelBen2013, DelBen2015, DelBen2015b, Rybkin2016, Wilhelm2016a, & - Wilhelm2016b, Wilhelm2017, Wilhelm2018 + Wilhelm2016b, Wilhelm2017, Wilhelm2018, Stein2022, Stein2024, Bussy2023 USE cp_eri_mme_interface, ONLY: create_eri_mme_section USE cp_output_handling, ONLY: add_last_numeric, & cp_print_key_section_create, & @@ -81,7 +81,8 @@ CONTAINS "RI-MP2, RI-SOS-MP2, RI-RPA and GW (inside RI-RPA). ", & n_keywords=4, n_subsections=7, repeats=.TRUE., & citations=(/DelBen2012, DelBen2013, DelBen2015, DelBen2015b, Rybkin2016, & - Wilhelm2016a, Wilhelm2016b, Wilhelm2017, Wilhelm2018/)) + Wilhelm2016a, Wilhelm2016b, Wilhelm2017, Wilhelm2018, Stein2022, & + Stein2024, Bussy2023/)) NULLIFY (keyword, subsection) @@ -442,7 +443,8 @@ CONTAINS CALL section_create(section, __LOCATION__, name="CANONICAL_GRADIENTS", & description="Parameters influencing gradient calculations of canonical RI methods. "// & "Ignored if the IM_TIME section is set.", & - n_keywords=3, n_subsections=1, repeats=.FALSE.) + n_keywords=3, n_subsections=1, repeats=.FALSE., & + citations=(/DelBen2015b, Rybkin2016, Stein2022, Stein2024/)) NULLIFY (subsection, keyword) CALL create_cphf(subsection) @@ -1484,7 +1486,8 @@ CONTAINS description="Cubic scaling RI-RPA, GW and Laplace-SOS-MP2 method using the imaginary time formalism. "// & "EPS_GRID in WFC_GPW section controls accuracy / req. memory for 3-center integrals. "// & "SORT_BASIS EXP should be specified in DFT section.", & - n_keywords=12, n_subsections=2, repeats=.FALSE.) + n_keywords=12, n_subsections=2, repeats=.FALSE., & + citations=(/Wilhelm2016b, Wilhelm2018, Bussy2023/)) NULLIFY (keyword) CALL keyword_create(keyword, __LOCATION__, name="_SECTION_PARAMETERS_", & diff --git a/src/mp2.F b/src/mp2.F index 158dde9aac..33b00b209a 100644 --- a/src/mp2.F +++ b/src/mp2.F @@ -17,9 +17,12 @@ MODULE mp2 admm_uncorrect_for_eigenvalues USE atomic_kind_types, ONLY: atomic_kind_type,& get_atomic_kind_set - USE bibliography, ONLY: DelBen2012,& + USE bibliography, ONLY: Bussy2023,& + DelBen2012,& DelBen2015b,& Rybkin2016,& + Stein2022,& + Stein2024,& cite_reference USE cp_blacs_env, ONLY: cp_blacs_env_type USE cp_control_types, ONLY: dft_control_type @@ -239,6 +242,9 @@ CONTAINS IF (calc_forces) THEN CALL cite_reference(DelBen2015b) CALL cite_reference(Rybkin2016) + CALL cite_reference(Stein2022) + CALL cite_reference(Bussy2023) + CALL cite_reference(Stein2024) !Gradients available for RI-MP2, and low-scaling Laplace MP2/RPA IF (.NOT. (mp2_env%method == ri_mp2_method_gpw .OR. mp2_env%method == ri_rpa_method_gpw .OR. mp2_env%method == ri_mp2_laplace)) THEN CPABORT("No forces/gradients for the selected method.")