The shortcomings of semi-local and hybrid functionals: what we can learn from surface science studies

A study of the adsorption of CO on late 4d and $5d$ transition metal (111) surfaces (Ru, Rh, Pd, Ag, Os, Ir, and Pt) considering atop and hollow site adsorption is presented. The applied functionals include the gradient corrected PBE and BLYP functional, and the corresponding hybrid Hartree-Fock density functionals HSE and B3LYP. We find that PBE based hybrid functionals (specifically HSE) yield, with the exception of Pt, the correct site order on all considered metals, but they also considerably overestimate the adsorption energies compared to experiment. On the other hand, the semi-local BLYP functional and the corresponding hybrid functional B3LYP yield very satisfactory adsorption energies and the correct adsorption site for all surfaces. We are thus faced with a Procrustean problem: the B3LYP and BLYP functionals seem to be the overall best choice for describing adsorption on metal surfaces, but they simultaneously fail to account well for the properties of the metal, vastly overestimating the equilibrium volume and underestimating the atomization energies. Setting out from these observations, general conclusions are drawn on the relative merits and drawbacks of various semi-local and hybrid functionals. The discussion includes a revised version of the PBE functional specifically optimized for bulk properties and surface energies (PBEsol), a revised version of the PBE functional specifically optimized to predict accurate adsorption energies (rPBE), as well as the aforementioned BLYP functional. We conclude that no semi-local functional is capable to describe all aspects properly, and including non-local exchange also only improves some, but worsens other properties.


I. INTRODUCTION
The adsorption of carbon monoxide has intrigued researchers in surface science for the last decades. This has obvious reasons. Carbon monoxide is a simple diatomic molecule, but despite its simplicity it shows a rich phase diagram on metal surfaces. [1,2] The CO dissociation reaction involves a very simple reaction mechanism, but nevertheless the reaction barriers can be greatly varied by alloying or roughening the surface. [2] It is fair to say that, together with hydrogen adsorption and dissociation, CO has become a classical benchmark for experimental and theoretical surface science techniques. This interest is also strongly driven by the importance of CO for many technologically relevant reactions. [3] In this light, it is unsatisfactory that state of the art density functional calculations fail to describe several aspects of the adsorption of CO on metal surfaces accurately. This concerns both the predicted adsorption site as well as the absolute magnitude of the adsorption energy. For the Cu, Rh and Pt (111) surface, the generalized gradient approximation (GGA) in the PBE version predicts that CO adsorbs at a high coordination site (typically the hollow site), whereas experiments unequivocally show that atop adsorption is preferred. For Ag and Au a degenerate site preference is predicted, which does not agree with experiments either. [4,5] Equally compelling is that adsorption energies are significantly overestimated using the PBE functional. [4,5] These discrepancies have first been noticed in a now classical paper by Feibelmann et al. and the name CO adsorption puzzle has been coined. [4] Since then there has been emerging evidence that the present local and semi-local functionals are not capable to describe the subtle balance between donation of charge to the substrate and back donation to the molecule correctly for many adsorption problems. However, it is still an open question whether this is a universal shortcoming of such functionals, or whether the functionals can be amended to improve the description of surface related properties without drastically worsening other important properties. [6] The common model to describe CO adsorption is the Blyholder model, [7] which invokes interactions of the two CO frontier orbitals, the 5σ HOMO (highest occupied molecular orbital) and the 2π * LUMO (lowest unoccupied molecular orbital), with the metal states. Due to the interaction with the metal states bonding 5σ-metal orbitals and anti-bonding 5σ-metal orbitals develop, and the latter are partly shifted above the Fermi-level of the metal, causing a net bonding interaction (donation). Likewise, bonding 2π * -metal hybrid states become populated (backdonation). Simple symmetry reasons tell that the highly directional 5σ-metal interaction is particularly strong for atop adsorption, whereas the 2π * interaction is dominating for hollow site adsorption. [2,8,9,10,11,12,13] Here we apply a set of now well established semi-local as well as hybrid functionals to the CO adsorption problem. Our focus is on the systematic variation across the periodic table, and on how hybrid functionals compare to the more traditional semi-local functionals. We show that none of the available functionals is capable to yield an equally satisfactory description of the metal and the CO adsorption.

II. THEORETICAL AND COMPUTATIONAL METHOD
The first principles density functional theory calculations and the hybrid functional calculations have been performed using the VASP code, PAW [14] potentials in the implementation of Kresse and Joubert [15] and a cutoff energy of 400 eV. The surfaces are modelled by a periodic four layer metal slab with c(2 × 4) symmetry, and the CO molecule is adsorbed on one side of the slab. We focus only on the atop and hollow fcc adsorption sites, since the difference in calculated adsorption energies between the hcp and fcc hollow adsorption site is generally small, on the order of 0.05 eV; see, e.g., Refs. 16,17,18,19. The Brillouin zone integration is performed on symmetry reduced Γ-centered 6×6×1 grids (i.e. roughly 12 × 12 × 1 in the primitive cell). Since we want to concentrate on trends, all metals are considered in the fcc structure, although Os and Ru crystallize in the hcp lattice structure.
The PBE [20] and BLYP [21,22] functionals are used for the GGA calculations. For the hybrid functional calculations, the HSE functional [23,24,25] has been applied. We use a variant of the HSE06 functional that observes the homogeneous electron gas limit and all important sum rules. Contrary to the conventional HSE06 functional, the screening parameter is set to ω = 0.300Å −1 in both the semi-local as well as non-local part of the exchange functional, whereas the recommended choice is ω = 0.207Å −1 . As has been shown previously, the specific choice of the screening parameter ω has very little influence on the total energies, but slightly affects band gaps. [25,26,27] Throughout this paper we will use the acronym HSE for this functional. As we have shown elsewhere, for Cu, Rh and Pt the HSE functional yields the same site order and almost the same energetics as the more popular hybrid PBEh (sometimes also termed PBE0 or PBE1PBE) functional [31,32]  Furthermore the widely adopted B3LYP [21,33] functional is applied in the present work. The difference to the previous functional is predominantly in the correlation energy which is approximated using the semi-local Lee-Yang-Parr correlation. Systematic tests for selected extended systems are described elsewhere. [26,27,28,29,30] Implementational details can be found in Ref. 26. For the PBE and HSE functionals, the bulk metal lattice constants were optimized and used in the periodic slab calculations. In order to save compute time, the BLYP and B3LYP calculations were performed at the HSE geometries, which also agree best with experiment.  Table I shows that the overall description is best on the HSE level, good on the PBE level, at best, modest using the B3LYP functional and quite bad using BLYP. The HSE functional always yields smaller lattice constants than the PBE functional, and the contraction is largest for half-filled d bands and smallest when the d band is entirely filled (Ag). Compared to experiment, this yields to a sizeable underestimation of the lattice constants for Ru, a sizeable overestimation for Ag and Pd, and fairly accurate values for the remaining elements. The B3LYP functional gives larger lattice constants with substantial errors for Ag, Pd and Pt (3.4 %, 2.6 % and 2.5 %). Finally, the BLYP functional yields the largest lattice constants with errors up to 4 % for Pd, Ag and Pt. Paier et al. have already shown in Ref. 35 that this is mainly related to significant errors in the LYP correlation functional for metallic systems.
In summary, ad-mixing non-local exchange decreases the lattice constants, whereas the LYP correlation gives larger lattice constants than the PBE correlation.

B. Metal d band parameters
In table II we show the calculated metal d band parameters using PBE and HSE: the center of the d band, ǫ d and the metal work function, Φ. The experimental work function is also reported. The d band center can not be defined unambigously within a plane-wave approach. As we are mainly interested in trends (i.e. PBE vs. HSE), we choose to calculate the centre of gravity of the occupied d band (integration up to the Fermi level), which is a lower bound for the true d band center (see Ref. 5 for details). In table II, we have also introduced a fractional d band filling (n d ), calculated as the d charge contained in a sphere centered at the metal atom. The sphere radius has been fixed in order to contain exactly (within numerical error estimated to be 0.02 e) the number of valence s and d electrons in the atomic configuration.
[38] Interestingly, it turns out that the volume of the sphere is almost equal to the volume per atom in the fcc lattice (Wigner Seitz sphere), providing a physically sound basis for our definition of n d .
We considered the bulk phase [40,42] at equilibrium geometry using the PBE functional for calculating n d . We will show below that this definition is helpful for discussing trends between and among the 4d and 5d metals. From table II, we see that ǫ d calculated using HSE (ǫ HSE d ) is always found at larger binding energies than ǫ d calculated using PBE (ǫ PBE d ). This is related to the tendency of HSE to give a larger d bandwidth and a slightly weaker self-interaction within the d shell than PBE. [18] The common picture of reactivity of metal surfaces suggests that the  lower in energy ǫ d is, the smaller the adsorption energies should be. We will show below that this is in contradiction with our results for HSE. Although ǫ HSE d < ǫ PBE d , the adsorption energies generally increase for HSE (see figure 1). The HSE work function is generally smaller than the one calculated using PBE. This can be traced back to the induced surface dipole changes due to the inclusion of the Fock exchange. Indeed, hybrid functionals cause a complex redistribution of the charge density in the surface layer: an analysis of the difference between the PBE and HSE charge density (not shown here) shows that HSE functional depletes charge from the d z 2 orbital in favor of d xz and d yz orbitals, thus decreasing the spill-out of electrons from the surface layer into the vacuum. This charge redistribution should also affect the initial steric repulsion, thus explaining why, for the top site, the HSE adsorption energies are larger than the PBE ones (see figures 1, 2 and 3). The net effect of this complex charge redistribution is to reduce slightly the surface dipole and the work function. [43] The overall agreement with the experimental work function is good for both functionals.

C. The gradient corrected PBE functional
The central results of the present work are summarized in figures 1 and 2, where we show the adsorption energies for the atop and fcc site, using PBE, HSE, BLYP and B3LYP. For B3LYP and BLYP, they are evaluated using the HSE equilibrium geometries. Experimental values are also shown with error bars. [44] We will first concentrate on the PBE results, which have been published in a similar form by Gajdos et al. [5] For atop adsorption, the adsorption energy progressively decreases towards the noble metals. The qualitative and to some extend quantitative behavior of the atop curve can be reproduced using the Hammer-Nørskov d band model, [39,45] in which the central parameter is the position of the d band with respect to the Fermi-level. It moves to lower energies as the d band becomes filled. A word of caution is in order here. Using our definition, the d band moves to higher energy as the d band filling increases (table II). This is not in contradiction to previous findings, [39] since the different behaviour is related to a different definition of the d band center. In Ref. 39 the full d band is considered, while here we focus only on the occupied portion. The same trends as found in table II have been reported in Ref. 5. As we are mainly interested in trends, our following discussions are not affected by the particular definition. We note, however, that models based on the weak interaction limit, which forms the basis for the derivation of the Hammer-Nørskov approach, possess only limited value in the case of the strong interaction present for CO adsorbed on a transition metal surface, in particular, at the hollow site. [2] In fact, the d band model describes well the overall trends in the periodic table, but it does not account for the site dependent geometric effects (i.e. top vs. hollow). A convenient way to present the results for the 4d and 5d transition metals is to show all values in a single graph and to choose as the abscissa the occupation of the d band (figure 3) according to the previous definition. Interestingly for the top site, all adsorption energies fall on a single curve, suggesting that for top site adsorption the filling of the d band and not the position of the d band is the most sensible parameter (descriptor) to characterize the interaction strength. [2] Figure 3 clearly shows a decreasing preference for atop binding in the sequence Os<Ru∼Ir<Rh<Pt<Pd. These findings agree with experimental evidence within experimental errors (see Ref. 44 for further details) and other theoretical calculations. [46] Furthermore, for many properties, for instance catalytic properties, Pt behaves similar to a fictitious element between Rh and Pd, as the d band filling is roughly in-between the latter two elements. According to conventional wisdom, the 5σ interaction dominates for atop coordination. Hence, the decrease of the adsorption energy with increasing n d filling can be explained in terms of occupation of antibonding 5σ − d metal states that destabilize the molecular surface bond.
Quantum chemical calculations have shown that the back donation mechanism is the most dominant for the hollow site (2π * −d interaction). [47,48] We find, in agreement with previous studies, [5] that the hollow site adsorption energy generally increases towards the late transition metals, with a final sharp decrease for the noble metals. This decrease for the noble metals is related to a fully occupied d band, located several eV below the Fermi-level. Therefore, for noble metals it becomes difficult to shift antibonding 5σ − d states above the Fermi-level, and additionally the interaction between the deep and contracted d orbitals and the 2π * frontier orbitals becomes small. It turns out that sensible descriptors are the work function and the d band filling. In fact, previous studies have shown that the back-donation interaction is very sensitive to the relative position of the 2π * orbitals and the d band, which also depends on the metal work function Φ. [46,49] Clearly, a larger Φ value at fixed or almost constant d filling (i.e. Os vs Ru; Ir vs Rh; Pt vs Pd) implies a downward shift of the Fermi energy, decreasing the back-donation into the 2π * − d bonding orbitals, thereby decreasing the fcc adsorption energy. As a result, the CO adsorption at the fcc site is less stable on Os, Ir and Pt than on Ru, Rh and Pd. When increasing the d filling, that is moving from left to right along the periodic table, the back-donation is also strengthened: since the bonding 2π * − d orbitals become increasingly occupied, the molecular surface bond is stabilized. Counteracting this trend is the change in the work function, which increases towards the late transition metals. Finally, steric repulsion, which is always stronger at the hollow site than at the atop site, also plays a role. It decreases when moving to the right of the periodic row. [50] Although we found it impossible to disentangle the aforementioned factors, the calculations show that the combined effect of these trends is that hollow site adsorption becomes stronger towards the late transition metals. Last, but not least, we want to note that the different qualitative trends for the atop and fcc adsorption energies as a function of the d filling reflect the different dominant interaction term at the atop and fcc site, that is donation and back-donation, respectively.
We now return to figures 1 and 2. It is clearly born out that the top site is preferred for early transition metals, but moving along the transition metal series the hollow site becomes progressively more stable. For the PBE functional, a change in site preference is observed between Ru and Rh for the 4d metals, and Ir and Pt for the 5d metals, earlier than observed in experiment. The noble metals behave somewhat differently. As already argued above, this is related to a very weak 5σ → d donation and weak d → 2π * back-donation. Here bonding is dominated by the s and p orbitals, but again the PBE functional fails to predict the right site order. A 0.1 eV shift of the fcc and top curves against each other would give the correct site preference for all metals, e.g. top for Rh and Pt, but maintaining hollow site adsorption for Pd

D. The hybrid HSE functional
The hybrid HSE functional shifts the 2π * orbitals towards the vacuum level and the 5σ orbital towards stronger binding energies. [18] In combination with a downshift of the d band center, caused by a reduced self-interaction within the d shell due to the inclusion of Fock-exchange, one might naively believe that the interaction between the 2π * orbitals and the metal d states should weaken for hybrid functionals. Although this conjecture is confirmed for the noble metals, for transition metals the weaker interaction is counteracted by the much increased d band width of the metal and a concomitant decrease of the work-function from PBE to the HSE functional (see table II).
It is clear from both, figures 1, 2 and 3, that the adsorption energy increases for most adsorption sites using HSE: (i) the increase is pronounced and roughly constant for the top site on the transition metals; (ii) for the hollow site a slight destabilization is visible for Os and Ru, and a slight stabilization is visible for Pd and Pt. It is remarkable how smooth the changes along the transition metal series are, supporting the high precision of the present calculations and excluding numerical artifacts. As a result of these shifts, the crossing points between the hollow and top site move towards the right, for both 4d and 5d metals, but not sufficiently so to favor atop adsorption on Pt(111). For the other critical cases, i.e. Rh, Ag (and Cu, not included here but discussed in Ref. 18), HSE is capable to recover the correct site order, but only just. The main reason why HSE fails to predict CO adsorption above the top Pt site is that the hollow site also becomes more stable for the late transition metals, in particular Pd and Pt, whereas the increase in stability is small at the hollow site for earlier transition metals. This increased stability may be related to the tendency of HSE to increase surface energies (see below) and the tendency to localize the charge at the metal atoms, causing charge depletion at the hollow sites. In figure 4 we show the difference between the PBE and HSE charge density, plotted for a (111) plane ∼ 0.7Å above the metal atoms of the clean surface, that is, roughly in between the CO molecule and the metal surface layer. The figure indeed indicates charge depletion above the hollow site (red regions), suggesting a decrease of the steric repulsion using HSE. Furthermore, figure 4 suggests that, going from PBE to HSE, the hcp site should increase more in stability than the fcc site. This is confirmed by our previous calculations for Rh [18] that gave a stabilization of ∆E PBE−HSE fcc = 0.007 eV for the fcc site, and ∆E PBE−HSE hcp = 0.027 eV for the hcp site.
The most unsatisfactory result is that the adsorption energies increase substantially for the transition metals from PBE to HSE, which leads to fairly strong disagreement with experiment on the HSE level-even more so than for the PBE functional. We note that the same site preference, but even larger adsorption energies are obtained using the standard PBEh hybrid functional (also termed PBE0 or PBE1PBE functional). This functional also includes 25 % non-local exchange but, contrary to the HSE functional, does not truncate the long range part of the non-local exchange. [18] Since the only difference between the PBE and PBEh functional is that the latter one replaces 25 % of the local approximation to the exchange by the Hartree-Fock exchange, the conclusion is that inclusion of 25 % non-local exchange, bare or long-range screened, hardly improves the description of CO adsorption on metal surfaces.
This confirms our previous study, [18,19] but contradicts a recent work, where it is claimed that the CO adsorption FIG. 4: Difference between the PBE and HSE charge density for Rh and Pt clean surfaces plotted in a (111) plane ∼ 0.7Å above the surface plane (i.e. between the C atom of CO adsorbed at the hollow site and the surface metal atoms). Blue (red) color corresponds to regions where the PBE charge density is larger (smaller) than the HSE one. Three isolines are drawn to guide the eyes. "fcc" and "hcp" indicates the location of the fcc and hcp hollow site, respectively.
puzzle has been solved for Pt using hybrid functionals. [51] In order to explain the contradiction, we have calculated the adsorption energies using norm-conserving and PAW potentials and we present the results in the Appendix. We show that the conclusions of Ref. 51 derive from an artifact of the pseudopotential approximation, and all-electron methods still yield the wrong site order for CO on Pt(111) using the PBEh functional.

E. The hybrid B3LYP functional
Our results also somewhat disagree with other reports, which have applied the B3LYP functional and found stronger preference for atop adsorption and weaker adsorption energies on Cu and Pt. [12,52,53,54] To resolve this issue, B3LYP values are reported in figure 2. In fact, the B3LYP functional shifts the crossover point even further to the right and lowers the adsorption energies significantly compared to the HSE and PBE functional. Since these "improvements" can not be related to the inclusion of non-local exchange (see above), the difference must stem from the semi-local part, and in all likelihood the Lee Yang Parr correlation energy (see e.g. Ref. 35). In agreement with this conjecture, the semi-local BLYP functional also predicts the right site order on all considered metals. In summary, improvements from the PBE to the B3LYP functional are mostly related to the LYP correlation functional, and only partly to the inclusion of non-local exchange, confirming the recent work by Alaei et al. [6] The seeming "improvements" for the BLYP and B3LYP functional over PBE and HSE, however, come at a considerable price best illustrated in table I. The equilibrium volume of Ag is overestimated by 12 %, that of Pt and Pd by about 9 %, an error that is by no means acceptable. Somewhat less troublesome but still worthwhile mentioning is that adsorption on the noble metals becomes slightly endothermic using B3LYP and BLYP (no geometry optimization was performed for the B3LYP case). We do not want to overemphasize this point, since presently DFT functionals generally lack weak van der Waals like dispersion forces, which might contribute to the CO bonding on noble metals.

IV. DISCUSSION
It has recently been suggested that the quest for improved semi-local functionals is not yet over, and that it might be possible to design semi-local functionals that better describe the adsorption on surfaces. [6] Unfortunately this is a misconception, more precisely, for semi-local functionals, it seems impossible to improve the description of adsorption without worsening other properties. To discuss this issue we have extended our study to three additional functionals, the rPBE functional by Hammer and Nørskov, [58] and the recently suggested PBEsol [60] and AM05 [59] functional. The last functional is based on the subsystem approach and designed to yield accurate surface energies. The other two functionals are modifications of the standard PBE functional. These modifications are in the opposite direction: the PBEsol functional follows similar ideas as the AM05 functional in order to improve surface energies, whereas the rPBE functional improves formation energies and adsorption energies.
Most gradient corrected functionals essentially modify the LDA exchange and correlation energy to account for areas of large gradients. This is usually done by increasing the exchange and correlation energies locally in those areas where the gradient is large. For the "canonical" PBE functional, the enhancement factors have been determined using sum-rules, and the result is a well balanced functional that describes molecular and solid state properties equally well. There are however failures, in particular small molecules are over-bound, whereas solids involving elements beyond the second row are usually under-bound. Resultantly, for the metals considered here the lattice constants are too large and the surface energies much too small (underbinding). To highlight this, figure 5 plots the surface energy for all considered functionals versus the theoretically predicted lattice constants of Rh and Pt. Similar results are to be expected for the other transition metals.
The AM05 functional has been designed to describe jellium and jellium surface energies accurately. The PBEsol functional mimics the AM05 functional by reducing the enhancement factor in the area of strongly varying gradients, thus approaching the LDA description. [60] As a rule of thumb, both increase the cost for surfaces, and hence systems involving surface areas become less favourable. Atoms are destabilized, the formation energies of molecules increase (molecules possess a smaller surface area per atom than atoms) and the atomization energies of bulk materials become much larger (no surface area). Furthermore bulk lattice constants shrink to reduce "internal" surfaces. The AM05 and PBEsol functionals shift the focus towards a proper description of solids, but in doing so small molecules are strongly overbound compared to atoms using the PBEsol and AM05 functional. [60] In figure 5 (a), the AM05 and PBEsol functionals yield smaller lattice constants and larger surface energies than the PBE functional; both are now in quite reasonable agreement with experiment.
The rPBE functional does exactly the opposite: it increases the enhancement factor in the area of strongly varying gradients, and it therefore moves further away from the LDA description. As a rule of thumb, the rPBE functional decreases the cost for surfaces. Therefore, molecules are destabilized compared to atoms, bulk materials become less stable than molecules, and bulk lattice constants increase. The decrease in the surface energy and increase in the lattice constant is again visible in figure 5 (a). The BLYP functional, characterized by an even larger enhancement factor (see e.g. Ref. 6), continues this trend.
We can now understand, why the adsorption energy decreases from AM05 and PBEsol, over PBE, rPBE to BLYP [ Figure 5 (b)]: when CO adsorbs on the surface, the total surface of the combined system CO + surface is effectively smaller than that of the isolated CO and the bare surface. [58] In AM05 and PBEsol, surfaces are expensive, and thus adsorption processes, which always lead to a reduction of the effective surface area, yield a larger energy gain. For rPBE and even more so for BLYP, surfaces are cheap, and less energy is gained during adsorption processes.
The change in the energy difference between the top and the hollow site [ figure 5 (c)] can be understood by reasoning along the same line. Atop adsorption reduces the surface area very little essentially along a single bond between the CO and the surface, whereas hollow site adsorption reduces the surface area along three bonds between the three neighbouring surface atoms and the CO. If surfaces are "expensive" (PBEsol), hollow site adsorption is preferred, since hollow site adsorption reduces the effective surface area most. If surfaces are "cheap", top site adsorption is preferred.
The binding energies, surface energies, lattice constants, and adsorption energies are thus all linearly related, and we could plot any quantity versus the other and we will always observe a linear relationship. In figure 5 we have chosen the lattice constant as independent variable, but the surface energy would work equally well (and arguably might be the better control parameter). Within the family of gradient corrected functionals, PBE seems to remain the best compromise. The AM05 and PBEsol functionals describe the lattice constants of the solids best and also yield improved surface energies, but unfortunately they overestimate the heats of formation, and important for the present study, the binding energies of molecules on surfaces. The BLYP functional arguably yields the best adsorption energies and the right site order, but for the price of crossly wrong bulk energies, much too large lattice constants, and much too small surface energies.
The hybrid functionals constitute a deviation-albeit only a subtle one-from this general trend. The inclusion of exact exchange also clearly increases the cost for surfaces. Resultantly, HSE leads to a smaller lattice constant than PBE, and B3LYP to a smaller lattice constant than BLYP. Even the average adsorption energies follow roughly the same trend as for semi-local functionals: in particular the adsorption energy for the top site lies on the same curve as for the semi-local functionals. The only major deviation is visible for the energy difference between the top and hollow site. This we believe to be related to the shift of the 2π * LUMO orbital towards the vacuum level, decreasing the capability of the 2π * orbital to accept charge from the substrate. As we have argued before this generally decreases the adsorption energy, with a stronger effect for the hollow site than the for top site. On passing, we note that on the basis of surface energy arguments alone, one would expect hybrid functionals to yield similar atomization energies for molecules and solids as a semi-local functional predicting similar lattice constants. This would imply that the atomization energies of molecules were significantly overestimated and that of solids slightly overestimated for the HSE functional. However, as we have reported previously, this is not the case, i.e. the atomization energies of molecules agree quite well with experiment, [34] whereas those of solids are underestimated. [26,35] The reason for this is that the inclusion of the exact exchange stabilizes spin polarized solutions compared to non-spin polarized solutions, e.g. spin polarized atoms are much more stable for hybrid functionals than for semi-local exchange correlation functionals. For sp elements and d metals this leads to a slight and sizeable underestimation of the atomization energies, respectively. [26] Neglecting these important spin-polarization effects, hybrid functionals do behave reasonably similar to semi-local functionals, with the tendency of increasing the cost for surfaces compared to the corresponding semi-local functionals.
Clearly, however, none of the functionals can predict the surface energies and the adsorption energies well at the same time. From our point of view, the good prediction of the B3LYP and BLYP for adsorption energies and the adsorption site is not related to a better description of the underlying physics but rather accidental and fortuitous. By reducing the surface energy to an, in fact, unrealistically low value-worsening the lattice constants-the right site order is eventually recovered, but there is little doubt, the reduction of the surface energy is not in accord with experiment. We believe that the main error of semi-local functionals is in the description of the strength of the back donation to the 2π * orbital which tends to be overestimated using any of the semi-local density functionals. Hybrid functionals improve on that aspect, but unfortunately simultaneously increase the d band width and concomitantly the overall adsorption energies.

V. SUMMARY AND CONCLUSIONS
Overall, the situation is unsatisfactory. Considering CO adsorption alone, the B3LYP and BLYP functionals seem to offer a very good description, but the price to pay is that the description of the metals, in particular d metals, is unsatisfactory. Errors of 10 % in the volume are unacceptable by today's standards, and worse, the atomization energies of d metals are wrong by up to 50 % for the B3LYP functional. [35] Previously published seemingly good B3LYP results [12,52,53,54] must be reconsidered, when the broader picture is taken into account. B3LYP and BLYP fail for metals and surface energies: we feel there is hardly a point to predict surface properties, if the underlying substrate and surface is not accurately and properly described.
The HSE functional performs significantly better for metals, yielding a better overall description of the equilibrium lattice constants than most semi-local functionals. But again, atomization energies for transition metals are wrong by typically 20 %, [26] and the overall description of CO on transition metal surfaces improves only little compared to the PBE functional. In particular, the increased adsorption energies are unsatisfactory. Overall, we are forced to conclude that for metals and metal surfaces, hybrid functionals are hardly a major step forward. In accord with our previous conclusions, [18] we believe that the inclusion of a sizable amount of non-local exchange is not capturing the physics properly in metals. The non-local exchange interaction is strongly screened in metals, both in the long range as well as medium range. In combination with a semi-local correlation functional, the hybrid HSE06 functional overestimates the exchange interactions at medium distances. [18,26] Finally, we have shown that trends from one to the other semi-local functional can be easily understood considering the cost for creating a surface. For the recently suggested AM05 and PBEsol functionals, surfaces are, in agreement with experiment, fairly expensive favouring bulk-like behaviour, whereas rPBE and BLYP underestimate surface energies and shift the stability towards smaller fragments (molecules and atoms). It is really puzzling that AM05 and PBEsol yield good bulk lattice constants and surface energies but overestimate adsorption energies, whereas rPBE and BLYP underestimate binding in solids and surface energies but describe adsorption energies well. A functional that could resolve both issues would be a tremendous step forward but we firmly believe that semi-local functionals can not achieve this goal.

VI. ACKNOWLEDGMENTS
This work was supported by the Austrian Fonds zur Förderung der wissenschaftlichen Forschung. A.S. thanks J. Paier for interesting and fruitful discussions. Figure 4 has been done using the XCrySDen package. [61]

VII. APPENDIX. EFFECT OF PSEUDOPOTENTIAL APPROXIMATION IN HYBRID-FUNCTIONAL CALCULATIONS
In a recent work [51] it was suggested that the proper site preference for CO on Pt(111) can be obtained using the PBEh hybrid functional. This observation contradicts our own finding, [18] and we show here that the previous result is an artefact related to the pseudopotential approximation applied in Ref. 51. To understand the discrepancy between our own work and that published in Ref. 51, we generated a norm conserving pseudopotential (NC-PP) mimicking the FHI NC-PP applied in Ref. 51. Figure 6 shows the density of states (DOS) evaluated using this normconserving pseudopotential and the PAW method. In our implementation, the PAW method has several advantages over NC-PPs. First, the exact shape of the all-electron orbitals is restored avoiding any shape approximations in the charge density, potentials and non-local exchange interaction. Second, although the PAW method applied here is using frozen core orbitals imported from a GGA calculation for the atom, the exchange and correlation energy is re-evaluated exactly and consistently using all electrons (core and valence), whereas in the pseudopotential approximation the hybrid functional is applied to the valence electrons only, and the core-valence interaction is implicitly calculated on the GGA level. Clearly, both approaches result in a different density of states for bulk Pt using the PBEh functional (see figure 6). We note that no difference is observed on the GGA level, for which the potentials were initially generated. Furthermore, we observe that the PAW results are entirely robust and do not change upon treating the 5s or 5p shell  as valence (not shown), whereas placing the 5s and 5p shell in the valence for the NC-PP (NC−5sp in figure 6) yields results that are practically identical to the PAW results. The implicit treatment of the core-valence interaction on the PBE level for the standard NC-PP causes fairly large errors.
In the second step we evaluated the CO adsorption energy at the fcc and top site for a √ 3 × √ 3 slab using these three potentials and 4 × 4 k-points (similar technical parameters were applied in Ref. 51). The results are summarized in table III. At the PBE level, all three potentials give very similar results to those of Ref. 51. Large differences are observed at the PBEh level. The standard NC potential yields the same site order as in Ref. 51 (but different adsorption energies), whereas the NC−5sp potential-treating the 5s and 5p electrons explicitly as valence -and the PAW method show that the top and fcc site are practically degenerate for this setup. For a lower coverage using a c(2 × 4) surface cell and using more k points, the fcc site is consistently preferred by 0.04 eV for the PBEh functional. We conclude: the PBEh functional does not predict the right site order on Pt(111) in the low coverage limit. Also the adsorption energy is significantly overestimated by the hybrid PBEh functional, and overall the description is thus not improved over the GGA case. [18] Furthermore, normconserving transition metal potentials generated using LDA/GGA functionals must be used with great care in hybrid functional calculations. We believe that this rule also applies to elements beyond the third row (4p, 5p and 6p elements), where the core valence overlap is often appreciable.