2.1. Chemical kinetics identification through a deep operator network

This section discusses the proposed framework for the physics-informed enforcement of mass and element conservation. As stated earlier, we integrate the framework with our recently developed React-DeepONet framework for chemistry acceleration. The framework is designed to predict the solution of a subset of the thermo-chemical scalars’ vectors between two prescribed time increments. Therefore, the framework is designed to advance the chemical source operator that can be coupled with transport in a CFD simulation. We illustrate the framework using data for a 0D, homogeneous mixture with prescribed pressure, composition (through an equivalence ratio), and initial temperature. The corresponding system of ODEs, which are solved as an initial value problem, takes the form:

(2.1) $$ \frac{d{\phi}_i}{dt}=f(\phi, \hskip.4em {p}_0),\hskip1em i=1,2,3,\dots, n $$

(2.2) $$ \phi (0)={\phi}_0 $$

where $ \phi $ is the thermo-chemical composition vector, which includes species mass fractions and temperature, the chemical source term $ f(\phi, \hskip.4em {p}_0) $ is a function of thermo-chemical scalars, $ \phi $ , and is typically modeled using the law of mass action. The chemical kinetics presented by the above system of stiff ODEs is learned through the React-DeepONet framework of Kumar and Echekki (Reference Kumar and Echekki2024), where a reduced set of representative species is considered instead of the full thermochemical solution vector. This choice yields significant savings in computational time and memory, but accounting for a partial solution of the chemical system must also address mass and element conservation.

The dimensionality of the reduced composition space ( $ {\phi}_i,i=1,2,\dots m\in {R}^m $ ) is far smaller than the dimensionality of the original solution vector ( $ {\phi}_i,i=1,2,\dots m,m+1,m+2,\dots n\in {R}^n $ ). In the Kumar and Echekki study (2024), the React-DeepONet that was used to capture the two-stage ignition process of n-dodecane at various initial conditions was implemented with seven representative species and temperature (i.e., $ m= $ 8), and the full mechanism required the integration of 348 species and temperature (i.e., $ n= $ 349).

It is important to note that for implementing the proposed approach in CFD calculations, the training data should be suitable for the targeted composition space in CFD. This data can be generated using 0D (Nguyen et al., Reference Nguyen, Domingo, Vervisch and Nguyen2021; Mehl and Aubagnac-Karkar, Reference Mehl and Aubagnac-Karkar2023) or 1D (e.g., 1D flamelets) reactors. Simpler 2D or 3D direct numerical simulations or a reduced scale or dimension represent alternatives to construct the targeted composition space (Owoyele and Echekki, Reference Owoyele and Echekki2017; Kumar et al., Reference Kumar, Rieth, Owoyele, Chen and Echekki2023).

The React-DeepONet architecture adopted by Kumar and Echekki (Reference Kumar and Echekki2024) is shown in Figure 1. It consists of two subnetworks, the Branch Net and the Trunk Net, which include as inputs the solution vector for the representative species and temperature at a time, $ t $ , and a prescribed time step, $ \delta t $ , on the Branch and Trunk nets, respectively. The outputs of the two subnetworks are combined in a concatenation layer, resulting in output for the solution vector for the representative species and temperature at a later time increment, $ t+\delta t $ . Based on its architecture, React-DeepONet maps the solution vector of a set of representative species and temperature at a time of the evolution of this solution across a time step $ \delta t $ , thus avoiding the costly direct chemistry integration and enabling the time step, $ \delta t $ , to be prescribed primarily by transport time scale requirements in reacting flow simulations. To train the React-DeepONet, the following loss function is evaluated, which compares the React-DeepONet predictions to the training data, which corresponds to the numerically integrated solution of the full chemical system:

(2.3) $$ {\mathcal{L}}_{\mathrm{DATA}}=\frac{1}{m}\sum \limits_{i=1}^m{({\phi}_i(t+\delta t)-{\phi}_i^{\ast }(t+\delta t))}^2 $$

where $ \delta t $ varies from $ 0 $ to a large fixed time step $ \Delta t $ and $ {\phi}_i^{\ast}\left(t+\delta t\right) $ is the corresponding ground-truth value. Therefore, React-DeepONet accommodates training for a range of values for $ \delta t $ within $ \Delta t $ . During prediction, any value of $ \delta t $ within this range can be used or adaptively changed during the simulation. The training prescribed above acts as a preconditioner for the model parameters, and additional refined training is implemented based on an autoregressive procedure (Kumar and Echekki, Reference Kumar and Echekki2024. Within the process, the refined model is trained to predict the solution vector at multiple time steps, $ {N}_l $ . Accordingly, the loss function of the refined training model combines the loss functions corresponding to steps from 1 to N $ {}_l $ .

Figure 1. React-DeepONet: ANN for predicting the evolution of representative species and temperature.

2.2. Mapping nonrepresentative species from the reduced set: Correlation Net

The objective of the proposed approach is twofold: (1) to constrain the React-DeepONet predictions to satisfy mass and element conservation and (2) to predict nonrepresentative species, which are not predicted using React-DeepONet or transported in CFD. These species, combined with the representative species, may account for the bulk of the mass in the mixture. Although the dimensionality of the full solution vector is very large, it resides on a low-dimensional manifold since most of the scalars are highly correlated (Pope, Reference Pope2013; Coussement et al., Reference Coussement, Isaac, Gicquel and Parente2016). The reduced set ( $ {\phi}_i,i=1,2,\dots m $ ), identified based on principal component analysis (PCA) of reaction rate data as proposed by Alqahtani and Echekki (Reference Alqahtani and Echekki2021), is representative of the thermo-chemical manifold. It strongly correlates with the nonrepresentative species ( $ {\phi}_i,i=m+1,m+2,\dots n $ ). However, other methods for selecting the representative species can be used including directed relations graphs (DRGs) methods (Lu and Law, Reference Lu and Law2005), the global pathway selection (GPS) algorithm (Gao et al., Reference Gao, Yang and Sun2016), or the more recent manifold-informed state vector subset (Zdybał et al., Reference Zdybał, Sutherland and Parente2023). Recently, we showed that the PCA-based technique and the other three techniques yield similar representative species for chemical systems (Gitushi and Echekki, Reference Gitushi and Echekki2024).

Leveraging their high correlation, nonrepresentative species can be mapped from the reduced set of representative species using an ANN denoted as Correlation Net as shown in Figure 2. The training of the ANNs is implemented using training data from simulations that capture the accessed composition space of the target simulations in CFD. Inputs and outputs to Correlation Net are the representative species and temperature and nonrepresentative species, respectively, as shown in Figure 2. The loss that is evaluated during training is expressed as

(2.4) $$ {\mathrm{\mathcal{L}}}_{\mathrm{C}-\mathrm{DATA}}=\frac{1}{\left(m-n\right)}\sum \limits_{i=m+1}^n{\left({\phi}_i-{\phi}_i^{\ast}\right)}^2 $$

Figure 2. Correlation Net: ANN for recovering nonrepresentative species.

Again, the loss compares predictions by Correlation Net with values obtained from the training data and resulting from the direct integration of the full thermo-chemical scalars’ solution vector. Since the full thermo-chemical space obeys certain physical laws, such as the sum of species mass fractions is one and elements remain conserved during the reaction as the solution vector evolves, Correlation Net is trained with additional losses corresponding to these physical constraints.

(2.5) $$ {\mathrm{\mathcal{L}}}_{\mathrm{C}-\mathrm{SUM}-\mathrm{MF}}=\left|\sum \limits_{i=2}^n{\phi}_i-1.0\right| $$

(2.6) $$ {\mathrm{\mathcal{L}}}_{\mathrm{C}-\mathrm{ELEM}}=\sum \limits_{K=1}^{N_{\mathrm{Element}}}\left|\sum \limits_{i=2}^n\frac{N_K^i{W}_K}{W_i}\left({\phi}_i-{\phi}_i^{\ast}\right)\right|, $$

where $ {N}_K^i $ is the Kth element contribution in ith species, $ {W}_K $ is the atomic weight of Kth element, and $ {W}_i $ is the molecular weight of the ith species. Also, $ {\phi}_1 $ is essentially temperature and is excluded from the conservation loss formulations. It should be noted that these two losses include the Correlation Net output and its input. The total loss function for Correlation Net takes this form:

(2.7) $$ {\mathrm{\mathcal{L}}}_C={W}_1\;{\mathrm{\mathcal{L}}}_{\mathrm{C}-\mathrm{DATA}}+{W}_2\;{\mathrm{\mathcal{L}}}_{\mathrm{C}-\mathrm{SUM}-\mathrm{MF}}+{W}_3\;{\mathrm{\mathcal{L}}}_{\mathrm{C}-\mathrm{ELEM}}, $$

which combines losses associated with correlations and mass and element conservation. In this expression, $ {W}_1 $ , $ {W}_2 $ , and $ {W}_3 $ are the weights for each loss term. The choice of the weights is dictated by the magnitudes of the various loss terms within Eq. (2.7) as discussed below. Although mass and element conservation may be redundant for reactions, they represent two useful constraints to enforce during the Correlation Net training. The results discussed below will investigate whether implementing both constraints improves the training of Correlation Net.

It is important to note that alternative approaches have been developed in the literature to reinforce mass and element conservation within the framework of ANN-based predictions of chemistry. Readshaw et al. (Reference Readshaw, Jones and Rigopoulos2023) implemented mass conservation in different steps. In the first step, the species are predicted by training an ANN without mass conservation constraints. The subsequent steps are implemented outside the ANN predictions. In these steps, residual species representing the main contributions to the elements are determined based on the ANN predictions. Then, the mass fractions of these species are recalculated to conserve their corresponding elements across the integration time step. An additional step is implemented if the new value of the residual species is negative by setting its value to zero and enforcing element conservation using two rectifying species.

Mehl and Aubagnac-Karkar (Reference Mehl and Aubagnac-Karkar2023) also adopted a correction scheme for element conservation on selected species. For H₂ oxidation, H and O elements are corrected by adding to or removing mass from H₂ and O₂ from the ANN predictions. They also implemented an element-based hybrid model to replace the ANN evaluation with direct integration predictions when an error proxy exceeded a threshold value. This error proxy is based on the relative errors achieved by the predictions of H and O elements’ mass fractions during the reaction.

Rohrhofer et al. (Reference Rohrhofer, Posch, Gößnitzer, García-Oliver and Geiger2024) implemented mass conservation within the network architecture by adding a softmax or a Sigmoid activation function to the last layer of the network to predict species from an input represented by a mixture fraction, a progress variable, and a dissipation rate using flamelet libraries. In contrast to the above studies, the present study implements the physical constraints of mass and element conservation through the loss function. While no additional constraints are imposed beyond the training process, the present results indicate that imposing these constraints can still significantly improve the predictions of our solution mapping approach. More importantly, it provides, within the same loss formulation, the option to weigh in the species predictions versus their mass/element constraints differently.

2.3. Physics-constrained chemical kinetics identification

As stated earlier, the prediction of the representative species and temperature is implemented through a DeepONet. The nonrepresentative solution vector can be reconstructed through the pre-trained Correlation Net for any reduced solution vector. The mapping from Correlation Net extends the reduced set to the full solution vector as shown in Figure 3. This figure illustrates the prediction of the entire thermo-chemical state solution vector using React-DeepONet for temperature and representative species and augmented with a pretrained correlation net to recover nonrepresentative species while enforcing conservation principles through loss evaluations. Utilizing the mapping provided by pretrained Correlation Net, the training process of React-DeepONet can now be enriched with the physical constraints of species and elemental mass conservation, as formulated in Section 2, with the corresponding losses defined:

(2.8) $$ {\mathrm{\mathcal{L}}}_{\mathrm{SUM}-\mathrm{MF}}=\left|\sum \limits_{i=2}^n{\phi}_i\left(t+\delta t\right)-1.0\right| $$

(2.9) $$ {\mathrm{\mathcal{L}}}_{\mathrm{ELEM}}=\sum \limits_{K=1}^{N_{\mathrm{Element}}}\left|\sum \limits_{i=2}^n\frac{N_K^i{W}_K}{W_i}\left({\phi}_i\left(t+\delta t\right)-{\phi}_i^{\ast}\left(t+\delta t\right)\right)\right| $$

in addition to the data loss that is defined as

(2.10) $$ {\mathrm{\mathcal{L}}}_{\mathrm{DATA}}=\frac{1}{m}\sum \limits_{i=1}^m{\left({\phi}_i\left(t+\delta t\right)-{\phi}_i^{\ast}\left(t+\delta t\right)\right)}^2 $$

Figure 3. Full ANN coupling of React-DeepONet and Correlation Net.

The total loss becomes

(2.11) $$ \mathrm{\mathcal{L}}={W}_1{\mathrm{\mathcal{L}}}_{\mathrm{DATA}}+{W}_2{\mathrm{\mathcal{L}}}_{\mathrm{SUM}-\mathrm{MF}}+{W}_3{\mathrm{\mathcal{L}}}_{\mathrm{ELEM}} $$

where $ {W}_1 $ , $ {W}_2 $ , and $ {W}_3 $ are the weights for each loss term. Again, element and mass conservation as in Eq. (2.7) may be redundant. Therefore, enforcing element conservation loss alone is sufficient to establish mass and element conservation, and either weight, $ {W}_1 $ or $ {W}_2 $ , can be set to zero. Although, preferably, $ {W}_2 $ should be kept and $ {W}_1 $ can be set to zero. We will again investigate whether enforcing this redundancy in mass and element conservation can improve the training in the Results section.

Finally, we note that the solution vector in the network shown in Figure 3 is subjected to both a transformation and a normalization. The main exception is where we implement the loss function evaluation, where we need to transform the species mass fractions to their actual values to impose mass and element conservation. The transformation and scaling procedure is implemented as follows:

• • A fifth-order transformation is implemented on the species profiles and temperature to accommodate their strong variations. For example, for the $ k $ th species mass fraction, $ {Y}_k $ , this transformation is implemented as follows:

  (2.12) $$ {y}_k={Y}_k^{1/5} $$

• • In a second step, we normalize the transformed variables, $ {y}_k $ , to yield values in the range of −1 to 1 among all the training data by using the minimum and maximum values of each transformed quantity.