Deep Learning-Assisted Investigation of Electric Field–Dipole Effects on Catalytic Ammonia Synthesis

External electric fields can modify binding energies of reactive surface species and enhance catalytic performance of heterogeneously catalyzed reactions. In this work, we used density functional theory (DFT) calculations—assisted and accelerated by a deep learning algorithm—to investigate the extent to which ruthenium-catalyzed ammonia synthesis would benefit from application of such external electric fields. This strategy allows us to determine which electronic properties control a molecule’s degree of interaction with external electric fields. Our results show that (1) field-dependent adsorption/reaction energies are closely correlated to the dipole moments of intermediates over the surface, (2) a positive field promotes ammonia synthesis by lowering the overall energetics and decreasing the activation barriers of the potential rate-limiting steps (e.g., NH2 hydrogenation) over Ru, (3) a positive field (>0.6 V/Å) favors the reaction mechanism by avoiding kinetically unfavorable N≡N bond dissociation over Ru(1013), and (4) local adsorption environments (i.e., dipole moments of the intermediates in the gas phase, surface defects, and surface coverage of intermediates) influence the resulting surface adsorbates’ dipole moments and further modify field-dependent reaction energetics. The deep learning algorithm developed here accelerates field-dependent energy predictions with acceptable accuracies by five orders of magnitudes compared to DFT alone and has the capacity of transferability, which can predict field-dependent energetics of other catalytic surfaces with high-quality performance using little training data.


Computational Setup
DFT calculations were conducted using the Vienna Ab-initio Simulation Package (VASP). These calculations were used to determine the ground state energy and the most favorable configuration of each species in the process of ammonia synthesis over Ru. There were two Ru surfaces tested, including the most thermodynamic favorable flat surface -Ru(0001), and the most active stepped surface for ammonia synthesis -Ru(1013). 1 In particular, the Ru(1013) surface is equivalent to a periodic "ladder" with B5 type adsorption sites ( Figure S1). 2 The Ru(0001) surface has four adsorption sites for which species can adsorb, including top, bridge, hcp, and fcc. The Ru (1013) surface has seven adsorption sites for the reaction-related species to adsorb, including hollow site For all calculations we used the Generalized Gradient Approximation with the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional. [3][4][5] The projector-augmented wave method was applied to solve the Kohn-Sham equations. 6 A first-order Methfessel-Paxton smearing method with width of 0.1 eV was used to account for partial occupancies around the Fermi level and to speed up the convergence. 7 We chose the k-point Monkhorst-Pack mesh (3 × 3 × 1) for p(4 × 4) Ru(0001) surface and Ru(1013) surface and the plane-wave energy cutoff (500 eV) for DFT calculations. For four-layer Ru(0001) and Ru(1013) slabs, we fixed the bottom two layers to their bulk position and relaxed the top two layers to electronically interact with the surface adsorbates.
The lattice constant of Ru was calculated as 2.726 Å, which is consistent as experimental work and other theoretical works. [8][9] Geometries were considered optimized when the energy had converged to 10 -5 eV and the forces were smaller than 0.03 eV/Å.
We have examined the energetics of the possible reaction paths for ammonia synthesis ( Figure S2) from dinitrogen reduction over Ru(0001) and Ru(1013) surfaces. Based on the thermodynamic energetics, the 1 st and 2 nd most favorable reaction paths are presented in Figure   S3. Thus, in our work, we only examine the electric field effects on the 1 st and 2 nd most favorable reaction paths.
To simulate the field-induced energetics through DFT calculations, we use the approach as proposed by Neugebauer and Scheffler, 10 wherein the uniform electric fields are generated without adding or removing any charge from the supercell. 11 Regarding the 'field emissions' effects, 12 our vacuum width is ~11 Å. This leads to a charge density in the vacuum, which is less than the magnitude of the Gibbs' oscillations (0.001 e/Å 3 ) associated with plane-wave cutoffs.
We first examined the possible adsorption configurations for the reaction intermediates during ammonia synthesis over Ru surfaces (Figures S4 -S19). To examine the electric field effects on the energetics (i.e., adsorption, reaction) of ammonia synthesis, we applied the negative electric fields from 0 to -1 V/Å with the interval of 0.2 V/Å and applied the positive electric fields from 0 to 1 V/Å with the interval of 0.2 V/Å to the most favorable adsorption geometries and optimized them. The adsorption energy ( !" ( )), in the presence and absence of electric fields, could also be calculated according to Eqn. (S1): where #$#!% ( ) represents the total energy of the species over the surface with fields, &'()!*+ ( ) represents the energy of the surface (either Ru(0001) or Ru(1013)) with fields, and !"&$(,!#+ represents the energy of the adsorbate in the absence of electric fields. Here, our assumption is that when an external electric field is present, only the energetics of adsorbates after adsorption and catalyst surfaces will be influenced by the fields. While the energetics of the gas phase molecules won't be influenced by the electric field effects due to the full degree of freedoms (rotation, translation, vibration) of the molecules in the gas phase. More negative value of the !" ( ) represents strong bonding between the adsorbates and the surface. The field-dependent adsorption energies of the ammonia synthesis-related reaction intermediates are presented in We then examined the co-adsorption configurations of the intermediates for the 1 st and 2 nd most favorable reaction mechanism during ammonia synthesis over Ru(0001) and Ru (1013) surfaces (Figures S22 -S30). The co-adsorption configuration setup of two species was based on the most favorable stable sites of each involved species, found previously in Figures S4 -S19.
The field-dependent reaction energy (∆ (-. ( )) was calculated using Eqn. (S2): where /0 ( ) and 10 ( ) are the total energies of the final state (FS) and the initial state (IS) of one elementary step during the ammonia synthesis process with and without electric fields.
Endothermic reactions are accompanied by positive values of ∆ (-. ( ) and vice versa.
Asterisk ( * ) represents the empty surface site.     The most favorable adsorption geometry is when N2 * sits at the hcp site on the Ru(0001) surface. The most favorable adsorption geometry is when N2 * sits at the S4 site on the Ru(1013) surface.   The most favorable adsorption geometry is when NH * sits at the hcp site on the Ru(0001) surface. When NH * adsorbs over the bridge site, it is not stable and it shifts to the nearby hcp site. The most favorable adsorption geometry is when NH * sits at the h2 site on the Ru(1013) surface. Figure S14. The possible adsorption configurations of NH2 * over the Ru(0001) surface and its coresponding energetics. The most favorable adsorption geometry is when NH2 * sits at the bridge site on the Ru(0001) surface.
When NH2 * sits at the fcc or hcp site, it is not stable and it shifts to the nearby bridge site.

Large Electric Field in Experimental Setups
In the research field of electric field-enhanced catalysis, to the best of our knowledge, there are three ways to generate a large electric field in the experimental works.
• Highly localized electric field within the catalyst structure. Large oriented electric fields can be created over molecular length scales in (metallo-)enzyme [26][27][28][29][30][31][32] and zeolite catalytic 33 active sites, on the order of >0.1 V/Å. These fields may decrease activation barriers in enzymatic and zeolite catalysis and in turn, govern the selectivity of the desired reaction path. [34][35] Especially, the Boxer group 27, 36 has identified the highly localized electric field,

Field-Dependent Reaction Barriers
The elementary reactions over Ru(1013) to calculate the field-dependent activation barriers are shown below: where the asterisk (⁎) represents the surface empty site. All field-dependent activation barrier results from DFT are shown in Figure S31. It is noticed that H2 * over Ru(1013) doesn't exist. This is because when H2 adsorbs over the surface, it automatically breaks into 2H * with the electric field of -1 V/Å, while H2 dissociation with the electric field of 0 and 1 V/Å is also a barrier-free process ( Figure S31a). The convergence criteria of TS are the same as that of the adsorption geometry optimization, i.e., the energy converged to 10 -5 eV and the forces smaller than 0.03 eV/Å.
For the elementary reaction of NH2 * + H * ↔ NH3 * + ⁎, we used CI-NEB method to get a better initial guess for the minimum energy pathway and used Dimer 79 method to find the true TS. The true TS was further verified by calculating their vibrational frequencies confirming that a single unique normal mode eigenvector was obtained which corresponds to the negative curvature at the saddle point. We then corrected the field-dependent activation barrier (∆ ǂ ( )) and field-dependent reaction energy ( ∆ 567 ( ) ) with zero-point energy (ZPE), entropy (S) and temperature (T) corrections as follows: Thus, only field-dependent vibrational entropy contribution ( =>? ( ) ) is considered for the temperature-entropy correction for surface reactions (Eqn. (S16)). The entropy for gas phase molecule (i.e., NH3) at a certain temperature (T) and standard pressure (0.1 MPa) can be referred to the JANAF thermochemical table. 80 We also fixed all Ru atoms and relaxed only the surface adsorbates to avoid vibrational frequency noise form the surface for vibrational frequency calculations. It is worth mentioning that the vibrational frequencies ( @ ( )) for the surface species were calculated with and without electric fields. : where G is Boltzmann constant and R is gas constant. We presented the ZPE and TS corrected field-dependent activation barriers in Table S1 and field-dependent reaction energies in Table S2.
Here, temperature conditions of 473 K 81 and 773 K 52 are to estimate low temperature and high temperature field-accelerated energetics of ammonia synthesis.

Charge Analysis
Differential charge density analysis and Bader charge analysis for various coverages of NH3 over Ru(0001) surfaces are shown in Figure S32. Figure S32. Differential Charge density analysis and Bader charge analysis for various coverages of NH3 over Ru(0001) surfaces. The isosurface of the differential charge density is 0.01 e/Bohr 3 .

Deep Learning Method
Let H denote the features and H stand for the similarities, where = ( , , ) is the set of adsorbates in the gas phase, catalytic surfaces, and adsorbate adsorption over the catalytic surfaces, respectively. Let H be the electric fields, and H be the field-dependent total energy, then our goal is to find a function : (  We designed a novel deep learning framework for predicting the field-dependent adsorption energies. Specially, we employed the Graph Neural Network (GNN) to capture the relationship among the geometries of adsorbate in the gas phase, pure catalytic surface, and adsorbate adsorption over the catalytic surface. The reason of applying GNN is a special deep learning architecture widely used to capture the information-exchange among the relational data.
Then we adopted Multi-Layer Perceptron (MLP) for regression analysis. Let = { J , K , ; } be the learnable parameters in the model, where J , K , ; are learnable parameters for GNN, hidden layers and prediction layers, respectively. Detailed objective function is as follows.
Graph Neural Network. We built a graph for each structure, where each node denotes an atom, and each edge is represented by the similarity between two atoms. The adjacency matrix of the graph is defined by: where 3 ( , ) is the squared distance between atom and atom . With this adjacency matrix, we then applied Graph Convolutional Network (GCN), 82 a concrete algorithm of GNN, to incorporate similarities of all atoms with their corresponding features, as follows: where V = + is the adjacency matrix with added self-connections, is the identity matrix, Then for each of the three parts, i.e., adsorbate in the gas phase, pure catalytic surface, and adsorbate adsorption over the catalytic surface, we got a hidden representation. Combining them with the energy of adsorbate and surface, we adopted another MLP as the prediction layer: where $ is the predicted the energy for adsorbate adsorption over the catalytic surface. Finally, we adopted Equation (S1) to calculate the adsorption energy of reactive intermediates under various electric fields.
Objective Function. We adopted absolute error as our objective function, which can be written as: We adopted Adam optimizer 84 to minimize the objective function.

Implementation Details
We adopted Adam optimizer 84 with a learning rate of 0.003, and run 8000 epochs to train the model. The experiments were run on a physical machine with Ubuntu 18.04, total memory of 64GB, AMD Ryzen Threadripper 2920X 12-Core Processor, and an NVIDIA GP102 GPU.

Training Structures
Our deep learning model for ammonia synthesis over Ru(1013) and Ru(0001) surfaces has three components, including the optimized geometries and the corresponding energies of (1) the radicals/molecules in the gas phase (19 species)  adsorption sites. We also include some species (e.g., NH, N-NH co-adsorption, etc.) at their less favorable adsorption sites. Considering that there is no public train/test split for this dataset, we adopted 10-fold cross validation, which splits all structures into 10 folds and perform the fitting procedure a total of ten times. Each fit is performed on a training set consisting of 9 folds of structures and the remaining fold as the test set. The training structures are shown in Tables   S3 and S4.  To broaden the impacts of our developed deep learning model, we show here that our deep learning model has the transferability to predict field-dependent adsorption energies, dipole moments, and polarizabilities of ammonia synthesis-related intermediates over other catalytic surfaces (e.g., Fe). The possible adsorption sites over Fe(110) and Fe(111) are shown in Figure   S34. We applied the Ru-trained deep learning model to predict the energetics of ammonia synthesis over Fe(110) and Fe (111)     H-H (hol-hol) NH-H (lb-hol) NH2 (lb) Note: hol represents hollow site and bri represents the bridge site. Table S6. The 20% training structures over Fe(110) and Fe (111) with no electric field and ± 0.6 V/Å.