High-Pressure Annular Seal Leakage and Rotordynamics with Application to Turbomachinery

The paper presents theoretical investigations that deal with the steady state mass flow rate, labyrinth seal-influence coefficients, and unbalance response of a cryogenic turbopump. A numerical method for calculating the leakage flow rates, stiffness, and damping coefficients of liquid hydrogen– lubricated seals is presented. Rotation and seal surface roughness effects are also included for leakage and dynamic force calculation. A finite element solution of a Reynolds equation is obtained for steady state as well as dynamic pressure distributions produced by a small amplitude whirl of the seal journal center (a first order perturbation solution). Film rupture is taken into account with the boundary to the ruptured film zone determined by an iterative procedure. A brief procedure is presented for modeling of the rotor/seals system. This work extends the previous theory for plain centered seals to large eccentricities using a perturbation analysis, which includes surface roughness effects. Investigations have been carried out for seal breadth-to-diameter ratio (BDR) of 0.25 to 0.5. The performance of a convergent and three-stepped seal is estimated. Numerical data for these seals for an eccentricity ratio up to 0.8 has been obtained.

employed in the space shuttle main engine (SSME) high-pressure fuel turbo pump (HPFTP) to improve stability.Seals for liquids with low kinematic viscosity, such as steam, hydraulic oil, and liquid hydrogen often have Reynolds number well above 2000.They tend to produce large seal stiffness and damping coefficients that can be beneficial to rotor stability.Hydrodynamic forces generated in these seals can drive a machine unstable or contribute strongly to the stability of another machine.Furthermore, the fluid flow dynamics in a seal can have dramatic effect on vibration of rotating machinery.Previous reviews of seal-rotordynamics indicate that large subsynchronous vibrations were encountered in tests of the space shuttle main cryogenic engine high-pressure fuel turbo pump (HPFTP) with labyrinth seals and the original shaft bearing supports (Black, 1969;Black and Jenssen, 1971;Allaire, Lee, and Gunter, 1978).Literature review shows several see-through and step labyrinth seal flow analyses involved simple correlation but no rotational effects on seal flow.Yamada (1962) experimentally measured the friction factor for axial Reynolds number up to 10 4 and gave the formula for laminar and turbulent flow.Those equations relied on Blasius-type friction factors.Koeing and Bowly (1972) introduced the computer as a design tool in their article.They provide a computer program based on Egli's (1935) method, which pinpoint the important design parameters for various labyrinth configurations.Further, the kinetic energy carry over factor (or equivalent) is typically not well known for new designs (Rhode and Hibbs, 1993).
On the other hand, seal-rotordynamic coefficients are required for correlation of high-speed test measurements.These coefficients depend on a host of non-dimensional parameters, like length-to-diameter ratio (b/d), pressure ratio (P u /P d ), aspect ratio, H g /c, and c/r, etc. Black, 1969;Black and Jenssen, 1971 first explained the influence of seal forces on the rotordynamic behavior of pumps.Allaire, Lee, and Gunter (1978) extended previous theories for short plain centered seals to large eccentricities using a perturbation analysis.However, the analysis was performed in a stationary reference frame that ignores rotational effects.Childs (1981) traces the development of 404 D. R. GOWDA AND B. S. PRABHU theoretical treatment of the dynamic behavior of high-pressure annular inter-stage seals.In 1982, Childs and Dressman tested a three-step interstage seal and reported lower stiffness coefficients by about 20%.On the basis of their results obtained for a three-step seal, and neglecting the added-mass and crosscoupled coefficients, they recommend for reduced modeling seals in rotordynamics analysis of pumps.Childs (1983Childs ( a, 1983b) ) predicted very good results for a finite length seal.They provided a comprehensive program for the analysis of turbulent annular seals based on the Hirs (1973) turbulence bulk-flow model and included inlet swirl velocity effects.Journal bearings are normally cavitated.They are quite nonlinear, yielding rotordynamic coefficients which are strongly dependent on the static eccentricity ratio.Due to the absence of cavitation in seals, the orbit is assumed to be circular, and they may remain linear out static eccentricity ratio of about 0.5 (Childs, 1993).Hence, the coefficients of equation of motion for small motion about a centered position are generally satisfactory for the analysis of pump seal-rotordynamics.In cryogenic turbopumps, seals are associated with high pressures and flow rates and rotor speed up to 45000 rpm.Navier Stokes momentum equations, the continuity equation, and the turbulence model (Simon and Frene, 1989) govern hybrid behavior of these seals.
In a real turbopump seal, the surface roughness varies from 0.3 to 0.4 microns and having a slope of 10 to 12 degrees.This is an important parameter to be considered.Any isotropic surface has a two-dimesional roughness, and widely found in finemachined surfaces.Christensen and Tonder (1973) and Ram Turaga, Sekhar, and Majumdar (1998) gave the details of such roughness structures.They showed the different trend in comparison with respect to longitudinal and transverse roughness.In the sealing area fluid film has been divided into triangular elements.The fluidity matrices for three-node are evaluated using the procedure given in Reddi (1969) and Booker and Huebner (1972).The finite element formulation used in the present study is similar to one used by Kilt and Lund (1986) for a smooth seal, with necessary modifications.Allaire, Kokur, and Nicholas (1984) developed a pressure parameter method for the finite-element solution for Reynolds types equation.They have considered two relative surface roughness values of 0.001 and 0.01.The computer program developed by Chen and Jackson (1985) is more generalized and the output are in good agreement.Their computer code can also be used to obtain results of a three-step HPFTP-seal design that was modeled by Allaire Lee, and Gunter (1978).Childs and Kim (1985) performed an experimental test for a damper seal with different stator roughness patterns.Chen and Jackson (1987) developed a generalized theory for eccentric and misalignment effects in high-pressure annular seals.Nelson andNguyen (1988a, 1988b) extended the perturbation technique to handle two-phase flows in seals and analyzed the rotordynamic coefficients of incompressible flows with a relative rough surface.Scharrer and Nuez (1989) presented an analysis for an incompressible flow in a centered annual seal with a deformed surface.A redesign of the SSME HPFTP inter-stage seal as separate pieces from the diffuser has brought into attention that, due to assembly and operating interferences, the resulting seal clearance is wavy in the axial direction instead of straight or tapered.Scharrer (1991) advanced a model for fixed-tilt annular seals and presented force coefficients for cryogenic seal application.San Andres (1991) has developed a CFD solution for plain seals using a bulk flow model with Moody-fictionfactor.He analyzed the dynamic forces of turbulent flow annular high-pressure wavy seals.San Andres solution compares well with the predictions by Nelson andNguyen (1988a, 1988b) over their calculated eccentricity range.As mentioned by Childs (1993), theoretical calculations with a code developed by San Andres (1991) show a fairly good agreement, except for C xy , which is underestimated at high eccentricity ratios.However, San Andres' code solves two-dimensional, bulk flow, incompressible, Navier Stokes equations to predict the static and rotor dynamic characteristics of bearings and seals, in centered and eccentric positions.Rame Gowda, Sekhar, and Prabhu (1999) studied physical parameters of a straight through seal by using a finite-element method (FEM).Their work consists of a combination of numerical and computational investigations.They compared labyrinth seal leakage flow parameter with measured values of (Prasad, Manavalan, and Nanjundarao 1997), and predicted using finite volume method (FLUENT).The results were in good agreement by about 10% with predicted values.The dynamic coefficients are calculated with the empirical relation given by Kramer (1993).
The present work extends the theoretical treatment for plain and labyrinth seals to larger eccentricities and includes surface roughness effects.The continuity equation and momentum equations are solved using a perturbation approach.An analysis of above cited turbopump is then used as an example to illustrate the application of the seal for rotating machinery stabilization.One can use the "Stochastic finite-element" model for isotropic surfaces of annular roughened seals.Using this technique, it is possible to study the spatially randomly varying quantities of interest, such as variance and covariance of the rough surfaces, in addition to the mean clearance values.A schematic view of a rotor/seal eccentric position of seal geometry is shown in Figure 1.

FIGURE 1
Coordinate system for seal analysis.

SEAL ANALYSIS
Two-dimensional Reynolds equation governs the pressure distribution of fluid in a seal can be written as, Equation ( 1) when non-dimensionalized, the governing equation for the pressure distribution of an incompressible fluid under dynamic condition can be written in θ − z coordinates in the form: For a seal journal shown above, the following transformations are made in Equations ( 1) and ( 2) as: x = r θ and U = r ω.The turbulence Reynolds-type equation that governs flow field in the clearance space of annular seal, which has been derived from the Equation ( 1) for the purpose of flow continuity as Allaire et al. (1984): where With reference to the straight-through seals, schematic (Figure 1), the terms k x and k y includes the turbulent correction factors, such that G θ and G z are similar but not identical.Above turbulent correction factors are derived from the Reynolds equation for an incompressible fluid.
The rotor rotates at constant angular speed ω about its axis.It is assumed that the rotor undergoes a steady whirl in an elliptic orbit with a frequency ω p about its mean steady state position.The governing equation for a seal gives the fluid pressure distribution, as a function of ε and φ.Assuming that the rotor whirls about its mean steady state position given by ε 0 and φ 0 , the first order perturbations (which are generally valid for small displacements) of the pressure and the local film thickness (with no misalignment) can be expressed as, where The steady state film thickness is given by h0 = 1 + ε 0 cos θ, and H θ = 1 + ε 0 cos θ + H g /c for grooved surfaces.Substituting Equation ( 4) to (5) in Equation ( 2) and ( 3) and collecting the zeroth and first order terms for ε 1 and ε 0 φ 1 , a 3-set of equations were obtained for a seal.Using the finite element method, these equations can be solved.Static and dynamic pressures P0 , P1 , and P2 are obtained, which are given in Appendix A.

FINITE-ELEMENT FORMULATION
The detailed finite element formulation for the stochastic roughness modeling is as follows: Stationary values are found if: ∂ I ( pi )/∂ pi = 0.In the interior of an element the pressure distribution p i (θ, z) is expressed in terms of known polynomials and unknown parametric values: where n is the number of element nodes and the index j indicates the node number and the equivalent for stationary nodal values are ∂ I ( pi )/∂ pij = 0. Substituting Equation ( 6) into Equation A7.1 to A7.3 and performing differentiation on these equations yields static pressure P0 , and perturbed (dynamic) pressures pi = p2 , which is the imaginary component.These three differential equations can be reduced to the form, The boundary conditions are as follow: P0 = 0, when θ = 0, 2π and other two sides the same holds good when P0 < 0, upstream and downstream pressures are specified.P1 and P2 are made equal to zero at θ and at nodes, when P0 was equal to zero.

Isotropic Roughness Parameters for a Seal
On perturbing equation (see Appendix A7.1 and A7.3), from mean-steady-state position and repeating the procedure as described, one can obtain another set of equations considering seal roughness both (stator and rotor) as: The above three sets of equations give P0 , P1 , and P2 , respectively, and these values are used for calculation of eight linearized dynamic coefficients.
For the roughened surfaces, total mean film thickness is calculated using the expression, where r denotes the radial clearance and ε the eccentricity ratio.In general, in Equation ( 11) the right hand term is a function of both space variables θ and z, its random character being indicated by the random variable ζ .Its precise form depends upon the structure of the surface roughness, the E ( ) operator is given by the relation, is the density function for the random variable h s .It has been found convenient to use a polynomial approximation rather than the Gaussian distribution itself.Such a function is This function terminates at c = ±3σ where σ is the standard deviation.When h ≈ 3σ (with in the hydrodynamic limit i.e., h > 3σ ), these may influence the seal performance significantly.
Once finite element formulation is done, the final expressions neglecting the squeeze term for the deterministic case using 3-node linear triangular elements can be written as,

Stochastic Finite-Element Formulation
Taking the expected values from Equations ( 13) to ( 16) for the rough surface seal, one can get the following equation, Q is the same as in the deterministic case.H is a stochastic variable and u is deterministic variable, which are approximated by interpolation functions.H has Gaussian or an approximation to Gaussian distribution.The pressure boundary conditions for solving Equation ( 17) are the same as those for deterministic finite element method and the flow is considered as linear normal outward flow.Seal leakage flow m can be calculated from

Unidirectional Load
With reference to Figure 1, the hydrodynamic forces and the steady-state characteristics are calculated from: Unidirectional applied load is assumed to be constant in both magnitude and direction and equal to the steady state load carrying capacity, i.e., W = W0 .The equations of motion under these conditions can be written as: Equations ( 21) and ( 22) are second-order differential equations in ε and φ.These are solved by using the fourth-order Runge-Kutta method with the initial boundary condition ε, φ = 0. From this we get the eccentricity ratio, attitude angle, and their derivatives.These are used in the solution of the Reynolds equation for the next time step to obtain the pressure distribution and, hence,

Stiffness and Damping Coefficients
The components of the dynamic seal load along the line of centers and perpendicular to the lines of centers corresponding to the perturbed pressure ε 1 p 1 e iωpt can be written as Again, since the seal journal executes small harmonic oscillation about its steady state position in an elliptic orbit, dynamic load carrying capacity can be expressed as a spring force and a isoviscous damping force as given below: Using Equation ( 24) in Equation ( 23), and integrating, the following non-dimensional components of stiffness and damping coefficients result: K rr , K rφ , Drr , Drφ .Similarly, considering dynamic displacements of the seal journal along the φ direction, the other coefficients K r φ , K φφ , Dφr , Dφφ can be written by analogy.In K i j the first suffix denotes the direction of the force and the second suffix denotes the direction of displacement and for Dij the first and second suffix denote force and velocity directions, respectively.As the dynamic pressure distributions P1 and P2 have been obtained by finite element method, the above eight coefficients can be obtained by integration.

SPACE SHUTTLE APPLICATION-CRYOGENIC PUMP
The seal geometry used in present analyses is taken from Allaire, Lee, and Gunter (1978) and Chen and Jackson (1987) but with a reduced length-to-diameter ratio of 2:1. Figure 2(

RESULTS AND DISCUSSION
The numerical computation has been carried out for four seals.A fine grid of 40 × 20 (1600 elements) was used with a three-node triangular element.The dynamic characteristics of the seal viz., stiffness, damping, and mass flow rates, have been obtained by using the finite-element perturbation method.The results presented in Table I shows the comparison of the seal-rotordynamic coefficients and mass flow rate for centered seal.The results from present work corroborate well with the results published by Allaire, Lee, and Gunter (1978), and Chen and Jackson (1987).The dynamic characteristics for a three-step separated labyrinth seal result are not satisfactory.But the other three labyrinth seals corroborate well.It has also been found that the combined labyrinth seal gives lower leakage when compared with other three seals.The combined seal has a lower flow rate, about 42% less as against the three-step labyrinth seal.
Figure 3 shows the variation in leakage versus pressure differences.As the pressure drop increases, the leakage also increases sharply.From the present results, it is also clear that smaller clearance provide less leakage rates, which is necessary from a labyrinth seal working point of view.Since the present study uses a positive clearance of 0.0556 mm, lower leakage rates are

FIGURE 6
Shows direct damping coefficients versus eccentricity ratio.

FIGURE 7
Shows cross-coupled damping coefficients versus eccentricity ratio.
predicted as compared to Allaire, Lee, and Gunter (1978), and Chen and Jackson (1987), which use slightly larger clearances of 0.07 mm and 0.65 mm.
Figure 4 shows variation of the direct stiffness coefficient versus eccentricity ratio.It is observed that the direct stiffness coefficient values of four seals do not vary much with an eccentricity ratio up to 0.6.In the case of present seal models, stiffness increases moderately with increase in eccentricity ratio.Figure 5 shows cross-coupled stiffness coefficients with eccentricity ratio.The coefficient values are well maintained up to 0.4 eccentricity ratio and thereafter there is a sharp increase.From Figure 6 it is observed that the variation of direct damping coefficients with respect to eccentricity ratio up to 0.5 is as expected.Beyond 0.5 it increases sharply.The three-step separated and combined seals show a sharp increase in damping coefficient with increase in eccentricity ratio.Figure 7 shows cross-coupled damping coefficients also increase with increase in eccentricity ratio from 0.2, for both separated and combined seals.The three-step seals of Allaire, Lee, and Gunter (1978), and Chen and Jackson (1987) show similar trends up to 0.7 eccentricity ratio.Figure 8 shows the variation of the whirl frequency ratio with static eccentricity ratio.Allaire's rough seal maintains whirl frequency rations almost equal to 0.99 up to eccentricity ratio of 0.6.This shows that it is operating in synchronized condition in that range.Separated seal maintains steadily up to 0.5 whirl frequency ratio

FIGURE 8
Shows whirl frequency ratio versus eccentricity ratio.

FIGURE 9
Shows unbalance response at different speeds without labyrinth seal.
with a negative sign, whereas combined seal sub-synchronizes up to 0.6 whirl frequency ratio and sharply decreases up to 0.8 eccentricity ratio.
The numerical example has been carried out for a two-stage turbine rotor with a rotor-labyrinth seal using ROTORDYNA, a code developed by second author.Using the combined seal coefficients, the unbalance response of a labyrinth seal for steadystate condition for various speeds is plotted in Figure 9 and Figure 10. Figure 9 shows the steady-state rotor response for different speeds without a seal.The response in both x and yaxes are shown for the above two cases.It is observed that the unbalance response amplitude is lower with a seal as shown in Figure 10.Using the above code, one can make a study of the total unbalance response of a turbopump model, which can have a number of labyrinth seals at different locations.

CONCLUSIONS
A combined numerical-computational method was developed to account for grooved-surface-roughness treatments in the rotor and stator of turbulent annular seals.A finite-element analysis has been carried out for labyrinth seals with low axial flow rates.The stiffness and damping of the seals have been cal-

FIGURE 10
Shows unbalance response at different speeds with labyrinth seal.
culated for wide range of eccentricity ratios.A combined seal has been used as a numerical model considering the effects of a uniformaly-distributed tooth height on rotor surface, which is differed in a step height by 1 mm.Four different seal configurations of the same length dynamic characteristics have been evaluated.This work extends the concentric seal theories to cover the seal leakage and dynamic coefficients with static eccentricity ratios up to 0.8.Leakage results are consistent with the model for both the centered and fully-eccentric operations for all seals.Contrary to the seal eccentricity effect, the seal leakage increases as the static eccentricity ratio increases.The seal rotor-dynamic characteristic changes significantly for a high eccentricity ratio, especially when the perturbation direction coincides with minimum clearance.For eccentric seals, the relationship between the direct damping coefficients has been preserved up to 0.3 static eccentricity ratio.Unbalance response shows that the amplitude is lower when the seal coefficients are included.

FIGURE 2
FIGURE 2 (a) and (b) show a sectional seal configuration.
FIGURE 4Shows direct stiffness coefficients versus eccentricity ratio.

TABLE I
Rotodynamic coefficients and mass flow rate of different seals The force components are then used in the solution of the equation of motion.The procedure of solving the Reynolds equation and the equation of motion is repeated until it completes the number of iterations.