Cavitation is a particular two-phase ﬂow with phase transition (vaporization/condensation) driven by pressure change without any heating. It can be interpreted as the rupture of the liquid continuum due to excessive stresses. Modeling of cavitating flows as a multi-fluid is a complex problem especially when the 3D consideration is adopted. However most recently research works are based on the mixture consideration of homogeneous fluid composed by two phases of liquid and vapor, which is also used by our model….
Minimizing the nuisance of cavitation is a great challenge in the design phase of a marine propeller. For efficiency reasons, the propeller usually needs to be operated in cavitating conditions but one still needs to avoid the effects of vibrations, noise and erosion. However, cavitation is a complex phenomenon not yet neither reliably assessable nor fully understood. Experimental observations can only give a part of the answer due to the obvious limitations in the measurement techniques; one example is measuring reentrant jets and internal flow, where flow features are hidden for optical measurement techniques by the cavity itself. Standard simulation tools used in design typically include potential flow solvers, lifting surface or boundary element approaches, with strict theoretical limits on cavitation modeling that only in the hands of an experienced designer may give satisfactory propeller designs. Adding to the challenge is a lack of theoretical knowledge of the physical mechanisms leading to harmful cavitation and thus how to modify a design if some form of nuisance is detected.
The transport equation models of cavitation suggested by Alajbegovic, Grogger et Philipp [ 1] and Yuan, et al., use the simplistic Rayleigh model. This model describes the limiting case of inertia-controlled growth of a spherical bubble in a liquid under a step variation in pressure of the surrounding liquid. However, this model cannot accurately describe bubble collapse and neglects a number of effects, which determine the behavior of cavitation bubbles.
Recently, progress has been made in the development of numerical models for calculation of cavitation flows. Though the models may differ in terms of realization (using the single-fluid or multi-fluid frame-work, the Eulerian- Eulerian or Eulerian- Lagrangian approaches), all of them are empirical to a certain level.
Modeling of cavitation flow as a multi-fluid is a complex problem which does not lead to satisfying results especially when the cavitation is modeled as 3D. Therefore in the engineering practice cavitation flow is often modeled as a single-fluid, where the cavitation area is handled as an area with the pressure lower then the vapour pressure. This approach always leads to the result, and the requirement of computer time is many times lower in comparison with multi-phase flow models. Moreover the steady solution of multiphase flow model may not be found at all due to the unsteady nature of cavitation flow.
Significant progress has been achieved recently in the development of homogeneous-mixture (single-fluid ) models for the simulation of three-dimensional transient cavitating flows (Chen and Heister  ; Kunz, Boger and Stinebring ; Ahuja, Hosangadi et Arunajatesan ; Yuan, et al. ; Singhal, et al.; Kubota, Kato et Yamaguchi ). These models allow single-fluid solvers to be applied to the conservation equations for the mixture, without increase in computational cost due to the increase in the number of conservation equations when applying the multi-fluid flow concept.
The present work is an investigation to develop a relevant physical model to simulate the cavitating flow. The main goal is the development of computational methodologies which can provide detailed description of the numerical set up for modeling and simulation the cavitation with the CFD code.
2. Mathematical formulation
To simulate cavitating flows, the two phases, liquid and vapour, need to be represented in the problem, as well as the phase transition mechanism between the two. Here, we consider a one fluid, single-fluid (mixture), introduced through the local vapour volume fraction and having the spatial and temporal variation of the vapour fraction described by a transport equation including source terms for the mass transfer rate between the phases. The numerical model solves the Reynolds averaged Navier-Stockes equations, coupled with a localized vapour transport model for predicting cavitation.
The fundamental equations governing the flow are taken with incompressible fluid case as given by the Navier-Stockes equations, which control the transport of momentum within the fluid (2) in addition to the mass conservation constraint (1):
The effective density and viscosity of the mixture are given (3) and (4) respectively :
Where α is the vapor fraction (α
Vapour fraction transport equation
The transport equation for the vapour scalar fraction α is given by:
This transport equation of volume fraction of vapour, with appropriate source terms to regulate the mass transfer between phases (liquid/vapor), is solved.
The proposed model formulation
This work deals with a numerical simulation of cavitation process around a hydrofoil. The numerical approach is based on predicting of the collapse process of pocket of vapor due to liquid compressibility. The model based is the Rayleigh Plesset (dynamics of spherical bubbles) integrated in the term source. Theinterface velocity investigated during the collapse process is given by:
finally, the proposed model can be express :
The value of is a calibrate parameter with experimental results; it is define the bubble density.
3. Numerical result
3.1. Numerical study
To validate the proposed model, a confrontation to experimental measures and to the numerical models (Yuan et al., Schnerr and Sauer and EOS). The application is a NACA0009 hydrofoil, truncated at 90% of the original chord length. It has the final dimensions of 100 mm of chord length.The 3D test section is modeled by a quasi 2D domain, with three rows of cells in spanwise direction for the numerical domain. The same mesh and numerical setup is used for the computations with the two models.
The domain is 9 blocks (See Fig. 1). The boundary conditions are set using a velocity inlet and a pressure at the outlet (the parameter which fixes the cavitation number).
The equations of the non-truncated hydrofoil are:
Numerical solutions of fluid flow and heat transfer problems are only approximate solutions. In addition to the errors that might be introduced in the course of the development of the solution algorithm, in programming or setting up the boundary conditions… The discretization approximations introduce errors which decrease as the grid is refined, and that the order of the approximation is a measure of accuracy. However, on a given grid, methods of the same order may produce solution errors which differ by as much as an order of magnitude. Theoretically, the errors in the solution related to the grid must disappear for an increasingly fine mesh [HYPERLINK \l "Fer96"9 ].
The pressure coefficient at σ=0.4 was taken as the parameter to evaluate six grids (See Fig. 2) and determine the influence of the mesh size on the solution. The selected convergence criteria were a maximum residual of 10−4. According to this figure, the grid with 53419 cells is considered to be sufficiently reliable to ensure mesh independence.
|Reference velocity@ Inlet||Cref||20 [m s-1]|
|Length of chord||0.11 [m]|
|Saturation pressure||3164 [Pa]|
|Liquiddensity||997 [kg m-3]|
|Vapourdensity||0.023 [kg m-3]|
|Reference time||0.0055 [s]|
|Pressure @ Outlet||Pout : calculate as a function of the value of|
Mass, momentum, turbulence and scalar transport equations were solved with Anys CFX which uses a coupled solver, which solves the hydrodynamic equations (for u, v, w, p) as a single system. This solution approach uses a fully implicit discretisation of the equations at any given time step. For steady state problems the time-step behaves like an ‘acceleration parameter’, to guide the approximate solutions in a physically based manner to a steady-state solution. This reduces the number of iterations required for convergence to a steady state, or to calculate the solution for each time step in a time dependent analysis. For the advection discretization, the High Resolution scheme is employed and the second order backward Euler scheme for the transient term, to assure accurate solution and to reduce the numerical diffusion for the solution .
In these calculations turbulence effects were considered using turbulence models, as the k-ε RNG models, with the modification of the turbulent viscosity. To model the flow close to the wall, standard wall-function approach was used. For this model, the used numerical scheme of the flow equations was the segregated implicit solver.
3.2. Validation steady flow
3.2.1. Turbulence model discussion
For steady state ﬂows, there is a priori no difference between the two equations formulations. The use of SST in the case of 2D steady hydrofoil computations instead ofk-ω or k-ε is justiﬁed in this way.
SST model works by solving a turbulence/frequency-based model (k–ω) at the wall and k-ε in the bulk flow. A blending function ensures a smooth transition between the two models. The SST model performance has been studied in a large number of cases. In a NASA Technical Memorandum, SST was rated the most accurate model for aerodynamic applications .
3.2.2. Pressure coefficient validation
A validation of the proposed model is described in this section. Firstly, we proceed with a distribution of the pressure coefficient around the surface of the hydrofoil and then we adjust with the profile velocity.
Fig. 3 presents the pressure coefficient distribution around the NACA 0009 for the proposed model and the Yuan et al. model with the experimental measurements. The proposed model shows a peak at the closing-pocket (condensation), this is due to the collapse velocity; it tends to compress the vapor pocket (collapse). The results from the proposed model can be satisfied compared to the experimental measurement.
We remark the good agreement of the numerical values compare to experimental data for R0=10-6 m and n0=1014 bubbles/m3. These constants are fixed for all this study.
Experimental data concerning the NACA0009 hydrofoil are reported for different cavitation numbers, and compared with the results of the computations in Fig 4.
These figures show satisfactory results of the proposed models in predicting the pressure distribution on the hydrofoil. As expected, the cavity becomes larger with decreasing cavitation number. However the models exhibit very different flow behavior at the cavitation detachment and closure regions. This is the result of dynamic bubble which is implementing in the proposed model. Clearly the model can perfectly predict the cavity of the vapor.
Figure 5 show satisfactory results of the proposed model in the pressure distribution on the hydrofoil. As expected, the cavity becomes larger with decreasing cavitation number. However the model exhibit very different flow behavior at the cavitation detachment and closure regions. This is the result of dynamic bubble which is taken account in the proposed model (interface velocity). Clearly the model can perfectly predict the cavity of the vapour.
3.2.3. Velocity distribution validation
Figure 6 show the predicted velocity profile of the main flow far from the wall is in good agreement with measurements. In the near-wall and in the wake region, the proposed model tends to reproduce the flow re-entrant jet.
3.3. Validation unsteady flow
3.3.1. 2D configuration(NACA0009)
The vapor cavity is characterized by a thick main cavity and an important shed cavity volume in a disorganized way, such as the cavity can be divided into many small cavities. For one period of cavity creation and collapse, one can remark two distinct life cycles highlighted by the lift and drag signals. Firstly, starting from the maximum cavity length, the closure region is in the small adverse pressure gradient and the reentrant jet is too thin, such as the cavity closure region is continuously broken into small vapor volumes. Secondly, as the reentrant jet reaches the region of the high pressure gradient, the reentrant jet is more important and the whole cavity is extracted and shed downstream. We consider the non-truncated variant of the 2D NACA0009 hydrofoil at high incidence angles. The flow parameters are i=5° and σ=1.2, and reflected by a high pressure gradient over the hydrofoil leading to an unsteady state flow behavior and shedding of large transient cavities. The domain is taken larger than the experimental test section to avoid numerical problem mainly due to reflections on the boundaries.
Cavitating flow are highly sensitive to turbulent fluctuations present in the flow  . The effects of modeling turbulence quantities have an enormous impact on the cavitation dynamics and the overall flow structure. The k-ε RNG turbulence model and LES model are adopted in this study.
3.3.2. Result and discussion
In this section we present the behavior of the proposed model with LES and the modified turbulence model. Most of characteristics of the unsteady flow is defined as a function of lift and drag coefficients, this coefficient can be expressed as:
For flows about 2D hydrofoils the dimensionless total vapour volume Vvap,2D is defined as :
where N is the total number of control volumes and α is the volume fraction of vapour in each control volume, with the volume Vi of the fluid in control volume.
The total vapour volume Vvap,2D, defined in equation (10), is a convenient parameter for understanding the transient evolution of the cavitating flow. The total vapour volume is calculated at each time step. After the start-up phase the growth and shedding of the vapour sheet and the collapse of the shed vapour cloud induce a self-oscillatory behavior, which is approximately periodic in time. The graphics of the analyzed variables (Fig.8) show a periodic signal, even if the identified periods are very different from each others. The use of Vvap,2D is the easiest way to identify the periodicity of the vapour formation and collapse during a life cycle. Furthermore, the time-history of the lift and drag coefficients are compared to those of the total vapour volume to correlate the occurring flow phenomena. The lift and drag coefficient has a cyclic time signal. Even if the periods are pretty clear, the fluctuations in a given period are very different from one to another. This phenomenon is mainly due to the dynamic of the shed cavities which are driven by the main flow-field downstream of the cavity closure. The growth and collapse mechanism of the cavity is driven by a cyclic phenomenon, whereas the dynamic of the cavity, when it is swept away can have very different non reproducible behavior. The vapour can be attached to the wall or far from it. This different ways of cavities shedding have an important influence on the pressure field at the hydrofoil wall, and there by on the lift and drag values.
In this section the cycles illustrated in figures 15 and 16 are considered. The solution for the volume fraction α above the hydrofoil is presented for a number of equidistant time-intervals during the cycle.
The vapour cavity is characterized by a thick main cavity and an important shed cavity volume in a disorganized way, such as the cavity can be divided into many small cavities. For one period of cavity creation and collapse, one can remark two distinct life cycles highlighted by the lift and drag signals.
Firstly, starting from the maximum cavity length, the closure region is in the small adverse pressure gradient and the reentrant jet is too thin, such as the cavity closure region is continuously broken into small vapour volumes. Secondly, as the reentrant jet reaches the region of the high pressure gradient, the reentrant jet is more important and the whole cavity is extracted and shed downstream.
3.3.3. 3D configuration(NACA66mod)
These excellent results obtained for 2-D cavitating flow conﬁrm the correct assumptions of the proposed numerical model in terms of vaporization and condensation processes, and to verify it’s performance with 3-D considerations, we present in the next some numerical results compared to experimental measurement  for cavitating flow around hydrofoil NACA66 (Fig 10).
In this section, a numerical simulation was performed to further study the performance of theproposed model and unsteady three-dimensional, the choice is justified because of the availability of experimental measurements. The profile used is NACA66 (mod) -312 a = 0.8 with a string of 0.15m, the width of the profile is 0.19m and is placed at 6 ° incidence relative to the upstream flow, not in an infinite medium -viscous. It has a thickness of approximately 12% and a camber on 2% to 45% and 50% of the leading edge along the rope. Theoretical points of this section have been interpolated by B-spline technique using the software mesh.
The study area is in 3D, it is made larger than the experimental test section to avoid numerical problems mainly due to reflections from the boundary conditions. The estate consists of 156 000 cells, the mesh topology is structured with a C-type conditions with wall on its border, the hydrofoil with the no-slip condition. (see Fig. 8). The steady state solution was used as an initial condition, the turbulence model used is k-ε RNG.
The proposed model was validated quantitatively and qualitatively compared in two-dimensional distribution of pressure coefficient and velocity profiles calculated with experimental measurements. We propose in the following one-dimensional confrontation, the application as we have described above is the NACA66 (mod), the advantage of selecting this application is the availability of measuring pressure along the profile for two cycles of detachment from the pocket of steam. We have reported three points of pressure measurement located respectively 50%, 70% and 90% chord.
The obtained numerical result presented by figure9 show that the model reproduce correctly the typical behaviour of partial cavity with development of re-entrant jet and the periodic shedding of cavitation clouds. Firstly, starting from the maximum cavity length, the closure region is in the small adverse pressure gradient and the re-entrant jet is too thin, such as the cavity closure region is continuously broken in to small vapour volumes. Secondly, as the re-entrant jet reaches the region of the high pressure gradient, the re-entrant jet is more important and the whole cavity is extracted and shed downstream.
This different ways of cavities shedding have an important influence on the pressure field at the hydrofoil wall, and there are characterized by the temporal evolution of the lift coefficient Cl compared to experimental result 
A comprehensive theoretical approach is done, and detailed formulations of the proposed model are presented. Influence of the numerical parameter has been widely studied. A comparative study is made for the steady flow. Good agreement with measurements was obtained for the proposed model. We have shown the importance of liquid compressibility on the pocket of vapour and how it is modeling in the proposed model. Computations with the proposed model is compared with experimental data and other numerical models, it shows the ability of the model to reproduce the steady-state developed cavitation flow fields. Finally, the unsteady behaviour of the cavitating flow depends strongly on the turbulence model, a k-ε RNG model is adopted in this steady.