Unsteady simulation of cavitating flows in Venturi

9th International Conference on Hydrodynamics October 11-15, 2010 Shanghai, China

711

2010, 22(5), supplement :753-758 DOI: 10.1016/S1001-6058(10)60026-1

Unsteady simulation of cavitating flows in Venturi Eric Goncalves*, Jean Decaix, Regiane Fortes Patella LEGI, Grenoble-INP, France E-mail: [email protected] ABSTRACT: A compressible, multiphase, one-fluid

RANS solver was developed to study turbulent cavitating flows. The interaction between turbulence and two-phase structures is complex and not well known. This constitutes a critical point to accurately simulate unsteady behaviours of cavity sheets. In the present study, different turbulence transport-equation models are investigated. Numerical results are given for a Venturi geometry and comparisons are made with experimental data. KEY WORDS: Cavitation, RANS simulations, Homogeneous Model, Turbulence Model

1 INTRODUCTION The simulation of cavitating flows is a challenging problem both in modelling of the physics and in developping robust numerical methodologies. Such flows are characterized by important variations of the local Mach number, compressibility effects on turbulence, and involve thermodynamic phase transition. For the simulation of these flows, the numerical method has to handle accurately any Mach number. Moreover, the turbulence modelling plays a major role in the capture of unsteady behaviours. Cavitation sheets that appear on solid bodies are characterized by a closure region which always fluctuates with the existence of a re-enrant jet. This one is mainly composed of liquid which flows upstream along the solid surface. Reynolds-Averaged Navier-Stokes (RANS) models are frequently used to simulate such unsteady cavitation flows. One fundamental problem with this approach is that turbulence models are tuned by steady-state mean flow data. Moreover, the standard eddy-viscosity models based on the Boussinesq relation are known to over-product eddy-viscosity, which reduces the development of the re-entrant jet and two-phase structure shedding [1]. The limitation of the turbulent viscosity is therefore a determinant point to capture realistic cavitation sheets.

Different ways have been investigated to limit or to correct standard turbulence models. An arbitrary modification was proposed by Reboud to reduce the tubulent viscosity [1], and has successfully been used by different authors [2, 3]. Other corrections are based on the modelling of compressibility effects of the vapour/liquid mixture in the turbulence model. Correction terms proposed by Wilcox [4] in the case of compressible flows have been tested for unsteady periodic cavitating flows [2]. A sensitivity analysis of constants Cε1 and Cε2, which directly influence the production and dissipation of turbulence kinetic energy, was conducted for a k-ε model and a cavitating hydrofoil case [5]. Finally, a filter-based method was investigated [6] by which the sub-filter stresses are constructed directly using the filter size and the k-ε turbulence closure. In this paper, an in-house finite-volume codes solving the RANS equations is described with an homogeneous approach. Cavitation phenomenon is modelled by a barotropic liquid-vapour mixture equation of state (EOS). Various transport-equation turbulence models are tested with different corrections and eddy viscosity limiters. The turbulence model influence is discussed, based on comparisons between experimental and numerical results. 2 THE NUMERICAL CODE The numerical simulations were carried out using an in-house CFD code solving the one-fluid compressible RANS system for multi-domain structured meshes. It is based on a cell-centered finite-volume discretization. For the mean flow, the convective flux density vector on a cell face is computed with the Jameson scheme [7] in which the dispersive error is cancelled. The viscous terms are discretized by a second-order space-centered scheme. For the turbulence transport equations, the upwind Roe scheme was used to obtain a more robust method. The temporal integration is performed through a low-storage matrix-free implicit method. Each pure

9th International Conference on Hydrodynamics October 11-15, 2010 Shanghai, China

712

phases follow the stiffened gas EOS. In the pure liquid, a preconditioning method based on the Turkel approach [8] is used. As the preconditioning method does not preserve the consistency in time, the dual time stepping approach is used to perform unsteady flows. The numerical treatment of the boundary conditions is based on the use of the preconditioned characteristic relations. More details concerning the code are given in [9]. 2.1 The barotropic model The cavitation model is based on a sinusoidal barotropic law proposed by Delannoy and Kueny [10]. This law is characterized by its maximum slope 2 . The quantity cmin is an adjustable parameter of 1/ cmin the model, which can be interpreted as the minimum speed of sound in the mixture.

ρ=

ρL + ρV ρL −ρV 2

+

2

⎛ p− pvap 2 ⎞ sin⎜ 2 ⎟ (1) ⎝ cmin ρL − ρV ⎠

The influence of cmin was studied in previous works [9]. The value of cmin is set to 0.472 m/s. The vapour pressure Pvap is fixed to 2340 Pa. 2.2 Turbulence models Various two-equation turbulence models were used in the present study: the Menter SST k-ω model (KWSST) [11], the high Reynolds version of the Jones-Launder k-ε model (KE) [12], and the one-equation Spalart-Allmaras model (SA) [13]. Turbulence models always leads to the generation of stable cavities, because very strong turbulent eddy viscosity μt inside the cavity avoids the re-entrant jet formation which plays the major role on the instability of partial sheet cavity. As a remedy to such problem, one can use a turbulent viscosity limiter in the mixture area. An empirical reduction was proposed by Reboud [1] for unsteady simulations. In the present study, we propose to test and to compare different eddy viscosity limiters: the Reboud formulation, the Menter SST (Shear Stress Tensor) correction [11] and Durbin realizability constraints [14]. For the modeling of the flow close to the wall, a wall law approach is used.

Reboud introduced a function f(ρ) in the computation of the turbulent viscosity for the k-ε model:

μ t = f (ρ ) ⋅ C μ

ε

n

(2)

(ρ −ρ ) L

V

where n is a parameter fixed to 10, Cμ=0.09 and α is the void ratio. 2.4 Shear stress tensor correction The Menter correction [11] is based on the empirical Bradshaw's assumption which binds the shear stress to the turbulent kinetic energy for two-dimensional boundary layer: k 0 .3 k ) ν t = Min( , (3) ω 2 Ω F2 ( y ) where F2 is a blending function that tends to zero outside the boundary layer, and Ω is the vorticity. 2.5 Durbin correction Based on the realizability principle, a minimal correction was derived for two-equation turbulence models and was shown to cure the stagnation-point anomaly [14]. A weakly non-linear model was thus obtained with a Cμ coefficient function of the dimensionless mean strain rate: c Cμ = min(0.09, ) (4) S where S is the stress tensor and c a parameter varying between 0 and 1. In the present study, c is fixed to 0.3. 3 COMPUTATIONAL RESULTS 3.1 Experimental conditions The proposed numerical models were applied in a two-dimensional Venturi geometry, characterized by a divergence angle of 4o (Fig. 1). All experiments were led at the CREMHyG (Centre d'Essais de Machines Hydrauliques de Grenoble). The selected operating point is characterized by reference[15]: Uinlet = 10.8 m/s : the inlet velocity Pinlet = 35000 Pa: the pressure in the inlet section P −P σ inlet = inlet vap2 = 0.55 : the cavitation parameter in 1 ρV 2 inlet the inlet section Lref = 0.252m : the reference lenght Re =

2.3 Reboud correction

k2

with: f (ρ )= ρV +(1−α )

6 U inlet Lref = 2.7 10 : the Reynolds number

ν

9th International Conference on Hydrodynamics October 11-15, 2010 Shanghai, China

713

stepping method is 100. The time of simulation is about 4s. 3.4 Global analyses

Fig. 1 Schematic view of the Venturi geometry

With these parameters, a cavity length L ranging between 70 mm and 85 mm has been obtained. The experimental visualizations showed for this geometry a quite stable cavity behavior. It is characterized by a length almost constant in time, although the closure region always fluctuates with the existence of a re-enrant jet and little vapor cloud shedding (Fig. 2).

Fig. 2 Photograph of the experimental cavitation sheet [15]

3.2 Mesh The grid is a H-type topology. It contains 251 nodes in the flow direction and 62 in the orthogonal direction. A special contraction of the mesh is applied in the main flow direction just after the throat to better simulate the two-phase flow area (Fig. 3). The y+ values of the adjacent cells to walls vary between 12 and 27 for a non cavitating computation.

Fig. 3 View of the mesh near the throat

Different calculations were performed by considering different turbulent models, summarized in Table 1. The goal was to obtain a quasi stable cavitation sheet whose length varied between 70 – 85 mm with a re-entrant jet. Table 1. Unsteady computations Lcav/mm Model σinlet SA 0.59 75 KE 0.63 85 SA Reboud 0.57 79.5 KE Reboud 0.56 77 KE Durbin 0.55 82 KW SST 0.62 80

Comments steady results weakly unsteady results aperiodic quasi stable sheet aperiodic quasi stable sheet aperiodic quasi stable sheet periodic quasi stable sheet

Firstly, the standard turbulence models led to steady or weakly unsteady solutions with a very small recirculating area. On the contray, all corrected models were able to capture a quasi stable cavity sheet in good agreement with the experimental vizualisation. The tested limiters: Reboud, SST and realizability constraints make possible the simulation of unsteady behaviours. Secondly, the SST Menter model captured a periodic self-oscillating cavity whereas all other models simulated an aperiodic cavity. The frequency, evaluated with a direct Fourier transformation of the vapour volume signal (Fig. 4), is around 10Hz. Moreover, the σinlet value is around 0.62, a larger value in comparison with the experimental value equal to 0.55.

Fig. 4 Evolution of the vapour volume, KWSST model

3.3 Numerical parameters

3.5 Velocity and void ratio profiles

For the non cavitating regime, computations are started from an uniform flow-field using a local time step. For the unsteady cavitating regime, computations are performed with the dual time stepping method and are started from the non cavitating numerical solution. Δt U inle = 0.01 The dimensionless time step Δt* = Lref and the number of sub-iterations of the dual time

Local analyses concern void ratio and velocity profile comparisons inside the cavity. The experimental void ratio and velocity profiles are obtained for five stations by a double optical probe (Fig. 1). The velocity is evaluated as the most probable value and the void ratio is obtained with a post-processing algorithm from the signal of the double optical probe. The numerical values were obtained by a time-averaged treatment.

714

9th International Conference on Hydrodynamics October 11-15, 2010 Shanghai, China

Fig. 5 Velocity profiles, from station 1 (top)

to 5 (bottom)

Fig. 6 Void ratio profiles, from station 1 (top) to 5 (bottom)

9th International Conference on Hydrodynamics October 11-15, 2010 Shanghai, China Fig. 5 shows the evolution of the longitudinal velocity for the experiments and the computations. The overall agreement seems good between the experimental data and the simulations. For stations 1 and 2, no re-entrant jet phenomena occured in the experiment. Yet, the Menter and the KE Reboud models simulated a re-entrant jet. One possible explanation is that the μt limiter is too strong, which allows the development of the jet. Associated with the Durbin realizability constraints, the k-ε showed a better behaviour. Further downstream, for stations 3, 4 and 5, experimental observation indicates a recirculating behaviour with a re-entrant jet extending roughly half the sheet thickness. According to experiments, this flow configuration is smoothly time fluctuating (Fig. 2). This recirculating behaviour with a re-entrant jet is well simulated by all computations. At station 3, the thickness of the recirculating area is correctly estimated by all calculations, except by the KE Reboud model. At station 5, large discrepancies appear between models. Fig. 6 illustrates experimental and numerical results concerning the void ratio. For the first station, close to the throat, the vaporization phenomenon is clearly represented. The void ratio value is almost equal to 0.9 near the wall. For all computations, the cavity thickness is very well estimated. Downstream, at the second station, the void ratio is higher (around 96%). The distribution is similar to that obtained for station 1, with a correct estimation of the sheet thickness. The maximum value of the void ratio is under-estimated with the Menter model. The void ratio profile computed with the KE Reboud model present a non monotonic profile. From the third station, the re-entrant jet becomes noticeable, as observed before in the velocity field analyses. The void ratio value at the wall is largely over-estimated by all computations. As noted above, the cavity thickness is well predicted by all calculations. At station 4, large discrepancies are observed for the cavity thickness, the maximum and the wall void ratio values. Only the KE Reboud model gives a wall value close to the experimental one. At the last station, all computations over-predicted the maximum void ratio value.

715

constant values Pvap in the cavity. Dowstream, large discrepancies are notable between models. The re-compression is under-estimated with the k-ε models, especially the realizable version.

Fig. 7 Time-averaged wall pressure

The RMS wall pressure fluctuations are plotted in Figs. 8-9 versus the distance x-xinlet. The RMS fluctuation P'rms is divided by the time-averaged pressure Pav. For all computations, the statistical treatment was performed on a simulation time of 1 s. Experimental data indicate an augmentation of pressure fluctuation at the end of the cavity sheet, with a peak located at the fifth station. The peak position varies between models. With the Spalart-Allmaras and the Menter models, a very intense peak is present downstream the fifth station, whereas the peak obtained with the k-ε Reboud models is upstream. The realizable k-ε model provides fluctuations in better agreement with experimental data. Yet, for the last stations, fluctuations are over-predicted.

3.6 Wall pressure and RMS fluctuations The wall pressure distribution (P-Pvap/Pvap) is plotted in Fig. 7 versus the distance x-xinlet. The first five data are located inside the cavity (where the void ratio and velocity profiles are measured). For all computations, the pressure remains at an almost

Fig. 8 RMS wall pressure fluctuations, Spalart-Allmaras and Menter SST models

716

9th International Conference on Hydrodynamics October 11-15, 2010 Shanghai, China REFERENCES

Fig. 9 RMS wall pressure fluctuations, k-ε models

4 CONCLUSIONS A study of an aperiodic quasi stable cavitation sheet in a 2D Venturi configuration was performed by numerical RANS simulations. Compared to previous works done in the same geometry, more precise and reliable results concerning void ratio and velocity profiles were obtained. Results show that the performed numerical simulations were able to describe qualitatively the unsteady behaviour of the cavitation sheet. Global and local analyses of flows were proposed to compare different transport-equation turbulence models. Simulations demonstrated the determinant role of the use of a viscosity limiter. The Reboud correction associated with the Spalart-Allmaras and the Jones-Launder models provided a quasi stable sheet with a significant re-entrant jet in close agreement with the experimental data. Similar results were obtained with the realizable k-ε models. The Menter SST model captured a different solution: a low-frequency periodic cavitation sheet with a re-entrant jet. Additional works are in progress to investigate other turbulence models and to pursue comparative analyses between numerical and experimental studies.

[1] Reboud J L, Stutz B, Coutier O. Two-phase flow structure of cavitation: experiment and modelling of unsteady effects [C]. Proceedings of the 3rd Int. Symp. on Cavitation, Grenoble, France, 1998. [2] Coutier-Delgosha O, Fortes-Patella R, Reboud J-L. Simulation of unsteady cavitation with a two-equation turbulence model including compressibility effects [J]. Journal of Turbulence, 2002,3(58). [3] Chen Y, Lu C J, Wu L. Modelling and computation of unsteady turbulent cavitation flows [J]. .Journal of Hydrodynamics, Ser. B, 2006,18(5): 559-566. [4] Wilcox D. Turbulence Modeling for CFD [M]. La Canada, CA: DCW Industries Inc, 1998. [5] Vaidyanathan R, Senocak I., Wu J. et al. Sensitivity evaluation of a transport-based turbulent cavitation model [J]. Journal of Fluids Engineering, 2003,125(5): 447-458. [6] Wu J, Wang G, Shyy W. Time-dependent turbulent cavitating flow computations with interfacial transport and filter-based models [J]. Int. Journal for Numerical Methods in Fluids,2005,49(7):739-761. [7] Jameson A, Schmidt W, Turkel E. Numerical simulation of the Euler equations by finite volume method using Runge-Kutta time stepping schemes [C]. AIAA Paper 81-1259; 14th Fluid and Plasma Dynamics Conference, Palo Alto, California; June 1981. [8] Turkel E. Preconditioned methods for solving the incompressible and low speed compressible equations [J]. Journal of Comp. Physics, 1987,172(2):.277-298. [9] Goncalves E, Fortes Patella R. Numerical Simulation of Cavitating Flows with Homogeneous Models [J]. Computers Fluids, 2009,38(9): 1682-1696. [10] Delannoy Y, Kueny J L. Two phase flow approach in unsteady cavitation modelling [C]. Cavitation and Multiphase Flow Forum, ASME-FED,1990,98: 153-158. [11] Menter F R. Two-equation eddy-viscosity turbulence models for engineering applications [J]. AIAA Journal, 1994,32(8):1598-1605. [12] Jones W P, Launder B E. The Prediction of Laminarization with a Two-Equation Model of Turbulence [J]. Int. Journal of Heat and Mass Transfer ,1972,15: 301-314. [13] Spalart P R, Allmaras S. A one-equation turbulence model for aerodynamic flows [C]. AIAA 92-0439. 30th Aerospace Sciences Meeting, Reno, Nevada, 1992. [14] Durbin P A. On the k-e stagnation point anomaly [J]. Int. J. of Heat and Fluid Flow,1996,17: 89-90. [15] Barre S, Rolland J, Boitel G, et al. Experiments and modelling of cavitating flows in Venturi: attached sheet cavitation [J]. European Journal of Mechanics, 2009,28: 444-464.