Numerical study of the pitching motions of supercavitating vehicles

951

2012,24(6):951-958 DOI: 10.1016/S1001-6058(11)60323-5

NUMERICAL STUDY OF THE PITCHING MOTIONS OF SUPERCAVITATING VEHICLES* YU Kai-ping, ZHANG Guang , ZHOU Jing-jun, ZOU Wang, LI Zhen-wang School of Aeronautics, Harbin Institute of Technology, Harbin 150001, China, E-mail: [email protected]

(Received March 14, 2012, Revised September 7, 2012) Abstract: The pitching motions of supercavitating vehicles could not be avoided due to the lost water buoyancy. In order to have some insight for the design of the supercavitating vehicles, the fixed frequency and free pitching motions are investigated. A numerical predicting method based on the relative motion principle and the non-inertia coordinate system is proposed to simulate the free pitching motions of supercavitating vehicles in the longitudinal plane. Homogeneous and two fluid multiphase models are used to predict the natural and the ventilated supercavitating flows. In the fixed frequency pitching motions, a variety of working conditions are considered, including the pitching angular velocities and the supercavity scales and the results are found to be consistent with the available experimental results in literature. The mesh deformation technology controlled by the moment of momentum equation is adopted to study the free pitching motions and finally to obtain the planing states proposed by Savchenko. The numerical method is validated for predicting the pitching motions of supercavitating vehicles and is found to enjoy better calculation efficiency as comparing with the mesh regeneration technology. Key words: supercavitating flow, two fluid multiphase flow model, homogeneous model, pitching motion, mesh deformation

Introduction  The natural supercavity is achieved when an underwater vehicle travels at a sufficient high speed, so that the fluid pressure drops to the saturated vapor pressure. When the vehicle is enveloped in the cavity, the viscous drag drops drastically thus the vehicle might move as in the air. Even for vehicles designed to travel at the natural supercavitating velocity, the drag must be firstly reduced by the ventilated supercavitating to enable the vehicles to accelerate to the conditions at which the natural supercavity can be sustained. It is well proven that the ventilated cavitation can achieve a significant drag reduction and, therefore it has attracted a growing research attentions, especially, by means of the CFD. While the ventilated supercavity can dramatically reduce the drag of the vehicle, the lost water buoyancy due to the cavity enclosure makes the hydrodynamics quite different from the whole wet state. Besides the

* Project support by the Major National Natural Science Founation of China (Grant No. 10832007). Biography: YU Kai-ping (1968-), Male, Ph. D., Professor

lift on the cavitator and the fins, another lift from the interaction between the vehicle and the cavity boundary can not generally be avoided[1-3] due to the moment balance, and the wet area is the main source of the stable restoring force[4,5]. Thus, a study of the restoring force of the aft body and the further pitching motion of the vehicle is of great significance to the design of supercavitating vehicles. In the published literature, only a few papers concern with the force due to the interaction between the aft section of the vehicle and the cavity boundary and most papers deal with the trajectory control based on the Hassan’s theory that describes the forces and the moments experienced by a cylindrical body steadily planing on flat and cylindrical free surfaces[6]. No papers are found to study the pitching motions with CFD predictions. In papers[7,8], the interaction of the aft body of the vehicle with the cavity wall pitches in the supercavity was discussed based on experiments in water tunnel. In papers[9,10], a series of experiments and numerical simulations were conducted to investigate the hydrodynamic drag, lift and moment characteristics at the aft section of underwater supercavitating bodies for different supercavity shapes. In papers[11,12] Kulkarn and Kirschner et al. studied the

952

impact dynamics of a supercavitating projectile with its nose being prescribed to move in a straight line. Their models predicted that the motion of the projectile consists of a series of impacts with the frequency increasing as the projectile slows down. In this paper, a numerical method based on the relative motion principle and the non-inertia coordinate system is proposed to study the fixed frequency and free pitching motions of supercavitating vehicles. The method is validated by comparing the calculated results with those in papers[7,8]. At the same time, the method is used to successfully predict the planing state proposed by Savchenko. Therefore, the method can be used in the design of supercavitating vehicles.

1. Numerical methods 1.1 Basic governing equations of two-fluid model In the multiphase system considered in this paper, the fluid is assumed to be isothermal and the densities of the fluid are functions of the pressure but not the temperature. Under this assumption, only the continuity and the momentum equations are used, but not the energy equations. To simulate the ventilated cavitation flows is to solve the standard 3-D Navier Stokes equations, the turbulence equations, the gas state equation and the volume fraction equation. The continuity equation for the single phase is

w (J D UD ) + ’ (J D UDVD ) = 0 , D = 1, 2 wt

(1)

The momentum equation for the single phase is w (J D UD VD ) + ’ [J D ( UD VD … VD )] = wt

^

`

’ J D PD ª¬’VD + (’VD )T º¼  J D ’pD +

M D + J D UD g , D = 1, 2

(2)

For a ventilation cavitating flow, D represents the gas phase or the liquid phase, for a natural cavitating flow D represents the vapor phase or the liquid phase. The volume fraction equation is N

¦J D = 1

(3)

D =1

where U is the density of the air, P is the local pressure, T is the temperature which keeps constant in the simulations. Formulas (1)-(4) are the whole governing equations. 1.2 Cavitation model In this paper, the Zwart-Gerber-Belamri model is used to predict the natural supercavitation. 2 pv  p , 3 Ul

m  = Fc

3D v U v RB

m + = Fe

3rnuc U v (1  D v ) 2 pv  p RB 3 Ul

(5)

where subscripts v , l respectively represent the vapor and liquid phases, rnuc is the volume fraction at the nucleation sites, RB represents the radius of the nucleation sites, while F is an empirical factor which may differ for condensation and vaporization, p is the pressure in the liquid surrounding the bubble, and pv is the saturated vapor pressure with the value of 3 170 Pa in this paper. 1.3 Turbulence model In this paper, with due consideration of the computational expense, a steady computation with the Shear Stress Transport (SST) turbulence model is carried out to give the initial data for the unsteady calculations in which the Detached Eddy Simulation (DES) turbulence model is utilized. A simple introduction of the two turbulence models is given below. 1.3.1 SST turbulence model The predominant disadvantage of the standard k  H turbulence model is that the model is designed for high Reynolds number flows thus it performs poorly for low Reynolds number cases and in the near wall regions. The Wilcox k  Z turbulence model was firstly introduced by Kolmogorov, some improvements were made with the most popular suggestions made by Wilcox, who formulated a low Reynolds number alternative to the standard k  H turbulence model. The deficiency of the k  Z turbulence model stems from the sensitivity of the model to Z . One possible solution is to use a combination of the k  Z model equations implemented near wall regions and the k  H model to be employed in the bulk flow region, which leads Menter to formulate the SST turbulence model. The details can be found in the paper[13].

The gas state equation is

U

U ( P, T )

(4)

1.3.2 DES turbulence model In order to improve the predictive capabilities of

953

turbulence models in highly separated regions, Spalar proposed a hybrid approach, which combines features of classical RANS formulations with elements of Large Eddy Simulations (LES) methods. The concept is termed the DES and it is based on the idea of covering the boundary layer by a RANS model and switching the model to a LES mode in detached regions. Ideally, the DES can predict the separation line from the underlying RANS model, but which means to capture the unsteady dynamics of the separated shear layer with the resolution of the developing turbulent structures. Compared to the classical LES methods, the DES saves the computing resources by orders of magnitude for high Reynolds number flows. Though this is due to the moderate costs of the RANS model in the boundary layer region, the DES still offers some of the advantages of an LES method in separated regions. 1.4 Numerical resolutions The simulations are based on 3-D calculations and a finite volume discretization of these equations is used. A solver of the coupled conservation equations of mass and momentum is adopted, with an implicit time scheme and the multigrid technology. The transient term is treated with a second-order implicit scheme.

boundary and extrapolated to the inlet boundaries. At the external boundary, the pressure and the volume fractions are extrapolated and the opening boundary condition is specified. The details are shown in Fig.1.

Fig.1 The boundary conditions

Fig.2 Two different models in simulations

Fig.3 The computational grids

w ( U) ) 1.5U n +1) n +1  2 U n) n + 0.5U n 1) n 1 = wt 't

(6)

where ) stands for u , v , or w , the three velocity components. The diffusive terms are calculated in a central difference manner. The advection schemes can be cast in the form

Models with and without fins are applied in the simulations as shown in Fig.2. The computational grids used in this paper are shown in Fig.3. Before carrying out the formal computations, different meshes (shown in Table 1) were used to predict the pressure coefficient on the cone and the cylinder sections of the body to test the influence of the grids on the computational results.

)ip = )up + E ’) 'r

Table 1 Different grid parameters

(7)

where )up is the value at the upwind node, r is the vector from the upwind node to ip and a high resolution scheme is adopted which uses the a special nonlinear recipe E at each node, being computed to be as close to 1 as possible without introducing new extrema. The recipe E is based on the boundedness principles.

2. Computational models and grids The velocity components, the volume fractions, the turbulence intensity and the length scale are specified at the velocity inlet boundary and extrapolated to the pressure outlet or the opening boundaries. The mass inlet boundary is defined at the blowhole. The pressure distribution is specified at the pressure outlet

Grid number

Grid size of the first layer on the model (m)

Axialuradialu circumferential

Grid-1

0.00025

210u48u68

Grid-2

0.00025

225u58u78

Grid-3

0.00025

245u70u88

Figure 4 shows the computational results of the pressure coefficient on the cone and the cylinder sections for model 1 with three kinds of grids. It can be seen that the differences between grid-2 and grid-3 are small, and with due consideration of the calculation accuracy and expense, the grid-2 is finally used for all subsequent calculations. 3. The computational results The following calculations mainly include two

954

parts: the fixed frequency pitching motions and the free pitching motions that are controlled by the moment of momentum equation in the longitudinal plane. The fixed frequency pitching motions are calculated preliminarily to determine the influencing factors on the cavity scale and the lift characteristics of the body and the computational results can be used to compare with the experimental results in paper[8] to validate the numerical method.

Figure 6 shows the pressure evolution at the measuring point as shown in Fig.5 during the fixed-frequency pitching motion. The simulation results show that the pressure on the point varies significantly, which is consistent with the experimental results in paper[8] and shows that the numerical method is valid. In Fig.6 t = t / T is a dimensionless variable. Now we consider the influencing factors in the pitching motions. The influences of the angular velocity on the cavity shape and the hydrodynamics are as follows. The pitching angle law: T = [1/ 30sin (Z t )] rad , Z = 4S, 6S, 8S , the free stream velocity is increased to 25 m/s, and the ambient pressure turns to 2.5 atm.

Fig.4 Comparison of predicted surface pressure coefficient distributions on the model 1 with different grid

Fig.5 The sketch map of measurement point

Fig.7 The change of the supercavity shape in the process of vibration

Fig.6 The change law of pressure at measuring point

Figure 7 shows the change of the cavity shape when Z = 4S , and the mass flow rate is 0.013456 kg/s. It can be seen that the supercavity is stable in the front part of the supercavitating vehicle, and large amounts of gas escape from its tail because of the periodical tail-clap of its vehicle.

3.1 The results of fixed frequency pitching motions The computational conditions in paper[8] are as follows: The upstream velocity Vf is 10.25 m/s, the upstream pressure Pf is 33 800 Pa, the total ventilation rate CQ = Q / Vf D 2 is 0.144, and the ventilated cavitation number is 0.0864. The pitching angle follows the law: T = T max sin (Z t ) , where T max = (S / 180) rad is the maximal pitching angle, Z = 10S The time step in the unsteady computation is 0.005 s. The initial state of all simulations in this paper is a steady calculated result with the pitching angle T = 0 rad and the vehicle is fully enveloped in the supercavity.

Fig.8 The final cavity shape with different pitching angular velocity Z = 4 S, 6S, 8S

Figure 8 shows the final cavity shape of model 1 for different angular velocities when the pitching

955

angular T = 0 rad after the cavity reaches a relatively stable state, and the three pitching motions have the same initial state. It can be seen that the pitching angular velocity has a great influence on the gas leakage. The pitching motions with higher frequencies make the gas leak more quickly.

Fig.9 The lift coefficient with different angular pitching angular velocity

cient are defined as Cl =

Fl Fd , Cd = 2 0.5U v S 0.5U v 2 S

(8)

3.2 The results of free pitching motions In paper[14] suggested four possible modes of vehicle motions according to the vehicle velocity and the supercavity development, among which two modes are demonstrated in what follows. The velocity of the first motion mode is about 70 m/s and that of the other is in the range of 50 m/s-200 m/s, as shown in Fig.11. In this paper, the velocity is different, although the motion modes are the same, as obtained by adjusting the ambient pressure. Two free velocities adopted in this paper are 100 m/s and 110 m/s when the vehicle is enveloped by a complete natural supercavity. The other two velocities are 25 m/s and 38 m/s.

Fig.11 Two possible modes of vehicle motion

Fig.10 The drag coefficient with different angular pitching angular velocity

Figure 9 and Fig.10 demonstrate the lift and drag variations of model 1 with different pitching angular velocities. It can be seen that the larger the pitching angular velocity, the larger the lift and the drag become, which shows that the impact velocity between the aft body and the cavity boundary would affect the lift and the drag. When the angular velocity Z = 4S , the lift and the drag vary regularly from beginning to end, because the aft body hardly affects the leakage of gas and on the other hand, when Z = 6S and Z = 8S , the lift and the drag change significantly in the process of the pitching motion from the initial state, however, finally, the maximal lift and drag come to a close value, which demonstrates that the cavity shape or the closure position of the cavity has a greater influence on the hydrodynamics of the body as comparing with that of the pitching angular velocity. On the other hand, the maximal and minimum lift coefficients show an asymmetric feature, that is why the asymmetric cavity shape affects the closure angle of the cavity on the body. The lift coefficient and the drag coeffi-

The body reference frame is built on the center of the cavitator, the dynamics equations in the longitudinal plane in the body coordinate are m(u  rv  r 2 x g ) = Fx

(9)

 g ) = Fy m( v + ru + rx

(10)

I z r + mxg v + mxg ru = M z

(11)

where u and v are the velocities along the directions of x and y in the body reference frame, x g is the coordinate of the center of gravity, r is the pitching angular velocity of the body. In this paper, the motion of the cavitator center is limited to the direction of X e in the inertial coordinate just like the axial restrained model experiments in paper[15] so that Eqs.(9) and (10) can not be solved due to the unknown forces Fx and Fy , only the moment of momentum Eq.(11) is solved and if Eqs.(9) and (10) are also solved then we will have the uncontrolled trajectory. The kinematics equations are

956

X e = u cos T  v sin T

(12)

Ye = u sin T + v cos T

(13)

lly keeps a planing state which can be seen in Fig.15 as it is shown in Fig.11(b) and the final pitching angle is smaller than 0.004 rad.

Where T is the pitching angle. It is assumed that X e = const and Ye = 0 , so we have u = X e cos T and v =  rX e cos T . Substitute u and v into Eq.(11), and the final equation is I z r  M z = 0

(14)

where I z = 1.34 kgm 2 , m = 8.8 kg , xg = 0.3 m . M z is calculated by numerical simulations. The initial state of the vehicle is enveloped by the supercavity with zero pitching angle and the steady numerical predicting results are used. The time step of the unsteady computation is 0.001 s.

Fig.14 The initial and final state of the pitching motion when

X e = 110 m / s

Fig.15 The change process of pitching angle when X e = Fig.12 The initial and final state of the pitching motion when

X e = 100 m / s

110 m / s

The above analysis shows that the closure position of the cavity would seriously affect the lift characteristics of the vehicle thus further affect the gesture of the vehicle and the similar conclusions are obtained in the fixed frequency pitching motions. In order to understand the free pitching motions with a low free stream velocity, we let the free velocity X e = 25 m / s and X e = 38 m / s , the mass flow rate be 0.0134 kg/s, and consider the free pitching motion with ventilated supercavity.

Fig.13 The change process of pitching angle when X e =

100 m / s

Figure 12 shows the initial state of the free pitching motion and the subsequent relatively stable state, with the vehicle being enveloped by the natural supercavity. It can be seen that the vehicle moves with a two-cavity state, and the initial state is obtained by adjusting the ambient pressure. Figure 13 shows the variations of the pitching angle with time, which shows that after about 0.6 s, the variations of the pitching angle become stable, in the range between 0.004 rad and 0.004 rad. However, with the initial state as shown in Fig.14 with the free velocity X e = 110 m / s , the vehicle fina-

Fig.16 The change process of pitching angle when X e = 25 m/s, 38 m/s

Figures 16 and 17 show that the vehicle finally enters the planing state, when the free velocity X e =

957

38 m / s , the vehicle reaches the stable state sooner than when X = 25 m / s , and the final stable pitching e

angle is smaller, which shows that the free stream velocity indeed affects the lift of the aft body.

Fig.17 The final state of the pitching motion when mass flow rate is 0.0134 kg/s

Fig.18 The motion state when X e = 25 m / s , mass flow rate is 0.0045 kg/s

Fig.19 The evolution of pitching angle when X e = 25 m / s , mass flow rate is 0.0045 kg/s

4. Conclusions In this paper, a two fluid multiphase model and the homogeneous model are used to simulate the ventilated supercavity and the natural supercavity, respectively. The SST turbulence model is used to carry out the steady calculation as the initial data of the subsequent unsteady simulations in which the DES turbulence is assumed. On this basis, the fixed frequency and free pitching motions under different working conditions are investigated. Firstly, the pressure measuring points on the vehicle are selected to obtain the variations of the pressure, which are found to be consistent with the conclusions in paper[8], so that the numerical method is validated. Secondly, the pitching angular velocity of the vehicle is found to have a great influence on the gas leakage. The pitching motions with higher frequencies make the gas leak more quickly and the final cavity scale becomes smaller. Finally, the free pitching motions are investigated by restricting the motion using the equations of the moment of momentum and the mesh deformation technology. The planing states are successfully predicted, including the states enclosed by the natural supercavity and the ventilated supercavity. The cavity scale and the closure position of the cavity are found to seriously affect the lift on the aft body and maybe lead to the trajectory instability. The numerical method in this paper can predict all kinds of gestures of the vehicle under different working conditions in the longitudinal plane, and can even predict the trajectory. The method has an obvious calculation advantage comparing with the mesh regeneration technology. The whole trajectory prediction in the longitudinal plane will be carried out in future. References

In order to know whether the vehicle might also reach a stable state when the cavity has a smaller scale and as a different closure way to verify the numerical method, we let the free velocity be 25 m/s, the mass flow rate be 0.0045kg/s, the initial and calculated termination states of the body are shown in Fig.18. The vehicle finally does not reach a steady state in the process of the vehicle motion. From Fig.16 and Fig.19, it can be seen that although the free velocity is the same, the vehicle has a different gesture due to the different cavity shape and closure way of the cavity on the body thus it has a different lift. Therefore, it can be seen that when the vehicle enters the planing state, the lift on the aft body is stable, so that the lift characteristics can be studied in water tunnel experiments.

[1]

[2]

[3]

[4]

[5]

KIRSCHNER I. N., KRING D. C. and STOKES A. W. et al. Control strategies for supercavitating vehicles[J]. Journal of Vibration and Control, 2002, 8(2): 219242. DENG Fei, ZHANG Yu-wen. On hydrodynamic design of underwater supercavitating high speed vehicle[J]. Journal of Northwestern Polytechnical University, 2004, 22(6): 806-810(in Chinese). PUTILIN S. I. Some features of a supercavitating model dynamics[J]. International Journal of Fluid Mechanics Research, 2001, 28(5): 631-643. ZHOU Jing-jun, YU Kai-ping and ZHANG Guang. Numerical simulation on the process of supercavity development and the planing state of supercavitating vehicle[J]. Journal of Ship Mechanics, 2011, 15(3): 200206(in Chinese). SAVCHENKO Y. N. Experimental investigation on supercavitating motion of bodies[C]. RTO AUT

958

[6]

[7]

[8]

[9]

[10]

Tecture Series on Supercavitating Flows. Brussels, Belgium: Von kámán Institute, 2001. HASSAN S. E. Analysis of hydrodynamic planing forces associated with cavity riding vehicles[C]. NUWCNOT Technical Memorandum 90085. 1999. LI Ji-tao, HU Tian-qun and HE You-sheng. Experimental investigation of supercavitating vehicle motion with oscillatory mode[C]. Proceedings of the Seventh China National Conference on Experimental Fluid Mechanics. Beidaihe, China, 2007(in Chinese). LEE Qi-tao, XUE Lei-ping and HE You-sheng. Experimental study of ventilated supercavities with a dynamic pitching model[J]. Journal of Hydrodynamics, 2008, 20(4): 456-460. JIANG Zeng-hui, YU Kai-ping and WANG Cong et al. Experimental research on hydrodynamic character of aft section of underwater supercavitating bodies[J]. Engineering Mechanics, 2008, 25(3): 26-30(in Chinese). YU Kai-ping, ZHOU Jing-jun and MIN Jing-xin et al. A contribution to study on the lift of ventilated supercavitating vehicle with low Froude number[J]. Journal of Fluids Engineering, 2010, 132(11): 1-7.

[11]

[12]

[13]

[14]

[15]

KULKARNI S., PRATAP R. Studies on the dynamics of a supercavitating projectile[J]. Applied Mathematical Modelling, 2004, 24(2): 113-129. KIRSCHNER I. N., ULMAN J. S. simplified dynamical system analysis of supercavitating high-speed bodies[C]. Fifth International Symposium on Cavitation. Osaka, Japan, 2003, No. Cav03-OS-7-005. ZHOU Jing-jun, YU Kai-ping and YANG Ming. Numerical simulation on lift of aft section of supercavitating bodies based on homogeneous multiphase model[J]. Chinese Journal of Hydrodynamics, 2010, 25(1): 113-118(in Chinese). FENG Guang, CHU Xue-sen and TAO Can-hui et al. Experimental study of supercaviting flow around body swinging at two degree of freedom[J]. Journal of Ship Mechanics, 2010, 14(6): 633-640(in Chinese). YI Shu-qun, ZHANG Ming-hui and ZHOU Jian-wei et al. Experimental research about the effects of attack angle on supercavitation of restrained model during axial accelerating[J]. Chinese Journal of Hydrodynamics, 2010, 25(3): 292-298(in Chinese).