WhereWulff: A Semiautonomous Workflow for Systematic Catalyst Surface Reactivity under Reaction Conditions
Rohan Yuri Sanspeur, Javier Heras-Domingo, John R. Kitchin and Zachary Ulissi, “WhereWulff: A Semiautonomous Workflow for Systematic Catalyst Surface Reactivity under Reaction Conditions”, Journal of Chemical Information and Modeling 63(8), 2427-2437 (2023). https://doi.org/10.1021/acs.jcim.3c00142
This summary was written by a language model from the paper and has not been checked line by line. Rely on the published paper for details.
Question
Computing the oxygen evolution reaction (OER) activity of a multicomponent oxide by density functional theory (DFT) involves many manual, expensive steps: optimizing the bulk and its magnetic state, cleaving many slabs, deciding which facets and which surface terminations matter under reaction conditions, and then placing and relaxing every intermediate. Large slabs of correlated oxides often exceed cluster wall-time limits. The paper asks whether this expert procedure can be encoded into an open, largely autonomous workflow that also decides which calculations are worth running, so that it can serve as a first-principles check on machine-learning screening pipelines (Sanspeur et al. 2023).
Approach
The authors built two FireWorks-based workflows on top of pymatgen, atomate and Custodian, with VASP as the DFT engine (Sanspeur et al. 2023). An auxiliary bulk workflow relaxes the input crystal, then fits equations of state for nonmagnetic, ferromagnetic and antiferromagnetic orderings, with initial moments set from oxidation states and crystal field theory. It picks the lowest-energy configuration and computes the energy above the convex hull and the Pourbaix decomposition free energy at a chosen pH and potential, as measures of synthesizability and electrochemical stability.
The main workflow, WhereWulff, takes that bulk, enumerates symmetric (nonpolar) slabs up to a user-set maximum Miller index, relaxes them, and computes surface energies, using a metal reference species and chemical potentials for nonstoichiometric slabs. It then builds the Wulff shape and keeps only facets that contribute above a user threshold (10% in this work) to the particle shape, which is the origin of the name. On each kept facet it places full monolayers of OH* and O* (with rotational sampling of OH*), builds a surface Pourbaix diagram over the clean, OH* and O* terminations, selects the most stable one at the requested conditions, and computes the OH*, O* and OOH* free energies on a chosen metal site, including OOH* rotations and two hydrogen orientations. The final step stores the four reaction free energies, the potential-determining step and the theoretical overpotential. A ContinueOptimizeFW task uses Custodian to stop VASP cleanly before the wall time runs out and inserts a continuation job through FireWorks detours, passing geometry, moments and wavefunction through UUID-tagged database records. An end-to-end regression test on IrO2 runs in GitHub Actions with simulated VASP, and the same mechanism allows “hybrid” workflows that restart from a midpoint. Calculations used PBE with Dudarev +U on Co 3d, a 500 eV cutoff and 0.05 eV/Å force convergence.
Key findings
As a benchmark, RuO2 and IrO2 (110) and (101) were both O*-covered under acidic OER conditions, with theoretical overpotentials of 0.52 and 0.53 V for RuO2 and 0.49 and 0.39 V for IrO2 on (110) and (101). The potential-determining step was O* to OOH* on (110) and OH* to O* on (101), consistent with the literature (Sanspeur et al. 2023).
The workflow was then applied to two double perovskites, BaSrCo2O6 and BaSnTi2O6, that an earlier machine-learning screening had proposed. For both, the (001), (101) and (110) facets passed the Wulff filter and six other terminations were dropped, saving at least 108 DFT calculations and cutting the total from about 240 to about 132. Under alkaline conditions, the BaSrCo2O6 surfaces were clean and the BaSnTi2O6 surfaces were OH*-covered. BaSrCo2O6-(001) had the lowest overpotential, 0.63 V, within the 0.5 to 0.75 V range previously reported, with OH* to O* as the limiting step on all three BaSrCo2O6 facets. The BaSnTi2O6 facets gave 0.73 to 0.86 V. The choice of termination mattered: BaSnTi2O6-(110) gave 2.06 V on an O*-covered surface but 0.73 V on the OH*-covered one, whereas BaSrCo2O6-(001) changed only from 0.63 to 0.62 V between clean and OH* coverage, which the authors attribute to its widely spaced Co sites. The workflow handled about 180 resubmissions without intervention, which on some hosts amounted to as much as about 75% of the root jobs under a 48 h limit.
One small inconsistency in the text: the machine-learning source of the two perovskites is described as Zheng et al., while the cited reference (ref. 30) is by Li, Achenie and Xin.
Significance
WhereWulff is presented as the first public, scalable encoding of the full OER modeling procedure for oxides, in which upstream stability information (surface energies and surface Pourbaix diagrams) decides which downstream activity calculations are run (Sanspeur et al. 2023). The authors see four uses: as a DFT source of truth for closed-loop discovery, as a data generator, as a teaching tool for experimentalists, and as a base for other reactions such as H2O2 production. They note that DFT cost remains the bottleneck and propose substituting OC22-trained surrogate models, and adding solvent and coverage effects to the surface energies before the Wulff step. The code is at github.com/ulissigroup/wherewulff and the data at Zenodo (10.5281/zenodo.7600476).