Journal of Mechanical Science and Technology 27 (3) (2013) 713~720 www.springerlink.com/content/1738-494x

DOI 10.1007/s12206-013-0120-0

Centrifugal compressor shape modification using a proposed inverse design method† Mahdi Nili-Ahmadabadi1 and Farzad Poursadegh2,* 1

Department of Mechanical Engineering, Isfahan University of Technology, Isfahan, Iran 2 School of Mechanical Engineering, Sharif University of Technology, Tehran, Iran

(Manuscript Received December 8, 2011; Revised July 29, 2012; Accepted September 12, 2012) ----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

Abstract This paper is concerned with a quasi-3D design method for the radial and axial diffusers of a centrifugal compressor on the meridional plane. The method integrates a novel inverse design algorithm, called ball-spine algorithm (BSA), and a quasi-3D analysis code. The Euler equation is solved on the meridional plane for a numerical domain, of which unknown boundaries (hub and shroud) are iteratively modified under the BSA until a prescribed pressure distribution is reached. In BSA, unknown walls are composed of a set of virtual balls that move freely along specified directions called spines. The difference between target and current pressure distributions causes the flexible boundary to deform at each modification step. In validating the quasi-3D analysis code, a full 3D Navier-Stokes code is used to analyze the existing and designed compressors numerically. Comparison of the quasi-3D analysis results with full 3D analysis results shows viable agreement. The 3D numerical analysis of the current compressor shows a huge total pressure loss on the 90° bend between the radial and axial diffusers. Geometric modification of the meridional plane causes the efficiency to improve by about 10%. Keywords: Inverse design; Centrifugal compressor; Meridional plane; Quasi-3D; 90° bend ----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

1. Introduction Inverse design in fluid flow problems includes finding a shape associated with a prescribed profile of wall pressure or velocity. Design issues of the border can be solved in two ways: non-iterative (coupled or direct) and iterative (noncoupled). In the direct solution method, a form of problem formulation can be used whereby coordination of the border is a dependent variable (implicitly or explicitly) in the governing differential equation. In the iterative method, flow and geometrical variables in the solution process are independent of one another. These methods commence with a starting initial presumption, then a corresponding analytical solution is solved and again the process is performed until convergence (which is generally the difference between the boundary pressure distribution and target pressure distribution) is achieved. Although iterative methods are general and powerful, they are often computationally costly and mathematically complex. These methods utilize analysis methods for the flow field solution as a black-box [1]. One of the iterative methods is based on the residualcorrection approach. In this method, the key problem is in relating the calculated differences between the current and target pressure profile (residual) to the required changes in the *

Corresponding author. Tel.: +98 937 110 1364, Fax.: +98 311 3912628 E-mail address: [email protected] † Recommended by Associate Editor Byeong Rog Shin © KSME & Springer 2013

geometry. Obviously, the approach for developing a residualcorrected method is to identify an optimum state between the computational effort (for determining the required geometry correction) and the number of iterations needed to obtain a converged solution. This geometry correction may be estimated by means of a simple correction rule that utilizes relations between geometry changes and pressure differences known from linearized flow theory [1]. Residual-correction decoupled solution methods attempt to utilize the analysis method as a black-box. Barger et al. [2] presented a streamlined curvature method that considers the possibility of relating a local change in surface curvature to a change in local velocity. A large number of methods have been developed following this concept [3-12]. Recently, Nili proposed flexible string algorithm (FSA), a two-dimensional inverse design algorithm for internal flow regimes. He developed his method for ideal compressible subsonic [13], supersonic [14], and viscous incompressible flows [15]. The role of quasi-3D flow analysis in the aerodynamic design process has been described in Ref. [16], relative to the detailed design of impellers. Computational speed is the key feature that keeps quasi-3D Euler codes firmly entrenched in virtually all design systems, even though more exact viscous computational fluid dynamics (CFD) codes are available. Quasi-3D flow analysis can generate a solution in a matter of seconds, even on a personal computer of modest capability. Currently, this is the only type of analysis method that can

714

M. N.-Ahmadabadi and F. Poursadegh / Journal of Mechanical Science and Technology 27 (3) (2013) 713~720

provide the type of immediate internal flow analysis a designer needs for efficient, detailed component design activity. For centrifugal compressor impellers, comparing the quasi3D flow analysis with experimental results shows close agreement, except near the trailing edge and leading edge, where the experimental data show higher velocities than those predicted by the analysis [16]. Ref. [17] gives the numerical techniques for quasi-3D flow analysis that are used to find the flow distribution in the meridional plane of a centrifugal compressor. This method solves velocity gradient equations with the assumption of a hub-toshroud mean stream surface. Instead of normal lines, a set of arbitrary straight lines from hub to shroud is used, called quasi-orthogonal lines, and these lines remain fixed regardless of any streamline changes. This work determines the velocities in the meridional plane of a backward-swept impeller, a radial impeller, and a vaned diffuser, as well as approximate blade surface velocities. The presence of a high loss region, as pointed out by Zangeneh et al. [18], is similar to the “jet/wake” phenomenon in centrifugal compressors. They indicated that the high loss region is generated by secondary flows that move low momentum fluids under the action of gradients under reduced static pressure toward the location of minimum reduced pressure at the shroud-suction corner near the trailing edge. In this research, an inverse design algorithm, called ballspine algorithm (BSA), is used to modify the current geometry based on modification of pressure distribution. The method is a scheme for inverse design rather than optimization. The BSA is a developed version of FSA [13-15]. The FSA had some restrictions for the design of axe-symmetric ducts. Consequently, this work developed a comprehensive one to cover diverse geometries. The BSA is of physical concept, whereas optimization approaches are based on mathematical concepts. In other words, the new feature of the method consists of considering the duct wall as a set of moving balls with mass. The BSA is categorized as a residual-correction method. Unlike other residualcorrection methods with mathematical basis using flow equation for inverse design problems, the BSA turns the inverse design problem into a fluid-solid interaction scheme with physical sense that uses pressure to deform the flexible wall. In this work, the BSA is applied for the aerodynamic design of a 90° bend between the radial and axial diffusers of a centrifugal compressor. For this purpose, a quasi-3D code that solves the Euler equation on the meridional plane of a centrifugal compressor is incorporated into the BSA in order to design the meridional plane. Decreasing the high loss region on the 90° bend between the radial and axial diffusers can be done by modifying the hub-shroud profiles and will improve compressor efficiency. Once the pressure distribution along the hub and shroud wall is modified, the walls are changed to satisfy the modified pressure distributions. The compressor with modified hub-shroud profile is numerically simulated and analyzed by a fully 3D Navier-Stokes code, and its per-

Fig. 1. Schematic of 90° bend duct with balls and spines.

formance obtained from numerical analysis is compared with the numerical performance of the existing compressor.

2. Fundamentals of ball-spine inverse design Fig. 1 illustrates a 2D flexible 90° bend. The bend duct walls are composed of imaginary balls that can be freely moved along specified directions. Passing fluid flow inside the bend duct applies a pressure profile (current) on the wet side. When a target pressure profile is applied to the outer side of the duct, the flexible walls begin to be deformed in a way that the inner side pressure profile (current) equals the outer pressure profile (target). In other words, the force from the difference between current and target pressure profiles anywhere along the wall leads the corresponding ball at any point to be displaced along the specified direction. As the target shape is achieved, the difference between two pressure profiles reaches zero, and the balls will stop moving. If the balls move in the direction of the applied force, adjacent balls may collide or move away from each other. Each ball in the process of deformation can move along a specified direction, called spine, to avoid this problem. Fig. 1 shows the spines for a 90° bend duct. The first points of each wall and the outlet plane must be fixed to be converged in the design procedure. The displacement of balls in each shape modification step is calculated from Eq. (1). Applying the calculated displacement to each ball obtains it’s a new position. ∆si = C ∆Pi .

(1)

In this equation, (∆s) is the displacement of each point along its spine, (∆P) is the difference between current and desired pressure profiles at each point, and (C) is an adjusting parameter for the convergence rate of BSA method. The lesser the amount of (C), the slower the convergence rate will be. Indeed, if the amount of (C) exceeds the limit, the design algorithm will diverge. In clarifying how kinematic relation of the abovementioned

M. N.-Ahmadabadi and F. Poursadegh / Journal of Mechanical Science and Technology 27 (3) (2013) 713~720

715

Residual

Fig. 2. Free body diagram of a ball.

flexible wall is driven, a uniform mass distribution along the wall is assumed. A free body diagram of a virtual ball on the wall is shown in Fig. 2. The relations can be derived as follows: FS = ∆P ⋅ A ⋅ cosθ = maS ∆S =

⇒

aS =

∆P ⋅ A ⋅ cos θ m

1 2 aS ( ∆t ) . 2

(2) (3)

If (ρ) is defined as wall surface density, the following relation is obtained:

( ∆t ) ∆P cosθ . 1 ∆P ⋅ A ⋅ cosθ 2 ( ∆t ) = 2 ρA 2ρ 2

∆S =

(4)

In the above, (∆t) is an interval for ball movement at each shape modification step. Based on the Newton’s second law, when a force is applied to a mass, displacing that mass takes some time. After determining the forces, which come from the computed pressure distributions, a virtual time step, ∆t, is needed to displace the flexible boundaries. Thus, the new position of each ball for a duct, such as the one shown in Fig. 1, is obtained from the following relation: (t + ∆t )

Si

(t )

= Si

( ∆t ) + 2ρ

2

(t )

∆Pi cosθi = Si + C ∆Pi .

(5)

Despite mathematical-based methods that require arbitrary choice of parameters, in the BSA method, the choice of ball mass (m) and interval for ball movement control the convergence and stability of the numerical solution. Similar to other regular design procedures, residual for the proposed design process is defined as follows: N

∑  P − ( P i

residual =

T arg et

i =1

N

(

)

) 

 PT arg et  ∑ i  i =1

i

.

(6)

Fig. 3. Inverse design flowchart.

In Eq. (6), (Ptarget) and (P) are the desired and current pressure profiles at each point, respectively. 10-2 is the convergence criteria whereby the shape modification process should be terminated and the obtained shape is the desired one. Similar to the case in other iterative inverse design methods, the flow field should be analyzed at each shape modification step. In this paper, the numerical techniques presented in Ref. [17] are used for quasi-3D flow analysis through the meridional plane of the radial diffuser, 90° bend, and deswirl vanes. The spines used for inverse design are the same quasiorthogonal lines used for quasi-3D analysis. Fig. 3 shows how BSA is typically incorporated into existing flow solution procedures. The computed pressure surfaces are normally obtained from partially converged numerical solutions of quasi3D flow equations.

3. Quasi-3D analysis procedure The meridional plane of the centrifugal compressor is analyzed by quasi-3D code. The input quantities for quasi-3D analysis include rotational speed, mass flux, number of blades, specific-heat ratio, inlet total temperature and density, gas constant, hub-to-shroud profile, mean blade shape, and a normal thickness distribution of blades. The detailed geometrical information and operational conditions of the compressor are described in Table 1.

4. Compressor numerical investigations The numerical simulation inside the compressor is accomplished to achieve precise insight into its flow field. The grid generation of the compressor, as the most time-consuming step of the numerical simulation, is done with structural elements in all domains. Then, the generated grids near the walls

Table 1. Centrifugal compressor parameters. Parameter (unit)

Value

Operational condition

Design point

Number of impeller blades

25

Impeller eye diameter (mm)

250

Impeller outlet diameter (mm)

340

Impeller outlet width (mm)

28

Radial diffuser length (mm)

47.5

Axial diffuser length (mm)

39

Number of radial diffuser vanes

24

Number of axial diffuser vanes

60

Specific-heat ratio

1.4

Gas constant (J/(kg·K))

287

Inlet flow angle (°)

0.0

Inlet total temperature (°K)

322.2

Inlet total density (kg/m3)

1.088

Design point mass flow rate (kg/s)

5.25

Design point Rotational speed (rpm)

16222

and blades are refined to capture the steep gradients. The Reynolds-averaged Navier-Stokes equations describe the conservation equations by finite volume method. The equations are discretized via a coupled implicit method. Moreover, shear stress transport model is used for turbulence modeling. The mixing plane interface model is applied for dividing the computational domain into stationary and moving zones. It utilizes relative motion among various zones to send calculated values among the zones. Furthermore, coriolis and centrifugal accelerations are added to the momentum equations in rotating zones. By using the mixing plane model and periodic boundary condition, modeling only one pitch of the impeller and diffusers is feasible. Mass flow and total temperature at the inlet and average static pressure at the outlet are used for this simulation. In studying the grid independency, the compressor pressure ratio and efficiency at the design point are considered as the parameters for examining three grid configurations (Fig. 4). The selected convergence criterion is the attainment of a constant value for drag, lift, and moment coefficients of walls. In Fig. 4, the calculated compressor pressure ratio and efficiency reach an asymptotic value as the number of elements increases. According to this figure, grid B (64,8538 elements) is considered to be sufficiently reliable to ensure mesh independency. The y+ value changes from 30 to 130 for grid B. Fig. 5 shows the generated grid for the impeller, radial, and axial diffusers. Three-D numerical analysis is performed at different rotational speeds. For each rotational speed, the values of mass flow rate and stagnation temperature at the inlet and static pressure at the outlet are applied as boundary conditions. Fig. 6 shows the contour of static pressure in the mid sec-

Total pressure ratio

M. N.-Ahmadabadi and F. Poursadegh / Journal of Mechanical Science and Technology 27 (3) (2013) 713~720 1.991 1.989 Grid A Grid B 648538

284914

Grid c 1477128

1.987 1.985 2

4

6

8

10

12

14

16

Number of elements*1e5 86.2 Efficiency

716

Grid A 284914

Grid B 648538

86.15

Grid c 1477128

86.1 2

4

6 8 10 12 Number of elements*1e5

14

16

Fig. 4. Effect of grid size on pressure ratio and efficiency of the compressor.

Fig. 5. Compressor geometry with its grid.

Fig. 6. Static pressure contour on the mid section.

tion between the hub and the shroud. The static pressure increases along the impeller and the diffuser, as expected. As shown in Fig. 7, the total pressure increases through the impeller due to the energy transfer to the fluid. In the radial diffuser, except in a small region at the tail end of the diffuser vanes, the stagnation pressure approximately remains constant. However, in the axial diffuser, because of the acute 90° bend, the stagnation pressure drastically decreases, which implies drastic loss in this part. Certainly, the stagnation pressure loss causes the compressor efficiency to decrease.

M. N.-Ahmadabadi and F. Poursadegh / Journal of Mechanical Science and Technology 27 (3) (2013) 713~720

717

Fig. 10. Geometry of current and modified diffusers. Fig. 7. Stagnation pressure on the mid section.

5. Modifications of 90° bend on the meridional plane

Fig. 8. Mach number contour on the mid section.

Fig. 9. Numerical efficiency calculated at the radial and axial diffuser outlets.

Fig. 8 shows the contour of absolute Mach number in the mid section between the hub and the shroud. The Mach number increases through the impeller due to the increase in the linear velocity of the impeller, approaching 1 at the impeller outlet, which represents choking of the compressor at higher rotational speeds. In the radial diffuser, the Mach number decreases due to diffusion. The high loss of momentum on the lower wall of the axial diffuser causes the mach number to be reduced and even reverse flows to be generated. In Fig. 9, the numerical efficiency calculated at the radial diffuser outlet is compared with that at the axial diffuser outlet at different rotational speeds. These differ by about 20%. This difference is related to the high loss region at the 90° bend.

As mentioned before, the acute 90° bend passage between the radial and axial diffusers hampers efficiency, which can be improved by geometric modification at the meridional plane. The meridional plane of two diffusers is modified by using the BSA design algorithm and modifying the current pressure distribution along the hub and shroud of the diffusers. In modifying the current pressure distribution on the hub and shroud of the 90° bend, the steep adverse pressure gradients are smoothed to cut down any probable separation due to these adverse gradients. Considering that the 90° bend causes dramatically adverse pressure gradients to be applied on the inner wall (hub), and given that the area section of flow increases from the radial diffuser outlet to the axial diffuser inlet through a very short length, the applied modification cannot remove the adverse pressure gradient thoroughly. The axial length of the entire compressor has a noticeable constraint. Fig. 10 shows the modified meridional plane of the diffusers. In this design, the additional adverse pressure gradient caused by the sudden deflection of the bend is removed. However, the inlet and outlet pressures are not changed. In other words, the area ratios of the diffusers are fixed. The 3D numerical simulation of the compressor with the modified 90° bend and with the same current impeller shows that the total efficiency is improved by about 10%. In Figs. 11 and 12, the flow field contours of the modified geometry are compared with those of the current geometry on the meridional plane; the efficiency is shown to improve. Although the entire separation is removed, it becomes negligible. As mentioned before, this phenomenon is due to the very short length of the 90° bend with high area ratio. The boundary conditions are based on mass flow inlet and average static pressure outlet. Therefore, the quantity of stagnation pressure on either side cannot be controlled. Indeed, the modified case meets more increase in total pressure from inlet to outlet, which is clearly seen in the modified case. Figs. 13-15 and 16-18, respectively, show the stagnation pressure contour for the current and modified axial diffusers, on the blade-to-blade plane, at three spans (hub, mid, and shroud). Comparison between two sets of three figures indicates the improvement of energy losses. In this work, the

718

M. N.-Ahmadabadi and F. Poursadegh / Journal of Mechanical Science and Technology 27 (3) (2013) 713~720

Fig. 11. Stagnation pressure contour on current diffusers.

Fig. 15. Stagnation pressure contour on the blade-to-blade plane of the current axial diffuser near the shroud.

Fig. 12. Stagnation pressure contour on modified diffusers.

Fig. 16. Stagnation pressure contour on the blade-to-blade plane of the modified axial diffuser near the hub.

Fig. 13. Stagnation pressure contour on the blade-to-blade plane of the current axial diffuser near the hub.

Fig. 17. Stagnation pressure contour on the blade-to-blade plane of the modified axial diffuser at mid section.

Fig. 14. Stagnation pressure contour on the blade-to-blade plane of the current axial diffuser at mid section.

Fig. 18. Stagnation pressure contour on the blade-to-blade plane of the modified axial diffuser near the shroud.

M. N.-Ahmadabadi and F. Poursadegh / Journal of Mechanical Science and Technology 27 (3) (2013) 713~720

719

Fig. 21. Efficiency improvement versus rotational speed.

Fig. 19. Stagnation pressure distribution at the current compressor outlet.

Fig. 22. Pressure ratio improvement versus rotational speed.

cial flow analysis software as a black-box. In this research, the BSA design procedure is incorporated into a quasi-3D analysis code for designing the 90° bend between the radial diffuser and deswirl vanes. Modifying the 90° bend increases the efficiency by 10%, which is a considerable improvement. Fig. 20. Stagnation pressure distribution at the modified compressor outlet.

geometry was modified only on the meridional plane not on the blade-to-blade one. Figs. 11-18 show that although decreased, the losses still exist. These losses can be removed by modifying the axial diffuser blade profile, which increases the efficiency. However, blade profile modification is not related to this research. Figs. 19 and 20 compare the stagnation pressure distribution at the compressor outlet, depicting total pressure improvement for the modified geometry. Also, the increment of pressure ratio and efficiency of the modified geometry at different rotational speeds is shown in Figs. 21 and 22. Besides, according to the fact that the modifications have been made on 90 degree bend, the results show no efficiency difference between the impellers in both cases.

6. Conclusions The BSA turns the inverse design problem into a fluid-solid interaction scheme that has a physical base. The BSA is a quick converging method and can efficiently utilize commer-

Nomenclature-----------------------------------------------------------------------a A BSA F FSA m n P P0 r TPD ∆ ∆t ∆Pi θ ρ ηtt

: Linear acceleration (m.s-2) : Local area of duct wall (m2) : Ball spine algorithm : Force vector (N) : Flexible string algorithm : Ball mass (kg) : Number of virtual balls on the wall : Static pressure (Pa) : Stagnation pressure (Pa) : Radial coordinate, radius (m) : Target pressure distribution : Difference : Time step(s) : Difference between TPD and CPD at each link : Angle between spine and force direction (deg) : Density of fluid (kg·m-3), surface density (kg·m-2) : Total to total efficiency

Subscripts i

: Balls index

720

S

M. N.-Ahmadabadi and F. Poursadegh / Journal of Mechanical Science and Technology 27 (3) (2013) 713~720

: Spine direction

Superscripts I.G. : Initial guess t+∆ t : Related to updated geometry t : Related to current geometry

References [1] F. Ghadak, A direct design method based on the laplace and euler equations with application to internal subsonic and supersonic flows, PhD Thesis, Sharif University of Technology, Aero Space Department, Iran (2005). [2] R. L. Barger and C. W. Brooks, A streamline curvature method for design of supercritical and subcritical airfoils, NASA TN D-7770 (1974). [3] P. Garabedian and G. McFadden, Design of supercritical swept wings, AIAA Journal, 30 (3) (1982) 444-446. [4] J. Malone, J. Vadyak and L. N. Sankar, Inverse aerodynamic design method for aircraft components, J. Aircraft, 24 (1) (1986) 8-9. [5] J. Malone, J. Vadyak and L. N. Sankar, A technique for the inverse aerodynamic design of nacelles and wing configurations, AIAA Paper 85-4096 (1985) . [6] R. L. Campbell and L. A. Smith, A hybrid algorithm for transonic airfoil and wing design, AIAA Paper 87-2552 (1987). [7] R. A. Bell and R. D. Cedar, An inverse method for the aerodynamic design of 3-D aircraft engine nacelles, Third International Conference on Inverse Design Concepts and Optimization in Engineering Sciences, ICIDES-III, G.S. Dulikravich, Washington, D.C., 23-25 (10) (1991) 405-417. [8] J. B. Malone, J. C. Narramore and L. N. Sankar, An efficient airfoil design method using the Navier-Stokes equations, AGARD (1989) Paper 5. [9] J. B. Malone, J. C. Narramore and L. N. Sankar, Airfoil design method using the Navier-Stokes equations, J. Aircraft, 28 (3) (1991) 216-224. [10] S. Takanashi, Iterative three-dimensional transonic wing design using integral equations, J. Aircraft, 22 (1985) 655660. [11] N. Hirose, S. Takanashi and N. Kawai, Transonic airfoil design procedure utilizing a Navier-Stokes analysis code, AIAA Journal, 25 (3) (1987) 353-359. [12] G. S. Dulikravich and D. P. Baker, Aerodynamic shape inverse design using a fourier series method, AIAA Paper

99-0185 (1999). [13] M. Nili-Ahmadabadi, M. Durali, A. Hajilouy and F. Ghadak, Inverse design of 2D subsonic ducts using flexible string algorithm, Journal of Inverse Problems in Science and Engineering, 17 (8) (2009) 1037-1057. [14] M. Nili-Ahmadabadi, A. Hajilouy, M. Durali and F. Ghadak, Duct design in subsonic & supersonic flow regimes with & without shock using flexible string algorithm, Proceedings of ASME Turbo Expo 2009, Florida, USA, GT2009-59744 (2009). [15] M. Nili-Ahmadabadi, A. Hajilouy, F. Ghadak and M. Durali, A novel 2-D incompressible viscous inverse design method for internal flows using flexible string algorithm, Journal of Fluids Engineering, ASME, 132, 031401-1-9 (2010). [16] R. H. Aungier, Centrifugal compressor, a strategy for aerodynamic design and analysis, ASME Press, New York, USA (2000). [17] M. R. Vanco, Fortran program for calculating velocities in the meridional plane of a turbomachine, NASA TN D-6701 (1972). [18] M. Zangeneh, W. N. Dawes and W. R. Hawthorne, Three dimensional flow in radial-inflow turbine, ASME paper No. 88-GT-103 (1988).

Mahdi Nili-Ahmadabadi is Assistant Professor and faculty member of the Department of Mechanical Engineering at Isfahan University of Technology. He received his Master and Ph.D. degrees from Sharif University of Technology in 2005 and 2010, respectively. His major research interests are in the fields of inverse design, turbomachinery, experimental aerodynamics, and PIV measurement. Farzad Poursadegh, the corresponding author of the paper, is Master of Science graduate from Sharif University of Technology. His major research interests are in the fields of turbomachinery and aerodynamics. He has recently focused on the development of inverse design methods in turbomachinery, especially in compressors. Moreover, he has been the chief specialist in thermo-fluids of oil turbo-compressor equipment since 2010.