Multiscale optimization of formic acid dehydrogenation process via linear model decision tree surrogates
Ethan M. Sunshine, Giovanna Bucci, Tanusree Chatterjee, Shyam Deo, Victoria M. Ehlinger, Wenqin Li, Thomas Moore, Corey Myers, Wenyu Sun, Bo-Xun Wang, Mengyao Yuan, John R. Kitchin, Carl D. Laird, Matthew J. McNenly, Sneha A. Akhade, “Multiscale optimization of formic acid dehydrogenation process via linear model decision tree surrogates”, Computers & Chemical Engineering 194, 108921 (2025). https://doi.org/10.1016/j.compchemeng.2024.108921
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
Simultaneously optimizing models of a process at several scales is known to find better designs than optimizing each scale on its own, but it requires every model to be written as algebraic constraints. Atomistic models, such as Kohn-Sham density functional theory (DFT) feeding a microkinetic model, are iterative calculations that cannot be dropped into an optimization problem directly. The paper asks whether data-driven surrogates with a functional form that converts exactly to a mixed-integer linear program (MIP) can stand in for such models, so that a catalyst, a reactor and a techno-economic cost model can be optimized together (Sunshine et al. 2025).
Approach
The demonstration is the release step of a liquid organic hydrogen carrier (LOHC) supply chain: dehydrogenation of formic acid to hydrogen and CO2 at or near a hydrogen customer such as a refueling station. Three underlying models, built by domain experts with different tools, are linked through shared cross-scale variables. The catalyst model predicts turnover frequency (TOF) for Pd, Cu or Pt nanoparticles from particle size, temperature and pressure, using DFT surface energies, a Wulff construction, site counting (including edges and corners) and microkinetic modeling. The reactor model is a Python CSTR with mixer, flash separator and recycle (converged with a SciPy tear-stream iteration) that returns hydrogen yield. The costing model is a bottom-up techno-economic analysis giving the levelized cost of hydrogen in $/kg H2.
Each underlying model is treated as a black box that generates data, and a linear model decision tree (LMDT), a regression tree with a linear fit in each leaf, is trained on it with the linear-tree package. The reactor and costing models each labeled 100,000 uniformly random points; the catalyst model was sampled on a grid of 20 temperatures (300 to 400 K) and particle sizes of 3, 4 and 5 nm, giving 3600 points; all used an 80/20 train/validation split. The trees were embedded in Pyomo via OMLT as generalized disjunctive programs, with the metal choice added as a further disjunction (spot metal prices of 38,000, 8 and 28,000 $/kg for Pd, Cu and Pt, and a fixed 0.3% volumetric loading). A convex hull transformation turned the problem into a MIP solved with GLPK in a few seconds. Because unbounded leaf regressions let the reactor surrogate predict yields above 100%, the yield passed to the cost model was constrained to be at most the surrogate output and at most 100%. Reactor volume was also limited to station capacity divided by 500.
Key findings
The surrogates had 439 leaves (catalyst, depth 12), 78 leaves (reactor, depth 13) and 245 leaves (costing, depth 11); the reactor surrogate was the least accurate, but the points it placed above 100% yield actually lay between 80 and 100%. A sequential “one-model” strategy fixed the catalyst for maximum activity (Pt, 3 nm, 400 K, 100 atm, TOF 37.0 Hz) and the reactor for maximum yield, then minimized cost, giving 23.32 $/kg H2. Fixing only the catalyst (two-model) gave 16.20 $/kg H2 by shrinking reactor volume to 0.39 m3 and lowering pressure. Optimizing all three models together selected a far less active but cheap Cu catalyst (5 nm, 300 K, 100 bar, TOF 0.0191 Hz) and reached 14.01 $/kg H2, about 40% below the sequential design. Station capacity went to its 1500 kg H2/day upper bound in every case, and all designs reached 100% yield, so pellet effectiveness factor and reaction order were not pushed to their bounds.
Sweeping station capacity with the metal fixed reproduced an economy of scale. Only Pt was feasible at small capacities because of the reactor volume limit and its high loading; Pd became feasible at larger capacity but was never cheapest, and Cu was the lowest-cost choice once it became feasible. At high capacity the cost fell linearly, an extrapolation artifact that would eventually predict negative costs.
Significance
The work shows a general recipe for multiscale design: wrap each model, whatever its software, in a MIP-representable surrogate and solve the coupled problem to global optimality, which here uncovered a counterintuitive catalyst choice. The authors list the limitations plainly: training data must be generated for every model (active learning could help), extrapolation is unreliable, optima from discontinuous LMDTs tend to sit on leaf boundaries where errors are largest, so uncertainty quantification and error propagation across linked surrogates are needed, and binary variables grow with tree number and depth (Sunshine et al. 2025). Code and data are in the LLNL systems2atoms repository.