3D numerical simulation of pump cavitating behavior

The quasi-steady cavitating behavior of three pumps was investigated by 3D unsteady viscous computations. The numerical model is based on the commercial code FINE/TURBO TM , which was adapted to take into account the cavitation phenomenon. The resolution resorts to a time-marching algorithm initially devoted to compressible flows. A low-speed preconditioner is applied to treat low Mach number flows. The vaporization and condensation processes are controlled by a barotropic state law that links the void ratio evolution to the pressure variations. A radial pump, a centrifugal pump, and a turbopump inducer were calculated and the cavitating behaviors obtained by the computations were compared to experimental measurements and visualizations. A reliable agreement is obtained for the two pumps concerning both the head drop charts and the extension of the vapor structures. A qualitative good agreement with experiments is also observed in the case of the turbopump inducer. The accuracy of the numerical model is discussed for the three geometries. These simulations are a first attempt to simulate the complete 3D cavitating flows in turbomachinery. Results are promising, since the quasi-steady behaviors of the pumps in cavitating condition are found quantitatively close to the experimental ones. A continuing effort is pursued to improve the prediction accuracy, and to simulate unsteady effects observed in experiments, as, for example, rotating cavitation


INTRODUCTION
To increase the performance of industrial hydraulic components, pumps should operate with high rotation speed and low inlet pressure.It results in cavitating conditions that can lead to severe disturbances, such as substantial performance losses, strong unsteady forces acting on the pump components, noise and erosion.
In rocket engine turbopumps, for example, the inducer stage is designed to allow a stable cavitation behavior at the lowest available inlet pressure.In centrifugal pumps, strong pressure gradients in the cavity wakes and associated unsteady cavitation behavior are involved in erosion phenomena, whose intensity directly acts on the machine life duration.
The general understanding of the unsteady behavior of cavitation is therefore of primary importance for the design and operation purposes.More particularly, the investigation of mechanisms that are responsible for a pump elevation drop in cavitating conditions is a necessary step to improve inducers design.These considerations have motivated the present study, which focuses on the quasi-steady effects of cavitation on the behavior of two centrifugal pumps and a turbopump inducer.
A 3D numerical model of unsteady cavitation is proposed with the final objective to predict cavitating flows in rocket engine turbopump inducers.Cavitation phenomenon was simulated by strong variations of the density of a single fluid.That model was implemented in the commercial code FINE/TURBO TM developed by NUMECA International.The numerical code solves the RANS equations, associated with a barotropic state law managing the fluid variable density.Models based on that simple law gave many promising results in 2D simulations of steady and unsteady cavitating flows ( [1] to [6]).Computations were applied in previous works to Venturi type sections designed to simulate an inducer blade suction side, and to inducer blade sections placed in a cavitation tunnel: the proposed model seemed to well describe the cavitation behavior in 2D geometries [7], [8].A first application to a centrifugal pump geometry also gave promising results [9].
The model is applied in the present paper to calculate 3D unsteady cavitating flows in a radial pump, in a centrifugal pump and in a rocket engine turbopump inducer.The influence of quasi-steady cavitation behavior on the pump characteristics and on the final head drop is simulated by transient calculations.Experimental and numerical results concerning pumps characteristics and performance breakdown are compared at different flow rates.Appearing types of cavitation and the spatial distribution of vapor structures inside the impellers are also analyzed for different operating points at nominal flow rate.
The investigation of these three different pump behaviors enhances the description of the performance of the numerical model, and points out its limitations at the present time.: τ/τ 0 τ 0 : reference cavitation parameter -

PHYSICAL AND NUMERICAL MODEL
The mean features of the applied physical and numerical models, as well as the calculation procedure are summarized in the present paper.More details are given by Coutier-Delgosha et al. [7].

1) Physical model
Cavitating flows are described by a single fluid model, based on previous numerical and physical work developed in LEGI by Delannoy and Kueny [1], Reboud et al. [4].This fluid is characterized by a density ρ that varies in the computational domain: when the density in a cell equals the liquid one (ρ l ), the whole cell is occupied by liquid, and if it equals the vapor one (ρ v ), the cell is full of vapor.Between these two extreme values, a liquid/vapor mixture, still considered as one single fluid, occupies the cell.The void fraction α = (ρ -ρ l ) / (ρ v -ρ l ) can thus be defined as the local ratio of vapor contained in this homogeneous mixture.
Velocities are assumed to be locally the same for liquid and for vapor.An empirical state law is used to manage the mass fluxes resulting from vaporization and condensation processes.The barotropic law links the density to the local static pressure ρ(P).When the pressure is higher or lower than vapor pressure, the fluid is supposed to be purely liquid or purely vapor, according to the incompressible state law ρ = ρ l or ρ = ρ v .The two fluid states are joined smoothly in the vaporpressure neighborhood.It results in the evolution law presented in Figure 1, characterized mainly by its maximum slope 1/c min 2 , where c 2 = ∂P/∂ρ.c min can thus be interpreted as the minimum speed of sound in the mixture.

2) Numerical model
The numerical model of cavitating flows based on this physical description is an extension of the commercial code FINE/TURBO TM developed by NUMECA International.FINE/TURBO TM solves the time dependant Reynolds-Averaged Navier-Stokes equations on three-dimensional structured meshes.Time accurate resolutions of the equations use the dual time stepping approach.Pseudo-time derivative terms are added to march the solution towards convergence at each physical time step.The range of application is extended to low-compressible or incompressible flows by introducing a preconditioning matrix [10].
The discretization is based on a finite volume approach.Convection terms are treated by a second order central scheme associated with artificial 2 nd order and 4 th order dissipation terms.The pseudo-time integration is made by a four-step Runge-Kutta procedure.The physical time-derivative terms are discretized with a second order backward difference scheme.The code resorts to a multigrid strategy to accelerate the convergence, associated with a local time stepping and an implicit residual smoothing.
The numerical model was adapted to treat the cavitation process by Coutier-Delgosha et al. [7].The key point of this adaptation is the modification of the fluid state law.Applied barotropic law implies the simultaneous treatment of two different cases: the fluid is highly compressible in the liquid/vapor mixture (the Mach number can be as high as 4 or 5) and is almost incompressible in the pure vapor or pure liquid areas.So the main difficulty consisted in managing these two different states of the fluid, without creating any spurious discontinuity in the flow field.Besides, cavitation consists in a very sharp and very rapid process.The density variations in time and space are smoothed to avoid numerical instabilities.A theoretical study of the influence of the preconditionner on the stability in the two-phase regions was performed.An adaptation of the preconditionning coefficients was found to be effective in 2D configurations, and work is in progress to apply it to calculations of rotating machineries.
We use for the simulations presented in this paper a Baldwin-Lomax turbulence model.More calculations considering κ-ε turbulence models are in progress to improve physical analyses.An investigation of the respective efficiency of different turbulence models associated with the single fluid modelization of cavitation is also proposed in [11].

3) Initial transient treatment
First of all, a steady step is carried out, with a pseudo vapor pressure low enough to ensure non-cavitating conditions in the whole computational domain.Then, the unsteady calculation is started: the NPSH is slowly reduced by increasing smoothly the pseudo vapor pressure at each new time step.Vapor structures spontaneously appear and grow during this process, in the regions of low static pressure.In the presented computations, the NPSH is decreased until the pump head drop, due to highly cavitating conditions, is reached.The physical time-step is high enough to only simulate the quasicstatic effects of the NPSH decrease.So the unsteady effects are not modeled in the computations presented hereafter.For a given pseudo vapor pressure value, the inlet pressure does not fluctuate so far the convergence is reached, and it is directly used to compute the NPSH values reported on the charts.

1) Geometry
This impeller geometry has been chosen by the Laboratory of Turbomachinery and Fluid Power of the Darmstadt University of Technology (TUD) [9] [12].Measuring techniques previously developed for a single hydrofoil in a test section were adapted to the pump [6].The impeller has 5 single-curved blades with two different radii at inlet and outlet.Hub and shroud are parallel to get an almost 2-dimensional blade-to-blade channel with constant width.
A radial-symmetric housing (Figure 2) is used to obtain almost constant conditions on the impeller outlet, so that comparable cavitation conditions are simulated in each channel, without any influence of a volute casing.
Nominal conditions are a rotational speed of 36 Hz and a flow-rate Q = 210 m 3 /h.Specific speed of the runner is n s = N Q 1/2 / H 3/4 = 20 (europ.value: N in rpm, Q in m 3 /s and H in m), the outer runner diameter equals 278 mm.Cavitation conditions are defined by NPSH value based on the upstream total pressure, water vapor pressure and density at ambient temperature.

2) Mesh and boundary conditions
The model was applied to the presented radial pump geometry.A 300.000 cells multi-blocks mesh of one single blade-to-blade channel was used (Figure 3).Boundary conditions applied for the simulations are the following: Velocity is imposed at the inlet of the suction pipe.Laws of the wall are imposed along solid boundaries.The relative motion between the inlet pipe walls and the runner is taken into account.On the other hand, the outlet housing shape is not described and the parallel walls are treated as hub and shroud extensions up to the outlet, at 1.5 times the runner outer radius, where a uniform static pressure is imposed.Periodic conditions are applied at the lateral boundaries of the blade-to-blade channels.
The c min value in the state law was calibrated by previous 2D studies reported in [4], and it is set to the reference value 2.25 m/s in this computation.

3) Calculation
Numerical simulation of the cavitation characteristics of the pump was performed at different flow rates.The shape of the cavitating structures is first compared to the experimental visualizations at nominal point of operation.When the NPSH decreases in the calculation, attached cavitation sheets grow both on the suction side and on the pressure side, as observed experimentally.Moreover, vapor structures appear at the inlet radius of the shroud.This cavitation caused by the local curvature of the streamlines is fully consistent with the observations reported by Hofmann et al. [9].The experimental visualization is compared to the result of the computation to enhance this reliable agreement (Figure 4).The cavitation number was adjusted to give the same global size of the cavitation structures than in the experiment: the numerical result then corresponds to a NPSH about 10% lower than the experimental one (7 m, instead of 8 m).The three cavitation areas are correctly simulated by the code: attached cavity on the suction side (A), extension of cavitation on the shroud along the blade (B), and cavitating flow on the inlet radius of the shroud (C).In the computation, attached sheet cavity (A) and its extension on the shroud (B) belong to the same vapor structure, while they look like two separated regions in the experiment.A more global comparison between experiments and computations is presented in Figure 5, by drawing the head drop curves for three different flow rates: namely 0.8Qn, Qn and 1.08Qn.In the presented calculations, the performance drop is represented first by a smooth decrease of the pump head, as a function of the NPSH.The final head-drop is only partially simulated because the computation rapidly becomes unstable and stops (this is more particularly the case at 1.08Qn).Our upstream boundary condition consists in imposing in a strict manner the mass flow rate passing in the pump.Because the coupling between the pump and the hydraulic loop is not taken into account, the effects of the cavitation blockage on the flow rate are neglected and the headdrop is less progressive than in the experiments.
Results obtained from this first simulation of the pump cavitation behavior are promising: the head drop is predicted with a good homogeneity with respect to the three flow rates.The NPSH values obtained for the 3% and 10% head drop are globally overestimated with respect to experimental values (of about 1m for the 10% head drop, and 1.5 to 3 m for the 3% head drop).

1) Geometry and mesh
The SHF pump is a shrouded seven blades centrifugal runner (figure 6).The scale and the calculation configuration are as close as possible to the experimental conditions reported by Combes and Archer [13], in order to simplify comparisons in cavitating conditions.The outlet radius is 330 mm, the inlet pipe length is 200 mm, the rotation speed equals 3000 rpm, and the nominal flow rate is 0.157 m 3 /s.Specific speed of the runner is n s ≈ 35.The experimental outlet has been slightly simplified for the numerical simulations: a rotating vaneless diffuser is used at the outlet, and the fluid leakage on the back of the shroud is not taken into account.Both hub and shroud are rotating, axial velocities are imposed at the inlet, and a static pressure is fixed at the diffuser outlet (figure 7).

Figure 7: Calculation conditions (meridional view)
A two-block mesh is applied on a single blade-to-blade channel.The influence of the mesh size was evaluated by noncavitating steady computations performed with meshes made respectively of 75000 cells, 250000 cells, 500000 cells and 10 6 cells.The effect on the total pressure head and on velocity profiles through a blade-to-blade channel was found negligible in the three last cases, and a 250000 cells mesh is applied for the non-cavitating and cavitating computations presented in this paper.Figure 8 presents a 3D view showing principally the hub treatment.For this geometry, the head drops are more completely simulated than for the previous one.Numerical instabilities that occur in highly cavitating conditions were mastered, so that 10% to 15% of head decrease could be obtained before the computation divergence.This improvement was reached mainly by eliminating the initial transient effects associated with the start of the unsteady computation, which probably interfere with the cavitation process and generate during the head drop some spurious numerical instability.The increase of the second order artificial dissipation level in the two-phase regions also removed some difficulties.
Results obtained from cavitating calculation of the pump are promising: the elevation decrease due to cavitation development, as well as the final head drop, are correctly predicted for both flow rates.This reliable agreement concerning the head drop prediction confirms the ability of the numerical model, already observed in the case of the purely radial pump presented in the previous section, to predict the quasi-steady behavior of centrifugal pumps in cavitating conditions.

3) Flow field visualizations at nominal flow rate
The head drop chart H (NPSH) obtained by the calculation at nominal flow rate is drawn on Figure 10. Figure 11 illustrates the general development of cavitation in the pump.The hub is colored in grey, the blades are represented in dark red, and the two-phase areas appear in yellow (the surface drawn corresponds to a 5% void ratio).It can be observed that cavitation mainly consists in a cavity that grows on the blades suction side, while no vapor is present on the pressure side until the three last configurations.

Figure 11: Cavity evolution on the blade suction side Vaporized areas in yellow, hub in grey, blades in dark red
The suction side cavity grows progressively: the cavitation sheet first appears near the hub, where its extension remains always larger than around the shroud (Conf. 1 and 2); then, the cavity length increases slowly until it reaches the throat (conf. 2 to 4) and after that it almost stabilizes (conf.5 and 6).In the same time, the detachment point of the cavitation sheet moves downstream on the blade.
Another cavity suddenly appears on the pressure side in configuration 4 (not visible here) in the vicinity of the shroud, and grows rapidly (conf.5 and 6) inside the blade-to-blade channels.This phenomenon is strongly correlated with the final head drop reported in Figure 10.
This behavior is consistent with experimental observations, which indicate at nominal flow rate a two-step cavitation evolution: as a matter of fact, vapor first appears near the hub on the blades suction side.It progressively grows towards the throat, but it does not enter the blade-to-blade channels.This first step only slightly affects the pump elevation.The second step consists in a sudden development of cavitation on the pressure side, which fills the blade-to-blade channels and induces the final head drop.
The good qualitative agreement obtained with the numerical simulation at nominal flow rate has to be associated with the accurate prediction of the head drop chart.
A comparison between the results of the computation and experimental visualizations performed at EDF R&D Division is proposed in Figure 12.The numerical results are presented with the same viewpoint as in experiments.The order of the figures corresponds to a decrease of the NPSH, from the vapor onset to the pump head drop.The experimental views are extracted from video, which leads to an uncertainty on the experimental NPSH of ± 0.5 m.
A good agreement with experiments can be observed concerning the shape of the suction side cavitation sheet, during the whole NPSH decrease.The cavity length and its position in the radial direction are systematically well predicted by the numerical model.
The main discrepancy appears at the leading edge, where the calculation shows a progressive shift of the cavity detachment point, when the NPSH decreases.This behavior is not observed in the experimental visualizations.

1) Geometry, mesh, and boundary conditions
The third geometry is a turbopump inducer, i.e. a fourblade axial runner especially designed to operate in cavitating conditions (Figure 13).

Figure 13: Inducer geometry
A 500.000 cells multi-blocks mesh of one single blade-toblade channel was used (Figure 14a).Only one channel is treated, and periodicity conditions are applied to the frontiers with the adjacent channels.Non-matching cells are used along these periodic boundaries, to avoid an excessive skewness of the grid (Figure 14b).A special attention is also paid to the mesh in the expected cavitating areas: grid is refined in the streamwise direction just after the blade leading edge, on the suction side, and also in the tip area, in the three directions.Boundary conditions applied for the simulations are the following: Velocity is imposed at the inlet of the suction pipe, about 2 diameters upstream from the runner, and the static pressure is imposed at the outlet.Laws of the wall are imposed along solid boundaries.The relative motion between the inlet pipe walls and the inducer is taken into account.The tip gap is also modelized.Periodic conditions are applied at the lateral boundaries of the blade-to-blade channels.
The calculation was performed in hydrogen.The minimum speed of sound c min in the two-phase mixture, which is necessary to determinate the maximum slope of the barotropic law, is calibrated to ensure the similarity between the calculations (in H 2 ) and the experiments (in water).The parameter c min /V ref is thus the same in both cases, which leads for the calculations to a value c min within the celerity range (30 m/s -100m/s) calculated with the Jakobsen's law [14].Figure 16 shows the development of cavitation corresponding to the eight operating points indicated on Figure 15.It illustrates the apparition and the growing of vapor/liquid areas on both faces of the blades, and finally the progressive filling of blade-to-blade channels by the vapor.Blades are colored in grey, the hub in blue, and the external shape of the two-phase areas (corresponding to a 5% void ratio) is represented in yellow.It can be observed that vapor first appears in the vicinity of the shroud, near the leading edge.This phenomenon slightly affects the inducer performance, which starts to decrease, even if only a little volume of vapor is present in the flow field.Vapor remains localized near the shroud and in the tip gap (config. 1 to 4) until τ' reaches 0.3.This behavior is consistent with the smooth decrease of the inducer head observed at the same time.
The development of cavitation on the blade suction side (config.5) directly induces a severe performance breakdown.This behavior can be explained by the important blockage generated immediately by the vapor in the upstream part of the blade-to-blade channels.It is a very sudden phenomenon (config.6 and 7), which leads to the presence of vapor structures characterized by a low void ratio, from hub to shroud and from pressure side to suction side at the channel inlet.Experimental visualisations are not yet available concerning this geometry, but the cavitating areas upstream from the inducer are consistent with behaviors observed previously in this type of machine.Calculations are in progress to investigate the influence of the inlet and outlet boundary conditions on the result.
3) Comparisons with experimental performance charts Numerical head drop charts are compared to the experimental ones at three flow rates, namely 0.8 Qn, Qn and 1.2 Qn (Figure 17).Details on the experiments can be found in [15].
The cavitation parameter τ' corresponding to the performance breakdown is in all cases over-predicted by the model: the gap between experimental and numerical results equals ∆τ' = 0.1.This discrepancy must be considered circumspectly.In the case of the previous centrifugal pumps, the order of magnitude of the cavitation parameter at head drop is notably higher, so that the same absolute error results in a much lower relative error (about 20% for the first pump, and less than 5% for the second one).Since an inducer has much better suction capacities (τ' at the performance breakdown is about 0.1), the order of magnitude of τ' is, in the present case, the same than the prediction error ∆τ'.As a matter of fact, in the case of an inducer, the simulation does not allow at the present time a sufficient accuracy for the useful prediction of the cavitation performance.The prediction of the head drop could be notably improved in the future by reducing the aspect ratio of the mesh in the two-phase areas, and by analyzing the influence of the artificial dissipation on the simulation of the performance breakdown.

CONCLUSION
The cavitating behaviors of a radial pump, a centrifugal pump, and an axial turbopump inducer were investigated with 3D unsteady calculations performed at various flow rates.The study of the quasi-steady development of cavitation was performed by drawing the head drop charts and visualizing the vapor areas for several operating points.The numerical simulations were compared for the three pumps to experimental results.
A good agreement was obtained in the case of the pure radial and the centrifugal pumps: head drops predictions are consistent with experimental measurements in most of the cases, and the two-phase structures observed during experiments are similar to the ones predicted by the numerical model.A good qualitative agreement was also obtained in the case of the turbopump inducer, but the very good suction capacity of this machine leads to a poor relative accuracy concerning the head drop location.More details and numerical tests concerning calculations performed with the radial 2D pump, the SHF pump, and the inducer can be found respectively in [9], [18], and [15].Work is in progress to improve the prediction for this type of machine, and to modelize the unsteady effects associated with cavitation, as already performed in 2D blade cascades [16], [17].

Figure 3 :
Figure 3: View of the mesh on shroud side of the pump.(The entire pump geometry is reconstructed by rotation of the single blade-to-blade channel)

Figure 5 :
Figure 5: Head-drop curves of the TUD pump Comparison at 0.8Qn, Qn and 1.08Qn

Figure 8 :
Figure 8: Reference mesh on hub side of the pump 2) Head drop charts Figure 9 presents the head drop curves H (NPSH) for two different flow rates: namely Qn and 1.3Qn.

Figure 9 :
Figure 9: Head drop charts at two flow rates

Figure 10 :6
Figure 10: Head drop chart at nominal flow rate Green line = onset of vapor in the flow field Squares = operating points represented on Fig. 11, 12.

Figure 12 :
Figure 12: Comparison experiment/calculation of the suction side cavity extension for several NPSH values SPATIAL TURBOPUMP INDUCER

Figure 14a :Figure 15 :
Figure 14a: View of the mesh on hub side of the inducer (The entire inducer geometry is reconstructed by rotation of a single blade-to-blade channel)

Figure 16 :
Figure 16: Development of the two-phase areas as the cavitation number decreases (corresponding to red squares reported on Figure 15)

Figure 17 :
Figure 17: Comparisons between numerical and experimental head drop charts