Yongjie Xia,
Liang Huangab,
Robert C. Forreyc and
Hansong Cheng*ab
aDepartment of Chemistry, National University of Singapore, 3 Science Drive 3, Singapore 117543
bSustainable Energy Laboratory, Faculty of Material Science and Chemistry, China University of Geosciences, Wuhan 430074, China. E-mail: chmch@nus.edu.sg
cDepartment of Physics, Penn State University, Berks Campus, USA
First published on 18th August 2014
Tungsten carbide has been proposed to be a potential inexpensive alternative to platinum for catalyzing hydrogenation reactions. Using density functional theory (DFT) calculations, we performed a systematic study on the geometries and electronic structures of small (WC)n (n = 1–10) clusters. A cubic-like growth of WC clusters was found to be the preferred pathway. The interactions between the WC clusters and hydrogen were compared to interactions of hydrogen with the WC(0001) surface. Both atomic clusters and the crystalline surface of WC can lead to facile dissociation of H2 while the activation energy on the former is even lower; the H diffusion barrier on WC clusters was found to be significantly higher than that on the WC(0001) surface, yet both of the calculated diffusion barriers are well above the values on the platinum counterpart. The lack of hydrogen mobility may exert a profound influence on the kinetics of hydrogen for participating in reactions on WC. The desorption strength of H on (WC)n (n = 2, 4, 6, 9) at various hydrogen coverages was explored and the hydrogen saturated (WC)9 cluster was used to rationalize the catalytic formation of hydrogen tungsten bronzes. The study provides insights which may be useful for efforts to replace Pt with WC in catalytic processes involving hydrogen.
Atomic clusters are amenable to theoretical treatment and have been used as a model for studying how material properties scale with an increasing size.11 Bare and supported Pt-group metal clusters have been systematically investigated to understand chemical reactivity with hydrogen.12–14 However, as a perceived “Pt-like” catalyst, the chemical properties of tungsten carbide clusters have rarely been addressed. To our knowledge, only (WC)n (n = 1, 2, 4, 6, 12, 15, 36) clusters have been constructed to understand the mechanical strength of tungsten carbide.15 In the present study, we conducted a systematic study using DFT on hydride formation (the hydride nature of (WC)nHx is verified by charge analysis, vide infra) of small (WC)n (n = 1–10) clusters with the objective to gain insight into the similarities and differences compared to Pt in catalyzing hydrogenation. The dissociative chemisorption of H2 on selected WC clusters and the diffusion of H atoms upon adsorption were examined and compared with the reaction on a WC surface described by a periodic slab model. In particular, key catalytic properties of (WC)n (n = 2, 4, 6, 9) including H2 dissociative chemisorption energies and the maximum capacity for accommodating H atoms at full hydrogen saturation are addressed with the purpose of providing chemical insight into realistic catalytic processes in which catalyst surfaces are fully covered by hydrogen at a finite pressure.16 This work is organized as follows: the computational methodology is briefly discussed in section 2; the geometries of (WC)n are presented in section 3.1, the dissociative adsorption and diffusion of hydrogen on a bare (WC)6 cluster is considered in section 3.2, the adsorption strengths of H on selected WC clusters with increasing hydrogen loading are presented in section 3.3. The adsorption, dissociation, and diffusion of hydrogen on the WC(0001) surface are discussed in section 3.4. The thermodynamics of H2 chemisorption on WC(0001) was investigated in 3.5. Finally, we summarized this work in section 4.
The chemical reactivity of hydrogen on the WC(0001) surface was computed using the periodic DFT method as implemented in the Vienna ab initio simulation package (VASP).24,25 Bulk WC adopts either a hexagonal α-WC lattice or a cubic β-WC lattice.15 The former is more stable below 2525 °C with the W-terminated (0001) surface being most stable, which is therefore selected in the present study.15 A (3 × 3) cell of the WC(0001) surface was modeled as periodically repeating slabs consisting of three WC bilayers. A 15 Å vacuum layer was inserted between adjacent slabs (Fig. 1a and b). The projector-augmented wave method and the PBE functional were used with an energy cutoff of 400 eV. The Brillouin zone integration was sampled within a 3 × 3 × 1 Monkhorst–Pack k-point mesh.26 Full structural optimizations with a conjugate-gradient method were performed until the forces became smaller than 0.05 eV Å−1. We used the climbing-image nudged elastic band (CI-NEB)27 method to calculate the minimum energy profiles along the prescribed reaction pathways.
To verify the stability of (WC)n clusters with full hydrogen saturation, we performed ab initio molecular dynamics (MD) simulations at 298 K for 6 ps with a time step of 1 fs in a NVE canonical ensemble. The MD simulations were conducted with VASP using a 20 × 20 × 20 Å3 box and Γ point for Brillouin zone sampling. The results of selected MD simulations are included in the ESI.†
![]() | ||
Fig. 2 The optimized structures of (WC)n (n = 1–10). For WC and (WC)2, the most stable structures are shown; for (WC)n (n = 3–10), the three lowest energy structures are illustrated, denoted as n–i for (WC)n with i = 1–3. Symmetry and relative energy (in eV) with respect to the most stable isomer are also provided. Gray sphere C, blue sphere W. Cartesian coordinates of the most stable configurations are provided in ESI.† |
The lowest energy structure of (WC)2 was found to be planar while the trimer (3-1) is a 3D structure. The tetramer (4-1) adopts a distorted cubic structure, which can be viewed as a building block of larger WC clusters. Indeed, the lowest energy structures of larger clusters, i.e. 6-1, 8-1 and 9-1, as well as the less stable isomers including 8-2, 10-2 and 10-3 all exhibit an irregular cubic pattern, which is distinctively different from the hexagonal pattern in bulk tungsten carbide. The cubic feature can also be discerned in less regular clusters including 5-1, 7-1 and 10-1. Our results suggest that small (WC)n clusters adopt a cubic structural growth pathway, which is consistent with the previous study.15
The calculated cluster binding energies of the minimum energy structures are presented in Fig. 3. As expected, the binding energy increases monotonically with the cluster size. We note that the calculated binding energy of (WC)10 is still well below the calculated bulk cohesive energy of WC (16.05 eV),28 and thus the properties of these small WC clusters may be distinctive from bulk WC.
We first investigate H2 dissociative adsorption and H diffusion at zero coverage. It has been shown that H2 can be adsorbed on a solid surface via quasi-molecular bonding by lengthening the H–H bond to 0.9 Å.31 By sampling several initial configurations for H2 to approach the (WC)6 cluster, we identified the adsorption mode with H2 on a vertex site Wver being energetically most stable (Fig. 4). The H–H bond is elongated to 0.923 Å (R1). The dissociation process is essentially barrierless (0.02 eV), similar to the case of H2 dissociation on Pt clusters for which a very low barrier (0.05 eV) was reported.16
![]() | ||
Fig. 4 The calculated energy diagram of H2 dissociative chemisorption on (WC)6 and the subsequent migrations of the H atoms on the cluster. The energies are relative to the sum of the energy of bare (WC)6 and H2. The hydrogen adsorption sites are also labeled, with Cedge, Wedge denoting edge sites and Cver, Wver denoting vertex sites respectively; H–H distances are labeled in parentheses, H–W (H–C) bond lengths are labeled in italic (in Å). Numbers of the imaginary frequencies for each transition state are given in italic. Gray C, blue W, green H (the same legend is used in Fig. 6, 7 and 9). Cartesian coordinates of all states are provided in ESI.† |
The quasi-dissociative state of the H2 molecule on the (WC)6 cluster is attributable to the favorable interaction between the lowest unoccupied molecular orbital (LUMO) of H2 and the highest occupied molecular orbital (HOMO) of the (WC)6 cluster, as illustrated in Fig. 5. The total energy is lowered to −1.11 eV, accompanied by lengthening of the H–H distance to 1.990 Å (P1). The H atoms can diffuse further to form the energetically most stable structure (P5). The overall process is highly exothermic and the total energy is lowered by 2.46 eV (R0 → P5). To examine the mobility of H on (WC)6, two H diffusion routes(P1 → P2 → P5 and P1 → P3 → P4 → P5) were sampled. As displayed in Fig. 4, the former route where an H atom diffuses from Wver to Wedge and further to another Wver has an activation barrier of 0.73 eV. The activation barrier for the later route is much higher and thus less favored. The barrier of 0.73 eV is still significantly higher than the H diffusion barrier on Pt (0.23 eV), Pd (0.11 eV) and Ni (0.19 eV) clusters.13 The results suggest that H diffusion on the selected WC cluster at zero coverage is moderately difficult, in contrast to the facile diffusion of H atoms on Pt clusters. The significantly lower H mobility is expected to substantially diminish the catalytic activity of the material.
(1). Each W site is capable of accommodating 1–3 H atoms, depending on the coordination environment. The W atom of the planar (WC)2 cluster can hold 3 H atoms while larger WC clusters with 3D structures can coordinate with two H atoms on a vertex site and one H atom on an edge site or a face site (the vertex site and edge site are depicted in Fig. 4 in (WC)6 and the face site is shown in Fig. 7a in the (WC)9 hydride).
![]() | ||
Fig. 7 (a) The fully saturated hydrides of the selected (WC)n (n = 2, 4, 6, 9) clusters, the face adsorption sites of (WC)9 are also presented. Cartesian coordinates of the four configurations are available in ESI.† (b) The calculated H–W bond length distributions of the (WC)n hydrides. g(r) was obtained by tabulating all the H–W distances at each step of the MD trajectories fitted with Gaussian functions. |
(2). Each C atom can accommodate only one H atom.
(3). H atoms generally prefer W atoms first; the third H atom on W for the case of (WC)2 is less preferred than on C; the preference of H on different W sites follows Wver > Wedge > Wface, the second H on Wver site has slightly lower adsorption strength than on C sites (detailed H adsorption strengths on various W and C sites of (WC)6 and (WC)9 are compared in Fig. S1†).
Structures with H atoms disobeying the rules usually display lower H2 adsorption strength, as exemplified with (WC)6 in Fig. S2.† The (WC)2, (WC)4, (WC)6 and (WC)9 clusters with full H-loading are displayed in Fig. 7a. The clusters can chemisorb up to 4, 6, 8 and 11 H2 molecules, respectively. Room temperature MD simulations were performed to ensure that no H–H recombination occurs, as illustrated in Fig. S3.† To confirm that a cluster is at full H loading, additional H atoms to those “fully saturated” WC hydrides would result in formation of a H2 molecule. The average H2 adsorption energies and sequential H desorption energies of different H coverage are defined as:
Eads = (2E(WC)n + mEH2 − 2E(WC)nHm)/m, | (1) |
Edes = EH − (E(WC)nHm − E(WC)nHm−2)/2, (n = 2, 4, 6, 9 and m = 2, 4, 6…), | (2) |
![]() | ||
Fig. 8 (a) The calculated average H2 adsorption energies; (b) the H desorption energies and (c) the average Hirshfeld charge gain per H atom for (WC)n (n = 2, 4, 6, 9) clusters vs. H coverage. |
We observe certain irregular behavior ((WC)9H4) in the Eads plot, which is largely due to the structure relaxation induced by hydrogen adsorption. However, the overall trend of Eads in the four cases decreases as H loading increases, as expected. In view of the fact that W atoms can accommodate more H than C at full H-loading and is energetically more favorable, the W–H bond is considered to be the main contribution to the adsorption strength. Therefore, the W–H bond length distribution g(r) is used to characterize H adsorption (Fig. 7b). The position of the peak of g(r) increases with cluster size. The bond length rationale is consistent with the calculated H2 dissociative chemisorption energies (Fig. 8a) which shows that hydrogen adsorption strength on (WC)n at full H-loading decreases with increasing cluster size within a narrow energy range (1.0–1.3 eV), similar to the narrow range of hydrogen adsorption strength for Ptn (0.9–1.1 eV), Pdn (0.6–0.8 eV) and Nin (0.8–1.1 eV).13 In parallel, the desorption energy of an H atom from fully saturated (WC)n (n = 2, 4, 6, 9) (Fig. 8b) also varies within a narrow range (2.15–2.49 eV). Again, the energy is in a similar range from 2.08 eV to 2.73 eV as reported previously for hydrogen adsorption on small Pt, Pd and Ni clusters.13 With an increase of H-loading, the bonding nature of the metal hydrides evolves from metallic to covalent, leading to localized chemical properties and the narrowness of chemisorption and the desorption energies at full H-loading.13 The capacity of accommodating hydrogen atoms for WC is smaller than Pt of the same size. For example, Pt9 can accommodate 34 H atoms, as compared with 22 H for (WC)9. In comparison, Pd clusters have comparable hydrogen uptake to that of WC clusters (22 H for Pd9), while the H capacity of Ni was found to be slightly lower than the value of both Pd and WC clusters (18 H for Ni9). We note here that the topological geometries of (WC)n clusters and their hydrides differ significantly from their counter-parts of nickel family clusters and the corresponding hydrides since the “super-atom” structure is not spherical.
As established in the previous section, the dissociative chemisorption of hydrogen on WC is facile. Another prerequisite for the WC catalyzed formation of HxWO3 is the facile H migration from WC to a substrate. Here, we consider the WO3(001) surface as the substrate. Thermodynamically, for this process to occur, the H adsorption strength on the WO3 surface must be stronger than the H desorption strength on WC. Indeed, the calculated Edes values for (WC)6 and (WC)9 at full H-loading are smaller than the reported H adsorption strength on the WO3(001) surface (2.84 eV),32 making the H migration process thermodynamically feasible. The slow kinetics of the WC catalyzed spillover might be attributable to the low hydrogen mobility (see Fig. 4). Obviously, if a (WC)n cluster is not fully saturated by H atoms, H-migration to the substrate surface may not occur. For example, the H desorption energy of (WC)9H14 was calculated to be 2.97 eV, which is larger than the H adsorption strength on the WO3(001) surface. For most practical hydrogenation processes, the hydrogen pressure would be large enough so that the WC cluster would be fully saturated. The “full-saturation” scenario may also be applicable to model other reactions where WC serves as a catalyst. The hydride nature of (WC)nHx is reflected by charge analysis (Fig. 8c), which shows that H is negatively charged at various hydrogen coverage. This is similar to small Pt, Pd and Ni clusters, which also form hydride upon hydrogen adsorption, indicating that the hydrogen adsorbed on WC would be nucleophilic. The average partial charge gain per H atom essentially decreases with increased H loading in all cases, consistent with the corresponding trends of Eads. Again, the discrepancy of a few points is attributed to the structure relaxation induced by hydrogen adsorption.
A previous study proposed a volcano relationship to correlate the log of exchange current density with hydrogen binding energy on surfaces of an electrode to explain the HER activity of various cathode materials.6 The d-band center theory has also been a widely-used indicator describing activity or selectivity of a catalyst.33 However, a correlation between kinetics of a catalytic hydrogenation reaction and H diffusion barrier has not been well understood. Empirically, the low H mobility is chiefly responsible for the mediocre performance of the WC catalyst, which suggests that H diffusion barriers should be an important parameter to be considered for design of new catalysts.
ΔG = G(2H/WC(0001)) − G(WC(0001)) − G(H2) | (3) |
ΔG = ΔH − T[S(2H/WC(0001)) − S(WC(0001)) − S(H2)] | (4) |
We then approximate the reaction enthalpy with the calculated reaction energies ΔE since the reaction enthalpy does not change with temperature significantly. The entropic change from WC(0001) and 2H/WC(0001) is much smaller than the gas phase entropy of molecular H2 and thus is neglected. Then eqn (3) becomes
ΔG ≈ ΔE + TS(H2) | (5) |
![]() | (6) |
The equilibrium pressure of hydrogen at a given temperature is plotted in Fig. 10. Remarkably, the equilibrium temperature for H2 at a partial pressure of 1 atm is ca. 1010 K, which suggests that WC can retain its catalytic activity up to high temperature.
Hydrogen adsorption strength of (WC)2,4,6,9 at various hydrogen coverage were obtained. Both the H2 adsorption energies (Eads) and H desorption energies (Edes) of fully saturated (WC)nHx clusters lie within a narrow range. The Edes of (WC)2,4,6,9 are comparable to the corresponding Ptn of the same size. We use the Edes of (WC)9 at full H-loading to account for the H migration from tungsten carbide to tungsten oxide, for which the full saturation model is highlighted. Thermodynamics analysis reveals that tungsten carbide can retain its catalytic activity at relatively high temperature.
The elucidated inert diffusion behavior of H on tungsten carbide appears to be the primary reason for the relatively low catalytic activity of tungsten carbide. We expect this study can provide useful insight into the catalytic behavior of tungsten carbide and shed a light for design of efficient, low cost hydrogenation catalysts.
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/c4ra06225f |
This journal is © The Royal Society of Chemistry 2014 |