Reaction Mechanism Reduction for Ozone-Enhanced CH4/Air Combustion by a Combination of Directed Relation Graph with Error Propagation, Sensitivity Analysis and Quasi-Steady State Assumption

In this study, an 18-steps, 22-species reduced global mechanism for ozone-enhanced CH4/air combustion processes was derived by coupling GRI-Mech 3.0 and a sub-mechanism for ozone decomposition. Three methods, namely, direct relation graphics with error propagation, (DRGRP), sensitivity analysis (SA), and quasi-steady-state assumption (QSSA), were used to downsize the detailed mechanism to the global mechanism. The verification of the accuracy of the skeletal mechanism in predicting the laminar flame speeds and distribution of the critical components showed that that the major species and the laminar flame speeds are well predicted by the skeletal mechanism. However, the pollutant NO was predicated inaccurately due to the precursors for generating NO were removed as redundant components. The laminar flame speeds calculated by the global mechanism fit the experimental data well. The comparisons of simulated results between the detailed mechanism and global mechanism were investigated and showed that the global mechanism could accurately predict the major and intermediate species and significantly reduced the time cost by 72%.


Introduction
With the development of advanced measurement technology and the research milestones on combustion processes, the size of detailed mechanisms has increased significantly.A mechanism for a hydrocarbon combustion process contains hundreds of species and thousands of reactions, like USC-Mech II.For instance, Qi et al. [1] developed a detailed mechanism of ethylbenzene and toluene including 176 species and 804 reactions by employing synchrotron vacuum ultraviolet photoionization mass spectrometry to detect the flame species.For the simplest hydrocarbon fuel methane, the well-known mechanism GRI-Mech 3.0 [2] consists of 53 species and 325 reactions.When simulating such a reaction process using computational fluid dynamics technique, an equation with convection, diffusion, and source terms must be introduced to describe the mass-fraction variance of each species [3,4].The use of a detailed mechanism rapidly increases the computational cost [5,6].Furthermore, the time scales of the elementary reactions extensively vary, thereby increasing the stiffness of the equation set.Consequently, employing a detailed mechanism for the simulation of a combustion process at a laboratory or pilot scale is considerably challenging for researchers.
Mechanism reduction must be employed to minimize the scale (i.e., the amount of species and reactions) of a detailed mechanism and overcome the aforementioned problem [7].The constructed reduced mechanism should maintain the main characteristics of the detailed one, that is, the reduced mechanism should be able to predict not only the mass-fraction distributions of major, minor, and even some trace species but also the important combustion characteristics, such as ignition delay time, and laminar flame speed, within acceptable limits of error.A decrease in scale can reduce computational cost and stiffness and enable researchers to apply the mechanism to simulation investigations in a laboratory, pilot, and even large-scale facility.For example, Lv et al. [8] constructed a reduced mechanism by sensitivity analysis and quasi-steady-state assumption for an selective non-catalytic reduction (SNCR) removal process and applied it to a 100 MW utility boiler.
The methods for mechanism reduction are divided into three kinds: (1) based on the removal of redundant components; (2) based on the removal of redundant reactions; (3) based on the decrease in the mechanism stiffness.The first kind of method evaluates the importance of species in the entire reaction process.With the removal of the less important species, a simple reduced mechanism can be constructed.Methods falling under the first kind include principal component analysis (PCA) [9], directed relation graph (DRG) method [10] and its improvement DRG with error propagation (DRGEP) [11], path flux analysis [12], and simulation error minimization connectivity method [13].The second kind of method handles reactions instead of species.All elementary reactions can be assessed according to their contributions to the entire reaction process.Once the contribution of an elementary reaction is under the threshold value, such a reaction will be removed from the reduced mechanism.Methods falling under the second kind include sensitivity analysis (SA) [8,14], eigenvalue analysis [15], and PCA of the matrix F with simulation error minimization [13].The third kind sorts the elementary reactions by stiffness and removes the species and reactions simultaneously if their characteristic time scales are minimal in relation to the entire reaction process.Methods falling under their kind include computational singular perturbation [16], intrinsic low-dimensional manifold method [17], and quasi-steady-state assumption (QSSA) [8,18].
Under certain circumstances, the aforementioned methods can significantly simplify the simulation processes.However, the feasibility of mechanism reduction is dependent on a number of issues, including the experience of researchers and their understanding of reaction process.Furthermore, the removal of species and reactions from the detailed mechanism will reduce the optimal operating conditions, which may be challenging for researchers and thus require in-depth investigation.A general reference method for reducing mechanisms is therefore needed to be developed in computational fluid dynamics (CFD) research field, e.g., a combination of different mechanism reduction methods.
Nowadays, engines like the jet engine at a high altitude with low temperature and oxygen concentration are exposed to more and more extreme environments.Therefore, there is a growing interest of researchers and practical engineers in the development of technologies aimed at the enhancement of ignition and combustion in internal combustion engines [19][20][21].Ozone is one of the most favorited species used for combustion enhancement due to its highly reactive characteristics [22].Tachibana et al. [23] investigated the effect of 1000 ppm ozone addition on the combustion of three different fuels in the diesel engine and they revealed that the ozone can shorten the ignition delay, increase the fuel cetane number and lower the compression ratio of ignition limit.Masurier et al. [24] studied the effect of ozone addition on the combustion of six different primary reference fuels (PRFs) in a homogeneous charge compression ignition (HCCI) engine.Experiments confirmed that injection of ozone moved forward the phasing of six fuels and significantly improved the combustion process and the effect is higher for fuels which have the greater octane number.Ombrello et al. [25] experimentally and numerically investigated the thermal and kinetic effects of O3 on the C3H8 flame propagation.They concluded that the addition of O3 can enhance the early oxidation of fuel in the pre-heat zone and provide additional enhancement of the flame propagation speed at the flame front.However, the current investigation of ozone enhancement on combustion is mostly based on experiments or detailed kinetics, which are very difficult to directly apply to CFD.
A review of the previous literature reveals that a general reference method for reducing mechanism and a reduced ozone enhanced combustion mechanism are needed to be developed.Within this context, the objective of the present study is twofold: (1) to propose a general method of mechanism reduction for hydrocarbon combustion processes; (2) to construct an accurate and reliable global mechanism for ozone-enhanced CH 4 /air combustion processes.

Methods for the Development of Reduced Global Mechanisms
A mechanism reduction flowchart based on the characteristics of every method is presented in Figure 1.DRGEP is first applied to a detailed mechanism; as a result, a first-stage skeletal mechanism is generated by species removal.Subsequently, SA is conducted to form the final-stage skeletal mechanism by the removing the elementary reactions.Finally, QSSA method is used to generate a global mechanism.The methods included in the flowchart are presented in detail in the following subsections.
A review of the previous literature reveals that a general reference method for reducing mechanism and a reduced ozone enhanced combustion mechanism are needed to be developed.Within this context, the objective of the present study is twofold: (1) to propose a general method of mechanism reduction for hydrocarbon combustion processes; (2) to construct an accurate and reliable global mechanism for ozone-enhanced CH4/air combustion processes.

Methods for the Development of Reduced Global Mechanisms
A mechanism reduction flowchart based on the characteristics of every method is presented in Figure 1.DRGEP is first applied to a detailed mechanism; as a result, a first-stage skeletal mechanism is generated by species removal.Subsequently, SA is conducted to form the final-stage skeletal mechanism by the removing the elementary reactions.Finally, QSSA method is used to generate a global mechanism.The methods included in the flowchart are presented in detail in the following subsections.

Preprocessing by DRGEP Method
As mentioned in Section 1, the computational time cost of mechanism reduction should be considered.Most detailed mechanisms for hydrocarbon combustion processes contain thousands of elementary reactions.Thus, analysis every elementary reaction is significantly time consuming.Furthermore, all the methods based on the analysis of sensitivity coefficients are not suitable at the beginning of mechanism reduction, because it has the same cost as reaction analysis.To overcome this problem, Lu et al. [26] proposed the DRG method, which is based on production rate analysis; this method is also suitable for the preprocessing of detailed mechanisms.A relationship diagram can be constructed using the species production rates.The relationship coefficient between two species is defined by Equation (1) [26]: where subscript i denotes the i-th elementary reaction, ωi is the production rate, νA,i is the stoichiometric coefficient of species A, and δBi is the determination coefficient of species B. If i-th elementary reaction involves species B, it is equivalent to 1; otherwise, it is 0. Evidently, if the normalized relationship coefficient is sufficiently large, the removal of species B from the skeletal mechanism is expected to induce a significant error on the production rate of species A. Consequently, if species A is retained in the skeletal mechanism, species B should also be retained considering the strong interdependence between the two species.When rAB is less than the threshold value ε, the interdependence between species A and B are negligible.Under this condition, the removal of species B is reasonable.However, in a reaction process, not all species have a direct relationship with the object.For example, in Figure 2a, which shows a DRG diagram of four species (A, B, C, and D), species A has no direct relationship with species D, but has a strong relationship with species B. If rAB < rCD < rAC and the proposed threshold satisfies rAB < ε < rCD < rAC, then, based on the criterion of DRG, species B should be removed from the skeletal mechanism first followed by species D. Given that coefficient rCD cannot represent the relationship between species A and D, it is irrational to remove species B first without previously comparing the relationships of species A with species B and D. Pepiot-Desjardins et al. [27] proposed a DRGEP method to solve this problem.When the relationship coefficient is redefined as an R-value, the effect of all the species to the object can be evaluated in the detailed mechanism.R-value is defined by Equation (2) [27]:

Preprocessing by DRGEP Method
As mentioned in Section 1, the computational time cost of mechanism reduction should be considered.Most detailed mechanisms for hydrocarbon combustion processes contain thousands of elementary reactions.Thus, analysis every elementary reaction is significantly time consuming.Furthermore, all the methods based on the analysis of sensitivity coefficients are not suitable at the beginning of mechanism reduction, because it has the same cost as reaction analysis.To overcome this problem, Lu et al. [26] proposed the DRG method, which is based on production rate analysis; this method is also suitable for the preprocessing of detailed mechanisms.A relationship diagram can be constructed using the species production rates.The relationship coefficient between two species is defined by Equation (1) [26]: where subscript i denotes the i-th elementary reaction, ω i is the production rate, ν A,i is the stoichiometric coefficient of species A, and δ Bi is the determination coefficient of species B. If i-th elementary reaction involves species B, it is equivalent to 1; otherwise, it is 0. Evidently, if the normalized relationship coefficient is sufficiently large, the removal of species B from the skeletal mechanism is expected to induce a significant error on the production rate of species A. Consequently, if species A is retained in the skeletal mechanism, species B should also be retained considering the strong interdependence between the two species.When r AB is less than the threshold value ε, the interdependence between species A and B are negligible.Under this condition, the removal of species B is reasonable.However, in a reaction process, not all species have a direct relationship with the object.For example, in Figure 2a, which shows a DRG diagram of four species (A, B, C, and D), species A has no direct relationship with species D, but has a strong relationship with species B. If r AB < r CD < r AC and the proposed threshold satisfies r AB < ε < r CD < r AC , then, based on the criterion of DRG, species B should be removed from the skeletal mechanism first followed by species D. Given that coefficient r CD cannot represent the relationship between species A and D, it is irrational to remove species B first without previously comparing the relationships of species A with species B and D. Pepiot-Desjardins et al. [27] proposed a DRGEP method to solve this problem.When the relationship coefficient is redefined as an R-value, the effect of all the species to the object can be evaluated in the detailed mechanism.R-value is defined by Equation (2) [27]: where S is the set of all reaction paths leading from species A to B and r ij is calculated by Equation ( 1).
According to this definition, the relationship between species A and D are calculated by R AD = r AC × r CD .As shown in Figure 2b, by comparing R AB , R AC , and R AD , the removal order could be easily confirmed.
Energies 2018, 11, x FOR PEER REVIEW 4 of 12 where S is the set of all reaction paths leading from species A to B and rij is calculated by Equation ( 1).
According to this definition, the relationship between species A and D are calculated by RAD = rAC × rCD.As shown in Figure 2b, by comparing RAB, RAC, and RAD, the removal order could be easily confirmed.
(a) (b) The removal of species would change the R-value in the DRGEP diagram.Thus, re-computing the R-value after the removal of a certain number of species is highly recommended.In this study, such a re-computation is performed after every removal of ~10 species.

Advanced Processing by SA
Sensitivity analysis is a method used to investigate the effect of parameter changes on the solution of mathematical models [14].In chemical kinetics, models are usually based on differential equations, which derive concentration-time or concentration-space distribution, reaction rates, and various kinetic features of the reaction.After selecting the important species and kinetic features, the importance of every reaction can be assessed.Given that SA in this study bears an impact on every reaction, the computational consumption is considerable if the mechanism contains a significant number of elementary reactions.Thus, the skeletal mechanism preprocessed by DRGEP should be employed.The normalized sensitivity coefficients are calculated using Equation (3) [14]: where Ai is the reaction rate of the i-th elementary reaction, Cj stands for the concentration of species j or another kinetic feature, such as laminar flame speed.Equation (3) represents the contribution of an elementary reaction to the production or consumption of concerned species or a kinetic feature.
If the normalized sensitivity coefficients of a reaction are considered trivial, such a reaction has minimal effect on the entire process and can be regarded as redundant.Thus, a large number of reactions will be removed from the skeletal mechanism by setting the threshold.
A high threshold will lead to the excessive removal of reactions; as a result, the skeletal mechanism cannot be accurately described.A low threshold will produce a skeletal mechanism that contains reactions of minimal influence, thereby deviating from the actual purpose of reduction.Furthermore, after every removal of reactions from the skeletal mechanism, the retained skeletal mechanism must be carefully checked to ensure that the reaction process can still be described within certain errors.

Global Mechanism Formation by QSSA Method
In the skeletal mechanism, species with significantly short time scales also exist.The derivative of concentration to time is also approximately zero, that is: The removal of species would change the R-value in the DRGEP diagram.Thus, re-computing the R-value after the removal of a certain number of species is highly recommended.In this study, such a re-computation is performed after every removal of ~10 species.

Advanced Processing by SA
Sensitivity analysis is a method used to investigate the effect of parameter changes on the solution of mathematical models [14].In chemical kinetics, models are usually based on differential equations, which derive concentration-time or concentration-space distribution, reaction rates, and various kinetic features of the reaction.After selecting the important species and kinetic features, the importance of every reaction can be assessed.Given that SA in this study bears an impact on every reaction, the computational consumption is considerable if the mechanism contains a significant number of elementary reactions.Thus, the skeletal mechanism preprocessed by DRGEP should be employed.The normalized sensitivity coefficients are calculated using Equation (3) [14]: where A i is the reaction rate of the i-th elementary reaction, C j stands for the concentration of species j or another kinetic feature, such as laminar flame speed.Equation (3) represents the contribution of an elementary reaction to the production or consumption of concerned species or a kinetic feature.If the normalized sensitivity coefficients of a reaction are considered trivial, such a reaction has minimal effect on the entire process and can be regarded as redundant.Thus, a large number of reactions will be removed from the skeletal mechanism by setting the threshold.
A high threshold will lead to the excessive removal of reactions; as a result, the skeletal mechanism cannot be accurately described.A low threshold will produce a skeletal mechanism that contains reactions of minimal influence, thereby deviating from the actual purpose of reduction.Furthermore, after every removal of reactions from the skeletal mechanism, the retained skeletal mechanism must be carefully checked to ensure that the reaction process can still be described within certain errors.

Global Mechanism Formation by QSSA Method
In the skeletal mechanism, species with significantly short time scales also exist.The derivative of concentration to time is also approximately zero, that is: By introducing QSSA, these species can be regarded to be in the steady state and their concentration can be described by an algebraic equation set.Thus, QSSA can further reduce the equation number of the differential equation set so as to decrease the computational costs.Giral et al. [28] proposed a method to identify the kind of these species; criterion used in this method is defined by Equation ( 5): where subscript i is the i-th species; η i stands for the criterion; and ω P i and ω C i re the production and consumption rates, respectively.If the criterion for a species is lower than the given threshold, this species can be treated as a quasi-steady-state (QSS) species.
Importantly, the net reaction rate of a QSS species is equal to zero, however, if the concentration of the species is not zero, the concentration can be calculated using an algebraic equation based on the concentrations of other species and the reaction rate coefficient.On the basis of linearly independent theory, for the mechanism containing E kinds of elements and S kinds of species, if Q kinds of QSS species are included, then the final global mechanism only requires that the global (S-Q-E) steps are given, and all these steps must be linearly independent.Global mechanism is a form of mechanism expression, wherein the reaction rates of every global step are calculated by the final stage skeletal mechanism.

Detailed Mechanism for Ozone-Enhanced CH 4 /Air Combustion Process
The decomposition of ozone in the pre-heating zone of the flame can initiate and accelerate chain-branching reactions [29]; as a result, the laminar flame speed is increased and the ignition delay time is decreased.Such an enhancement has already been supported both numerically and experimentally by many researchers [30][31][32][33].Wang et al. [31] and Halter et al. [30] proposed a sub-mechanism that could be used in conjunction with GRI-Mech 3.0 to predict such an enhancement.By taking into consideration the singlet oxygen transient generated by plasma, Konnov et al. [32] integrated the sub-mechanism and extended it to describe plasma-aided combustion.In this study, detailed mechanisms are generated by combining GRI-Mech 3.0 with the sub-mechanism proposed by Wang et al. [31].The sub-mechanism is presented in Table 1.

Skeletal Mechanism
The primary skeletal mechanism was established by DRGEP.The relationship diagram was generated based on the one-dimensional laminar flame speed calculation model for the combustion of CH 4 /Air under a standard condition (pressure of 1 atm and temperature of 298.15 K).Based on the cross-combination of equivalence ratios 0.6, 1.0, and 1.4 and added ozone concentrations 0 ppm, 3000 ppm, and 7000 ppm, nine cases were examined.The reactants and products, as well as OH, CH 2 O, and NO were chosen as the reference species.The laminar flame speed was selected as the reference parameter.In this study, the selection of the threshold of Equation ( 2) started with 0.001, and the threshold was increased using the iterative calculation method to overcome the dependence of experience.The appropriate threshold was obtained by the determination of the calculation error between the reference component and reference parameter.The criterion C was defined as follows: where f represents the reference concentration or reference parameter, the subscript o represents the original case (the case calculated adopting the GRI-Mech 3.0), the subscript t represents the test case without redundant reactions, ε rtol represents the relative tolerance, and ε atol represents absolute tolerance.Given this definition, the error of the test case was within the tolerance range when the final criterion is below 1; otherwise, the magnitude of the threshold might be extremely large, and the error of the simplified skeletal mechanism might exceed the tolerance range.The tolerances of the reference components and parameter are listed in Table 2, and the variations in the error at different thresholds for the first and second iterations of calculation are shown in Figure 3.After two iterative computations, the criterion in the concentration of OH exceeded 1 (1.584), with the threshold being 0.105.The simplified mechanism with an iterative threshold of 0.085 served as the primary skeletal mechanism, which had 36 species and 219 elementary reactions.

Skeletal Mechanism
The primary skeletal mechanism was established by DRGEP.The relationship diagram was generated based on the one-dimensional laminar flame speed calculation model for the combustion of CH4/Air under a standard condition (pressure of 1 atm and temperature of 298.15 K).Based on the cross-combination of equivalence ratios 0.6, 1.0, and 1.4 and added ozone concentrations 0 ppm, 3000 ppm, and 7000 ppm, nine cases were examined.The reactants and products, as well as OH, CH2O, and NO were chosen as the reference species.The laminar flame speed was selected as the reference parameter.In this study, the selection of the threshold of Equation ( 2) started with 0.001, and the threshold was increased using the iterative calculation method to overcome the dependence of experience.The appropriate threshold was obtained by the determination of the calculation error between the reference component and reference parameter.The criterion C was defined as follows: where f represents the reference concentration or reference parameter, the subscript o represents the original case (the case calculated adopting the GRI-Mech 3.0), the subscript t represents the test case without redundant reactions, εrtol represents the relative tolerance, and εatol represents absolute tolerance.Given this definition, the error of the test case was within the tolerance range when the final criterion is below 1; otherwise, the magnitude of the threshold might be extremely large, and the error of the simplified skeletal mechanism might exceed the tolerance range.The tolerances of the reference components and parameter are listed in Table 2, and the variations in the error at different thresholds for the first and second iterations of calculation are shown in Figure 3.After two iterative computations, the criterion in the concentration of OH exceeded 1 (1.584), with the threshold being 0.105.The simplified mechanism with an iterative threshold of 0.085 served as the primary skeletal mechanism, which had 36 species and 219 elementary reactions.SA was employed for the simplification of the skeletal mechanism.We started from the analysis of the elementary reactions and assumed the elementary reactions with low sensitivity coefficients in all the cases to be redundant reactions.As a result of the removal of elementary reactions, the skeletal mechanism was downsized to its minimum scale.SA was also employed in the nine cases to which DRGEP was applied.The chosen reference components, reference parameter, and selected threshold were same as those of the DRGEP method.The iterative computation method was also employed to obtain the threshold that allows C to have a maximum value of 1.After the SA treatment, a skeletal mechanism containing 25 species and 128 reactions was generated.The final skeletal mechanism is shown in Table S1 (Supplementary Materials).
The accuracy of the skeletal mechanism was determined by the prediction of the laminar flame speeds and distribution of the critical components.A comparison of the final skeletal mechanism and detailed mechanism in predicting the laminar flame speeds at different ozone additions and equivalence ratios is shown in Figure 4.The maximum relative error (8.91%) occurred in the case with 4000 ppm ozone addition of and equivalence ratio of 1.4.The maximum absolute error (2.48 cm/s) occurred in the case with the 6000 ppm ozone addition and equivalence ratio of 1.3.The cases with significant errors were mainly those with high ozone additions and rich fuel condition.This result can be explained by the removal of C 2 , C 3 , and some alcohols in the detailed mechanism.
SA was employed for the simplification of the skeletal mechanism.We started from the analysis of the elementary reactions and assumed the elementary reactions with low sensitivity coefficients in all the cases to be redundant reactions.As a result of the removal of elementary reactions, the skeletal mechanism was downsized to its minimum scale.SA was also employed in the nine cases to which DRGEP was applied.The chosen reference components, reference parameter, and selected threshold were same as those of the DRGEP method.The iterative computation method was also employed to obtain the threshold that allows C to have a maximum value of 1.After the SA treatment, a skeletal mechanism containing 25 species and 128 reactions was generated.The final skeletal mechanism is shown in Table S1 (Supplementary Materials).
The accuracy of the skeletal mechanism was determined by the prediction of the laminar flame speeds and distribution of the critical components.A comparison of the final skeletal mechanism and detailed mechanism in predicting the laminar flame speeds at different ozone additions and equivalence ratios is shown in Figure 4.The maximum relative error (8.91%) occurred in the case with 4000 ppm ozone addition of and equivalence ratio of 1.4.The maximum absolute error (2.48 cm/s) occurred in the case with the 6000 ppm ozone addition and equivalence ratio of 1.3.The cases with significant errors were mainly those with high ozone additions and rich fuel condition.This result can be explained by the removal of C2, C3, and some alcohols in the detailed mechanism.A one-dimensional laminar premixed flame model for 4000 ppm ozone addition was derived using the detailed mechanism and skeletal mechanism separately, and the results are compared in Figure 5.For four cases with different equivalence ratios, the temperatures predicted by the detailed mechanism and skeletal mechanism show no significant difference.As a result of the removal of C2, C3, and some alcohols, the case with the equivalence ratio of 1.3 demonstrated a small deviation from the predicted critical intermediate components.The deviations of the pollutant NO in all the cases were obvious, with the maximum error being 40 ppm.On the one hand, the precursors for generating NO were removed as redundant components; on the other hand, the concentration of NO (in ppm) was lower than those of the other components.Therefore, a slight change in the mechanism can significantly affect NO prediction.Furthermore, the NO predicted by the skeletal mechanism can only reflect the overall trend of distribution.A one-dimensional laminar premixed flame model for 4000 ppm ozone addition was derived using the detailed mechanism and skeletal mechanism separately, and the results are compared in Figure 5.For four cases with different equivalence ratios, the temperatures predicted by the detailed mechanism and skeletal mechanism show no significant difference.As a result of the removal of C 2 , C 3 , and some alcohols, the case with the equivalence ratio of 1.3 demonstrated a small deviation from the predicted critical intermediate components.The deviations of the pollutant NO in all the cases were obvious, with the maximum error being 40 ppm.On the one hand, the precursors for generating NO were removed as redundant components; on the other hand, the concentration of NO (in ppm) was lower than those of the other components.Therefore, a slight change in the mechanism can significantly affect NO prediction.Furthermore, the NO predicted by the skeletal mechanism can only reflect the overall trend of distribution.

Global Mechanism
After the successful construction of the skeletal mechanism, the QSS substance was evaluated using Equation (5) for the different cases.The same iterative calculation was used to obtain the optimal threshold, and the determinant variable under threshold was treated quasi-steadily.Although some components satisfy the criterion of the QSS substance, the substances were critical intermediate active components, e.g., O and H, of the combustion process; consequently, a significant deviation might exist in the calculation of the combustion process if these components were treated as QSS substances.In this study, the QSS substances were finally confirmed to be N, C2H5, and CH3O.According to the aforementioned linearly independent theory, owing to the involvement of C, H, O, and N elements in the reaction mechanism, the constructed global mechanism contained 22 species and 18 reactions.The overall reactions are listed in Table 3.
The reaction rates of the overall reactions and the substances' reaction rates based on the reaction rates of the overall reactions and the chemical reaction rates of the QSS substances are listed in Table S2 (in Supplementary Materials).
Three cases were selected to verify the accuracy of the global mechanism.A one-dimensional laminar premixed flame models for added ozone concentrations of 2333 ppm (case 1), 3733 ppm (case 2), and 5000 ppm (case 3) were separately derived by the detailed mechanism and global mechanism at different equivalence ratios.Cases 1 and 2 adopted the design of the experiments of Wang [31] and could verify accurately the laminar flame speeds.

Global Mechanism
After the successful construction of the skeletal mechanism, the QSS substance was evaluated using Equation (5) for the different cases.The same iterative calculation was used to obtain the optimal threshold, and the determinant variable under threshold was treated quasi-steadily.Although some components satisfy the criterion of the QSS substance, the substances were critical intermediate active components, e.g., O and H, of the combustion process; consequently, a significant deviation might exist in the calculation of the combustion process if these components were treated as QSS substances.In this study, the QSS substances were finally confirmed to be N, C 2 H 5 , and CH 3 O.According to the aforementioned linearly independent theory, owing to the involvement of C, H, O, and N elements in the reaction mechanism, the constructed global mechanism contained 22 species and 18 reactions.The overall reactions are listed in Table 3.
The reaction rates of the overall reactions and the substances' reaction rates based on the reaction rates of the overall reactions and the chemical reaction rates of the QSS substances are listed in Table S2 (in Supplementary Materials).
Three cases were selected to verify the accuracy of the global mechanism.A one-dimensional laminar premixed flame models for added ozone concentrations of 2333 ppm (case 1), 3733 ppm (case 2), and 5000 ppm (case 3) were separately derived by the detailed mechanism and global mechanism at different equivalence ratios.Cases 1 and 2 adopted the design of the experiments of Wang [31] and could verify accurately the laminar flame speeds.A comparison of the laminar flame speeds from the simulation results (using both a detailed mechanism and a global mechanism) and experimental data [31] at different ozone additions (2333 ppm and 3733 ppm) is shown in Figure 6.The laminar premixed speed predictions from the experimental data and simulation result obtained by the global mechanism exhibit no significant difference.The global mechanism is developed from the detailed mechanism by reducing the redundant species and reactions.The reaction paths for radicals are therefore decreased, which can results in the under prediction of flame speed [34].The maximum error of laminar flame speed prediction was 6.29% at the equivalence ratio of 0.6 with 2333 ppm ozone addition.A comparison of the mole fractions of the species and the temperatures predicted by the detailed mechanism and global mechanism with ozone addition of 5000 ppm is shown in Figure 7.The results predicted by the global mechanism are approximate to those predicted by the detailed mechanism.Likewise, the global mechanism predictions at equivalence ratios of 0.7 and 1.3 demonstrated minimal deviations in the prediction of critical intermediate components.The global mechanism prediction at an equivalence ratio of 1.0 fit well the detailed mechanism prediction.The comparisons showed that the final global mechanism retained the computational properties of the detailed mechanism for the important components, and the prediction accuracy of the global mechanism for the combustion process was within the acceptable range.The maximum errors of temperature prediction at the three equivalence ratios of 0.7, 1.0 and 1.3 were 6.78%, 1.14% and 7.8%, respectively.
The calculation costs of the skeletal mechanism and final global mechanism at different stages are listed in Figure 8.The calculation cost of the detailed mechanism is used as the benchmark.The final global mechanism reduced the time cost by 72% with respect to the detailed mechanism.A comparison of the laminar flame speeds from the simulation results (using both a detailed mechanism and a global mechanism) and experimental data [31] at different ozone additions (2333 ppm and 3733 ppm) is shown in Figure 6.The laminar premixed speed predictions from the experimental data and simulation result obtained by the global mechanism exhibit no significant difference.The global mechanism is developed from the detailed mechanism by reducing the redundant species and reactions.The reaction paths for radicals are therefore decreased, which can results in the under prediction of flame speed [34].The maximum error of laminar flame speed prediction was 6.29% at the equivalence ratio of 0.6 with 2333 ppm ozone addition.A comparison of the mole fractions of the species and the temperatures predicted by the detailed mechanism and global mechanism with ozone addition of 5000 ppm is shown in Figure 7.The results predicted by the global mechanism are approximate to those predicted by the detailed mechanism.Likewise, the global mechanism predictions at equivalence ratios of 0.7 and 1.3 demonstrated minimal deviations in the prediction of critical intermediate components.The global mechanism prediction at an equivalence ratio of 1.0 fit well the detailed mechanism prediction.The comparisons showed that the final global mechanism retained the computational properties of the detailed mechanism for the important components, and the prediction accuracy of the global mechanism for the combustion process was within the acceptable range.The maximum errors of temperature prediction at the three equivalence ratios of 0.7, 1.0 and 1.3 were 6.78%, 1.14% and 7.8%, respectively.
The calculation costs of the skeletal mechanism and final global mechanism at different stages are listed in Figure 8.The calculation cost of the detailed mechanism is used as the benchmark.The final global mechanism reduced the time cost by 72% with respect to the detailed mechanism.

Conclusions
From the detailed ozone-enhanced CH4/air combustion mechanism including 54 species and 341 elementary reactions, a simplified skeletal mechanism containing 25 species and 128 elementary reactions was constructed by DRGEP and SA.Based on the skeletal mechanism, a global mechanism including 18 reactions and 22 species was constructed by QSSA.
The generated reduced mechanism was evaluated in terms of removing redundant components and reactions.The accuracy of the final global mechanism in predicting ozone-enhanced CH4/air

Conclusions
From the detailed ozone-enhanced CH4/air combustion mechanism including 54 species and 341 elementary reactions, a simplified skeletal mechanism containing 25 species and 128 elementary reactions was constructed by DRGEP and SA.Based on the skeletal mechanism, a global mechanism including 18 reactions and 22 species was constructed by QSSA.
The generated reduced mechanism was evaluated in terms of removing redundant components and reactions.The accuracy of the final global mechanism in predicting ozone-enhanced CH4/air

Conclusions
From the detailed ozone-enhanced CH 4 /air combustion mechanism including 54 species and 341 elementary reactions, a simplified skeletal mechanism containing 25 species and 128 elementary reactions was constructed by DRGEP and SA.Based on the skeletal mechanism, a global mechanism including 18 reactions and 22 species was constructed by QSSA.
The generated reduced mechanism was evaluated in terms of removing redundant components and reactions.The accuracy of the final global mechanism in predicting ozone-enhanced CH 4 /air combustion was verified through test cases with different equivalence ratios and ozone additions.A comparison of the prediction by the global mechanism with the prediction by the detailed mechanism and experimental results showed that the global mechanism can predict the ozone-enhanced CH 4 /air combustion process accurately.The maximum error of laminar flame speed prediction was 6.29% at the equivalence ratio of 0.6 with 2333 ppm ozone addition.The characteristics of the detailed mechanism were retained even though some species and elementary reactions were removed, and the error of prediction was within the acceptable range.The maximum errors of temperature prediction at the three equivalence ratios of 0.7, 1.0 and 1.3 were 6.78%, 1.14% and 7.8%, respectively.A comparison of the calculation cost of the global mechanism with that of the original detailed mechanism showed that approximately 72% was reduced in the calculation time cost when the final global mechanism was used.However, the prediction of the pollutant NO by the global mechanism was inaccurate.If the accuracy of NO prediction was pursued, more reactions should be added to the global mechanism, and the calculation cost would consequently increase.
Author Contributions: Z.W. and K.C. contributed to build the mechanism reduction methods; Y.L., K.W., and L.L. designed and performed the kinetics simulations together; all authors contributed to writing and revising the paper.

Figure 3 .
Figure 3. Variations in criterions and components with respect to threshold.(a) First iteration; (b) Second iteration.

Figure 3 .
Figure 3. Variations in criterions and components with respect to threshold.(a) First iteration; (b) Second iteration.

Figure 4 .
Figure 4. Predicted laminar flame speeds of cases with different ozone additions (The speeds were calculated using the detailed mechanism and skeletal mechanism).

Figure 4 .
Figure 4. Predicted laminar flame speeds of cases with different ozone additions (The speeds were calculated using the detailed mechanism and skeletal mechanism).

Figure 6 .
Figure 6.Comparison of laminar flame speeds between simulation results and experimental data from Wang [31] at different added ozone concentrations: the scatter plot represents the experimental data, the solid line represents the simulation results by the detailed mechanism, and the dashed line represents the results by the global mechanism at added ozone concentrations of (a) 2333 ppm and (b) 3733 ppm.

Figure 6 .
Figure 6.Comparison of laminar flame speeds between simulation results and experimental data from Wang [31] at different added ozone concentrations: the scatter plot represents the experimental data, the solid line represents the simulation results by the detailed mechanism, and the dashed line represents the results by the global mechanism at added ozone concentrations of (a) 2333 ppm and (b) 3733 ppm.

Figure 8 .
Figure 8.Comparison of the calculation consumption of detailed, initial skeletal, final skeletal and global mechanisms.

Figure 8 .
Figure 8.Comparison of the calculation consumption of detailed, initial skeletal, final skeletal and global mechanisms.

Figure 8 .
Figure 8.Comparison of the calculation consumption of detailed, initial skeletal, final skeletal and global mechanisms.

Table 2 .
Tolerances of the reference components and parameter.

Table 2 .
Tolerances of the reference components and parameter.