Pourbaix Machine Learning Framework Identifies Acidic Water Oxidation Catalysts Exhibiting Suppressed Ruthenium Dissolution.

The demand for green hydrogen has raised concerns over the availability of iridium used in oxygen evolution reaction catalysts. We identify catalysts with the aid of a machine learning-aided computational pipeline trained on more than 36,000 mixed metal oxides. The pipeline accurately predicts Pourbaix decomposition energy (Gpbx) from unrelaxed structures with a mean absolute error of 77 meV per atom, enabling us to screen 2070 new metallic oxides with respect to their prospective stability under acidic conditions. The search identifies Ru0.6Cr0.2Ti0.2O2 as a candidate having the promise of increased durability: experimentally, we find that it provides an overpotential of 267 mV at 100 mA cm-2 and that it operates at this current density for over 200 h and exhibits a rate of overpotential increase of 25 µV h-1. Surface density functional theory calculations reveal that Ti increases metal-oxygen covalency, a potential route to increased stability, while Cr lowers the energy barrier of the HOO* formation rate-determining step, increasing activity compared to RuO2 and reducing overpotential by 40 mV at 100 mA cm-2 while maintaining stability. In situ X-ray absorption spectroscopy and ex situ ptychography-scanning transmission X-ray microscopy show the evolution of a metastable structure during the reaction, slowing Ru mass dissolution by 20× and suppressing lattice oxygen participation by >60% compared to RuO2.

Surface Segregation Studies in Ternary Noble Metal Alloys: Comparing DFT and Machine Learning with Experimental Data.

Surface segregation, whereby the surface composition of an alloy differs systematically from the bulk, has historically been hard to study, because it requires experimental and modeling methods that span alloy composition space. In this work, we study surface segregation in catalytically relevant noble and platinum-group metal alloys with a focus on three ternary systems: AgAuCu, AuCuPd, and CuPdPt. We develop a data set of 2478 fcc slabs with those compositions including all three low-index crystallographic orientations relaxed with Density Functional Theory using the PBEsol functional with D3 dispersion corrections. We fine-tune a machine learning model on this data and use the model in a series of 1800 Monte Carlo simulations spanning ternary composition space for each surface orientation and ternary chemical system. The results of these simulations are validated against prior experimental surface segregation data collected using composition spread alloy films for AgAuCu and AuCuPd. Our findings reveal that simulations conducted using the (110) orientation most closely match experimentally observed surface segregation trends, and while predicted trends qualitatively match observation, biases in the PBEsol functional limit numeric accuracy. This study advances understanding of surface segregation and the utility of computational studies and highlights the need for further improvements in simulation accuracy.

WhereWulff: A Semiautonomous Workflow for Systematic Catalyst Surface Reactivity under Reaction Conditions.

This paper introduces WhereWulff, a semiautonomous workflow for modeling the reactivity of catalyst surfaces. The workflow begins with a bulk optimization task that takes an initial bulk structure and returns the optimized bulk geometry and magnetic state, including stability under reaction conditions. The stable bulk structure is the input to a surface chemistry task that enumerates surfaces up to a user-specified maximum Miller index, computes relaxed surface energies for those surfaces, and then prioritizes those for subsequent adsorption energy calculations based on their contribution to the Wulff construction shape. The workflow handles computational resource constraints such as limited wall-time as well as automated job submission and analysis. We illustrate the workflow for oxygen evolution reaction (OER) intermediates on two double perovskites. WhereWulff nearly halved the number of Density Functional Theory (DFT) calculations from ∼240 to ∼132 by prioritizing terminations, up to a maximum Miller index of 1, based on surface stability. Additionally, it automatically handled the 180 additional resubmission jobs required to successfully converge 120+ atoms systems under a 48-h wall-time cluster constraint. There are four main use cases that we envision for WhereWulff: (1) as a first-principles source of truth to validate and update a closed-loop self-sustaining materials discovery pipeline, (2) as a data generation tool, (3) as an educational tool, allowing users (e.g., experimentalists) unfamiliar with OER modeling to probe materials they might be interested in before doing further in-domain analyses, (4) and finally, as a starting point for users to extend with reactions other than the OER, as part of a collaborative software community.

Applying Large Graph Neural Networks to Predict Transition Metal Complex Energies Using the tmQM_wB97MV Data Set.

Machine learning (ML) methods have shown promise for discovering novel catalysts but are often restricted to specific chemical domains. Generalizable ML models require large and diverse training data sets, which exist for heterogeneous catalysis but not for homogeneous catalysis. The tmQM data set, which contains properties of 86,665 transition metal complexes calculated at the TPSSh/def2-SVP level of density functional theory (DFT), provided a promising training data set for homogeneous catalyst systems. However, we find that ML models trained on tmQM consistently underpredict the energies of a chemically distinct subset of the data. To address this, we present the tmQM_wB97MV data set, which filters out several structures in tmQM found to be missing hydrogens and recomputes the energies of all other structures at the ωB97M-V/def2-SVPD level of DFT. ML models trained on tmQM_wB97MV show no pattern of consistently incorrect predictions and much lower errors than those trained on tmQM. The ML models tested on tmQM_wB97MV were, from best to worst, GemNet-T > PaiNN ≈ SpinConv > SchNet. Performance consistently improves when using only neutral structures instead of the entire data set. However, while models saturate with only neutral structures, more data continue to improve the models when including charged species, indicating the importance of accurately capturing a range of oxidation states in future data generation and model development. Furthermore, a fine-tuning approach in which weights were initialized from models trained on OC20 led to drastic improvements in model performance, indicating transferability between ML strategies of heterogeneous and homogeneous systems.

Semi-grand canonical Monte Carlo simulation of the acrolein induced surface segregation and aggregation of AgPd with machine learning surrogate models.

The single atom alloy of AgPd has been found to be a promising catalyst for the selective hydrogenation of acrolein. It is also known that the formation of Pd islands on the surface will greatly reduce the selectivity of the reaction. As a result, the surface segregation and aggregation of Pd on the AgPd surface under reaction conditions of selective hydrogenation of acrolein are of great interest. In this work, we lay out a workflow that can predict the surface segregation and aggregation of Pd on a FCC(111) AgPd surface with and without the presence of acrolein. We use machine learning surrogate models to predict the AgPd bulk energy, AgPd slab energy, and acrolein adsorption energy on AgPd slabs. Then, we use the semi-grand canonical Monte Carlo simulation to predict the surface segregation and aggregation under different bulk Pd concentrations. Under vacuum conditions, our method predicts that only trace amount of Pd will exist on the surface at Pd bulk concentrations less than 20%. However, with the presence of acrolein, Pd will start to aggregate as dimers on the surface at Pd bulk concentrations as low as 6.5%.

Machine-learning accelerated geometry optimization in molecular simulation.

Geometry optimization is an important part of both computational materials and surface science because it is the path to finding ground state atomic structures and reaction pathways. These properties are used in the estimation of thermodynamic and kinetic properties of molecular and crystal structures. This process is slow at the quantum level of theory because it involves an iterative calculation of forces using quantum chemical codes such as density functional theory (DFT), which are computationally expensive and which limit the speed of the optimization algorithms. It would be highly advantageous to accelerate this process because then one could do either the same amount of work in less time or more work in the same time. In this work, we provide a neural network (NN) ensemble based active learning method to accelerate the local geometry optimization for multiple configurations simultaneously. We illustrate the acceleration on several case studies including bare metal surfaces, surfaces with adsorbates, and nudged elastic band for two reactions. In all cases, the accelerated method requires fewer DFT calculations than the standard method. In addition, we provide an Atomic Simulation Environment (ASE)-optimizer Python package to make the usage of the NN ensemble active learning for geometry optimization easier.

Tuning oxide activity through modification of the crystal and electronic structure: from strain to potential polymorphs.

Discovering new materials with tailored chemical properties is vital for advancing key technologies in catalysis and energy conversion. One strategy is the modification of a material's crystal structure, and new methods allow for the synthesis and stabilization of potential materials in a range of crystal polymorph structures. We assess the potential reactivity of four metastable oxide polymorphs of MO2 (M = Ru, Rh, Pt, Ir) transition metal oxides. In spite of the similar local geometry and coordination between atoms in the metastable polymorphic and stable rutile structure, we find that polymorph reactivities cannot be explained by strain alone and offer tunable reactivity and increased stability. Atom-projected density of states reveals that the unique reactivity of polymorphs are caused by a redistribution of energy levels of the t2g-states. This structure-activity relationship is induced by slight distortions to the M-O bonds in polymorphic structures and is unattainable by strain. We predict columbite IrO2 to be more active than rutile IrO2 for oxygen evolution.

Relationships between the surface electronic and chemical properties of doped 4d and 5d late transition metal dioxides.

Density functional theory calculations were performed to elucidate the underlying physics describing the adsorption energies on doped late transition metal dioxide rutiles. Adsorption energies of atomic oxygen on doped rutiles M(D)-M(H)O2, where transition metal M(D) is doped into M(H)O2, were expressed in terms of a contribution from adsorption on the pure oxide of the dopant M(D) and perturbations to this adsorption energy caused by changing its neighboring metal cations and lattice parameters to that of the host oxide M(H)O2, which we call the ligand and strain effects, respectively. Our analysis of atom projected density of states revealed that the t2g-band center had the strongest correlation with adsorption energies. We show that charge transfer mediated shifts to the t2g-band center describe the ligand effect, and the radii of the atomic orbitals of metal cations can predict the magnitude and direction of this charge transfer. Strain produces systematic shifts to all features of the atom projected density of states, but correlations between the strain effect and the electronic structure were dependent on the chemical identity of the metal cation. The slope of these correlations can be related to the idealized d-band filling. This work elucidates the underlying physics describing adsorption on doped late transition metal oxides and establishes a foundation for models that use known chemical properties for the prediction of reactivity.

Accurate electronic and chemical properties of 3d transition metal oxides using a calculated linear response U and a DFT + U(V) method.

We validate the usage of the calculated, linear response Hubbard U for evaluating accurate electronic and chemical properties of bulk 3d transition metal oxides. We find calculated values of U lead to improved band gaps. For the evaluation of accurate reaction energies, we first identify and eliminate contributions to the reaction energies of bulk systems due only to changes in U and construct a thermodynamic cycle that references the total energies of unique U systems to a common point using a DFT + U(V) method, which we recast from a recently introduced DFT + U(R) method for molecular systems. We then introduce a semi-empirical method based on weighted DFT/DFT + U cohesive energies to calculate bulk oxidation energies of transition metal oxides using density functional theory and linear response calculated U values. We validate this method by calculating 14 reactions energies involving V, Cr, Mn, Fe, and Co oxides. We find up to an 85% reduction of the mean average error (MAE) compared to energies calculated with the Perdew-Burke-Ernzerhof functional. When our method is compared with DFT + U with empirically derived U values and the HSE06 hybrid functional, we find up to 65% and 39% reductions in the MAE, respectively.

Electrocatalytic oxygen evolution with an immobilized TAML activator.

Iron complexes of tetra-amido macrocyclic ligands are important members of the suite of oxidation catalysts known as TAML activators. TAML activators are known to be fast homogeneous water oxidation (WO) catalysts, producing oxygen in the presence of chemical oxidants, e.g., ceric ammonium nitrate. These homogeneous systems exhibited low turnover numbers (TONs). Here we demonstrate immobilization on glassy carbon and carbon paper in an ink composed of the prototype TAML activator, carbon black, and Nafion and the subsequent use of this composition in heterogeneous electrocatalytic WO. The immobilized TAML system is shown to readily produce O2 with much higher TONs than the homogeneous predecessors.

Structure Sensitive Reaction Kinetics of Chiral Molecules on Intrinsically Chiral Surfaces.

Enantiospecific heterogeneous catalysis utilizes chiral surfaces to resolve enantiomers via structure sensitive surface chemistry. The catalyst design challenge is the identification of chiral surface structures that maximize enantiospecificity. Herein, we develop data driven models for the enantiospecificity of tartaric acid reactions on chiral Cu(hkl)R&S surfaces. Measurements of enantiospecific rate constants were obtained by using curved Cu(hkl)R&S surfaces that enable kinetic measurements on hundreds of chiral surface orientations. One model uses feature vectors derived from generalized coordination numbers to capture the local structure around Cu atoms exposed by the Cu(hkl)R&S surfaces. The second model introduces the use of chiral cubic harmonic functions to capture the symmetry constraints of the face-centered cubic Cu structure. The model using 58 generalized coordination numbers has a fitting error similar to that of the model using only 5 cubic harmonic functions. The two models predict maxima in the enantiospecificity on surfaces with very similar surface orientations. The models developed in this work are applicable for any enantiospecific reaction happening on any chiral material with a cubic lattice structure, opening the way to understanding the surface structure sensitivity of the enantiospecific reaction kinetics.

Effects of strain, d-band filling, and oxidation state on the surface electronic structure and reactivity of 3d perovskite surfaces.

Trends in the dissociative oxygen adsorption energy and oxygen vacancy formation energy on cubic LaBO(3) and SrBO(3) perovskite (001) surfaces (where B = Sc, Ti, V, Cr, Mn, Fe, Co, Ni, and Cu) and their dependence on strain, d-band filling, and oxidation state were examined using density functional theory in the generalized gradient approximation. The effects of strain were found to be small compared to the effects of d-band filling and oxidations state. Electronic structure descriptors such as the d-band center of the B-atom were identified for trends in the dissociative oxygen adsorption energy and for the oxygen vacancy formation energy. A chemical correlation between these two reaction energies was also identified showing the trends in these reaction energies are not independent of each other.

Effects of strain, d-band filling, and oxidation state on the bulk electronic structure of cubic 3d perovskites.

The properties of the d-band structure of the transition metal atom in cubic LaBO(3) and SrBO(3) perovskites (where B = Sc, Ti, V, Cr, Mn, Fe, Co, Ni, and Cu) and their dependence on strain, d-band filling, and oxidation state were investigated using density functional theory calculations and atom-projected density of states. The strain dependence of the d-band width is shown to depend systematically on the size of the B atom. We show that the transition metal d-band width and center are linearly correlated with each other in agreement with a rectangular band model. A simple matrix element formalism based on the solid state table can readily predict the strain dependence of the d-band width.

Configurational correlations in the coverage dependent adsorption energies of oxygen atoms on late transition metal fcc(111) surfaces.

The coverage dependence of oxygen adsorption energies on the fcc(111) surfaces of seven different transition metals (Rh, Ir, Pd, Pt, Cu, Au, and Ag) is demonstrated through density functional theory calculations on 20 configurations ranging from one to five adsorption sites and coverages up to 1 ML. Atom projected densities of states are used to demonstrate that the d-band mediated adsorption mechanism is responsible for the coverage dependence of the adsorption energies. This common bonding mechanism results in a linear correlation that relates the adsorption energies of each adsorbate configuration across different metal surfaces to each other. The slope of this correlation is shown to be related to the characteristics of the valence d-orbitals and band structure of the surface metal atoms. Additionally, it is shown that geometric similarity of the configurations is essential to observe the configurational correlations.

Hydrogen dissociation and spillover on individual isolated palladium atoms.

Using a combination of low-temperature scanning tunneling microscopy and density functional theory it is demonstrated how the nature of an inert host metal of an alloy can affect the thermodynamics and kinetics of a reaction pathway in a much more profound way than simply a dilution, electronic, or geometric effect. This study reveals that individual, isolated Pd atoms can promote H2 dissociation and spillover onto a Cu(111) surface, but that the same mechanism is not observed for an identical array of Pd atoms in Au(111).

Step decoration of chiral metal surfaces.

Highly stepped metal surfaces can define intrinsically chiral structures and these chiral surfaces can potentially be used to separate chiral molecules. The decoration of steps on these surfaces with additional metal atoms is one potential avenue for improving the enantiospecificity of these surfaces. For a successful step decoration, the additional metal atoms should ideally remain at the kinked step sites on the surface. We performed density functional theory (DFT) calculations to identify pairs of metal adatoms and metal surfaces where this kind of step decoration could be thermodynamically stable. These calculations have identified multiple stable examples of step decoration. Using our DFT results, we developed a model to predict surface segregation on a wide range of stepped metal surfaces. With this model, we have estimated the stability of step decoration without further DFT calculations for surface segregation for all combinations of the 3d, 4d, and 5d metals.

First-Principles Investigation of the Epitaxial Stabilization of Oxide Polymorphs: TiO2 on (Sr,Ba)TiO3.

Metastable polymorphs, many of which have never been fabricated, have been predicted to exhibit interesting and technologically relevant properties. Epitaxial synthesis is a powerful structure-directing method that can produce metastable polymorphs but is typically done in a trial and error fashion. Unfortunately, the relevant thermodynamic terms governing epitaxial synthesis of new materials are unknown. Accurate calculation of the relevant thermodynamic terms and their incorporation into predictive models would accelerate the synthesis of metastable polymorphs by identifying thermodynamically favorable paths. Using density functional theory with three different functionals, we computed several relevant terms for TiO2 anatase (A) and rutile (R) film growth on low-index surfaces of SrTiO3 (STO) and BaTiO3 (BTO) cubic perovskites. After identifying potential coherent epitaxial interfaces based on experimental observations, the volumetric formation, volumetric strain, and areal substrate-film interface energies were calculated for (001)A∥(001)(S/B)TO, (102)A∥(011)(S/B)TO, (100)R∥(111)(S/B)TO, and (112)A∥(111)(S/B)TO coherent interfaces. These terms were integrated into a standard model of epitaxial nucleation, and the results yielded reasonable agreement between experimental observations and DFT predictions of the preferred epitaxial polymorph. Predicted trends in epitaxial stability were essentially independent of the three functionals used in the calculations. These results are discussed in light of their promise that DFT-informed epitaxial film growth can accelerate fabrication of new polymorphs. These results also validate the recently proposed 20 kJ/mol stability window for predicting which polymorphs could be epitaxially stabilized.

Alchemical Predictions for Computational Catalysis: Potential and Limitations.

Kohn-Sham density functional theory (DFT) is the workhorse method for calculating adsorbate binding energies relevant for catalysis. Unfortunately, this method is too computationally expensive to methodically and broadly search through catalyst candidate space. Here, we assess the promise of computational alchemy, a perturbation theory approach that allows for predictions of binding energies thousands of times faster than DFT. We first benchmark the binding energy predictions of oxygen reduction reaction intermediates on alloys of Pt, Pd, and Ni using alchemy against predictions from DFT. Far faster alchemical estimates yield binding energies within 0.1 eV of DFT values in many cases. We also identify distinct cases where alchemy performs significantly worse, indicating areas where modeling improvements are needed. Our results suggest that computational alchemy is a very promising tool that warrants further consideration for high-throughput screening of heterogeneous catalysts.

The atomic simulation environment-a Python library for working with atoms.

The atomic simulation environment (ASE) is a software package written in the Python programming language with the aim of setting up, steering, and analyzing atomistic simulations. In ASE, tasks are fully scripted in Python. The powerful syntax of Python combined with the NumPy array library make it possible to perform very complex simulation tasks. For example, a sequence of calculations may be performed with the use of a simple 'for-loop' construction. Calculations of energy, forces, stresses and other quantities are performed through interfaces to many external electronic structure codes or force fields using a uniform interface. On top of this calculator interface, ASE provides modules for performing many standard simulation tasks such as structure optimization, molecular dynamics, handling of constraints and performing nudged elastic band calculations.

Identifying potential BO2 oxide polymorphs for epitaxial growth candidates.

Transition metal dioxides (BO2) exhibit a number of polymorphic structures with distinct properties, but the isolation of different polymorphs for a given composition is carried out using trial and error experimentation. We present computational studies of the relative stabilities and equations of state for six polymorphs (anatase, brookite, rutile, columbite, pyrite, and fluorite) of five different BO2 dioxides (B = Ti, V, Ru, Ir, and Sn). These properties were computed in a consistent fashion using several exchange correlation functionals within the density functional theory formalism, and the effects of the different functionals are discussed relative to their impact on predictive synthesis. We compare the computational results to prior observations of high-pressure synthesis and epitaxial film growth and then use this discussion to predict new accessible polymorphs in the context of epitaxial stabilization using isostructural substrates. For example, the relative stabilities of the columbite polymorph for VO2 and RuO2 are similar to those of TiO2 and SnO2, the latter two of which have been previously stabilized as epitaxial films.