One-Dimensional Analysis of Double Annular Combustor for Reducing Harmful Emissions

: The number of aircraft ﬂights worldwide has increased steadily since the introduction of air transportation to the public. Accordingly, environmental issues caused by the exhaust gases of aircraft engines have emerged. In particular, international organizations have crafted emission regulations since gases exhausted during takeoff and landing have been identiﬁed as the direct cause of air pollution near airports. Nitrogen oxides (NOx) produced at high combustion temperatures and carbon monoxide (CO) due to incomplete combustion affect the performance of the combustion chamber. Therefore, annular combustors comprising two annular zones have been developed to reduce the emissions of these two substances, which occur under different conditions. Parameters that should be considered when modifying a conventional single annular (SAC) to a double annular combustor (DAC) are discussed herein. In this paper, an optimization algorithm for obtaining the main design parameters of the DAC is presented to minimize NOx and CO emissions and an operation solution for reducing carbon monoxide emission is identiﬁed. A thermodynamic model of a high-bypass turbofan engine (PS-90A) is used to establish the inlet and outlet conditions of the combustor. Analysis results show that NOx emissions can be effectively reduced by adjusting the design parameters and CO emissions can be signiﬁcantly decreased by partially turning off the fuel supply based on the engine cycle.


Introduction
The flight times of civilian aircraft equipped with gas turbine engines have increased steadily over the past few decades owing to the promotion of air transport. Hence, the effects of aircraft emissions on the atmosphere have increased. According to Onischik et al. [1] and Korinevskaya et al. [2], the eddy airflow caused by an aircraft's wake retains harmful substances in the area around the aerodrome, thereby adversely affecting humans. In addition, Zhu et al. [3] stated that fine dust at large international airports is proportional to the number of flights. Hence, the Committee on Aviation Environmental Protection (CAEP) of the International Civil Aviation Organization (ICAO) is continuing its effort in strengthening regulations pertaining to hazardous emissions (NOx, CO, THC, etc.) [4].
An important issue in the design of civil aviation gas turbine engines is compliance with international environmental standards. Various intermediate products and combustion products are discharged when the air-fuel mixture reacts in the engine combustor. Among them, emissions, except carbon dioxide and water, adversely affect the human body. Since the 1970s, NASA Glenn Laboratories, General Electric, and Pratt and Whitney have developed various types of combustors to solve this problem [4]. Premixing in the fuel-lean state can reduce the adiabatic flame temperature for suppressing NOx generation. Using this principle, NASA developed a heterogeneously catalyzed combustor [5] that reduces the possibility of combustion instability, a swirl-can-type combustor with a rapid quenching effect, and a combustor comprising several partial premixes of fuel and air [6].

Analysis Method
Kovner et al. [23] presented the shape of a general annular-type gas turbine engine combustor, as shown schematically in Figure 1. After the pressure and temperature are increased by the compressor, the air passes through the diffuser from the left side and enters the casing. In Figure 1, the middle diameter of the last-stage stator blade of the compressor is represented by the diameter Dc.mid. Some of the air passes through the swirler to form a recirculation zone, and fuel is sprayed at the center of the swirler and burned in the combustion zone of the flame tube. The combustion zone is defined as the space where the combustion process is completed. Its axial range is regarded as the space of the flame tube from the inlet of the swirler to the beginning of the dilution zone, and its length is indicated by lcz. The remaining mass fraction of air enters the flame tube in a stepwise manner from the annular space between the casing and the flame tube through the primary zone holes, secondary holes, dilution holes, and cooling slots to participate in the combustion, dilution, and cooling processes, respectively. Subsequently, the gas exits to the right to rotate the turbine. In Figure 1, the middle diameter of the turbine inlet stator is denoted as Dg.mid.  [23,24].
A sufficient volume must be ensured to realize a successful combustion process during combustor sizing. In the initial step of the flame tube sizing, a simple annular column having sufficient volume can be assumed, as shown in the blue dotted line in Figure 1. The geometric parameter hft represents the height of the initial flame tube. Among the factors affecting the combustor volume sizing, the engine restart condition during flight can be considered. It is known that the combustion efficiency should be approximately 75% or higher for a successful restart when the gas turbine engine is stopped in the air and changed into the autorotation state [23]. According to Doroshenko [25], the combustion efficiency depends on the volume of the combustion zone and the thermodynamic state of the air at the compressor outlet. Moreover, according to Kovner et al. [23], the volume of the combustion chamber can be related to the maximum restartable altitude Hign. At first, the volume of the combustion zone Vcz can be expressed by Equation (1) [24] using Doroshenko's equation [25]. Here, the subscript H1 represents the case in which the flight Mach number is Mf ≈ 1 at altitude H = Hign, and KV is a parameter of combustion augmentation similar to parameter θ suggested by Khandelwal et al. [22]. Onischik et al. [26] and Kovner et al. [23] suggested that excess air coefficient αcz.H1 and combustion augmentation parameter KV.cz.H1 in the combustion zone at the autorotation state are close to 1.7 and in the range of 0.45-0.55, respectively. If the volume Vcz is obtained, hft can be obtained as in Equation (2) [23]. When the average axial air velocity in the flame tube without burning Cft.x is used instead of Hign, the average cross-section area of the flame tube Fft and A sufficient volume must be ensured to realize a successful combustion process during combustor sizing. In the initial step of the flame tube sizing, a simple annular column having sufficient volume can be assumed, as shown in the blue dotted line in Figure 1. The geometric parameter h ft represents the height of the initial flame tube. Among the factors affecting the combustor volume sizing, the engine restart condition during flight can be considered. It is known that the combustion efficiency should be approximately 75% or higher for a successful restart when the gas turbine engine is stopped in the air and changed into the autorotation state [23]. According to Doroshenko [25], the combustion efficiency depends on the volume of the combustion zone and the thermodynamic state of the air at the compressor outlet. Moreover, according to Kovner et al. [23], the volume of the combustion chamber can be related to the maximum restartable altitude H ign . At first, the volume of the combustion zone V cz can be expressed by Equation (1) [24] using Doroshenko's equation [25]. Here, the subscript H1 represents the case in which the flight Mach number is M f ≈ 1 at altitude H = H ign , and K V is a parameter of combustion augmentation similar to parameter θ suggested by Khandelwal et al. [22]. Onischik et al. [26] and Kovner et al. [23] suggested that excess air coefficient α cz.H1 and combustion augmentation parameter K V.cz.H1 in the combustion zone at the autorotation state are close to 1.7 and in the range of 0.45-0.55, respectively. If the volume V cz is obtained, h ft can be obtained as in Equation (2) [23]. When the average axial air velocity in the flame tube without burning C ft.x is used instead of H ign , the average cross-section area of the flame tube F ft and h ft can be obtained using the continuous equation shown in Equations (3) and (4), respectively [23]. Then the flame tube volume is calculated by the semi-empirical Equation (5) [23]. (1) (2) Khandelwal et al. [22] stated that the design process of the DAC can be similar to that of the SAC, and that the concept can be presented as a flame tube composed of multiple annular combustors, as shown in Figure 2. To form two independent combustion spaces (pilot and main zones) while using a single compressor, air and fuel must be distributed to each annular combustor. Therefore, the airflow rate ratios (G a.pz = G a.ch.pz /G a.ch , G a.mz = G a.ch.mz /G a.ch ) and fuel flow ratios (q f .pz = q f .pz /q f .ch , q f .mz = q f .mz /q f .ch ) to the pilot and main zones are required for the design of the double annular combustor. (1) 0.12 (2) . . (3) Khandelwal et al. [22] stated that the design process of the DAC can be similar to that of the SAC, and that the concept can be presented as a flame tube composed of multiple annular combustors, as shown in Figure 2. To form two independent combustion spaces (pilot and main zones) while using a single compressor, air and fuel must be distributed to each annular combustor. Therefore, the airflow rate ratios ( ̅ .
.  To compare with the conventional single annular combustor, the geometric and thermodynamic parameters of the compressor outlet and turbine inlet must be fixed [27,28]. To compare with the conventional single annular combustor, the geometric and thermodynamic parameters of the compressor outlet and turbine inlet must be fixed [27,28]. In addition, the values at the design and off-design points were derived because the condition of the air entering the combustion chamber depends on the atmospheric conditions and engine operation modes. The method suggested by Agulnik et al. [21,29] was used to calculate a one-dimensional thermo-gas-dynamic turbofan engine model at the design point. The initial input values were determined by referring to the basic data of the ground stop state (M f = 0, H f = 0) of the PS-90A [24], as listed in Table 1.

Parameter Value Units
Maximum thrust, F oo 157 kN Bypass ratio, m 4.5 -Overall pressure ratio, π oo 30 -Total temperature of turbine, T g * 1640 K stoichiometric coefficient of fuel, L 0 14.8 kg air /kg fuel ICAO normalizes the aircraft landing and takeoff (LTO) cycle by categorizing it into four stages (takeoff, climb-out, approach, and taxi/idle), as shown in Figure 3 [4,30]. The climb-out and approach stages ranged up to 3000 ft (approximately 914 m) above the air, whereas the taxi/idle stage included all processes of the aircraft moving on the ground before LTO. Table 2 shows the values of the engine thrust mode (ratio to maximum thrust P oo ) and the mode duration of the subsonic aircraft for each stage, as standardized by ICAO. Based on the proposed LTO cycle, the off-design calculation was conducted at 100%, 85%, 30%, and 7% of the full thrust. The in-house code of the Moscow Aviation Institute was used to calculate the throttle curve [27].
Energies 2021, 14 In addition, the values at the design and off-design points were derived because the condition of the air entering the combustion chamber depends on the atmospheric conditions and engine operation modes. The method suggested by Agulnik et al. [21,29] was used to calculate a one-dimensional thermo-gas-dynamic turbofan engine model at the design point. The initial input values were determined by referring to the basic data of the ground stop state (Mf = 0, Hf = 0) of the PS-90A [24], as listed in Table 1.  Figure 3 [4,30]. The climb-out and approach stages ranged up to 3000 ft (approximately 914 m) above the air, whereas the taxi/idle stage included all processes of the aircraft moving on the ground before LTO. Table 2 shows the values of the engine thrust mode (ratio to maximum thrust Poo) and the mode duration of the subsonic aircraft for each stage, as standardized by ICAO. Based on the proposed LTO cycle, the off-design calculation was conducted at 100%, 85%, 30%, and 7% of the full thrust. The in-house code of the Moscow Aviation Institute was used to calculate the throttle curve [27].  The regulatory level (ωi, g/kN), which is the estimating factor in the optimization process, is defined as the mass (g) of emitted substance i per thrust (kN) at each stage of the cycle [30]. The emissions of nitrogen oxides (NOx) and carbon monoxide (CO) were obtained using Equation (6). In addition, the permissible regulatory level of NOx in engines of this class depends on the overall pressure ratio of the compressor, as shown in Equation (7) (CAEP/8-until 2023) [30]. The CO emission is limited to less than 118 g per unit thrust . 118 g/kN ) [30].  The regulatory level (ω i , g/kN), which is the estimating factor in the optimization process, is defined as the mass (g) of emitted substance i per thrust (kN) at each stage of the cycle [30]. The emissions of nitrogen oxides (NOx) and carbon monoxide (CO) were obtained using Equation (6). In addition, the permissible regulatory level of NOx in engines of this class depends on the overall pressure ratio of the compressor, as shown in Equation (7) (CAEP/8-until 2023) [30]. The CO emission is limited to less than 118 g per unit thrust (ω re f .CO = 118 (g/kN) [30].
ω re f .NO x = −9.88 + 2.0π oo = 50.12 g kN (7) where i is the type of emission material, j is each LTO step, and EI is the emission index determined by the mass (g) of the emission material per unit fuel mass (kg). In general, the amount of NOx produced increases with the combustion temperature. The excess air coefficient α in the combustion chamber, the gas residence time in the recirculation area (combustion area), and the temperature T c * of the inlet air affect the temperature in the combustion chamber. The penetration depth of air entering the flame tube affects the cooling and mixing processes. Therefore, the cross-sectional area (F f t / ∑ µF h f ) of the flame tube to the sum of the hole areas determines the temperature uniformity at the outlet of the combustion chamber. The mixing length and flame length decrease as the gas rotation intensity (W = W t /W x ) increases by the swirler. However, according to Sarkicov et al. [18], the possibility of leanburn occurrence in the combustion zone increases with the swirl numbers of the swirler and swirl injector, thereby resulting in a significant increase in the amount of carbon monoxide. In addition, the number of injectors per cross-sectional area of the flame tube (N inj /F f t 100) can affect the distribution uniformity of the fuel. This is related to the mixing efficiency in the combustion zone. As the flow rate of air entering the swirler increases, the excess air coefficient in the recirculation area increases, and the temperature distribution becomes uniform. This is associated with the ratio of the swirler area to the sum of the total hole area (F sw / ∑ µF h f ). Sarkicov et al. [17,18] presented the emission indices of NOx and CO (as shown in Equations (8) and (9), respectively), which comprised the engine regime parameter function and combustion chamber structure parameter function.
According to Mitrofanov [31] and Jeong [32], the coefficients A 1 and B 1 depend on the type of engine and fuel. Jeong [32] presented values such as A 1 = 1.2 × 10 −4 and B 1 = 1.26 × 10 9 based on a comparison with the combustion test results of the PS-90A engine using kerosene. Figure 4 [32] shows that the emission indexes of NOx and CO calculated using these coefficients are similar to the empirical correlation derived from the test results (R 2 > 0.9). blade angle of swirler φ was in the range of 45-55 degrees. The values of the above-mentioned parameters were determined by referring to actual engine data [24]. The optimization process for the double-annular combustor sizing is presented in Figure 5b using the subroutine. In this process, the Box-Wilson gradient method was primarily used to design a combustor that can minimize emissions for a specified engine model; subsequently, values similar to the actual optimal value were searched using the simplex optimization method.  Depending on the operation purpose of the DAC, the appropriate parameters must be selected for each annular combustion zone to reduce harmful emissions. The annular zone for reducing CO emissions is the pilot zone, whereas that for reducing NOx generation is the main zone. Some terms of the empirical equations shown in Equations (8) and (9) are parabolic functions with a minimum value; however, other monotonically decreasing functions exist. Therefore, in this study, the ranges of design and performance parameters from the hot-test data [17,18] of the combustors were considered to ensure the validity of the equations [24,32]. In particular, the parameter range was determined to satisfy the condition F sw / ∑ µF h f < 0.4 to avoid the deterioration of the starting characteristics and the risk of flow separation [17]. According to Kovner et al. [23], the cross-sectional area ratio F sw / ∑ µF h f can be considered to be approximately equal to the ratio of airflow rate, and it is within the range of 0.05~0.2. Table 3 lists the values of the main parameters selected for the design of the pilot and main zones. Mathematical simulations by Mazaheri and Shakeri [33] show that changes in combustor geometry can lead to a reduction in nitric oxide emissions. Furthermore, according to Mousemi et al. [19] and Saboohi et al. [20,21], it is possible to design an optimized combustor in terms of NOx and CO reduction by tuning the coupling between the design factors of combustor components. The basic sizing of the combustor and the calculation of NOx and CO emissions are shown in the flowchart of Figure 5a in the form of a subroutine. The values of l cz and l d are known to be in the range of 1.2~2.0 and 0.5~1.5, respectively [23]. The relative amount of air used for the turbine cooling δ cool can be cal-  (10) [29], where G a.cool is the mass flow rate of air used for turbine cooling.
The average air velocity in the holes of the flame tube C o is known to have values ranging from 30-100 m/s [23]. According to Startsev [34], the swirler diameter D sw and the relative pitch between swirlers t = D sw /t in a flame chamber with a single injector row have the ranges of 0.03-0.06 m and 0.2-0.5, respectively. Pchelkin [35] stated that the blade angle of swirler ϕ was in the range of 45-55 degrees. The values of the above-mentioned parameters were determined by referring to actual engine data [24].
The optimization process for the double-annular combustor sizing is presented in Figure 5b using the subroutine. In this process, the Box-Wilson gradient method was primarily used to design a combustor that can minimize emissions for a specified engine model; subsequently, values similar to the actual optimal value were searched using the simplex optimization method.

Optimization of DAC Design Parameters
As mentioned above, the number of parameters for the DAC design is sufficiently large, even when the positions of the compressor and turbine are fixed. In this regard, a gas turbine engine combustor design process that satisfies the design objectives has been obtained using a multivariate optimization method in a previous study [37]. In that study, the gradient method applying fractional 2 N designs was first applied; subsequently, the simplex strategy was performed in the previously derived range to obtain the parameters of the combustor that minimize NOx and CO emissions [32].

Box-Wilson Strategy (Gradient Method)
The response y affected by k factors can be expressed as , , … , . Optimization is the process of obtaining values x1, x2, …, xk at the minimum or maximum evaluation factor y. An optimization example can be considered with two factors (k = 2) affecting y, as shown in Figure 6. The univariate design method requires several subsearch processes from starting point A to attain the optimum value. In Figure 6, the red square represents the endpoint of each subsearch process. However, the optimal value can be identified more quickly and accurately when multiple factors are simultaneously varied along the gradient than by varying only one factor [38,39]. Therefore, a linear regression equation model with minimal factors must be obtained in advance such that it can be used in the gradient method [38].

Optimization of DAC Design Parameters
As mentioned above, the number of parameters for the DAC design is sufficiently large, even when the positions of the compressor and turbine are fixed. In this regard, a gas turbine engine combustor design process that satisfies the design objectives has been obtained using a multivariate optimization method in a previous study [37]. In that study, the gradient method applying fractional 2 N designs was first applied; subsequently, the simplex strategy was performed in the previously derived range to obtain the parameters of the combustor that minimize NOx and CO emissions [32].

Box-Wilson Strategy (Gradient Method)
The response y affected by k factors can be expressed as y = f (x 1 , x 2 , . . . , x k ). Optimization is the process of obtaining values x 1 , x 2 , . . . , x k at the minimum or maximum evaluation factor y. An optimization example can be considered with two factors (k = 2) affecting y, as shown in Figure 6. The univariate design method requires several subsearch processes from starting point A to attain the optimum value. In Figure 6, the red square represents the endpoint of each subsearch process. However, the optimal value can be identified more quickly and accurately when multiple factors are simultaneously varied along the gradient than by varying only one factor [38,39]. Therefore, a linear regression equation model with minimal factors must be obtained in advance such that it can be used in the gradient method [38].  When the engine shaft starts to rotate, air from the compressor outlet flows unconditionally into the two annular combustors. However, the fuel can be selectively supplied to each zone. A fractional factorial experiment with five design variables ( 5), i.e., the relative mass flow rate of air in the pilot zone ( ̅ . . / . ), dimensionless relative mass flow rate of fuel in each zone ( . , . ), average axial gas velocity in the pilot zone ( . . ), and maximum restart altitude ( ) were performed to obtain a mathematical model. The dimensionless relative mass flow rate of the fuel in the pilot and main zones can be expressed, as shown in Equations (11) and (12), respectively: . .
The maximum and minimum values of each factor were based on the recommended range for the modern aircraft gas turbine engine design [23,32], whereas the other values for establishing the design plan are shown in Table 4.
Considering the significant effect of the factor and the interaction between factors, was selected as the generator. Because of this generator, the terms , , , and were included in addition to the linear terms in the regression equation.  Table 5.  When the engine shaft starts to rotate, air from the compressor outlet flows unconditionally into the two annular combustors. However, the fuel can be selectively supplied to each zone. A fractional factorial experiment with five design variables (k = 5), i.e., the relative mass flow rate of air in the pilot zone (G a.pz = G a.pz /G a.ch ), dimensionless relative mass flow rate of fuel in each zone (q f .pz , q f .mz ), average axial gas velocity in the pilot zone (C f t.x.pz ), and maximum restart altitude (H ign ) were performed to obtain a mathematical model. The dimensionless relative mass flow rate of the fuel in the pilot and main zones can be expressed, as shown in Equations (11) and (12), respectively: The maximum and minimum values of each factor were based on the recommended range for the modern aircraft gas turbine engine design [23,32], whereas the other values for establishing the design plan are shown in Table 4. Considering the significant effect of the factor and the interaction between factors, x 5 = x 2 x 3 was selected as the generator. Because of this generator, the terms x 1 x 2 , x 1 x 3 , x 1 x 4 , and x 1 x 5 were included in addition to the linear terms in the regression equation. The number of design variables included the center point, excluding the number of generators (p) (N = 2 k − p + 1 = 17) [39]. The regulatory level of the emissions against the allowable value (ω NO x = ω NO x /ω re f . NO x , ω CO = ω CO /ω re f .CO ) and the multiplication of the relative Energies 2021, 14, 3930 11 of 21 regulatory levels of NOx and CO (ω NO x ω CO ) as responses of each experiment (the row of the plan). The plan of the fractional factorial experiment and the calculation results are presented in Table 5.  The mathematical models (regression equations) for the three responses (y 1 , y 2 , y 3 ) are expressed as shown in Equations (13)-(15), respectively. It was confirmed that the coefficient of determination (R 2 ) of the modelŷ for the calculated values y were higher than 0.97.
It is extremely difficult to visually confirm the optimization process of an equation comprising more than two variables. According to Marchukov et al. [38], the mathematical regression modelŷ obtained using Equation (16) can be used to determine the gradient ∇y to determine the optimal value, where b 1 , b 2 , . . . , b 5 are the coefficients of each term in the regression model.
As shown in Figure 6, if the increment ∆x i changes in proportion to these coefficients (b i ) starting from the initial position A and toward the optimum value, then the calculation point in space moves along the gradient direction. Furthermore, the actual value ∆ x i must be directly proportional to b i I i because the variation interval I is defined as In the combustor design algorithm of Kovner et al. [23], H ign and C f t.x.pz were considered as the initial input parameters for the combustion chamber design. Because the restart altitude H ign depends on the speed C f t.x.pz in the flame tube and for the fuel flow rate, the increment in C f t.x.pz [∆ x 2 = 1, d = ∆ x 2 /(b 2 I 2 )] was first provided to obtain the increments of the remaining factors, which were relatively easy to change. Using these incremental values ( x i ) and values at the center point ( x i ), the factor value at step n 5 f in the gradient direction can be calculated ( x in = x i + n 5 f ∆ x i ). For increment values of ∆ x 1 = −0.0187, ∆ x 3 = −0.1861, ∆ x 4 = −0.0015, and ∆ x 5 = −345.72 corresponding to each factor, ω NO x ω CO exhibited a minimum value at n 5 f = −1.39717, as shown in Figure 7. The values of the factors at this point were as follows: G a. pz = 0.3261, C f t.x.pz = 33.6028 m/s, q f . mz = 1.2600, q f . pz = 1.0021, and H ign = 10483.03 m. In the combustor design algorithm of Kovner et al. [23], and . .
were considered as the initial input parameters for the combustion chamber design. Because the restart altitude depends on the speed . .
in the flame tube and for the fuel flow rate, the increment in . .
[Δ 1, Δ / ] was first provided to obtain the increments of the remaining factors, which were relatively easy to change. Using these incremental values ( ) and values at the center point ( ̄ ), the factor value at step in the gradient direction can be calculated (  The mathematical models (Equations (13)-(15)) obtained using five factors showed that imposed a slight effect on the amount of emission substances. Therefore, the gradient method was performed again using factors excluding in the next step. Considering the significant effect of the factor and the interaction between factors, was selected as the generator. Table 6 shows a plan with four factors, and Equations (17) The five-factor gradient method was used, and the other factors had increments = 0.0098, = 0.1012, and = 0.0089 for the increment of the second factor ( = 1) in the equation, ̄ . Figure 8 shows the value of for each restart altitude ( 8000, 10,000, and 12,000 m). Consequently, it was confirmed that the minimum values of were obtained owing to the previous step. The results of the primary optimization process with five and four factors are shown in rows 2-5 of Table 7. The mathematical models (Equations (13)-(15)) obtained using five factors showed that H ign imposed a slight effect on the amount of emission substances. Therefore, the gradient method was performed again using factors excluding H ign in the next step. Considering the significant effect of the factor and the interaction between factors, x 4 = x 1 x 2 x 3 was selected as the generator. Table 6 shows a plan with four factors, and Equations (17)- (19) present the expressions for ω NO x , ω CO , and ω NO x ω CO , respectively.
The five-factor gradient method was used, and the other factors had increments ∆ x 1 = −0.0098, ∆ x 3 = −0.1012, and ∆ x 4 = −0.0089 for the increment of the second factor (∆ x 2 = 1) in the equation, x in = x i + n 4 f ∆ x i . Figure 8 shows the value of ω NO x ω CO for each restart altitude (H ign = 8000, 10,000 and 12,000 m). Consequently, it was confirmed that the minimum values of ω NO x ω CO were obtained owing to the previous step. The results of the primary optimization process with five and four factors are shown in rows 2-5 of Table 7. Table 6. Fractional factorial design and experimental results involving four factors.

Run (u) Factor
Response

Nelder-Mead Simplex Method
In the previous section, significant values of the main parameters were obtained using the gradient method. As the next step of the optimization process, the simplex method was used to obtain a value similar to the actual optimal value. As explained in Section 3.1, a function comprising two factors can be used to search for an optimal value by a group of trial points with the shape of continuous equilateral triangles (see Figure 9). However, it is impossible to visually confirm the optimal value of functions with three or more factors. The basic principle of the simplex method is to obtain the opposite vector ( → r N ) of the farthest vector ( → r worst ) from the optimal value in a k-dimensional space (k = number of factors).
Energies 2021, 14, x. https://doi.org/10.3390/xxxxx www.mdpi.com/journal/energies In the previous section, significant values of the main parameters were obtained using the gradient method. As the next step of the optimization process, the simplex method was used to obtain a value similar to the actual optimal value. As explained in Section 3.1, a function comprising two factors can be used to search for an optimal value by a group of trial points with the shape of continuous equilateral triangles (see Figure 9). However, it is impossible to visually confirm the optimal value of functions with three or more factors. The basic principle of the simplex method is to obtain the opposite vector ( ) of the farthest vector ( ) from the optimal value in a k-dimensional space (k = number of factors). The first simplex was formed with the average value of the results of the previous step, and the calculation was performed using the algorithm suggested by Marchukov et al. [38]. For a more detailed search, each factor was assigned a smaller variability interval (I) than that in the previous step ( 0.01; 1 m/s; I 0.01; 0.01; 100 m). The first simplex (Aij) is a type of lower triangular matrix, and the first simplex can be obtained using Equation (20) [38], where , is a vector of values obtained in the previous step (i is the number of factors, and j represents the number of experiments). , The matrices and vectors used in the calculation were as follows: The new position , constituting the (s + 1)-th simplex was obtained using Equation (21) with the farthest vector , from the optimal value in the s-th simplex comprising vectors , .
The results obtained through this process are listed in Table 7. The calculated values for the SAC before it is converted to the DAC are presented in the first row of Table 7 [32]. The first simplex was formed with the average value of the results of the previous step, and the calculation was performed using the algorithm suggested by Marchukov et al. [38]. For a more detailed search, each factor was assigned a smaller variability interval (I) than that in the previous step (I 1 = 0.01; I 2 = 1 m/s; I 3 = 0.01; I 4 = 0.01; I 5 = 100 m). The first simplex (A ij ) is a type of lower triangular matrix, and the first simplex x ij can be obtained using Equation (20) [38], where m 0,i is a vector of values obtained in the previous step (i is the number of factors, and j represents the number of experiments).
The matrices and vectors used in the calculation were as follows: The results obtained through this process are listed in Table 7. The calculated values for the SAC before it is converted to the DAC are presented in the first row of Table 7 [32].

Initiatives to Reduce CO Emission
The results obtained using the three-stage optimization algorithm show that the emission of nitrogen oxides decreased significantly, but the multiplication of the evaluation index (ω NO x ω CO ) exceeded the reference value. Based on Table 7, this result was caused by CO emission, which exceeded twice the reference value. As shown in Figure 10, the carbon Energies 2021, 14, 3930 15 of 21 monoxide emission per unit mass of fuel was significant in the taxi/idle and approach stages, indicating that the combustion efficiency was extremely low in these two stages. In particular, CO was generated in the main zone.
decreased from the value of 33.33 m/s. In the case of S, the values obtained using the simplex method remained constant for parameters other than . .
. The excess air coefficient of the main zone was selected as 10 considering the air-rich combustion instability range [29]. When the average gas velocity in the main zone decreased starting from the velocity value obtained by the simplex method ( . . 33.33 m/s), the CO emission decreased, whereas the NOx emissions increased. When combustion was terminated in one or more LTO regimes, the CO emissions remained below the reference value. However, to ensure the safe flight of an aircraft, supplying fuel to only one zone at the approach stage may not be a good option.   According to the basic principle of combustion theory, combustion efficiency can be increased by increasing the gas velocity in the main zone. Additionally, it has been reported that CO generation was suppressed in the air-rich state [29,40]. According to Gleason and Niedzwiecki [41], the fuel supply in the main zone for the two modes (approach and taxi/idle) can be turned off; however, this will yield poor results. Five cases (S, A, B, C, and D) were considered, as shown in Table 8, which reflect this method. Figure 11 shows that the relative regulatory levels varied as C f t.x.mz decreased from the value of 33.33 m/s. In the case of S, the values obtained using the simplex method remained constant for parameters other than C f t.x.mz . The excess air coefficient of the main zone α mz was selected as 10 considering the air-rich combustion instability range [29]. When the average gas velocity in the main zone decreased starting from the velocity value obtained by the simplex method (C f t.x.mz = 33.33 m/s), the CO emission decreased, whereas the NOx emissions increased. When combustion was terminated in one or more LTO regimes, the CO emissions remained below the reference value. However, to ensure the safe flight of an aircraft, supplying fuel to only one zone at the approach stage may not be a good option. An important aspect of gas turbine combustor design is the uniformity of the temperature distribution at the turbine inlet [23]. Hence, the relative fuel flow rate (q f ) of the two zones in the takeoff stage was assigned the same value in case B . In addition, case B" can be considered as a case involving a lower amount of air entering the pilot zone. According to Kovner et al. [23], less CO is emitted as the volume of the flame tube increases. When H ign was increased, the pilot zone volume increased while the airflow ratio (G a. pz ∼ = 0.33) remained fixed. According to the regression equations shown in Equations (9) and (10) in Section 3.1, as well as Figure 8, the effect of H ign on the emissions was slight. Therefore, H ign can be increased while maintaining the remaining values (case B-20 km in Table 8). Before reducing the velocity (at C f t.x.mz = 33.33 m/s), the evaluation index (ω NO x ω CO ) of B and D indicated minimum values, as shown in Figure 12. In cases B and B", the CO emissions reduced by almost one-half, but the NOx emissions exceeded the reference value. Both substances did not exceed the reference value in the case of B-20 km, but the amount of NOx increased slightly compared with that of case B. An important aspect of gas turbine combustor design is the uniformity of the temperature distribution at the turbine inlet [23]. Hence, the relative fuel flow rate ( ) of the two zones in the takeoff stage was assigned the same value in case B′. In addition, case B″ can be considered as a case involving a lower amount of air entering the pilot zone. According to Kovner et al. [23], less CO is emitted as the volume of the flame tube increases. When was increased, the pilot zone volume increased while the airflow ratio ( ̅ . ≅ 0.33) remained fixed. According to the regression equations shown in Equations (9) and (10) in Section 3.1, as well as Figure 8, the effect of on the emissions was slight. Therefore, can be increased while maintaining the remaining values (case B-20 km in Table 8). Before reducing the velocity (at . .

/ ), the evaluation index (
) of B and D indicated minimum values, as shown in Figure 12. In cases B′ and B″, the CO emissions reduced by almost one-half, but the NOx emissions exceeded the reference value. Both substances did not exceed the reference value in the case of B-20km, but the amount of NOx increased slightly compared with that of case B.
In the initial model (initial SAC in Figure 12) before conversion to the DAC, the amount of NOx emitted exceeded the CAEP/8 reference value. NOx emissions decreased considerably in the DAC, but CO emission increased in case B after applying the optimization process. As shown in Figure 12 and Table 9, a similar trend was observed in the actual SAC and DAC data of the CFM56 engine [42,43]. In the case of B-20km, the amount of NOx increased as the CO emission reduced, which is comparable to the data of GE-90. The relative NOx emission of the DAC, derived through several calculation steps, was similar to the actual engine data; however, it was confirmed that the CO emissions differed slightly.   In the initial model (initial SAC in Figure 12) before conversion to the DAC, the amount of NOx emitted exceeded the CAEP/8 reference value. NOx emissions decreased considerably in the DAC, but CO emission increased in case B after applying the optimiza-tion process. As shown in Figure 12 and Table 9, a similar trend was observed in the actual SAC and DAC data of the CFM56 engine [42,43]. In the case of B-20 km, the amount of NOx increased as the CO emission reduced, which is comparable to the data of GE-90. The relative NOx emission of the DAC, derived through several calculation steps, was similar to the actual engine data; however, it was confirmed that the CO emissions differed slightly.

Conclusions
Based on the thermodynamic theory of the turbofan engine and turbomachinery theory, a one-dimensional model of a SAC of the PS-90A engine was developed. The relationship between NOx and CO emissions and the design parameters of the combustor were derived by comparing the values obtained from the calculation of off-design points and the experimental data pertaining to four LTO modes of the ICAO. Empirical equations were included in the design algorithm to design the DAC. Because NOx and CO were generated under different conditions, the main and pilot zones constituting the combustor reduced the amount of each substance. Hence, optimization was applied to select the initial design values of the two zones.
In this study, the relative mass flow rate of air in the pilot zone (G a.pz ), dimensionless relative mass flow rate of fuel in each zone (q f .pz , q f .mz ), average axial gas velocity in the pilot zone (C f t.x.pz ), and maximum restart altitude (H ign ) were used as initial values. To perform the multivariate optimization algorithm using the gradient method and the simplex method in a stepwise manner, the experimental plan theory was applied in advance to obtain a regression model of the calculated values.
To comprehensively evaluate the amount of harmful emissions, the product of the relative regulatory levels of NOx and CO (ω NO x ω CO ) was introduced as an evaluation index. As a result of performing optimization using this evaluation index, the amount of NOx reduced significantly, whereas the combustion efficiency of the main zone in the taxi/idle regime decreased significantly. It was confirmed that terminating the mainzone fuel flow in the low-thrust regime can more effectively solve the problem than by sustaining combustion in this zone. Decreasing the gas velocity in the pilot zone can reduce CO emission; however, NOx emissions may increase. The values derived in this study did not differ significantly from the actual engine test data.
Hence, the possibility of designing a combustor that satisfies the environmental standards of CAEP/8 was indicated through the proposed DAC design method. It is thought that the optimal design parameters of the independent annular pilot and main zones can be obtained for reducing the emissions through one-dimensional design. Since the turbine inlet is shared, however, the pilot and main zones must be connected in the dilution zone. The CO emission derived from this study may be higher than the actual value as the combustion gas from the two independent combustion zones and the additionally introduced air are mixed in such a dilution zone. Compared with the GE90 engine test data, the case B-20 km showing the most reasonable results has a difference of more than 40% in regulatory levels of CO emission, while the difference is less than 5% for the NOx emission.