vdW-DF-ahcx: a range-separated van der Waals density functional hybrid

Hybrid density functionals replace a fraction of an underlying generalized-gradient approximation (GGA) exchange description with a Fock-exchange component. Range-separated hybrids (RSHs) also effectively screen the Fock-exchange component and thus open the door for characterizations of metals and adsorption at metal surfaces. The RSHs are traditionally based on a robust GGA, such as PBE [PRL $\textbf{77}$, 3865 (1996)], for example, as implemented in the HSE design [JPC $\textbf{118}$, 8207 (2003)]. Here we define a RSH extension to the van der Waals density functional (vdW-DF) method [ROPP $\textbf{78}$, 066501 (2015)], launching vdW-DF-ahcx. We use an analytical-hole (AH) framework [JCP $\textbf{128}$, 194105 (2008)] to characterize the GGA-type exchange in the vdW-DF-cx version [PRB $\textbf{89}$, 075148 (2014)], isolate the short-ranged exchange component, and define the new RSH. We find that the performance vdW-DF-ahcx compares favorably to (dispersion-corrected) HSE for descriptions of bulk (broad molecular) properties. We also find that it provides accurate descriptions of noble-metal surface properties, including CO adsorption.


I. INTRODUCTION
The search for a better computational-theory understanding of small-molecule/organics-substrate interfaces is directly motivated by technological and environmental challenges. There is, for example, a need for interface insight to improve catalysts, 1-4 batteries, [5][6][7][8][9] , gas adsorption 10,11 and photocurrent generation in organic solar cells. [12][13][14][15][16][17] However, molecular adsorption on metal and semiconductor surfaces challenge density functional theory (DFT). There is only one defining approximation, namely in the choice of the exchange-correlation (XC) energy functional. However, that choice determines accuracy, transferability, and hence usefulness of the DFT calculations. 18,19 The problem for DFT practitioners lies in finding an XC choice that is optimal for descriptions at both sides of the interfaces. 6,20,21 A van der Waals(vdW)-inclusive DFT approach is generally needed 21 and we have to go beyond the traditional generalized gradient approximations (GGAs). [23][24][25][26][27] Instead, one can use a ground-state energy functional with a dispersion correction, 4,28-41 a corresponding VV10based extension, [43][44][45][46][47] or move to the vdW density functional (vdW-DF) method. [48][49][50][51][52][53][54][55][56][57][58][59][60][61] The latter aims to track truly nonlocal correlation effects on the ground-state electron-density footing that DFT use to describe other types of interactions. There are many early examples of use for adsorption characterizations  because the method grew out of a fruitful feedback loop involving many-body-perturbation theory (MBPT) and surface science. 54, Meanwhile, use of a hybrid, i.e., mixing in a small fraction of Fock exchange, is also desirable. This is because hybrids are expected to improve the accounts of charge transfers. 19,[115][116][117][118] In this paper we introduce a new range-separated hybrid (RSH), termed vdW-DF-ahcx, based on the consistent-exchange CX version, and therefore having the same correlation description as in the first generalgeometry vdW-DF version. 51,52 Our new nonlocalcorrelation RSH is named to emphasize that it constitutes an analytical-hole 119 (AH) consistent-exchange (AHCX) formulation of in the vdW-DF family. It will be abbreviated AHCX below. The design starts with an analysis of the CX exchange design, termed cx13 or LV-rPW86, 56 that reflects the Lindhard screening logic and ensures current conservation. 57,60 The design is similar to that of HSE 116,117 and HSEsol 119,128 in that it focuses on isolating the short-range (SR) components, E CX,SR x and E SR FX , of cx13 and of Fock exchange, respectively. As in HSE, we rely on an error-function separation, 1 r 12 = erfc(γr 12 ) r 12 + erf(γr 12 ) of the Coulomb matrix elements for electrons separated by r 12 = |r 1 − r 2 |. However, the HSE 116,117 design was based on the EP model of the PBE exchange hole. 132 Here, we rely on the Henderson-Janesko-Scuseria (HJS) AH framework 119 for representing exchange in density functional approximations. This approach was also taken in the definition of HSEsol, below discussed as HJS-PBEsol. 119,128 There are also HJS-based rangeseparated screened hybrids that include the long-range Fock exchange. 43,130,131,133,134 We discuss AH formulations of PBEx 26 (in PBE), PBEsolx 27 (in PBEsol), and cx13 56 (in CX), to set the stage for defining AHCX. The error-function separation ensures that AHCX screens the long-ranged component of Fock exchange. The AHCX is therefore set up to even describe metals and organic-molecular binding at metal (and highdielectric-semiconductor) surfaces. The new RSH formulation aims to make the vdW-DF method available for bulk-and surface-application problems that require a more accurate description of charge transfer. The broad implementations of the adaptively compressed exchange (ACE) description 3,135 in DFT code packages, such as Quantum Espresso, 137,138 means that this strategy is also computationally feasible.
In fact, we find that AHCX is as fast as HSE for calculations of bulk and molecule properties. We document this for two bulk cases within, while we here provide an example comparison for medium-sized molecules, namely timing information extracted from our study of the C60ISO benchmark set on fullerene isomerizations, (Ref. 5 includes a presentation of the set). These are just 10 out of the roughly 2300 molecular problems (investigated in multiple functionals) that are part of this AHCX launching work, but they give an impression. We study these C 60 systems in cubic unit cells of length 18.6Å using the ONCV-SG15 1,2 pseudopotentials (PPs) at a 160 Ry wavefunction-energy cutoff. We average the total CPU-core-hour cost for completion across the ten geometries (although excluding a case where the Grimme-dispersion-corrected HSE-D3 5,39 took an exceptional long convergence path) finding these timing results: PBE at 6 core hours, CX at 13 core hours, HSE-D3 at 568 core hours, and AHCX at 434 core hours.
We note that the overall time consumption for hybrid studies of molecules is dominated by the Quantum Espresso ACE initialization. We also note that we have separately documented that the evaluation of the nonlocal correlation energy (used in AHCX and CX) scales well with cores and system size up to at least 10000 atoms, Ref. 141, and will not be a relevant bottle neck, at least not for a long time coming.
We motivate the AHCX design as a robust trulynonlocal correlation RSH through the quality of the AH exchange description that we supply for cx13, that is, for the GGA-type exchange in CX. We argue that we can port the AH exchange description 119 from PBEx 26 and PBEsolx 27 to cx13, which is constructed as a Padé interpolation 56 (LV-rPW86) of the Langreth-Vosko (LV) exchange 99,101 and of the revised PW86 exchange. 25,142 We note that related charge-conservation criteria are used in PBE/PBEsol and in the CX designs, and that they all reflect the MBPT analysis of exchange in the weakly perturbed electron gas. 23,54,99,101,[143][144][145] The rest of the paper is organized as follows. Section II summarizes the AH description of robust GGA exchange formulations, while Sec. III converts that insight into defining the AHCX. Section IV has computational details, a test of our AH-analysis approach, and a discussion of AHCX costs for bulk and other studies. Section V contains results and discussions, i.e., a documentation of performance for broad molecular properties, for bulk structure and cohesion, and for noble-metal surface properties, including CO adsorption. Sec. VI contains our summary and outlook. Finally, the appendix documents robustness of the adsorption-site-preference results with changes in the PP choice while the supplementary information (SI) material provides details of performance characterizations.

II. GGA EXCHANGE AND EXCHANGE HOLES
The local Fermi wavevector k F (r) = (2π 2 n(r)) 1/3 defines the (Slater) exchange energy density in the local density approximation (LDA), ε LDA x (r) = −(3/4π)k F (r). In a GGA, a local energy-per particle term also defines the exchange part of the XC energy functional, i.e., the exchange energy functional E x [n] = r n(r) ε GGA x (n(r); s(r)) . ( The GGA energy-per particle description ε GGA x (n(r); s(r)) depends on the electron density n(r) and the scaled density gradient s(r) = |∇n|/(2k F (r)n(r)). The ratio between the GGA and LDA energy-per-particle expressions defines the GGA exchange enhancement factor As indicated, this enhancement factor can exclusively depend on s, and the homogeneous electron gas (HEG) description is recouped by enforcing the limit value F GGA x (s → 0) = 1. The exchange hole n x (r; r ) represents the tendency for an electron at position r to inhibit occupation of an electron of the same spin at a neighboring point r . It is part of the total XC hole n xc (r; r ) which, in turn, is a full representation of all XC energy effects This XC hole is defined via the screened electrodynamical response in the electron gas. 23,51,54,60,96,99,101,[143][144][145][146][147] It can be computed from the adiabatic connection formula 145,146 (ACF) and one can obtain a full specification the density-density fluctuations for the homogeneous electron gas (HEG) by MBPT or from quantum Monte Carlo calculations. 146,[148][149][150][151][152] Such calculations can, in turn, be used to establish a model for the HEG XC hole, n HEG xc , that takes a weighted or modified Gaussian form 153,154 and defines LDA. 145,146,152 Through the HEG, one can also extract and separately analyze the LDA exchange hole using a Gaussian model form. 119,132,152 The exchange hole for general electron-density distributions,ñ can, in principle, be computed from inserting singleparticle orbitals in a one-particle density matrix Using thisñ x form in Eq. (4) gives the Fock-exchange result However, a direct use of this Fock exchange in combination with the correlation parts of either LDA, GGAs or vdW-DFs is not desirable, for example, because it gives divergences in the description of extended metallic systems. In fact, a use of Fock exchange is also inappropriate for molecules and insulators, because electrons respond to and therefore screen the Coulomb field; The impact can be substantial also for molecules. 127 Even a partial inclusion of the Fock-exchange description, Eq. (7), must (in general) be both compensated by correlation 18,23,24,51,123,131,145,146,[155][156][157][158] and described at an appropriate non-unity value of the dielectric constant. 127,[129][130][131]159 Starting with the LDA description, one seeks instead an exchange-hole description n x , and corresponding exchange energy functionals, Eq. (4), in which some XC cancellation is already built in. The LDA, GGA, and vdW-DF descriptions for n x (r; r ), and hence for E x , differ from the exchange description given by the Fock expression, Eqs. (5)- (7). In constraint-based GGA, the assumption of modified Gaussian hole form is adapted (from the LDA start). 119,132,160 This assumption also enters in a key role in the vdW-DF method by setting details of the underlying plasmon-response description. 51,[55][56][57][58]60 Since the GGA is given by the local value of the density n and of the scaled gradient s, the GGA exchange hole must also be approximated in those terms and we introduce a dimensionless representation 25,132,142,160 Here the separation |r − r | is scaled by the local value of the Fermi wavevector k F (r). The correspondingly scaled LDA exchange-hole form, J HEG x (s), must arise as the proper s → 0 limit of Eq. (8). For actual DFT calculations in the GGA, we need the resulting exchange enhancement factor. We get it by inserting the form of J GGA x (n(r), s(r)) into Eq. (4). The formal relation to the exchange energy E xc is given by the enhancement factor Perdew and co-workers furthermore defined a model for the GGA exchange hole, J GGA x (n(r), s(r)). The weighted-Gaussian-hole form closely resembles that of J HEG x , 152 but terms and exponents are now made functions of the scaled density gradient s. 25,27,142,160 Parameters are set by charge conservation, constraints, and physics arguments. Such GGA models for the exchange hole have been defined for PW86/rPW86 exchange and PBE and PBEsol exchange descriptions (denoted PBEx and PBEsolx). The details in setting J GGA x (n(r), s(r)) produce different scaled-gradient enhancement factors, Eq. (9), and hence different exchange energy functional descriptions, rPW86, PBEx, PBEsolx.
Connected with the PBE design, Ernzerhof and Perdew (EP) also designed an oscillation-free exchangehole form, 132 that simplified the definition of HSE as a PBE-based RSH. 116,117 Signatures of Friedel oscillations are present in, for example, the PBE exchange hole (in the s → 0 limit), 160 but they are not considered physical for descriptions of molecules. EP started by a slight modification of the LDA description to extract a nonoscillatory LDA-exchange hole form J LDA EP,x (s) and then repeated the constraint-based GGA-exchange design to craft the oscillation-free exchange hole form J GGA EP,x . It has been used to analyze both PBEx and PBEsolx. 132,161 This EP representation for PBEx was used with rangeseparation, Eq. (1), to establish HSE. 116,117 Characterizing the nature of the scaled exchange hole, Eq. (8) is equally relevant for the exchange components of vdW-DFs. This follows because the vdW-DF method presently relies on a GGA-type exchange, for example, cx13 (LV-rPW86) in the CX release.
A. HJS model of PBE and PBEsol exchange holes Figure 1 shows exchange hole n x descriptions for the semilocal PBE and PBEsol functionals, as well as for the truly nonlocal-correlation CX functional. The exchange part of the functionals are denoted PBEx, PBEsolx, and cx13. The hole shapes are here represented in the HJS or AH model framework for exchange 119 , using a form J GGA x,HJS (n(r), s(r)). However, the plots for PBEx and PBEsolx exchange can be directly compared with the EPbased analysis, J GGA x,EP (n(r), s(r)), 132 shown in Fig. 2 of Ref. 132 and Fig. 3 of Ref. 161; There are no discernible differences. There are also just small differences from the original PBEx and PBEsolx hole representations (apart from the removal of Friedel-oscillation signatures).
We work with the HJS framework for exchange holes because it gives several advantages. We call it an AH model of GGA exchange because it permits an explicit evaluation of Eq. (9) once J GGA x,HJS (n(r), s(r)) is inserted in Eq. (9). The evaluation is stated in Ref. 119 and it also permits a straightforward definition of functional derivatives. This is a clear advantage when the AH analysis is used to also define and code (as summarized in the following section) a RSH hybrid HJS-PBE 119 that mirrors HSE. The new design simply involves changing from the J GGA x,EP form to the J GGA x,HJS form. A key AH advantage is ease of generalization. It follows from having an analytical evaluation of Eq. (9). For example, HJS immediately established a plausible exchange-hole shape that reflects the PBEsol exchange enhancement factor. 119 In turn, this AH determination led to the definition of HJS-PBEsol, 119,128 that is, a RSH based on PBEsol. We shall use the AH-model for seeking generalization to the vdW-DF method in the following subsection.
To begin, we summarize the HJS AH framework for characterizing the PBE and PBEsol exchange nature. HJS revisits the mathematical framework for representing an LDA exchange hole TheF(s) function (affecting the y 2 term of the LDA That is, it is given by the shape of H(s). Meanwhile, the finalḠ(s) modification (off of the LDA description) is set by an overall charge conservation criteria, 23,25,[145][146][147] i.e., a limit on the integral of the exchange hole. The important observation is that all components of the AH exchange description, Eq. (12), are completely defined by the shape of H(s). The relation betweenḠ(s) and H(s) can be explicitly stated in the HJS AH exchange model, see Eq. (40) of Ref. 119. Overall, the gradient-corrections provide, effectively, an extra Gaussian damping (defined by H(s)), compared to the LDA description. The suppression arises at large values of the scaled distance y. This damping is, however, offset by an enhancement of the exchange hole at and around y = 1 (with the enhancement growing significantly with s). Nevertheless, all exchange hole components are completely set once we pick a plausible form for H(s). In practice, the HJS AH modeling, for a given exchange of the GGA family, proceeds in 3 steps. First, we assume a rational form of the Gaussian damping function and make a formal evaluation of Eq. (9): Here ζ = s 2 H(s), η =Ā+ζ, and λ = D+ζ. Second, since the GGA exchange is fully specified by the enhancement factor F x (s), we can invert Eq. (15) to get a numerical representation of H(s). Finally, we fit that numerical determination to the rational form, Eq. (14). This HJS procedure formally works for any GGA exchange -but it is an implicit assumption that we can still view the AH form J as reflecting an implicit, constrained maximum-entropy principle. 26,132,152,160 The Gaussian damping form, H(s), must reproduce a given GGA-type exchange enhancement factor F x (s) while using as little structure as possible; There is no physics rationale for keeping such structure. 132,160 In a HJS representation, Eq. (14), we must strive to use as few significant polynomial coefficients as possible. We must also check that we have, indeed, avoided fluctuations that could hide nonphysics input.
The left and middle panel of Fig. 1 contrast the dependence on scaled density gradient s values of our AH models for PBEx and PBEsolx, as refitted here by us.
The AH model for PBEsol exchange has a softer variation with the scaled separation y and with the scaled gradient s. In contrast, the AH for PBE exchange deepens more quickly. The GGA hole formulation is best motivated when the value of the scaled hole form J remains larger than −1 (as it reflects a suppression of the electron occupation). The PBE and PBEsol exchange hole descriptions comply with that criteria for s < 2.5, i.e., in regions that are expected to dominate the description of bonding in materials and often also among molecules. 20,26,27,56,60,160,162 The bottom panel of Fig. (2) shows the shape of the exponential suppression or damping factor H(s) for PBEx and PBEsolx (and for the cx13 exchange functional that is relevant for CX). For the PBEx analysis there is a close overlap with the exchange hole described in the original EP model, Ref. 132, as also observed in the HJS paper, Ref. 119. Similarly, the PBEsolx hole matches the EPmodel-based representation asserted in Ref. 161.
The AH hole characterizations should be seen as trusted models of the PBEx and PBEsolx hole. The HJS AH analysis of PBEx and PBEsolx also have nearly identical forms: Alone the scripted parameters denoted by an over bar in Eq. (12) differs slightly in their formal expression. There are difference in the precise formulation of H(s). However, the resulting AH descriptions of PBEx and PBEsolx can claim an important maximum entropy status because of the close similarity with the EP descriptions.
An important difference between the EP and the HJS exchange-hole descriptions lies in how these models are being used. The EP model 132 was introduced to discuss the enhancement factor PBEx in terms of a nonoscillatory exchange hole, and it led to HSE. 116,117 The HJS model 119 allows for a simpler PBE discussion, but it is also being used for reverse engineering, i.e., being used to assert a plausible exchange-hole shape from the exchange enhancement factor. This track gave rise to both HJS-PBEsol 119,128 and several recent long-ranged corrected hybrids. 134

B. Exchange-hole models for vdW-DFs
This paper formulates the AHCX RSH (below) from an expectation that the HJS procedure (for reverse engineering an exchange hole form) also remains valid for both rPW86 and the LV exchange descriptions. This is plausible because cx13 is formed as a Padé interpolation of the LV exchange 54,99,101 (at small to medium values of the scaled density gradient s) and of rPW86 25,142 (at large s values).
We first note that the PW86 exchange paper, Ref. 25, introduced a design strategy for setting the exchange hole within the electron gas tradition. The procedure involves four steps: 1) Establish the symmetry limits on how the density gradient modifies the exchange hole off of the well-understood LDA representations; 2) establish the small-s variation, for example, from MBPT results as in LV exchange; 3) impose an overall exchange-hole charge conservation criteria; and 4) extract the exchange enhancement factor from Eq. (9). Input beyond such physics analysis is always minimized in the electron-gas tradition that defines LDA, PW86/rPW86 exchange as well as both PBE, PBEsol and the vdW-DF method.
We focus our wider AH analysis on the exchange choices that are used in the Chalmers-Rutgers vdW-DF releases. The exchange forms are the revPBEx (used in vdW-DF1), the rPW86 (used in vdW-DF2), and cx13 (used in CX). One of these, rPW86, 142 is a refit of the PW86 that defined the strategy for making constraintbased GGA exchange design off of an model exchange hole. [25][26][27]116,119,128 We shall give a detailed motivation for trusting the HJS-type AH exchange characterizations for rPW86 and cx13, below.
We have fitted the parameters a 2−7 and b 0−9 of Eq. (14) for the AH representations of the cx13, revPBEx, and rPW86. This extends our PBEx and PBEsolx discussion in the previous subsection. Additional information is available in the SI material. Table I summarizes the most important such parameterization results: Our AH characterization for cx13 exchange and our AH refits of the HJS descriptions for PBEx and PBEsolx.
The right panel of Fig. 1 provides a practical motivation for trusting our AH analysis of the cx13 (LV-rPW96) exchange hole. The argument is given by noting similarities to the exchange-hole shapes of well-established exchange functionals. The cx13 shape begins (at low s) by reflecting a PBEsol nature and eventually it rolls over to the PBE-type (and rPW86-type) variation. Values greater than s = 2.5 lead to the deeper hole minima, as in PBE. The cross over is not surprising since the PBEsol builds on a MBPT analysis that is close to the LV description (with small exchange enhancements up to medium s values) while the rPW86 is known to have deep exchange holes. Figure 2 summarizes the exchange enhancements (top panel) and our AH analysis of corresponding exchange holes (bottom panels). The top panel confirms, in terms of exchange enhancements, that cx13 starts with a PBEsolx-(and PBEx-)like behavior at small s values but transforms to the rPW86 behavior at large s values.
A more formal argument for trusting the AH analysis of cx13 exchange can be stated as follows. We note that the CX leverages the formally exact electrodynamicresponse framework of the vdW-DF method, as far as it is possible. 60 The CX explicitly enforces current conservation up to the cross over in the cx13 (or LV-rPW86) Padé construction, namely at s ≈ 2.5. As such, for the lower-s LV end (up to 2-3), CX is an example of a consistent vdW-DF, 60 systematically relying on a plasmonbased response description that adheres to all known constraints and sum rules. 51,58 At the LV end, cx13 furthermore leverages the Lindhard screening logic to balance exchange and correlation. 60 There is, in fact, a strong formal connection to the PBE constraint-based design logic in that the CX emphasis on current conservation ensures an automatic compliance with the charge-conservation criteria. 60 The cx13 or LV-rPW86 is a well-motivated electron-gas construction at this low-to-mid-s LV end. Here the cx-13 also satisfies the PBE-type local implementation 26 of the Lieb-Oxford bound, 163 systematically having exchangeenhancement factors F x (s) < 1.804. For s < 2.5, the cx13 (and hence CX design) complies with the criteria on the exchange hole depth, J > −1.
At the high-s end, the cx13 exchange design rolls over into the rPW86. 142 It is there naturally within the realm of constraint-based exchange descriptions that started with PW86, Ref. 25. However, in this end, we must also discuss compliance with the actual, globally implemented, Lieb-Oxford bound. 163 This is because the PW86 was designed before the importance of the bound was understood. 26,27,160 For any given ground-state electron density n, we consider the ratio R of the total exchange to the total LDA exchange (in the unit cell). Rigorous bounds are R < 1.804 (as is hard wired at the local level in PBE 26 ) and R < 1.174 for two-electron systems (as is hard wired at the local level of SCAN 46,164 ). We note that there are unphysical (spherical-shell) systems where the density is constructed so that the scaled density gradient is constant and where the actual Lieb-Oxford bound can only be satisfied if implemented also at the local level (as in PBE and SCAN). 165 For such unphysical systems, the rPW86 and cx13 (and hence CX) will fail to comply with this exchange condition, 163,165 if we furthermore assume that the unphysical system has a large scaleddensity value, s > 3.
The important question remains, however, whether the cx13 (and hence CX) violates the actual Lieb-Oxford bound 163 in real systems, that is, as encountered in actual DFT studies. 166 The point is, that the Lieb-Oxford bound is formulated globally and must be checked on the unit-cell level in our periodic-system calculations. 163 Answering this question is straightforward for any given problem, as long as one saves the density after completing a Quantum Espresso calculation. Our now updated ppACF post-processing code (launched in Refs. 123 and 124 and committed to the Quantum Espresso package) gives a general coupling-constant analysis, an evaluation of kinetic-correlation energy components, as well as a per-system determination of actual exchange R ratio (for nonhybrids). For the 10 fullerene calculations discussed in the introduction, we find that the cx13/CX value for the actual exchange ratio R never exceeds 1.025. In practice, the CX exchange ratio remains far below even the stringent bound (1.174) that exists for two-electron systems. 164 The bottom panel of Fig. 2 provides further details of the cx13 design. It does so by comparing the Gaussian suppression factors H(s) that arise in our AH analysis of constrained exchange functionals. The shapes of the Gaussian suppression H(s) are similar in all cases. The s position of the maximum for cx13 (LV-rPW86) sits essentially on the PBEsolx position, slightly above the maximum position for rPW86 and PBEx. Meanwhile, the cx13 maximum value aligns with that of the PBE description. As expected, the asymptotic cx13 behavior coincides with that of the rPW86 while for low s values, the Gaussian suppression is almost identical to that which characterizes PBEsol. This is also expected as the PBEsolx was explicitly designed to move the exchange enhancement closer to the input that also defines LV exchange. That is, the PBEsolx is closer to the MBPT analysis of exchange in the weakly perturbed electron gas. 27,54,60,99,101 Overall, Fig. 2 confirms that there are no wild features in the AH characterizations of exchange in the vdW-DF releases, including CX. Moreover, the bottom panel confirms our assumption, that our AH model description of cx13 can be trusted. This follows because we find that it reflects a mixture of PBEsol-like, PBE-like, and rPW86-like behaviors for the Gaussian suppression factor of the hole. In fact, like PBE, the cx13 aims to serve as a compromise of staying close to MBPT results at low s (good for solids) and a more rapid enhancement rise (as in rPW86) for larger s values (good for descriptions of molecular binding energies 142,167 ).

III. RANGE-SEPARATED HYBRIDS
Hybrid functionals build on an underlying regular (density explicit) functional for XC energy. They simply replace some fraction α of the exchange description, E x [n], of that functional with a corresponding component extracted from the Fock-exchange term E FX [n]. The latter is evaluated from Kohn-Sham orbitals and calculations are carried to consistency in DFT. We note that the regular vdW-DFs all have a GGA-type exchange by design and the vdW-DF hybrid design can therefore be captured in the same overall discussion.
Simple, unscreened hybrids functional are described by the exchange component, 115,157,158,168,169 E hyb while the correlation term is kept unchanged. Here E FX [n] is evaluated as the Fock interaction term (7). Examples of such simple hybrids are PBE0 and the vdW-DF0 class. 122,123 The CX0p results when we also use a coupling-constant scaling analysis [155][156][157][158]169 to establish a plausible average value, α = 0.2, of the Fock-exchange mixing in the vdW-DF-cx0 design. 123 For general RSH designs we split both the Fock exchange and the functional exchange into short range (SR) and long range (LR) parts: We simply insert the error-function separation, Eq. (1), of the Coulomb matrix elements into Eq. (5) to extract, for example, E SR FX . Given a trusted AH exchange representation of the underlying functional 'DF', we can also extract the SR exchange component Use of the EP model 132 for a characterization of the PBE exchange hole, in combination with the range-separation Eq. (1), led directly to the design of the first RSH, namely the HSE. 116,117 More generally, for a trusted density functional 'DF' (that have a GGA-type exchange functional part E GGA x ,) we arrive at a HSE-type RSH or screened hybrid extension using where E DF c is the correlation part of 'DF'. This correlation can be semi-or truly nonlocal in nature. The simple recast brings out similarities with the design of unscreened hybrids. For PBE and PBEsol there are already trusted AH descriptions in the HJS-PBE 119 and HJS-PBEsol 119,128 formulations (besides the original HSE 116,117 obtained with the EP model-hole framework 132 ). To these we here add our own formulations of these analytical-hole screened hybrids, termed PBE-AH and PBEsol-AH. We do this to allow independent checks on our approach, as detailed in the following section.
The main results of this paper are the definition and launch of the AHCX RSH. It is built from our AH exchange description for CX, again, as summarized in Table  I and Figs. 1 and 2. The key observation is that we have trust in using the HJS type AH model of the CX exchange, so that we can rely on Eq. (20) in establishing an approximation for E SR,CX x . We set the Fock exchange fraction at α = 0.2 as in the CX0p. 123 Still, the AHCX differs from CX0p by the presence of the inverse screening length γ and thus a focus on correcting exclusively the short-range exchange description. We set γ = 0.106 (atomic units) as in the HSE06 formulation. 117 The AHCX is a vdW-DF RSH that is constructed in the electron gas tradition, being an all-from-ground-state density design. 19,57 That is, it uses one and the same plasmon pole model 56,60 for all but the Fock-exchange term. Even the Fock mixing can be seen as being set by the coupling-constant scaling of the consistentexchange CX version, and thus given from within the construction. 123,124 The AHCX is computationally more costly than meta-GGAs 4,46,164,170,171 and DFT+U. 172 These are other approaches that can also improve orbital descriptions and compensate for charge-transfer errors.
The SCAN functional 164 is perhaps the alternative that is closest in nature to CX and AHCX: they are all functionals of the electron-gas tradition and SCAN is constructed to retain some account of vdW interactions at intermediate distances 164 (thanks to input from a formal analysis of two-electron systems).
There are pros and cons of SCAN, CX, and AHCX use. SCAN and CX are certainly faster than AHCX. However, SCAN must be supplemented by a separately-defined semi-empirical addition, in SCAN+rVV10, 46 to capture nonlocal-correlation effects across separations. In contrast, the AHCX is set up (as a single XC functional design) to capture general interactions, for example, across the range of fragment separations at organics-metal interfaces. One could extend the AHCX framework for use with optical or MBPT-specified tuning, 127,[129][130][131]159 thus allowing for some motivated external parameters. However, the here-defined basic AHCX design is deliberately kept free of adjustable parameters.
In practical terms, the SCAN, CX or AHCX choice comes down to a discussion of the nature of the material system as well as to attention to accuracy needs. Below we simply exemplify the AHCX potential in a set of demonstrator challenges (including noble-metal adsorption) and we include comparisons with CX and with literature SCAN results 173,174 to illustrate differences in performance.

IV. COMPUTATIONAL DETAILS AND ANALYTICAL-HOLE MODEL VALIDATION
We have coded the AH exchange-hole model and AHCX (as well as HJS-PBE, PBE-AH, and PBEsol-AH) in an in-house version of Quantum Espresso. 3,137,138 It is a clear advantage of the HJS framework that it allows an analytical evaluation of the formal 'SR' enhancementfactor expression Eq. (20), for example, in terms of finding and coding functional derivative terms. The analytical evaluation is formally given in Eq. (43) of Ref. 119. For example, for our AHCX implementation we need simply to evaluate the terms for expressions and parameters specific to CX.
The implementation is fully parallel and can directly benefit from the computational acceleration that the ACE operator 3,135 provides for the Fock-exchange component F SR FX . This makes it possible to run efficient calculations of molecules, of bulk metals (requiring many k-points), and of surface slab systems on a standard highperformance-computer cluster.
For a demonstration of performance on bulk structure and cohesion, on broad molecular properties in the GMTKN55 suite, 5 and on CO adsorption, we use the electron-rich optimized normconserving Vanderbilt 1 (ONCV) PPs, in the ONCV-SG15 release 2 , at a 160 Ry wavefunction energy cut off.
Beyond the core documentation (bulk, molecule, and adsorption), we also include an AHCX demonstrator on workfunction and surface energy performance. These results involve computations of many slab geometries, all with surface relaxations, 60 and at times with cumbersome electronic convergence. Hybrid studies of noblemetal surface properties are expensive in the electron-rich ONCV-SG15 PPs. Meanwhile, use of norm-conserving PPs significantly helps stability of the ACE Fockexchange evaluation that enters all types of hybrid calculations, at least in the Quantum Espresso version where we placed our AHCX implementation. For the clean-noble-surface AHCX demonstrator work we therefore use the set of more electron-sparse AbInit normconserving PPs. 176 We argue for the validity of this approximation for noble metals and we track the likely impact by also providing workfunction and surface energy characterizations for PBE and CX in both the AbInit PPs and in the ONCV-SG15 PPs.
A. Analytical-hole model validation Table II reports a simple sanity check that we provide to argue robustness of our maximum-entropy description of the AH model and RSH constructions. We use the fact that Ref. 119 reports parameters for their AH characterization of PBE exchange, a form that is here termed HJS-PBE. We provide a modeling self-test using the ONCV-SG15 PPs at 160 Ry, noting that the details of the H(s) specification must not change the RSH results. Table II contrasts the results of such HSE-like descriptions for 6 atomization energies. For the HJS-PBE and PBE-AH we chose values of the Fock-exchange mixing α = 0.25 and for the screening γ = 0.106 (inverse Bohr's) that exist for default HSE runs. 117 The agreement with reference data is good. The alignment of HSE and HJS-PBE (or AH-PBE) descriptions is very good.
Most importantly, Table II illustrates that our code to extract AH descriptions for PBE fully aligns results of PBE-AH with those of HJS-PBE. We trust our AH model description of PBE exchange, as defined by parameters in Table. I. Interestingly, the details of our H specifications (Table  I) differ from those provided by HJS. 119 This is true even if there is no difference in the resulting RSH descriptions. This finding suggests robustness in the overall AH constructions of RSH, not only for the HJS-PBE but also for the here-defined vdW-DF-based RSH, the AHCX.

B. DFT calculations: molecules
For a survey of performance on molecular properties we have turned to the large parts of GMTKN55 that are immediately available for a planewave assessment using the ONCV-SG15 PP set. 1,2 The full GMTKN55 is organized into 6 groups of benchmarks: 1) small molecule properties, 2) large molecular properties, 3) barrier properties, 4) inter-molecular noncovalent interactions (NCIs), 5) intra-molecular NCIs, and 6) total NCIs. However, we exclude the G21EA and WATER27 benchmark sets of GMTKN55. 5 In effect, we focus on a subset, the remaining easily-accessible 'GMTKN53', in our discussion.
The G21EA benchmark is excluded because it contains negative charging of ions and small radicals (like OH − ). Such small negatively charged systems genuinely challenge our planewave assessment of broad molecular properties exactly because we seek to maintain a high precision. The fundamental problems arise by a combination of two factors. First, the self-interaction errors (SIEs) will, in these systems, push the highestoccupied molecular-orbital level up towards or above the vacuum floor, when it is done in a complete-basis set description. 121 Second, our planewave description rapidly approached that complete-basis set limit due to our efforts to also carefully control spurious inter cell electrostatics interactions down towards the 0.01 kcal/mol limit. 178 We conclude that a direct performance assessment on G21EA is meaningless in a planewave code and it is perhaps best left for a more advanced handling. 121 We note in passing that we can approximately assert the G21EA performance by using the planewave equivalent of the so-called-moderate-basis-set approach, 121 trapping the negatively charge ions in a small 6-to-8-Åcubed boxes. We can thus obtain an approximate comparison that shows that the variation in G21EA performances has no discernible impact on the statistical performance measures that we track in GMTKN53.
We furthermore removed WATER27 for a related SIEinduced 121 problem -the path to electronic-structure convergence for the OH − is exceedingly cumbersome once we seek to fully converge the WATER27 set with respect to unit-cell size. Of course, hybrids helps in the WATER27 (and in the G21EA) study -but then we can offer no fully converged comparisons. Again, using as separate (12Å) small-box handling of the OH − system makes it possible to complete an approximate WATER27 benchmarking; As will be detailed elsewhere, this assessment identifies vdW-DF2, CX, CX0p and AHCX as top performers.
Fortunately, having 53 (out of 55 sets of Ref. 5) still makes for ample statistics. Accordingly, we simply work with the easily-accessible GMTKN53 and compare with literature GMTKN55 results. 5 This is fair because, if anything, we have negatively offset our AHCX performance reporting (by omitting G21EA and WATER27).
We supplement the GMTKN53 characterization both by tracking performance differences (that arise for some individual benchmarks) and by analyzing the groups of benchmarks that mainly reflect NCIs. We essentially use the same computation setup as the pilot study included in a recent focused review of consistent vdW-DFs, Ref. 60. However, we have now moved off of the electron-sparse AbInit PPs and onto the electron-rich ONCV-SG15; 1,2 Moreover, we now make systematic use of large cubic cells and a Makov-Payne electrostatics decoupling; 179 In the pilot assessment, 60 the electrostatic decoupling was only done for charged systems.
We find that this upgrade of our benchmarking strategy accounts for some of the deviations that we previously incorrectly ascribed to limitations of the CX 56 and vdW-DF2-b86r 180 functionals. 60 The past limitations seem to, relatively speaking, especially impact our assessment of Group 4 (intermolecular NCI) performance, as represented in Refs. 4 and 60. This paper provides a more fair performance characterization of the performance of CX and vdW-DF2-b86r as well as of AHCX relative to other vdW-inclusive DFT versions, for example, dispersion-corrected metaGGA.

C. DFT calculations: bulk and surfaces
For bulk structure and cohesive energies, we used the Monkhorst-Pack scheme with an 8 × 8 × 8 k-point grid for the Brillouin-zone integration. Our demonstration survey includes 4 semiconductors, 3 ionic insulators, 1 simple metal, and 5 transition metals including Cu, Ag, and Au. Here we also used the electron-rich ONCV-SG15 1,2 PPs set at 160 Ry wavefunction energy cut-off. We tracked the total-energy variation with lattice constant and used a fourth-order polynomial fit 7 to identify optimal (DFT) structure, the cohesive energy, and bulk modulus.
Surface energies for Cu (111), Ag (111), and Au(111) were estimated from the energies calculated from 4 to 12 layers of the noble metals in unit-cell configurations with 15Å vacuum. Here we sought to offset some of the high computational cost by instead using the AbInit normconserving PPs 176 for extra demonstrations of AHCX usefulness. In the surface studies, we first determine the in-surface unit cell from the calculated bulk lattice constant. Next we allow the outermost atoms to relax in a slab-geometry description. The exception is for AHCX characterizations, where we instead relied on the results Scaling of computational cost for hybrid calculations of gold and diamond using the primitive cell. We compare central-processing-unit(cpu)-hour cost per Fockexchange evaluation step for calculations of the bulk cohesive energies as computed using a 8 × 8 × 8 k point mesh and 160 Ry wavefunction energy cut off for the ONCV-SG15 PPs. We track the dependence on using a so-called 'down-sampled' q mesh (setting k-point differences) for evaluating the Fockexchange contribution in the hybrid study. The inset tracks how the cohesive energies change with the choice of the q mesh.
of the CX structure characterizations.
To document AHCX accuracy for CO adsorption on noble-metal surfaces, we used again the electron-rich ONCV-SG15 PPs at 160 Ry. We consider the FCC and TOP sites in a 2x2 surface-unit-cell description and a six-layer slab geometry. We also track the adsorptioninduced relaxations on the top three layers and compute the CO adsorption energy as follows, Here E CO/metal (111) denotes the fully relaxed energy of the adsorbate system, while E metal (111) and E CO denote the total energy of the fully-relaxed (isolated) surface and molecule system, respectively. Figure 3 and Table S.II of the SI material detail the computation time required for completing one evaluation of the ACE Fock-exchange operator in hybrid calculations on a single 20-core node for gold and diamond, top and bottom panels, respectively. We include a comparison of HSE, 117 CX0p, 123 and AHCX for gold and diamond, using a 8x8x8 k point grid and the electron-rich ONCV-SG15 PPs. We track the variation with the socalled q mesh (of k point differences used in the Fock exchange evaluation).

D. Computational costs
There is no additional cost of doing AHCX over HSE (as in the molecular cases discussed in the introduction) and the AHCX cost is slightly lower than the cost of the corresponding unscreened hybrid CX0p. This is in part a statement of the computational efficiency of the new nonlocal-correlation RSH functional. It is also a statement that we have provided a robust implementation of the analytical hole description and that we have succeeded in leveraging the benefits of the Quantum Espresso ACE implementation. 3 An overview of scaling of AHCX (or HSE) hybrid calculations can be deduced from Ref. 3. The cost per Fock exchange evaluation scales with the FFT scope (size of unit cell) and with the number of bands (that is, number of electrons per unit cell) squared. 3 Also, for hybrid studies of bulk and surfaces one needs to consider multiple k points and a q mesh of k-point differences, causing an additional scaling factor. When there is symmetry (as in simple bulk) one can cut down this factor significantly by limiting the set of k and q points that enters in an actual Fock-exchange evaluation. Importantly, however, there are communication costs when jobs are spread across multiple nodes. 3 For example, we can push more hybrid calculations with electron-sparse PPs because they require fewer bands; By limiting bands we both directly accelerate the computation speed and we reduce the memory requirements, making our jobs fit on more types of computer resources.
The AHCX is more costly than non-hybrids, including metaGGAs, for example, as illustrated in the introduction. Also, as mentioned above, hybrid calculations are costly for surfaces where relaxations enviably produce low-symmetry geometries. However, the AHCX computational costs can be handled by choosing an appropriate parallelization strategy. We try to limit the number of nodes as we also try to accommodate the memory requirements.
For a specific illustration of AHCX costs we make the following observations from our studies. We report summary of nearly 5000 AHCX and HSE-type hybrid studies of molecular properties of the GMTKN55 (having up to 37-Å-cubed unit cells to carefully control electrostatic couplings); Those ONCV-SG15 jobs took us a total of 200000 core hours (5-6 times the cost of using the corresponding CX and PBE regular functionals). The cost for an individual bulk structure study in AHCX/HSE is (due to symmetry) relatively cheap (see Fig. 3); Here the 700+ computations of energy-versus structure variations that enters our structure optimization (of 11 simple cases and of 5 transition metals) took just over 100000 hours at the production stage. Finally, the cost for hybrid studies of noble-metal surface energies and adsorption are on a different scale. We were able to obtain the 6 frozengeometry AHCX CO-noble-metal adsorption energies in the electron-rich ONCV-SG15 PP setup at the costs of about 100000 core-hours each (on selected computer resources). For the additional noble-metal workfunction and surface energy demonstrators (involving many slabs) we cut the cost to about 1.5 million core hours by instead relying on the electron-sparse AbInit PPs for AHCX calculations.

V. RESULTS AND DISCUSSION
A. Broad molecular properties Figure 4 shows a performance comparison that we provide for 10 molecular benchmarks of the full GMTKN55 suite. The figure characterizes the performance in terms of mean absolute deviation (MAD) values. Use of hybrid AHCX provides significant improvement over CX in the case of self-interaction problems (in the SIE4x4 benchmark set), for atomization energies (in the W4-11) set, and for several barrier-type problems, including those reflected in the PX13 and WCPT18 benchmark sets. The improvements are perhaps largest in the case of C60 isomerizations, for large deformations having stretched bonds. Figure 5 summarizes our survey of performance across the entire NCI groups of the GMTKN55, although excluding the WATER27 benchmark set for reasons explained above. The figure characterizes the performance in terms of a total MAD, or TMAD, value for intermolecular NCIs (abscissa position) and of a TMAD value for intramolecular NCIs (ordinate position). These effective MAD values result as we first compute the MAD values for each NCI benchmark of the GMTKN55 and then average over these MADs over the number of intermolecular and intramolecular NCI benchmarks investigated in GMTKN55 Group 4 and Group 5, as in Ref. 4. The SI material provides details and lists the quantitative data that we provide for vdW-DF hybrids (AHCX and CX0p, circles) and for regular vdW-DFs (vdW-DF2, vdW-DF2-b86r and CX, upwards triangles), as well as for rVV10 and dispersion corrected revPBE-D3 and HSE-D3. The SI material also includes a performance characterization of vdW-DF1 51,52,54 and vdW-DF-ob86. 182 The NCI assessments in Fig. 5 can be compared to Fig. 3 of Ref. 4, a study that also summarized NCI results obtained in Refs. 5 and 183. The assessments differ in code nature (planewave versus orbital-based DFT) but we can use Fig. 5 for a broader discussion of performance. That is, our comparison for NCI performance also indicates how the new AHCX RSH vdW-DF hybrid fares in relation to the dispersion-corrected metaGGAs and dispersion-corrected hybrids discussed in Ref. 4. We find that the AHCX performance is better than those of the dispersion-corrected metaGGAs and of ωB97X-D3 (as reported in Ref. 5 and 183). However, our AHCX does not fully match the performance of DSD-BLYP-D3, 5 when it comes to the TMAD value for intermolecular NCI.
For completeness, Fig. 5 also shows the approximate assessments (downwards triangles) of CX and vdW-DF2-b86r performance that we obtained in a pilot benchmarking. 60 The differences in CX and vdW-DF2-b86r positions of upwards and downward filled triangles reflect the improvement that we have here made in benchmarking strategy relative to Ref. 60. The systematic use of decoupling of spurious dipolar inter-cell attractions is found to have a clear effect on performance for the group of intermolecular NCI benchmark sets.
We find that the use of vdW-DF hybrids leads to significant improvements in the description of problems characterized by strong NCIs. The performance of AHCX and CX0p are again nearly identical but significantly better than dispersion-corrected HSE-D3 and of the regular vdW-DFs. This is not surprising because the hy- Progress in performance from regular vdW-DFs (upward triangles) to hybrid vdW-DFs (circles) on benchmark sets reflecting noncovalent interactions (NCI), but excluding the WATER27 set. The SI material documents that vdW-DF1 can be found just outside the range. Our characterization can be directly compared to the broader survey reported in Fig. 3 of Ref. 4. As in that study, the functional performance is represented in terms of 'total mean-absolutedeviation' (TMAD) values (that average MAD values over benchmark sets, see text) obtained for intermolecular NCI sets (x-axis) and for intramolecular NCI sets (y-axis). The present characterizations for CX and vdW-DF2-b86r are more accurate than the assessment (shown by downward triangles) previously obtained in a pilot characterization. 60 For reference we also include characterizations of dispersion corrected revPBE-D3 and HSE-D3, as well as of rVV10.
brids are expected to have a more correct description of orbitals and reshaping orbitals (and hence the density variation) will directly affect the strength of the vdW attraction. 34,[184][185][186][187][188][189] For a summary of overall progress that the hybrid vdW-DFs may bring for general molecular properties, we rely on the 'easily available GMTKN53', defined above. This is the subset of the GMTKN55 suite that omits the WATER27 and G21EA benchmarks sets, for good reasons. 121 Fortunately, with GMTKN53 we remove just one benchmark set in two groups and can still provide a balanced account. It is still meaningful to discuss functional performance in terms of per-group results as in the GMTKN55 study. 5 For comparison of performance we therefore use the reference geometry and reference data of Ref. 5 to determine the deviation on each of the over 2000 molecular processes that are included. We evaluate the weighted total mean absolute deviation measures (WTMAD1 and WTMAD2) that Grimme and co-workers introduced to allow a comparison also among the GMTKA55 groups. 5 The adaptation of these measures to the present GMTKN53 survey simply involves adjusting the counting , and AHCX (a CX-based RSH). The comparison is given for a 53-benchmark subset of the full GMTKN55 suite, 5 namely those that are directly available to planewave benchmarking using the electron-rich ONCV-SG15 PPs. As indicated by the '*' superscripts, we have here omitted the G21EA benchmark set from the GMTKN55 Group 1 (small molecular properties) and the WATER27 benchmark set from Group 4 (intermolecular-noncovalent interactions) because they contain small negatively charged ions and radicals (see text); We include all benchmark sets from Group 2 (large molecular properties), Group 3 (barrier properties), and Group 5 (intramolecular properties), and thus essentially all sets in Group 6 (total noncovalent interactions). The performance is measured by the weighted total mean absolute deviation WTMAD1 measure introduced in Ref. 5. of participating benchmark sets for GMTKN55 groups 1 and 4, 5 see SI material for further details. Figure 6 shows our summary functional comparison based on this 'easily accessible GMTKN53', using WT-MAD1 measures. The SI material lists the data for both WTMAD1 and WTMAD2 comparisons. We find that the AHCX performance is significantly better than dispersion-corrected HSE-D3 for this assessment of broad molecular properties. The AHCX also performs systematically better than CX, which in turn performs at the level of dispersion-corrected revPBE-D3. The latter is reported to be one of the best performing dispersion corrected GGAs. 5 Importantly, the AHCX also performs at the same level and perhaps slightly better than the unscreened CX0p hybrid. Our primary motivation for defining AHCX is to  Table  III. In the violin plots, the white diamonds mark the mean deviation, the box identifies the range between the first and third quartile, and the black line indicates the median of the distribution. Moreover, whiskers indicate the range of data falling within 1.5*box-lengths of the box.
have a vdW-DF RSH that can also describe adsorption at metallic surfaces, which is not something that the unscreened CX0p is set up to do. The comparison in Fig. 6 shows that the AHCX is also good enough on the molecular side of interfaces. That is, it remains a candidate for also serving us for the interface problems, at least so far.
B. Bulk structure and cohesion Figure 7 contrasts the performance of hybrid functionals HSE and AHCX with those of the underlying PBE and CX functionals for bulk: 6 metals (elements Al, Cu, Rh, Ag, Pt, Au), 4 semiconductors (Si, C, SiC, and GaAs) and 3 ionic insulators (MgO, LiF, and NaCl).
Here we compare deviations from measurements of lattice constants a and cohesive energies ∆E (both back corrected for zero-point energy and thermal effects) as well as for bulk moduli B 0 . Symbols for the statistical measures are summarized in the figure caption. The PBE (AHCX) description for cohesive energies (for bulk moduli) have Ag and Au (Rh) as statistical outliers; The SI material contains a detailed presentation of the per-bulkmaterial performance.
The figure shows that AHCX gives the smallest median deviation and spread of data for the deviations for lattice parameters and bulk moduli; There are clear improvements with AHCX over both HSE and CX (which is the second-best performer overall). For cohesive energies, the CX has the smallest median deviation whereas AHCX has the smallest spread and both perform better than PBE and HSE for bulk properties.  Table III compares the results of PBE, CX, and AHCX descriptions of transition-metal structure and cohesion to back-corrected experiments values 6 (last column). The results are provided for the electron-rich ONCV-SG15 PPs (except where noted by an asterisk '*' or when taken from literature).
We first observe that use of ONCV-SG15 PP setup yields a state-of-the-art benchmarking on these transition metal systems. For PBE (and for CX) we can track deviations of our self-consistent ONCV results relative to the self-consistent PBE (non-selfconsistent CX) allelectron results, provided in Ref. 6, and to the fully self-consistent PAW-based assessments that one of us have previously provided. 166 For lattice constants the largest PBE (CX) deviation from self-consistent (nonselfconsistent) all-electron results 6 is 0.2% (0.1%). The ONCV characterization has but minute differences from (is spot on) from the PAW-based characterization of PBE (CX) performance. 166 Similarly, for the description of co-hesive energies, we find small PBE (CX) deviations, with a Rh maximum of 3.0% (1.4%), relative to all-electron results. 6 Here the ONCV characterizations are slightly less precise overall than the previous PAW-based characterizations (but have a smaller spread). Finally, for the bulk moduli, we find that our ONCV-SG15 benchmarking differs by at most 2.2% for Rh from the allelectron transition-metal results. Overall, we find that our ONCV-SG15 bulk-structure benchmarking is highly precise, for example, for characterizations of noble-metal properties.
Comparing next the calculated results against backcorrected experimental results, 6 Table III confirms the overall impression of AHCX promise (Fig. 7) also holds for the transition metals. The table shows that CX performs significantly better than PBE on the transition metal lattice constants (cohesive energies); This is also expected, perhaps especially for the noble metals. 6,60,90,166 A move from PBE to CX reduces the maximum absolute deviation (from experimental values) on transition-metal lattice constants (cohesive energies) from 2.4% (24.2%) to 0.8% (10.1%). Meanwhile, AHCX performance is on par with that of CX for lattice constants and a clear improvement over HSE (as detailed in the SI material). This pattern is repeated for the bulkmodulus assessments, see SI material.
Overall for AHCX we find the following deviations from experimental lattice-constant (cohesive-energy) values: 6 -0.3% (-4.7%) for Cu, 0.2% (-6.4%) for Ag, 0.8% (-10.3%) for Au, -0.2% (-5.8%) for Pt, and -0.7% (-9.3%) for Rh. The performance for cohesive energies is overall slightly worse than that for CX on transition metals, but significantly better than that of PBE and HSE, see SI material. While broader tests of bulk-transition-metal performance of AHCX, for example, extending Refs. 6 and 166 are desirable, it is also clear that AHCX is useful for descriptions of bulk properties in general and certainly for noble-metal studies.

C. Noble metal surfaces
Descriptions of (organic) molecule adsorption at metal surfaces are an important reason for defining and launching the AHCX. For metal problems, we must screen the long-range Fock-exchange component. The new RSH vdW-DF is set up to reflect the perfect electrostatic screening. 130 Below we begin the discussion by looking at the AHCX ability, as a RSH vdW-DF, to describe properties of noble-metal surfaces. We stay with the robust (but expensive) electron-rich ONCV-SG15 PPs (at 160 Ry) for the regular-functional descriptions and for the description of CO adsorption, next subsection; For a demonstration of AHCX's ability to describe workfunctions and surface energies, here, we shall use and discuss calculations obtained in the more electron-sparse (but still normconserving) Ab Init PPs (at 80 Ry).
We note that understanding surface energies of metal- lic nanoparticles is important for correctly setting up cost-effective catalysis usage in catalysis. This is because they define the Wulff shapes 194 of nanoparticles and hence the extent that they provide access to specific active sites. 162,195 As such this additional illustration could also be relevant for DFT practitioners.
Table IV compares our PBE, CX, and AHCX results for workfunctions and surface energies against experimental values and against literature SCAN values. 174 An asterisk '*' on a column indicates that those results are obtained using the AbInit PPs at 80 Ry, the rest are obtained by ONCV-SG15 PPs at 160 Ry. All results are based on a sequence of slab calculations and we first compute bulk lattice constants and then surface relaxations specific to the functional and PP choice. As indicated by the prime on the results in the AHCX column, however, we rely on the CX structure determination for our characterization of this CX-based RSH.
The table furthermore compares against experimental values for workfunctions (as summarized in Ref. 174), and surface energies for noble-metal (111) facets. 60,174 Experimental surface-energy values are available from observations of the liquid metal surface tension. 196,197 As such, they are defined by an average and it is natural to focus on the major facets 174 where we now introduce σ (other) = (σ (110) + σ (100) )/2. Previous works have compared measured surface energies σ exp (for various transition metals) and compared with per-facet surface-energy results -σ (111) , σ (110) , and σ (100) -computed in various functionals, including CX. 60,162,174,195 Here, we simplify the assessment task in two steps: a) we use our past CX results 60 to define a ratio x = σ (111) /σ (other) that reflects the relative weight of (111) in the major-facet averaging, 60,174 Eq. (25), and b) we extract the experimentally-guided estimate Table IV shows that moving between PBE, SCAN (as asserted in literature 174 ), CX, and AHCX generally increases the computed values for both workfunctions and surface energies for the noble-metal (111) facets; We note that the same trend also arise (among the PBE, CX, and AHCX descriptions) when we instead compute these with the Ab Init PPs at 80 Ry. It is alone for the Cu surface energy that we find an AHCX surface energy that is smaller than the CX results (AHCX gives the same Ag surface energy as CX).
Table IV also shows that switching between the electron sparse AbInit PPs and the electron-rich ONCV-SG15 PPs has essentially no impact on the workfunction and a 5% (10%) effect on the CX results for the surface energy of Cu and Au (of Ag). These shifts are directly correlated with the fact that the Cu and Au (Ag) lattice constant in the AbInit PP description is 1% (2%) larger than what results with the ONCV-SG15 PP choice. Table IV suggests that both CX and AHCX provide an accurate description of the noble-metal surface energies and workfunctions. We note that a significant uncertainty exists for surface energy reference values due to the need to interpret the input from high-temperature measurements. Unlike for PBE and SCAN, the CX and AHCX surface energies fall within the range of such observation-based estimates (For Au (111), the CX, but not the AHCX, result falls just outside the error bars). For workfunctions we find that CX and SCAN are slightly more accurate than AHCX for Cu (111), while AHCX is the best performer for the Ag (111) surface; For the To argue AHCX accuracy for clean-noble-metal surface properties, we make four observations. First, our use of AbInit PPs for AHCX likely impacts the comparison with the ONCV-SG15 based PBE and CX results for the surface energy σ. Second the trend from the PBE and CX results suggests that the impact on the PP choice on σ primarily arises by the inverse dependence on the surface-unit-cell area, an area-scaling that is discussed two paragraphs above. Third, Appendix A assesses the sensitivity of the CO-adsorption site-preference results to the PP choice, documenting good robustness for the noble metals; This robustness is expected by the fact that these systems have a nearly closed d band 118,198,199 and it makes it likely that the PP impact on the AHCX σ description will be dominated by the area-scaling effect. Fourth, adjusting our AbInit-PP based AHCX characterization upwards by a 5-10% area effect will generally push the AHCX description closer to the back-corrected experimental surface-energy values: The AHCX result for Ag will still reside on the edge of the range characterizing the Ag experimental data. In other words, Table  IV shows that AHCX is useful and that it may well be a strong performer for characterizations of clean-noblemetal surface properties.

D. Noble metal adsorption
Finally, we turn to discuss a case of actual metalsurface adsorption. We focus on CO adsorption as it allows us to use just a 2x2 in-surface unit-cell extension and because it may provide a hint of whether AHCX ultimately may help in understanding heterogeneous catalysis. The CO adsorption on late transition metals is a classic surface-science challenge. Figure 8 shows the CO adsorption geometry for 111 surfaces. The CO is standing upright with a Blyholdertype binding 200 through the carbon atoms (brown ball) which sits between the metal substrate (gold atoms) and the oxygen (red ball). There is no experimental ambiguity that CO molecule adsorbs on the TOP site (as illustrated in the right panel), rather than at the FCC or hollow site (left panel) at low coverage, Ref. 201. Table V summarizes the CO adsorption studies that we here provide using electron-rich ONCV-SG15 PPs. The most famous CO adsorption problem arises at the Pt(111) surface, Ref. 205-208, a problem that has repeatedly challenged DFT descriptions when used for XC energy approximations that also deliver an accurate description of the in-surface Pt lattice constant and surface energies, Refs. 173, 209-212. Here we focus on the noble-metal adsorption cases. The table compares our results for CO adsorption at FCC and TOP sites on Cu (111), Ag (111), and Au (111). As indicated by the prime, the AHCX characterization is again done frozen at the adsorption structure that results in the CX study. We find that the AHCX description stands out by predicting the correct site preference for CO on noble-metal (111) surfaces; None of the PBE, SCAN or CX functionals have that consistency. The AHCX prediction of the actual (TOP-site) adsorption energy is good in silver and gold; All functionals have some problems in accurately predicting the magnitude of the adsorption energy for CO/Cu (111).
The AHCX functional overestimates the binding in the case of Cu (111). Comparing with the CX description, we see that the AHCX still provides a reduction of the adsorption energy. This AHCX binding softening (compared to CX) is, in part, expected. For instance, use of either PBE0 and HSE reduce the Cu(111) adsorption energy compared with PBE, Ref. 209.
We also highlight an important point of the AHCX success: We find that AHCX simultaneously have accurate descriptions of molecules, bulk lattice constants, surface properties, and CO adsorption energies; There is no lucky hit for one quantity at the expense of others, compare Fig. 6 as well as Tables III, IV and V. The AHCX lattice constant for Cu has a 0.3% deviation from experiment, better than CX (see SI material). The AHCX descriptions of workfunctions are accurate for Cu. Ag, and Au, and the description of Cu, Ag, and Au surface energies are good. Last but not least, the AHCX performance for predicting the CO-adsorption site preference is excellent, and (except for Cu) accurate also on the predictions of the actual CO adsorption energies.

VI. SUMMARY AND CONCLUSION
We find that the new AHCX RSH performs clearly better than both CX and dispersion-corrected HSE for molecules and better than CX and HSE for extended matter. For molecules, our new AHCX RSH performs at the level of the corresponding unscreened CX0p hybrid.
Furthermore we find that the AHCX is accurate for the description of noble-metal surface properties, at least for the (111) facet. The AHCX improves both the PBE and CX descriptions in terms of being overall more accurate on surface workfunctions and surface energies for the (111) noble-metal facets. It furthermore stands out by predicting the right site preference for CO adsorption on the (111) facet across the noble metals. For Ag and Au the AHCX adsorption energies are accurate, but the binding strength is overestimated for Cu (111). However, at a geometry fixed by the CX structure characterization, the AHCX still improves the CX description of the CO/Cu(111) adsorption binding energy.
Overall the AHCX shows promise for tackling a longstanding problem of describing organics-metal interfaces. In the heterogeneous systems we often want a vdWinclusive hybrid for the molecule side, yet we cannot motivate the use of unscreened hybrid CX0p in cases with a metallic nature of conduction for the substrate. The new screened CX-based RSH, termed AHCX, shows overall, a significantly better performance on molecules than dispersion-corrected HSE. Our results also indicate a strong performance on semiconducting and metallic systems. It is a candidate for making DFT better at characterizing chemistry in heterogeneous systems.

SUPPLEMENTARY MATERIAL
See the supplementary material for the parameters of HJS analytical-hole, computational time scaling, molecular and bulk benchmark data. Hybrid calculations are most stable for normconserving PPs in the version of Quantum Espresso that we used for our AHCX implementation and subsequently for our documentaion work; A workaround was introduced recently in the Quantum Espresso code package, but it does not seem to always help and, in any case, it did not arive in time for us to benefit in this project. The stability concerns exist whether we use the existing HSE 116,117 or our coding of analytical-hole formulations. 119 Meanwhile, it is convenient to have a (nearly) complete suite for materials explorations. The AbInit PP 176 and the ONCV-SG15 1,2 PP releases are two such options.
We primarily use the electron-rich ONCV-SG15 because it included semi-core electrons as valence states to allow simpler discussion of AHCX performance examples. Use of the ONCV-SG15 is a costly option but also considered safe, for example, for hybrid studies.
The AbInit PPs 176 are also normconserving but electron sparse, relying instead on a nonlinear core correction to represent the remainder. One must expect a lower accuracy than what is available with ONCV-SG15 1,2 or from PAW-based structure descriptions, e.g., Ref. 166; Comparing the AHCX columns with and without an asterisk '*' in Table III gives an illustration. Also, relying on a nonlinear core-corrections, the use of the AbInit PPs can cause problems for hybrid descriptions of adsorption on open-d-shell systems. 118 However, since they need fewer band (lower the computation costs) and have memory requirements that fits more computer resources, the electron-sparse approach may sometimes be an interesting option. We used the AbInit PPs to provide an additional demonstration of the potential AHCX usefulness regarding noble-metal workfunctions and surface energies. We did that noting that we work with nearlyclosed-shell systems (where the importance of d-band rehybridization is generally reduced 198,199,209 ). Table VI reports an assessment of the impact of PP choices on CX and AHCX adsorption energy results for CO on noble-metal surfaces. Here and in a molecular survey (not shown), we find that using instead the electron-sparse AbInit PPs introduces no systematic changes in the conclusion that AHCX is accurate for molecules and for predicting the site preference for noble-metal adsorption. The ONCV-SG15 descriptions are more accurate in terms of giving better lattice constants and better absolute adsorption values across the noble metals. However, we also find working descriptions with the computationally cheaper AbInit PPs, for noble metals. This AHCX robustness on adsorption suggests that AHCX is, in fact, promising also on workfunction and surface energy results (the additional demonstrator included in Section V.C). A future ONCV-SG15 characterization of noblemetal surface energies and workfunctions may quantify the impact of this approximation, but some problems are simple enough that we can leverage the acceleration and memory savings, when desired.
Supplementary Materials for: vdW-DF-ahcx: a range-separated van der Waals density functional hybrid    Table S II shows the time for converging one complete Fock-exchange evaluation step in the Quantum ESPRESSO code with hybrid functionals HSE, vdW-DF-CX0p (CX0p,) and vdW-DF-AHCX (AHCX). The calculations of cohesive energies E coh were done at the (back-corrected) experimental lattice parameters of gold and diamond to allow the reader to also track the change in evaluation precision with the q mesh of k-point differences in the Fock exchange evaluation. We used an 8 × 8× k point sampling, the ONCV-sg15 1,2 set of normconserving pseudopotentials and a 160 Ry wavefunction energy cut off. We also used 2 extra bands for the semiconductor and 10 extra bands for the metal cases and note that (extra) bands cost in hybrid calculations. 3 Besides the resulting E coh description, we list the average central-processor-unit(cpu) core-hour cost per Fock-exchange evaluation step for various choices of q meshes. The timing results reflect calculations on 1 node (20 cores sharing 60 GB of memory) of a high-performing computer cluster. The data is plotted in Fig. 3 of the main text. The listed cohesive energies should not be considered as more than an indication of functional performance as they were not done at the lattice constant native to our functional. Instead, we refer to Table III, section V.B, and the appendix in the main text, as well as to SI Tables S VI and S V (below,) for a more complete assessment of AHCX performance on bulk properties.
Appendix C: Molecular benchmarks TABLE S III. Performance of vdW-DFs and of a few vdW-inclusive DFTs as asserted on the intermolecular and intramolecular NCI benchmark groups of the GMTKN55 suite. 5 As indicated by the asteriks '*', we exclude the WATER27 benchmark set in the intermolecular NCI group. We list so-called TMAD values in kcal/mol. The TMAD values are defined in the text and reflect mean absolute deviation (MAD) values averaged within the two NCI benchmark groups.  Table S III lists the data plotted in Fig. 5 in the main text. We track the performance across essentially the full set of noncovalent interaction (NCI) benchmarks of the GMTKN55 suite. We exclude the WATER27 set for reasons explained in the main text; This also implies that our vdW-DF assessments, plotted in Specifically, the data in NCI sets is resolved into and analyzed in terms of a total mean-average-deviation (MAD) values (denoted TMAD) for GMTKN55 groups 4 (intermolecular NCI) and for group 5 (intramolecular NCI). These measures arise by computing the MAD i (and MAD j ) values for every NCI benchmark set i in group 4 (and j in group 5), and then averaging over the number N = 11 (and M = 9) of sets in the group, for example, Table S IV summarizes the assessments we have performed for the 53-benchmark subset of the full GMTKN55 suite. 5 This data is organized into the benchmark groups of the GMTKN55 suite, omitting, however, the WATER27 set in group 4 and the G21EA set in Group 1, again see the main text. The data is plotted in Fig. 6.
We list the weighted total mean absolute deviation measures (WTMAD1 and WTMAD2), essentially as introduced by Grimme and co-workers in Ref. 5. However, we have adapted the measures to the slightly smaller range of benchmark sets, for example, Here, as in Ref. 5, we weight contributions of set i according to the reference value |∆E| i . The weighting w i is set to 10 (0.1) if |∆E| i is smaller than 7.5 kcal/mol (larger than 75 kcal/mol); The weighting w i is set to unity otherwise. For the alternative WTMAD2 scheme, one relies on a the overall suite average absolute deviation energy, |∆E| suite = 56.84 kcal/mol. It is computed from reference energies by considering all the reactions that enters in the GMTKN55 suite. Adapting to the reduced number (53) of sampled benchmark sets, we report and compare per-functional values  Fig. 6 of the paper. Our assessment is defined by the GMTKN55 suite 5 but we focus on the 53 sets that are easily accessible to our planewave benchmarking using the electron-rich ONCV-sg15 pseudopotentials at 160 Ry wavefunction energy cut off. Use of the WTMAD1 and WTMAD2 measures (adapted to averaging over 53 sets) permits a meaningful performance comparison both within and among the benchmark Groups 1-6 that are also defined in Ref. 5. As indicated by an astriks '*', we omit the G21EA sets in Group 1 and the WATER27 set in Group 4. For reference we also compare the performance of the CX against revPBE-D3, a strong performing dispersion-corrected GGA. 5 Best performance WTMAD1 numbers (as asserted in the 53 benchmark subset) for total-noncovalent interactions (Group 6) and overall are highlighted.

Functional
Measure Group 1 * Group 2 Group 3 Group 4 * where N i denotes the number of reactions in the individual benchmark set i. With a small further adjustment, we can also use Eq. (1) and (2) to provide a per-group assessment of functional performance, Table S IV and Fig. 6. For example, for WTMAD1, we simply a) limit the summations over i to the specific benchmark sets entering in each group and b) adjust the denominator to the number of benchmark sets that enters in the specific group.
While our WTMAD' forms, Eqs (1) and (2), are slightly adjusted, we still discuss them and compare them as estimates of full-suite WTMAD values. 5 This is done, for example, in Fig. 6 of the main text (while also using an asteriks '*' to emphasize the adjustment). Comparisons are motivated because we have merely removed one set out of 11 in Group 4, and 1 sets out of 18 in Group 1, omissions that can only make small changes in averaged assessments. In fact, we have validated that the weighted MAD values reported for AHCX in Table S IV undersells our new hybrid vdW-DF performance by 0.2 kcal/mol on the WTMAD1 measures and by 0.4 kcal/mol WTMAD2 measure when done at the level of electron-rich ONCV-sg15 PPs. 1,2 To that end, we obtained characterization in AHCX also for the WATER27 and G21EA sets that carefully controlled the impact of the the self-interaction errors, 8 as will be published elsewhere.