Journal of Mechanical Science and Technology 28 (3) (2014) 895~905 www.springerlink.com/content/1738-494x

DOI 10.1007/s12206-013-1157-9

Analysis and measurement of the impact of diffuser width on rotating stall in centrifugal compressors† Yong-Sang Yoon* and Seung Jin Song School of Mechanical and Aerospace Engineering, Seoul National University, Seoul, 151-742, Korea (Manuscript Received January 17, 2013; Revised September 8, 2013; Accepted October 6, 2013) ----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

Abstract In compression systems, instability has long been an important issue. However, compared to axial machines, relatively little work has been done on the stability of centrifugal machines. Especially, many analytical models of stabilities have been developed to predict and control rotating stall, using compressor characteristic. However, stability models for centrifugal compressors are not scarce. Much research on compressor stability has focused on stalling flow coefficient and rotating stall phenomenon at the stalling flow coefficient. Given this situation, this paper presents a stability analysis of centrifugal compressors to predict rotating stall inception as well as the speed and number of cells. This analysis involves the use of compressor geometries, a steady compressor characteristic, and threedimensional flow analysis in the diffuser. The flow field perturbations at the axial inlet duct, impeller, and radial exit duct are determined via an eigenvalue analysis. The predictions are validated against experimental results from compressors with three different diffuser widths. The model accurately predicts the rotating stall inception flow coefficient. As the compressor characteristic becomes less steep with increasing diffuser width, the stalling flow coefficient increases. Also, experiment validates the model prediction that, depending on the dominant mode of flow perturbation, the number of rotating stall cells can be changed from three to two cells in the tested configurations. Furthermore, the cell speed increases as the flow coefficient decreases for a given number of stall cells. However, when the stall cell number is reduced, the cell speed decreases. Keywords: Centrifugal compressor; Diffuser; Instability modeling; Rotating stall; Stall cell number; Stall cell speed ----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------

1. Introduction The operating range of compression systems is limited at low mass flow rates by the onset of fluid dynamic instability. As the flow rate decreases, the pressure rise across the compressor increases monotonously until flow becomes unstable. At this point, the steady, axisymmetric flow through the compression system becomes unsteady and non-uniform. Also, the annulus averaged pressure rise and mass flow through the compression system are reduced. There are two types of dynamic instability - surge and rotating stall. Surge is a onedimensional instability characterized by oscillations in the compressor annulus averaged mass flow and pressure rise that extend through the entire compression system. Rotating stall is a two- or three-dimensional disturbance in which regions of low, or reversed, mass flow (termed stall cells) rotate about the annulus of the compressor. To understand these instabilities, analytical models have been developed. Greitzer [1, 2] has identified the B-parameter *

Corresponding author. Tel.: +82 70 7147 4217, Fax.: +82 31 8018 3732 E-mail address: [email protected] † Recommended by Associate Editor Donghyun You © KSME & Springer 2014

which describes surge dynamics. Moore and Greitzer [3] suggested a two-dimensional model in time variation to explain rotating stall in axial compressors. Based on the two models, much research on instability and stability management has been done for axial compressors. Hynes and Greitzer [4] presented an axial compressor stability model to assess the effect of inlet flow distortion. Graf et al. [5] analytically and experimentally investigated the impact of nonuniform geometry on axial compressor stability. Gordon [6] extended Graf’s model to account for both full-span and part-span rotating stall modes of axial compressors, using a three-dimensional compressor model. According to Camp and Day [7], rotating stall can be initiated with two types of instability in axial compressors. One is a short length scale instability, ‘spike,’ and the other is a longlength scale ‘modal’ type instability. The rotating stall inception types in axial compressors are influenced by the rotor incidence angle and slope of compressor characteristic. Thus, as flow coefficient decreases, if the compressor characteristic curve has a peak before critical rotor incidence angle is reached, modal type stall inception occurs. On the other hand, if the peak occurs after the critical rotor incidence angle, spike type stall is initiated.

896

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

However, these models have all been focused on axial compressor instabilities. For centrifugal compression system instability, Spakovszky [8, 9] renewed both axial and radial compressor dynamic system modeling, introducing ‘transmission matrix’ which can be linked to any other system component. This approach has an advantage that each component has its own transmission matrix and these matrices can be coupled to other connected components. This enables modeling complex systems such as centrifugal compressors with efficiency and flexibility. Instead of compressor characteristic, this model requires loss coefficient vs. flow angle as input parameters. These stability models enable control of rotating instability to enhance the compressor operating range. Based on this kind of two-dimensional rotating stall model, dynamic control of rotating stall has been implemented by Paduano [10], using the control of IGV stagger angle. Gysling [11] has also managed the rotating stall at the initiation point by air injection. However, such research has only focused on stall inception, so it is difficult to understand how the rotating stall develops as flow coefficient decreases. To understand the rotating stall phenomena, some research has been conducted. Kaemmer and Rautenberg [12] have found rotating stall with a cell number varying between 1 and 2 in a centrifugal compressor. Haupt et al. [13] experimentally found several stall patterns of two and three cells at one operating point in a centrifugal compressor. However, the mechanism of how rotating stall is developed with specified cell number is not clear. Despite such preceding work, research about rotating instability in a centrifugal compressor is still relatively scarce. Therefore, this paper aims to analyze and measure rotating stall in “modal” centrifugal compressors. Specific objectives to be addressed are as follows. (1) To develop a centrifugal compressor rotating stall model with compressor characteristic as input, which is a relatively easy parameter to measure. (2) To predict flow coefficient where rotating stall begins. (3) To predict and understand the relationship between stall cell number, speed, and flow coefficient. Both analytical and experimental methods are used in this study.

2. Centrifugal stability modeling 2.1 Model description The new analytical model for centrifugal compressor stability adopts the transmission matrix method used by Spakovszky [8, 9]. Spakovszky has developed a novel centrifugal compressor modeling approach in which the main components (e.g., impeller, diffuser) are modeled as separate modules. Flows are considered to be two-dimensional in this model. Furthermore, common compressor parameters (number of impeller blades, slip factor, impeller blade exit angle, e.g.) are taken into account to make the model generally applicable to different compressor configurations. However, the values of some model parameters are difficult to obtain with

x=∞ Axial duct x=0 Impeller

r2

r3

Diffuser (Radial duct)

Fig. 1. Schematic of centrifugal compressors.

reasonable accuracy. In particular, the total pressure rise in the impeller is assumed to be isentropic and compensated with loss coefficient with respect to flow angle. To simplify this input parameter, in this paper, the impeller pressure rise is directly obtained from the compressor characteristic. Therefore, the model’s flow domain consists of the inlet and exit ducts connected by a compressor impeller characteristic curve. Centrifugal compressors consist of an axial duct, impeller, and radial diffuser as shown in Fig. 1. The total system model produces an eigenvalue problem by boundary conditions. The flow parameters in the compressor, such as stream function, velocity, and pressure, are assumed to be composed of mean and unsteady perturbation parts. For example, velocity is expressed as: V (t ,q ) = V + d v(t ,q ) .

(1)

All of the variables are nondimensionalized by density r , the impeller tip radius R, and the impeller inlet tip speed U t . The impeller frequency derived from these values, W = is used. Lengths are nondimensionalized as x = as v =

Ut , R

X , velocities R

P V , pressure as p = , and time as t = T W . Ut rU t 2

2.2 Model assumptions The flow in this system is assumed to be two-dimensional, inviscid, and incompressible. The impeller is modeled as semi-actuator disk to connect the two flows in the axial duct and radial ducts (Fig. 1). Since the test centrifugal compressor has a vaneless diffuser, the diffuser flow is assumed to be inviscid and two-dimensional. Also, the mean flow in the diffuser is assumed to be a free vortex flow. 2.3 Axial duct flow The flow is assumed to be two-dimensional in the axial ( x ) and circumferential ( q ) directions. The axial and circumferential velocity perturbations and pressure perturbation are obtained from the vorticity and the stream function equations. The vorticity equation is written as:

897

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

¶dw ¶dw ¶dw + vx + vq =0 ¶t ¶x ¶q

(2)

and vorticity can be described in terms of the stream function as ¶ 2dY ¶ 2dY + = -dw ¶x 2 ¶q 2

(3)

æl v ö é ù - çç n + jn q ÷÷ x v vx ø ê ú jne è x ê ú æ ö l v ê - çç n + jn q ÷÷ x ú æ ö l v v vx ø ú Tax , n (:,3) = êç n + jn q ÷ e è x . êè v x ú vx ø ê ú 0 ê ú ê ú êë úû

(8c)

2.4 Radial duct flow where vx and vq are the mean flow velocities. Assuming a periodic solution in q , the vorticity and stream function perturbations can be decomposed into spatial harmonics. ¥

dw =

åw% ( x, t ) × e

jnq

(4a)

n

n=0

In centrifugal compressors, the flow enters the impeller axially and exits radially. Thus, the flow in a diffuser is also assumed to be two-dimensional in the radial ( r ) and circumferential ( q ) direction. The procedure to obtain the transmission matrix of the flow in the radial duct is similar to that for the axial duct. The vorticity equation can be written as:

¥

dY =

å Y% ( x, t ) × e

jnq

n

.

(4b)

n=0

In addition, w% n ( x, t ) can also be decomposed into functions of x and t as w% n ( x, t ) = z%n (t ) × c ( x) .

(5)

Thus, the vorticity equation, Eq. (2), can be represented as follows: z&%n v c% ¢ + jnvq c% n =- x n = ln c% n z%n

(6)

where ln is the eigenvalue which determines stability. Next, by using the tangential momentum equation, pressure perturbation can be obtained. Thus, the flow field perturbation solution can be denoted with the following transmission matrix and is written as: é d Vx ù êd V ú = ê qú êë d P úû

ì é A%u , n ù ü ï ê ú ï l t + jnq íTax , n × ê B%u , n ú ý × e n n =1 ï êC% ú ï ë u,n û þ î

(7)

where Q and G are the mass flow rate and circulation which are defined as r vr (r ) and r vq (r ) , respectively. In addition, the vorticity and stream function give the follwing equation: ¶ 2dY 1 ¶dY 1 ¶ 2dY + + 2 = -dw . r ¶r r ¶q 2 ¶r 2

(10)

Also, from the tangential momentum equation, the pressure perturbation in diffuser can be obtained. Thus, the flow perturbation field is given by the radial transmission matrix as follows: éd Vr ù êd V ú = ê qú êë d P úû

ì é D% u , n ù ü ï ê ú ï l t + jnq íTrad , n × ê E% u , n ú ý × e n n =1 ï ê F%u , n ú ï ë ûþ î ¥

å

é jnr n -1 ê -nr n -1 Trad , n (:,1) = ê ê ê ( - jnQ + nG ) r n - 2 - jln r n ë

{

where é ù jne nx ê ú nx Tax , n (:,1) = ê -ne ú ê nx ú êë( - jln + nvq - jnvx )e úû é ù jne - nx ê ú - nx Tax , n (:, 2 ) = ê ne ú ê - nx ú êë( jln - nvq - jnvx )e úû

(9)

(11)

where Trad , n is a 3 by 3 matrix

¥

å

¶dw Q ¶dw G ¶dw + + 2 =0 r ¶r ¶t r ¶q

(8a)

(8b)

}

ù ú ú ú ú û

é jnr - n -1 ê Trad , n (:, 2) = ê nr - n -1 ê ê ( - jnQ - nG ) r - n - 2 + jln r - n ë

{

(12a)

}

ù ú ú ú ú û

é ù R ê ú jn n ê ú r ê ú ¶Rn ú Trad , n (:,3) = ê ¶r ê ú êì ú 2 ü æ ö l d R r dR Q Q G ï ï n êí +ç + + n ÷ n ýú 2 êë îï jn dr jn ø dr þïúû è r jnr

(12b)

(12c)

898

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

angle, respectively. Eqs. (15)-(17) are linearized and cast into a transmission matrix resulting in

and

Rn =

ò

l æ G ö r - ç jn log x + n x 2 ÷ 2Q ø è Q

e

r0

(r x

n - n +1

)

- r - nx n +1 dx .

(12d)

2.5 Impeller actuator A semi-actuator disk model is presented to model the flow in the impeller. P2 - Pt1 = Y ts - limp

¶vx ¶v - limp x ¶q ¶t

(

ARimp log ARimp ARimp - 1

)×s

imp

.

(14)

r 2 A2 is the impeller area-density ratio and r1 A1 is the mean passage length in the streamwise direction.

Here, ARimp = simp

Eq. (13) shows that the unsteady local pressure difference across the compressor (the exit static minus the inlet total) is given by the steady flow compressor pumping characteristic plus a correlation to account for the fact that the flow in the impeller is non-uniform in the circumferential direction. Linearizing Eq. (13) yields Eq. (15). Though first derived by Hynes and Greitzer [4], Eq. (15) has been applied to centrifugal compressor stability analysis for the first time. d Y ts ¶d vx ¶d vx . d P2 - d Pt1 = d vx - limp - limp dj ¶q ¶t

In Eq. (15),

(15)

d Y ts is the derivative of the pressure rise dj

characteristic with respect to the flow coefficient, or the slope of the impeller characteristic. Assuming a spatial harmonic solution to Eq. (15), the relation between the pressure perturbations of the impeller exit and inlet can be obtained. The radial velocity at the exit of impeller can be obtained from continuity as Vr 2 =

1 Vx1 . ARimp

(16)

The circumferential velocity is obtained from the assumption that the flow follows the blades perfectly. Blockage and deviation are not accounted for in this model. Thus, Vq 2 = U t + Vr 2 tan b 2 = Vr 2 tan a 2

ì é v%x1, n ù ü ï ê úï l í Bimp , n êv%q 1, n ú ý e n ï n =0 ê P%1, n ú ï ë ûþ î ¥

å

jnq

.

(18)

Thus, the n-th harmonic transmission matrix for a radial impeller is defined as:

(13)

where limp presents the inertia of the fluid in the impeller passage and is defined as:

limp =

éd Vr 2 ù êd V ú = ê q2ú êë d P2 úû

(17)

where b 2 and a 2 are a relative exit angle and absolute exit

Bimp , n

é ù 1 0 0ú ê ARimp ê ú ê ú tan b 2 =ê 0 0ú . ê ú ARimp ê ú êæ ¶Y ts ú ö + + jn V V 1 l l ( ) êç n imp x1 ÷ ú q1 ¶ j ø ëè û

(19)

In this modeling, three different diffuser widths are analyzed. The different diffuser width changes the ARimp value so diffuser width is reflected in ARimp . 2.6 Transmission matrix stacking As shown in Fig. 1, the perturbation relationship between the impeller inlet and exit flows connects the axial duct flow to the radial duct flow. From Eqs. (7), (11) and (18), the transmission matrices are stacked into the following equation: ì é A%u , n ù ü é D% u , n ù ï ê úï ê% ú Trad , n ( r2 ) × ê Eu , n ú = Bimp , n × íTax , n × ê B%u , n ú ý . ï êC% ú ï ê F%u , n ú ë û ë u,n û þ î

Thus, the spatial harmonic coefficients D% u , n , % Fu , n in the radial diffuser are connected to A%u , n , C%u , n in the inlet duct as é A%u , n ù é D% u , n ù ê% ú ê% ú -1 ê Eu , n ú = Trad , n (r2 ) Bimp , nTax , n × ê Bu , n ú . êC% ú ê F%u , n ú ë û ë u,n û

(20)

E% u , n , and B%u , n , and

(21)

2.7 Boundary conditions To connect upstream and downstream flows, boundary conditions are needed. In the axial duct, the flow is irrotational. Therefore, C%u , n = 0 in the upstream region. Far upstream of compressor, the potential perturbation should be zero since the compressor is the only energy source. Thus, B%u , n = 0 . Downstream of the compressor, the volute dimensions are much larger than those of the diffuser. Therefore, the static pressure perturbation, d Ps ( r = r3 ) , is assumed to be zero at the inlet of

899

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

Table 1. Impeller and diffuser specifications.

the plenum Eq. (22c). C%u , n = 0

(22a)

B%u , n = 0

(22b)

é D% u , n ( s ) ù ê ú éëT3rad , n (3,1) T3rad , n (3, 2) T3rad , n (3,3) ùû × ê E% u , n ( s ) ú = 0 . ê F%u , n ( s ) ú ë û

(22c)

2.8 Eigenvalue determination The six unknown spatial harmonic coefficients, A%u , n , B%u , n , % Cu , n , D% u , n , E% u , n , and F%u , n , can then be solved for by applying the three boundary conditions (Eqs. (22a)-(22c)) to Eq. (21). é A%u , n ù ê ú -1 T3rad , n (3,:) × Trad , n (r2 ) × Bimp , n × Tax , n × ê 0 ú = 0 . ê 0 ú ë û

(23)

Impeller exit diameter (2R2t)

110 mm

Inducer tip diameter (2R1t)

63.4 mm

Inducer hub diameter (2R1h)

20.4 mm

Backsweep angle (from radial)

-35 degrees

No. of impeller blades

18

Test rpm

30,000 rpm (500 Hz)

Diffuser type

Vaneless

Diffuser width (b2)

5.6 mm

Diffuser exit radius

143 mm

Table 2. Test configurations. Configurations

mm, (%, width/R2t)

Case #1

5.6 mm (10.2 %)

Case #2

4.6 mm (8.4 %)

Case #3

3.6 mm (6.5 %)

Air Supply

A%u , n should have a non-trivial solution. Therefore, Eq. (23)

RPM Control Valve

can be rewritten as: T3rad , n (3,:) × Trad , n (r2 ) -1 × Bimp , n × Tax , n (:,1) = 0 .

(24)

Thus, for the n-th mode, the eigenvalue, ln can be found. In the ln = s n + jwn satisfying Eq. (24), s n indicates the growth rate and wn the rotation rate. If s n has a positive value, the flow becomes unstable and develops into a rotating stall with a speed of wn .

3. Description of experiment To validate the new analytical model, rotating stall in a centrifugal compressor with three different diffuser widths has been measured, because this centrifugal compressor has a diffuser inlet rotating stall first from the results by Yoon et al. [14]. As reported in Yoon et al. [14], three different diffuser width cases show different rotating stall phenomena and they can be good validation of this stability modeling. Thus, this centrifugal compressor with a vaneless diffuser at Seoul National University has been used. This facility consists of an air supply system, a centrifugal compressor driven by a radial turbine, a flowmeter, and control valves. The schematic of this facility is shown in Fig. 2. The turbine is driven by the flow from a diesel engine-driven air supply system; the turbine inlet flow has a working pressure of 7 bar and a flow rate of 24.1 m3/min. The valve between the air supply and the turbine is controlled to set the rotational speed of the compressor. Opening the valve increases the rotational speed of the compressor and the turbine. For throttling the compressor, a second valve downstream of the centrifugal compressor is controlled separately. Upstream of the centrifugal compressor, the flow passes through a flowmeter which consists of an 80-

Chamber and Flowmeter

Centrifugal Compressor

Bearing

Turbine

Throttle Valve

Fig. 2. Schematic of test facility.

mm diameter nozzle and four 60-mm diameter nozzles in a chamber with dimensions of 800 mm × 800 mm × 1100 mm (height × width × length). Two screens are located in front of this nozzle bank to block dust and make the flow uniform. Subsequently, the flow enters the centrifugal compressor and exits to atmosphere through a throttle valve. The specifications of the centrifugal compressor are indicated in Table 1 and geometry parameters are shown in Fig. 3. The impeller has 18 impeller blades without splitter blades. In Table 2, test configurations are listed. For these configurations, shims are inserted in the diffuser to reduce the diffuser width as shown in Fig. 4. For the accurate model validation, the rotating stall phenomena are investigated in these 3 configurations. The compressor test conditions are 30000 rpm and 0.12 kg/s at which the relative impeller exit Mach number is 0.17 and the relative impeller exit Reynolds number is about 4×105. Relative impeller exit Reynolds number is based on the properties at the impeller exit, the relative velocity at the impeller exit, and the meridional length of the impeller blade. For each configuration, the compressor characteristic has been obtained. The compressor operating characteristic is

900

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

2 Case #1 Case #2 Case #3

1.5

y

ts

1 0.5 0 -0.5 -1 0

0.2

0.4 0.6 Flow coefficient

0.8

1

Fig. 5. Total-to-static pressure coefficient across whole compressor vs. flow coefficient for 3 diffuser width cases. Fig. 3. Centrifugal compressor geometries.

4. Results and discussion 4.1 Impeller and compressor characteristic The compressor characteristic (an input parameter in this model) is denoted by the flow coefficient and the pressure rise across the compressor. The flow coefficient is the ratio of the axial flow velocity at the impeller inlet and the impeller inlet tip velocity (Eq. (25)). j=

Vx1 . Ut

(25)

The pressure coefficient is the pressure rise in the compressor normalized by impeller inlet tip dynamic pressure as in Eq. (26). Fig. 4. Diffuser width configurations.

y ts =

determined by the flow coefficient and the pressure coefficient. Mass flow rate is measured with the nozzle flowmeter and temperature is measured at the impeller inlet and the diffuser exit with K-type thermocouples. Static pressure has been measured at the impeller inlet, diffuser inlet, and diffuser exit. Each pressure tap is connected to a Scanivalve with 48 channels, and the pressure in each tap is measured with an MKS pressure transducer. To detect stall, six high frequency response pressure transducers, Kulite XCQ-080, have been installed circumferentially at the diffuser inlet. Since the full scale output of this sensor is 100 mV, TEAC SA-59 amplifiers with a maximum gain of 2000 have been used to amplify the signal. The pressure signals from the amplifier are saved in a digital data recording oscilloscope (YOKOGAWA DL750) with a maximum sampling rate of 10 MS/s and 12 bit resolution. The data are then transferred to the computer through an SCSI card.

DPts

rU t 2

.

(26)

Fig. 5 shows the overall pressure coefficient (impeller and diffuser) plotted versus the flow coefficient for three different diffuser width configurations. Also, instability onset point is denoted for each case. As the diffuser width is decreased, the instability is initiated at a lower flow coefficient. This diffuser width affects stability as well as performance. Fig. 6 shows the total-to-static pressure coefficient across only the impeller. d Y ts in Eq. (15) is obtained from these data. The slope dj Even in different diffuser widths, the impeller characteristics for those configurations have similar curves. It means impeller characteristic is not affected by diffuser width. However, the instability in this compressor is more sensitive to the diffuser than the impeller. This is because this instability is initiated at the inlet of the diffuser. This result was also shown in Yoon et al. [14].

901

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905 phi=0.315

20

2 Case #1 Case #2 Case #3

1.5

0

phase

ts

1

y

mode 1

-20 0 10

0.5

10

30

10

20

30

50

-0.5 -1 0

40

50

mode 2

0

-10 0 100

0

20

40

50

mode 3

0 0

10

20

30

40

50

impeller revolution

0.2

0.4 0.6 Flow coefficient

0.8

1

Fig. 6. Total-to-static pressure coefficient across only the impeller for 3 diffuser width configurations. 20

Fig. 8. Spatial Fourier Coefficient at f = 0.315 for Case #1.

divides the signal with a periodic function so that the mode of the signal can be found as

phi=0.315 10

ìï

å a (t ) × e îï

0

magnitude

¥

d P (q n , t ) = Re í

ikq

k

k =0

100

üï ý þï

(27)

phi=0.287 50

where d P (q , t ) is the pressure perturbation, q is the circumferential location, ak (t ) is the spatial Fourier coefficient of the pressure perturbation. If a set of pressure perturbations at N circumferential locations is measured, the spatial Fourier coefficients are derived as follows:

0 100

phi=0.246 50 0 0

0.1

0.2

0.3

0.4

0.5

0.6

normalized frequency by shaft frequency ))

Fig. 7. FFT analysis at each operating point for Case #1.

4.2 Rotating stall inception point Garnier et al. [15] has defined the stall inception point by using power spectrum density (PSD) of fast Fourier transform (FFT) from unsteady pressure data. Thus, the stall onset point in this research is defined to be where the magnitude of the FFT result reaches 30% of the FFT magnitude in the fully developed rotating stall. In addition, precursor is defined as an instability with a magnitude smaller than 30% of the fully developed stall FFT magnitude. Garnier [15] has found the PSD increases as rotating stall developings. In Fig. 7, FFT amplitudes are plotted vs. frequency for each flow coefficient in Case #1. The frequency is normalized by shaft frequency. Here, the rotating stall is fully developed at f = 0.246 and the magnitude is about 100. Thus, the stall flow coefficient is defined at f = 0.287 because the magnitude at f = 0.287 is over 30 and less 30 at f = 0.315 . From this way, rotating stall inception flow coefficient can be experimentally determined. However, at f = 0.315 , a small peak at 0.27 normalized frequency exists. To analyze this peak, spatial Fourier transform (SFT) is also adopted as Garnier et al. [15], Tryfonidis et al. [16], and Kang [17] used. As shown in Eq. (27), this SFT

ak (t ) =

1 N

N

åd P(q , t ) × e n

ikq

(28)

n =1

where q n are the sensor location, k is the spatial mode (stall cell number), and N is the number of location where the pressure is measured. The spatial Fourier coefficient (SFC), ak (t ) , is a complex function of time. When small-amplitude perturbations grow, the phase of the SFC increases linearly. This linear increase region is the traditional stall precursor (Tryfonidis et al. [16]). Fig. 8 shows the SFC at j = 0.315 for Case #1. As the phase of Mode 3 increases linearly, the small peak on 0.27 of normalized frequency at f = 0.315 in Fig. 7 is due to this third mode stall precursor. After that, in Fig. 9, SFC at f = 0.287 has a linear increase in Mode 2 when rotating stall cell is developing. The other stall flow coefficients and precursors for Case #2 and #3 are experimentally determined by this method. As a result, the precursor with Mode 2 at Case #2 appears at f = 0.234 , corresponding to the small bump near 0.14 normalized frequency in Fig. 10 because the SFC in Fig. 11 shows linear increase at 2nd mode. After that, stall cell is developed with Mode 2 at f = 0.210 in the same manner. As flow coefficient decreases more from this stalling flow coefficient, the stall cell number is kept to be two. In addition, as seen in Figs. 12 and 13, Case #3 also has a Mode 2 precursor, first at f = 0.173 and two rotating stall cells after stalling flow coefficient of f = 0.152 .

902

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905 phi=0.287

20

50 0

0

10

20

30

40

0

50

magnitude

phase

-50 0 10

mode 2

-10 0 40

10

20

30

40

50

100

0 100

mode 3

10

20

30

40

phi=0.152

50

phi=0.135

50

20 0 0

phi=0.173

10

mode 1

0 0

50

0.1 0.2 0.3 0.4 0.5 normalized frequency by shaft frequency

impeller revolution

0.6

Fig. 12. FFT analysis at each operating point for Case #3.

Fig. 9. Spatial Fourier coefficient at f = 0.287 for Case #1.

50

20

25

phi=0.234

10

phase

magnitude

100

phi=0.210

50 0

phi=0.192

0 0

0.1 0.2 0.3 0.4 0.5 normalized frequency by shaft frequency

10

20

30

40

100

50

60

mode 2

0 -100 0 200

100 50

mode 1

0 -25 0 200

0

10

20

30

40

50

60

mode 3

100 0.6

0 0

10

20 30 40 impeller revolution

Fig. 10. FFT analysis at each operating point for Case #2.

50

60

Fig. 13. Spatial Fourier coefficient at f = 0.173 for Case #3. 100 50

10

20

30

40

60

mode 2

0 10

20

30

40

50

50

60

0 -0.1 -0.2

mode 3

0 -50 0

50

0.1

100 -100 0 100

mode 1 mode 2 mode 3 mode 4

0.2

Growth rate

phase

-50 0 200

0.3

mode 1

0

-0.3

10

20 30 40 impeller revolution

50

60

Fig. 11. Spatial Fourier coefficient at f = 0.234 for Case #2.

-0.4 0.2

0.25

0.3 Flow coefficient

0.35

0.4

Fig. 14. Predicted growth rate for Case #1.

The stalling flow coefficient predicted by the new model can be obtained from the growth rate of the eigenvalue. Figs. 14-16 show the growth rate vs. flow coefficient for cases #1, #2, and #3, respectively. As previously explained, the positive growth rate (growth rate > 0) indicates instability in the centrifugal compressor because the velocity and pressure in Eq. (7) have divergent values in time. For Case #1 (Fig. 14) which has the biggest diffuser width, the growth rate becomes zero first in Mode 3 at f = 0.298. Thus, the predicted stalling flow coefficient for Case #1 is 0.298. Similarly, the predicted stalling flow coefficients for Case #2 and #3 are 0.247 and

0.198, respectively. In addition, in both Case #2 and #3, Mode 2 reaches zero first rather than the other modes. Thus, the predicted and measured stalling flow coefficients are plotted vs. diffuser width in Fig. 17. The predicted and measured values show similar dependence on diffuser width. 4.3 Stall cell number and speed Using the impeller characteristic slope from experiment as an input, the model can predict not only the stall inception

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

903

0.2 mode 1 mode 2 mode 3 mode 4

0.1

Growth rate

0 -0.1 -0.2 -0.3 -0.4 -0.5

0.2

0.25 0.3 Flow coefficient

0.35

0.4

Fig. 15. Predicted growth rate for Case #2.

Fig. 18. Predicted and measured rotation rate for the Case #1.

0.2 mode 1 mode 2 mode 3 mode 4

0.1

Growth rate

0 -0.1 -0.2 -0.3 -0.4 -0.5 -0.6 0.1

0.15

0.2 0.25 Flow coefficient

0.3

0.35

Fig. 16. Predicted growth rate for Case #3.

Fig. 19. Predicted and measured rotation rate for the Case #2.

stalling flow coefficient

0.3

0.25

0.2

modeling experiment

0.15 60

70

80 diffuser width

90

100

Fig. 17. Predicted and measured stalling flow coefficients for 3 diffuser width cases.

flow coefficient but also the number and speed of the stall cells. From the time lag among signals from different transducers, the rotating stall cells' number, propagation direction, and speed can be found. As the stalling flow coefficient was predicted in Fig. 14, the third mode’s growth rate among other modes first becomes positive at f = 0.298 in Case #1. Thus, the predicted rotating stall onset occurs at f = 0.298 with three stall cells. However, the dominant unstable mode changes into the second mode at f = 0.277 (Fig. 14). Thus, the

model predicts a reduction in cell numbers from three to two. When it comes to other cases, Fig. 15 shows that Case #2 has a dominant mode of Mode 2 and instability at. In addition, the dominant mode of the highest growth rate at a given flow coefficient less than stalling flow coefficient doesn’t change in Fig. 15. As flow coefficient decreases, Mode 2 first becomes zero earlier at f = 0.247 and then maintains a greater growth rate than any other modes. Thus, this indicates the instability in this case has two stall cells during the rotating stall development process. In Fig. 16, Case #3 shows same trend with Case #2. The instability starts first at f = 0.198 in Mode 2 and the dominant mode does not change, either. The speeds of rotating stall cell at each case are also shown in Figs. 18-20 for Case #1, #2, and #3, respectively. Each mode of instability can have its own rotation rate corresponding to the speed of rotating stall cell. However, as explained before, the main mode in the growth rate dominates the main instability mode, so the rotating stall cell speed is indicated with a thick line on the rotation rate of each mode. Since the instability at Case #1 is initiated with Mode 3 at f = 0.298 and then the main mode is changed into Mode 2 at f = 0.277 in Fig. 14, the speed of real instability is shown with thick lines in Fig. 18. Here, the cell speed of Mode 3 is predicted to be increased from f = 0.298 to f = 0.277 . How-

904

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

Fig. 20. Predicted and measured rotation rate for the Case #3.

ever, when the number of stall cell changes, the speed can be predicted to be dramatically reduced at f = 0.277 and then increased again from f = 0.277 as flow coefficient decreases. In Fig. 18, measured stall cell speeds have also been marked with circles and cross markers as well as cell number. The circle marker indicates rotating stall and the cross marker indicates a stall precursor. Experimentally, a three-cell stall precursor is visible at f = 0.315 , and a two-cell stall occurs at f = 0.287. The predicted stall cell speeds increase slightly with decreasing flow coefficient, and this trend is also shown in this experimental result. Consequently, in both analysis and experiment, the stall cell number changes from three to two. Furthermore, in modeling, the cell speed is increased at same stall cell number as flow coefficient decreases. However, once cell number is reduced, it can be predicted that the speed is decreased dramatically. Such change in stall cell number is not uncommon. Kaemmer and Rautenberg [12] experimentally found rotating stall with a cell number varying between one and two in a centrifugal compressor. Haupt et al. [13] also reported, by experiment, several stall patterns of two and three cells at one operating point in a centrifugal compressor. Such research found stall cell number changes in the throttling process but no mechanism was explained. This paper is the first to present the cell number change by the aforementioned modeling. The results of rotating stall cell speed in the other cases are shown in Figs. 19 and 20. As flow coefficient decreases, Case #2 and #3 have increasing rotating stall cell speed at Mode 2 because the main mode of instability is Mode 2 in growth rate in Figs. 15 and 16. Experimental results also show similar trend to this result. Thus, the experimental results can validate this model. From the analyitical results, the diffuser width affects the stall cell phenomena and the decreased diffuser width enhances stability, reducing the stalling flow coefficient and reducing stall cell number up to a certain point.

5. Conclusions The conclusions from this study are as follows.

(1) To understand rotating instability in centrifugal compressors, a new 2-D stability model using compressor characteristics has been developed. Inlet flow in the axial duct and exit flow in the radial duct are connected with impeller pressure rise. This relation produces eigen-equation and enables prediction of the stalling flow coefficient, stall cell number, and speed. These predictions have been compared to experimental results from compressors with three different diffuser widths. (2) The new model can accurately predict centrifugal compressor instability inception, rotating stall cell number, and speed from the steady impeller characteristic data. (3) The rotating stall cell number can change from three cells to two cells in some configurations. In the stability system, a dominant mode exists and it determines the stall cell number. In this model, it is found that the dominant mode can be changed as flow coefficient decreases at a given geometry. (4) The stall cell speed increases for a given same stall cell number as the flow coefficient decreases. However, the stall cell speed can be reduced when the cell number decreases. (5) Diffuser width changes only the impeller exit flow, and decreased diffuser width decreases the stalling flow coefficient. In addition, small diffuser width can reduce the stall cell number as Case #1 with large width has three cells but Case #2 and #3 have two cells at stall inception.

Acknowledgment This research has been supported by BK21 from the Korean government.

Nomenclature-----------------------------------------------------------------------A AR B C D E F P Q R T U V X p r t v x

G W a

: Coefficient in Fourier series [-] : Area ratio [-] : Coefficient in Fourier series [-] : Coefficient in Fourier series [-] : Coefficient in Fourier series [-] : Coefficient in Fourier series [-] : Coefficient in Fourier series [-] : Pressure [N/m2] : Normalized mass flow rate [-] : Mean radius [m] : Transmission matrix [-] : Speed [m/s] : Velocity [m/s] : Axial position [m] : Normalized pressure [-] : Normalized radius [-] : Time [sec] : Velocity normalized by tip velocity [-] : Normalized axial length [-] : Normalized circulation [-] : Impeller frequency [1/sec] : Absolute exit angle [rad]

Y.-S. Yoon and S. J. Song / Journal of Mechanical Science and Technology 28 (3) (2014) 895~905

b c d e f g l q

r t w y z

: Relative exit angle [rad] : Vorticity perturbation in axial position [-] : Perturbation [-] : Eccentricity [m] : Flow coefficient [-] : Stagger angle [rad] : Rotor/impeller inertia [-] : Circumferential position [rad] : Density [kg/m3] : Normalized time [-] : Vorticity [-] : Pressure rise coefficient, stream function [-] : Vorticity perturbation in time [-]

Subscripts and superscripts ax imp n r rad t ts u x

q 1 2 3

: Axial duct : Impeller : n-th harmonic : Radial direction : Radial duct : Impeller inlet tip : Total-to-static : Unsteady, unknown value : Axial direction : Circumferential direction : Impeller inlet : Impeller exit : Diffuser inlet

Operators

d() 0 0%

: Perturbation quantity : Mean value : Fourier spatial coefficient

References [1] E. M. Greitzer, Surge and rotating stall in axial flow compressors - 1. Theoretical compression system model, Journal of Engineering Power, Transactions ASME, 98 (2) (1976) 190-198. [2] E. M. Greitzer, Surge and rotating stall in axial flow compressors – 2. Experimental results and comparison with theory, Journal of Engineering Power, Transactions ASME, 98 (2) (1976) 199-217. [3] F. K. Moore and E. M. Greitzer, A theory of post-stall transients in axial compressors: Part I – Development of the equations, Journal of Engineering for Gas Turbines and Power, Transactions ASME, 108 (1) (1986) 68-76. [4] T. P. Hynes and E. M. Greitzer, A method for assessing effects of circumferential flow distortion on compressor stability, Journal of Turbomachinery, Transactions ASME, 109 (3) (1987) 371-379. [5] M. B. Graf, T. S. Wong, E. M. Greitzer, C. S. Tan, H.-W.

905

Shin and D. C. Wisler, Effects of nonaxisymmetric tip clearance on axial compressor performance and stability, Journal of Turbomachinery, Transactions ASME, 120 (4) (1998) 648-661. [6] K. Gordon, Three-dimensional rotating stall inception and effects of rotating tip clearance asymmetry in axial compressors, Ph D thesis, Department of Aeronautics and Astronautics, MIT. (1999). [7] T. R. Camp and I. J. Day, A study of spike and modal stall phenomena in a low-speed axial compressor, Journal of Turbomachinery, Transactions ASME, 120 (3) (1998) 393-401. [8] Z. S. Spakovszky, Applications of axial and radial compressor dynamic system modeling, Ph.D. thesis, Department of Aeronautics and Astronautics, MIT. (2001). [9] Z. S. Spakovszky, Backward traveling rotating stall waves in centrifugal compressors, Journal of Turbomachinery, Transactions ASME, 126 (1) (2004) 1-12. [10] J. D. Paduano, A. H. Epstein, L. Valavani, J. P. Longley, E. M. Greitzer and G. R. Guenette, Active control of rotating stall in a low-speed axial compressor, Journal of Turbomachinery, Transactions ASME, 115 (1) (1993) 48-56. [11] D. L. Gysling, Dynamic control of rotating stall in axial flow compressors using aeromechanical feedback, Ph.D. thesis, Department of Aeronautics and Astronautics, MIT. (1993). [12] N. Kaemmer and M. Rautenberg, A distinction between different types of stall in a centrifugal compressor stage, Journal of Engineering for Gas Turbines and Power, 108 (1) (1986) 83-92. [13] U. Haupt, A. N. Abdelhamid, N. Kaemmer and M. Rautenberg, Excitation of blade vibration by flow instability in centrifugal compressors, ASME Paper No. 86-GT-283. [14] Y. S. Yoon, Shin. H. Kang and S. J. Song, Clearance effects on centrifugal compressor characteristic and stability, ASME Paper No. 2007-GT-27566. [15] V. H. Garnier, A. H. Epstein and E. M. Greitzer, Rotating waves as a stall inception indication in axial compressors, Journal of Turbomachinery, Transactions ASME, 113 (2) (1991) 290-302. [16] M. Tryfonidis, O. Etchevers, J. D. Paduano, A. H. Epstein and G. J. Hendricks, Prestall behavior of several high-speed compressors, Journal of Turbomachinery, Transactions ASME, 117 (1) (1995) 62-80. [17] J. Kang and S. Kang, Stall inception in a high-speed centrifugal compressor, ASME Paper No. 2001-GT-0301 (2001).

Yong-Sang Yoon received his Ph.D at Mechanical and aerospace engineering from Seoul National University in 2006. He had worked at LG Electronics(from 2006 to 2011) and Whittle Laboratory at University of Cambridge(from 2011 to 2013). He is now a senior research engineer in Samsung Techwin, Korea.