Optimal Design and Analysis of a High-Load Supersonic Compressor Based on a Surrogate Model

Abstract

To explore the internal flow mechanism and improve the performance of a supersonic compressor, an efficient global optimization design method was developed for an axial flow compressor and applied in the optimization design of a prototype supersonic compressor. Based on the multiple circular arc (MCA) blade parameters, the method can be used to parameterize the elementary stage of the blade. The optimized solution is obtained by changing the elementary stage and stacking lines of the blade during the optimization process. It has the advantages of fewer optimization variables, strong physical intuition, and a smooth surface. The optimization results show that a change in the rotor blade shape parameters has a significant effect on the compressor efficiency under design conditions, while a change in the skewed-swept parameters of the stator is the main factor that improves the compressor’s performance under near-stall conditions. Further numerical results show that the optimized rotor changes the form of the shock, weakens the degree of shock boundary layer interference, inhibits the radial migration flow of the supersonic rotor, reduces the loss of the rotor blade top, and improves the performance of the compressor under design conditions. The stator’s optimization restrains the generation of a concentrated shedding vortex at the root of the blades and greatly improves the stall margin of the compressor. Finally, the total pressure ratio and flow rate are less than 1% of the values based on the prototype operating conditions, the design mass flow of the optimized high-load supersonic compressor is increased by 0.25%, the isentropic efficiency is increased by 1.05%, and the stall margin is enhanced by 3.5%, thus verifying the effectiveness of the optimization method.

1. Introduction

As a key component of an aero-engine, the compressor accounts for up to 30% of the engine’s weight. To meet the needs of aero-engine performance improvement, the design of compressors is developing in the direction of a high-pressure ratio, high efficiency, and high stall margin [1]. However, with the increases in the compressor load and blade tip tangential speed, the significantly enhanced reverse pressure gradient and shock intensity inevitably restrict the compressor efficiency and stall margin, resulting in great challenges for the aerodynamic design of compressors.

A supersonic compressor is considered as any compressor that utilizes supersonic flow relative to one or more of the blade rows. As a rule, this method is applied only to designs that entail supersonic velocities along the entire blade span. Designs that utilize Mach numbers ranging from supersonic to subsonic along the span are generally classified as transonic. Research on supersonic compressors began in the early 20th century. Before 1950, traditional axial flow compressor designs were limited by tip-relative Mach numbers below 0.7, which resulted in severe limitations on the mass flow, speed, and stage pressure ratio of compressors. If the restriction on the relative Mach number of the blade tip is removed, the potential of the compressor for performance improvement will be significantly increased. The use of shocks in the blade row as a means of increasing the static pressure was first proposed in 1937 [2]. In this paper, it is proposed that for a supersonic flow with a relatively low Mach number, theoretically, the lower total pressure loss can result in the sufficient efficiency of the compressor; for example, the shock formed by an incoming flow of Mach 1.4 can yield a static pressure rise of 2.12, while maintaining an isentropic efficiency of 94%. In 1945, Kantrowitz and Donaldson [3] conducted a preliminary exploration of the supersonic compressor. In 1946, Kantrowitz [4] formed some basic aerodynamic assumptions on the single-stage supersonic compressor, and Wennerstorm [5] designed a supersonic compressor with swept stator blades at the hub and tip, with a total pressure ratio of 3. However, in early supersonic compressor research, the flow mechanisms of the shock and boundary layer interaction were not completely understood; thus, the separation loss caused by the strong shock decreased the efficiency of the supersonic compressor to a level much lower than that of the conventional compressor design. To improve the efficiency of the supersonic compressor, some progress was made in controlling the interference between the shocks and boundary layer using tandem blades, vortex generators, and boundary layer suction [6,7,8]. However, due to the fact that the testing of supersonic compressors relies on Freon as the working medium, experimental research based on supersonic compressors was limited by the need for environmental protection. After 1990, with the development of computational fluid dynamics technology, the numerical simulation of a three-dimensional viscous flow field could be used to accurately simulate the flow field of supersonic compressors, which effectively mitigates the impacts of shocks on supersonic compressors. Kuters and Schreiber [9] simulated the flow of supersonic cascades and proved that the results of their numerical simulation and experimental results were consistent. Severin [10] conducted a detailed numerical study on the supersonic compressor UEET, designed by NASA, and determined the flow instability mechanism of the UEET compressor and its relationship with the shock channel. Since then, numerical simulation methods have been used in the design of supersonic compressors.

With the rapid development of computer science, an increasing number of automatic optimization methods have been used in the fine design of compressors. Benini [11] developed a multi-objective optimization algorithm for the redesign of the Rotor37. Lian and Lion [12] developed a highly effective method for multi-objective optimization problems with high computational requirements. Their method integrated an experimental design, the establishment of a response surface model, and CFD calculation. Pan et al. [13] established a multi-objective optimization platform based on the Kriging surrogate model and used the improved multi-objective particle swarm optimization algorithm to optimize the supersonic tandem rotor. These methods combine numerical optimization and numerical calculation. They replace the method of artificial adjustment of the blade parameters based on experience, thus greatly reducing the design cycle. Venturelli and Benini [14] used a multi-objective optimization method to optimize a two-dimensional supersonic cascade, and the optimization results showed that the total pressure loss of the cascade was reduced by 25%. Cui [15] used a genetic algorithm to optimize the two-dimensional blade shape and three-dimensional blade shape of a supersonic rotor and studied the influences of the characteristic parameters of the blade on the rotor margin, obtaining an ultrasonic rotor design with a 2.8 pressure ratio, isentropic efficiency of 90.2%, and stall margin of 18.11%. Based on the Kriging surrogate model, Liu [16] adjusted the solidity and aspect ratio of a transonic rotor, considerably widening the stall margin of the inlet stage, using the transonic rotor of an F-class heavy-duty gas turbine compressor. Yi [17] optimized a high-pressure ratio centrifugal compressor, increasing its efficiency by 1.24%, using 41 variables to parameterize the impeller of the centrifugal compressor with SVR and random forest as surrogate models. Ma [18] achieved multi-objective global optimization by combining SVR and the NSGA-III multi-objective optimization method in the impeller optimization of a centrifugal compressor. Wang [19] proposed a similar framework for the aerodynamic optimization of turbomachinery by coupling the well-known multi-objective genetic algorithm NSGA-II and back propagation neural network. As a result, the centrifugal impeller’s efficiency increased throughout the whole speed range. These approaches are integrated with techniques for numerical computation and optimization. They significantly shorten the design period and act as a replacement for the manual adjustment of blade settings based on experience.

Although the aforementioned research has greatly advanced the aerodynamic optimization of compressors, the parameterization techniques used have a large number of variables and cannot be applied to the optimization of multi-row blades. Thus, previous research has primarily focused on rotors’ aerodynamic optimization. In this paper, a novel optimization approach is proposed to address this issue. The MCA (multiple circular arc) blades are chosen for this optimization method to parameterize the two-dimensional blade profile of the compressor. In order to maintain the continuity of the blade shape, the third-order Bessel curve is used to control the spanwise distribution of the blade profile parameter. The Bessel curve is also used to regulate the stacking line of the blades. This parameterization method greatly reduces the number of optimization variables; thus, it can be used for the optimization of the supersonic compressor stage. At the same time, the Kriging model is used as a surrogate model, as its effectiveness has been verified in previous studies [13,16]. Additionally, the best parameters for the surrogate model are identified using the genetic algorithm, which reduces the amount of time needed for optimization. The significance of the optimization variables is then examined based on the trained Kriging model, and their impacts on the performance parameters are tracked throughout the optimization process. Furthermore, by contrasting the internal flow fields of the supersonic compressors, the flow mechanism of the performance change in the supersonic compressor before and after optimization is elucidated.

2. Research Object and Numerical Simulation Method

2.1. Research Object

According to the design requirements, the preliminary design of a high-load supersonic compressor is carried out. The flow calculation is based on the Navier–Stokes method with circumferential average and is solved using the time-advanced finite volume method, which greatly reduces the dependence on the empirical model. The shock position can be captured more accurately, the separated flow of the blocked flow and meridian flow field can be predicted, and the influences of the ring wall boundary layer and radial mixing can also be taken into account [20,21]. The modified Koch–Smith model [22] is used to calculate the blade loss, as well as the blockage and transonic loss of the hub and casing; the NASA-SP-36 model, checked using a large number of experimental data points, is used as the lag angle model; the Roberts loss model [23] is used to modify the secondary flow loss distribution and secondary flow lag angle in the compressor radial direction; and the Gallimore model [24] is used to calculate the mixing loss. Given the total inlet temperature, total pressure, air flow angle, and the average static pressure at the outlet, the Spalart–Allmaras turbulence model is used to calculate these parameters [25,26]. Parameters such as the inlet and outlet angles of the blade row geometry are reasonably selected through the flow design.

Table 1. Key geometric and aerodynamic parameters of the compressor.

Parameters	Value
Loading factor	0.35
Design speed/(rpm)	25,000
Design mass flow/(kg/s)	9.8
Design total pressure ratio	2.35
Design isentropic efficiency	0.87
Blade number	65, 34, 67
Diameter of inlet casing (mm)	190
Inlet hub tip ratio	0.74
Rotor tip clearance (mm)	0.2
Rotor tip speed (m/s)	497
Rotor average aspect ratio	0.74
Rotor average solidity	1.88
Stator average aspect ratio	0.95
Stator average solidity	1.86
Rotor average diffusion factor	0.54
Stator average diffusion factor	0.48

lists the key geometry and aerodynamic parameters of the compressor.

2.2. Numerical Method

The Euranus solver from the NUMECA FINE/Turbo software package is used for all the three-dimensional numerical calculations. This solver uses Jameson’s finite volume difference scheme and the Spalart–Allmaras turbulence model, which has shown good results and has been widely used in research on supersonic compressors [13,15] to solve the three-dimensional Reynolds-averaged Navier–Stokes equation in the relative coordinate system. The spatial discretization scheme adopts the second-order precision central difference scheme, and the time advance adopts the explicit fourth-order Runge–Kutta method. The non-reflecting 1D model is selected as the interface model. To improve the computational efficiency, the multi-grid method, local time step, and residual smoothing are used to accelerate the convergence. The boundary conditions are as follows: total inlet temperature, total pressure, flow angle, and average static pressure at the outlet. The wall uses an adiabatic no-slip boundary condition, and the hub wall and the blade wall connected to the rotor blades rotate, while the casing and other parts of the hub wall are defined as stationary.

The fully structured grid is used to generate the compressor grid through the IGG/Autogrid5 module using the NUMECA FINE/Turbo package. The O4H topology is used for the blade surface and the O topology is used for the tip clearance. The thickness of the first grid layer near the wall is Y+ < 2. To verify the validity of the numerical method, a 1.5-stage compressor is checked using the numerical method. Firstly, Figure 1 shows that when the total number of compressor grid nodes reaches 2.3 million, the peak efficiency and corresponding total pressure ratio are unchanged when the number of grids is increased, which indicates that the influence of the number of grids on the calculation results can be ignored. Therefore, to ensure the accuracy of the calculation and save calculation costs, a 2.3-million-node grid was selected for the final study reported in this paper. Then, the calculation results were compared with the experimental results. As shown in Figure 2, the calculated total pressure ratio and efficiency are in agreement with the test data for the whole dimensionless flow range, meeting the requirement of engineering calculation accuracy. Therefore, the numerical calculation method selected for this paper can be used to predict the aerodynamic performance of the 1.5-stage compressor and to analyze the flow mechanism.

3. Optimum Design Method for a High-Load Supersonic Compressor

3.1. Parameterization Method

For the parameterization method used in blade shape optimization, researchers mostly use the control points of Basel surface or free-form surface transformation as control variables [27,28]. The advantages of this method are the simple parameterization and easily adaptable control parameters. However, such parameterization methods need to artificially limit the interaction between different data points to prevent the singularity of the majority of lobes identified by the algorithm from meeting the actual engineering needs. In addition, from a data point of view, the optimization results are only related to the coordinate points, which makes it impossible to quantify the effects of independent parameters on the compressor’s performance and limits the designer’s awareness of the real physical flow field. Based on these problems, important design parameters of the MCA blade are extracted as control variables. As is well-known, a 3D compressor blade can be composed of several 2D profiles stacked in the radial direction according to a certain stacking line. The parameters of the MCA blade are shown in Figure 3. Based on the above parameters, the middle arc of the 2D profile can be obtained, the leading edge and trailing edge of the blade are determined using a fixed oval without changes in the optimization process, and the length of the short axis of the ellipse is equal to the diameter of the circle of the blade shape, so as to ensure the continuity of the connection between the front and tail edges, on the one hand, and the suction surface and pressure surface, on the other. Thus, we can obtain a 2D profile. In this study, the 3D blade is composed of 11 sections following a certain stacking line, which controls the sweep and lean of the blade. Figure 4 shows this process clearly.

The inlet and outlet airflow angles are determined in the initial stage of the design and are no longer taken as optimization variables. At the same time, to ensure that the blade strength meets the engineering requirements, the maximum thickness of the 2D profile is not changed. Therefore, only the front camber angle and the maximum thickness of the 2D profile are selected as optimization variables in this paper. Based on the above parameter changes, the 3D compressor blade is composed of several 2D profiles stacked in the radial direction according to a certain stacking rule. To control the continuity of the blade surface and reduce the number of design variables, a third-order Basel curve is used to control the distribution of the 2D profile design parameters along the blade height. At the same time, two Bessel curves are used to define the stack line of the blades in order to control the blade, specifically its lean and sweep. If the control point of the stack line changes in the direction of the gas flow, the blade will sweep, and if it changes in the direction of rotation, the blade will lean.

3.2. Numerical Optimization Method

The optimization strategy discussed in this paper utilizes an adaptive surrogate model. The key concept is depicted in Figure 5 and is effective in striking a balance between the global search capability and local exploration capability of the solution, resulting in accurate and efficient optimization. The initial sample set is generated using the Latin hypercube sampling (LHS) method, and the aerodynamic performance of each sample is obtained using the high-precision numerical method. Based on the database obtained from the experimental design stage, the surrogate model is constructed. The surrogate model used in this paper is the Kriging model, whose feasibility for aerodynamic optimization was verified in [14]. The surrogate model establishes a relationship between the design parameters and the results of simulations, and it is used to evaluate the objective and constraint functions and, ultimately, determine the optimal value. In this paper, our focus is placed on the optimization of the isentropic efficiency under both design conditions and near-stall conditions. To address multiple sub-objectives, the linear weighting method is employed to convert the optimization problem into a single-objective problem, using Formula (1), and the evolutionary algorithm is used for optimization.

$$\min {\displaystyle \sum _{i=1}^r{\omega }_i\times f_i(X)}$$

(1)

In this formula, ${\omega }_i$ is the weight coefficient and satisfies the following equation:

$${\displaystyle \sum _{i=1}^r{\omega }_i=1}$$

(2)

The accurate simulation of the optimization results is used to evaluate the true objective function and constraints of the optimal value in order to verify the accuracy of the surrogate under the optimal solution. New simulation results (and possibly those from other points) are added to the database, and new design points are continuously improved to produce an increased number of accurate surrogate models throughout the design process. One repeats this design cycle until the specified maximum number of design iterations has been reached. If no better solution is obtained after the set number of iterations has been reached, then the number of iterations can be continued or the parameter selection of the optimization variable can be re-evaluated. Therefore, the design iteration process includes the following steps:

(1)

Apply the optimization algorithm to the surrogate model.

(2)

After the optimization process, verify whether it is the best solution through accurate model calculation of the optimization results.

(3)

Add a new set of samples to the database, including the optimization solutions of the previous generation, to enrich the model and build a better surrogate model. Then, verify the results of the optimization solutions using the high-fidelity numerical method in order to ensure the accuracy of the new surrogate model.

Repeated retraining of the surrogate model ensures that it will continue to represent the newly defined search area. Typically, the quality of the model increases with the number of design cycles, since new points are added to the database. If the optimization results meet the set goals, the optimization can be terminated.

3.3. Optimization Variables

The preliminary design determines the distribution of the metal angle of the blades. Thus, the airflow angle is no longer selected as the optimization variable in the subsequent optimization process, and the blade shape parameters are selected for optimization. In this paper, the front camber angle and merge chord, which corresponds to the maximum thickness position mentioned above, of the blade type are selected as the optimization variables of the 2D profile. The physical meanings of the design parameters can be observed in Figure 3. The rotor and stator blades are assembled in 11 segments to form the blades. The blade shape optimization variables in the extension direction are governed by a third-order Bessel curve, as demonstrated in Figure 6. The stack line of the blades is optimized utilizing two Bezier curves to control the blades’ sweep and lean. The ranges of all the optimization variables can be observed in

Table 2. Initial values and variables range for optimization variables.

Parameters	Value
RF_Point1(deg)	−13.96 [−20, −8]
RF_Point2(deg)	−4.7 [−15, 0]
RF_Point3(deg)	−0.05 [−7, 3]
RF_Point4(deg)	−0.04 [−5, 5]
SF_Point1(deg)	32 [25, 40]
SF_Point2(deg)	36.74 [30, 40]
SF_Point3(deg)	39.38 [33, 45]
SF_Point4(deg)	39.97 [35, 50]
RM_Point1(%)	55 [50, 60]
RM_Point2(%)	56.83 [50, 63]
RM_Point3(%)	63.5 [58, 70]
RM_Point4(%)	65 [60, 75]
SM_Point1(%)	60.05 [55, 70]
SM_Point2(%)	54.78 [50, 65]
SM_Point3(%)	50.94 [45, 60]
SM_Point4(%)	48.53 [40, 55]
Rotor_Lean_Point1(m)	0 [−0.005, 0.001]
Rotor_Lean_Point2(m)	0 [−0.002, 0.002]
Rotor_Lean_Point3(m)	0 [−0.002, 0.003]
Rotor_Lean_Point4(m)	0 [−0.001, 0.003]
Rotor_Sweep_Point1(m)	0 [−0.002, 0.002]
Rotor_Sweep_Point2(m)	0 [−0.002, 0.002]
Rotor_Sweep_Point3(m)	0 [−0.002, 0.002]
Rotor_Sweep_Point4(m)	0 [−0.002, 0.002]
Stator_Lean_Point1(m)	0 [−0.005, 0.001]
Stator_Lean_Point2(m)	0 [−0.004, 0.001]
Stator_Lean_Point3(m)	0 [−0.002, 0.002]
Stator_Sweep_Point1(m)	0 [−0.002, 0.004]
Stator_Sweep_Point2(m)	0 [−0.002, 0.004]
Stator_Sweep_Point3(m)	0 [−0.002, 0.003]

.

3.4. Optimization Objective Function and Constraints

In this study, the rotor and stator of the supersonic compressor were optimized. The objective of this optimization is to maximize the isentropic efficiency and the stall margin of the compressor under the condition in which the design pressure ratio remains unchanged. Therefore, the outlet back pressure corresponding to the compressor prototype design conditions and the near-stall conditions are taken as the outlet boundary conditions of different optimization design sample points. The objective function is to maximize the sum of the isentropic efficiency F under two operating conditions. The restraint condition determines that the relative changes in the total pressure ratio and flow rate under the two operating conditions do not exceed 1% of the values for the respective operating conditions of the prototype, as shown in Formulas (3)–(5).

$$\mathrm{F}=\max \left({\eta }_{PE}+{\eta }_{NS}\right)$$

(3)

$$\left|m_{100\%\mathrm{n}}^{opt}-m_{100\%\mathrm{n}}^{ori}\right|\le 1\% m_{100\%\mathrm{n}}^{ori}$$

(4)

$$\left|{\pi }_{100\%\mathrm{n}}^{opt}-{\pi }_{100\%\mathrm{n}}^{ori}\right|\le 1\% {\pi }_{100\%\mathrm{n}}^{ori}$$

(5)

The penalties are constructed as follows:

$$P=W\cdot {\left(\frac{Q_{ref}-Q}{Q_{ref}}\right)}^k$$

(6)

In the above equation, $Q$ is the computed quantity, and $Q_{ref}$ is a reference value used to non-dimensionalize the penalty term, which is usually selected as the calculated value of the prototype. k is the exponent of the penalty term, which is usually 2. The variable W is a weight factor that allows one to scale up or down the influence of a penalty term on the global objective function.

4. Optimization Results and Analysis

4.1. Analysis of Optimization Algorithm Model Results

The starting sample points used for training the model were 520 in number. The new points were stored in the surrogate model during the optimization process for each generation and, at the end of the optimization, two points were selected from among all the stored points. After the two points were evaluated using the high-fidelity numerical method, this iteration was finished, and the two points were added into the database, which was used to train the new surrogate model. Throughout 18 iterations, the final sample count in the database rose to 556. To evaluate the surrogate model’s accuracy in the final iteration, for this paper, we employed the leave-one-out method. Figure 7a shows the prediction results for the near-peak efficiency condition, with R² = 0.759. The results for the near-stall condition can be observed in Figure 7b, where R² = 0.847. Most of the sample points are close to the 45° regression line, with only a few deviations observed. As the iteration progressed, the optimization space contracted, and points were added to the optimal area, leading to a high fitting accuracy of the surrogate model near the optimal solution. However, the accuracy of the model in terms of the non-optimal areas was not as high, which resulted in a lower global regression coefficient R² value.

To study the influences of the optimization variables on the compressor’s performance, we used variance decomposition analysis to conduct a global sensitivity analysis based on the surrogate model. The main idea is that according to the surrogate model, a functional relationship between an n-dimensional input variable and the response can be obtained as follows:

$$y=f\left(x_1,x_2,\dots ,x_n\right)$$

(7)

The total variance of the function can be expressed as follows:

$$D=V_{X\in R^n}(f(X))={\displaystyle {\int }_{X\in R^n}{f}^2}\left(x_1,x_2,\dots ,x_n\right)d{x}_{1}d{x}_2\dots d{x}_n-f_0^2$$

(8)

The multi-dimensional integral of this last equation can be computed using Monte Carlo integration. Similarly, the following equation yields $f_i\left(x_i\right)$.

$$f_i\left(x_i\right)={\displaystyle {\iint }_1\dots }{\displaystyle {\int }_{n-1}f}\left(x_1,x_2,\dots ,x_n\right)d{x}_{1}d{x}_2\dots d{x}_{i-1}d{x}_{i+1}\dots d{x}_n-f_0$$

(9)

Again, the integration is performed using the Monte Carlo method, but this time,$x_i$ is held as constant. By performing this step multiple times, we can conduct the next integration to obtain $V\left[f_i\left(x_i\right)\right]$.

$$D_i=V_{x_i}\left(f_i\left(x_i\right)\right)={\displaystyle {\int }_{x_i}{f}_i^2}\left(x_i\right)d{x}_i$$

(10)

Then, the first-order sensitivity indices are given as follows:

$$S_i=\frac{D_i}{D}$$

(11)

The virtue of this method is that it takes all the non-linearities and interactions into account, as compared to one variable at a time by way of sensitivity analysis. Moreover, even if we conduct a sensitivity analysis only to the first order, we will always gain an idea of what remains unexplained by calculating the following equation:

$$S_{other\_order}=1-{\displaystyle \sum _{i=1}^n{S}_i}$$

(12)

As Figure 8 illustrates, the alteration in the rotor parameters under design conditions significantly affects the compressor efficiency. The front camber angle of the rotors, controlled by four control points, can increase the efficiency by over 3%. The front camber angle of the lower half of the blade has a larger impact, suggesting that the parameter selection at this point is not optimal during the design stage. The impacts of rotor sweep and lean and the merge chord change the compressor efficiency by less than 2%, and this is not shown in detail in the pie chart. In near-stall conditions, the stator parameter changes have a more significant effect on the compressor efficiency. Unlike blade shape changes, which have a larger impact on the compressor efficiency under design conditions, the stacking line of the stator is a crucial factor affecting efficiency in the near-stall conditions, accounting for about 45.7%. Parameters with an influence of less than 2% were not taken into consideration. The rotor blade shape parameters also play a role in the efficiency of the compressor in near-stall conditions, coming second only to the stator sweep parameters. Meanwhile, the stator blade shape parameters have a relatively consistent impact on the efficiency, accounting for over 2%. Additionally, the blade tip parameters have a more significant impact on the efficiency in near-stall conditions than the blade root parameters.

Based on the findings of the sensitivity analysis in Figure 8, the most influential optimization variables for each operating condition were selected. In this case, two parameters, the front camber angle of the rotors and the sweeping shape of the stators, were chosen for an analysis of their impacts on the compressor efficiency. The results are displayed in Figure 9. RF_Point1 and RF_Point3 failed to reach the optimal values due to limitations on the parameter range in this optimization. A more appropriate range should be considered for future optimizations. However, RF_Point2 and RF_Point4 revealed that the isentropic efficiency of the compressor changed from an increasing to a decreasing trend with the changing values, allowing for the determination of the optimal values of these parameters. The sweeping shape of the stator also had a limited parameter range. In further optimizations, the value range of the stator’s sweep point could be increased appropriately based on the consideration of the axial distance of the blade and the difficulty of blade processing.

4.2. Analysis of Aerodynamic Performance of Compressors

4.2.1. Overall Performance Comparison

A geometric comparison of the blade before and after optimization is shown in Figure 10 and Figure 11, and the distribution of the design parameter changes in the span direction is shown in Figure 12. The front camber angle near the root of the rotor decreases, and the load on the blade moves backward. This trend gradually decreases with the increase in the blade height. It can be observed that the front camber angle of the blade increases slightly at the 70–90% spanwise position, but due to the smaller front camber angle of the blade near the blade tip, the slight fluctuation in the value also causes a significant change in the blade shape, and at the same time, the merge chord of the rotor blade moves forward as a whole. Due to the local supersonic incoming flow at the root of the stator, the front camber angle near the root of the stator is further reduced. With the increase in the spanwise height, the Mach number of the stator incoming flow decreases, and the front camber angle of the optimized blade shape increases. At the same time, the merge chord of the stator changes slightly near the end wall of the blade after optimization, but the position of the maximum thickness in the middle part of the blade moves backward as a whole.

To depict the respective influences of changes in the 2D profiles and stacking lines on the compressor performance, Figure 13 presents the performance lines of the prototype, 2D_Opt, which only considers 2D profile optimization, and 3D_Opt, which takes into account both the 2D profile and stacking line optimizations. Optimization of the 2D profile alone leads to a 0.8% increase in peak efficiency, a 0.7% increase in the maximum pressure ratio, a 0.53% increase in the stability margin (relative to 4.33%), and a negligible 0.07% reduction in the clogging flow. With the addition of stacking line optimization, the peak efficiency rises by 1.05%, the maximum pressure ratio increases by 2.8%, the stability margin increases by 3.5%, and the maximum flow rate increases by 0.25%. The performance curve of 2D_Opt moves upwards and aligns with that of the prototype. The 3D_Opt still shows a noticeable drop in efficiency after the clogging condition is removed, but the decline is much smaller compared to that of both the prototype and 2D_Opt. Further analysis of the compressor’s flow field before and after optimization will be conducted to clarify the mechanism of the improvement in the compressor performance.

4.2.2. Impact of Optimization on Design Conditions

Figure 14 presents the static pressure distribution on the rotor blade surface according to various blade heights under the design conditions. As depicted in Figure 14, the shock at the root of the prototype rotor blade is largely removed, resulting in the absence of shock on the pressure surface. The reduction in the front camber angle of the root of the 2D_Opt rotor, however, reduces the blade load, causing the airflow to remain unaltered until it reaches 30% of the chord length. Past this point, the airflow begins to decelerate due to a series of compression waves, leading to an increase in static pressure on the blade surface. The static pressure curve of the blade surface of the 3D_Opt is the same as that of the 2D_Opt. In all three cases involving the blade root’s suction surface, the airflow undergoes deceleration due to a shock and then continues to accelerate, ultimately forming a weak shock on the suction surface. The shock structure in the channel at 50% of the blade height changes, with the pressure surface now having two weak shocks instead of the original strong shock, due to the impact of blade shape alteration. Meanwhile, the shape of the two shocks on the suction surface remains unchanged. Despite this fact, the static pressure rise distribution reveals that the 2D profile optimization strengthens the first shock on the suction surface and reduces the second shock’s intensity. Although the shape of the shock does not change dramatically compared to the 3D_Opt, the position of the first shock is shifted slightly forward. At 95% of the blade height, the shock structure in the blade channel undergoes only a minor change, with the optimization of the stacking line reducing the load on the leading edge of the blade and weakening the strength of the first shock on the suction surface.

The relative Mach number contour for each blade height of the rotor is shown in Figure 15. At 5% of the blade height, the front end of the blade is relatively straight, causing the airflow to decelerate after encountering an inlet oblique shock and then continuing to accelerate in the passage to form a second shock on the suction surface of the blade. After 2D profile optimization, the front camber angle decreases, the straight section of the front end of the blade is elongated, and the acceleration of airflow in the passage is more pronounced, resulting in the formation of a stronger shock. The shock structure at the blade roots remains unchanged with the optimization of the stacking line. The intensity of the second shock is comparable between the prototype and the 2D_Opt. At 50% of the blade height, there is a significant change in the shock structure in both the prototype and optimized blade passages. The prototype channel exhibits two oblique shocks and one normal shock, causing a noticeable separation of the suction surface boundary layer. After the 2D profile optimization, the channel shock structure remains as two oblique shocks and one normal shock. However, due to the decrease in the front camber angle and the increase in the merge chord, the second oblique shock is detached from the normal shock, and the position of the normal shock is closer to the blade trailing edge. This change causes the pressure surface to form two shock waves, and the strength of the positive shock wave in the channel is weakened, which has the benefit of decreasing the loss caused by the shock. The shock structure in the channel becomes consistent with that of the 2D profile optimization when combined with stacking line optimization. In the upper area, the optimized configuration is consistent with the original shock structure, and the shock position remains unchanged after 2D profile optimization, as the blade shape change in the upper area is not significant. The normal shock in the upper area moves closer to the trailing edge with the optimization of the stacking line.

The limited streamline on the suction surface of the rotor blade was analyzed further, and it was found that 2D profile optimization greatly inhibits the intensity of radial migration flows on the blade surface. This is due to the acceleration between shocks and the backward movement of the shock position, which reduces the influence of radial migration flow on the mainstream. As shown in Figure 16, stacking line optimization results in an exaggerated backward position of the normal shock in the upper half of the rotor, reducing the impact of the shock boundary layer interference on the blade passage. However, variation in the stacking line of the blades enhances the radial migration flow in the lower half of the rotor.

Figure 17 illustrates the pressure distribution on the stator blade surface at the design point, providing insight into the flow changes in the stator passage before and after optimization. The suction surface of the prototype stator blade root initially experiences a shock, leading to airflow separation on the blade surface and further deceleration and expansion. The pressure surface, however, lacks any shocks. On the other hand, the stator blade root flow in the 2D_Opt forms a shock before reaching 20% of the chord length, resulting in a rise in static pressure on the blade surface. Upon combination with the stacking line optimization, the suction surface of the stator blade root undergoes a more pronounced shock formed before 30% of the chord length, causing a sharp increase in the static pressure on the blade surface, and yet no boundary layer separation takes place on the blade surface. No significant changes are observed in the shock structure in the passage at 50% of the blade height. The airflow acceleration on the suction surface of the 3D_Opt is slightly greater than that of the other two configurations, but the static pressure distribution at 95% of the blade height remains unchanged.

Figure 18 compares the relative Mach numbers at various heights on the stator blade and reveals an increase in the intensity of the shock at the root of the 2D_Opt stator, which is particularly pronounced with the addition of stacking line optimization. The prototype stator’s first shock at the root separates from the blade surface, compressing the flow channel and causing the main flow area to continue to accelerate so as to form a second, weaker shock. After 2D profile optimization, the front camber angle of the stator root decreases, which extends the shock intensity but weakens the separation on the suction surface. Combined with stacking line optimization, this eliminates the separation zone at the stator root and creates a single detached shock in the blade passage. The optimization also reduces the wake area at the tip of the stator, while slightly increasing the wake area in the blade, which suggests a transfer of low-energy fluid from the tip to the blade.

To better demonstrate the impact of stacking line optimization on the flow in the stator passage, Figure 19 displays the limited streamline of the suction surface of the stator and the contour of axial velocity at the stator outlet. It can be observed that the shock causes the separation of the stator roots in the prototype conditions, leading to a vortex on the suction surface of the hub and blade. Additionally, the high load on the top of the blade results in a small recirculation area on the suction surface. The 2D profile optimization reduces the separation area at the roots of the blades. The contour of axial velocity shows a significant decrease in low-speed areas at the roots of the blades, but the recirculation area at the top of the blades remains unchanged. After the stacking line optimization, the shock-induced separation zones at the blade roots are eliminated, and the transfer of low-energy fluid from the top to the middle of the blade eliminates the recirculation zones at the top of the blade. The $Vz$ contour reveals a wider distribution of low-speed zones along the span, suggesting that the sweep and lean of the stators improve the flow in the end zone of the stator by directing low-energy fluid to the middle of the blade.

In conclusion, the improvement in efficiency observed in the optimized configuration under design conditions is due to the combination of 2D profile optimization and stacking line optimization. The 2D profile optimization changes the shock structure of the rotor, reduces the loss caused by the interaction between the shock and the rotor’s boundary layer, and weakens the intensity of the radial migration flow. On the other hand, the stacking line optimization brings the shock closer to the trailing edge in the upper half of the rotor, increasing the intensity of the radial migration flow. Additionally, the stacking line optimization eliminates the separation zone at the root of the stator blade and facilitates the migration of low-energy fluid from the blade tip to the blade, thus reducing the reflux zone at the tip.

4.2.3. Effect of Optimization on Near-Stall Conditions

To understand the cause of the stalled condition affecting the prototype supersonic compressor, the compressor was analyzed both with and without the stator. In Figure 20, the legend “Base” represents the prototype compressor, and the legend “IGV+Rotor” represents the condition where the stator is removed from the prototype compressor. As observed in Figure 20, the compressor was found to work at a lower flow coefficient without the stator, suggesting that the stall originates from the stator. Figure 21 illustrates the pressure distribution on the stator blade surface at the design point and near-stall point, providing insight into the flow changes in the stator passage. It can be observed that the shock waves in the channel are enhanced under near-stall conditions due to the change in the attack angle. The relative Mach number distribution was then analyzed at various blade heights in the baseline design under both design conditions and near-stall conditions, as depicted in Figure 22. It was observed that the flow field in the stator blade root area was severely degraded, and the low-speed region expanded to encompass the entire blade passage. The shocks became stronger and moved forward in the mid-blade region, while the Mach number at the leading edge of the blade tip increased, and the low-speed zone at the trailing edge decreased. Figure 23 compares the limited streamline distribution on the stator suction surface of the baseline design and reveals that under near-stall conditions, the small-scale, low-speed zone near the blade root caused by the shock developed into a full-blown three-dimensional stall, resulting in numerical divergence. This increased the entrainment of concentrated shedding vortices [29,30] and also led to an increase in the axial velocity of the tip fluid, explaining the decrease in the low-speed zone at the trailing edge of the blade tip.

The 3D optimization of the supersonic compressor resulted in a 3.5% increase in the stall margin. To understand the physical cause of this improvement, the prototype and the optimized configuration (3D_Opt) were compared under near-stall conditions, chosen based on the principle of equal stall margins. As observed in Figure 24, the optimization led to the elimination of large-scale separation of the stator roots, resulting in an improved flow compared to the prototype. The leading edge of the blade at 50% of the height had a longer acceleration period and weaker shock intensity. At 95% of the blade height, optimization caused the stator load to shift forward, leading to separation, due to an excessive reverse pressure gradient in the middle of the stator suction surface.

To illustrate the flow changes in the stator passage before and after optimization from multiple perspectives, Figure 25 compares the limited streamlines of the stator suction surface before and after optimization under near-stall conditions, with the $Vz$ contour at the stator outlet. It can be observed that compared with the prototype, after 3D optimization, the range of low-energy fluid accumulation on the stator suction surface is reduced, and the focus near the trailing edge of the blade root suction surface disappears, which means that the centralized shedding vortex disappears, the range of flow separation is reduced, and the degree of blockage is diminished. However, due to the increase in suction surface separation on the leading edge of the blade top and the absence of entrainment of the concentrated shedding vortex, the velocity in the passage of the blade top decreases, and the migration of low-energy fluid to the lower half of the height increases significantly under the action of stacking line optimization. It can be observed more clearly from the contour that the blade passage of the prototype stator is blocked due to the separation near the hub area, while the low-speed return area of the upper half of the stator is reduced due to the entrainment of concentrated shedding vortices. The blockage zone at the blade root of the 3D-optimized configuration is eliminated, and the distribution of the low-speed zone is more uniform in the span, thus avoiding premature stalling and improving the stall margin.

5. Conclusions

In this study, a high-load supersonic compressor was optimized using an automatic optimization platform. The optimization focused on the 2D profile and stacking line of the compressor. The key findings of this research are summarized as follows:

• The MCA blade was used to identify key design parameters, which were then utilized to create a 2D profile. The blade parameters were regulated throughout the span through the use of a third-order Bessel curve, providing a smooth and continuous shape with minimal optimization variables and the ability to incorporate design knowledge. The optimization process was carried out through a multi-objective approach, using the Kriging model as a surrogate model. As a result, the flow rate at the design point improved by 0.25%, the isentropic efficiency increased by 1.05%, and the stall margin expanded by 3.5%, leading to the significantly improved performance of the compressor.
• The results of the sensitivity analysis, carried out based on the optimized parameters derived from the trained Kriging model, indicate that the compressor stator optimally operates under the design conditions, where the variations in the rotor blade profile have a greater impact on the performance of the compressor in this stage. Under near-stall conditions, significant changes in the compressor stator sweep have a pronounced effect on the compressor’s performance in this stage.
• The optimization of the rotors’ 2D profile changes the shock structure within the supersonic compressor channel, decreases the boundary layer separation on the rotor surface, and reduces the strength of radial flow migration. Therefore, the performance of the supersonic rotor tip is improved. By optimizing the stacking line of the rotors, the flow rate of the rotor is increased, the blocking flow of the compressor is improved, and the position of the shock in the blade channel is shifted backward. However, there is a slight increase in radial flow migration on the rotor surface as a result of the stacking line optimization.
• The concentrated shedding vortex at the blade roots is weakened, but not eliminated, by the optimization of the stators’ 2D profile. With the increase in the outlet static pressure, the concentrated shedding vortex at the prototype stator roots gradually develops to produce an angular stall with a strong three-dimensional spatial flow, which results in the divergence of the calculation results of the prototype supersonic compressor. The stacking line optimization enables low-energy fluid to converge in the middle of the blade, resulting in a more uniform distribution of the low-speed zone along the span of the blade, eliminating the concentrated shedding vortex at the blade root caused by the shock and greatly improving the comprehensive margin of the compressor.