**Multi-Condition Optimization of Cavitation Performance on a Double-Suction Centrifugal Pump Based on ANN and NSGA-II**

**Wenjie Wang 1 , Yanpin Li 2 , Majeed Koranteng Osman 1,3,\* , Shouqi Yuan 1 , Benying Zhang 1 and Jun Liu 4**

Received: 25 July 2020; Accepted: 4 September 2020; Published: 10 September 2020

**Abstract:** Double-suction centrifugal pumps form an integral part of power plant systems in maintaining operational stability. However, there has been a common problem of achieving a better cavitation performance over a wider operating range because the traditional approach for impeller design often leads to the design effect not meeting the operational needs at off-design conditions. In addressing the problem, an optimization scheme was designed with the hub and shroud inlet angles of the double-suction impeller to minimize the suction performance at non-design flow conditions. A practical approach that speeds up the cavitation simulation process was applied to solve the experimental design, and a multi-layer feed forward artificial neural network (ANN) was combined with the non-dominated sorting genetic algorithm II to solve the multi-objective problem into three-dimensional (3D) Pareto optimal solutions that meet the optimization objective. At the design point, the suction performance was improved by 6.9%. At non-design flow conditions, the cavitation performance was improved by 3.5% at 1.2*Qd* overload condition, 4% at 0.8*Qd*, and 5% at 0.6*Qd*. Additionally, there was significant reduction in the attached cavity distribution in the impeller and suction domains when the optimized model was compared to the original model at off-design points. Finally, the optimization established a faster method for a three-objective optimization of cavitation performance using ANN and 3D Pareto solutions.

**Keywords:** multi-condition optimization; cavitation performance; artificial neural networks (ANN); net positive suction head (NPSH); double suction

#### **1. Introduction**

In power plants, cooling water circulation pumps play an important role of distributing water to various parts of the plant to maintain operational stability. These pumps, mostly centrifugal pumps with double-suction impellers are usually operated continuously for long running hours for systems with large capacity demands due to their ability to deliver large flow rates at a constant head [1,2]. There is however a common problem of achieving a better cavitation performance over a wide operating range because traditional approach for impeller design often leads to the design effect not meeting the requirements at off-design conditions. This puts the downstream components at risk of being damaged since in most cases the flow passage gets blocked by the attached cavities [3,4]. In view of

this, the cooling water pumps are usually operated at rated design conditions to minimize the risk of cavitation which can lead to costly damage.

There has been an extensive literature on cavitation in centrifugal pump impellers over the past years, and these are readily available in open literature. Quite a number of these investigations predicted that there is some correlation between cavitation instabilities and changes in flow angle, the pressure gradient at cavity closure, as well as the flow rate and cavitation number [5,6]. More so, the unsteadiness in the pressure field vis-à-vis the rapid head-drop phenomenon that occurs during cavitation has been periodically analyzed in a centrifugal pump, and the reason for the head drop was attributed to the vortex generation in the cavitating region [7,8].

Pei [9] applied the L9 (33 ) orthogonal design of experiments (DOE) and computational methods to improve the required net positive suction head (NPSHr) by 0.63 m. The relationship between the hub and shroud inlet angles have also been studied and established that a shroud blade angle of 30◦ improved cavitation performance [10]. The inverse design methodology has also been applied to parameterize the blade profile to improve suction performance and efficiency on a mixed flow pump impeller and diffuser [11]. Xu [12] performed a multi-objective optimization study using the Taguchi method. In his study, five design parameters of a centrifugal pump were sampled into 16 different impeller designs and optimized to improve design efficiency and cavitation performance.

Recently, there has been a shift from the single objective optimization methods to the application of genetic algorithms (GA) and surrogate models to find optimum global geometrical parameter combinations that can solve multi-objective optimization problems [13,14]. Jin [15] extensively reviewed the various surrogate models using their prediction accuracy, efficiency, and robustness, and concluded that for higher-order nonlinear problems, neural networks should be used. The application of artificial neural networks has gained much recognition in turbomachinery optimization. In centrifugal pumps, the radial bias neural network (RBNN) had the best prediction accuracy when compared with kriging and the response surface approximation model [16]. Moreover, Pei [13] carried out a multi-objective optimization on the inlet pipe shape of a vertical inline pump using artificial neural network (ANN) and multi-objective genetic algorithm (MOGA) to increase the efficiency over a wider range. Furthermore, the multi-layer neural network has been established to have a better prediction accuracy compared to the single-layer neural network [17].

Despite several optimization works on cavitation performance in centrifugal pumps, most of the studies have been focused on optimizing the inlet blade angle without focusing much on other impeller parameters such as the hub and shroud angles whose effect on cavitation performance has not been well researched. In addition, there is a common problem of achieving a better cavitation performance over a wider operating range in centrifugal pumps [18] since cavitation optimization studies have usually been centered on design and near-design flow conditions.

Therefore, in order to adapt to its multi-operating condition and broaden its non-cavitation operating range, an orthogonal scheme was designed to execute a multi-flow condition optimization for cavitation performance using a multi-layer feedforward neural network and the non-dominated sorting genetic algorithm II (NSGA-II). Since the traditional approach to predicting NPSHr is time consuming, a faster approach [19,20] was applied to reduce the simulation time during the optimization process. For this paper, a numerical anti-cavitation optimization design of the double-suction impeller was carried out at part load, the design point and at an overload condition to establish a set of optimum parameters that would best improve the suction performance at off-design flow conditions.

#### **2. Optimization Procedure**

Figure 1 shows the optimization procedure. The initial process was to sample the input bound variables based on an orthogonal design of experiment. Secondly, series of 3D model double-suction impellers were designed by bladegen in ANSYS 2019 R3 and CREO 5.0. ANSYS was adopted for the numerical simulations to obtain the objective function values, which was used to train the surrogate models. The multi-objective problem was then solved to obtain 3D Pareto solutions, which were then verified by numerical simulation to determine the reliability of the optimization.

**Figure 1.** Flowchart of optimization process.

#### *2.1. Objective Functions*

During the optimization process of centrifugal pumps, the goals and targets are considered as significant indexes for performance evaluation. For this study, the NPSHr at three flow conditions were selected as the optimization targets. The objective function is to minimize the NPSHr at part loads (0.8*Qd*), the nominal flow condition (1.0*Qd*), and overload (1.2*Qd*) flow conditions. These conditions are obtained by numerical calculations. The mathematical relation for NPSH is

$$\text{NPSH} = \frac{p - p\_v}{\rho g} + \frac{V\_{\text{in}}}{2g} \tag{1}$$

*p* and *Vin* represents static pressure and velocity at the pump inlet respectively, *pv* is the vapor pressure, ρ is density, and *g* is the gravitational acceleration [21].

#### *2.2. Design of Experiment*

In deciding the optimization variables, the inlet diameter of the impeller, *D*1, the impeller diameter at the outlet, *D*2, the diameter of the hub *Dh*, and the blade width at outlet *b*2 were held constant. This was to maintain the shape of the impeller due to space constrains of the suction and volute casing of the double-suction centrifugal pump. The choice of decision variables was limited to the blade inlet profile at hub and shroud. Three geometrical parameters, namely, hub inlet angle β1 hub, middle inlet angle β2 middle, and the shroud inlet angle β3 shroud were selected as optimization variables. Each parameter was given a set of five values. Table 1 shows the decision variables and their set of values used for the parameterization. Orthogonal design of experiment was applied here to design the experimental scheme. From Table 1, there were 3 factors and 5 levels for the design. An orthogonal scheme of L25 (35 ) was designed according to Equation (2). Twenty-five impeller models were generated based on the sampling results from the DOE using ANSYS bladegen. Table 2 is the

orthogonal scheme and the meridional shape of the impeller is shown in Figure 2. The position of the inlet profile was adjusted by a Bézier curve with three points.

$$L\_m(n^p) \tag{2}$$

**Table 1.** Range of design variables.

**Table 2.** Orthogonal scheme.

**Figure 2.** Meridional shape.

#### *2.3. Surrogate Training and Optimum Solution Solving*

Artificial neural network (ANN) was used for the surrogate modeling. The concept of ANN and MOGA has been explained in the previous works [13,14]. The Levenberg–Marquardt algorithm was adopted as the training algorithm since it has the fastest convergence rate. Tangent hyperbolic activation function (*tanh*) was applied because it has a much better recognition accuracy for multi-layer neural networks. For this study, a dual-layer feed-forward ANN with 10 hidden neurons was adopted. Figure 3 is the multi-layer feedforward ANN structures used for training the NPSHr objectives.

**Figure 3.** Artificial neural network (ANN) structure for the required net positive suction head (NPSHr) objective.

The mathematical relation for the ANN function is written as Equation (3), the activation function, *tanh* is written as Equation (4) and the linear function as Equation (5).

$$y = g\left(\sum\_{j=1}^{n} w\_j^2 \times f\left(\sum\_{k=1}^{m} w\_{k,j}^1 \mathbf{x}\_k + b\_n^1\right)\right) + b^2 \tag{3}$$

$$f(\mathbf{x}) = \left[\frac{2}{(1 + e^{-2\mathbf{x}})}\right] - 1\tag{4}$$

$$g(\mathbf{x}) = a\mathbf{x} + b \tag{5}$$

To obtain the global Pareto frontier for the three-objective functions, the multi-objective genetic algorithm was applied to solve the problem. The problem for the three-objective optimization can be described as follows.

$$\begin{cases} \begin{array}{c} \text{minimize NPSHr } 0.8Q\_d = f\_1\\ \text{INPSHr } 1.0Q\_d = f\_2\\ \text{NPSHr } 1.2Q\_d = f\_3 \end{array} \\ \begin{array}{c} \text{subject to} \\ 17^\circ \le A \ge 25^\circ\\ 13^\circ \le B \ge 23^\circ\\ 11^\circ \le C \ge 21^\circ \end{array} \end{cases} \tag{6}$$

Since the problem was a three-objective problem, a 3D Pareto optimal solution was required. In construction of the Pareto solutions, the following input parameters were used. Population size was 200. A 0.8 Pareto-front population and 0.85 crossover fraction were applied with 1000 generations at a function tolerance of 10.

#### **3. Tested Pump and Computational Domain**

#### *3.1. Description of Test Pump*

The test object is a 250GS40 double-suction centrifugal pump that has a specific speed of 89.5. The flow domain has been divided into the semi-spiral suction domain, a shrouded impeller which has six twisted blades and a volute which also serves as the outlet domain. The test object is shown in Figure 4. Figure 4a shows the computational domain of the volute and suction unit and Figure 4b shows that of the impeller. The design parameters of the pump are presented in Table 3. The specific speed was derived as follows:

$$m\_s = 3.65 \times \frac{N \times \sqrt{Q\_d}}{60 \times H^{\frac{3}{4}}} \tag{7}$$

where *ns* is specific speed, *N* (rpm) represents rotating speed, *Qd* (m3 /h) is the nominal flow rate, and *H* (m)*,* is the design point head.

**Figure 4.** Computational domain of impeller, suction and volute.

**Table 3.** Design specifications of model pump.

#### *3.2. NPSHr Prediction Procedure*

A novel NPSHr prediction method originally introduced by Ding et al. [19] was applied to speed up the calculation process due to the large sample size of the optimization cases. This new method has been compared with the traditional method of NPSHr prediction using the multi reference frame (MRF) model [20], and further applied to cavitation studies in the double-suction centrifugal pump by Pei [22] and Wang [14]. There are three main prediction steps, and a convergence method, which speeds up to the desired accuracy where convergence is not reached after the third step.

Step 1: The pump head is calculated by introducing a new boundary pair, which is flow rate at the inlet and static pressure at the outlet. The static pressure at the outlet *pout* is estimated as follows.

$$P\_{\rm out}(1) = H(0) + P\_{\rm tot-in}(0) - H\_D(0) = P\_{\rm out-in}(0) \tag{8}$$

*ptot*−*in* here is the total pressure at the inlet.

Step 2: In step 2 the static pressure at the outlet is calculated to correspond to a 3% drop in the pump head. Here the static pressure is estimate as

$$P\_{\rm out}(2) = 0.97 H\_{100} - H\_D(1) \tag{9}$$

Here, *H*100, the head at 0% drop in head and *HD* is the dynamic head.

Step 3: The third step is to adjust the errors to converge at the 3% head drop. The outlet static pressure for this step is calculated as

$$P\_{\rm out}(\mathfrak{Z}) = \ 0.97 H\_{100} + P\_{\rm tot-in}(\mathfrak{Z}) - H\_D(\mathfrak{Z}) \tag{10}$$

Step (*n* + 1): This step is applied when convergence is not achieved after the third step. This step has been explained in detail by Ding [19] and Pei [20]. Outlet static pressure is calculated as

$$P\_{\rm out}(n+1) = 0.97H\_{100} + n\_3 - H\_{\rm D}(n) \tag{11}$$

#### *3.3. Governing Equations*

The time dependent Navier–Stokes equation is derived from the continuity equation [23] is given as

$$\frac{\partial \rho\_m}{\partial t} + \frac{\partial}{\partial \mathbf{x}\_j} (\rho\_m u\_j) = 0 \tag{12}$$

$$\frac{\partial(\rho\_{\text{m}}\mu\_{i}\mu\_{j})}{\partial\mathbf{x}\_{j}} + \frac{\partial}{\partial\_{\text{l}}}(\rho\_{\text{m}}\mu\_{i}) = \frac{\partial P}{\partial\mathbf{x}\_{i}} + \frac{\partial}{\partial\mathbf{x}\_{j}} \left[ (\mu + \mu\_{l}) \left( \frac{\partial u\_{i}}{\partial\mathbf{x}\_{j}} + \frac{\partial u\_{j}}{\partial\mathbf{x}\_{i}} + \frac{2}{3} \frac{\partial u\_{k}}{\partial\mathbf{x}\_{k}} \delta\_{ij} \right) \right] \tag{13}$$

Shear stress transport (SST *k*-ω) was applied since it has the advantages of both the *k-*ω and *k-*ǫ turbulence models [24,25]. To simulate cavitating flows, the equation of mass transport is built from the Rayleigh–Plesset equation, and is expressed as

$$\frac{\partial(\rho\_{\upsilon}\alpha\_{\upsilon})}{\partial\_{t}} + \frac{\partial}{\partial \mathbf{x}\_{j}}(\rho\_{\upsilon}\alpha\_{\upsilon}\boldsymbol{u}\_{j}) = \dot{\boldsymbol{m}} = \dot{\boldsymbol{m}}^{+} - \dot{\boldsymbol{m}}^{-} \tag{14}$$

$$m^{+} = \mathcal{C}\_{\text{cap}} \frac{3r\_{\mathcal{S}}(1 - \alpha\_{\upsilon})\rho\_{\upsilon}}{R\_{b}} \sqrt{\frac{2}{3} \frac{\max(P\_{\upsilon} - P\_{\prime} 0)}{\rho\_{l}}} \tag{15}$$

$$m^{-} = \mathcal{C}\_{\text{cond}} \frac{3\alpha\_{v}\rho\_{v}}{R\_{b}} \sqrt{\frac{2}{3} \frac{\max(P - P\_{v\nu}0)}{\rho\_{l}}} \tag{16}$$

From literature, the standardized values are *Cvap* = 50, *Ccond* = 0.01, *rg* = 5 × 10−4 , *Rb* = 10−6 m, ρ*v* = 0.554 kg/m3 , ρ*l* = 1000 kg/m3 , and *pv* = 3169 Pa [26,27].

#### *3.4. Test of Grid Independence*

For maximum simulation accuracy, the test object was meshed using high quality structural hexahedral mesh with ANSYS ICEM. To attain higher precision and boundary motion features, the grids were refined with large numbers and concentrated near the walls. Figure 5 presents an overview of the mesh topology. The volute tongue can be seen in Figure 5a, the impeller in Figure 5b, and suction tongue in Figure 5c. After building the mesh, a test of grid independence was carried out by building 5 independent grids to speed up the calculation time and still maintain accuracy. This was done at the nominal flow condition. It has been established that performance parameters such as head and efficiency are not sufficient enough during the test for grid independence [22,28]; hence, in this research, the influence of the different grid numbers on the volute wall pressure was considered. Monitor points were set in the volute for investigation as shown in Figure 6 and presented in Table 4. The effect of the mesh density on the pump head, efficiency and pressure was least felt at monitor points V1 and V7 as the total mesh number reached 4,266,423. Mesh III was therefore adapted for the numerical simulations to reduce computation load and time. The final mesh has an average Y-plus less than 50, and the distribution can be seen in Figure 7.

**Figure 5.** Mesh of calculation domain: (**a**) volute tongue, (**b**) impeller, and (**c**) suction tongue.

**Figure 7.** Yplus distribution on blade surface.

**Table 4.** Grid cells of the selected mesh.

#### *3.5. Numerical Calculation Setup*

The computational domain was calculated using ANSYS 2019 R3 to solve the Reynolds-averaged Navier–Stokes equations. Water was used as the working fluid while assuming an isothermal heat transfer of 25 ◦C. The SST *(k-*ω) turbulence model was applied since its combined advantages can be extended by the automatic wall treatment to ensure the accuracy of the pressure gradient regardless of the distance to the nearest wall [24,25]. The direction of the flow was set normal to the boundary condition, and a smooth non-slip wall was applied. The default turbulence intensity level of 5% was applied at inlet based on similar works on centrifugal pumps with similar characteristics [14,29]. A high-resolution upwind scheme was applied to ensure accuracy and convergence consistency. At the domain interfaces, a frozen rotor with pitch angle of 360◦ was set between the rotor and the stator. The boundary conditions were specified as flow rate at inlet and static pressure at the outlet. The Zwart–Gerber–Belamri (ZGB) model was used for cavitation simulations. To determine the convergence criteria, maximum residual values of 10−5 were maintained, and iterations had to be periodic stable. Performance indicators were calculated as follows [21].

The Head:

$$H = \frac{P\_{\text{tot-out}} - P\_{\text{tot-in}}}{\rho \text{g}} \tag{17}$$

*ptot-out* and *ptot-in* are the total pressures at outlet and inlet.

Shaft power:

$$P\_s = T\omega \tag{18}$$

*T* denotes shaft torque while ω, is the angular speed. Efficiency:

$$\eta = \frac{\rho g H \mathcal{Q}}{P\_s} \tag{19}$$

#### *3.6. Description of Test System*

The experimental test was performed in an open test rig system. Test schematics can be found in Figures 8 and 9 is the experimental setup system. The suction pipe has a diameter of 250 mm while the delivery pipe has a diameter of 200 mm. High precision WT200 pressure transmitters were installed directly in the pipes at the suction and discharge points. The pressure transmitter at the suction end had a measuring range of −0.1 to 0.1 MPa, whereas the outlet transmitter had a measuring range of 0–1.6 MPa. Their measurement uncertainty was 0.5%. The flow rate was measured with LWGY-200A electromagnetic flowmeter with an uncertainty of 0.07%. Hydraulic and cavitation tests were performed at design and non-design flow conditions. To determine the experimental uncertainty, the Type B evaluation was used, and the systematic uncertainty of the experiment was ±0.411%. Further details of the experiments have been provided in our previous studies [20].

**Figure 8.** Schematics of the test setup. 1: Inlet pipe, 2(8): Valve, 3(6): Pressure transducer, 4: Tested pump, 5: Driven motor, 7: Magnetic flow meter, and 9: Outlet pipe.

**Figure 9.** Experimental setup.

#### **4. Discussion of Results**

#### *4.1. Validation of Numerical Method*

To prove the reliability of the numerical approach, the test results were compared with the computation results and presented in Figure 10. The efficiency curve rose gently from deep part loads to the maximum point of 86.63% at 41.49 m head, and started a steeper decline. The trend was same for both the experiments and numerical results with deviations of 3.01% and 2.03% between the efficiency and head respectively. In the numerical simulations, a smooth wall was assumed. Practically, it is very expensive and difficult to achieve a very smooth surface, and this could account for the minimal deviations between the tests and simulation results. The numerical approach was therefore said to reliable and suitable for the optimization studies since it agreed with the experiments.

**Figure 10.** Validation with test results.

#### *4.2. Cavitation Model Validation*

The cavitation model was also validated with the test results to further confirm the reliability of the numerical method for cavitation flow simulations. At the nominal flow condition, (Figure 11a), the NPSHr was calculated in only three simulation attempts, which resulted in 3.12% drop in head at NPSHr = 2.628 m. At non-design flow conditions (Figure 11b), the suction performance for the

experimental test was well predicted by the numerical results and the maximum deviation was 4.05%. For all flow conditions, the NPSHr did not exceed the design NPSHr of the pump. At 0.6*Qd* the computational results were higher than the test values but lower at design and overload conditions. This could also be attributed to the assumption of a smooth non-slip wall during the calculations. The minimum deviations confirm that the cavitation model can be reliable for NPSHr prediction for the optimization studies.

**Figure 11.** Comparison of suction performance with test results at design and off design points.

#### *4.3. Results from Optimization Studies*

#### 4.3.1. Orthogonal Test Results

The 25 impellers designed from Table 2 were simulated for the three flow conditions. Three of the parameter combinations were invalid. Table 5 gives the upper and lower bounds for the decision variables used for surrogate training, and the results from the orthogonal test is presented in Table 6. From direct analysis, it is observed that, majority of the schemes can satisfy the individual objectives; however, deciding the best scheme from the orthogonal results is difficult. It becomes necessary to apply surrogate models to search for the optimum parameters that would meet the objectives.

**Table 5.** Boundaries and variables for ANN.

#### 4.3.2. Training of Surrogate Model

Artificial neural networks were adopted as the metamodel to build the relationship between the objective functions and decision variables. ANN was then trained to study the relationship between the input and output parameters. The strength of the relationship between the linear model and the dependent variable was determined through R2 analysis. This determines the suitability of the surrogates for further optimization. Figure 12 shows the R2 of the ANN models of the objective functions calculated as 0.9877 for 0.8*Qd* (Figure 12a), 0.9805 for 1.0*Qd* (Figure 12b), and 0.9759 for 1.2*Qd* (Figure 12c), indicating that the ANN models for could be applied to the three-objective optimization since the prediction accuracy is high enough. Validation of a surrogate is a requirement hence a comparison of the ANN prediction and the CFD simulation results from Table 6 is drawn in Figures 13–15. From the graphs, the predicted ANN model values are in agreement with the CFD simulation values. The three conditions can all be used as objective functions in the three-objective optimization.

**Table 6.** Orthogonal scheme results.

**Figure 12.** R 2 analysis of NPSHr at 0.8*Qd* , 1.0*Qd* and 1.2*Qd* .

**Figure 13.** Validation of ANN NPSHr prediction with CFD at 0.8*Qd* .

**Figure 14.** Validation of ANN prediction with CFD at 1.0 *Qd* .

**Figure 15.** Validation of ANN prediction with CFD at 1.2 *Qd* .

#### 4.3.3. Solutions to the Three-Objective Problem

The Pareto frontiers from ANN for the three objective functions have been presented in Figure 16. The optimum solutions were presented as blue points in the Pareto optimum front, and each point represents a model impeller. The Pareto solutions presented a set of 200 optimized impellers schemes that satisfied the optimization objective. Three best impeller schemes were selected and built according to the optimum decision variables in Table 7 and calculated by CFD.

**Figure 16.** Pareto-frontiers from ANN.

**Table 7.** Variables for optimum cases.

Table 8 presents the CFD results from the optimized cases. The CFD predictions were very close to the 3D Pareto predictions rendering the results valid. At part load and design conditions, optimized case 1 performed better than the other cases. At overload condition, optimized case 2 performed better that optimized case 1. For all three cases, optimized case 3 had the worst cavitation performance although the optimization objective was achieved in all cases. Moreover, the head at the design condition was compared for the three flow conditions.  Comparing the three cases, the best single case for the optimization objective is case 1. The head comparison at design point is in Table 9.

**Table 8.** Comparison of original case and best optimal cases.

#### 4.3.4. Comparison of Suction Performance—Optimized and Original Design

The suction performance of the three optimum cases were compared in Figure 17. In each optimized condition, the results were good. All the three optimum cases performed better, however optimized case 1 performed best at all the flow conditions. At 0.8*Qd,* the effect of the optimization was less for all flow conditions. The selected optimum case was compared with the original model in Figure 18. It was seen that for the original pump, reducing the hub inlet angle from 17◦ to 16.77◦ and increasing the hub inlet angle from 15◦ to 19.49◦ gave a remarkable improvement of 6.9% in suction performance at 1.0*Qd*. Cavitation performance was improved by 3.5% at overload condition, 4% 0.8*Qd* and 5% 0.6*Qd* at non-design conditions. Comparing with the test results the optimization achieved its objective of improving the cavitation performance at non-design points.

**Figure 17.** Suction performance of optimized cases and experiments.

**Figure 18.** Verification of optimized case 1 with original model.

#### 4.3.5. Internal Flow Analysis

The internal flow of the optimized impeller (case 1) was compared with the original model to reveal the optimization effect on the flow structure. Figure 19 compares the variation of pressure on the blade surface for the three flow conditions that were optimized. At 0.8*Qd* flow rate, the pressure contour for the original and optimized blades were almost uniform. Looking at the pressure variation the difference was very minimal. At the design condition, the very low pressure areas at the blade leading regions had reduced after optimization. Compared to the original model the low pressure regions covered a wider area than in the optimized model which is an indication that cavitation performance could improve. At overload condition, the low pressure regions within the optimized blade had narrowed as well. The difference was quite clear as compared to the design condition and part load. The extremely low pressure contours at the leading regions had reduced too. For the three flow conditions, the low pressure distribution on the blade leading regions were reduced after optimization and this can translate to an improvement in suction performance.

**Figure 19.** Pressure distribution on blade surface at 0% head drop.

In Figure 20, the pressure fluctuation on the blade surface at critical cavitation conditions were compared. The comparison was done at the critical cavitation point of the original impeller (at NPSH = 2.16 m, NPSH = 2.35 m, and NPSH = 3.16 m), for 0.8*Qd*, 1.0*Qd*, and 1.2*Qd*, respectively. For all the flow conditions, the variation of pressure along the blade surface was improved as compared to the original model. The original model had lower blade surface pressures than the optimized case, indication that suction performance would improve for all the three conditions investigated. The streamline distribution in the impeller flow channel was also analyzed for the various flow conditions in Figure 21 at no cavitation condition. At 0.8*Qd* part load condition, the flow was distorted in two impeller channels of the original impeller. In one of the channels, some vortex was observed between the leading region and the middle of the impeller. Towards the trailing region the flow separated before normalizing. For the optimized impeller, the flow was uniform, and there was no vortex or separation. The streamline in the channels were uniform in both the original and optimized models at 1.0*Qd* and 1.2*Qd*, respectively.

#### 4.3.6. Attached Cavity Distribution in the Flow Domain

To clearly see the difference in attached cavity distribution, the critical cavitation point in the original model for each flow rate was used as the reference condition. The cavity distribution in the impeller, suction, and volute of the optimized impeller was compared to the original impeller for that point. Figure 22 is a comparison of the attached cavity distribution in the impeller for different flow conditions. At 0.8*Qd*, the attached cavities in the impeller reduced slightly at NPSH = 2.176, which is the critical cavitation point of the original model. At the design condition, a clear difference was observed. The bubbles covered the entire leading regions in the original impeller. In the optimized model, however, the cavity distribution had reduced, and portions of the leading edge were free from cavitation bubbles at NPSH = 2.53 m. The overload condition was similar to the design condition. Some parts of the leading edge were free from attached cavities; however, the difference in cavity distribution between the two impellers at overload condition was less.

**Figure 21.** Streamline pattern in impeller channel.

In the suction domain, cavitation occurred around the suction rings. At part load, the attached cavity distribution was not symmetric in both designs. The suction ring at impeller inlet 1 (Figure 4) had cavities made up of bubbles forming around it. It was, however, observed after the optimization that the attached cavities in the suction unit with the optimized impeller had reduced at the critical cavitation point (NPSH = 2.176 m) of the original model. Cavity distribution in the suction at the design point was compared at NPSH = 2.532 m. This was similar to part load condition. Bubbles appeared on suction ring at impeller inlet 1, and there was a reduction in the attached cavities around the suction ring. For 1.2*Qd*, the original impeller had bubbles on both rings. In the optimized model, the cavities had reduced to a single suction ring when the comparison was done at NPSH = 3.39 m. Only the suction ring at impeller inlet 2 had attached cavity distributions. Generally, there was reduction in

cavitation bubbles for all flow conditions. This meant that the suction performance had improved for all three conditions and the optimization target has been achieved. Figure 23 shows the bubble distribution in the suction domain.

**Figure 23.** Attached cavity distribution in the suction.

#### **5. Conclusions**

Cavitation has been a treat to cooling water circulation pumps leading to reduced efficiency and damage of pump impellers. Hence, to improve the anti-cavitation performance and ensure longer life and reliability of cooling water circulation pumps, a multi-parameter and a multi-condition optimization approach was designed. In the optimization design, a dual-layer feedforward artificial neural network and genetic algorithm was then applied to improve the cavitation performance of the double-suction impeller over a wider range of operating conditions. The results from the study are as follows:

**Author Contributions:** Conceptualization, W.W. and Y.L.; methodology, M.K.O. and S.Y.; software, M.K.O.; validation, J.L. and B.Z.; formal analysis, W.W., Y.L., M.K.O., and B.Z.; writing—original draft preparation, M.K.O.; writing—review and editing, M.K.O. and W.W.; and project administration, S.Y. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work is supported by the Natural Science Foundation of Jiangsu Province (Grant No. BK20190851), Primary Research & Development Plan of Shandong Province (Grant No. 2019TSLH0304), Natural Science Foundation of China (Grant No. 51879121, 51579104) and Primary Research & Development Plan of Jiangsu Province (Grant No. BE2019009-1).

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Research on the Application of Partial Similarity for a 1-1**/**2 Axial Compressor**

#### **Hong Xie , Moru Song and Bo Yang \***

School of Mechanical Engineering, Shanghai Jiaotong University, Shanghai 200240, China<EMAIL_ADDRESS>(H.X<EMAIL_ADDRESS>(M.S.)

**\*** Correspondence<EMAIL_ADDRESS>Tel.: +86-21-3420-6871

Received: 18 August 2020; Accepted: 7 September 2020; Published: 9 September 2020

**Abstract:** In this paper, a method based on the partial similarity principle is presented to improve the aerodynamic design with low cost and high accuracy for a 1-1/2 axial compressor. By means of this method, during the process of a similar design, the machine Mach number and flowrate coefficient are maintained. The flow similarity between the prototype and its large-scaled alternative was observed, according to a detailed analysis of flow fields of rotor and stator. As well, the relative discrepancies of isentropic efficiency and pressure ratio between two models are 1.25% and 0.4% at design point, respectively. Besides, their performance curves agreed very well in the whole operating range. Moreover, it was also found that the flow similarity between the two models can be maintained under unsteady working conditions. Thereafter, in order to investigate the impact of stability optimization method on the similarity principle, casing treatment with single circumferential groove was applied to these two models. The flow similarity was still maintained and the flowrate near the stall was reduced about 1.1% with negligible deterioration of the overall performance.

**Keywords:** partial similarity principle; flow similarity; stability improvement

#### **1. Introduction**

A high-pressure axial compressor is one of the three major parts in the state-of-the-art aero-engine, and this kind of compressor is one of the most complicated products among all compressor technology [1]. With the increasing performance and stability demands, it is vital to obtain the in-depth understanding of its internal intricate flow regimes [2]. Under the current experimental conditions, it is too difficult and hazardous to conduct experiments on the high-pressure and high-speed compressors. Consequently, low-speed and large-scale model testing developed into an achievable and reliable method, due to the benefits of low cost, better accuracy and relatively low risk [3]. The prototype can be scaled-up and its rotation speed is reduced accordingly; then thorough measurements can be realized by easily to investigate flow interactions.

The concept of low-speed model testing was first presented by Wisler in the study of exit stages in the core compressor [4]. Then, this method was applied for loss reduction and performance improvement with structural optimization [5]. Robinson improved the end-wall flow pattern in a four-stage, low-speed axial compressor with end-bent blading technology [6]. Lyes optimized the high-pressure compressor through research on blade bending and sweeping based on a low-speed model [7], which was also used for sweep and dihedral blading studies by RR Company and Cambridge University [8,9]. Besides, Boos and Lange also carried out a high loaded blading investigation on this large-scale, low-speed compressor [10,11]. Zhang to reached performance optimization of a low-speed axial rotor by three-dimensional (3D) optimization with numerical computation, and then the experiments were carried out for validation of the numerical results [12,13]. Zhang summarized four procedures with which to achieve the whole process of compressor scaling and performance optimization. The procedures included establishing a low-speed model that behaves similarly in a flow field to the high-speed compressor; detailed and accurate flow field measurement in the low-speed model; revising the design for the performance improvement; and the verification of new design methodology for the high-speed compressor at last [14]. The research above all employed the same modeling principles as that of Wisler [3]. These principles focused on some geometric parameters and aerodynamic parameters, such as the blade surface pressure and velocity distribution, and ignored the Mach number and Reynolds number.

However, the model principles mentioned above brought about some shortages. At first, the flow field is complex in 3D space, so it is not adequate to just keep the surface aerodynamic parameters the same. Then, the low-speed four-repeating-stage, studied by Wisler [3] and Zhang [14], required achieving a repeating condition for the similar inlet and outlet aerodynamic parameters. However, only the third stage can be used for the investigation, leading to wasted time and money. Additionally, this low-speed model can only be used for the investigation of a single stage, and cannot be used for multi-stage investigations. In addition, this model can be only used for the investigation at design point, but investigations at off-design points need to be similar. Besides, this kind of scaling process, based on the inverse method, is time-consuming and complicated. Finally, there is an unavoidable shortcoming of this model principle. The principle cannot be used for the transonic compressor or stage, because the shock waves' effects cannot be evaluated. Consequently, it is necessary to propose a simple, convenient, wider-in-applicability and relatively less disadvantages methodology for high-to-low-speed compressor transformation.

In recent decades, the scaling method for the low-speed model had been adopted by many researchers. The effects of Re on performance similarity were theoretically studied by Ma [15], and it has been proven that it is very important to maintain Ma in some cases [16,17]. In order to fulfill complete flow similarity, Ma and Re need to be identical for a high-speed compressor and low-speed model [18,19]. Nevertheless, it is almost impossible to maintain Re and Ma simultaneously, resulting in a performance difference between the prototype and the low-speed model. In practical applications, effects caused by changed Re can be compensated by some modifications of structural parameters [20,21]. Then, the Re could be ignored and the scaling process could be simplified accordingly, and this scaling laws can be summarized as a partial similarity principle. In this study, the partial similarity principle is innovatively used for the high-to-low transformation in the axial compressor.

The optimization of the compressor is the ultimate objective after successful establishment of the low-speed model. The stable operating range is one of the most vital objectives, and the delay of stall is an efficient method. Spike disturbance is often occurred in a low-speed compressor [22], and there are two criteria necessary for the emergence of spike disturbances, both of which are related to the tip clearance flow. One of the criteria is that the interface between the tip leakage flow and main flow aligns with the leading-edge plane [23]. Casing treatment technology has been proven to be a successful method for enhancing the stability margin of an axial compressor [24]. The interface can be delayed downstream under the impact of the casing treatment, resulting in better stall margin.

Two kinds of casing treatment are most commonly investigated. The first employs a lot of slots which are cut into the end-wall [25]. The stall margin of a compressor with slots can be improved largely, but with a large efficiency penalty [26]. The other one employs a single or several grooves which are cut into the end-wall circumferentially [27]. The stall margin improvements of a compressor with grooves are smaller than those of the first one, but the efficiency penalty is decreased [28]. Zhao concluded that the axial location had the most important effect on the stall margin improvement [29]. Houghton and Day obtained two optimum axial locations by applying a sequence of single-groove casings in two low-speed compressors—one was near the leading edge and the other near the mid-chord [30]. Du depicted that the optimal groove location in the low-speed compressor is 57% along the axial chord nearing the mid-chord [31], which coincided with the experiment results of Bailey [32].

In this study, a low-speed and large-scale axial compressor was fabricated from a high-pressure and high-speed prototype based on the partial similarity principle, the similar flow fields of which were validated by numerical computation. Then, in order to investigate whether the flow similarity can be contained when using the similarity principle, casing treatment with single circumferential groove was applied to these two models. Meanwhile, in this study, performance experiments were carried out to validate the reliability of the numerical method. The aim of this paper is to propose a method which can be used conveniently and quickly for a type of compressor with small dimensions and high rotation speed, rather than a certain compressor.

#### **2. Scaling Laws**

The design of low-speed model is different from general compressor design. The low-speed model is modeled from a high-speed prototype based on the similarity scaling laws. In general, the scaling method can be divided into exact similarity method and partial similarity method. The key difference between the two methods is whether Re remains constant or not. In this study, the Re of the prototype is calculated as Equation (1).

$$\text{Re} = \frac{\text{U} \text{b}\_2}{\nu} \tag{1}$$

where ν is the inlet kinematic viscosity.

In the practical application of the partial similarity principle, there are two dimensionless parameters that need to be kept constant [33]: the mass flow coefficient φ and the machine Mach number *MaU*. The two parameters can be expressed as follows.

$$\phi = \frac{\mathcal{Q}}{\pi/4 \cdot \mathcal{U} \cdot \mathcal{D}\_2^{-2}} \tag{2}$$

$$
\Delta q\_{II} = \frac{\mathcal{U}}{\sqrt{\kappa RT\_0}}\tag{3}
$$

The complete dimensional analysis process, based on the Buckingham–PI theorem, is described as follows. Generally, there are 7 variables overall when applying the partial similarity principle in this case, which are ρ, *u*2, *D*2, *P*, *T*,*t*, µ. The length dimension *L*, mass dimension *M* and time dimension *t* are chosen to be the base dimensions. Simultaneously, ρ, *u*2, *D*2 are selected to be the base variables. Then, the flow characteristics and the dimensional analysis can be depicted as follows.

$$f(\rho, u\_2, D\_2, P, T, t, \mu) = 0,\tag{4}$$

where µ is the dynamic viscosity. The deduction of similarity criterion can be described as follows.

$$\begin{cases} \pi\_1 = \rho^{a\_1} u\_2^{b\_1} D\_2^{\mathbb{C}\_1} P \\ \pi\_2 = \rho^{a\_2} u\_2^{b\_2} D\_2^{\mathbb{C}\_2} T \\ \pi\_3 = \rho^{a\_3} u\_2^{b\_3} D\_2^{\mathbb{C}\_3} t \\ \pi\_4 = \rho^{a\_4} u\_2^{b\_4} D\_2^{\mathbb{C}\_4} \mu \end{cases} \tag{5}$$

$$\begin{aligned} \left[\pi\_1\right] &= M^0 L^0 t^0 = \left(ML^{-3}\right)^{a\_1} \left(Lt^{-1}\right)^{b\_1} \left(L\right)^{c\_1} \left(ML^{-1}t^{-2}\right) \\ \left[\pi\_2\right] &= M^0 L^0 t^0 = \left(ML^{-3}\right)^{a\_2} \left(Lt^{-1}\right)^{b\_2} \left(L\right)^{c\_2} \left(L^2t^{-2}\right) \\ \left[\pi\_3\right] &= M^0 L^0 t^0 = \left(ML^{-3}\right)^{a\_3} \left(Lt^{-1}\right)^{b\_3} \left(L\right)^{c\_3} t \\ \left[\pi\_4\right] &= M^0 L^0 t^0 = \left(ML^{-3}\right)^{a\_4} \left(Lt^{-1}\right)^{b\_4} \left(L\right)^{c\_4} \left(ML^{-1}t^{-1}\right) \end{aligned} \tag{6}$$

$$\begin{cases} a\_1 = -1, b\_1 = -2, c\_1 = 0 \\ a\_2 = 0, b\_2 = -2, c\_2 = 0 \\ a\_3 = 0, b\_3 = 1, c\_3 = -1 \\ a\_4 = -1, b\_4 = -1, c\_4 = -1 \end{cases} \tag{7}$$

$$\begin{cases} \pi\_1 = \rho^{-1} u\_2^{-2} D\_2^0 P = \frac{P}{\rho u\_2^2} = \frac{\rho \text{KT}}{\rho u\_2^2} = \frac{\text{KT}}{u\_2^2} = \frac{1}{\kappa} \frac{\kappa \text{RT}}{u\_2^2} = \frac{1}{\kappa} \left(\frac{1}{Ma}\right)^2\\ \pi\_2 = \rho^0 u\_2^{-2} D\_2^0 T = \frac{T}{u\_2^2} = \frac{1}{\kappa \text{R}} \frac{\kappa \text{RT}}{u\_2^2} = \frac{1}{\kappa \text{R}} \left(\frac{1}{Ma}\right)^2\\ \pi\_3 = \rho^0 u\_2^1 D\_2^{-1} t = \frac{u\_2 t}{D\_2} = \frac{\frac{\pi}{4} D\_2^2 u\_2}{\frac{\pi}{4} D\_2^2 \frac{u\_2}{2}} = \frac{\frac{\pi}{4} D\_2^2 u\_2}{Q} = \frac{1}{\phi} \\ \pi\_4 = \rho^{-1} u\_2^{-1} D\_2^{-1} \mu = \frac{\mu}{\rho u\_2 D\_2} = \frac{1}{\text{Re}} \end{cases} \tag{8}$$

Finally, the flow characteristics are summarized as Equation (8).

$$F\left(\frac{1}{\kappa}\left(\frac{1}{Ma}\right)^2, \frac{1}{\kappa R}\left(\frac{1}{\mathbf{Ma}}\right)^2, \frac{1}{\phi'}, \frac{1}{\mathbf{Re}}\right) = 0,\tag{9}$$

Consequently, in this case, if the mass flow coefficient φ and machine Mach number *Mau*2 are constant and the Re is ignored, the flow similarity can be achieved in partial similarity principle.

The aerodynamic parameters at the design point of the prototype were chosen as the modeled conditions. In this study, the scaling factor is 10, so the speed of scaled-up model is low enough according to the scaling laws. The scaled-up compressor possesses completely the same relative geometric parameters as the prototype, including solidity, aspect ratio, hub-to-tip ratio, axial-space-chord ratio, airfoil stacking, relative radial clearance, etc. In the following study, the scaling process is qualified by not only the performance curve, but also the detailed flow field characteristics of the prototype and scaled-up model. The procedure, including scaling-up and operating range optimization, is shown in Figure 1. The whole process can be divided into the scaling-up process and the improving process. Firstly, the prototype is scaled-up to be a large-scale and low-speed model, which achieves enough flow similarity with the prototype. Then, the same casing treatment is used in the prototype and large model for the performance optimization. The following sections introduce those steps one by one.

**Figure 1.** Scheme of performance analysis and improvement.

#### **3. Numerical Method**

Numerical simulations were accomplished with commercial code EURANUS, which contains steady and unsteady solvers by calculating the conservative Reynolds-averaged N-S (Navier-Stokes) equations. It discretized the N-S equations using a cell-centered finite volume formulation in space, and applied Gauss's theorem and central difference method to determine viscous flux. In order to obtain the flow regimes accurately in the flow passage and the detailed flow pattern near the casing, a second-order central spatial discretization scheme was selected to estimate the inviscid fluxes. Frozen rotor technique was employed to deal with the interface between the inlet guide vane (IGV)

and rotor and between the rotor and stator. According to published research [34], Spalart-Allmaras's model can give fairly good results in the numerical simulation of turbomachinery. Thus, in this case, Spalart-Allmaras's model was used for turbulence simulation to evaluate the eddy viscosity [35].

Upstream and downstream, the numerical computational domain was extended, among which mesh was generated with NUMECA Auto Grid 5. The minimum grid spacing of the first layer was set to be 0.001 mm. According to Numeca FineTurbo user guide [36], the y+ is recommended to be below 10 when the Spalart-Allmaras model is used. Besides, the value of y+ is recommended to be below 5 in some published research [35,37]. In this paper, the y+ is about 1, which satisfies the need of the Spalart-Allmaras model. In order to analyze flow field in the tip gap, the number of total layers within tip clearance was set to be 17. A single-passage model was used in the steady computation to validate the consistency between numerical simulations and experiments. However, a single-passage model is not adequate to capture some flow characteristics of the stall inception process, such as the precursors, occurrence, evolution and initiation of the stall and the cross-passage flow. Moreover, the disturbances would propagate circumferentially when approaching the stall point, which plays a significant role in the simulation near the stall. Thus, after validating the reliability of numerical simulations, all the following numerical simulations employed the multi-passage model. As the sliding grid scheme requires two connected regions to have the same circumferential length, the same circumferential length is required for unsteady calculations. The quantities of blades can be changed based on the research of Rai [38], and the numbers of IGVs, rotors and stators were changed to 30, 30 and 30, respectively, so there were 30 passages in the full-annulus. Considering the limitations of the computing resources, a five-passage model was selected. Naturally, the performance of a five-passage model was different from that of the prototype because of the revised blade numbers, but it was considered to be enough to predict the general flow structure for the comparison among and optimization of prototype and low-speed models. The mesh numbers of single-passage and five-passage are about 3 million and 15 million, respectively, which can be seen in Figure 2. The mesh density is appropriate from the validation result of the grid independence. A three level multi-grid was adopted to accelerate the calculation.

**Figure 2.** Mesh of computation domain: (**a**) single-passage model; (**b**) five-passage model.

No-slip and adiabatic conditions were set at all the solid boundaries. Velocity direction, total pressure and total temperature were imposed at the inlet and average static pressure at the outlet as the boundary conditions. The back pressure was increased gradually to obtain the global stable performance maps. The last converged point was inferred to the near-stall point in this study. The steady computation results were used for the performance comparison and analysis, and the unsteady computation was also conducted to capture the unstable characteristics of the flow field. The dual-time step (DTS) method was employed to improve the time marching in the unsteady simulation. The number of angular positions was 10 in one single passage and 20 inner iterations were set for two successive computations, which is enough to get detailed unsteady information.

#### **4. Results and Analysis**

#### *4.1. Evaluation of the Numerical Method*

For the validation of the numerical method, the test was carried out in the 1-1/2 stage axial subsonic compressor in Shanghai Jiao Tong University, as shown in Figure 3 [39]. The blade numbers were 32, and 29 and 37 for the inlet guide vane (IGV), rotor and stator respectively. The section profile of all blades is NACA0012. Hub/tip ratio (ratio of the diameter of blade hub to the diameter of blade tip) was 0.76. The rotational speed was 12,000 rpm, with which the tip Mach number was less than 0.9 at design point. The height of tip clearance was 1% chord length at the blade tip. To acquire the aerodynamic parameters of inlet and outlet, temperature sensors and pressure sensors were arranged at the upstream (plane A) and downstream (plane D) areas, as shown in Figure 3a. On plane A, two static pressure probes and two five-hole total pressure probes (accuracy 0.02%) were installed to obtain inlet static pressure and total pressure, respectively. One five-point temperature probe (accuracy 0.2%) was set to measure the inlet total temperature. The same types sensors were equipped on the plane D to obtain outlet pressure and temperature. In this case, the absolute measurement errors of total pressure and temperature were about 0.02 Kpa and 0.5 K, respectively. The measurement error of prototype-stage isentropic efficiency was about 1.5%. The overall performance of the prototype was acquired with the adjustment of throttle area at outlet. A flowmeter was placed upstream of IGV to measure the flow rate. The geometric parameters of prototype are listed in Table 1.

**Figure 3.** Test facility: (**a**) meridional plane of prototype; (**b**) 3D model of prototype; (**c**) test rig.

1

1

1

1 0 1 0

2

**Table 1.** Design parameters of the prototype.

The experiment and numerical performance versus normalized flow rate are shown in Figure 4. Squares, triangles and circles were applied to make a distinction among the performance curves of experiment results and numerical results of single-passage and five-passage, respectively. The upper lines depict the isentropic efficiency characteristic and the lines below depict the pressure ratio characteristic. In this case, the isentropic efficiency η*s* is calculated as Equation (10).

$$\eta\_{\rm s} = \frac{\left(P\_1/P\_0\right)^{\frac{\kappa - 1}{\kappa}} - 1}{T\_1/T\_0 - 1} \tag{10}$$

where κ is the Adiabatic exponent, *P*0 the inlet total pressure, *P*1 the outlet total pressure, *T*0 the inlet total temperature and *T*1 the outlet total temperature.

**Figure 4.** Overall performances of prototype stages.

Three points are labeled in Figure 4 for depicting the near-stall points conveniently. NSE represents the last point of the characteristic at minimum normalized flow rate Q/Qref=0.77 in experiment, NS1 and NS5 represent the solution limit under the steady computation of single-passage and five-passage, respectively. The computation results of single-passage and experiment results correspond well with each other. The discrepancy at design point is less than 1.5%, and the tendencies of the curves agree well. The major difference at the off-design points is that the computation results overestimate the limit flow rate about 4.2%; that is because of the limitation of the steady simulation. However, the NS1 is very close to the NSE, which means that the steady simulation can exactly capture the point near the stall. The errors exist between the curves of single-passage and five-passage because of the changed

2

blade number. The five-passage model shows better performance in a relatively large flow rate, and the single-passage model shows a better stall margin. In generally, the steady simulation was adequately accurate to obtain the performance and the point near the stall. The five-passage model was used for the follow-up study to capture the circumferential flow characteristics, and the five-passage high-speed model was set to be the prototype accordingly.

#### *4.2. Performance Analysis of the Scaled-Up Compressors*

#### 4.2.1. Comparison of the Prototype and Scaled-Up, Low-Speed Model A at Design Point

The low-speed model A was scaled-up with the scaling factor 1:10 from the prototype geometrically, and the five-passage model was too. Model A maintained the machine Mach number and mass flow coefficient and kept the inlet and outlet parameters. The geometric and aerodynamic parameters of the prototype and model A are listed in Table 2.

**Table 2.** Design parameters of prototype and compressor A.

The distributions of averaged circumferential static pressure coefficient, *Cp*, at five spanwise points of the rotor and stator for the prototype and model A at design point are shown in Figures 5 and 6. *Cp* is defined as

$$\mathcal{C}\_p = \frac{P}{P\_{\text{total,in}}} \tag{11}$$

where *P* and *Ptotal*,*in* denote the surface static pressure and the total pressure at inlet. The static pressure distributions of numerical results of prototype and model A are shown in Figure 5. Overall, the distributions of static pressure of rotor between prototype and model A are almost exactly the same in the whole flow passage. However, there is still small a difference at mid-span, which will be analyzed next. The flow field near the hub and the casing reached good flow similarity, which was not realized by the former studies. The distributions of static pressure of stators also agree well, as shown in Figure 6. Nevertheless, the flow similarity of the stator is not as good as that of the rotors. Moreover, the similarity of the stator over 50% spanwise was worse than that of rotor. Besides, it can be seen that the static pressure rise of stator of model A was larger than that of prototype at all spanwise values.

The comparisons of averaged circumferential values of aerodynamic parameters along the blade height of rotor and stator are shown in Figures 7 and 8. In general, the tendencies of inlet and outlet relative flow angles of the rotor match well. The distributions of flow turning angle near the end-wall are nearly the same. However, the discrepancy of occurs from 5% to 95% blade height, for which the maximum is about 1.5 degrees. The good agreement also appears at the axial velocity distributions at inlet and outlet of rotor. However, the outlet axial velocity of model A is a little smaller than that of the prototype at the same region as the flow turning angle distribution. It can be inferred that the decreased axial velocity attributed the decreasing turning angle, which may be related to the variation of Re. The flow angle and axial velocity distributions of stators still agree well for the prototype and model A, though the error of outlet flow angle is bigger than that of rotor, which may not only be affected by the enlarged Re but the upstream effects of the rotor.

Averaged circumferential isentropic efficiency of the rotor and loss coefficient distributions of rotor and stator along blade height are shown in Figures 9 and 10, respectively. The total pressure loss coefficients of rotor ̟*R* and stator ̟*S* are expressed as:

$$
\Delta \sigma\_{\mathbb{R}} = (P\_{0w} - P\_{1w}) / (P\_{0w} - P\_0) \tag{12}
$$

$$
\rho\_S = \left(P\_{1t} - P\_{2t}\right) / \left(P\_{1t} - P\_1\right) \tag{13}
$$

where *P*0 is the inlet static pressure of the rotor, *P*0*w* is the average inlet relative total pressure of the rotor,*P*1*w* is the average outlet relative total pressure of the rotor, *P*1 is the inlet static pressure of the stator,*P*2*t* is the average outlet absolute total pressure of the stator and *P*1*t* is the average inlet absolute total pressure of the stator. The efficiency of model A is improved from 10% to 85% blade height, and the loss is reduced slightly accordingly. The efficiency and loss curves are approaching. There are two regions with obvious error that can be observed in the loss curves of stators; one is near the hub end-wall and one is from 70% to 100% blade height. The discrepancy may have resulted from the different impact on the downstream stator, caused by the wake flow, leakage flow and main flow. This impact varied with the varied wake flow, leakage flow and main flow along the radial direction, because of the different absolute values of blade height and axial gap between the rotor and stator. Thus, the difference in the upper parts of the downstream stators was brought about. Besides, the difference near the hub end-wall may be relative to the different friction loss, caused by the different boundary effects because of different Re and blade height. , ,

**Figure 5.** *Cont.*

**Figure 5.** Static pressure distributions of the rotor.

**Figure 6.** *Cont.*

**Figure 6.** Static pressure distributions of the stator.

**Figure 7.** *Cont.*

**Figure 7.** Averaged circumferential values of aerodynamic parameters of the rotor.

**Figure 8.** *Cont.*

**Figure 8.** Averaged circumferential values of aerodynamic parameters of the stator.

1 1

2 1

**Figure 10.** Loss coefficient distributions of the rotor and stator: (**a**) loss coefficient of rotor; (**b**) loss coefficient of stator.

According to the comparison and analysis of this section, the low-speed model was successfully fabricated based on the scaling laws. Exact flow similarity was achieved in the rotor and stator, although there was a little error. It is very convenient to realize the scaling process without a complicated iterative process based on the complex inverse method, and the flow similarity of the rotor and stator is reached at the same time. After the successful process at design point, the study on the off-design followed, which was for discussing the flow similarity in the whole operating range.

#### 4.2.2. Comparison of the Prototype and Scaled-Up, Low-Speed Model A in Operating Conditions

The performance comparison between the prototype and model A is shown as a function of normalized mass flow rate in Figure 11. The tendency of performance curves agrees well for the prototype and model A under the operating flow rate. The discrepancies of efficiency and pressure were less than 1.5% and 0.5%, and they were only 1.25% and 0.4% at design point, respectively. It was proven that the flow similarity is not only achieved at design point but for the whole operating flow rate based on the partial similarity principle. However, the flow rate near the stall of model A (the normalized flowrate at this point is namely, NSA) was about 4.2% less than that of prototype (the normalized flowrate at this point is namely, NSA), meaning that the stability is worsened when conducting the scaling-up process.

**Figure 11.** Performance in the whole operating range.

The static pressure and entropy distributions near the casing (99.5% spanwise) of model A at NSA, for prototype at NSP, and for the prototype at NSA, are compared in Figure 12. The red dotted line is static pressure and it is used to detect the trajectory of tip clearance vortex [40]. One of the characteristics for spike disturbance initiating rotating stall is that the interface, where the oncoming and tip leakage flows meet, is aligned to the leading-edge plane. In other words, the tip leakage flow spilled below the blade tip [32]. As shown in Figure 12a, the tip leakage flow of model A at NSA nearly spills below the blade tip, and this flow phenomenon can be also observed in the prototype at NSP. Furthermore, the flow field between model and prototype is similar, so it can be inferred that the flow similarity is maintained near the stall point, and the type of stall inception is also the same. The tip leakage vortex of prototype migrates downstream at NSA, near the stall point of model A. It illustrates that the scaling-up model promotes the tip leakage vortex moving upstream, which results in rotating stalling before the prototype. The entropy distribution can be applied for distinguishing the interface between main flow and tip leakage flow [31]. The black dashed line, where the entropy increases severely, describes the location of the interface. The black line is successive and parallel to the leading-edge plane at NSA in model A and at NSP in the prototype, meaning the inception of stalling. The interface of prototype at NSA is still in the passage, corresponding with the analysis on the static pressure distribution.

Based on the above sections, the high-speed compressor is successfully transformed to the low-speed, large-scale compressor in the whole operating range, although there are some tiny errors. These tiny errors may be caused for several reasons, such as the changed Re, the variations of absolute values of geometric parameters and so on. The further research about diminishing the errors was the content of the next study to achieve further flow similarity. It is worth noting that almost complete flow similarity was achieved under steady conditions, as shown in Figures 5, 7, 9 and 10. Then, it was necessary to study whether the flow similarity was contained in the unsteady state when using the scaling laws, because some important flow characteristics are closely related to the unsteady condition.

**Figure 12.** Comparisons of static pressure and entropy distributions at 99.5% spanwise: (**a**) static pressure distribution; (**b**) entropy distribution.

4.2.3. Comparison of the Prototype and Scaled-Up, Low-Speed Model A under Unsteady Condition

The flow similarity between the high-speed prototype and low-speed large-scale model A is realized under steady condition. The transient normalized relative axial velocity contours at 99.5% and 97% spanwise of prototype and model A at design point are shown in Figure 13. It can be observed that the flow field within the tip clearance and near the casing is stable, and there is not migration of tip leakage vortex. The velocity at 99.5% spanwise of model A is a little larger than that of prototype along the trajectory of leakage vortex, while it is opposite at 97% spanwise. Yet, on the whole, the evolution of the flow field for prototype is similar to that of and model, implying that the flow similarity is preserved under the unsteady condition at design point. The transient characteristics of velocity at 99.5% and 97% spanwise between prototype and model A at the near-stall point are compared in Figure 13.

The instantaneous velocity contours near the stall are shown in Figure 14. The fluctuation of the red dotted curve, detecting the trajectory of tip leakage vortex, indicates that the tip leakage vortex is unstable. The vortex is generated at the leading edge of the blade, and then migrates circumferentially, which can be deduced from the region depicted with A and B. When the vortex impinges the adjacent blades, its intensity is decreased, and the secondary leakage even multi-leakage appears. The red ellipse also presents the movement of the vortex. Again, the unstable characteristics of the prototype and model near the stall are similar, demonstrating that the unsteady flow features are maintained in the design and off-design operating conditions.

**Figure 13.** Transient axial relative velocity of the prototype and model A at 99.5% and 97% spanwise at design point.

**Figure 14.** Transient axial relative velocity of the prototype and model A at 99.5% and 97% spanwise near the stall.

#### *4.3. Stability Improvement through Casing Treatment*

In this study, the single circumferential groove was applied as the casing treatment technology. The location and geometric parameters, and the mesh of the flow passage with groove of model A are shown in Figure 15, and the mesh number of the groove was about 5 million, which is fine enough to provide grid independence.

**Figure 15.** Sketch map of groove geometry and mesh.

The characteristics of the smooth solid wall and the grooved casing treatment are shown in Figure 16. The model A with a groove is named model A CT, and NSACT means the near-stall point of model A CT. The efficiency and pressure ratio of model A were decreased about 0.4% and 0.3% at design point. The reduction of performance decreased along the direction of small flowrate and increased towards the large flowrate.

**Figure 16.** Pressure ratio and isentropic efficiency of a smooth solid wall and grooved casing treatment.

The reason for the decrease can be explained in Figure 17, which shows the apparent increasing entropy around and downstream the region existing groove. The two black lines added in this figure represent the upstream and downstream edges of groove. Besides, it can be observed that the impact of groove not only exists in the tip clearance but in the region near the casing. However, the solution limit near the stall is obviously moved toward the small flowrate, elucidating that the stability improvement of model A was achieved. The reduction of the flowrate near the stall of the prototype was about 1.2%.

**Figure 17.** Entropy distributions at 99.5% and 98% spanwise.

The stability improvement of the prototype was the ultimate object of this research, so the single circumferential groove was also applied to the prototype. The location of the groove was the same as for model A CT, and the geometric parameters were scaled down by the scaling factor 1:10 accordingly. The prototype applying the groove was named prototype CT, and NSPCT means the near-stall point of model A CT. The simulation results for the prototype and prototype CT are shown in Figure 18. The efficiency and pressure ratio of prototype decreased about 0.54% and 0.35% at design point—a little more than for model A, and the tendencies of the curves are similar to those of model A. The reduction of the flowrate of the prototype was 1.1%, which is very close to that of model A.

**Figure 18.** Pressure ratio and isentropic efficiency of the smooth solid wall and grooved casing treatment.

As shown in Figure 19, the entropy also increased apparently around and downstream the region existing groove. The higher entropy distribution around groove of prototype explain the larger performance decrease of prototype. Furthermore, the entropy distribution of prototype CT is similar to that of model A CT, indicating that the effect of groove on the two models is similar and the flow similarity is maintained. The mechanism of stability improvement applying the circumferential groove has been elucidated clearly in many publications [27–29], which is outside the scope of this study, so it will not be covered in this paper. The flow similarity between prototype and model A prove that the research of the casing treatment on the model A can be used for the prototype.

**Figure 19.** Entropy distributions at 99.5% and 98% spanwise.

#### **5. Discussion**

In this study, the flow similarity at the design and off-design points has been achieved well, except that there is a little distinction in the solution limit flowrate near the stall. Thus, further flow similarities will be focused on in the following research, such as the geometric corrections, including the value of tip clearance, the axial gap between rotor and stator, etc. Besides, the casing treatment is one of the effective technologies for stability and performance improvement; there are many other technologies can be used. Furthermore, the detailed measurement of aerodynamic parameters through establishing low-speed large-scale test facility will be carried out after the adequate flow similarity is achieved.

Some comments on the effect of Reynolds number on aerodynamic performance are significant. When scaling-up the prototype, the Re increases as the absolute dimensions of the prototype are enlarged, and the relationship of Re between the prototype and scaled model has been presented in the research of Perter [41]. The friction coefficient value of scaled-up model is decreased with increasing Re, which leads to better efficiency. Nevertheless, the decrease of friction coefficient is finite. The friction coefficient is determined not only by Re but also by the relative height of the roughness of the surface, which is limited. However, CFD computations in this paper were carried out for the cases of smooth surfaces and the roughness of the surface was neglected in this paper. Besides, the displacement thickness of the boundary layer is also effected by the changed Re [42]. Therefore, further investigation is necessary and will be discussed in future research.

#### **6. Conclusions**

In this paper, a method based on the partial similarity principle is proposed. In practical applications, this method can be used conveniently and quickly for the transformation of a high-speed, small compressor. In our case, a 1-1/2 axial compressor has been scaled up to be a low-speed, large-scale model, based on the high-to-low-speed transforming methodology—the partial similarity principle. Thereafter, casing treatment was used for the stability improvement and the flow similarity still remained. Numerical simulations were used for the analysis of the scaled-up process and the stability optimization process. In addition, experimental tests were carried out for the validation of simulation results. On the basis of studies, the following conclusions were drawn.

**Author Contributions:** Conceptualization, H.X. and B.Y.; methodology, H.X.; software, M.S.; validation, H.X., M.S.; formal analysis, H.X.; investigation, H.X.; resources, H.X., M.S.; data curation, H.X.; writing—original draft preparation, H.X.; writing—review and editing, B.Y.; visualization, H.X.; supervision, B.Y.; project administration, H.X.; funding acquisition, B.Y. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the National Science and Technology Major Project (2017-V-0012-0064).

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

*Article*


*1.2.6. Air distribution manifolds*

The supply pipe work inter-connecting the five pressure vessels is of 6 inch inside pipe di‐ ameter pipe. There is a 6 inch branch to the supersonic tunnel. Gas dynamics rigs in the Aer‐ odynamics Laboratory are supplied from two 4 inch pipe manifolds, one wall mounted and the other suspended from the ceiling. A four inch line, reducing to 3 inch, supply air to the Hydraulics Laboratory.

Maximum intermittent flow rates are about 10 lb/sec through the 6 inch branch supplying the supersonic tunnel and 4 lb/sec through the 4 inch manifolds. At these flow rates, the pressure drop between reservoirs and manifold outlets does not exceed about 3 psi. The hy‐ draulics Laboratory supply system permits an intermittent flow rate of about 5 psi.

A 2 inch dump line is provided, together with a control valve and attenuating duct silencer to empty the storage vessel contents or to permit an adjustable air bled for stabilization purposes.


Radio Frequency Identification from System to Applications

*Edited by Mamun Bin Ibne Reaz*

**RADIO FREQUENCY**

Edited by **Mamun Bin Ibne Reaz**

**SYSTEM TO**

**APPLICATIONS**

**IDENTIFICATION FROM**


**1. Introduction**

The current level of global climate change has been unprecedented in the past decades or even nearly a thousand years. Almost the entire world is experiencing a warming process, which is mainly characterized by rising temperatures, rising sea levels, retreating glaciers, and so on [1,2]. Global warming has exacerbated the global water cycle over the past century, causing a significant increase in the number of extreme weather events, such as storms, heat waves, floods, and droughts [3–9]. Droughts are one of the most threatening natural disasters in the world. They are caused by the

below-average level of precipitation over a long period of time and are generally characterized by their high frequency, long duration, and wide range [10,11]. In the context of global warming, the frequency and intensity of droughts have increased significantly [12,13], which has seriously a ffected agricultural production, water resources, and ecosystems, and have led to economic losses, famines, epidemics, and desertification [14–21]. Because of the complexity of drought variability, it is challenging to objectively quantify the intensity, duration, and spatial extent of droughts [22–24]. Thus, numerous studies have attempted to improve drought detection and monitoring; a few objective indices have been developed on the basis of readily available climate data, such as the Palmer Drought Severity Index (PDSI) [25], Standard Precipitation Index (SPI) [26], and Standard Precipitation Evapotranspiration Index (SPEI) [27]. Among these indices, the PDSI, which is based on the supply and demand in the water balance, is one of the most widely used drought indices in the world. However, the PDSI has several deficiencies, including the strong influence of the calibration period, the limitation in spatial comparability, and the subjectivity in relating drought conditions to the index values [28]. The SPI can e ffectively represent the multiscale characteristics of droughts; however only the precipitation variability is considered in the SPI calculation, and the role of temperature is ignored. The e ffect of temperature is evident in initiating droughts, although droughts are primarily caused by a below-average level of precipitation [29]. Therefore, the SPEI was developed by Vicente-Serrano et al. [27], which not only considered the effects of temperature on drought severity but also considered the multiscale characteristics that were incorporated in the SPI. Since it was proposed in 2010, the SPEI has been widely used to monitor and assess the dry-wet conditions around the world [30–33].

Vegetation, linking the atmosphere, hydrosphere, and biosphere [34], is an important component in the terrestrial ecological system and has an obvious relationship with climate change through the physiological responses of plants, such as plant photosynthesis, respiration, and evapotranspiration. The dynamic changes of vegetation play a predictive role in regional climate change [35,36]. Climate change can also a ffect the spatial-temporal pattern of vegetation. Drought is one of the most frequent natural disasters and the response of vegetation to drought is a considerable scientific problem [37]. In general, an increased frequency of extreme drought was associated with decreased vegetation growth [19,38,39]. For example, Symeonakis et al. (2004) pointed out that drought was the main factor resulting in vegetation and soil degradation in sub-Saharan Africa [39]. Ahmadi et al. (2019) indicated that drought could a ffect the e fficiency of water use in the ecosystem, subsequently disturbing the composition and functionality of terrestrial ecosystems [19]. In Northern China, drought-induced water stress caused a reduction in the terrestrial gross primary production [40]. Studies on the Qinghai-Tibet Plateau (QTP) and the Loess Plateau have revealed that there is a remarkable correlation between vegetation cover reduction and climate change [41,42]. Nevertheless, the magnitude of the response of vegetation to dry-wet conditions remains uncertain due to the complexity of the dry-wet transition and intrinsic drought sensitivity among vegetation types [21,43]. The intensity, duration, and timing of drought partly determine the e ffect of drought on vegetation productivity, where moderate drought with higher temperatures increases the net primary production (NPP), while severe drought causes a delayed response of NPP to precipitation [22,44]. As a satisfactory indicator of vegetation activities, the Normalized Di fference Vegetation Index (NDVI), which was used in this study, has been widely and successfully used to detect the vegetation variations [36–44].

The Yarlung Zangbo River (YZR) basin located in the southeast QTP, is the most important river to understand the water cycle in the QTP because it is not only the largest river system in QTP with the largest mean annual flow (56% flows from the QTP), but also an important moisture transportation channel from the Indian Ocean to the inner region of QTP. Owing to great spatial heterogeneities of climatic conditions and enormous biological diversity, the YZR basin has always been the crucial area of global diversity and ecological protection [45]. A considerable amount of research on the impact of climate change on vegetation variations in the YZR basin has been conducted, revealing that vegetation and precipitation are positively correlated in the whole basin and that the vegetation cover change is restrained by the dry-wet regime, terrain, and other factors simultaneously [46–48]. Song et al. (2011) and Li et al. (2015) have demonstrated that the warming rate of the YZR basin was significantly higher than that of the global average, and the duration and magnitude of the drought have gradually aggravated [49,50]. The water resources problem would become more severe under the impact of precipitation and temperature due to the significant warming of the YZR basin in the future [51]. The rising temperature drives a basin-wide vegetation cover improvement, however it seems that decreasing precipitation does not inverse the overall vegetation greening trend [51]. Similar results occurred in Nepal, bounded by the Tibetan highland and the Himalaya, indicating that the correlation between NDVI and temperature was significantly positive, while NDVI exhibited a negative relationship with precipitation [52]. To sum up, the dry-wet transition is of great importance for the vegetation dynamics, whereas the response mechanism of vegetation to dynamic variations is still unclear. In addition, the special terrain, vulnerable ecological environment, and sensitivity to climate change make the YZR basin one of the hotspot regions for the studies of water-ecology-environment sustainable development under global warming. Therefore, it is crucial to investigate the transition characteristics of the dry-wet regime in the YZR basin and quantify its e ffect on the vegetation dynamic variations, which could provide a scientific reference for the sustainable development of the environment and ecosystem in the Qinghai-Tibet Plateau. The primary objectives of this study are: (1) to detect the changes of the dry-wet regime in the YZR basin; (2) to quantify the spatial-temporal variations of vegetation from long-term satellite-based NDVI data; and (3) to investigate the dynamic responses of vegetation to possible drivers of the dry-wet transition in the YZR basin.

#### **2. Materials and Methods**

## *2.1. Study Area*

The Yarlung Zangbo River originates from the Gyama Langdzom Glacier is one of the highest rivers around the world. It is mainly composed of five tributaries, namely the Duoxiong Zangbo River, the Nianchu River, the Lhasa River, the Nyang River, and the Parlung Tsangpo River. The YZR basin, with a latitude of 28◦00 'N–31◦16' N and longitude of 82◦00 'E–97◦07' E, is located in the southeast QTP (Figure 1), with an area of about 24 km2 and average altitude of more than 4000 m (ranging from 132 m to 7258 m). The climate of the YZR basin is characterized as cold plateau mountain climate with intense solar radiation and low air temperatures. The amount of precipitation gradually increases from northwest to southeast in the basin, which is mainly a ffected by the warm and humid airflow from the Bay of Bengal and the Indian Ocean. The mean annual precipitation in the basin is 300–500 mm, and the rate increases with elevation by 10–30 mm/100 m [52]. The total rainfall from June to September accounts for 60%–90% of the mean annual precipitation in the whole basin, indicating the precipitation has an uneven distribution within a year. Due to the complex topographical features and high altitude, the vegetation cover within the area exhibits distinct vertical zonality and varies from mountain forest, mountain broad-leaved forest, mountain coniferous forest, and subalpine shrub meadow to alpine meadow along with rising elevation [53].

**Figure 1.** Location of the Yarlung Zangbo River basin and distribution of the hydro-meteorological stations.