Numerical analysis for supercavitating flows around axisymmetric cavitators
Numerical analysis for supercavitating flows around axisymmetric cavitators
International Journal of Naval Architecture and Ocean Engineering. 2013. Sep, 5(3): 325-332
Copyright ©2013, The Society of Naval Architects of Korea
  • Published : September 30, 2013
Export by style
Cited by
About the Authors
Young Kyun Kwack
Sung Ho Ko

Diffuse interface model for numerical analysis was used to compute supercavitating flows around various cavitators. The ambient pressures of 2 atm permitted cavitation studies in a range of cavitation number, σ=0.1 to 1.0 on selected conical and disk-headed cavitors. The computed results were compared with relation by Reichardt. Drag coefficient obtained from pressure forces acting on the cavitator also compared well with those obtained from analytical relations.
The movement of underwater vehicles with a suitable combination of body speed and water depth will produce local cavitation on the body. The extent and stability of the vapor cavities so formed may vary from small, discrete, transient bubbles to large and essentially stable bubbles enveloping the rear of the body. When a vapor- or gas-filled cavity grows until it is very long compared with the body dimensions, it is named a supercavity. Supercavities are formed by the growth of a fixed cavity or by the displacement of the liquid in a hydrodynamic wake, either because of a spread of vaporous cavitation, or because of the admission of gas into the low-pressure zones of the wake.
Compressible multiphase flow arises in many natural and industrial situations including bubble dynamics, shock wave interaction with material discontinuities, detonation of high energetic materials, hypervelocity impacts, cavitating flows, combustion systems to name only a few. The motivation of the present study is the accurate and computationally efficient resolution of interface problems in extreme flow conditions, as well as the computation of dynamic appearance of interfaces, that occur in cavitating flows. These interfaces are often separating pure media but also mixtures of materials in which wave dynamics is also important. Godunov type schemes and variants have reached a level of maturity to solve single phase flows in the presence of discontinuities. However, the presence of large discontinuities of thermodynamic variables and equations of state at material interfaces result in numerical instabilities, oscillations and computational failure (Abgrall, 1996; Karni, 1994).
To resolve these difficulties, two classes of methods have been developed. They are Sharp Interface Methods (SIM) and Diffuse Interface Methods (DIM).
In the Lagrangian class of SIM (Hirt et al., 1974; Farhat and Roux, 1991), the computational mesh moves and distorts with the material interface. However, when dealing with fluid flows, deformations are unbounded and resulting mesh distortions can make the Lagrangian approach unpractical (Scheffer and Zukas, 2000). Eulerian methods use a fixed mesh with an additional equation for tracking or reconstructing the material interface. In the volume of fluid (VOF) approach (Hirt and Nichols, 1981) each computational cell is assumed to possibly contain a mixture of both fluids and the volume occupied by each fluid is represented by the volume fraction, transported with the flow. This method is widely used for incompressible flows as there is no special thermodynamics to compute in mixture cells (Gueyffier et al., 1999).
The second type of methods (DIM) considers interfaces as numerically diffused zones, like contact discontinuities in gas dynamics. Diffuse interfaces correspond to artificial mixtures created by numerical diffusion. A pioneering work in this direction was performed by Abgrall (1996). Determination of thermodynamic flow variables in these zones is achieved on the basis of multiphase flow theory (Saurel and Abgrall, 1999; Abgrall and Saurel, 2003; Saurel et al., 2003; Murrone and Guillard, 2005; Abrall and Perrier, 2006; Saurel et al., 2007; Petitpas et al., 2007). The challenge is to derive physically, mathematically, and numerically consistent thermodynamic laws for the artificial mixture. The important issue is to fulfill interface conditions within this artificial mixture.
In this study, transient axisymmetric simulations of cavitating flows are performed using a diffuse interface model. The numerical method was used for cavitators with different geometries in a wide range of cavitation numbers and the results are compared with those of analytical relations.
- Numerical model
The diffused interface model, proposed in Saurel et al. (2009), is used in the computations of multiphase flow with an arbitrary number of fluids. The model is able to compute interfaces and fluid mixtures evolving under a unique pressure and unique velocity. It is thus particularly well adapted to cavitation studies with dynamic appearances of supercavity. The model is as follows:
Lager Image
where αk , ρk and
Lager Image
are the volume fraction, the density, the velocity vector of the phase k.
Lager Image
p and e are the mixture total energy, the mixture pressure and the mixture internal energy.
The mixture internal energy is defined as follows:
Lager Image
where Yk and ρ are the mass fraction and mixture density as follows:
Lager Image
Each fluid is governed by its own convex equation of state (EOS), e k = e k k , p), that allows the determination of the phases’ sound speed, c k = c k k , p). In the particular case of fluids governed by the stiffened gas EOS,
Lager Image
the resulting mixture EOS reads,
Lager Image
System (1) is hyperbolic with the three wave’s speeds u, u + c w and u−c w where c w is the multiphase extension of the non-monotonic Wood sound speed (Wood, 1930):
Lager Image
System (1) is also thermodynamically consistent as it agrees with the second law of thermodynamics. But system (1) presents inherent difficulties for its numerical solution. The difficulties are as follows:
  • (1) Volume fraction positivity: How to treat the non-conservative term in the volume fraction equation when shocks or strong rarefaction waves are present?
  • (2) The volume fraction varies across acoustic waves: Riemann solver is difficult to construct.
  • (3) The mixture sound speed has a non-monotonic behavior.
A pressure non-equilibrium model is built to circumvent these difficulties (Saurel et al., 2009). It is solved by a 4-step method:
  • (1) At each cell boundary solve the Riemann problem with the HLLC solver (Toro et al., 1994; Toro, 1997).
  • (2) The N-phase pressure non-equilibrium flow model is solved without relaxation effects with a Godunov-type scheme.
  • (3) A pressure relaxation solver is used to reach pressure equilibrium.
  • (4) Internal energies rest is achieved with the help of the mixture total energy to guarantee shock waves transmission through interfaces and improve robustness.
- Computational domain and flow conditions
In this study, transient axisymmetric simulations for 2D supercavitating flows shown in Fig. 1 are performed using a diffuse interface model. Calculations for supercavitating flows around various cavitators (30°, 45° cones and disks) were performed. The left boundary of the computational domain is velocity inlet boundary condition with given velocity (20, 25, 30, 40, 50 and 60 m/s ). The right boundary is an outflow condition with imposed pressure. The ambient pressure is taken equal to 2 atm . The top boundary is wall condition and the bottom boundary is axisymmetric condition.
Grids generated by the cartesian grid were concentrated around a head and end of the body, as shown in Fig. 2 . The grid consisted of approximately 30,000 cells. Courant-Friedrichs-Lewy condition (CFL) number was set to 0.4. Using Xeon (2.4 GHz ) for parallel computing, the computation time of 15 hours per case was spent.
Lager Image
Schematic for model and boundary conditions.
Lager Image
Grid generation.
When the pressure in a liquid flow falls below the corresponding vapor pressure, the liquid evaporates. This phenomenon, named cavitation, has a crucial effect on the performance of hydraulic systems. Pumps, nozzles, turbine blades, and hydrofoils are just a few examples where the occurrence of cavitation is an important design issue. Cavitation is categorised by a non-dimension number called cavitation number as follows:
Lager Image
In this, pv is the vapor pressure, ρ is the liquid density, and P , V are the main flow pressure and velocity. When a liquid flows over a solid object, as the fluid velocity increases five different cavitation regimes are observed in the flow over the body: incipient-, shear-, cloud-, partial- and supercavitation. Partial and cloud cavitation refer to the situation where vapor region extends over a subsection of the body. On the other hand, supercavitation refers to a long steady cavity that extends more than the cavitator length and closes in the liquid. In all cavitation regimes, there is a constant movement of reentrant liquid jet in the cavity closure region. In the cloud cavitation regime, this inward movement of liquid results in a detachment of large vapor structures from the main body. In the partial cavitation, however, there is no vapor detachment and the cavity length remains constant (Franc and Michel, 2004). Supercavitation is the main focus of this paper.
Shape of supercavity and drag coefficients are important parameters on the study of supercavitating flows. The relations which are indicated as theory are based on Reichardt’s analysis (Reichardt, 1946; Self and Ripken, 1955). The relations are as follows:
Lager Image
Lager Image
Lager Image
Lager Image
In these, dmax is the maximum cavity diameter, lc is the cavity length, D is the test cavitator diameter, CD is the drag coefficient at cavitation number σ and CD0 is the drag coefficient at σ = 0. Values of CD0 were extrapolated from experimental data; thus the equations are actually semi-empirical. Reichardt found CD0 for disks to be 0.79. The CD0 values for the disk and the 45- degree cone were given to 0.80 and 0.26, respectively (Self and Ripken, 1955). Since Reichardt’s studies were limited to cavitation indexes below 0.12, the preceding equations may not be applicable to larger values. And the CD0 values for the 30- degree cone were given 0.205 by Plesset and Shaffer (1948).
Disk drag coefficient by the approximate theory of Pleseet and Shaffer can be represented for σ up to 1.5 by an interpolation formula as follows:
Lager Image
This result was obtained by rotating the two-dimensional flat-plate pressure distribution about the axis of symmetry.
The cavity dimensions obtained from numerical results for 30°, 45° cones and disks are shown in Fig. 3 . It is notable that for cavitation numbers above approximately 0.44, the cavitators with larger angle β tended to produce larger cavities. Decreasing σ results in an increase in the length of the supercavity, with the dramatic increase in low cavitation number.
In Figs. 4 to 6 , the computed results were compared with analytical relations. It can be seen that the numerical results compare well with analytical relations. In low cavitation number, the discrepancy between the relations and numerical results may be attributed due to the effects of the wall. According to Franc and Michel (2004), in a flow confined by solid walls, the cavity becomes infinite for a value of the cavity pressure Pc smaller than the pressure at infinity P for a non-zero value of the cavitation number. This phenomenon corresponds to blockage of the flow.
Lager Image
Cavity length vs. cavitation number for various cavitators.
Lager Image
Cavity length vs. cavitation number for a 30° conical cavitator.
Lager Image
Cavity length vs. cavitation number for a 45° conical cavitator.
Lager Image
Cavity length vs. cavitation number for a disk-headed cavitator.
Lager Image
Comparison of the cavity length non-dimensionalized by chord length for various cones.
Lager Image
Resulting pressure force acting on the cavitator versus time (σ = 0.99, v = 20 m/s, d = 15 mm, angle = 30°).
Also, predicted cavity length was compared with an analytic solution derived by Newman (1977) in Fig. 6 . The relation shows the relationship between the cavitation number (σ) and the cavity length ( lc ) for a two dimensional symmetric body as follows:
Lager Image
where y 0 ( t ) is the half diameter of the body. Present results are in good agreement with analytic solution as the cavitation numer decreases.
The evolution of the resulting pressure force on the cavitator is shown in the Fig. 8 . During the first times, the force is very big as results of initial conditions stiffness that produce a shock wave. When the flow reaches quasi steady state, it is convergence.
Lager Image
Drag coefficient vs. σ for various cavitators.
Drag coefficient for each cavitor was obtained when considering the definition as shown below:
Lager Image
where A = (1/4)πD 2 .
The drag coefficient against cavitation number is plotted in Fig. 9 . As shown in the figure, a good agreement is observed between the results of the model with analytical relations. Important point is that as cavitation number decreases, in other words velocity increases, the drag coefficient decreases.
Diffuse interface model for numerical analysis was used to predict the cavitation that occurs behind axisymmetric cavitator in a liquid flow. The computed results were compared with relation by Reichardt. Drag coefficient obtained from pressure forces acting on the cavitator also compared well with those obtained from analytical relations. The size of cavity is associated with a certain head form. Decreasing cavitation number yields an increase in the length of the supercavity, with the dramatic increase in low cavitation number.
This work was carried out with the support of Agency for Defense Development (ADD) of Korea under Grand No. 09-01-05-23.
Abgrall R. 1996 How to prevent pressure oscillations in multicomponent flow calculations: a quasi conservative approach. Journal of Computational Physics 125 (1) 150 - 160    DOI : 10.1006/jcph.1996.0085
Abgrall R , Saurel R. 2003 Discrete equations for physical and numerical compressible multiphase mixtures. Journal of Computational Physics 186 (2) 361 - 396    DOI : 10.1016/S0021-9991(03)00011-1
Abgrall R. , Perrier V. 2006 Asymptotic expansion of a multiscale numerical scheme for compressible multiphaseflow. SIAM Journal of Multiscale Modeling and Simulation 5 (1) 84 - 115    DOI : 10.1137/050623851
Farhat C. , Roux F.X. 1991 A method for finite element tearing and interconnecting and its parallel solution algorithm. International Journal for Numerical Methods in Engineering 32 (6) 1205 - 1227    DOI : 10.1002/nme.1620320604
Franc J.P. , Michel J.M. 2004 Fundamentals of cavitation. Springer
Gueyffier D. , Li J. , Nadim A. , Scardovelli R. , Zaleski S. 1999 Volume-of-fluid interface tracking with smoothedsurface stress methods for three-dimensional flows. Journal of Computational Physics 152 (2) 423 - 456    DOI : 10.1006/jcph.1998.6168
Hirt C.W. , Nichols B.D. 1981 Volume of fluid (VOF) method for the dynamics of free boundaries. Journal of Computational Physics 39 201 - 225    DOI : 10.1016/0021-9991(81)90145-5
Hirt C.W. , Amsden A.A. , Cook J.L. 1974 An arbitrary Lagrangian Eulerian computing method for all flow speeds. Journal of Computational Physics 14 (3) 227 - 253    DOI : 10.1016/0021-9991(74)90051-5
Karni S. 1994 Multicomponent flow calculations by a consistent primitive algorithm. Journal of Computational Physics 112 (1) 31 - 43    DOI : 10.1006/jcph.1994.1080
Murrone A. , Guillard H. 2005 A five equation reduced model for compressible two phase flow problems. Journal of Computational Physics 202 (2) 664 - 698    DOI : 10.1016/
Newman J.N. 1977 Marine hydrodynamics. The MIT Press Cambridge
Petitpas F. , Franquet E. , Saurel R. , Le Metayer O. 2007 A relaxation-projection method for compressible flows. PartII. The artificial heat exchange for multiphase shocks. Journal of Computational Physics 225 (2) 2214 - 2248    DOI : 10.1016/
Plesset M.S. , Shaffer P.A. 1948 Cavity drag in two and three dimensions. Journal of Applied Physics 19 (10) 934 - 939    DOI : 10.1063/1.1697899
Reichardt H. 1946 The laws of cavitation bubbles at axially symmetric bodies in a flow. Ministry of Aircraft Production (Britain) Report and Translations No. 766
Saurel R. , Abgrall R. 1999 A multiphase Godunov method for compressible multifluid and multiphase flows. Journal of Computational Physics 150 425 - 467    DOI : 10.1006/jcph.1999.6187
Saurel R. , Gavrilyuk S , Renaud F 2003 A multiphase model with internal degrees of freedom: application to shockbubbleinteraction. Journal of Fluid Mechanics 495 283 - 321    DOI : 10.1017/S002211200300630X
Saurel R. , Franquet E. , Daniel E. , Le Metayer O. 2007 A relaxation-projection method for compressible flows. PartI. The numerical equation of state for the Euler equations. Journal of Computational Physics 223 (2) 822 - 845
Saurel R. , Petitpas F. , Berry R. 2009 Simple and efficient relaxation methods for interfaces separating compressiblefluids, cavitating flows and shocks in multiphase mixtures. Journal of Computational Physics 228 (5) 1678 - 1712    DOI : 10.1016/
Scheffer D.R. , Zukas J.A. 2000 Practical aspects of numerical simulation of dynamic events: material interfaces. International Journal of Impact Engineering 24 (8) 821 - 842    DOI : 10.1016/S0734-743X(00)00003-8
Self M.W. , Ripken J.F. 1955 Steady-state cavity studies in a free-jet water tunnel. St. Anthony Falls Hydraulic Laboratory Report No. 47
Toro E.F. , Spruce M. , Speares W. 1994 Restoration of the contact surface in the HLL Riemann solver. Shock Waves 4 25 - 34    DOI : 10.1007/BF01414629
Toro E.F. 1997 Riemann solvers and numerical methods for fluids dynamics : a practical introduction. Springer New York
Wood A.B. 1930 A textbook of sound. Bell and Sons Ltd London