Linear scaling relations have led to an understanding of trends in catalytic activity and selectivity of many reactions in heterogeneous and electro-catalysis. However, linear scaling between the chemisorption energies of any two small molecule adsorbates is not guaranteed. A prominent example is the lack of scaling between the chemisorption energies of carbon and oxygen on transition metal surfaces. In this work, we show that this lack of scaling originates from different re-normalized adsorbate valence energies of lower-lying oxygen vs higher-lying carbon. We develop a model for chemisorption of small molecule adsorbates within the d-band model by combining a modified form of the Newns–Anderson hybridization energy with an effective orthogonalization term. We develop a general descriptor to a priori determine if two adsorbates are likely to scale with each other.

Understanding the formation of a bond between a metallic surface and an adsorbate is critical to the field of heterogeneous catalysis.1–10 In principle, the bond (free) energies of all intermediates and transition states of a process fully determines the rate, but, in most cases, the number of such bond energies is so large that it becomes very difficult to single out why a given catalyst surface is the optimal catalyst for the reaction studied.11 Scaling relations have provided a method to define a few bond energies as descriptors of the catalytic rate.11–14 It has been found quite generally that, for a given surface structure, the bond energies of intermediates and transition states scale linearly with a few bond energies, and once these relations have been found, the full kinetics can be described as a function of the descriptors. These scaling relations result in the “volcano” relationships between rates/selectivities and descriptors that are used extensively in heterogeneous catalysis. It provides a theoretical basis for the heuristic Sabatier principle,15 which suggests that the optimal catalyst is one with optimal bond energies for relevant intermediate(s). The scaling relations show which intermediate(s) are relevant and the optimal value(s) for the binding energy of these intermediate(s).

Linear or piece-wise linear scaling relationships are observed when two adsorbates bind to the metal through the same atom, such as C*, CH*, CH2*, and CH3*,16,17 where * denotes the surface. However, it is not always observed when the binding atom is different, such as in the case of C* and O*.11 There is currently no a priori measure to suggest if the chemisorption energies of two arbitrary adsorbates scale with one another.

In this work, we develop a model that allows for an a priori estimation of whether the chemisorption energies of two adsorbates scale with each other. We use the Hammer–Nørskov model, which suggests that the trends in adsorption energies of transition metals are related to the bonding between the adsorbate valence states and metal d-bands.5 The hybridization energy is described using a modified form of the Newns–Anderson model,18 combined with an effective expression for the orthogonalization energy. The model demonstrates that the root cause of lack of scaling between two adsorbed species lies in differences of their re-normalized adsorbate valence energies. For instance, O – p states are considerably lower in energy than C – p states. Through this analysis, we identify the conditions under which carbon and oxygen do scale with each other.

In this section, we discuss the lack of linear scaling relations between the chemisorption energies of C* and O*. We use density functional theory (DFT) calculations to obtain these chemisorption energies (see Sec. VI for computational details). In general, for any adsorbate A, the chemisorption energy ΔEA is defined as

ΔEA=EA*E*EA,

where A* represents A adsorbed on the surface, * is the free surface, and A is the adsorbate in vacuum. Figure 1(a) shows the chemisorption energies of O* plotted against those of C*. It is evident that there is a significant amount of scatter and no linear scaling relations are observed.

FIG. 1.

(a) Scatter plot of DFT chemisorption energies of O, ΔEO and C, ΔEC; the color refers to the row in the Periodic Table, 3d metals are red, 4d metals are blue, and 5d metals are yellow. Note that there is no scaling between the adsorbates (R2 = 0.56); scaling of (b) ΔEO and (c) ΔEC with the d-band center. Same points as in part (a), splitting the energies into (d) 3d, (e) 4d, and (f) 5d metals. The adsorbate is always placed at the most stable site, which is at hollow for all transition metals.

FIG. 1.

(a) Scatter plot of DFT chemisorption energies of O, ΔEO and C, ΔEC; the color refers to the row in the Periodic Table, 3d metals are red, 4d metals are blue, and 5d metals are yellow. Note that there is no scaling between the adsorbates (R2 = 0.56); scaling of (b) ΔEO and (c) ΔEC with the d-band center. Same points as in part (a), splitting the energies into (d) 3d, (e) 4d, and (f) 5d metals. The adsorbate is always placed at the most stable site, which is at hollow for all transition metals.

Close modal

Despite the lack of linear scaling, there is broadly an increasing trend for both C* and O* energies in Fig. 1(a). Metals that bind C* strongly (for example, Fe) also bind O* strongly and those that bind C* weakly (Au, Ag) also bind O* weakly. However, due to the noticeable amount of scatter, this general “rule” does not always hold. Consider the example of platinum: For a constant ΔEC, ΔEO is about 2 eV weaker than what linear scaling would predict. This weaker than expected binding affects the catalytic activity of Pt, where it is considered weak binding for oxygen based reactions19,20 but has a lot of interesting carbon chemistry.12,17,21–28 Furthermore, the noble metals (Au, Ag, and Cu) appear to tail-off, having significantly weaker O* and C* chemisorption energies.

Figures 1(b) and 1(c) show the chemisorption energies of O* and C* plotted against the d-band center, a commonly used descriptor in heterogeneous catalysis.5 The tailing-off behavior of the chemisorption energies at large negative ϵd values is similar to that seen in Fig. 1(a).

Scatter in Fig. 1(a) is markedly reduced when each row of transition metals is plotted separately as shown in Figs. 1(d)1(f). Mid-row transition metals for each series display linear behavior (Fe, Co, and Ni for 3d; Ru, Rh, and Pd for 4d; and Os, Ir, and Pt for 5d), with the final elements of each row (Cu, Ag, and Au) falling off a linear fit.

These observations raise two fundamental questions for understanding scaling relations:

  1. Why do the noble metals (Cu, Ag, and Au) deviate from the linear scaling relationships constructed along a given row of transition metals [such as in Figs. 1(d)1(f)]?

  2. What is the source of the scatter in scaling relations across rows of the Periodic Table [such as in Fig. 1(a), deviation in red, blue, and yellow points]?

In the following section, we present a model of chemisorption, which is able to explain these observations.

In this section, we develop a model for the chemisorption of mono-atomic adsorbates that includes hybridization and orthogonalization contributions. Following the Hammer–Nørskov d-band model,5 we split the total chemisorption energy into the aforementioned components as

Echem=EhybHybridization+EorthoOrthogonalization+constant.
(1)

The hybridization component of the energy, Ehyb, is incorporated through the Newns–Anderson model

Ehyb=2πϵf=0arctanΔd+Δ0ϵϵaΛdϵ2ϵa,
(2)

where Δd represents a continuous, ϵ-dependent coupling element, Vak, between the metal d-density of states and the adsorbate state. In line with Ref. 18, it is taken to be a semi-ellipse centered around ϵd with a width wd. Δ0 is a constant contribution applied to approximate the broadening of the adsorbate states into resonances upon interaction with sp-states of the metal, Λ is the Hilbert transform18 of Δd + Δ0, and ϵa is the renormalized energy level of the adsorbate, after interaction with the sp states of the metal [see the supplementary material (Sec. S2) for further details] and ϵf is the Fermi level (set to 0). The factor of 2 in both terms represents degeneracy of the spin up and down states considered in the model. The adsorbate projected density of states is given by

ρaa(ϵ)=π1Δd+Δ0ϵϵaΛ2+(Δd+Δ0)2,
(3)

where the ϵ dependence of Δ and Λ has been omitted for clarity. The occupancy of the adsorbate state is given simply by na=ϵf=0ρaa(ϵ)dϵ.

We obtain the orthogonalization penalty to the chemisorption energy by comparing the eigenenergies from the narrow-band limit of the Newns–Anderson model with a two-level problem that includes a constant coupling term S,5,29

Eortho=2na+fVakS,
(4)

where f is the filling of the metal states, given by

f=0Δd+Δ0dϵΔd+Δ0dϵ,s.t.0<f<1.

Note that spin-polarization can be considered in Echem by performing a sum over individual spins in Eq. (2), as in Ref. 18, along with different Eortho values for each spin. In this work, we will instead deal with the spin-paired scenario, as our goal is to determine the factors that cause lack of scaling in Fig. 1.

With the above functional form for Ehyb and Eortho, we now parameterize the sum of the two, i.e., Echem [through Eq. (1)]. We choose the O – p and C – p states to have ϵa = −5 and −1 eV, respectively.30 The projected densities of states resulting from Eq. (3) with the chosen ϵa for the O* and C* also qualitatively agree with our DFT calculated projected density of states shown in Fig. 2(a) (see Sec. S12), where the localized O* states (shown in orange) are consistently lower-lying in energy as compared to those of C* (shown in purple). Note that the adsorbate projected density of states of both O* and C* show similar features to experimental x-ray photoemission spectra, such as localized states on either side of the d-band from Ref. 31 for Cu and Ni, which further validates our approach of using a single renormalized energy for O* and C*.

FIG. 2.

(a) DFT projected density of states for the d-states of the 4d transition metal surfaces (in black, analogous to Δ) and p-states of O (orange) and C (purple); chemisorption energy from the model, Echem, plotted against the DFT energies based on the parameters in Table I for (b) O* and (c) C*. Colors indicate the row of the element in the Periodic Table; red for the 3d, blue for 4d, and yellow for 5d transition metals. Chemisorption energies against ϵd for (d) O* and (e) C*. (f) Scaling relations between O* and C* chemisorption energies serve as a model analog to the DFT computed energies in Fig. 1(a).

FIG. 2.

(a) DFT projected density of states for the d-states of the 4d transition metal surfaces (in black, analogous to Δ) and p-states of O (orange) and C (purple); chemisorption energy from the model, Echem, plotted against the DFT energies based on the parameters in Table I for (b) O* and (c) C*. Colors indicate the row of the element in the Periodic Table; red for the 3d, blue for 4d, and yellow for 5d transition metals. Chemisorption energies against ϵd for (d) O* and (e) C*. (f) Scaling relations between O* and C* chemisorption energies serve as a model analog to the DFT computed energies in Fig. 1(a).

Close modal

As in Ref. 5, we assume that the coupling elements are proportional to those from the linear muffin tin orbital (LMTO) framework (denoted as Vsd),32 i.e., Vak2=βVsd2, where β only depends on the adsorbate species. Analogous to Ref. 33, we write Vsd as

Vsd=ηMaMdrla+ld+1,
(5)

where η are structure factors, which do not depend on bond-length between adsorbate and metal, Ma depends on only the adsorbate states.34,Md=s2ld+1Δl12, where ld and la are the quantum numbers of the metal d-states (ld = 2) and adsorbate p states (la = 2), respectively, s is the Wigner–Seitz radii, and Δl is the Anderson-width (taken from bulk calculations).35 In this work, we set r = s to incorporate the varying bond length between the adsorbate and metal into the coupling elements (see Sec. S3 for parameters).

We also assume that the overlap matrix element is proportional to the coupling element, S = −αVak, where α is another constant that is identical for all metals.

The only free parameters in the model for a given adsorbate are α and β, which are independent of the metallic surface species. We obtain these parameters by fitting Echem to the DFT energies in Fig. 1. Table I shows the least-squares errors fitted α and β for C* and O*. Note that Echem is multiplied by 3 to account for three bonds due to the preferred adsorption site being the hollow site (see Sec. S4 for ontop site fitting). Figures 2(b) and 2(c) show parity plots for fitted chemisorption energies from the model vs computed DFT energies. Fitting parameters for CO* based on this model, and a comparison with Ref. 5 is shown in Sec. S5.

TABLE I.

Fitting parameters for C* and O* from a least-squares error minimization routine.

Adsorbateα (eV−1)β (eV2)Δ0 (eV)
C* (ϵa = −1 eV) 0.094 1.44 0.1 
O* (ϵa = −5 eV) 0.062 4.55 0.1 
Adsorbateα (eV−1)β (eV2)Δ0 (eV)
C* (ϵa = −1 eV) 0.094 1.44 0.1 
O* (ϵa = −5 eV) 0.062 4.55 0.1 

Note that the simultaneous addition of hybridization and orthogonalization components to Echem has already been incorporated through the Newns–Anderson–Grimley model of chemisorption36 and, more recently, has been applied to study chemisorption and catalytic reactions.10,37 Within this framework, Echem is obtained from the Anderson Hamiltonian using a non-orthonormal basis, i.e., the obtained energies from the model implicitly include the effects of overlap.

In this work, we instead choose a simpler representation of Echem through Eq. (1) and subsequent parameterization of α and β to describe trends in chemisorption. As mentioned in Ref. 36, the same Echem may be obtained with and without explicit consideration of overlap, if Vak is changed accordingly. By parameterizing Eq. (1), as done in this section, we effectively alter Vak for each adsorbate (through β) while accounting for the repulsive interactions through Eortho. This parameterization allows us to model trends with a simpler model which functions in the same way as the Newns–Anderson–Grimley model. Equivalency between the two models is demonstrated in Sec. S6, where the parameterization procedure described above is performed with the Newns–Anderson–Grimley model instead of Eq. (1), yielding qualitatively similar behavior in terms of trends in variation of Echem, Ehyb, and Eortho with ϵd.

Furthermore, while the orthogonalization model presented in this section is simple and parameterizable, in Secs. S6 and S7, we show that other orthogonalization models yield near identical behavior of Echem, with qualitatively similar behavior of variation along ϵd for Ehyb and Eortho.

Having fit α and β to chemisorption energies from DFT, we compute the variation of Echem for O* and C* with ϵd for 3d, 4d, and 5d series of transition metals shown in Figs. 2(d) and 2(e), respectively. Note that, for a fixed row, the chemisorption energy is a function of only ϵd, as Vak2 and wd for each row of the Periodic Table are functions of ϵd.

The model reproduces two important features present in the DFT calculated energies. First, the variation of Echem for O* and C* with ϵd in Figs. 2(d) and 2(e) shows nearly linear behavior up to ϵd ≈ −3 eV. At more negative ϵd values, Echem decreases, consistent with the DFT computed trends in Figs. 1(b) and 1(c). Second, the lack of scaling between different rows of transition metals is also seen in the model chemisorption energies. Figure 2(f) shows Echem for O* and C*, plotted against each other. Similar to the trends from DFT calculations in Figs. 1(b) and 1(c), taken together the points do not scale, but they do scale when each transition metal row is plotted separately, up to the weak binding noble metals.

Having developed and parameterized the model, we now use it to explain the variation of the DFT energies with ϵd in Sec. IV.

The lack of scaling between Echem of O* and C* arises from differences in how Ehyb and Eortho of the two adsorbates vary with respect to ϵd. Previous studies have demonstrated that either hybridization5 or orthogonalization38 dominate for certain adsorbates and metal surfaces. In this section, we generalize these observations through the model in Sec. III and describe (1) why noble metals Cu, Ag, and Au fall off otherwise linear scaling relationships, (2) why scaling of transition metals is better performed across a fixed row. These two observations come directly from computed energies in Sec. II. We trace the origin of both these effects to differences in ϵa values of the two adsorbates.

Figures 3(a) and 3(b) show Ehyb of O* and C* vs ϵd, which has three distinct regions. At positive ϵd values, Ehyb increases approximately linearly with ϵd. At intermediate ϵd values (−2.5 < ϵd < 0 eV for 4d and 5d metals and −1.6 < ϵd < 0.9 eV for 3d metals), Ehyb decreases with ϵd. Finally, at large, negative ϵd values (−3.5 eV for 4d and 5d metals and −1.9 eV for 3d metals), Ehyb → 0, i.e., it saturates completely for O*, while just starting to saturate for C*.

FIG. 3.

Ehyb for (a) O* and (b) C*; there is a clear saturation, i.e., Ehyb → 0 in the case of O*, but only starting to saturate for C*. Ehyb is shown as dashed lines; star denotes ϵs. (c) Scaling between Ehyb between O* and C* where more opaque points indicate more negative ϵd values and the direction of the annotated arrow indicates increasing ϵd; Eortho for (d) O* and (e) C*; dashed lines represent na + f (f) scaling of Eortho between O* and C* where more opaque points indicate more negative ϵd values and the direction of the annotated arrow indicates increasing ϵd; (g) adsorbate projected density of states for (orange) O* and (purple) C* with Δ (black); all quantities are normalized by their maximum value; na for both adsorbates for Rh are annotated.

FIG. 3.

Ehyb for (a) O* and (b) C*; there is a clear saturation, i.e., Ehyb → 0 in the case of O*, but only starting to saturate for C*. Ehyb is shown as dashed lines; star denotes ϵs. (c) Scaling between Ehyb between O* and C* where more opaque points indicate more negative ϵd values and the direction of the annotated arrow indicates increasing ϵd; Eortho for (d) O* and (e) C*; dashed lines represent na + f (f) scaling of Eortho between O* and C* where more opaque points indicate more negative ϵd values and the direction of the annotated arrow indicates increasing ϵd; (g) adsorbate projected density of states for (orange) O* and (purple) C* with Δ (black); all quantities are normalized by their maximum value; na for both adsorbates for Rh are annotated.

Close modal

The difference in the variation of Ehyb with ϵd at large, negative ϵd (−3.5 eV for 4d and 5d metals and −1.9 eV for 3d metals) gives rise to scatter in linear scaling relations between O* and C* along a fixed row of transition metals. For O*, the variation of Echem in this large, negative ϵd region is determined solely by Eortho as Ehyb → 0. For C*, all regions of ϵd are jointly determined by the contributions of Ehyb and Eortho. This mismatch of energy contributions, Eortho alone for O* and Ehyb + Eortho for C*, explains why metals such as Cu, Ag, and Au in the large, negative ϵd region fall off linear scaling relationships.

To quantify the ϵd value at which Ehyb saturates, we compute the derivative of Ehyb with respect to ϵd, Ehyb, and determine the largestϵd value at which it tends to zero, which we denote as ϵs,

ϵs=argmax(ϵd)forEhyb(ϵd)0Ehyb0,
(6)

where Ehyb is given by

Ehyb(ϵd)=dEhybdϵd(ϵd)=2πϵf=0ΔdϵϵaΛ+ΛΔd+Δ0ϵϵaΛ2+Δd+Δ02dϵ,
(7)

where Δd(ϵ,ϵd)=Δd/ϵd and Λ′(ϵ, ϵd) = Λ/∂ϵd (see Sec. S8).

Figures 3(a) and 3(b) (dashed lines) show Ehyb for both O* and C*, respectively. Consistent with the aforementioned trends in Ehyb, there are three regions of Ehyb. At positive ϵd values, Ehyb is negative, in line with linear increasing variation of Ehyb. At intermediate ϵd values, Ehyb starts to increase, eventually vanishing at large, negative ϵd. ϵs is indicated by a starred point for each row (in red, blue, and yellow for 3d, 4d, and 5d, respectively). A larger, more positive value of ϵs indicates that saturation of Ehyb is expected to occur at more positive values of ϵd, as seen in O* as opposed to more negative values for C*.

Ehyb alone cannot explain the row dependence (3d, 4d, and 5d) of the scaling relations seen in Fig. 1. Figure 3(c) shows the scaling between Ehyb for O* and C*. For almost the entire range, there is an overlap between the three rows of transition metals (with the exception of the early transition metals), implying that this quantity does not cause the row dependence of scaling seen in the computed energies of Fig. 1(a). To explain this row dependence, we study the other component of Echem, i.e., Eortho.

Figures 3(d) and 3(e) (solid lines) show Eortho as a function of ϵd and dashed lines show the sum of na and f, which are multiplicative factors to Eortho [Eq. (4)]. At positive ϵd values, Eortho increases with decreasing ϵd due to the increase in Vsd2, na, and f. For ϵd < 0 eV, Eortho decreases with decreasing ϵd due to a decrease in Vsd2, though na + f still tend upward.

Differences in the behavior of na + f with ϵd between −3 < ϵd < 0 for O* and C* cause row-dependence of Eortho and, hence, Echem. Figure 3(f) shows the scaling between Eortho of  O* and C*; noticeable deviations are seen in linear scaling of Eortho but also in the overlap between the rows of transition metals.

Furthermore, differing behavior of na + f with ϵd is a consequence of different ρaa(ϵ) [determined from Eq. (3)] for different adsorbates, i.e., different ϵa. As an example, Fig. 3(g) shows ρaa(ϵ) for C* (purple) and O* (orange) on Rh, a mid-row transition metal. Since ϵa for C* is more positive than O*, the adsorbate projected density of states broadens and splits differently. In the case of C*, ρaa moves toward the upper edge of the d-band, resulting in a larger fraction of its peak lying above the Fermi level as compared to O*. This disparate ρaa behavior results in different na, causing differences in Eortho [from Eq. (4)] and, hence, row-dependent scaling.

The variations in Eortho at ϵd ≈ −2 eV are also caused by changes in ρaa, as the peaks start to cross the Fermi level. We note that this energy variation is less prominent when using the Newns–Anderson–Grimley model, as seen in Sec. S6. However, these variations in energy are of the order of 0.1–0.2 eV and do not change the findings of Fig. 3(f), where scaling of Eortho is different by ≈1 eV.

We now generalize our findings beyond C* vs O* scaling by proposing a descriptor to determine the ϵd value up to which scaling is expected for a fixed row, i.e., regions where the chemisorption energies are dictated by both hybridization and repulsive interactions.

In general, if two adsorbates have similar ϵs values, they are likely to scale. Figure 4(a) shows ϵs for a range of re-normalized energy levels ϵa. Since each ϵs value corresponds to a specific ϵa value for each row of transition metals, it suffices to say that two adsorbates with similar ϵa are likely to scale for a given row of transition metal elements. This condition is satisfied for scaling relations between adsorbates binding to the metal with the same atom, such as C* vs CHx* (x = 1, 2, 3),16 COOH* vs CO*,39 or OOH* vs OH*.20 It does not hold for adsorbates with very different ϵs and corresponding ϵa, such as C* and O*.

FIG. 4.

(a) ϵs, a measure of the saturation in Ehyb in Figs. 3(a) and 3(b) for a range of adsorbate energy levels. (b) Eortho for different adsorbates, C* (ϵa = −1 eV), N* (ϵa = −3 eV), S* (ϵa = −4 eV), and O* (ϵa = −5 eV) plotted against Eortho for O* at ϵa = −5 eV.

FIG. 4.

(a) ϵs, a measure of the saturation in Ehyb in Figs. 3(a) and 3(b) for a range of adsorbate energy levels. (b) Eortho for different adsorbates, C* (ϵa = −1 eV), N* (ϵa = −3 eV), S* (ϵa = −4 eV), and O* (ϵa = −5 eV) plotted against Eortho for O* at ϵa = −5 eV.

Close modal

Finally, Fig. 4(b) shows the effect that different ϵa has on the row-dependent scaling behavior. By plotting Eortho for O* (for ϵa = −5) against Eortho for S* (ϵa = −4 eV), N* (ϵa = −3 eV), and C* (ϵa = −1 eV), we determine the range of ϵa for which row dependent scaling is to be expected. For example, for N* (ϵa = −3 eV), Eortho for 3d, 4d, and 5d elements show less mismatch than for the case of ϵa = −1 eV (see Sec. S10 for ϵa for adsorbates), indicating that adsorbates with similar ϵa values are likely to not show row-dependence in their scaling behavior. The range of allowable ϵa values depends on the variation of Eortho across and along rows of the Periodic Table. Note that for adsorbates with near-identical ϵa (such as hydrogenated versions of the above monoatomic adsorbates), Fig. 4(b) essentially becomes a parity plot, leading to there being no row-dependence to scaling, as has been observed in a previous work.16 

We note in passing that the value of ϵs is dependent on how scaling is constructed. We expect different ϵs and scaling behavior when a scaling line is constructed using just a single metal with varied facets [Vak2 is a constant for a fixed metal, rather than =f(ϵd), see Sec. S9 for further discussion]. Since scaling lines are, in general, constructed across metals, we emphasize this type of scaling in this work.

Our analysis of C* and O* scaling, in this section, produces two main insights for the use of scaling lines in constructing activity volcanoes.

  1. When constructing scaling lines for elements of the same row, metals with large, negative ϵd values, such as noble metals may be left out. The reason is that they are likely to fall off linear scaling lines, as their energies are largely governed by orthogonalization and not hybridization. A more rigorous analysis may include using ϵs to gauge which metals linear scaling should be fit to.

  2. Linear scaling is likely to be much improved when all metals belong to the same row of transition metals. The suitability of using multiple rows can be judged based on the Eortho scaling between the two adsorbates.

Note that, in practice, there could be several other reasons for lack of scaling between two adsorbates. Relaxation of the surface and site dependence of the adsorbate might cause scatter on simple transition metal surfaces. For that reason, one should always consider scaling for a fixed surface structure. Adsorbates other than those that are monoatomic could have alternative effects that predominate due to greater degrees of freedom. However, the adsorbates and surfaces we describe in this work are transition metals with the simplest of mono-atomic adsorbates, with chemisorption energies in the dilute coverage limit. The fact that these simple systems do not scale in an obvious manner shows the complexity associated with surface chemical bond formation. Care must be taken while fitting linear scaling relations for transition metals for more complex reactions and intermediates.

In this work, we develop a model to understand why chemisorption energies of O* and C* do not follow a linear scaling relation. Following the d-band model of chemisorption, our model combines a modified version of the Newns–Anderson hybridization energy with effective terms for the orthogonalization repulsion. We show that lack of linear scaling along a row of the Periodic Table is a consequence of the saturation of the hybridization energy of O* for metals with lower lying d-band centers, while the hybridization energy for C* shows little to no saturation at relevant ranges of ϵd. We suggest a measure, ϵs, for the highest ϵd value at which saturation is expected. We quantify the row-dependence of scaling by determining the variation of the orthogonalization energy with the d-band center and show that adsorbates with similar ϵa are likely across rows. Finally, we show that linear scaling relations can be recovered in the case of C* and O* by considering metals of the same row of the transition metal series, bar noble metals.

Density functional theory calculations were carried out using Quantum ESPRESSO.40 Core electrons were treated using pseudopotentials from the standard solid-state pseudopotentials (SSSP) Precision pseudopotential database.41 Valence electrons were described using plane-waves with kinetic energy and density cutoffs of at least 80 and 600 Ry, respectively. In cases where the recommended cutoffs are greater that these values, the larger of the two is used. Gaussian smearing of width 0.0075 Ry (≈0.1 eV) was used for all calculations. A k-point mesh of (4, 4, 1) was used for all surface calculations. The PBE42 functional was used to describe exchange and correlation effects. Calculations were considered converged if the maximum force on each atom was less than 10−4 Ry bohr−1 and the energy difference between successive SCF cycles is less than 2 · 10−5 Ry. All calculations were performed without spin polarization. Structures were prepared using the atomic simulation environment.43 The adsorbate was initially placed at symmetrically inequivalent sites at the start of the relaxation. It is found that adsorbates at the hollow site for all metals have the lowest energy. Provenance of each calculation in this work was stored using AiiDA44 and performed using the AiiDA quantumespresso plugin.

The Newns–Anderson expressions were implemented using trigonometric functions from the Python numpy45 library. Multi-precision C-library arb46 was used to integrate the projected density of states using the acb_calc_integrate function. Decimal precision of 50 was used throughout this work. Least-squares error fitting were performed using the Python scipy47 library.

See the supplementary material for the list of symbols, model development, model parameterization details, and extended discussion on parameterization of the Newns–Anderson–Grimley model.

This research received funding from the European Union’s Horizon 2020 Research and Innovation Program under Grant Agreement No. 851441, SELECTCO2, and from The VILLUM Center for the Science of Sustainable Fuels and Chemicals (9455) from VILLUM FONDEN. Sudarshan Vijay and Karen Chan acknowledge computational resources from PRACE (Project No. 2020235596) and the Juelich Supercomputing Center.

The authors have no conflicts to disclose.

Python codes used to generate all figures in this manuscript are available at https://github.com/sudarshanv01/newns-anderson. The Newns–Anderson equations and those used in the text for the chemisorption energies and their derivatives are implemented in Python and can be accessed through the following package: https://github.com/sudarshanv01/CatChemi. AiiDA archive files are available at 10.24435/materialscloud:m2-fp.

1.
J. K.
Nørskov
,
T.
Bligaard
,
J.
Rossmeisl
, and
C. H.
Christensen
, “
Towards the computational design of solid catalysts
,”
Nat. Chem.
1
,
37
46
(
2009
).
2.
A.
Vojvodic
,
J. K.
Nørskov
, and
F.
Abild-Pedersen
, “
Electronic structure effects in transition metal surface chemistry
,”
Top. Catal.
57
,
25
32
(
2014
).
3.
C. T.
Campbell
and
J. R. V.
Sellers
, “
Enthalpies and entropies of adsorption on well-defined oxide surfaces: Experimental measurements
,”
Chem. Rev.
113
,
4106
4135
(
2013
).
4.
C. T.
Campbell
, “
The degree of rate control: A powerful tool for catalysis research
,”
ACS Catal.
7
,
2770
2779
(
2017
).
5.
B.
Hammer
,
Y.
Morikawa
, and
J. K.
Nørskov
, “
CO chemisorption at metal surfaces and overlayers
,”
Phys. Rev. Lett.
76
,
2141
2144
(
1996
).
6.
B.
Hammer
,
L. B.
Hansen
, and
J. K.
Nørskov
, “
Improved adsorption energetics within density-functional theory using revised Perdew–Burke–Ernzerhof functionals
,”
Phys. Rev. B: Condens. Matter Mater. Phys.
59
,
7413
7421
(
1999
).
7.
L. G. M.
Pettersson
and
A.
Nilsson
, “
A molecular perspective on the d-band model: Synergy between experiment and theory
,”
Top. Catal.
57
,
2
13
(
2014
).
8.
S.
Bhattacharjee
,
U. V.
Waghmare
, and
S.-C.
Lee
, “
An improved d-band model of the catalytic activity of magnetic transition metal surfaces
,”
Sci. Rep.
6
,
35916
(
2016
); arXiv:1610.01746.
9.
M.
Andersen
,
S. V.
Levchenko
,
M.
Scheffler
, and
K.
Reuter
, “
Beyond scaling relations for the description of catalytic materials
,”
ACS Catal.
9
,
2752
2759
(
2019
); arXiv:1902.07495.
10.
M. T.
Greiner
,
T. E.
Jones
,
S.
Beeg
,
L.
Zwiener
,
M.
Scherzer
,
F.
Girgsdies
,
S.
Piccinin
,
M.
Armbrüster
,
A.
Knop-Gericke
, and
R.
Schlögl
, “
Free-atom-like d states in single-atom alloy catalysts
,”
Nat. Chem.
10
,
1008
1015
(
2018
).
11.
A. J.
Medford
,
A.
Vojvodic
,
J. S.
Hummelshøj
,
J.
Voss
,
F.
Abild-Pedersen
,
F.
Studt
,
T.
Bligaard
,
A.
Nilsson
, and
J. K.
Nørskov
, “
From the Sabatier principle to a predictive theory of transition-metal heterogeneous catalysis
,”
J. Catal.
328
,
36
42
(
2015
).
12.
A. J.
Medford
,
J.
Wellendorff
,
A.
Vojvodic
,
F.
Studt
,
F.
Abild-Pedersen
,
K. W.
Jacobsen
,
T.
Bligaard
, and
J. K.
Nørskov
, “
Assessing the reliability of calculated catalytic ammonia synthesis rates
,”
Science
345
,
197
200
(
2014
).
13.
A. J.
Medford
,
C.
Shi
,
M. J.
Hoffmann
,
A. C.
Lausche
,
S. R.
Fitzgibbon
,
T.
Bligaard
, and
J. K.
Nørskov
, “
CatMAP: A software package for descriptor-based microkinetic mapping of catalytic trends
,”
Catal. Lett.
145
,
794
807
(
2015
).
14.
Z. W.
Seh
,
J.
Kibsgaard
,
C. F.
Dickens
,
I.
Chorkendorff
,
J. K.
Nørskov
, and
T. F.
Jaramillo
, “
Combining theory and experiment in electrocatalysis: Insights into materials design
,”
Science
,
355
379
380
(
2017
).
15.
P.
Sabatier
,
La Catalyse en Chimie Organique
(
Librairie Polytechnique
,
1920
).
16.
F.
Abild-Pedersen
,
J.
Greeley
,
F.
Studt
,
J.
Rossmeisl
,
T. R.
Munter
,
P. G.
Moses
,
E.
Skúlason
,
T.
Bligaard
, and
J. K.
Nørskov
, “
Scaling properties of adsorption energies for hydrogen-containing molecules on transition-metal surfaces
,”
Phys. Rev. Lett.
99
,
016105
(
2007
).
17.
F.
Studt
,
F.
Abild-Pedersen
,
T.
Bligaard
,
R. Z.
Sørensen
,
C. H.
Christensen
, and
J. K.
Nørskov
, “
Identification of non-precious metal alloy catalysts for selective hydrogenation of acetylene
,”
Science
320
,
1320
1322
(
2008
).
18.
D. M.
Newns
, “
Self-consistent model of hydrogen chemisorption
,”
Phys. Rev.
178
,
1123
1135
(
1969
).
19.
I. E. L.
Stephens
,
A. S.
Bondarenko
,
U.
Grønbjerg
,
J.
Rossmeisl
, and
I.
Chorkendorff
, “
Understanding the electrocatalysis of oxygen reduction on platinum and its alloys
,”
Energy Environ. Sci.
5
,
6744
(
2012
).
20.
J. K.
Nørskov
,
J.
Rossmeisl
,
A.
Logadottir
,
L.
Lindqvist
,
J. R.
Kitchin
,
T.
Bligaard
, and
H.
Jónsson
, “
Origin of the overpotential for oxygen reduction at a fuel-cell cathode
,”
J. Phys. Chem. B
108
,
17886
17892
(
2004
).
21.
R. D.
Cortright
,
J. M.
Hill
, and
J. A.
Dumesic
, “
Selective dehydrogenation of isobutane over supported Pt/Sn catalysts
,”
Catal. Today
55
,
213
223
(
2000
).
22.
S.
Bhandari
,
S.
Rangarajan
,
C. T.
Maravelias
,
J. A.
Dumesic
, and
M.
Mavrikakis
, “
Reaction mechanism of vapor-phase formic acid decomposition over platinum catalysts: DFT, reaction kinetics experiments, and microkinetic modeling
,”
ACS Catal.
10
,
4112
4126
(
2020
).
23.
A. D.
Allian
,
K.
Takanabe
,
K. L.
Fujdala
,
X.
Hao
,
T. J.
Truex
,
J.
Cai
,
C.
Buda
,
M.
Neurock
, and
E.
Iglesia
, “
Chemisorption of CO and mechanism of CO oxidation on supported platinum nanoclusters
,”
J. Am. Chem. Soc.
133
,
4498
4517
(
2011
).
24.
C.
Pedrero
,
T.
Waku
, and
E.
Iglesia
, “
Oxidation of CO in H2-CO mixtures catalyzed by platinum: Alkali effects on rates and selectivity
,”
J. Catal.
233
,
242
255
(
2005
).
25.
E.
Becker
,
P.
Carlsson
,
H.
Grönbeck
, and
M.
Skoglundh
, “
Methane oxidation over alumina supported platinum investigated by time-resolved in situ XANES spectroscopy
,”
J. Catal.
252
,
11
17
(
2007
).
26.
K.
Fujimoto
,
M.
Kameyama
, and
T.
Kunugi
, “
Hydrogenation of adsorbed carbon monoxide on supported platinum group metals
,”
J. Catal.
61
,
7
14
(
1980
).
27.
G. A.
Somorjai
, “
The surface structure of and catalysis by platinum single crystal surfaces
,”
Catal. Rev.
7
,
87
120
(
1972
).
28.
W.
Rachmady
and
M.
Vannice
, “
Acetic acid hydrogenation over supported platinum catalysts
,”
J. Catal.
192
,
322
334
(
2000
).
29.
S.
Wang
,
H. S.
Pillai
, and
H.
Xin
, “
Bayesian learning of chemisorption for bridging the complexity of electronic descriptors
,”
Nat. Commun.
11
,
6132
(
2020
).
30.
B.
Hammer
and
J. K.
Nørskov
, “
Theory of adsorption and surface reactions
,” in
Chemisorption and Reactivity on Supported Clusters and Thin Films
(
Springer
,
1997
).
31.
A.
Nilsson
and
L. G. M.
Pettersson
, “
Adsorbate electronic structure and bonding on metal surfaces
,” in
Chemical Bonding at Surfaces and Interfaces
(
Elsevier
,
2008
), Technical Report.
32.
J. K.
Nørskov
,
F.
Studt
,
F.
Abild-Pedersen
, and
T.
Bligaard
,
Fundamental Concepts in Heterogeneous Catalysis
(
Wiley Blackwell
,
2014
), ISBN: 978-1-118-88895-7, pp.
1
196
.
33.
A.
Ruban
,
B.
Hammer
,
P.
Stoltze
,
H. L.
Skriver
, and
J. K.
Nørskov
, “
Surface electronic structure and reactivity of transition and noble metals
,”
J. Mol. Catal. A: Chem.
,
115
,
421
429
(
1997
).
34.
J. K.
Nørskov
, “
Effective medium potentials for molecule-surface interactions: H2 on Cu and Ni surfaces
,”
J. Chem. Phys.
90
,
7461
7471
(
1989
).
35.
O.
Andersen
,
O.
Jepsen
, and
D.
Glötzel
, “
Canonical description of the band structures of metals
,” in
Proc. Int. Sch. Physics, Course LXXXIX
(
North-Holland
,
1985
),
59
.
36.
T. B.
Grimley
, “
Overlap effects in the theory of adsorption using Anderson’s Hamiltonian
,”
J. Phys. C: Solid State Phys.
3
,
1934
1942
(
1970
).
37.
T. D.
Spivey
and
A.
Holewinski
, “
Selective interactions between free-atom-like d-states in single-atom alloy catalysts and near-frontier molecular orbitals
,”
J. Am. Chem. Soc.
143
,
11897
11902
(
2021
).
38.
H.
Xin
and
S.
Linic
, “
Communications: Exceptions to the d-band model of chemisorption on metal surfaces: The dominant role of repulsion between adsorbate states and metal d-states
,”
J. Chem. Phys.
132
,
221101
(
2010
).
39.
H. A.
Hansen
,
V.
Viswanathan
, and
J. K.
Nørskov
, “
Unifying kinetic and thermodynamic analysis of 2 e- and 4 e-reduction of oxygen on metal surfaces
,”
J. Phys. Chem. C
118
,
6706
6718
(
2014
).
40.
P.
Giannozzi
,
S.
Baroni
,
N.
Bonini
,
M.
Calandra
,
R.
Car
,
C.
Cavazzoni
,
D.
Ceresoli
,
G. L.
Chiarotti
,
M.
Cococcioni
,
I.
Dabo
,
A.
Dal Corso
,
S.
De Gironcoli
,
S.
Fabris
,
G.
Fratesi
,
R.
Gebauer
,
U.
Gerstmann
,
C.
Gougoussis
,
A.
Kokalj
,
M.
Lazzeri
,
L.
Martin-Samos
,
N.
Marzari
,
F.
Mauri
,
R.
Mazzarello
,
S.
Paolini
,
A.
Pasquarello
,
L.
Paulatto
,
C.
Sbraccia
,
S.
Scandolo
,
G.
Sclauzero
,
A. P.
Seitsonen
,
A.
Smogunov
,
P.
Umari
, and
R. M.
Wentzcovitch
, “
QUANTUM ESPRESSO: A modular and open-source software project for quantum simulations of materials
,”
J. Phys.: Condens. Matter
21
,
395502
(
2009
); arXiv:0906.2569.
41.
G.
Prandini
,
A.
Marrazzo
,
I. E.
Castelli
,
N.
Mounet
, and
N.
Marzari
, “
Precision and efficiency in solid-state pseudopotential calculations
,”
npj Comput. Mater.
4
,
72
(
2018
); arXiv:1806.05609.
42.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
, “
Generalized gradient approximation made simple
,”
Phys. Rev. Lett.
77
,
3865
3868
(
1996
).
43.
A.
Hjorth Larsen
,
J.
Jørgen Mortensen
,
J.
Blomqvist
,
I. E.
Castelli
,
R.
Christensen
,
M.
Dułak
,
J.
Friis
,
M. N.
Groves
,
B.
Hammer
,
C.
Hargus
,
E. D.
Hermes
,
P. C.
Jennings
,
P.
Bjerre Jensen
,
J.
Kermode
,
J. R.
Kitchin
,
E.
Leonhard Kolsbjerg
,
J.
Kubal
,
K.
Kaasbjerg
,
S.
Lysgaard
,
J.
Bergmann Maronsson
,
T.
Maxson
,
T.
Olsen
,
L.
Pastewka
,
A.
Peterson
,
C.
Rostgaard
,
J.
Schiøtz
,
O.
Schütt
,
M.
Strange
,
K. S.
Thygesen
,
T.
Vegge
,
L.
Vilhelmsen
,
M.
Walter
,
Z.
Zeng
, and
K. W.
Jacobsen
, “
The atomic simulation environment—A Python library for working with atoms
,”
J. Phys.: Condens. Matter
29
,
273002
(
2017
).
44.
S. P.
Huber
,
S.
Zoupanos
,
M.
Uhrin
,
L.
Talirz
,
L.
Kahle
,
R.
Häuselmann
,
D.
Gresch
,
T.
Müller
,
A. V.
Yakutovich
,
C. W.
Andersen
,
F. F.
Ramirez
,
C. S.
Adorf
,
F.
Gargiulo
,
S.
Kumbhar
,
E.
Passaro
,
C.
Johnston
,
A.
Merkys
,
A.
Cepellotti
,
N.
Mounet
,
N.
Marzari
,
B.
Kozinsky
, and
G.
Pizzi
, “
AiiDA 1.0, a scalable computational infrastructure for automated reproducible workflows and data provenance
,”
Sci. Data
7
,
300
(
2020
); arXiv:2003.12476.
45.
S.
Van Der Walt
,
S. C.
Colbert
, and
G.
Varoquaux
, “
The NumPy array: A structure for efficient numerical computation
,”
Comput. Sci. Eng.
13
,
22
30
(
2011
); arXiv:1102.1523.
46.
F.
Johansson
, “
Arb: Efficient arbitrary-precision midpoint-radius interval arithmetic
,”
IEEE Trans. Comput.
66
,
1281
1292
(
2017
).
47.
P.
Virtanen
,
R.
Gommers
,
T. E.
Oliphant
,
M.
Haberland
,
T.
Reddy
,
D.
Cournapeau
,
E.
Burovski
,
P.
Peterson
,
W.
Weckesser
,
J.
Bright
,
S. J.
van der Walt
,
M.
Brett
,
J.
Wilson
,
K. J.
Millman
,
N.
Mayorov
,
A. R. J.
Nelson
,
E.
Jones
,
R.
Kern
,
E.
Larson
,
C. J.
Carey
,
İ.
Polat
,
Y.
Feng
,
E. W.
Moore
,
J.
VanderPlas
,
D.
Laxalde
,
J.
Perktold
,
R.
Cimrman
,
I.
Henriksen
,
E. A.
Quintero
,
C. R.
Harris
,
A. M.
Archibald
,
A. H.
Ribeiro
,
F.
Pedregosa
,
P.
van Mulbregt
,
A.
Vijaykumar
,
A. P.
Bardelli
,
A.
Rothberg
,
A.
Hilboll
,
A.
Kloeckner
,
A.
Scopatz
,
A.
Lee
,
A.
Rokem
,
C. N.
Woods
,
C.
Fulton
,
C.
Masson
,
C.
Häggström
,
C.
Fitzgerald
,
D. A.
Nicholson
,
D. R.
Hagen
,
D. V.
Pasechnik
,
E.
Olivetti
,
E.
Martin
,
E.
Wieser
,
F.
Silva
,
F.
Lenders
,
F.
Wilhelm
,
G.
Young
,
G. A.
Price
,
G.-L.
Ingold
,
G. E.
Allen
,
G. R.
Lee
,
H.
Audren
,
I.
Probst
,
J. P.
Dietrich
,
J.
Silterra
,
J. T.
Webber
,
J.
Slavič
,
J.
Nothman
,
J.
Buchner
,
J.
Kulick
,
J. L.
Schönberger
,
J. V.
de Miranda Cardoso
,
J.
Reimer
,
J.
Harrington
,
J. L. C.
Rodríguez
,
J.
Nunez-Iglesias
,
J.
Kuczynski
,
K.
Tritz
,
M.
Thoma
,
M.
Newville
,
M.
Kümmerer
,
M.
Bolingbroke
,
M.
Tartre
,
M.
Pak
,
N. J.
Smith
,
N.
Nowaczyk
,
N.
Shebanov
,
O.
Pavlyk
,
P. A.
Brodtkorb
,
P.
Lee
,
R. T.
McGibbon
,
R.
Feldbauer
,
S.
Lewis
,
S.
Tygier
,
S.
Sievert
,
S.
Vigna
,
S.
Peterson
,
S.
More
,
T.
Pudlik
,
T.
Oshima
,
T. J.
Pingel
,
T. P.
Robitaille
,
T.
Spura
,
T. R.
Jones
,
T.
Cera
,
T.
Leslie
,
T.
Zito
,
T.
Krauss
,
U.
Upadhyay
,
Y. O.
Halchenko
, and
Y.
Vázquez-Baeza
, “
SciPy 1.0: Fundamental algorithms for scientific computing in Python
,”
Nat. Methods
17
,
261
272
(
2020
); arXiv:1907.10121.

Supplementary Material