On the chemistry and mobility of hydrogen in the interstitial space of layered crystals h -BN, MoS 2 , and graphite

: Recent experiments have demonstrated transport and separation of hydrogen isotopes in layered materials, such as hexagonal boron nitride and molybdenum disulphide. Here, based on first-principles calculations combined with well-tempered metadynamics simulations, we report the chemical interactions and mobility of protons (H + ) and protium (H atoms) in the interstitial space of these layered materials. We show that both H + and H can be transported between the layers of h -BN and MoS 2 with low free energy barriers, while they are immobilized in graphite, in accordance with the experimental observations. In h -BN and MoS 2 , the transport mechanism involves a hopping process between the adjacent layers, which is assisted by the low-energy phonon shear modes. Defects present in MoS 2 suppress the transport and act as traps for H species. , Methods : Well-tempered metadynamics 4, 5 (WTMetaD) simulations were performed using 3.0

ABSTRACT: Recent experiments have demonstrated transport and separation of hydrogen isotopes in layered materials, such as hexagonal boron nitride and molybdenum disulphide. Here, based on first-principles calculations combined with well-tempered metadynamics simulations, we report the chemical interactions and mobility of protons (H + ) and protium (H atoms) in the interstitial space of these layered materials. We show that both H + and H can be transported between the layers of h-BN and MoS 2 with low free energy barriers, while they are immobilized in graphite, in accordance with the experimental observations. In h-BN and MoS 2 , the transport mechanism involves a hopping process between the adjacent layers, which is assisted by the lowenergy phonon shear modes. Defects present in MoS 2 suppress the transport and act as traps for H species.
Introduction. Recently, Hu et al. 1 reported hydrogen isotope separation by sieving through the interstitial space of layered materials, namely, hexagonal boron nitride (h-BN) and molybdenum disulfide (MoS 2 ). While the difference in the entry barrier of the hydrogen isotopes into the interstitial space of the layered materials was found to be the driving factor for the isotope separation, there is no conclusive picture yet on the chemical interaction of interstitial hydrogen with the layered materials, and on the transport mechanism through them. Spatial confinement has a strong impact on the chemical properties of molecules. [2][3][4] The layers of h-BN, MoS 2 , and graphene are chemically very stable, if not inert, and the interlayer interactions impose pressure on the intercalated species. [5][6][7] So, the chemistry in the interstitial space of a layered material is expected to be significantly different from the surface-adsorbed counterpart.
The experimental setup in the work of Hu et al. 1 does not allow to conclude if protons (H + ) or protiums (H atoms) enter the layered material: For protons, immediate neutralization is expected, as the proton's electron affinity of 13.6 eV is significantly higher than the work function of the host materials, even for wide band gap insulators such as h-BN. On the other hand, any two H atoms encountering each other are expected to recombine immediately due to the high formation energy of H 2 (4.52 eV in gas phase). 8 The observed transport process of hydrogen through h-BN and MoS 2 requires, therefore, further investigation.
Here, we explore, on grounds of density-functional theory (DFT -for details see Methods section and Supporting Information (SI)), the chemical interaction and mobility of both H and H + species in the layered materials of h-BN, MoS 2 , and graphite. We show that both H + and H can be accommodated in h-BN and MoS 2 , and that diffusion is possible via a low-energy transport mechanism that is assisted by the rigid-layer shear modes, 9 which have low frequency due to weak interlayer-interactions 10 (these modes have been observed in the Raman spectra at 52.5 cm -1 and 33.7 cm -1 for h-BN 9 and MoS 2 , 11 respectively). For both H and H + transport follows a zigzag hopping motion between adjacent layers. Diffusion in MoS 2 has a lower barrier than in h-BN, however, sulfur vacancies -a typical defect in MoS 2 12 -act as hydrogen traps and suppress the transport. In graphite, protons are immediately neutralized to H atoms, which themselves are almost immobile once bound to the carbon atoms. 13,14 Results and Discussions

Structure and formation of proton-and protium-intercalated h-BN, MoS 2 , and graphite .
We first study the static, fully optimized structures with single proton (H + ) or protium (H atom) confined in the interstitial space of layered crystals of h-BN, MoS 2 , and graphite. Details of the simulations are shown in Methods section and SI (Figures S1-S5). A proton located in the interstitial space of layered h-BN binds to a N atom, forming a N-H bond with a bond length of 1.05 Å, and locally slightly decreases interlayer separation. On the contrary, a H atom forms a B-H bond of 1.32 Å length and slightly increased interlayer separation. Details are shown in Figure   1 and an electron (~0.9 eboth in h-BN and MoS 2 ) and is bound to the crystal lattice. This charge is transferred mainly from the atoms that are in closest vicinity. Thus, the proton itself is almost neutralized, while the N or S atoms, to which the proton binds, become positively charged. Also, the second neighbors become slightly more positive, but the effect is much weaker than for the binding atom. This effect is also shown by the electron density isosurfaces (see Figure S6 in SI).
Even though the protons get immediately neutralized by the surrounding atoms of the crystal, in the remainder, we will still use the term proton in order to acknowledge the local charge of the system. On the other hand, no charge transfer between the h-BN, MoS 2 or graphite layers and the bound H atom is observed.
Diffusion process of proton and protium inside h-BN, MoS 2 , and graphite .
Next, we have simulated the diffusion process of the protons and H atoms through the interstitial space between the layers of 2D materials, again focusing first on h-BN. As both H and H + bind to atoms of the crystal lattice, we investigated the transport process between individual binding sites. As the system has many degrees of freedom, including low-energy vibrational modes, with rather large amplitudes, we employed well-tempered metadynamics 15 (WTMetaD) to determine the free-energy surface (FES) and to analyze the diffusion barriers of H and H + travelling through the interstitial space of the layered materials (for details of WTMetD and DFT calculations see SI). The FES and the diffusion barriers for H + and H atom are discussed below in detail. Note that WTMetaD diffusion is governed by Brownian motion, thus, the directionality of the diffusion is not defined. However, below, we will refer to this as "path". The information that is the most important from these simulations is the lowest free energy barrier for the transfer of hydrogen species in between the layers. At the end, we will also comment on the directionality of the diffusion for both species, assuming the conditions present in the experiment.
We BN, exhibit low-energy shearing modes. 9, 10 Once the shearing brings a N atom of the adjacent layer to the vicinity of the proton, there is a high probability that proton jumps to the closest N atom (N2, red) in the adjacent layer. The next jump (N2 to N3, green) then goes back to the first layer, and so forth, resulting in a zigzag transport path between layers. Thus, after each two successive jumps, the proton can be displaced by 2.5 Å, or the shortest intralayer N-N distance.
The diffusion process of H atom follows a similar shear-assisted process, with the difference that the binding sites are now the B atoms, see Figure 2 c and d, and that the FES is shallower due to the lower binding strength of H atom to the lattice, when compared to protons. The trajectory for this process is shown in Figure S3 in SI. It is important to note that proton carries its charge during the hopping process: even though the proton is nearly neutralized while being bound to N atom, when dislocating to the next N adsorption site, it leaves the electron behind. The atomic charge analysis shows that the electronic configuration of the layer recovers to the pristine state after the proton has left.
Close inspection of the WTMetD pathway suggests that the shear mode, where the individual layers slide against each other, significantly shortens the distance between the two next lowenergy hopping sites (the next N or B atom in the neighboring layer), and thus, has a strong effect on the hopping barrier (and even more so on the tunneling probability, which is neglected in this work), as depicted in Figure S4. Indeed, the barrier for the interlayer shear in h-BN has been found to be rather small, as low as 8 meV and 2 meV, for the shearing between the most stable AA' stacking to the A'B (N over N) and to AB' (B over B), respectively. 16,17 We estimate that the transition barrier between the layers in static calculations, as shown in Figure S4, substantially lowers the hopping energy barrier by ~0.5 eV, which is expected to lead to a very small free energy barrier of the overall motion. In conclusion, our simulations suggest that the transport of protons, and also of H atoms, through the interstitial space, follows a shear-mode assisted zigzag hopping pathway. The hypothesis of a kinetic barrier was tested by placing two H atoms into the simulation box.
Three scenarios, as shown in Figure 3,   Table S2 in SI.
We now turn to the case of MoS 2 layered crystals, for which we performed the same type of WTMetaD simulations of H + and H inside the interstitials. Both protium and proton inside MoS 2 bind to S atoms with very similar binding strengths (see SI for details). The transport mechanism is very similar to that in h-BN: Initially, H + or H bind to sulfur atom S1 (see Figure 4a; black).
As in h-BN, the diffusion is assisted by the shear modes, which allow the transfer from S1 to S2  Finally, we have investigated H atom inside graphite (proton gets immediately neutralized inside this system). Here, however, we have not observed any H atom transfers with free energy barriers similar to the other materials, meaning that this process would require much higher activation energies and would be improbable in the experimental conditions (see Figure S9).
This finding is consistent with the experimental results, 1,18 which show that graphite does not conduct hydrogen.
In conclusion, both H + and H atom can diffuse easily through the interstitial space of layered materials h-BN and MoS 2 , while they are immobile in graphite. The overall lattice of the layered materials remains intact. Locally, protons bind to N (S for MoS 2 ), pick an electron from the lattices and leave a charged site that extends over a few atoms. At low concentration, this charged state will be maintained. Protons could be transported between the layers of such 2D materials due to the Coulomb repulsion between the protons and the electric field gradient between the electrodes. On the other hand, protiums will bind to B (S for MoS 2 ) atoms. If two protiums are getting close, they will form proton-hydride pairs, and at high concentration, they will recombine to H 2 . Thus, for low concentration, both pathways are in principle possible. In both carrier materials, the transport process is supported by the shear modes of the material and follows a hopping process between neighboring layers. Our results should be of interest for hydrogen transport in other layered van der Waals materials, which could be important for applications, such as hydrogen isotope separation or proton exchange membranes.   Figure   S1). For all the systems, the interaction of proton (H + ) and hydrogen atom (H) with the layers was investigated. We used the following supercells to describe the interactions: the 4×4×2 supercell was used for both h-BN and graphite, while 4×4×1 supercell for MoS 2 , resulting in 128 atoms in h-BN and graphite, and 96 atoms in MoS 2 . The fully optimized lattice parameters of each system are shown in Table S1 and are in a fairly good agreement with the experimental values. 1-3 Figure S1. Top and side views of layered materials considered in the present work. Table S1. Calculated lattice parameters (a, b, and c) and the interlayer distances (d) in perfect h-BN, MoS 2 , and graphite, shown with and without H species bound to the layers. Available experimental data are given in parenthesis. [1][2][3] Bond distances between H + or H atom and the binding atoms are also given.  For the CVs defined in each system, the following reference H-X (X = N, C, B, or S) distances (R 0 ) were used: 2.3 bohr for H + @hBN, 2.5 bohr for H@graphite, 2.9 bohr for H@hBN, 3.0 bohr for H@MoS 2 and H + @ MoS 2 . These parameters are slighltly larger (up to 10%) than the average value in the fluctuation of the corresponding H-X distances during the BOMD simulations. Thus, when the H or H + is close to the transition state region, the coordination number of H to X is close to 0. Note that in the present models, WTMetaD setup, and the definition of the CVs, the hopping between the layers can be discussed, however, the direction of the transport is not defined and is governed by Brownian motion. To understand the latter, more CVs would need to be defined, e.g., the so-called path CVs, in several different directions, in order to sample the whole free-energy landscape (this is beyond the scope of this work). We have also investigated the case with three CVs specified, where the third one describes the hopping of H species within the same layer. However, such a process was not observed with the employed number of hills, thus in all cases, the free energy barrier for transport within the same layer was much higher than between the layers. Thus, it is a safe choice to reduce the number of CVs to just the two discussed above.

Proton transport in defective MoS 2 : The WTMetaD simulations of H + and H atom in a defective
MoS 2 system (with one sulfur atom vacancy) is shown in Figure S7. The transfer barriers in the case of proton are still similar to the case of a perfect material, however, the protons are strongly attracted by the S atoms in the vicinity of the vacancy, hindering the movement of the species.
An even stronger effect is observed for the H atom, when different free-energy barriers are found. This comes from the fact that the S-H binding strength is different when H species are close or far from the vacancy (see Figure S8). The energy difference between H + bound close to the vacancy or far from it is only 2.4 kJ/mol, while for H atom, this energy difference is as large as 12.1 kJ/mol, leading to the asymmetric FES. We can conclude here that the defects present in MoS 2 will suppress the overall transport between the layers.

Spin flip energy:
When two H atoms are in the interstitial space of a layered material, different spin states can be considered: singlet, triplet or antiferromagnetic singlet. The spin-flip energy for the latter two spin states of two H atoms were calculated with Crystal 17 at the PBE/TZVP level (see Table S2). These energies are very small, therefore, results obtained for the triplet or antiferromagnetic singlet are nearly the same.