Abstract

A new design optimization method is proposed for the problem of high-precision aerodynamic design of multistage axial compressors. The method mainly contains three aspects: full-blade surface parametrization can significantly reduce the number of control variables per blade row and increase the degrees of freedom of the leading edge blade angle compared with the traditional semiblade parametric method; secondly, the artificial bee colony algorithm improved initialization and food source exploration and exploitation mechanism to enhance the global optimization ability and convergence speed, and a distributed optimization system is built on the supercomputing platform based on this method; finally, a phased optimization strategy based on the “synchronous change in multirow blades” is proposed, and expert experience is introduced to achieve a better balance between exploration and exploitation. The optimization method is applied to the AL-31F four-stage low-pressure compressor. As a result, the adiabatic efficiency is improved by 0.67% and the surge margin is improved by 3.1% under the premise that the total pressure ratio and mass flow rate satisfy the constraints, which verifies the effectiveness and engineering practicality of the proposed optimization method in the field of multistage axial flow compressor aerodynamic optimization.

1. Introduction

Aeroengines and gas turbines are needed for strategic national development, and one of the core components of aeroengines and gas turbines is the multistage axial flow compressor. For the design of multistage axial flow compressors, an aerodynamic optimization design method can effectively reduce the dependence on manual experience and greatly improve the design capability. It is well known that the aerodynamic optimization design of a multistage axial flow compressor has the characteristics of typical high-dimensionality, expensive cost, and black box (HEB) [1] problem. With an increase in optimization control variables, the design space will grow exponentially and fall into the “curse-of-dimensionality.” Therefore, obtaining the optimized solution within the acceptable time range of engineering is a frontier problem in the field of aerodynamic optimization of multistage axial flow compressors, which has a very important engineering application value.

Before the year 2000, limited by the computing power at that time, aerodynamic optimization of the multistage axial flow compressor was mainly limited to one dimension and two dimensions [2–5]. With the great improvement in computing power, three-dimensional steady aerodynamic optimization of compressors has been developed rapidly, and a few three-dimensional unsteady optimization methods have even been proposed [6]. However, the cost of unsteady calculation is too large, much simplification is required, the actual calculation accuracy is limited, and there is no engineering practicability. At present, research on multistage axial flow compressors still mainly focuses on three-dimensional steady optimization. With the development of optimization regarding multidiscipline and high-fidelity problems, the HEB problem mentioned above has increasingly become the key to restrict the design of multistage axial flow compressors.

Traditional methods for solving the HEB problem of aerodynamic optimization of multistage axial flow compressors can be divided into six categories: (1) screening variables [7, 8], (2) stage-by-stage optimization [9], (3) spatial decoupling [10], (4) surrogate model method [11–15], (5) adjoint algorithm [16–18], and (6) parametric dimensionality reduction [19, 20]. Among them, the screening method is only suitable for examples with obvious bad flow fields but is not universal. The stage-by-stage optimization process is complex, and the coupling effect of the flow field between the front and back stages is ignored. The spatial decoupling method proposed in literature [10] is only suitable for the case of a weak radial flow but not for a blade with a low aspect ratio. The surrogate model method is difficult to train in high-dimensional problems, and the number of samples needed increases exponentially with the dimensionality, so the sample training cost is difficult to handle. Although the adjoint algorithm has the advantage that optimization is independent of the dimensionality of variables, it easily falls into local extrema [21] and has defects in handling multiconstraint multiobjective problems.

The literature [22, 23] notes that one of the most promising strategies for solving the HEB problem is parametric mapping dimensionality reduction. For aerodynamic optimization of multistage axial flow compressors, the aim is to find a parametric method to reduce the design variables and achieve dimensionality reduction while retaining as much of the optimal solution from the original design space as possible. The traditional parametric method has remarkably high-dimensional characteristics because each section of the blade is modeled separately and independently of each other [24]. An effective method to constrain the control parameters of each section is the surface parametric method, which was first applied in outflow aerodynamic optimization and is an effective dimensionality reduction method [25–28]. Lee et al. compared B-spline and CST surface methods and concluded that the two parametric methods have their own advantages and disadvantages and are suitable for different situations. In 2003, the semiblade surface parametric method (blade suction surface modification by a Bezier surface) was applied in the internal flow field of compressors by Burguburu and le Pape [29], which was further developed by Wang and Zhou [30] and Cheng et al. [1]. However, the defect of this surface parametric method is that only half of the blade surface is modified, and the influence of the other half surface on the flow field is ignored. In addition, this method cannot change the blade angle at the leading edge; this defect plays an important role in the development of the flow field as it results in limiting the improvement space of the aerodynamic performance for blade geometry. In this paper, a full-blade surface parametric method is presented that regards the suction surface and the pressure surface as a whole surface and superimposes the perturbed surface on it. This method reduces the control variables and realizes an effective dimensionality reduction under the condition of ensuring the degree of freedom of the leading edge, trailing edge, and blade body.

Another important part of the optimization process is the optimization algorithm. To avoid falling into a local extremum, this paper adopts an improved artificial bee colony (IABC) algorithm with good global optimization ability [31], which does not depend on the initial solution. In addition, due to the improvements in initialization, food source exploration, and exploitation mechanisms from the standard artificial bee colony [32] (ABC) algorithm, its convergence speed is improved. Because the evolutionary algorithm has the advantage that it does not depend on the gradient information of the function, it can perform a simultaneous search of multiple individuals in the design space with the potential for multitask concurrency. Combined with a commercial supercomputing platform with powerful parallel capabilities, the IABC algorithm can greatly reduce the optimization time.

However, because evolutionary algorithms essentially require many high-fidelity performance evaluations, the above strategy is still not enough to reduce the computational cost to an acceptable range. This paper proposes a phased optimization strategy based on the “synchronous change in multirow blades,” which can greatly reduce the time cost of optimization. The core of this strategy is to draw on the idea of “exploration first, exploitation later” in the optimization process and divide the entire optimization process into two phases to effectively reduce the number of iterations required for the spatial search, thereby reducing the optimization time.

The structure of this paper is as follows: Section 2 shows the optimization method; Section 3 performs the optimization case verification and analysis; Section 4 draws a few conclusions.

2. Optimization Method Based on Full-Blade Surface Parameterization

2.1. Parametric Method

The traditional parametric method shapes each radial section of the blade separately; on the one hand, it makes the blade profile difficult to smooth; on the other hand, it makes the blade control parameters numerous (a single row of blades can have hundreds of parameters). To better promote the solution of the HEB problem of multistage axial compressor aerodynamic optimization, the surface parametric method with significant low-dimensional characteristics was selected. This method understands the blade as a combination of the suction surface and pressure surface rather than the accumulation of two-dimensional sections. This paper proposes a full-blade surface parametric method that is an improvement of the semiblade surface parametric method [29]. Figure 1 shows the construction process of the method. The blade is opened radially from the leading edge to form an unfolded surface, and the Bezier perturbed surface is superimposed to form a new unfolded surface. Then, the new unfolded surface is stitched at the leading edge in the radial direction to generate a new blade. The parametric method can be divided into the following four steps: (1)The blade opens radially at the leading edge, forming an unfolded surface(2)Chord length parameterization is applied. Because the three-dimensional surface of the physical space and the unit Bezier surface are to be superimposed, the coordinates of each point of the blade in the physical space need to be converted to the unit plane in the calculation space, and and are the axes of the computational space. The formula for chord length parameterization is as follows:where ; and refer to the number of data points of a section and the total number of sections, respectively; refers to the sum of the chord lengths of each segment of the th section in the direction; refers to the sum of the chord lengths of the th section in the direction; is the th chord length of the th section in the direction; and is the th chord length of the th section in the direction (3)Bezier perturbed surface generation is applied. The perturbed surface is defined by the following formula:where refers to the change value of each point on the perturbed surface; is the control points of the perturbed surface; and refer to the Bernstein basis functions calculated by equation (4); and refer to the independent coordinate variables of the perturbed surface, and its variation range is ; and is the combination number calculated by equation (5) (4)The value of each point of the obtained perturbation surface is superimposed on the circumferential direction of the corresponding point of the unfolded surface, thereby generating a new blade

Figure 1 

The construction process of parametric method.

The full-blade surface parametric method is a new surface parametric method that retains the high-order smoothness characteristic while improving the existing semiblade surface parametric method, which is mainly manifested in the following two points: (1)Since the full-blade surface parametric method combines the suction surface and the pressure surface as a whole, it can take into account the influence of the suction surface and pressure on the aerodynamic performance while not significantly increasing the number of control variables, and the low-dimensional characteristics are maintained well. The control parameters of a single blade row can be reduced to 16(2)The full-blade surface parametric method increases the degree of freedom of the geometric variation at the leading edge. The semiblade surface parametric method lacks the freedom of the blade angle at the leading edge, so it has certain restrictions on the geometric exploration function. Figure 2 shows a schematic diagram of the point distribution of the full-blade parametric method. The yellow point is a free active point, and the blue point is a limited active point; that is, all the blue points need to keep changing synchronously so that the leading edge of the suction pressure can be smoothly connected

Figure 2 

Schematic diagram of the point distribution of the full-blade parametric method.

2.2. Optimization Algorithm

For the purpose of searching for the global optimal solution of the multistage axial flow compressor, an evolutionary algorithm with global search ability is adopted. Since evolutionary algorithms do not require the gradient information of function values, multiple individuals can be concurrent, which can greatly reduce the optimization time. Among various evolutionary algorithms, the ABC [32] algorithm has a better optimization accuracy and convergence speed. In this paper, the ABC algorithm is developed through an initialization mechanism, a food source exploration, and an exploitation mechanism, and its schematics are shown in Figure 3.

Figure 3 

Schematic diagram of the IABC algorithm.

The IABC [31] algorithm involves employed bees, onlooker bees, and detector bees. At the beginning, the employed and onlooker bees constitute the entire bee colony, and the number of bees is divided equally between these types. The food source in the schematic refers to a set of feasible solutions, and the food source concentration refers to the fitness of the feasible solutions. The detailed steps of the IABC algorithm are as follows: (1)Food source initialization. To obtain as much design space information as possible from as few sample points as possible, instead of the random initialization mechanism in the standard ABC algorithm, the optimal Latin square experiment method [33] is utilized to initialize the distribution of the food source, and the feasible solution of -dimensional vectors is generated(2)The concentration of the initial food source is calculated, the best value is recorded, and the concentration is ranked according to the food source. The first half of the food source’s fitness is taken as the target of the employed bees, and the rest is left to the onlooker bees(3)According to formula (1), all bees explore the vicinity of the initial source and calculate the concentration of the food source. When the concentration of the new food source is greater than that of the original one, the original food source will be replaced by the better one. If the opposite is true, the bees will continue to explorewhere refers to the th component of the food source location of the th bee, are random numbers between and is the best food source location in the general information of the food source (4)Onlooker bees explore the food source. The onlooker bees use the Russian roulette rule to select food sources, as shown in formula (2); that is, the food source is selected with a probability proportional to the concentration of the food sources determined by the employed bees and is explored in the vicinity according to formula (3). If the concentration of the obtained food source is higher than the concentration of the original food source, the original food source is replaced, and if not, the exploration is continuedwhere the value of is , referring to the number of employed bees. refers to the probability of selecting a certain food source; refers to the concentration of the th and refers to the th food sources; refers to the food source position that has the largest concentration in the adjacent area; and refers to the location of the new food source of onlooker bees. The Chebyshev distance formula is used to calculate the best point in the neighborhood, which is shown in formula (4): where and is the Chebyshev distance between and .

The neighborhood of is determined by formula (5), which means that if the Chebyshev distance between and is less than or equal to the product of (the neighborhood radius) and (the average Chebyshev distance), then locates in the neighborhood of ; otherwise, it does not: where is the average Chebyshev distance between and the entire onlooker bee population, is the neighborhood of , and is the radius of the neighborhood. When is 1, the algorithm has the best convergence.

is calculated by where refers to the individual fitness of each bee.

For the purpose of testing the performance of the above algorithm, the ABC algorithm and the genetic algorithm [34] (GA) are used for comparative analysis.

The algorithm parameters of IABC and ABC are set the same, and the total number of bee colonies is 200, among which the number of employed bees is 100, the number of onlooker bees is 100, the maximum number of explorations is 100, and the maximum number of iterations is 200. The GA algorithm has 200 chromosomes and 200 iterations and a crossover probability of 0.8 and mutation probability of 0.005. Since these algorithms are probabilistic algorithms, the result of each calculation is uncertain, so the average minimum value of 30 independent cycles is used for comparison. The comparison results are shown in Table 1.

Table 1 shows that the IABC algorithm has obvious advantages for the first three benchmark functions. The minimum value obtained is approximately 10 orders of magnitude smaller than that for the other two algorithms. For the Rosenbrock function, the ABC algorithm is optimal, but the IABC algorithm is at the same level, and the disadvantage is not obvious. For the Zakharov function, the GA is optimal, but the IABC algorithm is suboptimal by less than an order of magnitude; hence, the disadvantage is acceptable.

Figure 4 shows the convergence history of these benchmark functions. It can be seen that IABC has a better and faster convergence performance in the first four benchmark function tests than the other algorithms.

(a) Griewank

(b) Sphere

(c) Rastrigin

(d) Rosenbrock

(e) Zakharov

(a) Griewank
(b) Sphere
(c) Rastrigin
(d) Rosenbrock
(e) Zakharov

Figure 4 

Comparison of the convergence history of the three algorithms based on the 5 benchmark functions.

2.3. Optimization Strategy

Engineering optimization often does not require an optimal solution, but only an optimized solution that is better than the original one. In this case, a certain optimization strategy is needed to shorten the optimization cost. The advanced parametric methods and optimization algorithms introduced in Sections 2 and 3 reduce the dimensionality of the optimization task for the multistage axial compressor, thus reducing the optimization calculation cost to a certain extent. However, considering the high cost of the high-fidelity numerical calculation of multistage axial compressors, it is necessary to further adopt an appropriate optimization strategy to complete the optimization task within the acceptable time range of the project. This article uses a combination of the following three optimization strategies.

2.3.1. Coarse Grid Optimization and Fine Grid Verification

The literature [36, 37] notes that the coarse mesh has been able to capture the near wall viscosity, so coarse and fine meshes can predict the same trend of performance variation during optimization. Although the computing power of computer is improving, it is still insufficient for the high-precision design optimization of multistage axial compressor. Therefore, it is still necessary to adopt certain optimization strategies to reduce the time consumption of single flow field calculation. “Coarse mesh optimization, fine mesh verification” is a common and effective strategy to reduce the time consumption of single flow field calculation. The flow field calculation in the optimization process uses a coarser grid, and after obtaining the optimization solution, a fine grid is used to check the calculation to obtain the optimized solution under the high-fidelity flow field calculation.

2.3.2. Multitask Concurrent Parallel Computing Based on the IABC Algorithm and Supercomputing Platform

Since evolutionary algorithm optimization does not depend on the gradient information of the function, it naturally possesses multitask concurrent computing capabilities. Using China’s advanced “Taihu Lake Supernatural Light” supercomputing platform, a multitask concurrent optimization system based on the IABC algorithm is built, which can reduce the optimization time considerably. Figure 5 shows a schematic diagram of the multitask concurrent system based on the IABC algorithm and the supercomputing platform. On the user’s computer, the initial blade is used as the input of the program to start the optimization algorithm; initialize the food source distribution; iteratively conduct employed bee, onlooker bee, and detector bee optimization; and record the location of the best food source for each bee colony generation (i.e., the optimal solution for each generation). If the loop exit condition has not been met, then the loop is continued. In any loop iteration, the initial food source, employed bees, onlooker bees, and detector bees all adopt a concurrent mechanism, so the login node needs to be started as an intermediary. Since there is an upper limit for supercomputing concurrence, the number of concurrent bee colonies must be limited to Max_Run; that is, Max_Run bees are concurrently processed. After counting one bee, the next new bee is started to join the parallel. The process of concurrent bees is shown on the user’s computer as Max_Run. New geometry files are generated through the parametric model and then placed in the corresponding Max_Run new engineering task folder. At the login node, multitask files are generated through scripts, and then, multitask files are assigned. The function of the multitask file is to start a remote computing node to perform concurrent computing. On the computing node, each task runs independently and must complete three steps, namely, grid generation, flow field calculation, and postprocessing analysis. On the login node, the script responsible for regularly detecting the status of the computing node runs at all times. If an engineering task in the computing node ends, performance results will be written in donelist.txt, and a line of donelist.txt will be added. At this time, by regularly synchronizing donelist.txt on the user’s computer and the login node, the operation status on the local user’s computer is known, and a decision on whether to dispatch a new bee can be made. Eventually, the entire optimization cycle reaches the set number of iteration steps and exits the cycle to obtain the optimized blade for this case.

Figure 5 

Schematic diagram of the multitask concurrent system based on the IABC algorithm and supercomputing platform.

2.3.3. Phased Parameterization Strategy

To find the closest possible global optimal solution of the design space within the limited number of iterations, it is necessary to further study the exploration and exploitation strategy of the design space. The IABC algorithm mentioned in the third section uses a Latin square distribution instead of a random distribution, which can improve the efficiency of space exploration in the initialization phase, but it has not yet been able to affect the efficiency of the exploration and exploitation of the design space in the optimization process. To solve this problem, drawing on the sample training in the surrogate model and the idea of “exploration first, exploitation later” of the evolutionary algorithm, this paper proposes a phased optimization strategy based on the “synchronous change in multirow blades”; that is, the aerodynamic optimization process of the multistage axial compressor is divided into two phases, as Figure 6 shows.

Figure 6 

Optimization strategy diagram.

The first phase: the “synchronous change in multirow blades” strategy is adopted; that is, for the purpose of exploring the design space at a large step, each blade row of the multistage axial compressor is modified in the same way and by the same amount. After the optimization reaches a certain number of steps, which is determined according to the designer’s experience, the decision of when to end the first phase of the optimization exploration is made, and the optimized blade of the first phase is obtained. The first phase focuses on the exploration of the design space.

After the first phase, the optimization process comes to a temporary end and the existing optimization results are analyzed to identify the most promising blade rows and corresponding blade regions (largest loss in this area). These blade rows and blade regions are selected for the second phase of automatic optimization. In this analysis and selection process, the experience of the designer is incorporated.

The second phase: it is the independent optimization of the above selected potential blade rows. The aim is to build on the results of the first phase of exploration to enable further exploitation of the corresponding design space to obtain better performing blade shapes.

The phased optimization strategy based on the “synchronous change in multirow blades” can play a practical role in the aerodynamic optimization of multistage axial compressors for the following two reasons: (1)The optimized solution only exists in a few small local areas in the huge design space. In the evolutionary algorithm, the change in each blade row in each optimization iteration is very small, and it is difficult to quickly explore the design space. If only the number of variables of each optimization in the algorithm is increased, then it is easy to increase the possibility of the mutual cancellation of the effects of different variables, but it is also difficult to quickly explore the space. The “synchronous change in multirow blades” strategy can increase the pace of exploration and greatly increase the probability of finding the optimized solution neighborhood within a limited number of iterations(2)After the first phase, the designer’s experience is utilized to analyze the flow field so that fewer optimization variables need to be used for targeted optimization exploitation in the second phase, avoiding invalid exploration of the design space and improving the optimization efficiency

3. Optimization Case Verification and Analysis

3.1. Optimization Object

The optimization example used in this paper to verify the effectiveness of the above optimization method is the AL-31F four-stage axial low-pressure compressor, with a total of ten rows of blades, as shown in Figure 7. The inlet is a guide vane, and the outlet is tandem. The reason for using a tandem blade in the last stage of the stator blade is to ensure that the outlet airflow is along the axial direction and at the same time to avoid large separation. If only one stator blade is used, the blade-camber angle of that stator blade will be too large, and thus, the airflow will be prone to large separation at the back of the blade, causing great losses and even stalling and surging. Compared with the performance of international advanced low-pressure compressors, as shown in Table 2, it can be seen that the average stage pressure ratio, adiabatic efficiency, and surge margin of the AL-31F four-stage low-pressure compressors have a certain gap compared with other advanced low-pressure compressors. Therefore, there is some room for improvement in the AL-31F’s aerodynamic optimization design.

(a) Mass flow-pressure ratio

(b) Mass flow-efficiency

(a) Mass flow-pressure ratio
(b) Mass flow-efficiency

Figure 7 

Comparison of the experimental and calculated values of the compressor performance at the design speed.

3.2. Numerical Method Verification

A CFD numerical simulation uses the Fine module of NUMECA to solve the three-dimensional Reynolds average N-S equation in finite volume form. The Autogrid5 module of NUMECA is used to generate the grid, and the grid topology of the flow passage and blade is established automatically. In this test case, the model without tip clearance and the Spalart-Allmaras one-equation turbulence model are adopted. The time discretization is solved by the fourth-order explicit Runge-Kutta method, while the false oscillation near the shock discontinuity and other small oscillations are controlled by the central difference scheme and artificial viscosity. To accelerate calculation convergence, the local time step, the implicit residual fairing, and multigrid technique are utilized.

The calculation domain and boundary condition setting for four-stage low-pressure compressor are shown in Figure 8. The total temperature at the inlet boundary is 293.15 K, the total pressure is 101,325 Pa, and the inlet air flow direction is axial. The boundary of the solid wall is adiabatic and slip-free. Given the back pressure at the outlet boundary, by gradually adjusting the outlet back pressure during the calculation process, it will advance towards the near stall point and the blocking point. The near stall point is the condition before calculating divergence that arises from increasing the back pressure of the design point, and the blocking point is the condition before calculating divergence that arises from reducing the back pressure of the design point.

Figure 8 

Calculation domain and boundary condition setting for four-stage low-pressure compressor.

Figure 7 shows the comparison between the characteristic curve of the low-pressure compressor at the design speed obtained by numerical simulation and the experimental data. At the design speed, the maximum pressure ratio measured in the experiment is 4.05, the maximum pressure ratio obtained by simulation is 4.18, and the simulated value is 3.2% larger than the experimental value. Figure 7(b) illustrates that the simulated value of the maximum efficiency is 1.5% larger than the experimental value. At the design speed, the calculated pressure ratio range and efficiency range are larger than the experimental values, and the flow margin is small. The simulated design speed margin is 10.2%, and the experimentally measured margin is 14.9%. However, the calculated characteristic curve of the low-pressure compressor is consistent with the trend of the experimental curve, which proves the validity of the numerical simulation method.

To make the verification results more credible, grid independence verification is performed on the optimized objects. Figure 9(a) shows the grid of the three-dimensional blade, Figure 9(b) shows the B2B grids of hub section, and Figure 9(c) shows the calculation performance of four sets of grids under a backpressure of 300,000 Pa. The mesh quality satisfies the near wall condition . As seen from Figure 9(c), when the number of single channel grids in each row of blades reaches approximately 1 million, the changes in efficiency and pressure ratio are very small, so it is reasonable to believe that the third set of grids (on average, 930,000 per single passage of a blade) has met the grid independence requirements. However, by considering further shortening the optimization time, the “coarse grid optimization, fine grid verification” strategy is adopted, in which the first set of grids (average 250,000 per single passage of a blade) is used in the optimization process and the third set of grids is used for the flow field calculation of the optimization results.

(a) 3D blade grids

(b) Grids of B2B

(c) Verification of grid independence

(a) 3D blade grids
(b) Grids of B2B
(c) Verification of grid independence

Figure 9 

Grid independence verification.

3.3. Optimization Framework Construction

Figure 10 shows the optimization process. First, the flow field of the initial blade is calculated. If the performance does not meet the requirements, it enters into the optimization cycle. The IABC algorithm is adopted, which includes four steps: initialization, exploration of employed bees, exploration of onlooker bees, and establishment of detector bees. Each step will produce a corresponding number of new food sources, and each food source position represents a new set of control parameters. A new blade is generated under the full-blade surface parametric control method, and then, mesh generation, flow field calculation, and performance analysis are carried out for the new blade. If the requirements are met, the optimization will exit the cycle, and the optimized blade will be obtained; otherwise, the iterative optimization of the cycle will be continued. Usually, the decision regarding when to exit the cycle should depend on the designer’s experience.

Figure 10 

Schematic diagram of the optimization process.

3.4. Optimization Settings and Result Analysis

The three optimization strategies described in Section 4 are combined to conduct the optimization task of the four-stage low-pressure compressor in two phases.

3.4.1. The First Phase: Optimization Based on the “Synchronous Change in Multirow Blades”

(1) Optimization Goals and Constraints in the First Phase. From the comparison of the above numerical simulation and experimental results, it can be seen that the efficiency and surge margin of the design speed of the AL-31F four-stage low-pressure compressor still have room for improvement. Therefore, the optimization objective of the first phase is mainly to maximize the efficiency of the design point, taking into account the surge margin improvement. According to the author’s optimization experience, when 340,000 Pa of back pressure near the stall point is selected as the optimization condition and the adiabatic efficiency under the back pressure is taken as the objective function of the optimization process, an optimization solution that can improve the design point efficiency and surge margin can be obtained under the condition that only one flow field calculation is carried out in each optimization iteration. The flow rate constraint changes within 1%, and the total pressure ratio does not decrease, and the setting of the objective function and constraint conditions in the actual optimization process is shown in

where refers to the adiabatic efficiency under optimized conditions; and refer to the total pressure ratio in the optimization process and the total pressure ratio of the original compressor, respectively; and refer to the flow rate in the optimization process and the flow rate of the original compressor, respectively; minus is the minimum value given by experts; and refers to the optimization variable, whose lower and upper limits are and , respectively. The solutions that do not satisfy the constraints are directly excluded from the IABC program, and the next generation of bees will not converge towards this food source, or the probability of convergence towards this source is almost zero. In this way, the function of constraint can be realized in an unconstrained algorithm.

The formulas for , , and SM are as follows: where and are the total pressure ratios at the outlet and inlet, respectively; is the specific heat ratio; and are the total temperature at the outlet and inlet, respectively; SM refers to the comprehensive surge margin; and and refer to the total pressure ratios at the surge point and design point, respectively.

(2) Optimization Settings in the First Phase. The full-blade surface parametric method is utilized, and all the blade rows are changed synchronously, with 16 design variables per row. The four parameters controlling the leading edge are within [-5.0, 5.0] mm, and the 12 parameters controlling the blade body geometry are within [-6.0, 6.0] mm. In the IABC optimization algorithm, the population size is 100, and the number of iterations is 8. The maximum number of mines is 3, which means that the maximum number of times the next generation of bees can explore at a particular previous generation of food source cannot exceed three times. That is, if the calculated fitness is smaller than the original one three times in a row, a new food source needs to be generated randomly again. On the supercomputing platform, approximately 850 calculations were carried out, and among them, the effective point was approximately 800 (the invalid points that did not meet the constraints were removed). Each optimization task needs 27 cores in parallel, and up to 24 tasks can be performed concurrently (when the task is not queued after supercomputing). The optimization time of each optimization task is 30 min; this process takes approximately 30 hours in total.

(3) Optimization Results and Comparative Analysis in the First Phase. Table 3 shows the performance of the four-stage low-pressure compressor at the design point with the design speed. It can be seen that under the premise of the “synchronous change in multirow blades” strategy in the first phase, the optimized blade meets the constraints of the flow and pressure ratio, and the efficiency increases by 0.48% and the surge margin 4.8%.

As shown in Figure 11, the static pressure and Mach number distribution in the hub, middle, and tip sections of the optimized solution in the first phase are presented. The losses in the R3 and R4 flow fields are not shown because they are not significant. It is illustrated that the positive incidence angles of airflow of R1 and R2 at the hub and middle sections are too large; the airflow at the leading edge is accelerated, thus increasing the shock wave intensity; and the maximum Mach number at the middle reaches 1.4, leading to a corresponding large loss at the tip section. At the same time, it can also be noted that there are two other relatively large low-speed areas: one is the tip trailing edge of S1 and the other is the outlet of S5 at the tip, thus causing corresponding loss. Therefore, these areas should be the target of the exploitation optimization in the second phase.

(a) Hub

(b) Middle

(c) Tip

(a) Hub
(b) Middle
(c) Tip

Figure 11 

Static pressure and Mach number distribution at the hub, middle, and tip sections in the first phase.

3.4.2. The Second Phase: Optimization Based on Exploitation Strategy

(1) Optimization Objectives and Constraints in the Second Phase. According to the flow field analysis of the optimized blade in the first phase, the control variables of the leading edge at the hub and the middle of R1 and R2, the leading edge at the tip of R3, the leading edge at the tip of S1, and the blade body at the middle and tip of S5 are selected as the design variables of the second phase. The rest of the blade row is fixed to the optimized blade obtained in the first phase.

The optimization objectives, conditions, and constraints of the second phase are the same as those in the first phase.

(2) Optimization Parameter Settings in the Second Phase. After the flow field calculation and analysis of the optimized geometry obtained in the first phase, a total of five rows of blade geometry of R1, R2, R3, S1, and S5 are selected as optimized in the second phase. To reduce the optimization variables for R1, R2, R3, and S1, the leading edge control parameters are optimized to reduce the positive incidence angle. To reduce the optimization variables for S1 and S5, the blade body is controlled in addition to the leading edge variables. Refer to Figure 12 for the specified control variables of each row: all the green points are leading edge control points, which are variable in the radial direction, but the four points of the same row must keep pace to avoid distortion at the leading edge; all the red points are blade-independent control points, and all the black points are fixed points.

(a) R1, R2

(b) R3

(c) S1

(d) S5

(a) R1, R2
(b) R3
(c) S1
(d) S5

Figure 12 

Distribution diagram of the control points of each selected blade in the second phase.

In the second phase, the optimization parameters are set as follows: the full-blade surface parametric method is still adopted, with a total of 25 optimization variables. The 3 control parameters of the leading edge of R1 and R2 and the 2 control parameters of the leading edge at the hub of R3 all change within the range [-6.0, 0.0] mm. The variation range of the 8 control parameters of S1 and 9 control parameters of S5 is within [-6.0, 6.0] mm. The number of iterations is 10, the size of the bee colony in the IABC optimization algorithm is 100, and the maximum number of exploitations is 3. On the supercomputing platform, a total of approximately 1,000 calculations are performed, and the number of effective points is approximately 980 (without the invalid points violating the constraints). Each optimization task requires 27 cores in parallel, and a maximum of 24 tasks can be performed concurrently (under the condition that the tasks are not queued up). The optimization time of each optimization task is 30 min, and the optimization takes approximately 48 hours in the second phase.

The time consumption of the optimization method proposed in this paper is compared with the conventional method in Table 4. From Table 4, it can be seen that the proposed optimization method can greatly reduce the number of control variables in the optimization process, thus greatly reducing the times of flow field calculations required. At the same time, the optimization time can be reduced by an order of magnitude due to the construction of a supercomputing-based concurrent system. Although the optimization time of the DOE (design of experiment) method is not too long, it only has the function of global exploration and does not have the ability of local exploitation. The proposed method in this paper has the shortest optimization consumption time among several methods.

(3) Comparative Analysis of the Optimization Results in the Second Phase. Since a single-objective optimization (adiabatic efficiency) is adopted in the second phase, it cannot find a solution that is superior to the optimized solution found in the first phase in terms of adiabatic efficiency and the surge margin in the results. According to the designer’s experience and needs, the point where the adiabatic efficiency is most improved and the surge margin is slightly decreased is selected as our final optimized solution. Table 5 shows a performance comparison of the optimization design point before and after optimization in the second phase. It can be seen from the table that after two phases of optimization, the mass flow rate and pressure ratio at the design point of the optimized blade are almost unchanged, while the adiabatic efficiency increases by 0.67% and the surge margin 3.1%. Comparing the optimization results of the first phase and the second phase, it can be seen that the adiabatic efficiency of the design point in the second phase increases by 0.19%, while the surge margin decreases by 1.7%.

Figure 13 shows a performance comparison between the optimized blade of the second phase and the original one. The change in the flow rate of the four-stage low-pressure compressor after optimization is almost 0. This is a common situation in multistage axial flow compressors. This situation indicates that although the back pressure is changing, there is a flow rate restriction in one or more stages in the multistage compressor, probably due to the fact that it has stalled there. We tend to expand the flow surge margin by adjusting the angle of multiple rows of stator blades in the normal operating condition of a multistage axial compressor. Our optimization process does not consider the adjustment of the angle of the stator blades, so it appears as a straight line on the pressurizer characteristic graph. The optimized blade can maintain a high flow rate even under higher back pressure and thus improve its surge margin.

(a) Mass flow-total pressure ratio

(b) Mass flow-efficiency

(a) Mass flow-total pressure ratio
(b) Mass flow-efficiency

Figure 13 

Comparison of the performance between the optimized blade and original blade in the second phase.

Each blade row changed before and after optimization. Figure 14 shows the geometric comparison of R1, R2, R3, R4, S1, S2, and S5 before and after optimization. Among them, R4 and S2 underwent only one optimization, and the rest of the blade deformations that underwent only one optimization were similar to them. R1, R2, R3, S1, and S5 underwent two optimizations.

(a) R1

(b) R2

(c) R3

(d) R4

(e) S1

(f) S2

(g) S5

(a) R1
(b) R2
(c) R3
(d) R4
(e) S1
(f) S2
(g) S5

Figure 14 

Geometric changes before and after optimization.

From Figures 14(d) and 14(f), it can be seen that the hub and tip shape of the optimized R4 becomes thinner and the blade-camber angle increases slightly; the leading edge part of the middle airfoil is basically unchanged, while the blade angle at the trailing edge increases slightly; thus, the blade-camber angle increases slightly; while the optimized S2 is the opposite of R4, after optimization, the hub and tip airfoils of S2 become thicker and the blade-camber angle becomes smaller, while the blade-camber angle of the middle airfoil decreases slightly.

From Figure 14(a), it can be seen that the hub blade-camber angle of R1 slightly increases, and the middle and tip airfoils do not change much with only a slight normal shift. From Figures 14(b) and 14(c), it can be seen that the hub blade-camber angle of the optimized R2 and R3 increases significantly and the leading edge blade angle decreases. The middle airfoil is similar to the hub, except that the degree is reduced.

For S1, the optimized blade-camber angle at the hub, middle, and tip regions was reduced, most significantly at the tip regions, as shown in Figure 14(e). For S5, the optimized blade-camber angles at the hub and the tip were significantly reduced, most significantly at the tip, and the blade-camber angle at the middle increased, as shown in Figure 14(g).

(4) Comparative Analysis in the Second Phase. The flow field before and after optimization should be compared and analyzed. However, it is difficult to highlight the details of the local flow field if all of them are displayed indiscriminately. Therefore, in order to better reveal the reason for the performance improvement of the optimized blade, typical blade rows are chosen for flow field analysis. The flow field changes in the four rotor blades are similar, and the flow field details of the typical blade row regions are selected to show and analyze the changes before and after optimization. Figure 15 shows a comparison of the relative Mach number distribution at the hub, midspan, and tip of R2, R3, and S5. Figure 16 illustrates the static pressure distribution of each row of blades at the hub, midspan, and tip. Combined with Figures 14–16, it can be seen that although the blade-camber angle of the hub and mid airfoil shapes of the optimized R2 and R3 becomes larger, the inverse pressure gradient near the trailing edge is reduced due to the change in load and slightly thinner thickness, resulting in a significant reduction in the airflow separation area at this location and a reduction in the airflow loss there. The situation of R3 is basically the same as that of R2. As seen from Figures 15(g) and 15(h), as the blade-camber angle of the optimized tandem S5 tip is reduced, the back pressure gradient near the trailing edge is reduced, and hence, the low-speed region caused by airflow mixing at the trailing edge of S5 is significantly reduced. As a result, the flow in the channel is smoother, thereby reducing the airflow loss in the S5 outlet area.

Figure 15 

Comparison of the relative Mach number at the hub of R2, R3 and at the tip of S5 before and after optimization.

(a) Hub

(b) Middle

(c) Tip

(a) Hub
(b) Middle
(c) Tip

Figure 16 

Comparison of the static pressure distribution before and after optimization.

The limit streamlines on the suction side of R1, R2, R3, and R4 are shown in Figure 17. From Figures 17(a)–17(d), it can be seen that for R1 and R2, the changes before and after optimization are not significant; for R3, there are large changes before and after optimization. The closed separation area at the tip of the optimized R3 is significantly reduced after optimization, thus reducing the corresponding loss. In addition, the separation line of the optimized R3 is shortened in the radial direction and thus reduce the radial separation area. At the same time, it can be observed that the optimized R3 will be attached again after separation, and then, a second separation line will be formed, reducing the separation area along the flow direction and thus reducing the corresponding separation loss. For R4, the optimized R4 has an additional reflux region at the tip, which increases the loss at the tip. For R4, original R4 has an additional reflux region at the tip, which increases the loss at the tip.

Figure 17 

Comparison of limit streamlines on the suction surfaces before and after optimization.

Figure 18 shows the entropy distribution at midspan and tip section before and after optimization, from which it can be seen that the entropy value of the airflow is significantly lower after optimization (such as regions A, B, C, D, E, and F), which indicates a smoother airflow and a corresponding reduction in losses.

(a) Midspan

(b) Tip

(a) Midspan
(b) Tip

Figure 18 

Comparison of entropy distribution at the midspan and tip section before and after optimization.

Based on the above analysis, it can be seen that after two phases of optimization, the R2 and R3 blade has a reduced reverse pressure gradient near the trailing edge at the hub and middle, resulting in a significant reduction in the separation area after the shock wave, and a separation bubble for reattachment is formed under these conditions (shown as optimized R3), thus reducing the corresponding loss and improving the efficiency. The reverse pressure gradient near the trailing edge of the S5 blade is also reduced due to the reduction of the blade-camber angle at the tip, narrowing the low-speed area near the trailing edge and downstream and reducing the corresponding airflow loss.

4. Conclusion

This paper proposes a new optimization method for high-fidelity global optimization of multistage axial compressors based on the full-blade surface parameterization. The new method is applied to the optimization of the AL-31F four-stage low-pressure compressor aerodynamics. The following conclusions are drawn: (1)Compared with the traditional parametric method, the full-blade surface parametric method can effectively reduce the dimensionality and promote the solution of the HEB problem for a multistage axial flow compressor. At the same time, compared with the semiblade surface parametric method, the full-blade surface parametric method increases the degree of freedom of the blade angle change, therefore enhancing the probability of the existence of an optimal solution(2)After verifying on the benchmark function, the IABC algorithm is shown to have better optimization performance than the GA and ABC algorithms. The initialization method and exploration mechanism are applicable to the multipeak problem of multistage axial compressor optimization(3)The phased optimization strategy based on the “synchronous change in multirow blades” achieves a good balance between global exploration and local development. By introducing expert experience, the exploration of an invalid design space by optimization algorithms is avoided, and the optimization efficiency is greatly improved(4)The new optimization method based on the full-blade parameterization can be successfully applied to the AL-31F four-stage low-pressure compressor. When the flow and pressure ratio meet the constraints, the adiabatic efficiency and surge margin are increased by 0.67% and 3.1%, respectively

Abbreviations

Data Availability

The geometry and performance data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This study was supported by the National Major Science and Technology Projects of China (CN) (2017-II-0010-0024).