(371b) Forcing Signal Optimization and Simultaneous Simulation Approaches Applied to Dynamic/Programmable Catalysis

The search for the best catalyst for a specific chemical reaction has for long relied on the use of the Sabatier Volcano plot. This curve describes the relation between the adsorption or binding energy (BE) of a substrate on a surface and a descriptor of overall reaction rate, usually the turnover frequency (TOF). Generally, the plot shows the transition of the kinetically relevant step from one elementary step to another as the BE increases. At an intermediate BE value, an optimal point (highest TOF) arises, where there exists a balance between individual step’s rates.

To move beyond the Sabatier Volcano, researchers have proposed and investigated the use of dynamic or programmable catalysis. By applying this technique, the binding energy is set (or programmed) to oscillate with a predetermined shape, leading to the alternating promotion of the elementary steps. By doing so, the overall limitations associated with having a rate limiting step are avoided. The effectiveness of the idea has been demonstrated computationally^1-5, and experimentally^6,7, leading to interest in exploring the novel technique. In this work, we describe a methodology to efficiently simulate the cyclic steady state (CSS) during a dynamic catalytic system operation and the optimization of parameters of the forcing signal to maximize the time-averaged turnover frequency (avTOF) achieved. We discuss how conclusions and insights from this work can be projected into what we call the Extended Sabatier Volcano, in which the dynamic behavior is visualized.

As in previous studies, we used a model unimolecular reaction catalyzed on a surface as the catalytic system for the proceedings. We considered three reversible steps: adsorption, surface reaction and desorption. The rate constants of the reaction and desorption steps were considered functions of the binding energy, and are, therefore, the steps that are periodically favored with the simulation of a modular behavior for the binding energy. The system evolution was described by four ordinary differential equations (ODEs), describing gas phase concentrations and coverage fractions changes over time, and an algebraic equation for the site balance. Previous works have performed the simulations by applying forward integration of the ODEs until CSS is reached^{1, 8, 9}, a process that can be expensive, as it is highly dependent on the stiffness of the system.

In this work, we formulate the problem as a Boundary Value Problem (BVP) with limit cycle conditions to directly obtain results to the CSS. We modeled the system with Pyomo¹⁰, a Python-based, open-source optimization modeling language, using Pyomo.DAE¹¹ to support automatic transcription (discretization) of the differential equations. The system is discretized and the derivatives are approximated in each time point by using an implicit Euler method. The solver IPOPT¹² solves the square system of nonlinear algebraic equations by employing a modified Newton method with a line search algorithm. Therefore, the solutions for concentrations of gas phase compounds and fractions of coverage are simultaneously found for all discrete points in the time horizon. A previous work also framed the problem as a BVP¹³, and the basic ideas of this limit-cycle conditions formulation was described in previous studies, either for dynamic catalysis¹⁴ or other applications¹⁵, but we implemented it in a flexible, easy to deploy and interpret manner.

The results with this simultaneous approach are obtained in around one second, independently of the stiffness of the system. With the sequential approach, forward integrating the ODEs in time, run times vary broadly with system stiffness and conditions, ranging between many orders of magnitude, going up to around 9500 seconds within our tests. The computational efficiency of the simultaneous approach allowed the implementation of derivative-free optimization methods wrapped around the simulation to obtain optimal parameters (within bounds) for the forcing signal that describes BE oscillations.

We adopted the square wave as the main working waveform, but also deployed the methodologies to the sine and triangle/sawtooth waves. Additionally, we used a smoothed version of the square wave, in which we could adjust the asymmetry and smoothness, and a Gaussian pulse function in order to somewhat approximate the square wave with a continuous function. For the continuous forcing functions, it was possible to implement gradient-based optimization methods and also to modify the Pyomo/IPOPT framework to turn the wave parameters into Pyomo decision variables (degrees of freedom) rather than parameters to solve the optimization problem of maximizing the avTOF using IPOPT.

For the square wave, an increase in four orders of magnitude on the avTOF when compared to the peak of the Sabatier of the static system was verified. This result demonstrated not only the potential of using dynamic catalysis, but the value of using optimization techniques to identify sets of parameters and achieve the best results possible within bounds. This work serves as proof of concept for the implementation of the methodologies discussed here in real, complex reactions on the surfaces of real materials and with input data and bounds related to the specific dynamic catalytic system and stimulus source. Future computational work will focus on adding physical constraints to the model to improve the accuracy between models and experimental results and on methods for acquiring arbitrary waveforms, in opposition to calculating the ideal parameters of fixed waveforms.

[1] M. A. Ardagh, O. A. Abdelrahman, and P. J. Dauenhauer, “Principles of Dynamic Heterogeneous Catalysis: Surface Resonance and Turnover Frequency Response,” ACS Catal., vol. 9, no. 8, pp. 6929–6937, Aug. 2019, doi: 10.1021/acscatal.9b01606.

[2] G. R. Wittreich, S. Liu, P. J. Dauenhauer, and D. G. Vlachos, “Catalytic resonance of ammonia synthesis by simulated dynamic ruthenium crystal strain,” Sci. Adv., vol. 8, no. 4, p. eabl6576, Jan. 2022, doi: 10.1126/sciadv.abl6576.

[3] M. A. Murphy, S. R. Gathmann, C. J. Bartel, O. A. Abdelrahman, and P. J. Dauenhauer, “Catalytic resonance theory: Circumfluence of programmable catalytic loops,” J. Catal., vol. 430, p. 115343, Feb. 2024, doi: 10.1016/j.jcat.2024.115343.

[4] S. Gathmann, C. Bartel, L. Grabow, O. Abdelrahman, C. D. Frisbie, and P. Dauenhauer, “Dynamic Promotion of the Oxygen Evolution Reaction via Programmable Metal Oxides,” Chemistry, preprint, Jan. 2024. doi: 10.26434/chemrxiv-2024-gs6zn.

[5] V. V. R. Vempatti, S. Wang, O. A. Abdelrahman, P. J. Dauenhauer, and L. C. Grabow, “Catalytic Resonance of Methane Steam Reforming by Dynamically Applied Charges,” Chemistry, preprint, Jan. 2024. doi: 10.26434/chemrxiv-2024-8sxpx.

[6] J. Gopeesingh, M. A. Ardagh, M. Shetty, S. T. Burke, P. J. Dauenhauer, and O. A. Abdelrahman, “Resonance-Promoted Formic Acid Oxidation via Dynamic Electrocatalytic Modulation,” ACS Catal., vol. 10, no. 17, pp. 9932–9942, Sep. 2020, doi: 10.1021/acscatal.0c02201.

[7] J. Qi et al., “Dynamic Control of Elementary Step Energetics via Pulsed Illumination Enhances Photocatalysis on Metal Nanoparticles,” ACS Energy Lett., vol. 5, no. 11, pp. 3518–3525, Nov. 2020, doi: 10.1021/acsenergylett.0c01978.

[8] M. A. Ardagh, T. Birol, Q. Zhang, O. A. Abdelrahman, and P. J. Dauenhauer, “Catalytic resonance theory: superVolcanoes, catalytic molecular pumps, and oscillatory steady state,” Catal. Sci. Technol., vol. 9, no. 18, pp. 5058–5076, 2019, doi: 10.1039/C9CY01543D.

[9] S. R. Gathmann, M. A. Ardagh, and P. J. Dauenhauer, “Catalytic resonance theory: Negative dynamic surfaces for programmable catalysts,” Chem Catal., vol. 2, no. 1, pp. 140–163, Jan. 2022, doi: 10.1016/j.checat.2021.12.006.

[10] M. L. Bynum et al., Pyomo — Optimization Modeling in Python, vol. 67. in Springer Optimization and Its Applications, vol. 67. Cham: Springer International Publishing, 2021. doi: 10.1007/978-3-030-68928-5.

[11]B. Nicholson, J. D. Siirola, J.-P. Watson, V. M. Zavala, and L. T. Biegler, “pyomo.dae: a modeling and automatic discretization framework for optimization with differential and algebraic equations,” Math. Program. Comput., vol. 10, no. 2, pp. 187–223, Jun. 2018, doi: 10.1007/s12532-017-0127-0.

[12] A. Wächter and L. T. Biegler, “On the implementation of an interior-point filter line-search algorithm for large-scale nonlinear programming,” Math. Program., vol. 106, no. 1, pp. 25–57, Mar. 2006, doi: 10.1007/s10107-004-0559-y.

[13] Y. M. Psarellis, M. Kavousanakis, P. J. Dauenhauer, and I. G. Kevrekidis, “Writing the Programs of Programmable Catalysis,” Feb. 2023, doi: 10.26434/chemrxiv-2023-99wqg.

[14] B. Foley and N. Razdan, “Dynamic Catalysis Fundamentals: I. Fast calculation of limit cycles in dynamic catalysis,” Chemistry, preprint, Nov. 2021. doi: 10.26434/chemrxiv-2021-10hk4.

[15] Y. Kawajiri and L. T. Biegler, “Optimization strategies for simulated moving bed and PowerFeed processes,” AIChE J., vol. 52, no. 4, pp. 1343–1350, Apr. 2006, doi: 10.1002/aic.10736.