Activating Nitrogen for Electrochemical Ammonia Synthesis via an Electrified Transition-Metal Dichalcogenide Catalyst

The complex interplay between local chemistry, the solvent microenvironment, and electrified interfaces frequently present in electrocatalytic reactions has motivated the development of quantum chemical methods that can accurately model these effects. Here, we predict the thermodynamics of the nitrogen reduction reaction (NRR) at sulfur vacancies in 1T′-phase MoS2 and highlight how the realistic treatment of potential within grand canonical density functional theory (GC-DFT) seamlessly captures the multiple competing effects of applied potential on a catalyst interface interacting with solvated molecules. In the canonical approach, the computational hydrogen electrode is widely used and predicts that adsorbed N2 structure properties are potential-independent. In contrast, GC-DFT calculations show that reductive potentials activate N2 toward electroreduction by controlling its back-bonding strength and lengthening the N–N triple bond while decreasing its bond order. Similar trends are observed for another classic back-bonding ligand in CO, suggesting that this mechanism may be broadly relevant to other electrochemistries involving back-bonded adsorbates. Furthermore, reductive potentials are required to make the subsequent N2 hydrogenation steps favorable but simultaneously destabilizes the N2 adsorbed structure resulting in a trade-off between the favorability of N2 adsorption and the subsequent reaction steps. We show that GC-DFT facilitates modeling all these phenomena and that together they can have important implications in predicting electrocatalyst selectivity for the NRR and potentially other reactions.


S1 Structural relaxation and adsorption site sampling S1.1 Vacancy Relaxation
The 4x4 supercell of 1T-MoS 2 was fully relaxed upon the introduction of a sulfur vacancy, which caused Mo-atom dimerization along one in-plane lattice vector as expected for a 1T to 1T ′ -phase change and is illustrated in Figure S1.

S1.2 Selection of Relaxation Radius
A relaxation radius (radius from center of vacancy in wich atoms were allowed to relax) was selected after observing that if the adsorbed intermediate structures were allowed to fully relax, unphysical distortions of the MoS 2 lattice would occur causing unrealistic deviations in energy along the pathway.Examples of these distortions are illustrated for two selected intermediates in Figure S2.lattice when all atoms are allowed to relax for two selected intermediates.
To verify the effect of the size of the relaxation radius on the adsorption energy, we calculated the N 2 adsorption energy for different sized radii.For N 2 adsorption, we show that the adsorption energy fluctuates within 0.04 eV and the main effect of the relaxation radius is to prevent lattice distortions.

S7 Density of States and Crystal Orbital Hamilton Population analysis
In order to gain a deeper understanding of the orbital interactions involved in bonding, we analyzed the orbital projected density of states (pDOS) of the 1T ′ -MoS 2 -Sv structure, the N 2 , CO, and H at Mo adsorption structures shown in main text Figure 6, as well as free N 2 and CO.We also performed chemical-bonding analysis to obtain the projected Crystal Orbital Hamilton Population (pCOHP) 1,2 for the N 2 bound structure.The pDOS calculations were performed in JDFTx 3 using a Gaussian broadening (E σ ) of 0.001 Ha and additional smoothing was applied using a Gaussian filter (σ = 5).For all pDOS figures, the Fermi level was shifted to 0 eV.The pCOHP curves were obtained using the Local Orbital Basis Suite Towards Electronic-Structure Reconstruction (LOBSTER) v 4.1.0package. 4,5LOBSTER enables projection of plane wave basis sets onto localized orbitals and thus can give the bonding and antibonding contributions to the pDOS from adjacent atoms.Since LOBSTER calculations cannot currently be run on JDFTx outputs, we ran single point (fixed geometry) calculations using the Vienna ab initio simulation package (VASP) [6][7][8] on the optimized structures and number of electrons from JDFTx to obtain inputs for LOBSTER.Calculation settings in VASP were chosen to closely match our JDFTx calculations considering the use of projected augemented wave (PAW) pseudopotentials. 9We performed spinpolarized density functional theory calculations using the PAW method, cutoff energy of 600 eV, and Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional.A 3 × 3 × 1 k-point mesh was used with convergence threshold set to 1 × 10 −5 eV.Solvation effects were treated by using the implicit solvation model implemented in VASPsol. 10,11he pDOS plots for the neat 1T ′ -MoS 2 -Sv and 1T ′ -MoS 2 -Sv+H@Mo are shown in Figure S13 while the pDOS plots for 1T ′ -MoS 2 -Sv+N 2 and 1T ′ -MoS 2 -Sv+CO and their corresponding free molecules are shown in Figure S14.We find that the molecular 2π* orbital of free N 2 broadens and hybridizes with the Mo d-states due to adsorption.The adsorbed N 2 2σ*, 2π, and 3σ orbitals also hybridize with Mo d-states, but still remain relatively sharp.Similar electronic structure observations for N 2 adsorption to single atom catalysts have been made. 12,13In contrast, this interaction is not observed for 1T ′ -MoS 2 -Sv+H@Mo, while no significant change is observed in the Mo 4d and molecular 2p state overlap in the pDOS of CO bound to 1T ′ -MoS 2 -Sv, which remains strongly adsorbed across potentials.To understand bonding within the 1T ′ -MoS 2 -Sv+N 2 , we also performed COHP analysis of the Mo-*N A and *N A -N B interactions in the 1T ′ -MoS 2 -Sv+N 2 structure (Figure S15).The COHP analysis shows that the region in the DOS where the 2σ*, 2π, and 3σ orbitals N 2p states appear below the Fermi level corresponds to bonding interactions between Mo and N while the broad states above the Fermi level have primarily antibonding character.Additionally, the 3σ region flips to antibonding when the interaction between the nitrogen atoms are considered, indicating the electron density is likely being donated to the Mo in this bond.In the highlighted region of the DOS in S14b, we note the overlap between the broad Mo 4d and the N 3σ states correlates with the favorable to unfavorable adsorption of nitrogen with potential.That is, there is more overlap in these states when adsorption is most favorable (under no applied potential), and less overlap when adsorption is least favorable (-0.8 V applied potential).Finally, since GGA level DFT tends to underestimate energy gaps, we validated these findings by running a single point hybrid calculation on each of the PBE structures using PBE0.Despite the increase in spread of the density of states, the overlap trends for the Mo 4d and the N 3σ remained the same as shown in Figure S15.We conjecture that the large spread in energy of Mo 4d states is key to their ability to interact with bonding and antibonding states of backbonding molecules well above and below the Fermi energy.

Figure S1 :
Figure S1: Top down view of 4x4 supercells of MoS 2 showing structural relaxations after introduction of a sulfur vacancy (sulfur atoms are yellow and molybdenum atoms are purple).(a) 4x4 supercell before a vacancy is introduced in the 1T-phase (b) 4x4 supercell after a vacancy is introduced showing dimerization along the b-axis.(c) Overlayed images of 1T-MoS 2 with (Mo-atoms in darker purple) and without a vacancy (Mo-atoms in lighter purple) with arrows showing the direction of some Mo atoms upon relaxation.

Figure S2 :
Figure S2: Top down view of 4x4 supercells of MoS 2 showing structural distortions of MoS 2 lattice when all atoms are allowed to relax for two selected intermediates.

FigureFigure S5 :Figure S7 :Figure S8 :
Figure S4: Mo-N bond length in the 1T ′ -MoS 2 -Sv+N 2 structure under solvation with no potential in the canonical ensemble (NP) and applied potentials using GC-DFT

Figure S12 :
Figure S12: CO adsorption structures in top down view (a) and side on view (b,c) where (c) shows the charge density difference at the ±0.0001 e/Å 3 isosurface level for a −0.2 V step towards negative potential (0.0 to −0.2 V shown, however, similar across all potential steps).Red indicates regions of charge loss, blue indicates region of charge gain.

Figure
Figure S15: -pCOHP curves for Mo-N A and N A -N B atom interactions showing bonding (positive) and antibonding (negative) states.N A represents the nitrogen bound to the Mo site and N B the outer nitrogen as defined in the main text.

Table S1 :
N 2 adsorption energies for different relaxation radius sizes.

S1.3 Adsorption site sampling 6 7 8 10 11 Figure
S3: Top down view of 1T ′ -MoS 2 -Sv structure.Red circles indicate N 2 binding sites that were sampled.E ads values given in TableS2with a note on final adsorption geometries.

Table S2 :
Table of initial N 2 adsorption geometries explored where N