(118e) Stability-Guaranteed Inference of Reduced-Order Models with Control of Power-to-X Processes

Introduction
Industry 4.0 and sustainability demands have shifted industrial operations from steady-state to dynamic control systems. Differential equation-based models have traditionally been used to understand dynamic processes. However, they face challenges when dealing with uncertain parameters, large state-spaces, and nonlinearities. These limitations hinder outer-loop computations such as uncertainty estimation, optimization, and control (Biegler et al., 2014). To overcome these challenges and maintain computational efficiency and accuracy, it is necessary to develop compact dynamic system models (Benner et al., 2021).

Model order reduction (MOR) has emerged as a strategy to obtain reduced-order models (ROMs) that simplify complex models (Bremer et al., 2017). Operator inference (OpInf) (Peherstorfer and Willcox, 2016), a non-intrusive method, constructs ROMs solely from data. These methods are particularly beneficial in scenarios where the original model is either inaccessible and only available as a black box, or when there is a wealth of low-noise data from the process itself. However, these inferred models are not guaranteed to be stable. Hence, we focus on the framework proposed by Goyal et al. (2023) and Pontes et al. (2024), which allows the inference of dynamical models that are guaranteed to be stable.

Power-to-X (PtX) exemplifies a dynamic process that converts renewable energy into green chemicals. The economic viability of PtX operations depends on their ability to respond dynamically to fluctuations in the availability of renewable energy sources. Within the PtX process, our primary focus is on the CO₂ methanation. Previous studies (Peterson et al., 2023, Gosea et al., 2024) have shown that OpInf can accurately represent single trajectories describing the CO₂ methanation. For these ROMs to be used effectively in control scenarios, they must be able to predict trajectories as a function of control parameters (Bremer et al., 2021). By integrating control into OpInf, we can develop ROMs that are tailored for operation under varying conditions.

Model Order Reduction via Operator Inference
OpInf fits structured models in reduced coordinates using snapshot data. The approach involves data acquisition, low-dimensional data representation, and model derivation.

During data acquisition, we collect discretized state data, time derivatives and input data (control values) for the full order model (FOM). These data are organized column by column into three matrices: the state snapshot matrix 𝑸∈ℝ^{𝑛×𝑘}, the derivative matrix 𝑸̇∈ℝ^{𝑛×𝑘} and the input matrix 𝑼∈ℝ^{𝑚×𝑘}. Here, n is the dimension of the state, m is the dimension of the input and k is the number of datapoints. If the time derivatives cannot be computed directly, they are inferred from state snapshots using numerical schemes.

We use Principal Component Analysis (PCA) via Singular Value Decomposition (SVD) on the state snapshot matrix 𝑸 to decompose it into left and right singular vectors, and singular values. Using the dominant left singular vectors, we obtain a low-dimensional representation of the data. Precisely, the basis matrix 𝑽_r∈ℝ^{𝑛×𝑟}, formed by selecting the first r columns of V (the left singular vectors), represents the state as a linear combination of r orthonormal vectors. Using 𝑽_r, we construct the reduced states 𝒒̂ and their derivatives 𝒒̂̇ as per 𝒒̂(𝑡)≈𝑽_r^𝑇x(t). In essence, the reduced states are approximations obtained by projecting the full states onto the subspace defined by the column space of 𝑽_r.

Using the reduced states, we first determine the ODE form and then derive the ROM. By choosing vector-valued functions, called operators, that have polynomial properties with respect to the reduced state and the input, we typically choose a quadratic polynomial due to the prevalence of quadratic systems in physics. The ROM is formulated using the reduced operators 𝐜̂∈ℝ^𝑟 (constant operator), 𝐀̂∈ℝ^{𝑟×𝑟} (linear operator), 𝐇̂∈ℝ^{𝑟×𝑟**2} (quadratic operator), and 𝐁̂∈ℝ^{𝑟×m} (input operator) as follows:

𝒒̂̇(𝑡)=𝒄̂+𝐀̂ 𝒒̂(𝑡)+𝐇̂(𝒒̂(𝑡)⊗𝒒̂(𝑡))+𝐁̂𝑢, 𝒒̂(0)=𝒒̂₀ (1)

Next, we learn the ROM operators, which is theoretically possible using the least squares method. Since the problem is often ill-conditioned and the stability of the ROM is not guaranteed, we follow the parameterization guidelines by (Goyal et al., 2023, Pontes et al., 2024).

Application of Interest
This methanation reaction, which is exothermic, produces heat and leads to the formation of hotspots, which in turn have a significant impact on the efficiency of the reactor. Therefore, management of these hotspots becomes critical, especially in the context of fluctuating inputs.

We use a one-dimensional mechanistic reactor model based on Zimmermann et al. (2022). This model consists of a system of coupled PDEs representing an energy balance (expressed in terms of temperature, 𝑇) and a mass balance (expressed in terms of CO₂ conversion, 𝑋) incorporating constants, parameters, state equations, and boundary conditions. The model’s parameters are based on an industrial fixed bed reactor given by Bremer (2020), with the coolant temperature (𝑇_cool) as a controlled variable. The governing equations for 𝑋 and 𝑇 over the axial coordinate (𝑧) and time (𝑡) are given in Figure 1a. In our study of 𝑇_cool variations, both equations are first solved for 𝑇_cool=600 K. The resulting steady-state profile serves as the initial condition for subsequent trajectory calculations.

Numerical Results
The mechanistic model’s equations are discretized into 1000 control volumes using Python’s diffrax library and the finite volume method. A total of 12 trajectories (with ten for training and two for testing) are generated by varying 𝑇_cool between 500 K and 750 K. These trajectories are then compiled into the state snapshot matrix 𝑸. All states are shifted by the mean steady-state of the training set.

PCA via SVD transforms the state snapshots into a lower-dimensional representation. To capture the dynamics of the system, we derive ROMs over r dimensions with the goal of capturing between 99.9 % and 99.99 % of the total energy present in the data. Figure 1b illustrates the total energy among the singular values.

The reduced operators 𝐀̂, 𝐇̂ and 𝐁̂ are derived using the Adam optimization algorithm from PyTorch, with a CyclicLR scheduler. We apply stability parameterization and add a regularization factor of 10^-4 to 𝐇̂. The inferred ROMs are simulated numerically, compared to the FOM, and the model with the lowest Frobenius norm is selected.

In our study, the ROM capturing 99.5 % of the signal energy fits the training data best. This model, with a rank of 𝑟=10 (𝑟_X=5 and 𝑟_T=5), has an average relative Frobenius norm of 2.45 % during training. We validate the model’s ability to reproduce the system dynamics using two unseen trajectories, achieving relative Frobenius norms of 2.34 % and 3.78 %. Figure 1c and 1d show three-dimensional plots of the original data for the second test trajectory and the OpInf results, respectively. The computational complexity of solving the ROM is proportional to the number of degrees of freedom in the low-dimensional representation of the state. As a result, the ROM solves the initial value problem in 24% of the computational time required by the FOM.

Conclusion
Our study shows that a stable ROM with quadratic nonlinearity accurately captures system dynamics and predicts trajectories for unknown control parameters. Future works will potentially cover the implementation of the ‘lift-and-learn’ approach and parameterize OpInf to improve the efficiency and adaptability of predicting dynamic systems. With these advances, OpInf will develop ROMs that can be used in computationally intensive computations, thereby increasing system efficiency and adaptability under varying conditions.

References

Benner, P., Breiten, T., Faßbender, H., Hinze, M., Stykel, T., & Zimmermann, R. (Eds.). (2021). Model reduction of complex dynamical systems. Springer International Publishing.
Biegler, L. T., Lang, Y. D., & Lin, W. (2014). Multi-scale optimization for process systems engineering. Computers & Chemical Engineering, 60, 17-30.
Bremer, J., Goyal, P., Feng, L., Benner, P., & Sundmacher, K. (2017). POD-DEIM for efficient reduction of a dynamic 2D catalytic reactor model. Computers & Chemical Engineering, 106, 777-784.
Bremer, J. (2020). Advanced Operating Strategies for Non-Isothermal Fixed-Bed Reactors Exemplified for CO₂ Methanation (Doctoral dissertation, Magdeburg: Otto-von-Guericke-Universität).
Bremer, J., Heiland, J., Benner, P., & Sundmacher, K. (2021). Non-intrusive Time-POD for Optimal Control of a Fixed-Bed Reactor for CO₂ Methanation. IFAC-PapersOnLine, 54(3), 122-127.
Duff, I. P., Goyal, P., & Benner, P. (2024). Stability-Certified Learning of Control Systems with Quadratic Nonlinearities. arXiv preprint arXiv:2403.00646.
Gosea, I. V., Peterson, L., Goyal, P., Bremer, J., Sundmacher, K., & Benner, P. (2024). Learning reduced-order Quadratic-Linear models in Process Engineering using Operator Inference. arXiv preprint arXiv:2402.17698.
Goyal, P., Duff, I. P., & Benner, P. (2023). Guaranteed stable quadratic models and their applications in SINDy and operator inference. arXiv preprint arXiv:2308.13819.
Peherstorfer, B., & Willcox, K. (2016). Data-driven operator inference for nonintrusive projection-based model reduction. Computer Methods in Applied Mechanics and Engineering, 306, 196-215.
Peterson, L., Goyal, P., Gosea, I. V., Bremer, J., Benner, P. & Sundmacher, K. (2023). Learning reduced-order models for dynamic CO₂ methanation using operator inference. Max Planck Institute Magdeburg, Preprint, accepted at the Proceedings of the 34th ESCAPE.
Zimmermann, R. T., Bremer, J., & Sundmacher, K. (2022). Load-flexible fixed-bed reactors by multi-period design optimization. Chemical Engineering Journal, 428, 130771.