One Decade of Computational Studies on Single-Atom Alloys: Is In Silico Design within Reach?

ConspectusSingle-Atom alloys (SAAs) are an emerging class of materials consisting of a coinage metal (Cu, Ag, and Au) doped, at the single-atom limit, with another metal. As catalysts, coinage metals are rarely very active on their own, but when they are, they exhibit high selectivity. On the other hand, transition metals are usually very active but not as selective. Incorporating transition metals (guest elements) into coinage metals (host material) is therefore appealing for combining the activity and selectivity of each constituent in a balanced way. Additionally, first-principles calculations have shown that single atoms embedded in the surface of a coinage metal can exhibit emergent properties. Here, we describe how computational studies based on density functional theory (DFT) and kinetic Monte Carlo (KMC) simulations, often undertaken in close collaboration with experimental research groups, have shaped, over the past decade, the way we understand SAA catalysis.This Account reviews our contributions in elucidating the stability of SAAs, their electronic structure, and the way adsorbates interact and react on SAA catalytic surfaces. By studying in detail the processes that affect the stability of the SAA phase, we have shown that out of several bimetallic combinations of coinage metals with prominent Pt-group metals only PtCu and PdCu are stable surface alloys under vacuum. However, more surface alloy structures are possible in the presence of adsorbates because the latter can stabilize, via strong binding, dopants in the surface of the material. More interestingly, a large number of these surface alloys are resistant to the aggregation of dopant atoms into clusters, thereby favoring the SAA structure. These major results from DFT calculations serve as a guide for experimentalists to explore new SAA catalysts. Further analysis has shown that SAAs have a unique electronic structure with a very sharp d-band feature close to the Fermi level, analogous to the electronic structure of molecular entities. This is one of the reasons that SAAs are particularly sought after: although they are metallic nanoparticles, they have properties akin to those of homogeneous catalysts. In this context, we have contributed extensive screening studies, focusing on molecular fragments of catalytic relevance on a range of SAAs, which have driven the identification of new catalysts. We have also explored the rich chemistry of two-adsorbate systems via kinetic modeling, demonstrating how a spectator species with greater affinity for the dopant can modulate the reactivity of the catalyst via the so-called (punctured) molecular cork effect.Since the first experimental characterization of SAAs about a decade ago, theoretical models have been able to support and explain various experimental observations. These models have served as benchmarks for assessing the predictive capability of the underlying theoretical methods. In turn, the predictions that have been delivered have guided and continue to guide the experimental research efforts in the field. These advancements show that the in silico design of new SAA catalysts is now within reach.

Exact distributed kinetic Monte Carlo simulations for on-lattice chemical kinetics: lessons learnt from medium- and large-scale benchmarks.

Kinetic Monte Carlo (KMC) simulations have been instrumental in multiscale catalysis studies, enabling the elucidation of the complex dynamics of heterogeneous catalysts and the prediction of macroscopic performance metrics, such as activity and selectivity. However, the accessible length- and time-scales have been a limiting factor in such simulations. For instance, handling lattices containing millions of sites with 'traditional' sequential KMC implementations is prohibitive owing to large memory requirements and long simulation times. We have recently established an approach for exact, distributed, lattice-based simulations of catalytic kinetics which couples the Time-Warp algorithm with the Graph-Theoretical KMC framework, enabling the handling of complex adsorbate lateral interactions and reaction events within large lattices. In this work, we develop a lattice-based variant of the Brusselator system, a prototype chemical oscillator pioneered by Prigogine and Lefever in the late 60s, to benchmark and demonstrate our approach. This system can form spiral wave patterns, which would be computationally intractable with sequential KMC, while our distributed KMC approach can simulate such patterns 15 and 36 times faster with 625 and 1600 processors, respectively. The medium- and large-scale benchmarks thus conducted, demonstrate the robustness of the approach, and reveal computational bottlenecks that could be targeted in further development efforts. This article is part of a discussion meeting issue 'Supercomputing simulations of advanced materials'.

Large-scale benchmarks of the time-warp/graph-theoretical kinetic Monte Carlo approach for distributed on-lattice simulations of catalytic kinetics.

Motivated by the need to perform large-scale kinetic Monte Carlo (KMC) simulations, in the context of unravelling complex phenomena such as catalyst reconstruction and pattern formation, we extend the work of Ravipati et al. [S. Ravipati, G. D. Savva, I.-A. Christidi, R. Guichard, J. Nielsen, R. Réocreux and M. Stamatakis, Comput. Phys. Commun., 2022, 270, 108148] in benchmarking the performance of a distributed-computing, on-lattice KMC approach. The latter, implemented in our software package Zacros, combines the graph-theoretical KMC framework with the Time-Warp algorithm for parallel discrete event simulations, and entails dividing the lattice into subdomains, each assigned to a processor. The cornerstone of the Time-Warp algorithm is the state queue, to which snapshots of the simulation state are saved regularly, enabling historical KMC information to be corrected when conflicts occur at subdomain boundaries. Focusing on three model systems, we highlight the key Time-Warp parameters that can be tuned to optimise performance. The frequency of state saving, controlled by the state saving interval, δsnap, is shown to have the largest effect on performance, which favours balancing the overhead of re-simulating KMC history with that of writing state snapshots to memory. Also important is the global virtual time (GVT) computation interval, ΔτGVT, which has little direct effect on the progress of the simulation but controls how often the state queue memory can be freed up. We also find that pre-allocating memory for the state queue data structure favours performance. These findings will guide users in maximising the efficiency of Zacros or other distributed KMC software, which is a vital step towards realising accurate, meso-scale simulations of heterogeneous catalysis.

Speeding up the Detection of Adsorbate Lateral Interactions in Graph-Theoretical Kinetic Monte Carlo Simulations.

Kinetic Monte Carlo (KMC) has become an indispensable tool in heterogeneous catalyst discovery, but realistic simulations remain computationally demanding on account of the need to capture complex and long-range lateral interactions between adsorbates. The Zacros software package (https://zacros.org) adopts a graph-theoretical cluster expansion (CE) framework that allows such interactions to be computed with a high degree of generality and fidelity. This involves solving a series of subgraph isomorphism problems in order to identify relevant interaction patterns in the lattice. In an effort to reduce the computational burden, we have adapted two well-known subgraph isomorphism algorithms, namely, VF2 and RI, for use in KMC simulations and implemented them in Zacros. To benchmark their performance, we simulate a previously established model of catalytic NO oxidation, treating the O* lateral interactions with a series of progressively larger CEs. For CEs with long-range interactions, VF2 and RI are found to provide impressive speedups relative to simpler algorithms. RI performs best, giving speedups reaching more than 150× when combined with OpenMP parallelization. We also simulate a recently developed methane cracking model, showing that RI offers significant improvements in performance at high surface coverages.

Single-Atom Alloys for the Electrochemical Oxygen Reduction Reaction.

Single-atom alloys (SAAs) consisting of isolated transition-metal atoms doped in the surface of coinage metal hosts exhibit unique catalytic properties, harnessing the high activity of the dopant metals with the selectivity of the coinage metal hosts. Here we use density functional theory (DFT) to study SAAs comprised of Ni, Pd, Pt, Co and Rh doped into Ag and Au hosts, as candidate electrocatalysts for the oxygen reduction reaction (ORR) in proton-exchange membrane (PEM) fuel-cells. Our calculations reveal that the PdAu SAA exhibits a slightly lower theoretical overpotential, enhanced selectivity for 4-e- ORR, and tolerance to CO-poisoning compared to Pt(111). While the number of active sites of PdAu SAA is lower than that of Pt(111), the aforementioned desirable properties could bring the overall catalytic performance thereof close to that of Pt/C, indicating that the PdAu SAA could be a viable material for electrocatalytic ORR in PEM fuel-cells.

DFT benchmark studies on representative species and poisons of methane steam reforming on Ni(111).

Ni catalysts used in methane steam reforming (MSR) are highly susceptible to poisoning by carbon-based species, which poses a major impediment to the productivity of industrial operations. These species encompass graphitic carbon-like formations that are typically modelled as graphene. First principles-based approaches, such as density functional theory (DFT) calculations, can provide valuable insight into the mechanism of graphene growth in the MSR reaction. It is, however, critical that a DFT model of this reaction can accurately describe the interactions of Ni(111) with the MSR intermediates as well as graphene. In this work, a systematic benchmark study has been carried out to identify a suitable DFT functional for the graphene-MSR system. The binding energies of graphene and important MSR species, as well as the reaction energies of methane dissociation and carbon oxidation, were computed on Ni(111) using GGA functionals, DFT-D and van der Waals density functionals (vdW-DF). It is well-established that the GGA functionals are inadequate for describing graphene-Ni(111) interactions. In the case of vdW-DF, the optPBE-vdW functional predicts the binding energies of graphene and several important MSR species with reasonable accuracy; however, it provides poor estimates of CO and O binding energies. Among the DFT-D functionals, PBE-D3 and PBE-dDsC have been found to exhibit acceptable accuracy for graphene and most MSR species (excluding adsorbed CO), and therefore, both functionals are promising for elucidating carbon-based catalytic poisoning in the MSR reaction. Overall, no single DFT functional could estimate the binding energies of all the species with equally high accuracy.

Mechanistic insights into carbon-carbon coupling on NiAu and PdAu single-atom alloys.

Carbon-carbon coupling is an important step in many catalytic reactions, and performing sp3-sp3 carbon-carbon coupling heterogeneously is particularly challenging. It has been reported that PdAu single-atom alloy (SAA) model catalytic surfaces are able to selectively couple methyl groups, producing ethane from methyl iodide. Herein, we extend this study to NiAu SAAs and find that Ni atoms in Au are active for C-I cleavage and selective sp3-sp3 carbon-carbon coupling to produce ethane. Furthermore, we perform ab initio kinetic Monte Carlo simulations that include the effect of the iodine atom, which was previously considered a bystander species. We find that model NiAu surfaces exhibit a similar chemistry to PdAu, but the reason for the similarity is due to the role the iodine atoms play in terms of blocking the Ni atom active sites. Specifically, on NiAu SAAs, the iodine atoms outcompete the methyl groups for occupancy of the Ni sites leaving the Me groups on Au, while on PdAu SAAs, the binding strengths of methyl groups and iodine atoms at the Pd atom active site are more similar. These simulations shed light on the mechanism of this important sp3-sp3 carbon-carbon coupling chemistry on SAAs. Furthermore, we discuss the effect of the iodine atoms on the reaction energetics and make an analogy between the effect of iodine as an active site blocker on this model heterogeneous catalyst and homogeneous catalysts in which ligands must detach in order for the active site to be accessed by the reactants.

A DFT and KMC based study on the mechanism of the water gas shift reaction on the Pd(100) surface.

We present a combined density functional theory (DFT) and Kinetic Monte Carlo (KMC) study of the water gas shift (WGS) reaction on the Pd(100) surface. We propose a mechanism comprising both the redox and the associative pathways for the WGS within a single framework, which consists of seven core elementary steps, which in turn involve splitting of a water molecule followed by the production of an H-atom and an OH-species on the Pd(100) surface. In the following steps, these intermediates then recombine with each other and with CO leading to the evolution of CO2, and H2. Seven other elementary steps, involving the diffusion and adsorption of the surface intermediate species are also considered for a complete description of the mechanism. The geometrical and electronic properties of each of the reactants, products, and the transition states of the core elementary steps are presented. We also discuss the analysis of Bader charges and spin densities for the reactants, transition states and the products of these elementary steps. Our study indicates that the WGS reaction progresses simultaneously via the direct oxidation and the carboxyl paths on the Pd(100) surface.

Comparison of Queueing Data-Structures for Kinetic Monte Carlo Simulations of Heterogeneous Catalysts.

On-lattice kinetic Monte Carlo (KMC) is a computational method used to simulate (among others) physicochemical processes on catalytic surfaces. The KMC algorithm propagates the system through discrete configurations by selecting (with the use of random numbers) the next elementary process to be simulated, e.g., adsorption, desorption, diffusion, or reaction. An implementation of such a selection procedure is the first-reaction method in which all realizable elementary processes are identified and assigned a random occurrence time based on their rate constant. The next event to be executed will then be the one with the minimum interarrival time. Thus, a fast and efficient algorithm for selecting the most imminent process and performing all of the necessary updates on the list of realizable processes post execution is of great importance. In the current work, we implement five data-structures to handle the elementary process queue during a KMC run: an unsorted list, a binary heap, a pairing heap, a one-way skip list, and finally, a novel two-way skip list with a mapping array specialized for KMC simulations. We also investigate the effect of compiler optimizations on the performance of these data-structures on three benchmark models, capturing CO oxidation, a simplified water gas shift mechanism, and a temperature-programmed desorption run. Excluding the least efficient and impractical for large-problems unsorted list, we observe a 3× speedup of the binary or pairing heaps (most efficient) compared to the one-way skip list (least efficient). Compiler optimizations deliver a speedup of up to 1.8×. These benchmarks provide valuable insight into the importance of, often-overlooked, implementation-related aspects of KMC simulations, such as the queueing data-structures. Our results could be particularly useful in guiding the choice of data-structures and algorithms that would minimize the computational cost of large-scale simulations.

A Caching Scheme To Accelerate Kinetic Monte Carlo Simulations of Catalytic Reactions.

Kinetic Monte Carlo (KMC) simulations have been instrumental in advancing our fundamental understanding of heterogeneously catalyzed reactions, with particular emphasis on structure sensitivity, ensemble effects, and the interplay between adlayer structure and adsorbate-adsorbate lateral interactions in shaping the observed kinetics. Yet, the computational cost of KMC remains high, thereby motivating the development of acceleration schemes that would improve the simulation efficiency. We present an exact such scheme, which implements a caching algorithm along with shared-memory parallelization to improve the computational performance of simulations incorporating long-range adsorbate-adsorbate lateral interactions. This scheme is based on caching information about the energetic interaction patterns associated with the products of each possible lattice process (adsorption, desorption, reaction etc.). Thus, every time a reaction occurs ("ongoing reaction"), it enables fast updates of the rate constants of "affected reactions", i.e., possible reactions in the region of influence of the "ongoing reaction". Benchmarks on KMC simulations of NOx oxidation/reduction, yielded acceleration factors of up to 20, when comparing single-thread runs without caching to runs on 16 threads with caching, for simulations with a cluster expansion Hamiltonian that incorporates up to 8th-nearest-neighbor interactions.

Surface facet dependence of competing alloying mechanisms.

Metal alloys are ubiquitous in many branches of heterogeneous catalysis, and it is now fairly well established that the local atomic structure of an alloy can have a profound influence on its chemical reactivity. While these effects can be difficult to probe in nanoparticle catalysts, model studies using well defined single crystal surfaces alloyed with dopants enable these structure-function correlations to be drawn. The first step in this approach involves understanding the alloying mechanism and the type of ensembles formed. In this study, we examined the atomic structure of RhCu single-atom alloys formed on Cu(111), Cu(100), and Cu(110) surfaces. Our results show a striking difference between Rh atoms alloying in Cu(111) vs the more open Cu(100) and Cu(110) surface facets. Unlike Cu(111) on which Rh atoms preferentially place-exchange with Cu atoms in the local regions above step edges leaving the majority of the Cu surface free of Rh, highly dispersed, homogeneous alloys are formed on the Cu(100) and (110) surfaces. These dramatically different alloying mechanisms are understood by quantifying the energetic barriers for atomic hopping, exchange, swapping, and vacancy filling events for Rh atoms on different Cu surfaces through theoretical calculations. Density functional theory results indicate that the observed differences in the alloying mechanism can be attributed to a faster hopping rate, relatively high atomic exchange barriers, and stronger binding of Rh atoms in the vicinity of step edges on Cu(111) compared to Cu(110) and Cu(100). These model systems will serve as useful platforms for examining structure sensitive chemistry on single-atom alloys.

Adlayer structure and lattice size effects on catalytic rates predicted from KMC simulations: NO oxidation on Pt(111).

Repulsive and/or attractive interactions between surface adsorbates have an important effect on the structure of the adsorbate layer and consequently on the rate of heterogeneous catalytic reactions. Thus, developing reaction models that take into account adsorbate-adsorbate interactions is crucial for making accurate predictions of the catalytic rate and surface coverage during reaction. In the present work, we employ kinetic Monte Carlo simulation to model the catalytic NO oxidation on Pt (111), adopting a cluster expansion (CE) Hamiltonian approach for treating the aforementioned interactions. We investigate CEs of increasing complexity, ranging from pairwise 1st nearest neighbor to long-range and many-body terms. We show that energetic models incorporating solely short-range interactions result in ordered adlayer structures, which are disrupted by anti-phase boundaries and defective regions when the size of the periodic lattice is non-commensurate to the structure of the stable adlayer. We find that O2 dissociates on sites located in these defective regions, which are predominantly responsible for the activity, and the predicted catalytic rate is strongly depended on the lattice size. Such effects are absent when employing non-periodic lattices, whereon the catalytic activity appears more intense on edges/corner sites. Finally, inclusion of long-range interactions in the model Hamiltonian induces relative disorder in the adsorbate layer, which is ascribed to the "softening" of the repulsive interactions between adspecies. Under these circumstances, the distribution of activation energies for O2 dissociation is broader as compared to short-range interaction models and on this basis we explain the disparate catalytic rate predictions when using different CEs.

Non-Equilibrium Thermodynamics and Stochastic Dynamics of a Bistable Catalytic Surface Reaction.

Catalytic surface reaction networks exhibit nonlinear dissipative phenomena, such as bistability. Macroscopic rate law descriptions predict that the reaction system resides on one of the two steady-state branches of the bistable region for an indefinite period of time. However, the smaller the catalytic surface, the greater the influence of coverage fluctuations, given that their amplitude normally scales as the square root of the system size. Thus, one can observe fluctuation-induced transitions between the steady-states. In this work, a model for the bistable catalytic CO oxidation on small surfaces is studied. After a brief introduction of the average stochastic modelling framework and its corresponding deterministic limit, we discuss the non-equilibrium conditions necessary for bistability. The entropy production rate, an important thermodynamic quantity measuring dissipation in a system, is compared across the two approaches. We conclude that, in our catalytic model, the most favorable non-equilibrium steady state is not necessary the state with the maximum or minimum entropy production rate.

A machine learning approach to graph-theoretical cluster expansions of the energy of adsorbate layers.

The accurate description of the energy of adsorbate layers is crucial for the understanding of chemistry at interfaces. For heterogeneous catalysis, not only the interaction of the adsorbate with the surface but also the adsorbate-adsorbate lateral interactions significantly affect the activation energies of reactions. Modeling the interactions of the adsorbates with the catalyst surface and with each other can be efficiently achieved in the cluster expansion Hamiltonian formalism, which has recently been implemented in a graph-theoretical kinetic Monte Carlo (kMC) scheme to describe multi-dentate species. Automating the development of the cluster expansion Hamiltonians for catalytic systems is challenging and requires the mapping of adsorbate configurations for extended adsorbates onto a graphical lattice. The current work adopts machine learning methods to reach this goal. Clusters are automatically detected based on formalized, but intuitive chemical concepts. The corresponding energy coefficients for the cluster expansion are calculated by an inversion scheme. The potential of this method is demonstrated for the example of ethylene adsorption on Pd(111), for which we propose several expansions, depending on the graphical lattice. It turns out that for this system, the best description is obtained as a combination of single molecule patterns and a few coupling terms accounting for lateral interactions.

From Biomass-Derived Furans to Aromatics with Ethanol over Zeolite.

We report a novel catalytic conversion of biomass-derived furans and alcohols to aromatics over zeolite catalysts. Aromatics are formed via Diels-Alder cycloaddition with ethylene, which is produced in situ from ethanol dehydration. The use of liquid ethanol instead of gaseous ethylene, as the source of dienophile in this one-pot synthesis, makes the aromatics production much simpler and renewable, circumventing the use of ethylene at high pressure. More importantly, both our experiments and theoretical studies demonstrate that the use of ethanol instead of ethylene, results in significantly higher rates and higher selectivity to aromatics, due to lower activation barriers over the solid acid sites. Synchrotron-diffraction experiments and proton-affinity calculations clearly suggest that a preferred protonation of ethanol over the furan is a key step facilitating the Diels-Alder and dehydration reactions in the acid sites of the zeolite.

Controlling a spillover pathway with the molecular cork effect.

Spillover of reactants from one active site to another is important in heterogeneous catalysis and has recently been shown to enhance hydrogen storage in a variety of materials. The spillover of hydrogen is notoriously hard to detect or control. We report herein that the hydrogen spillover pathway on a Pd/Cu alloy can be controlled by reversible adsorption of a spectator molecule. Pd atoms in the Cu surface serve as hydrogen dissociation sites from which H atoms can spillover onto surrounding Cu regions. Selective adsorption of CO at these atomic Pd sites is shown to either prevent the uptake of hydrogen on, or inhibit its desorption from, the surface. In this way, the hydrogen coverage on the whole surface can be controlled by molecular adsorption at a minority site, which we term a 'molecular cork' effect. We show that the molecular cork effect is present during a surface catalysed hydrogenation reaction and illustrate how it can be used as a method for controlling uptake and release of hydrogen in a model storage system.

Transition Metal Carbides as Supports for Catalytic Metal Particles: Recent Progress and Opportunities.

Transition metal carbides (TMCs) constitute excellent alternatives to traditional oxide-based supports for small metal particles, leading to strong metal-support interactions, which drastically modify the catalytic properties of the supported metal atoms. Moreover, they possess extremely high melting points and good resistance to carbon deposition and sulfur poisoning, and the catalytic activities of some TMCs per se have been shown to be similar to those of Pt-group metals for a considerable number of reactions. Therefore, the use of TMCs as supports can give rise to bifunctional catalysts with multiple active sites. However, at present, only TiC and MoxC have been tested experimentally as supports for metal particles, and it is largely unclear which combinations may best catalyze which chemical reactions. In this Perspective, we review the most significant works on the use of TMCs as supports for catalytic applications, assess the current status of the field, and identify key advances being made and challenges, with an eye to the future.

Ten-electron count rule for the binding of adsorbates on single-atom alloy catalysts.

Single-atom alloys have recently emerged as highly active and selective alloy catalysts. Unlike pure metals, single-atom alloys escape the well-established conceptual framework developed nearly three decades ago for predicting catalytic performance. Although this offers the opportunity to explore so far unattainable chemistries, this leaves us without a simple guide for the design of single-atom alloys able to catalyse targeted reactions. Here, based on thousands of density functional theory calculations, we reveal a 10-electron count rule for the binding of adsorbates on the dopant atoms, usually the active sites, of single-atom alloy surfaces. A simple molecular orbital approach rationalizes this rule and the nature of the adsorbate-dopant interaction. In addition, our intuitive model can accelerate the rational design of single-atom alloy catalysts. Indeed, we illustrate how the unique insights provided by the electron count rule help identify the most promising dopant for an industrially relevant hydrogenation reaction, thereby reducing the number of potential materials by more than one order of magnitude.