Multiscale Perturbation Methods for Dynamic/Programmable Catalysis
Carolina Colombo Tedesco, Carl D. Laird, John R. Kitchin, Aditya S. Khair, “Multiscale Perturbation Methods for Dynamic/Programmable Catalysis”, Industrial & Engineering Chemistry Research 64(45), 21438–21448 (2025). https://doi.org/10.1021/acs.iecr.5c03023
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.
Note on the source text: the extracted text is the open-access (CC-BY 4.0) article. Prose is readable, but the two-column layout interleaves sentences and most equations are garbled, so the mathematical forms below are described in words.
Question
In dynamic or programmable catalysis, the binding energies of adsorbates are oscillated periodically so that different elementary steps are favored at different moments, with the aim of exceeding the turnover frequency allowed by the Sabatier principle for any static catalyst. Simulating these systems is expensive because the time-dependent rate constants give stiff ODEs that must either be integrated over many periods to reach a cyclic steady state, or solved as a boundary value problem (BVP), as in the authors’ earlier Pyomo-based work. This paper asks whether perturbation theory can instead give closed-form expressions for the time-averaged surface coverages and for the limit cycle around them (Colombo Tedesco et al. 2025).
Approach
The authors use a multiple time-scale expansion. Each rate constant is written as an amplitude times a waveform plus an offset. A slow time is defined by scaling with the smallest rate-constant amplitude, and a fast time by dividing by the forcing period; their ratio is a small parameter equal to that smallest amplitude times the period. Coverages are expanded as a leading-order term plus a small correction proportional to this parameter. At leading order the coverage is constant over the fast time, so it equals the cycle average. The averages follow from averaging the rate equations over one period and setting the time derivative to zero, using the approximation that the average of a product of rate constant and coverage equals the product of their averages (that is, neglecting their covariance). With constant gas-phase activities the averaged equations are linear algebraic equations. Subtracting the averaged equation from the first-order equation gives the fast fluctuation, which is obtained by direct integration over the period, with the integration constant fixed by requiring the fluctuation to have zero mean. For square waves all integrals are analytic and piecewise linear in the fast time. The results were compared with the authors’ BVP simulations, and the code is provided as Colab notebooks on Zenodo.
Key findings
Two linear schemes taken from Foley and Razdan were analyzed. Scheme I has one adsorbed intermediate whose formation and consumption rate constants alternate as symmetric square waves (one between 1 and 0 per second, the other between 0 and 1000 per second). The average coverage reduces to the averaged forward rate constant divided by the sum of averaged forward and consumption rate constants, giving 0.001, matching Foley and Razdan without needing the limit cycle. The first-order correction gives a piecewise linear limit cycle that agrees closely with the BVP result, apart from a small horizontal offset.
Scheme II has two adsorbates, A and B, in a three-step sequence with forward and reverse rate constants computed from a microkinetic model under a square-wave binding-energy oscillation (amplitude 0.1 eV, offset 1.4 eV, duty cycle 0.5, 1000 Hz). Solving the two averaged linear equations gave average coverages of 0.3828 for A and 0.6172 for B, and the analytical piecewise limit cycles again agreed closely with simulation. The analytical result was essentially instantaneous, compared with about 1 s for the BVP approach and about 40 s for forward initial-value integration.
The method requires clear time-scale separation. Taking the largest rate constant as the characteristic rate (1000 Hz for Scheme I and 20 Hz for Scheme II), agreement was good for forcing frequencies above that threshold. At 3 Hz the simulated system follows the stepwise regime, tracking the input with curved transients, while the analytical solution stays quasi-static and its oscillations become exaggerated, possibly violating the site balance. A larger amplitude of 0.3 eV at 1000 Hz also broke the limit-cycle prediction (underestimating swings in A, overshooting in B, with coverages not summing to one because site balance is not enforced), but the average coverages still differed from simulation by at most 0.04. Computed covariances between rate constants and coverages grow as frequency falls or amplitude rises, providing a diagnostic for when the method fails. The authors argue that small amplitudes correspond to points near the peak of their extended Sabatier volcano and are the easier case experimentally. For nonlinear systems, such as time-varying gas activities, the averaged equations become nonlinear algebraic equations that could be solved numerically within the same framework.
Significance
The work gives a fast, analytical alternative to numerical cyclic steady-state simulation for programmable catalysis in the high-frequency, small-amplitude regime, useful for rapid screening, preliminary design, or initializing numerical solvers. The authors also draw analogies to RC circuits, molecular ratchets and other periodically forced processes, suggesting the ideas may extend beyond catalysis (Colombo Tedesco et al. 2025).