Reaction Mechanism with Thermodynamic Structural Screening for

Jul 8, 2018 - In this work, DFT calculations were carried out to obtain the Pourbaix ... (39−41) A cutoff energy of 500 eV and a two-dimensional ...
1 downloads 0 Views 2MB Size
Subscriber access provided by - Access paid by the | UCSB Libraries

Article

Reaction Mechanism with Thermodynamic Structural Screening for Electrochemical Hydrogen Evolution on Monolayer 1T' phase MoS 2

Shiqi Chen, Xiaobo Chen, Guangjin Wang, Lu Liu, Qiaoqiao He, Xi-Bo Li, and Ni Cui Chem. Mater., Just Accepted Manuscript • DOI: 10.1021/acs.chemmater.8b02236 • Publication Date (Web): 08 Jul 2018 Downloaded from http://pubs.acs.org on July 8, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Chemistry of Materials

Reaction Mechanism with Thermodynamic Structural Screening for Electrochemical Hydrogen Evolution on Monolayer 1T’ phase MoS2 Shiqi Chen1, Xiaobo Chen*,1, Guangjin Wang2,3, Lu Liu1, Qiaoqiao He1, Xi-Bo Li1, and Ni Cui1 1

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. 2 Center for Functional Nanomaterials, Brookhaven National Laboratory, Upton, New York 11973, United States. 3 College of Chemistry and Materials Science, Hubei Engineering University, Xiaogan 432000, China.

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 diffusion-assisted Tafel manner, rather than the Heyrovsky manner assumed before. 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 towards HER should be done on the premise that the lowest chemical potential of protons and electrons is obtained.

1. INTRODUCTION 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 attentions, 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-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).2529 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 (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 used structural models. Hydrogen adsorption models fail to consid-

1

ACS Paragon Plus Environment

Chemistry of Materials 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 10

rized 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 before.

er aqueous environment and interfacial electric field, which might result in 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 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 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 of 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 summa-

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 labelled 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.

2. MODELS AND METHODS 2.1. Model Systems. The reaction mechanism was studied by constructing a solid/liquid interfacial model, which contained a monolayer of 1T’-MoS2 covered by two layers of water, as shown in Figure 1 (a) and (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 be 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 configuration 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 Abinitio 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 twodimensional 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 method (CINEB)42 was used to locate transition states. Free energies were obtained by correcting zero-point energies and entropies, as elaborated in Sec. 2 of the Supporting Information (SI). 2.3. Interfacial Thermodynamics. The electrode potential USHE, which is set up by interfacial dipoles, is calculated by

2

ACS Paragon Plus Environment

Page 3 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Chemistry of Materials

referencing the system’s work function  to that of the Standard Hydrogen Electrode ( ).      / (1) 43  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 ( 

 ) in an interfacial system could be calculated with 

      0 

, , 



 



  e  2.3!"pH. (3) -eUSHE describes the contribution of USHE to the chemical potential of electrons, while -2.3kT·pH 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 

      0 ,32 



      

      0      . (4) 

 , (2)

By comparing 

     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, 

      0 was calculated first, and then 

     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.

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.  is the free energy of an H2 molecule at 1 atm and 298 K, while ,  and , 0 refer to the free energies of an interfacial structure with  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

Figure 2. Calculated 

     vs.  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  . The shaded area represents the typical reaction condition of HER. Different proton concentrations are differentiated by the color bar, while the characteristics on 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 co-

3

ACS Paragon Plus Environment

Chemistry of Materials 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 10

existence 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 (↓).

3. RESULTS AND DISCUSSION 3.1. Free Energies and Pourbaix Diagrams. Figure 2 (a) and (b) show the calculated free energies, 

     , 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 since they represent the lowest 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. Since 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 2 (d). When pH = 0, five intersections indicated by the dashed lines divide the potential range into 6 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 2 (c). The Pourbaix diagram exhibits 6 phase regions separated by 5 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 impact. In overall, lowering of USHE leads to increase of proton concentration and shift from hydroxyl-covered surfaces to proton-covered ones. This is consistent with the experimental observation that catalytic performance always enhances when the electrode potential decreases. For a given proton concentration, lowering of USHE induces proton transfer from the catalyst surface to the solution and re-orientation 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 the 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 gradually transforms to red right triangles when pH increases. 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 water environment to obtain the lowest electrochemical potential, the surface coverage is proved to be much higher than that anticipated before. 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 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 region IV and VI and thus it contains the interfacial types of both the 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 they 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 to be 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 Sec. 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. Since 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 

      0 for adsorbed protons on site I and II, which are 0.05 and 0 eV, respectively. This suggests that S-II sites are more favorable for proton adsorption. Figure 3 (a) shows the minimum energy path (MEP) of the Volmer reaction at USHE = 0.15 V. The structural diagrams on 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

4

ACS Paragon Plus Environment

Page 5 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Chemistry of Materials

cate the work function change along the path. (b) Linearly fitted Volmer reaction barriers for different potentials.

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 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 3 (b). The barrier increases with the electrode potential in 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 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 4 (a). Since 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 4 (b) and (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. Sec. 4 of SI for the detail), the room-temperature free energy barriers (∆GRT) is determined to be 0.95 eV. Obviously, this twostep process is more favorable than the direct recombination manner. Particularly, the moderate barrier of 0.95 eV demonstrates 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 down 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 potential close to the equilibrium are not considered since 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 5 (a). The corrected barrier is as high as 1.94 eV, indicating that the Heyrovsky reaction on S-II sites is rather difficult. The reason could be ascribed to the strong proton adsorption on SII 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 5 (b) and (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 takes also a small barrier of 0.59

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 charge-extrapolation corrected result (red dots). The only imaginary frequency of the saddle point is shown. Blue crosses indi-

5

ACS Paragon Plus Environment

Chemistry of Materials 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

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.

Page 6 of 10

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).

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. Since no charge transfer across the interface, all the energies are not corrected for potential deviation. Imaginary frequencies are labelled nearby 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 S-II 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 since there is no charge transfer across the interface. Imaginary frequencies are labelled nearby the corresponding saddle points.

6

ACS Paragon Plus Environment

Page 7 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Chemistry of Materials

which corresponds to a TOF of ca. 0.02 s-1. Using the rate expression of the transition state theory,

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 protons desorb from the catalyst surface through the diffusionassisted Tafel mechanism, rather than the Heyrovsky mechanism assumed before. 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.

k =

 − ∆GRT k BT exp h  k BT

 (5)  

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 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 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 towards 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, since 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.

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

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

ASSOCIATED CONTENT Supporting Information. 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 bar-

7

ACS Paragon Plus Environment

Chemistry of Materials 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 10

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, 10971103. (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. (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, 49534961. (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, 56055609. (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.

riers. This material is available free of charge via the Internet at http://pubs.acs.org.

AUTHOR INFORMATION Corresponding Author *[email protected]

Notes The authors declare no competing financial interest.

ACKNOWLEDGMENT Mr. Xing Lu, the high-performance computing platform of Jinan University and the Guangzhou Key Research Center of 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, 21317050, 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. 11617330), and the Natural Science Foundation of Hubei Province (NO. 2018CFB669).

REFERENCES (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

8

ACS Paragon Plus Environment

Page 9 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Chemistry of Materials

(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 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 1996, 54, 11169-11186. (37) Blöchl, P. E., Projector augmented-wave method. Phys. Rev. B 1994, 50, 17953-17979. (38) Kresse, G.; Joubert, D., From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B 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. 2008, 103, 096102. (41) Klimeš, J.; Bowler, D. R.; Michaelides, A., Van der Waals density functionals applied to solids. Phys. Rev. B 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.

9

ACS Paragon Plus Environment

Chemistry of Materials 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Paragon Plus Environment

Page 10 of 10

10