Meta-GGA Performance in Solids at Almost GGA Cost
aa r X i v : . [ phy s i c s . c h e m - ph ] A ug Meta-GGA Performance in Solids at Almost GGA Cost
Daniel Mej´ıa-Rodr´ıguez ∗ and S.B. Trickey † Center for Molecular Magnetic Quantum Materials, Quantum Theory Project,Department of Physics, University of Florida, Gainesville, FL 32611 (Dated: 27 Aug. 2020)A recent modification, r SCAN, of the SCAN (strongly constrained and appropriately normed)meta-GGA exchange-correlation functional mostly eliminates numerical instabilities and attendantintegration grid sensitivities exhibited by SCAN. Here we show that the successful deorbitalizationof SCAN to SCAN-L (SCAN with density Laplacian dependence) carries over directly to yieldr SCAN-L. A major benefit is that the high iteration counts that hindered use of SCAN-L areeliminated in r SCAN-L. It therefore is a computationally much faster meta-GGA than its orbital-dependent antecedent. Validation data for molecular heats of formation, bond lengths, and vibrationfrequencies (G3/99X, T96-R, T82-F test sets respectively) and on lattice constants, and cohesiveenergies (for 55 solids) and bulk moduli (for 40 solids) are provided. In addition, we show that theover-magnetization of bcc Fe from SCAN persists in r SCAN but does not appear in r SCAN-L,just as with SCAN-L.
Setting and Motivation - Recognition of chemicallysignificant electron density inhomogeneities by use ofan indicator function (usually denoted α ; see below) isthe critical mechanism by which a meta-GGA exchange-correlation (XC) functional improves upon a general-ized gradient approximation (GGA). The most success-ful meta-GGA so far (see [1] and references therein),is SCAN, the strongly constrained and appropriatelynormed functional [2, 3]. Its success is attributed to en-forcement upon it of all known rigorous constraints whicha meta-GGA can meet, together with calibration to theenergies of certain primitive physical systems (the “ap-propriate norms”; see Supplemental Material to Ref. [2]).Despite its successes, SCAN introduced a numericalproblem and exacerbated a methodological challenge.We defer discussion of the methodological issue briefly.The numerical problem has two elements. SCAN exhibitshigh sensitivity to numerical integration grid density.Handling that requires extremely dense, hence costly,grids. The other element is instability of self-consistentfield convergence that is hard to foresee, hence control,for a given system (especially in periodic solids).Those two numerical issues with SCAN were addressedby Bart´ok and Yates [4] by a simple renormalization ofthe denominator of α , a rescaling, and replacement of theSCAN switching function f ( α ) with a smoother seventh-degree polynomial for α < .
5. The resulting revisedSCAN (rSCAN) is far better behaved computationallythan SCAN. rSCAN preserves both the good molecularbond lengths and vibrational frequencies performance ofSCAN. Unfortunately, rSCAN does not preserve the goodperformance of SCAN for benchmark molecular heats offormation [5]. In periodic solids, SCAN and rSCAN areabout the same for lattice constants and cohesive energies[5] on a 55 solid test set [6] and for bulk moduli on a 44solid set [7].Very recently Furness et al. [8] have shown that theshortcomings of rSCAN stem from the fact that its regu- larization resulted in violation of several constraints satis-fied by SCAN. They adopted the smooth switching func-tion strategy of rSCAN combined with modifications torestore compliance with all but one of the constraintssatisfied in SCAN. The result, their regularized-restoredSCAN functional or r SCAN, recovers the strong perfor-mance trends of SCAN relative to molecular and soliddata sets while maintaining the numerical stability ofrSCAN.The r SCAN combination of accuracy and stabilityopens an opportunity for equally improved response tothe methodological challenge. That comes from the reg-ularized chemical region detector α ( r ) := τ s − τ W τ T F + ητ W . (1)Here τ s = (1 / P f j |∇ ϕ j ( r ) | is the positive-definiteKohn-Sham (KS) kinetic energy density in terms of theKS orbitals ϕ j and occupations f j , τ W and τ T F are thevon Weizs¨acker and Thomas-Fermi kinetic energy densi-ties respectively, and η = 10 − is a small regularizationparameter. The original α has η = 0. Obviously theorbital dependence of the XC energy introduced by α disqualifies SCAN or r SCAN from being used in orbital-free DFT. Worse, that orbital-dependence makes an ordi-nary KS calculation almost prohibitively costly becauseit necessitates an optimized effective potential calculation[9–12] at every SCF step. Usual practice to evade thatcost is to do generalized-KS (gKS) calculations. The gKSEuler equation follows from variation of the density func-tional approximation with respect to the orbitals, not thedensity. For meta-GGA and higher-rung functionals theKS and gKS equations are not equivalent [13, 14].We addressed this challenge by deorbitalization [15–17]. The deorbitalized version of SCAN, SCAN-L, dif-fers from SCAN only in using an approximate orbital-independent approximation for τ s in the original α to give α [ n, ∇ n, ∇ n ] (with n the electron number density). De-orbitalization restores use of the KS equation. Further-more, a SCAN-L calculation should be much faster thanSCAN. In practice that advantage oft-times went unreal-ized because numerical instabilities caused very slow SCFconvergence. Experience[18] suggested that the problemmight be rooted in the ∇ n dependence. By deorbitaliz-ing r SCAN into r SCAN-L, we show here that much ofthe problem actually was inherited from SCAN.
Procedure and Results -The deorbitalization of r SCAN used precisely thesame deorbitalized τ s form and parametrization as wasused for SCAN-L in Refs. [15, 16]. Molecular calculationswere done with a locally modified developers’ version ofthe NWChem code [19]. Similarly, the periodic solid cal-culations were done with a locally modified version ofVASP 5.4.4 [20, 21]. Note the remarks about coding inthe VASP meta-GGA trunk in Ref. [16]. As in that ref-erence, we did calculations with coding implemented inthat trunk (to check unambiguously against the VASPresults of Ref. [8]) and coding in the GGA trunk (to as-certain optimal speed-up from the deorbitalization). Ba-sis sets, cutoffs, and other matters of technique were asdocumented in [15, 16] with one exception, the PAWs,documented below.Table I summarizes the results for the molecular testsets in the form of mean absolute errors (MAEs) with re-spect to experiment for heats of formation (G3/99X set[22]), bond lengths (T96-R set [23]), and harmonic vi-brational frequencies (T82-F set [23]) obtained with theNWChem HUGE grid. For lower-density grids, Table Ishows the mean absolute deviation (MAD) and maxi-mum absolute deviation (MAX) with respect to the
HUGE grid results.Table I confirms the necessity of very dense numericalintegration grids for both SCAN and SCAN-L. Even the
XFINE preset grid yields deviations above 1 kcal/mol and30 cm − (bond lengths are less problematic). In contrast,r SCAN and r SCAN-L results are well-converged withthe
MEDIUM and
FINE grid presets. This is a major im-provement both in reliability as well as in performance,since numerical integration easily can become a computa-tional bottleneck. See below regarding calculation of ∇ n on the integration grid in the context of a Gaussian-typeorbital (GTO) basis.It is important to note the disentanglement of insta-bilities due to the functional form versus the Laplacian-dependence. It now is clear that SCAN-L exhibitsroughly the same stability difficulties as SCAN becausetheir common structure. However, the highly stabler SCAN-L shows thermochemical deviations about threetimes larger than those for r SCAN on a given presetgrid. This residual grid sensitivity seems directly at-tributable to the Laplacian-dependence of r SCAN-L.Different from the setup described in [16], the peri-odic solid calculations in VASP used the ”no-suffix” PAWdatasets instead of the GW ones. (Discussion of this choice is below.) Since the PAW datasets used hereare softer than the GW sets, we lowered the kineticenergy cutoff to 600 eV. The k-point sampling densitywas increased to match that reported in [8] by using the
KSPACING=0.1 keyword.Table II shows MAEs with respect to experimental re-sults for three crystalline test sets [16]: 55 equilibriumlattice constants (with cubic or hexagonal symmetry),40 bulk moduli (cubic symmetry), and 55 cohesive ener-gies. Zero-point effects were removed from experimentallattice parameters and bulk moduli.Similar to the molecular case, Table II shows that deor-bitalization of r SCAN is achieved with the same successas with SCAN. (Note that the values in that table are notdirectly comparable with our previous reports [5, 16] be-cause of the PAW dataset change.) In both cases, latticeparameters are well-treated. Predictions of bulk moduliand cohesive energies with both deorbitalized function-als show large percentage deviations, but the deviationmagnitudes are nonetheless quite small (large percent-age error in a small quantity). Some physics also is in-volved. Part of the cohesive energy MAE difference be-tween r SCAN-L and r SCAN comes from the differentelectronic configurations found for the W atom. r SCAN-L, SCAN, and SCAN-L all find an 6s valence whenthe configuration is unconstrained, while r SCAN findsthe correct 6s configuration.The total time needed to converge the 660 single-pointcalculations (12 per each solid) was used as a surrogatemeasure of the speed and stability of each functional.Consistent with expectations, the total times relative tothe SCAN benchmark were 0.924 for r SCAN, 0.438 forSCAN-L, 0.272 for r SCAN-L, and 0.227 for PBE [24].In other words, r SCAN-L computational cost in a plane-wave basis is almost as inexpensive as a standard GGAfunctional, even though numerical demands associatedwith the Laplacian-dependence remain.There is an important caveat. The SCF stability of allthe approximate functionals, as measured by the numberof iterations needed for convergence, can be degraded byuse of the GW PAW datasets. The Laplacian-dependentfunctionals are significantly worse in this regard thanthe orbital-dependent ones. What one sees with some,but not all, of the GW datasets is erratic SCF conver-gence. In those cases, near-SCF-convergence from itera-tion steps N − N often is followed by drastic wors-ening at step N + 1. Exploration of the origins of thisbehavior is outside the scope of the present work.Despite the fact that SCAN-L inherits many proper-ties from SCAN, SCAN-L avoids the over-magnetizationof simple elemental solids such as bcc Fe [25], We there-fore tested r SCAN and r SCAN-L in bcc Fe at the ex-perimental lattice constant (2.86 ˚A). As with SCAN andrSCAN, r SCAN predicts the magnetic moment of bccFe to be 2.63 µ B /atom. r SCAN-L lowers that to 2.27 µ B /atom, in line with other GGA functionals and SCAN- TABLE I. Performance of SCAN, r SCAN, and r SCAN-L for heats of formation (kcal/mol), bond lengths (˚A), and vibrationalfrequencies (cm − ) at various grid densities. Mean absolute errors with respect to experiment from the NWChem HUGE gridcalculations are in boldface . For the lower-density grids, mean absolute deviation and maximum absolute deviation (inparenthesis), with respect to those
HUGE results are shown.SCAN SCAN-L r SCAN r SCAN-LHeats of formation coarse medium fine xfine coarse medium fine xfine coarse medium fine xfine SCAN, and r SCAN-L XC functionals for the solid statetest sets. SCAN SCAN-L r SCAN r SCAN-LLattice parameters [˚A] 0.034 0.038 0.037 0.039Bulk moduli [GPa] 7.4 8.8 6.0 8.9Cohesive energies [eV/atom] 0.21 0.21 0.20 0.33 L. Summary - The slow (sometimes extremely so) SCFconvergence and numerical sensitivities of SCAN-L orig-inate mostly in structural characteristics of SCAN. Theremoval of those provided by r SCAN leads to a simi-larly well-behaved deorbitalized version, r SCAN-L. Ex-cept for the elemental 3d solid magnetization discrep-ancy, r SCAN-L replicates the behavior of r SCAN onmajor molecular and solid benchmarks. r SCAN-L ad-ditionally provides a pure KS treatment (hence enablesband-gap and optical excitation comparison with gKS re-sults from r SCAN), and should, in most cases supportsignificantly faster solid calculations than r SCAN, onthe time scale of an ordinary GGA. Further speedup ofmolecular calculations from r SCAN-L in a GTO basiswill require addressing the remaining computational bot-tleneck, calculation of ∇ n from GTO second-derivativeson the integration grid (rather than directly as in a plane-wave code.) Nonetheless calculations with r SCAN-Lwith the
MEDIUM
NWChem grid are as fast as those withSCAN with the
XFINE grid.This work was supported by U.S. Dept. of EnergyEnergy Frontier Research Center grant DE-SC 0019330.
SUPPLEMENTAL MATERIAL
Detailed molecular test set results obtained withr SCAN and r SCAN-L, as well as the detailed resultsfor the periodic solid test sets for all functionals, can befound in the supplementary material [26]. ∗ dmejiarodriguez@ufl.edu † [email protected]fl.edu[1] E. B. Isaacs and C. Wolverton, Phys. Rev. Mat. , 063801(2018).[2] J. Sun, A. Ruzsinszky, and J.P. Perdew, Phys. Rev. Lett. , 036402 (2015).[3] J. Sun, R.C. Remsing, Y. Zhang, Z. Sun, A. Ruzsinszky,H. Peng, Z. Yang, A. Paul, U. Waghmare, X. Wu, M.L.Klein, and J.P. Perdew, Nature Chemistry , 831 (2016).[4] A. P. Bart´ok and J. R. Yates, J. Chem. Phys. , 161101(2019).[5] D. Mej´ıa-Rodr´ıguez and S.B. Trickey, J. Chem. Phys. , 207101 (2019)].[6] H. Peng, Z.-H. Yang, J.P. Perdew and J. Sun, Phys. Rev.X , 041005 (2016). [7] F. Tran, J. Stelzl and P. Blaha, J. Chem. Phys. ,204120 (2016).[8] James W. Furness, Aaron D. Kaplan, Jinliang Ning, JohnP. Perdew, and Jianwei Sun, arXiv 2008.03374v1, 11Aug. 2020.[9] M. St¨adele, J.A. Majewski, P. Vogl and A. G¨orling, Phys.Rev. Lett. , 2089 (1997)[10] T. Grabo, T. Kreibich and E.K.U. Gross, Molec. Eng. ,27 (1997).[11] T. Grabo, T. Kreibach, S. Kurth, and E.K.U. Gross,in Strong Coulomb Correlations in Electronic Structure:Beyond the Local Density Approximation , V.I. Anisimoved. (Gordon and Breach, Tokyo, 2000) 203.[12] A. Heßelmann and A. G¨orling, Chem. Phys. Lett. ,110 (2008) and refs. therein.[13] Z.H. Yang, H. Peng, J. Sun, and J.P. Perdew, Phys. Rev.B , 205205 (2016)[14] J.P. Perdew, W. Yang, K. Burke, Z. Yang, E.K.U. Gross,M. Scheffler, G.E. Scuseria, T.M. Henderson, I.Y. Zhang,A. Ruzsinszky, H. Peng, J. Sun, E. Trushin, and A.G¨orling, Proc. Nat. Acad. Sci. (USA) , 2801 (2017).[15] D. Mej´ıa-Rodr´ıguez and S.B. Trickey, Phys. Rev. A ,052512 (2017).[16] D. Mej´ıa-Rodr´ıguez and S.B. Trickey, Phys. Rev. B ,115161 (2018).[17] F. Tran, P. Kovacs, L. Kalantari, G.K.H. Madsen, andP. Blaha, J. Chem. Phys. , 144105 (2018).[18] “Laplacian-based generalized gradient approximationsfor the exchange energy”, Antonio C. Cancio and ChrisE. Wagner, arXiv 1308.3744[19] E. Apr`a, E.J. Bylaska, W.A. de Jong, N. Govind, K.Kowalski, T.P. Straatsma, M. Valiev, H.J.J. van Dam, Y.Alexeev, J. Anchell, V. Anisimov, F.W. Aquino, R. Atta-Fynn, J. Autschbach, N.P. Bauman, J.C. Becca, D.E.Bernholdt, K. Bhaskaran-Nair, S. Bogatko, P. Borowski,J. Boschen, J. Brabec, A. Bruner, E. Cau¨et, Y. Chen,G.N. Chuev, C.J. Cramer, J. Daily, M.J.O. Deegan, T.H.Dunning Jr., M. Dupuis, K. G. Dyall, G.I. Fann, S.A. Fischer, A. Fonari, H. Fr´uchtl, L. Gagliardi, J. Garza,N. Gawande, S. Ghosh, K. Glaesemann, A. W. G¨otz,J. Hammond, V. Helms, E.D. Hermes, K. Hirao, S. Hi-rata, M. Jacquelin, L. Jensen, B.G. Johnson, H. J´onsson,R.A. Kendall, M. Klemm, R. Kobayashi, V. Konkov,S. Krishnamoorthy, M. Krishnan, Z. Lin, R.D. Lins,R.J. Littlefield, A.J. Logsdail, K. Lopata, W. Ma, A.V.Marenich, J. Mart´ın del Campo, D. Mej´ıa Rodr´ıguez,J.E. Moore, J.M. Mullin, T. Nakajima, D.R. Nascimento,J.A. Nichols, P.J. Nichols, J. Nieplocha, A. Otero-de-la-Roza, B. Palmer, A. Panyala, T. Pirojsirikul, B. Peng,R. Peverati, J. Pittner, L. Pollack, R.M. Richard, P. Sa-dayappan, G.C. Schatz, W.A. Shelton, D.W. Silverstein,D.M.A. Smith, T.A. Soares, D. Song, M. Swart, H.L.Taylor, G. S. Thomas, V. Tipparaju, D.G. Truhlar, K.Tsemekhman, T. Van Voorhis, ´A. V´azquez-Mayagoitia,P. Verma, O. Villa, A. Vishnu, K.D. Vogiatzis, D. Wang,J.H. Weare, M.J. Williamson, T.L. Windus, K. Wolin-ski, A.T. Wong, Q. Wu, C. Yang, Q. Yu, M. Zacharias,Z. Zhang, Y. Zhao, and R.J. Harrison, J. Chem. Phys. , 184102 (2020)[20] G. Kresse and J. Furthm¨uller, Phys. Rev. B , 11169(1996).[21] G. Kresse and D. Joubert, Phys. Rev. B , 1758 (1999).[22] L.A. Curtiss, K. Raghavachari, P.C. Redfern, andJ.A. Pople, J. Chem. Phys. , 1063 (1997); L.A. Cur-tiss, P.C. Redfern, K. Raghavachari, and J.A. Pople, J.Chem. Phys. , 108 (2001).[23] V.N. Staroverov, G.E. Scuseria, J. Tao, and J.P Perdew,J. Chem. Phys. , 3865 (1996); erratum ibid. , 1396 (1997).[25] D. Mej´ıa-Rodr´ıguez, and S.B. Trickey, Phys. Rev. B100