Reaction Mechanism with Thermodynamic Structural Screening for

Jul 8, 2018 - Reaction Mechanism with Thermodynamic Structural Screening for Electrochemical Hydrogen Evolution on Monolayer 1T' phase MoS2...
0 downloads 0 Views 4MB Size
Article Cite This: Chem. Mater. 2018, 30, 5404−5411

pubs.acs.org/cm

Reaction Mechanism with Thermodynamic Structural Screening for Electrochemical Hydrogen Evolution on Monolayer 1T′ Phase MoS2 Shiqi Chen,† Xiaobo Chen,*,† Guangjin Wang,‡,§ Lu Liu,† Qiaoqiao He,† Xi-Bo Li,† and Ni Cui† †

Downloaded via KAOHSIUNG MEDICAL UNIV on October 8, 2018 at 20:27:53 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

Siyuan Laboratory, Guangzhou Key Laboratory of Vacuum Coating Technologies and New Energy Materials, Guangdong Provincial Engineering Technology Research Center of Vacuum Coating Technologies and New Energy Materials, Department of Physics, Jinan University, Guangzhou, Guangdong 510632, China ‡ Center for Functional Nanomaterials, Brookhaven National Laboratory, Upton, New York 11973, United States § College of Chemistry and Materials Science, Hubei Engineering University, Xiaogan 432000, China S Supporting Information *

ABSTRACT: We report on a density functional theory (DFT) calculation to determine the reaction pathway and barrier of the hydrogen evolution reaction (HER) at the interface of monolayer 1T′ phase MoS2 and water. By screening the interfacial structures with the lowest chemical potential of protons and electrons, key structural characteristics that are important for prediction of reaction mechanism are identified. Under typical reaction conditions of HER, the catalyst surface features a high proton coverage of ca. 37% while the aqueous solution has a relatively low hydronium concentration of no more than ca. 1.8%. This contrast leads to proton desorption from the catalyst surface through a diffusionassisted Tafel manner, rather than the Heyrovsky manner assumed previously. The result is supported by the agreement of the calculated reaction barrier and surface coverage with those of experimental estimate. In prediction of catalytic activity, hydrogen adsorption energies of reaction intermediates are widely used as the thermodynamic descriptor, while reaction barriers usually serve as kinetic parameters. We suggest that both thermodynamic and kinetic description toward HER should be performed on the premise that the lowest chemical potential of protons and electrons is obtained.

1. INTRODUCTION

effective catalysts and establishing of heterogeneous catalytic theory for prediction of reactivity. Having a very high hydrogen coverage (ca. 100%), a Pt surface is known to catalyze hydrogen generation through the Volmer−Tafel mechanism, and the desorption step is the ratelimiting step.22,23 In contrast, the low hydrogen coverage on the edge sites of MoS2 makes the Volmer−Heyrovsky mechanism dominant for edges’ reactivity, which has been confirmed theoretically and experimentally.11,24 Using a hydrogen adsorption model without water addition, the hydrogen coverage on the basal plane of 1T′ phase MoS2 was calculated to be ca. 10−20% based on density functional theory (DFT) calculations.17,25 This leads to the prediction and wide adoption of the Volmer−Heyrovsky pathway as the mechanism of the HER on MoS2 and other transition metal disulfides (TMDs).25−29 Although this prediction is supported by a recent volcano plot combining experimental turnover frequencies (TOFs) and DFT-calculated ΔGH,30 experimental mechanistic studies per se are inconclusive because of the discrepancy of the wide range of Tafel slopes reported.19−21,30 Furthermore, the experimentally estimated active site density

The hydrogen evolution reaction (HER) has attracted a great deal of attention due to the ongoing pursuit of hydrogen as a future energy carrier.1−3 Platinum is the most efficient inorganic HER catalyst, but its application in large scale hydrogen generation is limited by the high cost.3 Among all the potential alternatives, MoS2 has gained the most attention, not only because of its nearly thermoneutral hydrogen adsorption energy (ΔGH) but also due to its diversity in both crystal phases and nanoscale morphology, which provides great opportunities for engineering.4−9 The HER activity arising from the edge sites of MoS2 was first confirmed,7 which subsequently stimulated extensive effort to synthesize MoS2 catalysts with rich edge sites.6,10−14 However, the edges of 2H phase MoS2 suffer from the problem of low active site proportion of only ca. 20%,4,7 and the catalytic performance is limited by the poor conductivity due to the semiconducting nature.11,15,16 In contrast, 1T′ phase MoS2 is intrinsically metallic, which facilitates electron conduction.17,18 Experiments show that both MoS2 and WS2 nanosheets of 1T′ phase present superior HER activity that outperforms the best records of their corresponding 2H counterparts.19−21 Understanding of the detailed reaction mechanism of the HER has become a subject of general interest in both exploring of cost© 2018 American Chemical Society

Received: May 28, 2018 Revised: July 7, 2018 Published: July 8, 2018 5404

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411

Article

Chemistry of Materials

Figure 1. (a) Side and (b) top view of an electrochemical double layer (EDL) interfacial model consisting of monolayer 1T′ phase MoS2 and two layers of water molecules. Two categories of sulfur atoms labeled by S-I and S-II are identified on the catalyst surface, which differs from each other in the normal height. The rectangular 2 × 2 super cell lined out by the black box in (b) is used for DFT calculations, which consists of 8 MoS2 formula units and 14 water molecules.

(ca. 4.5 × 1014 sites·cm−2) is much higher than that predicted by calculations.20,21 To conclude the reaction mechanism, one needs to match the calculated reaction barrier and hydrogen coverage with those of experiments. The main problem in computational aspects relates to the structural models used. Hydrogen adsorption models fail to consider the aqueous environment and interfacial electric field, which might result in an incorrect estimate of hydrogen coverage and reaction mechanism. Even with water added, a solid/liquid interfacial model should also be selected with great care, because the interfacial configuration is strongly dependent on pH and electrode potential.31,32 Even under a given condition, aqueous structure is amorphous and changeable.33 Selection of different interfacial structures may result in variation of calculated reaction barriers and influence evaluation of the reaction mechanism. Now a question may arise on what kind of interfacial structure should be selected for calculations. The general thermodynamic criterion is that the structure should have the lowest chemical potential of protons and electrons (i.e. the electrochemical potential), so that the system could be maximally stabilized against forming molecular hydrogen. This criterion requires that both pH and potential enter into the thermodynamic screening process. The obtained structural characteristics at different pH and potential could be summarized in a Pourbaix diagram,31,32 which provide important reference for evaluating of the mechanism of heterogeneous catalytic reactions in general and the HER in particular. In this work, DFT calculations were carried out to obtain the Pourbaix diagram for the electrochemical double layer (EDL) system consisting of a monolayer of 1T′ phase MoS2 and two layers of water. Atomic-level structural characteristics including proton concentration, surface proton coverage, and water orientation were obtained for wide ranges of pH and potential. Under typical reaction conditions of HER, where the electrode potential is close to zero and pH lies in between 0 and 1, the interfacial system is evaluated to have a high hydrogen coverage of ca. 37% and a relatively low hydronium concentration of no more than 1.8%, which facilitates proton desorption through the Tafel manner. Kinetic barrier calculations confirmed that protons desorbed from the catalyst surface via a diffusion-assisted Tafel process, rather than the Heyrovsky process assumed previously.

monolayer of 1T′-MoS2 covered by two layers of water, as shown in Figure 1a,b. A rectangular super cell consisting of 8 MoS2 formula units (f.u.) and 14 water molecules was used. In the monolayer of MoS2, zigzag-type Mo plane distortion could be seen, which is related to the charge density wave (CDW) induced by strong Coulomb correlation.34,35 As a result, two categories of sulfur atoms are identified on the surface. One is anchored on the zigzag line of Mo ions and is higher in the normal direction. The other lies in between two zigzag lines and hence is lower than the first one. The former is marked with site I, while the latter is denoted with site II. Interfacial structures with different water configurations were obtained by combining molecular dynamics simulations and geometry optimization. Proton concentration could be changed by adding hydrogen atoms into or removing them from the system. H atoms either adsorb on the catalyst surface or get solvated to form protons. Such interfacial models include physical water structures, which simulate the real scenario of a Helmholtz double layer (HDL) in an acid electrolyte. Under HER conditions, electrons on the electrode surface cause electrostatic repulsion of acid radicals. In contrast, protons come up to balance the negative surface charges. A charge balance between the injected electrons in the electrode and protons in the solution is then set up. 2.2. Computational Methods. Density functional theory (DFT) calculations were carried out using the Vienna Ab-initio Simulation Package (VASP).36 Interactions of electrons with ion cores were represented by the projector augmented wave (PAW) potential.37,38 The optB86b exchange-correlation functional was used to describe the weak van der Waals (vdW) interactions.39−41 A cutoff energy of 500 eV and a two-dimensional Monkhorst−Pack special mesh of 4 × 8 were used for the supercell. Full structural relaxations were performed until all the forces acting on all the ions were smaller than 0.02 eV/Å. The Climbing-Image Nudged Elastic Band (CI-NEB)42 method was used to locate transition states. Free energies were obtained by correcting zero-point energies and entropies, as elaborated in section 2 of Supporting Information (SI). 2.3. Interfacial Thermodynamics. The electrode potential USHE, which is set up by interfacial dipoles, is calculated by referencing the system’s work function Φ to that of the Standard Hydrogen Electrode (ΦSHE).

USHE = (Φ − ΦSHE)/e

(1)

ΦSHE is experimentally measured to be 4.44 eV. Under the SHE condition, protons and electrons are in equilibrium with H2 molecules. The total free energy of protons and electrons (GH++e−) in an interfacial system could be calculated with 43

G H++ e−(μH++ e− = 0) =

G(N , n) − G(N , 0) n − GH N 2N 2

(2)

where N is the number of water molecules, and n is the number of protons, including both hydroniums and adsorbed protons. A negative n indicates the number of hydroxyls in the system. GH2 is the free energy of an H2 molecule at 1 atm and 298 K, while G(N,n) and

2. MODELS AND METHODS 2.1. Model Systems. The reaction mechanism was studied by constructing a solid/liquid interfacial model, which contained a 5405

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411

Article

Chemistry of Materials

Figure 2. Calculated GH++e−(μH++e−) vs USHE for the 1T′-MoS2/H2O interface at (a) pH = 0 and (b) pH = 14. (c) Simulated Pourbaix diagram for the 1T′-MoS2/solution interface, which exhibits critical structural characteristics as a function of pH and USHE. The shaded area represents the typical reaction condition of HER. Different proton concentrations are differentiated by the color bar, while the characteristics of proton distribution and water orientation are tabulated in d. Psol and Pad represents hydroniums and adsorbed protons, respectively, and mPsol + nPad denotes the distribution with m hydroniums in solution and n adsorbed protons on the catalyst surface. Upper triangles indicate possible adsorbed species, i.e. either adsorbed protons of different numbers or one adsorbed hydroxyl. Pentacles represent coexistence of four different proton distribution situations. Water orientation is indicated by symbols’ size, which is scaled by the number ratio of H-up waters (↑) to H-down ones (↓). G(N,0) refer to the free energies of an interfacial structure with n protons and the reference system without protons, respectively. For conditions away from the SHE, the total chemical potential of protons and electrons varies with USHE and pH through31 μH++ e− = − eUSHE − 2.3kT pH

free energies for each proton concentration. Six different proton concentrations are considered, which are indicated by distinct colors shown in the color bar. Data for equal proton concentration align linearly. For each of the proton concentrations, the interfacial structures with different proton distributions across the interface are further differentiated by symbol shapes. Generally, protons tend to be solvated at lower potentials but prefer to be adsorbed on the catalyst surface at higher potentials. Because dipole directions of water structures also influence the electrode potential, all data dots are scaled in size according to the number ratio of H-up waters (↑) to Hdown ones (↓), as indicated by the legends in Figure 2d. When pH = 0, five intersections indicated by the dashed lines divide the potential range into six regions, each of which has unique structural characteristics in terms of proton concentration, proton distribution, and water orientation for the most stable interfacial configuration. As pH increases, these regions move to the negative side with varied width. By linear interpolation of the structural characteristics of each region from pH = 0 to pH = 14, the Pourbaix diagram is obtained and shown in Figure 2c. The Pourbaix diagram exhibits six phase regions separated by five dashed lines. The slopes of these dashed lines deviate more or less from the ideal value of −0.059 (−2.3kT), indicating that the thermodynamics of this system does not conform to a linear URHE (i.e. the reversible hydrogen electrode) dependence. Instead, water environment, interface electric field, and the surface condition may have an impact. Overall, lowering of USHE leads to an increase of proton

(3)

-eUSHE describes the contribution of USHE to the chemical potential of electrons, while −2.3kTpH quantifies the contribution of pH to the chemical potential of protons. The total free energy of electrons and protons can then be obtained by a linear extrapolation of GH++e−(μH++e− = 0),32

G H++ e−(μH++ e− ) = G H++ e−(μH++ e− = 0) −

n μ + − N H +e

(4)

By comparing GH++e−(μH++e−) of different interfacial structures, the most stable configuration with the lowest electrochemical potential at specified potential and pH could be obtained. In this work, we have considered a large set of interfacial structures (>130), which differ from each other in terms of proton concentration, proton distribution across the interface, and water orientation. For each of the interfacial structures, GH++e−(μH++e− = 0) was calculated first, and then GH++e−(μH++e−) for different potential and pH was obtained using eq 4. Kinetic barrier calculations at specified potential and pH were performed for the most stable structure.

3. RESULTS AND DISCUSSION 3.1. Free Energies and Pourbaix Diagrams. Figure 2a,b shows the calculated free energies, GH++e−(μH++e−), as a function of USHE for pH = 0 and 14, respectively. The scenario for pH = 7 is displayed in Figure S1. These data points are selected from the full data set of Figure S2 because they represent the lowest 5406

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411

Article

Chemistry of Materials

This suggests that S-II sites are more favorable for proton adsorption. Figure 3a shows the minimum energy path (MEP) of the Volmer reaction at USHE = 0.15 V. The structural diagrams on

concentration and shift from hydroxyl-covered surfaces to proton-covered ones. This is consistent with the experimental observation that catalytic performance is always enhanced when the electrode potential is decreased. For a given proton concentration, lowering of USHE induces proton transfer from the catalyst surface to the solution and reorientation of waters from the H-down configuration to the H-up one. We now focus on the pH range from 0 to 1, where the HER activity is usually measured. Within this range, phase region II has positive potential and contains no excess protons in the system. This explains why no HER activity is experimentally observed in this region. With the decrease of USHE, the number of hydroniums increases. When entering region III, it reaches the molar concentration of 1/14 (one proton per supercell). The typical interfacial configurations of region III are represented by olive lower triangles. Further lowering of USHE induces proton adsorption onto the catalyst surface. In phase region IV, all protons in solution have been transferred onto the surface, resulting in a high surface coverage of 37.5% and simultaneously a “depletion” of hydroniums in solution. The structural characteristic of this region is represented by red upper triangles, which is gradually transformed to red right triangles when pH is increased. For the latter case, the hydronium concentration is 1/14 and the surface coverage becomes 25%. Interestingly, the predicted surface coverage using models without water layers is only ca. 12%.17,25,44 By including a water environment to obtain the lowest electrochemical potential, the surface coverage is proved to be much higher than that anticipated previously. This high coverage is stabilized by interfacial hydrogen bonds and charge transfer, which facilitates proton desorption through the Tafel process. On the contrary, the Heyrovsky process should be difficult due to the relative scarcity of hydroniums. In fact, HER is usually operated in the pH range of 0−1 and in the potential range of −0.3−0 V, which just corresponds to the shaded area on the upper left of region IV (cf. also Figure S3). Within this area, the concentration of hydroniums is only ca. 0.18%−1.8%, far less than the surface adsorbed protons. Phase region V is the transition zone between regions IV and VI and thus contains the interfacial types of both latter two regions. Except for the interfacial structures represented by red upper triangles and right triangles, two new interfacial types (blue squares and upper triangles) corresponding to a higher proton concentration of 0.29 enter this region. These interfacial configurations have quite similar free energies and therefore are all represented by pentacles. For the configuration of blue squares, the molar concentration of hydroniums is 1/14, while the surface proton coverage remains at 37.5%. For that of blue upper triangles, there are no hydroniums in the solution. However, the surface proton coverage achieves 50%, i.e. all the surface S-II sites are covered by protons. For the latter case, the Heyrovsky manner for desorption is probably hindered. In section 3 of SI, we have introduced a short discussion on the scenario that larger supercells are used or more water molecules are added. 3.2. The Volmer Reaction. The Volmer reaction is the first step of the HER, in which a hydronium transfers from the solution to the surface of a catalyst. Because there are two types of S sites on the surface, one needs to determine which one is more favorable for proton adsorption. We use the reference system to calculate GH++e−(μH++e− = 0) for adsorbed protons on sites I and II, which are 0.05 and 0 eV, respectively.

Figure 3. (a) Top: initial state (IS), transition state (TS), and final state (FS) structures for a Volmer reaction at USHE = 0.15 V. Bottom: DFT-calculated minimum energy path (black dots) and the chargeextrapolation corrected result (red dots). The only imaginary frequency of the saddle point is shown. Blue crosses indicate the work function change along the path. (b) Linearly fitted Volmer reaction barriers for different potentials.

the upper panel depict the initial state (IS), transition state (TS), and final state (FS) of the reaction, from which the process of proton transfer is clear to see. The IS structure is selected from the Pourbaix diagram, which has one hydronium in solution but no adsorbed protons on the catalyst surface. The black curve represents the DFT calculated energies along the reaction pathway, while the blue one shows the synchronal change of the work function. The DFT calculated barrier (Ea) is 0.43 eV, while the reaction energy Er = EFS − EIS is calculated to be 0.23 eV. Due to the effect of finite supercell size, the work function fluctuates by several hundred meV along the reaction pathway, which nullifies the intention to get accurate barriers and reaction energies at constant potential. We use the charge extrapolation correction scheme proposed recently by Norskov et al.45 to extrapolate all the energies 5407

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411

Article

Chemistry of Materials

Figure 4. (a) Calculated MEP for the Tafel reaction between two neighboring protons adsorbed on S-II sites at USHE = 0.09 V. The IS structure contains no hydroniums but has three adsorbed protons on the surface S-II sites. (b) MEP for proton migration from an S-II site to a neighboring S-I site. (c) MEP for the Tafel reaction between the migrated proton and an adsorbed proton at neighboring S-II site. Structural changes of these processes are shown on the upper panel. Because no charge transfers across the interface, all the energies are uncorrected for potential deviation. Imaginary frequencies are labeled near the corresponding saddle points.

Figure 5. (a) Calculated MEP for the Heyrovsky reaction between a proton adsorbed on an S-II site and a hydronium in solution at USHE = −0.63 V. The IS structure has a hydronium in the solution and three adsorbed protons on the surface S-II sites. (b) MEP for proton migration from an SII site to a neighboring S-I site. (c) MEP for the Heyrovsky reaction between the migrated proton and a hydronium. Black dots are DFT calculated energies, while red dots are charge-extrapolation corrected results. Blue crosses indicate the work function change along the path. Structural changes of these processes are shown on the upper panel. The proton migration path in b is not corrected because there is no charge transfer across the interface. Imaginary frequencies are labeled near the corresponding saddle points.

the Tafel manner. For verification, we select an interfacial structure at USHE = 0.09 V and pH = 0 to simulate the Tafel reaction. This structure has a potential close to the experimental range, which facilitates direct understanding of the Tafel process under HER conditions. The MEP of direct recombination of two neighboring protons adsorbed on S-II sites is shown in Figure 4a. Because there is no charge transfer across the interface, all the energies along the path are not corrected for potential deviation. The barrier of direct recombination is calculated to be 1.29 eV. As an alternative, it might also be possible for an adsorbed proton to migrate to a neighboring S-I site first and then to combine with an S-II proton to form a H2 molecule. This twostep process is shown in Figure 4b,c. Interestingly, the proton migration process takes only a small barrier of 0.79 eV. After that, the barrier of the Tafel process drops to 1.06 eV. Upon corrections of zero-point energies and entropies (cf. section 4 of SI for details), the room-temperature free energy barrier

along the MEP to those of the initial potential. The results are shown by red dots. After correction, the barrier is reduced to 0.35 eV and the reaction energy becomes −0.06 eV, indicating the exothermic nature of this reaction. The Volmer reaction at other potentials is generally similar. By selecting the most stable configurations as the IS structures, we obtain the Volmer reaction barriers at several different potentials, as shown in Figure 3b. The barrier increases with the electrode potential overall. By a linear extrapolation, Ea at the equilibrium potential is determined to be 0.30 eV. Such a small barrier suggests that the Volmer reaction is quite easy, in agreement with previous theoretical and experimental reports.19−21,24,25 3.3. The Tafel Reaction. According to the Pourbaix diagram, under HER conditions a 1T′-MoS2/H2O interface has a high surface coverage of 37.5% and a relatively low hydronium concentration of no more than 1.8%. This structural characteristic facilitates proton desorption through 5408

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411

Article

Chemistry of Materials (ΔGRT) is determined to be 0.95 eV. Obviously, this two-step process is more favorable than the direct recombination manner. Particularly, the moderate barrier of 0.95 eV demonstrates a great possibility for the Tafel mechanism to dominate the desorption process. For potentials lower than zero, the scenario is quite similar, which is exemplified by Figure S4. 3.4. The Heyrovsky Reaction. We also consider the possibility of the Heyrovsky reaction, in which a hydronium in solution combines with an adsorbed proton on the catalyst surface. This process could be slowed if either hydroniums or adsorbed protons are relatively scarce. As aforementioned, the reaction conditions of HER just belong to this case. For confirmation from the kinetic viewpoint, we use the most stable interfacial structure at USHE = −0.63 V and pH = 0 to study the Heyrovsky reaction. This structure has one hydronium in solution and three adsorbed protons on the surface S-II sites, which represents one of the typical interfacial types of region V. Structures with potentials close to the equilibrium are not considered because no hydroniums are available in supercells of the current size. The DFT calculated reaction MEP and the charge-extrapolation corrected result are both shown in Figure 5a. The corrected barrier is as high as 1.94 eV, indicating that the Heyrovsky reaction on S-II sites is rather difficult. The reason can be ascribed to the strong proton adsorption on S-II sites. As an alternative manner, we also consider the possibility for an adsorbed proton to migrate to a neighboring S-I site, after which it reacts with a hydronium in the solution. This two-step pathway is shown in Figure 5b,c, where all energies have been corrected by being extrapolated to those of USHE = −0.63 V using the charge extrapolation scheme. Interestingly, the proton migration process also takes a small barrier of 0.59 eV, confirming again that proton diffusion on the catalyst surface is easy. After the migration, the Heyrovsky process between the migrated proton and a hydronium reduces to 1.45 eV, indicating that the Heyrovsky reaction also takes place in a diffusion-assisted manner. This finding is also confirmed for another potential of −0.79 V (cf. Figure S5). After corrections of zero-point energies and entropies, we obtain a ΔGRT of 1.20 eV. For potentials close to zero, one could expect a much higher barrier when being extrapolated, which makes the Heyrovsky process less favorable than the Tafel one. All barriers here have been obtained using the interfacial structures with the lowest electrochemical potential. This ensures that the calculated barrier reflect the realistic kinetic process at given potential and pH. As a contrast, we recalculate the Heyrovsky barrier around −0.79 V by using an interfacial structure with a higher energy. This structure has a proton concentration of 0.21 and does not represent the most stable configuration at this potential. The calculated barrier (cf. Figure S6) based on this structure is only 0.89 eV, much smaller than the 1.20 eV for the scenario of the lowest electrochemical potential (Figure S5). 3.5. The Overall Reaction. To clarify the rate-limiting step, we calculate ΔGRT for both the Tafel and the Heyrovsky reactions at several different potentials for pH = 0. The results are plotted in Figure 6. It is clear to see the positive correlation of ΔGRT with potential for both the reactions. By linear fitting, the barrier of either of the reactions at specified potential could be obtained. Interestingly, within the considered potential range the diffusion-assisted Tafel reaction exhibits much lower ΔGRT than that of the Heyrovsky reaction, confirming that

Figure 6. Linearly fitted room-temperature free energy barriers (ΔGRT) for the Heyrovsky (black dots) and the Tafel reactions (red dots) at pH = 0.

protons desorb from the catalyst surface through the diffusionassisted Tafel mechanism, rather than the Heyrovsky mechanism assumed previously. The main reason responsible for this discrepancy is ascribed to the used models, which have different surface coverage and represent the lowest electrochemical potential under given condition. We can now explain why basal planes of 1T′ phase MoS2 have good catalytic activity. It is the proton distribution across the EDL interface that plays a critical role. For one hand, the relative scarcity of hydroniums in solution might be unfavorable for the Heyrovsky process. For another, the high surface coverage facilitates fast migration and facile desorption of protons through the Tafel manner. Interestingly, this interfacial characteristic is also supported by a recent experiment of Voiry et al.20 Using the copper underpotential deposition method, they estimated the active site density of strained 1T′ phase WS2 (with ca. 2% dilatation) to be ca. 4.5 × 1014 sites·cm−2. This data corresponds to a surface coverage of ca. 45%, or equivalently 3.5 S sites per supercell of 8 WS2 f.u. According to a previous calculation, a tensile strain of 2% may induce an increase of coverage by ca. 10%.44 It is therefore reasonable to estimate that the unstrained sample may have a surface coverage of ca. 35%. Considering the similarity between WS2 and MoS2, the predicted coverage of 37.5% is in good agreement with the experimental one. The calculated ΔGRT can also be compared with the estimated reaction barriers from experiments. At the condition of USHE = 0 V and pH = 0.3, Voiry estimated a turnover frequency (TOF) of ca. 0.04 s−1 per surface S site for hydrogen evolution on 1T′ phase WS2.20 For the same condition, Yin et al. gave an exchange current density of 12.6 μA·cm−2 for 1T′ phase MoS2,21 which corresponds to a TOF of ca. 0.02 s−1. Using the rate expression of the transition state theory, kBT ji −ΔG RT zyz expjjj z j kBT zz h (5) k { we estimate that the experimental barriers of 1T′-WS2 and 1T′-MoS2 are 0.85 and 0.87 eV, respectively. After adjusting the chemical potential of protons to pH = 0.3 and that of electrons to USHE = 0 V, the corrected Tafel barrier becomes 0.91 eV, very close to the experimental results. Considering that Voiry’s experiment is based on a strained sample with a higher coverage, it is reasonable to estimate that the unstrained k=

5409

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411

Article

Chemistry of Materials

Public Opinion and Big Data are thanked for providing computational support. This work is supported by the National Natural Science Foundation of China (no. 21303237, 21703081, 11404142, and 11647131), the Natural Science Foundation of Guangdong Province (no. 2018A030313386), the Special Program for Applied Research on Super Computation of the NSFC-Guangdong Joint Fund (no. U1501501), the Fundamental Research Foundation for the Central Universities (no. 21617330), and the Natural Science Foundation of Hubei Province (no. 2018CFB669).

sample may have a slightly larger barrier that is closer to our calculated result.

4. CONCLUSIONS By screening the 1T′-MoS2/H2O interfacial structures with the lowest electrochemical potential, we identified the key structural characteristics that are important for prediction of the reaction mechanism. Under the experimental condition of HER, the interfacial system has a high surface coverage of 37.5% and a relatively low hydronium concentration of no more than 1.8%, consistent with the estimate from the experimental density of active sites. This contrast of proton number declines the possibility of the Heyrovsky process as the rate-limiting step. Instead, hydrogen evolution goes through the Volmer−Tafel mechanism, in which proton migration on the catalyst surface assists the desorption process. The calculated barrier for this mechanism is 0.91 eV, in good agreement with the 0.87 eV estimated from experiments. Hydrogen adsorption energies of reaction intermediates have long been viewed as the thermodynamic descriptor toward HER, which justifies the use of a simple hydrogen adsorption model without water addition to describe the reactivity of catalysts. However, we find that the reaction mechanism of HER could not be determined properly using this model, because it fails to give proper description of protonic chemical potential and then proton distribution. A reasonable model system should contain the catalyst/water interface and meanwhile has the lowest chemical potential of protons and electrons, which ensure that both the thermodynamic descriptors and kinetic calculations are valid. In conclusion, to design a more active catalyst based on 1T′ phase MoS2, one should focus on lowering the barrier of the Tafel process and accelerating the diffusion of protons on the catalyst’s surface. This work provides new and comprehensive insights for exploring of the reaction mechanism of heterogeneous catalysis, which invites further theoretical and experimental verification for other solid/liquid interfacial systems.





(1) Turner, J. A. Sustainable Hydrogen Production. Science 2004, 305, 972−974. (2) Hou, Y.; Abrams, B. L.; Vesborg, P. C. K.; Björketun, M. E.; Herbst, K.; Bech, L.; Setti, A. M.; Damsgaard, C. D.; Pedersen, T.; Hansen, O.; Rossmeisl, J.; Dahl, S.; Nørskov, J. K.; Chorkendorff, I. Bioinspired molecular co-catalysts bonded to a silicon photocathode for solar hydrogen evolution. Nat. Mater. 2011, 10, 434−438. (3) Faber, M. S.; Jin, S. Earth-abundant inorganic electrocatalysts and their nanostructures for energy conversion applications. Energy Environ. Sci. 2014, 7, 3519−3542. (4) Hinnemann, B.; Moses, P. G.; Bonde, J.; Jørgensen, K. P.; Nielsen, J. H.; Horch, S.; Chorkendorff, I.; Nørskov, J. K. Biomimetic Hydrogen Evolution: MoS2 Nanoparticles as Catalyst for Hydrogen Evolution. J. Am. Chem. Soc. 2005, 127, 5308−5309. (5) Laursen, A. B.; Kegnaes, S.; Dahl, S.; Chorkendorff, I. Molybdenum sulfides-efficient and viable materials for electro - and photoelectrocatalytic hydrogen evolution. Energy Environ. Sci. 2012, 5, 5577−5591. (6) Karunadasa, H. I.; Montalvo, E.; Sun, Y.; Majda, M.; Long, J. R.; Chang, C. J. A Molecular MoS2 Edge Site Mimic for Catalytic Hydrogen Generation. Science 2012, 335, 698−702. (7) Jaramillo, T. F.; Jørgensen, K. P.; Bonde, J.; Nielsen, J. H.; Horch, S.; Chorkendorff, I. Identification of Active Edge Sites for Electrochemical H2 Evolution from MoS2 Nanocatalysts. Science 2007, 317, 100−102. (8) Merki, D.; Hu, X. Recent developments of molybdenum and tungsten sulfides as hydrogen evolution catalysts. Energy Environ. Sci. 2011, 4, 3878−3888. (9) Wilson, J. A.; Yoffe, A. D. The transition metal dichalcogenides discussion and interpretation of the observed optical, electrical and structural properties. Adv. Phys. 1969, 18, 193−335. (10) Kong, D.; Wang, H.; Cha, J. J.; Pasta, M.; Koski, K. J.; Yao, J.; Cui, Y. Synthesis of MoS2 and MoSe2 Films with Vertically Aligned Layers. Nano Lett. 2013, 13, 1341−1347. (11) Li, Y.; Wang, H.; Xie, L.; Liang, Y.; Hong, G.; Dai, H. MoS2 Nanoparticles Grown on Graphene: An Advanced Catalyst for the Hydrogen Evolution Reaction. J. Am. Chem. Soc. 2011, 133, 7296− 7299. (12) Wang, T.; Liu, L.; Zhu, Z.; Papakonstantinou, P.; Hu, J.; Liu, H.; Li, M. Enhanced electrocatalytic activity for hydrogen evolution reaction from self-assembled monodispersed molybdenum sulfide nanoparticles on an Au electrode. Energy Environ. Sci. 2013, 6, 625− 633. (13) Benck, J. D.; Hellstern, T. R.; Kibsgaard, J.; Chakthranont, P.; Jaramillo, T. F. Catalyzing the Hydrogen Evolution Reaction (HER) with Molybdenum Sulfide Nanomaterials. ACS Catal. 2014, 4, 3957− 3971. (14) Ye, G.; Gong, Y.; Lin, J.; Li, B.; He, Y.; Pantelides, S. T.; Zhou, W.; Vajtai, R.; Ajayan, P. M. Defects Engineered Monolayer MoS2 for Improved Hydrogen Evolution Reaction. Nano Lett. 2016, 16, 1097− 1103. (15) Tsai, C.; Abild-Pedersen, F.; Nørskov, J. K. Tuning the MoS2 Edge-Site Activity for Hydrogen Evolution via Support Interactions. Nano Lett. 2014, 14, 1381−1387.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.chemmater.8b02236. Supplementary figures and tables, calculation details for free energies, Pourbaix diagram for lager supercells with more water molecules, and calculations of room temperature free energy barriers (PDF)



REFERENCES

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Xiaobo Chen: 0000-0002-8042-5927 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS Mr. Xing Lu, the high-performance computing platform of Jinan University, and the Guangzhou Key Research Center of 5410

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411

Article

Chemistry of Materials (16) Yu, Y.; Huang, S.-Y.; Li, Y.; Steinmann, S. N.; Yang, W.; Cao, L. Layer-Dependent Electrocatalysis of MoS2 for Hydrogen Evolution. Nano Lett. 2014, 14, 553−558. (17) Chen, X.; Gu, Y.; Tao, G.; Pei, Y.; Wang, G.; Cui, N. Origin of hydrogen evolution activity on MS2 (M = Mo or Nb) monolayers. J. Mater. Chem. A 2015, 3, 18898−18905. (18) Gao, G.; Jiao, Y.; Ma, F.; Jiao, Y.; Waclawik, E.; Du, A. Charge Mediated Semiconducting-to-Metallic Phase Transition in Molybdenum Disulfide Monolayer and Hydrogen Evolution Reaction in New 1T′ Phase. J. Phys. Chem. C 2015, 119, 13124−13128. (19) Lukowski, M. A.; Daniel, A. S.; Meng, F.; Forticaux, A.; Li, L.; Jin, S. Enhanced Hydrogen Evolution Catalysis from Chemically Exfoliated Metallic MoS2 Nanosheets. J. Am. Chem. Soc. 2013, 135, 10274−10277. (20) Voiry, D.; Yamaguchi, H.; Li, J.; Silva, R.; Alves, D. C. B.; Fujita, T.; Chen, M.; Asefa, T.; Shenoy, V. B.; Eda, G.; Chhowalla, M. Enhanced catalytic activity in strained chemically exfoliated WS2 nanosheets for hydrogen evolution. Nat. Mater. 2013, 12, 850−855. (21) Yin, Y.; Han, J.; Zhang, Y.; Zhang, X.; Xu, P.; Yuan, Q.; Samad, L.; Wang, X.; Wang, Y.; Zhang, Z.; Zhang, P.; Cao, X.; Song, B.; Jin, S. Contributions of Phase, Sulfur Vacancies, and Edges to the Hydrogen Evolution Reaction Catalytic Activity of Porous Molybdenum Disulfide Nanosheets. J. Am. Chem. Soc. 2016, 138, 7965−7972. (22) Conway, B. E.; Tilak, B. V. Interfacial processes involving electrocatalytic evolution and oxidation of H2, and the role of chemisorbed H. Electrochim. Acta 2002, 47, 3571−3594. (23) Skúlason, E.; Tripkovic, V.; Björketun, M. E.; Gudmundsdóttir, S.; Karlberg, G.; Rossmeisl, J.; Bligaard, T.; Jónsson, H.; Nørskov, J. K. Modeling the Electrochemical Hydrogen Oxidation and Evolution Reactions on the Basis of Density Functional Theory Calculations. J. Phys. Chem. C 2010, 114, 18182−18197. (24) Huang, Y.; Nielsen, R. J.; Goddard, W. A.; Soriaga, M. P. The Reaction Mechanism with Free Energy Barriers for Electrochemical Dihydrogen Evolution on MoS2. J. Am. Chem. Soc. 2015, 137, 6692− 6698. (25) Tang, Q.; Jiang, D. Mechanism of Hydrogen Evolution Reaction on 1T-MoS2 from First Principles. ACS Catal. 2016, 6, 4953−4961. (26) Tsai, C.; Chan, K.; Nørskov, J. K.; Abild-Pedersen, F. Theoretical insights into the hydrogen evolution activity of layered transition metal dichalcogenides. Surf. Sci. 2015, 640, 133−140. (27) Zhang, Y.; Chen, X.; Huang, Y.; Zhang, C.; Li, F.; Shu, H. The Role of Intrinsic Defects in Electrocatalytic Activity of Monolayer VS2 Basal Planes for the Hydrogen Evolution Reaction. J. Phys. Chem. C 2017, 121, 1530−1536. (28) Yuan, J.; Wu, J.; Hardy, W. J.; Loya, P.; Lou, M.; Yang, Y.; Najmaei, S.; Jiang, M.; Qin, F.; Keyshar, K.; Ji, H.; Gao, W.; Bao, J.; Kono, J.; Natelson, D.; Ajayan, P. M.; Lou, J. Facile Synthesis of Single Crystal Vanadium Disulfide Nanosheets by Chemical Vapor Deposition for Efficient Hydrogen Evolution Reaction. Adv. Mater. 2015, 27, 5605−5609. (29) Liang, H.; Shi, H.; Zhang, D.; Ming, F.; Wang, R.; Zhuo, J.; Wang, Z. Solution Growth of Vertical VS2 Nanoplate Arrays for Electrocatalytic Hydrogen Evolution. Chem. Mater. 2016, 28, 5587− 5591. (30) Zhang, J.; Wu, J.; Guo, H.; Chen, W.; Yuan, J.; Martinez, U.; Gupta, G.; Mohite, A.; Ajayan, P. M.; Lou, J. Unveiling Active Sites for the Hydrogen Evolution Reaction on Monolayer MoS2. Adv. Mater. 2017, 29, 1701955. (31) Hansen, H. A.; Rossmeisl, J.; Norskov, J. K. Surface Pourbaix diagrams and oxygen reduction activity of Pt, Ag and Ni(111) surfaces studied by DFT. Phys. Chem. Chem. Phys. 2008, 10, 3722−3730. (32) Rossmeisl, J.; Chan, K.; Ahmed, R.; Tripkovic, V.; Bjorketun, M. E. pH in atomic scale simulations of electrochemical interfaces. Phys. Chem. Chem. Phys. 2013, 15, 10321−10325. (33) Björneholm, O.; Hansen, M. H.; Hodgson, A.; Liu, L.-M.; Limmer, D. T.; Michaelides, A.; Pedevilla, P.; Rossmeisl, J.; Shen, H.; Tocci, G.; Tyrode, E.; Walz, M.-M.; Werner, J.; Bluhm, H. Water at Interfaces. Chem. Rev. 2016, 116, 7698−7726.

(34) Shen, D. W.; Xie, B. P.; Zhao, J. F.; Yang, L. X.; Fang, L.; Shi, J.; He, R. H.; Lu, D. H.; Wen, H. H.; Feng, D. L. Novel Mechanism of a Charge Density Wave in a Transition Metal Dichalcogenide. Phys. Rev. Lett. 2007, 99, 216404. (35) Zhu, Z.; Cheng, Y.; Schwingenschlögl, U. Origin of the charge density wave in 1T-TiSe2. Phys. Rev. B: Condens. Matter Mater. Phys. 2012, 85, 245133. (36) Kresse, G.; Furthmüller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 11169−11186. (37) Blöchl, P. E. Projector augmented-wave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 50, 17953−17979. (38) Kresse, G.; Joubert, D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 1758−1775. (39) Dion, M.; Rydberg, H.; Schröder, E.; Langreth, D. C.; Lundqvist, B. I. Van der Waals Density Functional for General Geometries. Phys. Rev. Lett. 2004, 92, 246401. (40) Román-Pérez, G.; Soler, J. M. Efficient implementation of a van der Waals density functional: application to double-wall carbon nanotubes. Phys. Rev. Lett. 2009, 103, 096102. (41) Klimeš, J.; Bowler, D. R.; Michaelides, A. Van der Waals density functionals applied to solids. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 83, 195131. (42) Henkelman, G. A climbing image nudged elastic band method for finding saddle points and minimum energy paths. J. Chem. Phys. 2000, 113, 9901−9904. (43) Randles, J. E. B. The real hydration energies of ions. Trans. Faraday Soc. 1956, 52, 1573−1581. (44) Chen, X.; Wang, G. Tuning the hydrogen evolution activity of MS2 (M = Mo or Nb) monolayers by strain engineering. Phys. Chem. Chem. Phys. 2016, 18, 9388−9395. (45) Chan, K.; Nørskov, J. K. Electrochemical Barriers Made Simple. J. Phys. Chem. Lett. 2015, 6, 2663−2668.

5411

DOI: 10.1021/acs.chemmater.8b02236 Chem. Mater. 2018, 30, 5404−5411