D. L.

Rhode
Assistant Professor.
Prediction of Incompressible Flow
J. A. Demko
Graduate Student. in Labyrinth Seals
A new approach was developed and tested for alleviating the substantial con-
U. K. Traegner vergence difficulty which results from implementation of the QUICK differencing
Graduate Student. scheme into a TEACH-type computer code. It is relatively simple, and the resulting
CPU time and number of numerical iterations required to obtain a solution compare
G. L. Morrison favorably with a previously recommended method. This approach has been
employed in developing a computer code for calculating the pressure drop for a
Associate Professor.
specified incompressible flow leakage rate in a labyrinth seal. The numerical model
is widely applicable and does not require an estimate of the kinetic energy carry-over
S. R. Sobolik coefficient for example, whose value is often uncertain. Good agreement with
Graduate Student. measurements is demonstrated for both straight-through and stepped labyrinths.
These new detailed results are examined, and several suggestions are offered for the
Mechanical Engineering Department, advancement of simple analytical leakage as well as rotordynamic stability models.
Texas A&M University,
College Station, TX 77843

Introduction
Labyrinth seals play a vital role in gas and steam turbines, perimental data correlations to purely analytical expressions
compressors, and high-capacity pumps. In pumps the sealing involving friction factor for example. One example of an
objective is generally to minimize the leakage flow around the analytical model is that developed by Bilgen and Akgungor
impellers, whereas for turbines it is often to maintain a [3]. They assumed turbulent couette and Poiseuille flows in
minimal leakage flow for cooling hot parts. The seal designer the tangential and axial directions, respectively. In estimating
accomplishes this by providing a highly dissipative flow path the shear stress along the free shear layer, Han [4] recently
between high and low pressure regions. The flow passage considered the leakage flow energy consumed in driving the
through a labyrinth is illustrated in Fig. 1 for a simple straight- cavity vortex while estimating the wall shear stresses using
through seal configuration. Each tooth converts a portion of Prandtl's flat plate boundary layer solutions. Another exam-
the available pressure head into mean flow kinetic energy, ple is the work of Nikitin and Ipatov [5], who approximated
some of which is dissipated within the cavity immediately frictional losses via analysis and data correlation, showing
downstream. how the friction coefficient varies with tooth width, etc.
Numerical simulation of labyrinth seal leakage employing
the governing partial differential equations is relatively com-
plicated and expensive. A recent literature search indicated Objective
that only Stoff [1] and Rhode, et al. [2] have reported such an
investigation. Stoff's solution, which utilized the upwind Depending on flowfield details, labyrinth seal solutions us-
Hybrid scheme, corresponded to his large scale experimental ing the finite difference procedure of Stoff [1] may suffer
water facility consisting of a straight-through series of generic from false diffusion numerical error, which is discussed in a
cavities as illustrated in Fig. 1. He compared a single radial later section. Simple analytical models often predict leakage
profile of predicted mean swirl velocity, rms swirl velocity and rates only within about a factor of two, even for simple
turbulence dissipation rate with his corresponding straight-through seal geometries such as that of Fig. 1.
measurements for an axial station midway between adjacent However, these modellers have done well considering the
teeth. Rhode, et al. developed a finite difference procedure to flowfield complexity and the limited detailed information
compute the compressible flow cavity pressure drop at various available regarding rotating cavity flows. Thus there is a need
leakage Mach numbers. These values were subsequently used for a leakage model of improved accuracy which is widely ap-
in a Fanno analysis to predict the leakage rate and cavity-to- plicable and does not involve the uncertainty of an assumed
cavity stator wall pressure distribution. kinetic energy carry-over coefficient for example. The present
For incompressible flow, several investigators have threefold objective is: to present a new approach for
developed simple algebraic leakage models ranging from ex- alleviating the convergence difficulty that occurs upon im-
plementation of the QUICK scheme; to demonstrate the wide-
ly applicable capability for obtaining realistic labyrinth seal
Contributed by the Fluids Engineering Division of THE AMERICAN SOCIETY OF flowfield solutions; and to provide additional flowfield details
MECHANICAL ENGINEERS and presented at the Energy-Sources Technology Con-
ference and Exhibit, New Orleans, La., February 12-16, 1984. Manuscript to allow more realistic approximations for the refinement of
received the Fluids Engineering Division, May 4, 1984. simple analytical models.

Journal of Fluids Engineering MARCH 1986, Vol. 1 0 8 / 1 9

Copyright © 1986 by ASME
Downloaded From: http://fluidsengineering.asmedigitalcollection.asme.org/ on 01/20/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use
 y-Siator Housing 1- 11 2 3... HI

'X - :: : " ^
IrT GUT

O O o r—" i™~ x
1
1
i
_i
i3-
.: ~rn
1
1
Rotating Surface- 1
Fig. 1 Labyrinth seal with generic cavities illustrating the expected I 3-
1 L_
streamline pattern

Numerical Procedure
| h.
General Methodology. The numerical procedure em- r.v
ployed is based on that of the TEACH computer program of e,w
Gosman and Pun [6]. The recent QUICK convective differenc-
-e-
ing scheme of Leonard, which is described below, was im- Fig. 2 Computational domain for a typical labyrinth seal cavity with a
plemented in a version of the code which had previously been coarse grid illustrating the stairstep approximation
extended by Lilley and Rhode [7] to include swirl momentum.
The computational domain consists of a single cavity as shown
in Fig. 2. streamwise periodic flow occurring downstream of the first
The governing equations for this axisymmetric flow may be two or three cavities.
expressed in the general form The user specifies a desired leakage flow rate, and the com-
puter code calculates the associated cavity bulk pressure drop.
.9
Observe that by making several such computations, each with
T[ dx
(pur<t>) + — (pvr<j>)
dr dx (<v£) a different leakage rate, one may plot the leakage-pressure
drop characteristic which is of considerable interest to the
dr dr
-)]-* designer.
There are several techniques for simulating sloping wall
where <j> represents any of the dependent variables, and the boundaries such as that of the typical cavity geometry shown
equations differ primarily in their final source terms S*[7]. It in Fig. 2. Sophisticated approaches include utilization of a
is the standard high Reynolds number version of the two- body-fitted coordinate system or a coordinate transformation
equation k-e turbulence model [8] which has been employed for mapping the physical plane onto an idealized computa-
thus far. The presence of swirl momentum directly influences tional domain. These methods require extensive re-
the values of k and e, and in turn, the values of turbulent formulation of the governing differential equations, boundary
viscosity through the swirl velocity derivatives in the expres- conditions and formulae for length, area and volume.
sion for turbulence energy production. Moreover, expected numerical instabilities using the new coor-
Boundary values at the inlet of a particular cavity are dinates along with the added complexity of the QUICK
naturally very important, but unfortunately were unknown, as scheme render this a difficult task. Thus, in view of computer
necessary quantitative flow measurements were unavailable budget and manpower limitations the simpler stairstep ap-
even at a single operating condition. Due to this lack of inlet proximation illustrated in Fig. 2 was adopted for the typical
boundary data, for each numerical iteration the inlet values of cavity configuration.
each variable (except pressure) were assumed equal to the Locally, wall shear stress and heat flux on such a simulated
latest corresponding outlet values. This practice assumes the wall will not generally agree closely with measurements.
presence of a series of geometrically identical cavities and a However, bulk as well as local velocities, etc. a slight distance

Nomenclature

A = coefficient of finite difference equation

B = auxiliary coefficient
c = tooth clearance r = turbulent exchange coefficient
d = radial distance from cavity base to stator A = internodal mesh spacing
wall e = turbulence dissipation rate
HDS = hybrid upwind/central differencing scheme <j> = generalized dependent variable
k = turbulence kinetic energy p = time-mean density
L = axial length of seal cavity v = kinematic viscosity
Pe = grid Peclet number r = shear stress
R = stator wall radius Q = shaft speed
rsh = shaft radius to the base of a cavity
Re = Reynolds number Subscripts
Sp,Su = components of linearized source term n,s,e,w = north, south, east, west faces of a control
S = source term volume
Ta = Taylor number N,S,E,W = first, i.e., immediately adjacent neighboring
U = bulk axial velocity at flow inlet grid points to the North, South, East, West
V=(u,v,w) = time-mean velocity (in x, r, 8 direction) NN,SS, . . . = second neighboring grid points to the North,
W = bulk swirl velocity South, East, West
x,r,8 = axial, radial, azimuthal cylindrical polar eff = effective
coordinates lam = laminar

20/Vol. 108, MARCH 1986 Transactions of the ASME

Downloaded From: http://fluidsengineering.asmedigitalcollection.asme.org/ on 01/20/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use

 NN 1.0
fN
.75
<p*
y/a .50 HDS
grTTTTTTnil QUICK 3p,
QUICK 5p.
WW w E EE .25

.6 .8 1.2
ss u/U
Fig. 4 Predicted axial velocity values at x/a = 1.0 using: (a) Hybrid; (6)
Fig. 3 Illustration of quadratic upwind interpolation for <t>w on the west three-point QUICK; and (c) five-point QUICK schemes for laminar flow
face of a control volume with flow from left to right between parallel plates separated by distance 2a
away from the wall have been found to be surprisingly ac- where the last four neighbors were not involved previously us-
curate if careful attention is devoted to stairstep implementa- ing the Hybrid scheme.
tion details. For example Rhode, et al. [9] and Syed [10] found The method of solving the equations is the same as for the
close agreement with measurements using this approach for previous Hybrid scheme formulation. For the solution of the
flow past a backstep with a sloping wall. Also, Sindir [11] difference equations for a vertical column of cells, the
employed nearly the same implementation details and ob- neighbor values along adjacent vertical columns are assumed
tained close agreement with experiment for backstep flows to be temporarily known. Thus the equation for each cell in-
with a deflected wall. Further, the flow along the simulated volves only three unknowns, and the equation is arranged so
sloping wall exhibits extremely low velocity and is far removed that these appear on the left-hand side as
from the all-important leakage flow region.
Implementation of QUICK Differencing Scheme. False The set of these equations for the column of cells forms a
diffusion is a second-order truncation numerical error. The tridiagonal coefficient matrix which is solved using an
upwind portion of the previous Hybrid upwind/central dif- iterative, line-by-line application of the well known tri-
ferencing scheme can result in an erroneous solution which is diagonal matrix algorithm.
overly diffusive. This difficulty occurs where convection Numerical convergence may be impaired as the Af coeffi-
dominates (i.e., where the grid Peclet number P e = WlA/T^ cients may become negative when using QUICK. No false
exceeds 2.0) in the presence of both streamline-to-grid source transient terms are employed as a remedy as in
skewness and diffusive transport normal to the flow direction. reference [13]. Instead, for each iteration the difference equa-
Recall that the upwind scheme is based on the uniform <j> tion for certain cells is re-arranged as needed to promote
distribution assumption in evaluating <j> at any face of a con- diagonal dominance of the tridiagonal coefficient matrix. This
trol volume. The recent QUICK (Quadratic Upstream Inter- re-arrangement is applied to every cell for which any coeffi-
polation for Convective Kinematics) scheme of Leonard [12] cient, for example A%E, is negative. In effect, this ad hoc pro-
greatly reduces false diffusion by utilizing either a five-point cedure subtracts fi*y4|E0P from both sides of the difference
or a three-point, upwind-shifted quadratic interpolation for- equation in a very economical fashion. Observe that if
mula in evaluating <j>. Figure 3 graphically illustrates the three- 5* = 1.0, A$E no longer makes a negative contribution to the
point formula for a west face for example. The three-point in- diagonal coefficient Ap, and thus the diagonal dominance of
terpolation expressions for this face, using a uniform grid for the coefficient matrix is enhanced.
simplicity, are The source term Sfj becomes

ir(<l>p + <t>w)—^(^E -2«p + < A st--••Sb-B*A&<I> +

where <j>^ is the current "in-storevalue of (f>P. Values of 1.0
if u„ is positive and have routinely been used without difficulty for B",B°, and Bw
1 1 even though this value allows the possibility that /l£ = 0.0. If
-(<t>P + </>w) - y ( < ^ - 2 < £w + 0ww) necessary this can easily be precluded by specifying values near
0.8 for example.
if uw is negative. The first term in these expressions represents At present the k and e equations employ the Hybrid scheme
the centered difference formula, and the second is the crucial for computational economy reasons. Leschziner and Rodi [14]
stabilizing upstream-weighted curvature contribution. have supported this practice in explaining that solutions for k
The QUICK expressions for an arbitrary non-uniform grid and e were found to be unaffected by the convected differenc-
were incorporated by expressing each convected quantity <> / in ing scheme, espeically in the shear layer bordering a recircula-
the convection terms of the difference equations with the ap- tion zone. They attribute this to the fact that the source terms
propriate interpolation formula. The standard expressions of these equations are dominant.
were used for all other quantities in the equation. Upon col-
lecting terms the conventional difference equation of the form Numerical Testing
The implementation of QUICK was extensively checked by
AUp=EAt*j+st comparison with corresponding Hybrid scheme results and
measurements where available. The test problems selected in-
where
clude elliptic and parabolic problems with and without tur-
bulence and with and without swirl momentum. The simplest
Af^Af-Sf. of these concerns laminar flow in the boundary layer develop-
ment region between parallel plates separated by a distance of
was obtained. However, using the three-point expressions the 2a. The flow enters with a uniform velocity profile at
summation occurs over Re = 150. This problem exhibits very slight streamline-to-grid
j = E, W, N, S, EE, WW, NN, SS skewness, and thus negligible false diffusion, so that the solu-

Joumal of Fluids Engineering MARCH 1986, Vol. 108/21

Downloaded From: http://fluidsengineering.asmedigitalcollection.asme.org/ on 01/20/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use

 x/D
Table 1 Convergence data for flow between parallel plates
Scheme T
/THDS N/NHDS (T/N)/(T H D S /N H D S )
HDS 1.00 1.00 1.00
QUICK 3 p. 1.34 0.99 1.34 '
r/D QUICK 5 p. 1.76 1.05 1.66
QUICKE 0.84 0.60 1.41
QUICKER 1.85 1.32 1.40
T-computational time
N-number of iterations

Fig. 5 Axial velocity of swirling flow through a pipe expansion showing

predictions using (a) the Hybrid { ) and (b) the three-point QUICK
Table 2 Convergence data for flow through pipe expansion
scheme (-••) along with previous measurements [16] ( + + + ) Scheme T/T N/NHDS (T/N)/(T H D S /N H D S )
HDS
x/D HDS 1.00 1.00 1.00
.275 .55 QUICK 3 p. 1.31 1.00 1.31
QUICK 5 p. 1.62 0.99 1.64
T - computational time
N-number of iterations
Laminar flow through an abrupt pipe expansion ot
diameter ratio 2.0 constitutes a more complicated case in
which a large recirculation eddy exists behind the step. A fully
developed laminar velocity profile was specified at the inlet
where Re = 200. Once again, there was little false diffusion in
the Hybrid scheme solution resulting in no significant devia-
tions among the three solutions employing a 34x22 grid.
Table 2 shows the relative convergence data to be very similar
to that for the case of parallel plates. Note that the three-point
scheme is naturally more economical. Since Leonard also
found little loss of accuracy in using the three-point as op-
posed to the five-point version of QUICK, the five-point for-
Fig. 6 Axial velocity of coswirling annular jets in a pipe showing
predictions using (a) the Hybrid ( ) and (b) the three-point QUICK mula was discontinued.
scheme ( ) along with previous measurements [17] ( + + + ) The next case is that measured by Rhode, et al. [16]. Here
again the flow passes through an abrupt pipe expansion of
77 7 V • diameter ratio 2.0. However, this problem is much more com-
0.470 . 8V\ plicated, involving turbulent swirling flow with an inlet
77
0
\ Re= 1.1 x 105 and upstream swirl vane angles oriented at 38
degrees from the downstream direction. A nonuniform 40 x 25
0.460 -
grid was utilized. Figure 5 shows the improvement of the
V 7
QUICK results over that of the Hybrid scheme for the annular
E 0.450 - 07
I peak value of axial velocity at x/D = 0.5 and /•/£)= 0.4. This
l_
\ improvement was anticipated because substantial false diffu-
0.440 _ 7*
sion arises at x/D = 0.25 from the slope of the streamlines
which is near 45 deg. However, closer agreement with
0.430 .
0.425 .
I 1 1
?<y
v
V
1 1
7 7
measurements was expected. One possible reason for this is er-
ror in the measurement of inlet boundary values due to
02 0.4 0.6 0.8 adverse conditions there, a location immediately downstream
wA
sh<> of the swirler where a non-axisymmetric region of low velocity
Fig. 7 Comparison of present prediction ( ) of dimensionless occurs.
swirl velocity with corresponding LDA measurements (o o o , v v v) of The turbulent co-swirling concentric jet flow in a pipe
Stoff [1]atx/L = 0.5 measured by Vu and Gouldin [17] constitutes another test
tion from both the three-point and five-point QUICK for- case. In this flow, however, measurements were available for
mulations should agree with that of the Hybrid approach. use as inlet boundary values for every quantity including k and
Figure 4 indeed shows a maximum deviation among these e. Based on chamber diameter and bulk velocity,
solutions of only 1.6 percent using a 48 x 20 grid. The value is Re = 2.96 x 105 and the swirl parameter for the inner and outer
only half this much farther downstream. jets was 0.58 and 0.54, respectively. A highly nonuniform
This flow problem, grid distribution, boundary condition 35 x 30 grid was employed. Streamline skewness is very slight,
data, etc. is precisely identical to that used by Pollard and Siu and as shown in Fig. 6, the two solutions are almost identical.
[15] in evaluating their important implementation approaches Also, observe that both solutions are in good agreement with
of the three-point QUICK scheme. Their formulation named measurements.
QUICKE would not always converge to a solution, and thus
they recommend their QUICKER approach. By normalizing Results and Discussion
the required CPU time and number of numerical iterations by
the corresponding Hybrid scheme values, somewhat mean- Experimental Verification for Labyrinth Seals. A
ingful comparisons between the previous QUICKER and the literature search indicated that only the paper by Stoff [1]
present formulations are possible. These values are found in gives detailed quantitative velocity measurements within a
Table 1. Observe that the present three-point implementation labyrinth seal cavity. He utilized a large scale test section
using diagonal dominance enhancement seems to be somewhat model of a straight-through seal with rectangular cavities
more desirable. through which water flowed. The Reynolds number

22/Vol. 108, MARCH 1986 Transactions of the ASME

Downloaded From: http://fluidsengineering.asmedigitalcollection.asme.org/ on 01/20/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use

 • x
65 .85 1.
STREAMWISE
8.0 •
X PERIODIC
REGION
o .SS5
Q-
x
-SC
6.0
Cvl
'O
X
UJ r/R
4.0 •
rr
~>
<n
(j)
I>I
a: Z.O •
a. X

0.0
IN 2 3 4 5 6 OUT
0
TOOTH LOCATION
' k/U 2 x20
Fig. 8 Measured tooth-to-tooth pressure distribution ( x x x ) the (a) GENERIC CAVITY
streamwise periodic prediction ( - ) plotted to match the measured
pressure at the third tooth

u/U r/R
1
(a) GENERIC CAVITY
i ,55 .5 .65

u/U
k/U 2 x20
(b) TYPICAL CAVITY (b) TYPICAL CAVITY
Fig. 9 Axial velocity prediction for (a) the generic and (b) the typical Fig. 10 Predicted turbulence kinetic energy for (a) the generic and (b)
seal cavity employing 33 x 31 QUICK ( ), 53 x 53 Hybrid ( ) and the typical seal cavity using the QUICK scheme
33 x 31 Hybrid ( ) solutions

Re = 2 Uc/v « 3 .4 X 1 0 2 a n d t h e T a y l o r n u m b e r
Ta=(f^rf/j/)(c?//- s / 1 ) w «l.lxl0 4 . Figure 7 shows that swirl tions are: shaft speed fl = 35,410 rpm, cavity axial length
velocity predictions are within 6.0 percent of Stoff's cor- L = 1.113 mm, stator wall radius J? = 42.89 mm, tooth radial
responding measurements. clearance c = 0.216 mm, radial distance from cavity base to
The measurements of Morrison, et al. [18] for the tooth-to- stator wall d= 1.105 mm, and kinematic viscosity
tooth pressure distribution of stepped labyrinth seals flowing j - = 1 . 4 6 2 x l 0 - 7 m 2 /s.
water were employed in a series of further prediction tests. In Not included here are the streamline plots, which show that
the experimental facility, a single static pressure tap was the dividing streamline exhibits a reattachment stagnation
located at each tooth on the stator wall. The dimensionless point which is very slightly below the peripheral corner of the
flow parameters were Re = 4.95 X 103 and Ta = 6.7 x 103. The downstream tooth. This occurs for both cavity geometries. In
result for one case with six cavities is given in Fig. 8. The both cases the slope of a straight line approximation for the
predicted value of bulk pressure drop per cavity for the dividing streamline is near 2 degrees below horizontal. This
streamwise periodic cavities was 105 kPa. This was translated serves to aid the developers of simple analytical models, as an
into the predicted streamwise periodic pressure distribution, assumed value near 7 degrees has been employed by Jerie [19]
which was plotted by graphical construction using this slope for example.
while matching the pressure measured at the third tooth. Note Axial velocities from several solutions for the generic case
that streamwise periodicity commences at approximately the are shown in Fig. 9(a) which reveals that the 33 x 31 QUICK
third tooth. The measured pressure drop per cavity averaged and 53x53 Hybrid solutions are essentially identical. The
123 kPa in this region, yielding a discrepancy of under 15.0 33 X 31 Hybrid and 53 x 53 Hybrid solutions differ very slight-
percent. ly. As with the generic cavity, Fig. 9(b) shows that the 33 x 31
QUICK and the 53 x 53 Hybrid solutions of axial velocity for
Seal Flowfields Considered. The operating conditions of the typical cavity are nearly identical. Note that in this case the
the two cavity geometries investigated are identical. Liquid 33x31 Hybrid solution is not grid independent in the free
hydrogen at 42 °R enters a cavity with tooth-clearance bulk shear layer. The axial velocity values of this solution deviate
velocity [7=338 m/s. The primary dimensionless flow from those of the other two by as much as 20 percent.
Parameters are Re= l.Ox 106 and T a = 1.3 x 105. Other condi- Observe at this point that the QUICK approach exhibits an

Journal of Fluids Engineering MARCH 1986, Vol. 108/23

Downloaded From: http://fluidsengineering.asmedigitalcollection.asme.org/ on 01/20/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use

 x/L
x/L
0. .15 .35 .5 .65 .85
0. .15 .35 .5 .65 .85 1.0

r/R
r/R

0W 5x10 3
x 10 (a) GENERIC CAVITY
0.5/.U x/L
(a)GENERIC CAVITY 0. .15 .35 .5 .65 .85 1.0

r/R

r/R

5xl03
(b) TYPICAL CAVITY
Fig. 12 Predicted effective viscosity from the k-t model for (a) a generic
and (b) a typical seal cavity using the QUICK scheme

0.5/>U in these figures that the typical cavity is expected to achieve a

higher pressure drop due to higher levels of turbulence energy.
(b) TYPICAL CAVITY The distribution of pressure relative to the cavity inlet stator
Fig. 11 Predicted relative pressure for (a) the generic and (b) the typical wall pressure P o w is shown in Figs. 11(a) and 11(6). In the
seal cavity using the QUICK scheme leakage flow region of the generic cavity, the pressure varia-
tion is partly a reflection of the Bernoulli effect of pressure
decrease at constrictions in flow area. A sharp positive
economic advantage over the Hybrid approach. In obtaining a pressure peak, which occurs in an annular fashion near
grid independent solution, utilization of the Hybrid method x/L = 0.85, results from the flow stagnation on the
requires 53 X 53 grid points whereas utilization of QUICK re- downstream tooth. The overall pressure drop from inlet to
quires only 33 X 31 points. The execution time required for the outlet of this cavity is predicted as approximately 6 . 0 x l 0 2
QUICK solution is only 44% of that required by the Hybrid kPa.
version. This is mostly attributed to the well known apparent The pressure distribution within the typical cavity is seen in
sluggish response to grid refinement of the upwind portion of Fig. 11(6) to exhibit similar trends. However, leakage flow
the Hybrid scheme. As demonstrated by Han, et al. [13], this pressure varies in the streamwise direction much more
lack of response is actually a false diffusion effect. gradually for this geometry upstream of x/L = 0.65.
Figures 10(a) and 10(b) illustrate the fundamental physical Downstream of this location it decreases more rapidly than for
mechanism of labyrinth seal performance. The large value of the generic case as the flow area constricts in passing over the
du/dr in the free shear layer yields a high generation rate of downstream tooth. The predicted bulk pressure drop for this
turbulence kinetic energy in this region. A portion of this tur- cavity is approximately 8.0 X 102 kPa, an increase of 33 per-
bulence energy is transported into the recirculation region cent over the generic geometry.
which is clearly seen to act as an effective region of viscous Swirl velocity distributions are not included here. They are
dissipation. This turbulence generation and dissipation within quite uniform with values of 0.65 and 0.69 for the generic and
each cavity is the basic means by which the upstream total typical configurations, respectively, when nondimensionalized
mechanical energy is irreversibly converted into heat. Observe by the circumferential velocity of the cavity base. Simple

24/Vol. 108, MARCH 1986 Transactions of the ASME

Downloaded From: http://fluidsengineering.asmedigitalcollection.asme.org/ on 01/20/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use

 space limitations, the solutions presented here are only a sam-
ple of those computed in evaluating this approach.
The labyrinth seal computer code is seen to provide realistic
results for both straight-through and stepped labyrinths.
There is no worrisome uncertainty of assuming a value for the
mean flow kinetic energy carry-over coefficient for example.
The present solutions allow the advancement of simple
analytical leakage models via refinement of various important
re re sh approximations. Results of particular interest in this regard
for generally similar cavities are:
(i) the dimensionless swirl velocity is accurately
represented by a uniform value near 0.65;
A D T* =T /3Pr U2
(a) FREE SHEAR LAYER rx rx (ii) the slope of the dividing streamline may be approx-
STRESSES imated as 2 degrees below the horizontal;
T * =T /.5/.U 2
xr xr (iii) the predicted shear stress distributions on planes where
2.0 they are often required by analytical methods are presented
10.0 for modelling approximations as appropriate.
1.6
8.0 Acknowledgment
•°o 1-2 e The developments and results presented here are a by-
% product of two studies for which the financial support of
KQ8 AFOSR and NASA Marshall Space Flight Center are grate-
4.0*-
fully acknowledged. Also, the expert typing skills of Ms.
0.4 Bich-N. Ho are greatly appreciated.
2.0
References
0.0
0.0 1 Stoff, H., "Incompressible Flow in a Labyrinth Seal," Journal of Fluid
A B C D
(b) DISSIPATION CHAMBER WALL STRESSES Mechanics, Vol. 100, 1980, pp. 817-829.
2 Rhode, D. L., and Sobolik, S. R., "Simulation of Subsonic Flow in
Fig. 13 Predicted shear stress components along various planes (see Labyrinth Seals," accepted for publication in ASME Journal of Engineering for
Fig. 10) of the generic cavity Gas Turbines and Power, Vol. and page nos. unavailable.
3 Bilgen, E., and Akgungor, A. C , "Leakage and Frictional Characteristics
models may employ uniform values in this range for similar of Large Diameter Labyrinth Seals," National Conference on Fluid Power,
configurations. 1972, pp. 85-101.
The isotropic turbulent viscosity values resulting from the 4 Han, J . T . , "A Fluid Mechanics Model to Estimate the Leakage of Incom-
pressible Fluids through Labyrinth Seals," ASME Paper 79-FE-4, presented at
k-e turbulence model via the Joint ASME/CSME Applied Mechanics, Fluids Engineering and
Bioengineering Conference, Niagra Falls, June 18-20, 1979.
/*eff = C „ p £ 2 / e + /*lam 5 Nikitin, G. A., and Ipatov, A. M., "Design of Labyrinth Seals in
are given in Fig. 12(a) and 12(d). In accordance with tur- Hydraulic Equipment," Russian Engineering Journal, Vol. 53, No. 10, 1973,
bulence energy values, the typical cavity exhibits the greater pp. 26-30.
6 Oosman, A. D., and Pun, W. M., "Calculation of Recirculating Flows,"
viscosity values. Viscosity profiles are much more uniform Rept. No. HTS/74/2, 1974, Dept. of Mech. Engr., Imperial College, London,
than turbulence energy profiles because, even though ,ueff England.
varies with turbulence energy squared, the dissipation values 7 Lilley, D. G., and Rhode, D. L., " A Computer Code for Swirling Tur-
in the free shear region are more than a hundred times larger bulent Axisymmetric Recirculating Flows in Practical Isothermal Combustor
Geometries," NASA CR-3442, Feb. 1982.
than in the recirculation zone. 8 Launder, B. E., and Spalding, D. B., "The Numerical Computation of
Figure 13(a) gives the magnitude of both components of Turbulent Flows," Comp. Methods in Appl. Mech. andEngrg., Vol. 3, 1974,
shear stress coefficient for the generic cavity along plane AD pp. 269-289.
9 Rhode, D. L., Lilley, D. G., and McLaughlin, D. K., "On the Prediction
of the free shear layer as labeled in Fig. 10(a). Observe that the of Swirling Flowfields Found in Axisymmetric Combustor Geometries," ASME
axial component is an order of magnitude larger than the cir- Journal of Fluids Engineering, Vol. 104, 1982, pp. 378-384.
cumferential component. Further, it drastically increases near 10 Syed, S. A., personal communication, 1984.
the reattachment point, as the flow undergoes large gradients 11 Sindir, M., "Calculation of Deflected-Walled Backward-Facing Step
Flows: Effects of Angle of Deflection on the Performance of Four Models of
in this region of maximum effective viscosity. Turbulence," ASME paper No. 83-FE-16.
The shear stress components along each of the three walls of 12 Leonard, B. P., "A Stable and Accurate Convective Modelling Procedure
the dissipation chamber have been evaluated in simple Based on Quadratic Upstream Interpolation," Comp. Meths. Appl. Mech.
analytical models [4] using Prandtl's flat plate turbulent boun- Eng., Vol. 19, 1979, pp. 59-98.
13 Han, T., Humphrey, J. A. C , and Launder, B. E., " A Comparison of
dary layer solution for example. The upstream corner of each Hybrid and Quadratic-Upstream Differencing in High Reynolds Number Ellip-
wall was taken as the leading edge. Figure 13(d) shows the pre- tic Flows," Comp. Meths. Appl. Mech. Eng., Vol. 29, 1981, pp. 81-95.
sent predictions which indicate different behavior. The peaks 14 Leschziner, M. A., and Rodi, W., "Calculation of Annular and Twin
arising midway along CB and BA are largely attributed to the Parallel Jets Using Various Discretization Schemes and Turbulence-Model
Variations," ASME JOURNAL OF FLUIDS ENGINEERING, Vol. 103, 1981, pp.
significant pressure gradients which have not been accounted 352-360.
for in the boundary layer solution. The abrupt variation near 15 Pollard, A., and Siu, A. L., "The Calculation of Some Laminar Flows
point D results from the high viscosity near the stagnation Using Various Discretization Schemes," Comp. Meths. Appl. Mech. Eng., Vol.
point. 35, 1982, pp. 293-313.
16 Rhode, D. L., Lilley, D. G., and McLaughlin, D. K., "Mean Flowfields in
Axisymmetric Combustor Geometries with Swirl," AIAA Journal, Vol. 21, No.
Conclusion 4, 1983, pp. 593-600.
17 Vu, Bach T., and Gouldin, F. C , "Flow Measurements in a Model Swirl
A new and relatively simple method was developed for Combustor," AIAA Journal, Vol. 20, No. 5, 1982, pp. 642-651.
alleviating the convergence difficulty which arises upon im- 18 Morrison, G. L., Rhode, D. L., Cogan, K. C , Chi, D., "Labyrinth Seals
plementation of the QUICK scheme. The required CPU time for Incompressible Flow," Final Report for NASA Contract NAS8-34536,
November 1983.
and number of numerical iterations, as a percentage of that re- 19 Jerie, J., "Flow Through Straight-Through Labyrinth Seals," Pro-
quired by the previous Hybrid scheme, are shown to compare ceedings of the Seventh International Congress of Applied Mechanics, 1948,
favorably with a previously recommended approach. Due to Vol. 2, pp. 70-82.

Journal of Fluids Engineering MARCH 1986, Vol. 108/25

Downloaded From: http://fluidsengineering.asmedigitalcollection.asme.org/ on 01/20/2016 Terms of Use: http://www.asme.org/about-asme/terms-of-use