|
Items in cart [0] |
Questions? Call 888-624-8373 |
| Copyright © 2009. National Academy of Sciences. All rights reserved. Terms of Use and Privacy Statement |
Below are the first 10 and last 10 pages of uncorrected machine-read text (when available) of this chapter, followed by the top 30 algorithmically extracted key phrases from the chapter as a whole.
Intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text on the opening pages of each chapter.
Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.
Do not use for reproduction, copying, pasting, or reading; exclusively for search engines.
OCR for page 667
Finite Difference Analysis of Unsteady Cavitation
on a Two-Dimensional Hydrofoil
A. Kubota, H. Kato and H. Yamaguchi
University of Tokyo
Tokyo, Japan
Abstract
The authors present a new cavity
model, which is named a Bubble Two-phase
Flow (BTF) model, to explain the
interaction between viscous fluid and
bubble dynamics. This BTF cavity model
treats the inside and outside of a
cavity as one continuum by regarding the
cavity as a compressible viscous fluid
whose density varies widely. Navier-
Stokes equations including cavitation
bubble clusters are solved in finite
difference form by a time-marching
scheme. The growth and collapse of a
bubble cluster is given by a modified
Rayleigh's equation. Computation was
made on a two-dimensional flow field
around a hydrofoil NACA0015 at an angle
of attack of 8.0 deg. The Reynolds
number was 3x105. The computational
results showed the occurrence of the
attached cavity at the foil leading
edge. Furthermore, the present results
predicted the generation of cavitation
clouds with large-scale vortices. The
newly proposed BTF cavity model is very
flexible and promising.
1. Introduction
When the pressure of a liquid is
reduced at constant temperature by
either static or dynamic means, a state
is reached at which vapor or vapor-
filled bubbles exist. This phenomenon is
called cavitation [1]. Cavitation
results in performance deterioration of
hydraulic machinery, generation of
noise, vibration and erosion. This is
because cavitation is a dynamic
phenomenon, as it concerns itself with
the growth and collapse of cavities. In
spite of many excellent studies, the
actual structure of cavitation is not
667
yet fully understood. In order to
accurately predict the generation of
noise or onset of erosion, it is
particularly important to illuminate the
unsteady structure of cavitation.
Vortex cavitation is often observed
in the flow downstream of attached
cavitation. It is caused by vorticity
shed into the flow field just downstream
of the cavity. Such vortex cavitation
generates a large cavitation cloud under
certain conditions. The vortex
cavitation impinges on the body and its
subsequent collapse results in erosion
[2]. In previous work, the authors
performed an experimental investigation
of the unsteady structure (velocity
distribution) of cloud cavitation on a
stationary two-dimensional hydrofoil
using a conditional sampling technique
[3]. It was found that the cloud
cavitation observed in the experiment
was a large-scale vortex with many small
cavitation bubbles. Consequently, the
importance of the interaction between
large-scale coherent vortices in the
flow field and cavitation bubbles was
recognized.
Much theoretical work also has been
done in order to obtain a better
understanding of the physics of
cavitation phenomenon. Researchers have
continually developed new models of
cavitating flow based mainly on the
assumption that the flow is irrotational
(inviscid). For example, Tulin [4][5]
proposed small-perturbation (linearized)
theory for the case of a supercavitating
hydrofoil. He has treated the cavity as
a single vapor film and assumed that the
pressure inside the cavity is constant.
This is a sort of macroscopic analysis
of cavitation. Using this model, many
researchers have gradually improved the
OCR for page 668
calculation methods [6][7][8][9][10][11].
The single vapor film model is now well
established. It can also predict
macroscopic cavity characteristics
fairly well. However, they all dealt
with only steady cavitation. Furness and
Hutton [12] treated the case of an
unsteady attached cavity on a stationary
two-dimensional body by the singularity
method. This calculation result showed
unsteadiness of the cavity surface and a
reentrant jet [1]. Their methods,
however, could not predict the behavior
of a detached cavity after the attached
cavitation splits into two parts. Tulin
and Hsu [13] and van Houten [14] have
solved the unsteady cavity problem on a
periodically oscillating hydrofoil.
Their method also could not predict the
generation of detached cavitation
clouds. This is because of a limitation
of the cavity model, which treats the
cavity as a single vapor film where
pressure is constant. Therefore, a new
model of cavitation is required to
theoretically study the breakoffs of
attached cavitation and cavitation
clouds.
How to model the cavity trailing
edge, where the cavity collapses, is the
most difficult problem for the above-
mentioned single film, constant pressure
cavity model. When one observes actual
cavitation, it is found that the sheet
cavity splits into minute bubbles with
vortices in the end region. Then the
bubbles collapse. One can often observe
many vortex cavities in this region even
if the sheet cavity is stable. Van
Wijngaaden [15][16], March [17], Chahine
and Lie [18], d'Agostino and Brennen
[19] and others have studied the
dynamics of bubble clusters. The bubble
cluster is a kind of microscopic
modeling of cavitation. However, these
studies treated only a cluster of
collapsing bubbles under given
conditions. Hence, they could not answer
how the unsteady attached cavity sheds
cavitation clouds.
Highly vertical fluid motion such as
a cavitation cloud is often observed
downstream of a cavity. Experimental
observation shows a close relationship
between large-scale coherent vortex and
cavitation [3]. Therefore, it is
necessary to develop a theoretical
elucidation of the mechanism that
generates the large-scale vortex
structure. There is a strong interaction
between the large-scale vortex and
cavitation. The occurrence of attached
cavitation yields boundary layer
separation. The separated shear layer
rolls up, thus turning into a large-
scale vortex [20]. On the contrary, the
large-scale vortex yields a low pressure
region at its center. In the low
pressure region, bubbles grow and
remain. Existing cavity models are
powerless to explain the nonlinear
vortex dynamics of the above problem
since they assume inviscid flow.
In this paper, the authors will
propose a new cavity model that can
explain the interactions between
vortices and bubbles. The authors call
this new model a Bubble Two-phase Flow
(BTF) cavity model. In a macroscopic
(coarse-grained) view, this model treats
the cavity flow field phenomenologically
as a compressible viscous fluid whose
density varies greatly. It can treat the
inside and outside of the cavitation as
a single continuum and hence it can
express the detached cavitation clouds.
In a microscopic view, this model treats
cavitation structurally as bubble
clusters. By coupling these two views,
the BTF cavity model can clarify the
nonlinear interaction between
macroscopic vortex motion and
microscopic bubble dynamics. The authors
have developed a program code SACT-III
(Solution Algorithm for Cavitation and
Turbulence, version III) to solve the
BTF model equations using the finite
difference method. They will apply this
program to cavity flow around a
stationary two-dimensional hydrofoil in
order to verify its ability.
2. Formulation
An attached cavity, which is formed
at the leading edge of a cavitating
hydrofoil, collapses at the mid-portion
of the hydrofoil. It is also well known
that the attached cavity oscillates
cyclically within a certain range of
cavitation number. This unsteady cavity
sheds a cavitation cloud in each cycle
[3]. The front part of the attached
cavity is a film of vapor where pressure
is constant. At its rear part, the vapor
film splits up into tiny bubbles. A
large-scale vortex caused by the cavity
rolls up the bubbles, thus generating a
cavitation cloud. There are two cavity
types in a microscopic view. As shown in
Figure 1, these are the vapor film and
the bubble cluster (with vortices).
668
OCR for page 669
/ / . ;.; ~.\ CAVITY SURFACE
// N[~
( ( ~ ° ~ o of ~ ) ~
At_
-
I;igure 1 Two microscopic views of cavitation,
vapor f i lm types (above) and
bubble cluster types (below).
There is, however, a third type of
cavitation. This is traveling cavitation
or bubble cavitation [1]. This
phenomenon is not treated in this paper.
2.1 Macroscopic Modeling
In the macroscopic view, the Bubble
Two-phase Flow (BTF) model treats the
inside and outside of the cavitation as
one continuum. This is because it
regards the cavity flow field as a
compressible viscous fluid whose density
varies greatly. According to this
phenomenological modeling, contour lines
of void fraction (volume fraction of
cavities) express the shape of the
cavity as shown in Figure 2. Governing
equations of the macroscopic flow field
are as shown below.
The equation of continuity is
50e+ Ad pv)=o , (~) p=(l-f g) AL
669
~ I CROSCOP I C V I EW\
At\
1~ ,
\ LOCAL HOMOGENEOUS MODEL /
\ (LHM) ~
CAVITATION CLOUD
_
MACROSCOP I C V I EW
COMPRESSIBLE VISCOUS FLOW
l~igure 2 Modelling concept of Bubble Two-phase
Flow (BTF) cavi ty model.
where t, pv, v and p are time, mass flux
vector(pu, pv, pw), velocity vector
(u,v,w) and density of the mixture,
respectively.
The conservation equation for
momentum: Navier-Stokes equation is
8(G~+V(pvv)=-VP+Rlep{v2v~lv(v·v)} ,(2)
where P is the pressure in the mixture,
~ is the viscosity of the mixture, and
Re is the Reynolds number. This equation
is in conservation form [21]. The
nondimensionalized quantities based on
the uniform flow velocity and a
reference length have been employed in
the above two equations.
The BTF cavity model assumes that a
fluid of variable density replaces the
water-vapor mixture. The density of the
water containing bubbles (mixture) is
defined as follows:
(3)
OCR for page 670
where PL is the water density and fg is
the local void fraction. The mass and
momentum of the vapor are ignored, since
they are small compared with those of
the liquid. The actual ratio comparing
the density of vapor to that of liquid
is of the order 10-4. The change of
liquid mass due to the phase change is
also ignored.
The viscosity of the mixture is
assumed to be as follows:
P=(1-fg)PL+fgI`G ' (4)
where ILL is the water viscosity and FIG
is the vapor viscosity.
2.2 Microsconin Mod ~ 1 i n ~
- Local Homogeneous Model -
To calculate the macroscopic flow
field, it is necessary to know the local
void fraction function fg(t,x,y~z). The
greatest problem is to develop a model
that gives the relationship between the
flow field condition and the void
fraction. The present BTF model treats
cavitation microscopically as bubble
clusters. This is because bubbles play
important roles in the cavity inception
[22]. Furthermore, one of the main
purposes of the BTF model is to study
the mechanism of the cavity collapse.
The sheet-type cavity splits up into
tiny bubbles there. This structural
microscopic model cannot be basically
applied to the vapor film type
cavitation. However, the type of the
microscopic model has less effect on the
macroscopic model when the void fraction
is high. The present BTF model
introduces a Local Homogeneous Model
(Lam), which is a sort of Mean Field
Approximation (MFA), for simplicity. It
treats the cavity as a local homogeneous
cluster of spherical bubbles as shown in
Figure 2. Bubble number density and a
typical radius are assumed locally. This
typical bubble radius is obtained from
the growth-collapse equation of the
bubble cluster. This equation is derived
from the growth-collapse equation of one
spherical bubble (Rayleigh's equation).
The LHM gives the local void fraction
fg by coupling the bubble density and
the typical bubble radius as follows:
~ =n.4~R3 ,
where n is the bubble number density and
R is the typical bubble radius (see
Figure 3). In this paper, the bubble
density is assumed to be constant all
over the computational domain though
real cavity flows have a distributed
bubble density. This is because it is
too difficult to formulate coalescence
and fragmentation of the bubbles.
, ~ — if. _
~ I.. ..~... .N
... . ..., . ..
.. ....... ..... ...... ...
,) age:-.. \ //, ''.:
f = n 4lrR3
p = ( l~fg ) PI,
Figure 3 Local homogeneous microscopic cavity
model (LHM).
The LHM assumes that cavities form in
the shape of spherical bubbles. The
bubbles remain separate and distant
enough from each other so that their
shapes remain spherical. Interaction
between bubbles occurs through the local
pressure that develops in the liquid as
bubbles grow.
Lord Rayleigh originally derived the
equation of radial motion (growth and
collapse) of an isolated spherical
bubble in a homogeneous infinite medium
[23]. This equation is widely known as
Rayleigh's equation. It takes the
following form, neglecting the effect of
surface tension and viscous damping,
Rd2R+3(dR)2=Pv P (6)
where Pv is the vapor pressure. In this
study, the vapor pressure is assumed
constant. This is because the behavior
of the bubble is nearly isothermal and
gas inside the bubble is also ignored.
The finite-difference method was
employed in SACT-III. In this method, a
continuous domain is discretized into
finite grid points. Hence, an
(5) interaction between individual bubbles
within the grid spacing must be
considered to eliminate the effect of
the computational grid. This effect is a
sub-grid-scale (SOS) bubble interaction.
670
OCR for page 671
OCR for page 673
OCR for page 674
OCR for page 675
OCR for page 676
OCR for page 677
OCR for page 678
OCR for page 679
OCR for page 680
OCR for page 681
OCR for page 682
OCR for page 683
OCR for page 684
Representative terms from entire chapter:
cavity model
Next, let us consider the SGS bubble
interaction of the LHM, deriving
analytically the equation of motion of
the bubble cluster.
~ MA r
~~ r
Figure 4 Sub-grid-scale (SOS) bubble interaction
mode 1.
As shown in Figure 4, we consider the
influence of the other bubbles which
exist inside of the distance Or (=grid
spacing). The total velocity potential
due to the other bubbles
at the origin O is
i 1 dri 2
>: ri dtRi
(7)
when Ri<
dimensional rectangular cell version).
The results computed by SACT-II
explained cavity formation caused by
large-scale coherent vortices behind a
rectangular obstacle on a wall [24]. The
present SACT-III employs the finite
difference method in the body fitted
coordinates to solve the governing
partial differential equations given in
the preceding section. This section
explains the computational procedure and
the finite difference scheme of SACT-III.
The computational procedure is basically
the same as the Marker-and-Cell(MAC)
method [25] except for the use of a
regular mesh system, instead of the
staggerd mesh system.
3.1 Quasi-Poisson Equation for Pressure
By taking the divergence of the
Navier-Stokes equation(2), the following
Poisson equation for pressure is given:
where
~(~v,v)=-V[V(pvv)+Ilep{V2v+1V(V v)}] .(17)
Substituting the continuity equation(1)
into Equation(17),
2 2
V P= ~ 2p+~(pv,v) . (
at
From Equations(3) and (5), which
represent the basic assumptions of the
LEM,
v2p=_V(~ )+~(~v,v)
=-~t{V(~)}~(~v,v) ,
(16)
18)
~=(l-n.3~R3)~L . (19)
To simplify the LHM, n=constant is
assumed in this study as mentioned in
section 2.2. By differentiating
Equation(19) twice with respect to t, we
have the following equation:
^=-~L4n~{R a 2+2R(~t) } . (20)
and
From the LHM's equation of motion(13),
~ R=s7(pv,v,~R,R)+~(P) , (21)
where
( pv,v, aR,R)
=-{2(v.V)3R+(v.V)(v.V)R}
_+4~Ar2nR BR 2 (22)
2 2 { t+(v.V)R}
(1+2~^r nR)R ~
_ 2~Ar2R2 {8n+(v~v)n}{8nt+(v~~)R} ,
(l+2~r2nR)R ~t ~
672
P. -P
*(P)= 2 (23)
(1+2~Ar nR)RPL
If n=constant, the last term of the
right-hand side becomes zero in
Equation(22).
Substituting Equations(20) and (21)
into (18), we obtain the quasi-Poisson
equation for pressure including the
motion of the bubble cluster as follows:
V2P+~'(P)=~'(pv,v,~R,R)+~(pv,v) ,(24)
where
~'(P)=-pL4n~R ·~(P) , (25)
a~ (pv,v,~3R,R)
PL4n7tR { g! ( PV, V, a8R, R ) + 2 ( BR ) 2 }
26)
The left-hand side in Equation(24) is
approximated by the second-order finite
differencing scheme. As a consequence,
simultaneous equations of pressure P are
obtained if the right-hand side is given
in Equation(24). SACT-III solves the
simultaneous equations derived from
Equation(24) with a point successive
relaxation method. Equation(24) is
equivalent to the normal MAC method's
Poisson equation of incompressible flow
[25] when
~' =97' =0
They are always set to zero for non-
cavitating conditions. If the mixture is
filled with liquid, Equation(13) cannot
be solved since the bubble radius R
becomes zero. If the mixture is filled
with vapor, Equation(2) cannot be solved
since density of the mixture becomes
zero. Hence, when the void fraction fg
isless than fgmin(>0 °) or more than
f~max (<1.0), the bubble radius becomes
f~xed. Then, ~ and a'' are also set to
zero. In the following computations,
fgmin and fgmax are set as follows:
i =n 4~Ro3 '
f =0.95 ,
where Ro is the initial bubble radius.
3.2 Numerical Methods
Equations(2) and (13) are time-
integrated with and Euler explicit
scheme using the value of pressure P
obtained by solving Equation(24). To
solve a high-Reynolds-number flow, each
nonlinear term in Equation(2), for
example ~3x(~uv), was approximated with
the fourth-order centered finite-
differencing scheme with the fourth-
derivative term:
4
ul. a4(pv)^x3 . (27)
ax
The fourth-derivative term plays an
important role in stabilizing the
calculation. Physically, fourth-
derivative term means shorter-range
diffusion compared with the second-
derivative viscous term [26]. The
fourth-derivative term consequently
stabilizes the computation without
decreasing the Reynolds number and
introducing any turbulent models. The
universal availability of this fourth-
Downst ream
derivative term nas not been ascertained
yet for the turbulent flow calculation. Boundary
However, SACT-III introduces no
turbulent model since there exists no
established one for the two-phase
(cavity) flow at present.
All the other space differential
terms in Equations(2) and (13) are
approximated with the second-order
centered differencing scheme.
Equation(13), however, has no spatial
diffusive term for bubble radius R and
its time-derivative. The second-
derivative term is accordingly added in
equation (13) to eliminate the
instability of the nonlinear terms. For
example,
lu1 ·B Rex
ax
is added to u8R . This term means
diffusion.
To compute the high-Reynolds-number
flow around a body of arbitrary shape,
it is convenient to use body-fitted
coordinates through coordinate
transformation. Figure 5 shows the grid
system for the present problem of flow
around a two-dimensional hydrofoil. This
system is called a C-type grid [27]. The
connected physical (x,y,z) domain around
the hydrofoil is mapped onto the
rectangular computational (\,q,;)
domain. Here the pair of planes forming
the branch cut are both on the same
plane of the transformed region. The
surface of the body is also mapped on
the same plane with the branch cut.
A regular mesh system is employed.
Velocities, pressure and bubble radius
are given on the grid points. As shown
in Figure 5, the uniform flow boundary
PLOW ~
~ =
PHYS I CAL DOMA I N
Outer-flow Boundary Side Boundary
Downstream
Boundary
... . j
~ Wall Boundary ~
L ~ _ _ _ ~
Branch Cut
TRANSFORMED DOMA I N
Fi Sure 5 C-grid system around a two-dimensional
hydrofoi 1.
conditions are imposed at the outer-flow
boundary. These are:
u=l, v=w=O,
Bu_8v_8w_O
B~ B9~ all- ,
R=R , ~R=0 P=0 .
O at '
(29)
28) At the downstream and side boundaries,
the boundary conditions of the zero-th
order (zero-gradient) extrapolation are
imposed. At the branch cut boundaries,
the periodic boundary conditions are
imposed. At the wall boundary, the
following boundary conditions are
imposed.
u=v=w= 0,
BP=BR=0
(30)
First order (linear) extrapolation in
(t,Q,4) domain is used for the
velocities to calculate the nonlinear
terms of Equation(2). The terms are
approximated by the fourth-order
centered differencing scheme with the
fourth-derivative term.
4. Computational Results and Discussion
4.1 Condition of Computation
673
A hydrofoil section with a simple
mathematical configuration, NACA0015
[28], was chosen for the computation.
The computation was performed at angles
of attack ~ of 0.0, 8.0 and 20.0 degrees
[29]. However, this paper discusses only
the results at ~ =0.0 and 8.0 degrees.
The Reynolds number Re was 3x105 in all
the computations. It is based on the
uniform flow velocity and chord length
of the hydrofoil The viscosity ratio
G/pL was 9.12x10-3. The computation at
~ =0.0 deg. was performed only for non-
cavitating conditions to evaluate
numerical accuracy. Experimental cP
observation at ~ =8.0 deg. shows laminar
separation without bursting near the
leading edge. For cavitating conditions,
an attached type cavity accordingly
occurs near the foil leading edge.
The hydrofoil was accelerated from
u=0 to the steady speed of 1 for T=O-1.
The cavitation number was also decreased
gradually for T=2-3 so as to compute
stably for cavitating conditions. The
time increment t was determined in each
computational step to keep the Courant
number less than 0.25. The relaxation
factor was set to be 0.8 when
Equation(24) was solved with the
successive relaxation method. However,
it was reduced to 0.3 for cavitating
conditions. The convergence condition of
the pressure computation is as follows:
max(lAPI) < 0.001
(for non-cavitating conditions)
max(lAPI) < 0.01
(for cavitating conditions),
where UP is the residue of pressure in
an iterative calculation.
Figure 6 C-grid system around NACA0015 hydrofoil
at a=O. Odeg, lOl(~)x31~ ~ x3~.
Preceeding the computations for
cavitating conditions, it is necessary
to evaluate numerical accuracy of
SACT-III for non-cavitating conditions.
The angle of attack was 0.0 degrees.
Figure 6 shows the c-grid system. The
grid was uniform in the spanwise
section. The number of grid points was
lOl(~)x31(~)x3(~). The distance between
the trailing edge and the upper or lower
boundaries was 1.2. The distance between
674
the trailing edge and the downstream
boundary was 3.0. The shape of the front
part of the outer-flow boundary was oval
as shown in Figure 6. A direct numerical
method [30] was used for the grid
generation. The minimum grid spacings at
the leading and trailing edges were
0.70xlO 4 and 1.54x10-4, respectively .
The minimum spacing was about one
twelfth of Re~0 5 at the trailing edge.
-2.0 - _
-1.0- _
· Experiment (Re=6 x 10 5)
06 Present Calculation (Re=3x 105)
Hess-Smi th Method Hi th
Boundary Layer Calculation
(Re=3 x lOs)
~o57
1
Am_
I', , , 1
`49E ,^~ ,______ ~
Figure 7 Foil surface and wake pressure,
NACl\OO 1 S; ~ =0. Odeg. ; T=2. O.
Figure 7 shows the pressure
distributions on the foil surface at
T=2. The circles show the distribution
of Cp on the foil surface. The triangles
show the distribution of Cp in the wake,
i.e., on the branch cut shown in Figure
6. The computational results are
compared with the computation using the
Hess-Smith method [31] and the
measurement at Re=6xlO5 by Izumida [32].
The Hess-Smith method is a sort of
numerical solution method of potential
flow based on the boundary element
method. The foil shape was modified by
adding the computed displacement
thickness of the foil surface boundary
layer [33]. The laminar boundary layer
was computed using Thwaite's method,
with the empirical constants derived by
Curle and Skan [34]. The length of the
separation bubble was 150 times as long
as the momentum thickness at the laminar
separation point. The Head's entrainment
method modified by Cebeci [35] predicted
the turbulent boundary layer
development. The turbulent separation
was predicted to occur when the form
factor H12 exceeds 2.1. The computed
pressure distribution agrees very well
with the others as shown in Figure 7.
Furthermore, the pressure coefficient is
0.982, which is almost equal to 1.0, at
the front stagnation point. As a
consequence, the present numerical
method has good accuracy when the grid
system is fine enough. The computed
pressure distributions disagree with the
experimental result when the coarser
grid systems are used [29].
Convective Term
_~
~ ~ DMAX=14.4
Fourth-derivative Numerical Dissipation Term
Pressure Term
Truncation Error of Pressure Term
Figure 8 Distribution of convective, fourth-
derivative, pressure and its truncation O
error terms in x-momentum equation, AL.
NACA0015; a=O. Odeg. ;Re=3. Ox 105;T=2. O. ---
Figure 8 shows the distributions of
the convective term, the fourth-
derivative term, the pressure terms and
its truncation error in the x-momentum
equation. The contour interval is 1.0.
DMAX shows the maximum. The maximum of
the fourth-derivative term is not
negligible compared with the convective
and pressure terms. However, it is
distributed only near the foil surface
around the leading and trailing edges.
This result shows that the effect of the
fourth-derivative term is local but
important to stabilize the computation.
The truncation error of the pressure
term is sufficiently small compared with
the main differences terms.
4.2 Unsteady Attached Cavity
Figure 9 shows a close-up of the grid
system around the NACA0015 hydrofoil.
The angle of attack was 8.0 deg. The
number of grid points was lOlx31x3. It
was the same as that at ~ =0.0 deg. The
minimum spacing was 0.70x10-4 at the
Figure 9 Close-up of the grid system around NACA0015
hydrofoil, a=8.0deg., 101~) X31~7~) x3~.
trailing edge. It was almost the same as
that of the grid system at ~ =0.0 deg.
Figure 10 shows the time-averaged
velocity vectors from T=2 to 4 for non-
cavitating conditions. The computation
was performed stably. The boundary layer
separates at X=0.74 on the back side.
Instantaneous velocity vectors, however,
show unsteady vortex shedding from the
foil trailing edge region. In the other
region, the flow is almost steady. No
separation occurs near the leading edge.
SeparationPoint (X=0.74)
~_~
_~` _
675
Figure 10 Time-averaged (T=2~4) velocity vectors
around NACA0015 hydrofoil, a=8.0deg.;
Re=3. 0 x 10 s.
Figure 11 shows the chordwise
distribution of boundary layer
displacement thickness on the foil
surface. The present result agrees well
with the boundary layer calculation at
the front part of the back side.
However, it cannot predicte a laminar
separation bubble because of the
insufficiency of the grid points in the
}-coordinate. The separation point is
further upstream than in the boundary
layer calculation. On the contrary, the
agreement of the separation point is
good on the face side. However, the
boundary layer is thicker near the
leading edge.
Figure 12 shows an example of the
velocity profiles using wall variables.
The present profile closely follows the
low of the wall. The present method can
consequently express the turbulent
boundary layer without any turbulent
models.
~ 1
a -
ut
20:
10
o
BACK SIDE
Present Cal.
--— Boundary Layer Ca 1.
·: Separat ion Po i nt
t: Separ at i on Bubb l e
> T~rh~ 1 An t ~ 7
FACE SIDE
However, the shape of the pressure
distribution agrees well with the
experiments. Hence, the disagreement of
pressure hardly affects the nature of
the unsteady cavitation.
////// ~
Pv- P
p
x
_ = -O. 1
GO\
L 1` ~ A: non ~ r3 = to
B: n3'T ~ r3 = 1
1 C: n-n ~ r3 = 10
*1 C ~
|B | D: n3'T ~ r3 = 100
Present Cal.
- - - Boundary Layer Ca 1. / P v- P T i ~ e t
p
·: Separation Point ~
Gin , x
0 0.5
+0. 1
3.2
|.0
t.8 _
2.'
Figure 11 Boundary layer displacement thickness .= 2.2
distribution on NAChOOlS hydrofoil. ~
ct=8.0deg. ;2e=3.0X105. ~ 2.0
_ I.~
~ l..
nut ye
'/Z/571ny++5.5
. :Present Cal. (I=57)
Of = 5. 91x 10-3
v*O= 8. 07x 10-2
Figure 12 Velocity profile in boundary layer on
NACA0015 hydrofoil using wall variables,
~ =8.0deg.;Re=3.0X 105;I=57(X=0.0827~.
While the velocity profiles are
predicted well, the pressure is not
predicted quite so well. The present
lift coefficient CL is only about 58%
that obtained by the experiment.
~ A 4 ~~ 3
3
B : n 3 7r ~ r 3 = 1
- C: non ~ r:
~ D 4 a 3
-
= 100
0 1 2 3 ~ 5 e 7 8 ~ 10
Tise t
Figure 13 Influence of the SGS bubble interaction
model on bubble collapse (above) and
growth (below).
Before analyzing the cavitating flow
field, we should check the effect of the
SGS bubble interaction model. Figure 13
shows the effect of the SGS bubble
interaction when the ambient pressure P
has changed stepwise with the amplitude
of +(Pv~P) /PL=0 .1. n(4/3)= ~r3 is the
number of bubbles within the grid scale
a r. The convective terms were neglected
in Equation(13). The Runge-Kutta-Gill
method was used to solve the
of differential Equation(13). The bubble
radius is nondimensionalized using its
676
initial value. The case A is the
calculation of isolated bubbles, i.e.,
with no-bubble interaction effects. As
shown in this figure, presence of other
bubbles delays growth and collapse of
the bubbles. Fujikawa et al. [36]
carried out a theoretical analysis on
the interaction between two bubbles.
Their result shows that the collapse is
delayed when the two bubbles have the
same radius. The present result shows
the same tendency.
In the following computation, the
grid scale of the SGS bubble interaction
model is assumed as follows:
r=(~)2
where
g' =X.,y~-x~y em
(31)
(two-dimensional Jacobian).
This is because the flow structure is
two-dimensional. Accordingly, the SGS
bubble interaction effect is independent
of the grid spacing in the spanwise
direction.
/; Vertical Exaggeration of 4:1
Figure 14 Void fraction contours around NACAOOlS
hydrofo i l, cr = 1. 5; ~ =8. Odeg. ; Re=3. 0 x 10 5 -
the contour interval is 0.1 except for
the outermost line.
Figure 14 shows contour lines of void
fraction at ~ =1.5. The initial bubble
radius Ro and the bubble number density
n are 4x10-4 and lx1O6, respectively. In
this paper, we define a cavity as a
region where void fraction is more than
0.1. The bold contour lines are those of
fg=O.1 in this figure. Contour lines are
drawn at intervals of O.1 except for the
most outer line of fg=O.01. As seen in
this figure, thin steady attached cavity
incepts smoothly.
Figure 15 shows the time series
contour lines of void fraction at
0.1 0.01
W\_
~F~
\
T = ~ OF ~
I; igure IS Void fraction contours around NACA0015
hydrofoi l, ~ =1. 2; a =8. Odeg. ; Re=3. 0 x 10 5;
the contour interval is 0.1 except for
the outermost line.
~ =1.2. The void fraction is more than
O.9 at the center of the cavity. The
rear portion of the cavity oscillates
cyclically. Not only does the cavity
length change but the cavity itself
rises up at its rear part. The unsteady
characteristics computed here agree with
the experimental observations of sheet-
type cavitation [32]. It is therefore
concluded that the BTF cavity model can
express the features of sheet-type
cavitation beyond its microscopic model,
which is essentially suitable to
the bubble cluster flow.
Figure 16 shows the time-averaged
pressure distribution on the foil
677
T I ME-AVERAGED
CP-DISTRIBUTION
CP (T=3. 0~ 5. 0)
NIiCt0015; a =8O, ~ =1. 2
O
.
_1 n
n
1.0
CL = 0. 5472
CD = 0. 0693
........
it_ ~
L. Ed H.C -A—-
~ - ' — ~
0.0
~_0.
· P.~v i t v A rea
~ - ~ ~` ,J,,,_~ ~
Figure 16 Time-averaged pressure distribution and
void fraction contours around NACA0015 I;
hydrofoil, a=1.2;a=8.0deg. ;Re=3.0X105,
the contour interval is 0.1 except for
the outermost line.
surface and void fraction contour lines.
The foil surface pressure where cavity
exists is almost constant and equal to
the vapor pressure (Cp=-1.2). However, a
small pressure peak exists at the front
part of the cavity. The pressure
distribution is similar to the
calculated result by a nonlinear free-
stream line theory [10]. Furthermore,
the time-averaged cavity shape is
similar to the experimental observation
of sheet-type cavity. These facts also
suggest that the present BTF model can
be applied to sheet type cavitation
beyond its micro structure limitation.
Vertical Exaggeration of 4:1
Contour Interval of 0.1
~\ _ _ _
o
°L
i gu re 17 T i me-averaged pres sure coef f i c i ent contours
and velocity vectors around NI\CAOOlS
hydrofoi l, a =1. 2; a =8. Odeg. ;Re=3. 0 x 105,
the contour interval is 0.1, the bold broken
lines are void fraction contours of 0.1.
Y r
o.l r o.lo
~ I MU = U
O t.0 0 1.0
OISPLRCEt1ENT Tl~lct~rlEss = 0.0288q DlSPLPCEMENr THICKNESS = 0.02112
hlOMft`TUb TFilCht~fSS ~ O.00q97 MOMfNTUM THICKNESS ~ 0.00557
FOnM FnCTORll112) = S.8O FORM FRCrOlltH12} = 3.79
o. 10
~ r
O 1=67 1 0.10 1=69
X=O. 4812 X=O. 5832 1
o 1.0 o 1.0
Figure 17 shows the time-averaged
pressure contours and flow velocity
vectors. The contour line of fg=O.1 (the
bold broken line) agrees approximately
with that of Cp=-1.2. This is the reason
why the cavity is defined as a region
where void fraction is more than 0.1.
The reverse flow is observed clearly in a=1.2
the neighborhood of the end of the ~ - No-cavitation
Bold lines in Figure 18 show time-
averaged velocity profiles in the cavity
wake boundary layer along the
1-coordinate. Fine ones are profiles
DISPtnCEh1Et1T li1ICKt1ESS = 0.01~! OlSPLnCEMENT TFfICKNE'SS = 0.01507
MOl1ENTUM Tt1ICK'lESS · 0.00771 MOMENTUM THICKNESS ~ 0.00813
FORt1 fnCTOnt'll2, = t.f39 FORH F0CTOR{~12) = I.f3S
I;igure 18 Time-averaged velocity profiles in cavity
wake region, a =1. 2; a =8. Odeg. ;Re=3. 0 x 105.
fine lines show the boundary layer velocity
profiles for non-cavitating conditions.
678
~ 1
r=1.2
O. 03
O . O
O. O
to
No-cavitation
////
/
/
/ I
/1
/ ,{
l
Figure 19 Comparison of boundary layer displacement
thickness distribution on the back side on
4.3 Cloud Cavitation
Figure 20 shows void fraction contour
lines at ~ =1.0 with a time increment of
0.2. The initial bubble radius Ro and
the bubble number density n are 1x10-3
and 1x106, respectively. For this
condition, unsteady cavities
continuously grow and collapse. The
highly distorted attached cavity sheds
cavitation clouds cyclically (T=5.9 and
7.1), which soon collapse. This
phenomenon agrees well with many
experimental observations
[32][39][40][41][42][43]. Figure 21
shows an example of the photographs of
cloud cavitation [29]. The position of
the cavity break off point agrees well
with the computational result.
Figure 22 shows velocity vectors
around the foil. Overlaid bold broken
lines are void fraction contours of 0.1.
As mentioned before, they show
instantaneous cavity shapes. As shown in
this figure, the unsteady attached
cavity sheds not only cavitation clouds
but also vortices (see marks A, B and
C). The experimental result has
confirmed such vortex shedding phenomena
[3]. The position of cavitation clouds,
however, does not agree very well with
That of shedding vortices. There are two
~ - - - - - - - 5 possible explanations for this
NACA0015 hydrofoil. a =8.0deg.;Re=3.0X10 .
for non-cavitating conditions. As shown
in the section of 1=61 (X=0.2115), the
flow inside the cavity is quite slow
except near the foil surface. At the
section of 1=64 (X=0.3375), strong
reverse flow occurs. At 1=67 and 69
(X=0.4812 and 0.5832), the flow
reattaches and inflection points exist
in the velocity profiles. The
measurement result downstream of a
stable sheet cavity shows similar
inflection points in the velocity
profile [37][38]. As shown in this
figure, the generation of cavity causes
an increase in the boundary layer
thickness behind it. Figure 19 shows the
comparison of the chordwise
distributions of the boundary layer
displacement thickness for cavitating
and non-cavitating conditions. For
cavitating conditions, the displacement
thickness decreases at the cavity
collapsing region. It becomes minimum
near the mid-chord, then it begins to
increase again. This is also the same
tendency that the experimental results
show [38].
discrepancy. One is the phase delay of
the cavity growth. The other is that the
present cavity model does not consider
nuclei convection towards the center of
vortices.
Figure 23 shows a close-up of the
velocity vectors around the cavity. This
figure elucidates the mechanism of
cavitation cloud shedding. At T=5.5, a
new separation vortex occurs at the
cavity leading edge. Then it induces the
flow toward the foil surface (T=5.7).
Fluid density and pressure on the foil
surface increase due to the impinging
flow. It causes the cavity to break and
tear off (T=5.9, separation of the
cavitation cloud). The impinging flow
turns into a jet along the foil surface.
The jet sweeps away the cavitation cloud
(T=6.1). This is the scenario of the
generation of a cavitation cloud.
5. Concluding Remarks
In this study, the authors presented
a new modeling concept of cavitation
called BTF (Bubble Two-phase Flow). In a
679
~-
T = 5.
T = 5. 5 \
it_
-I=
~~\
_~-U.01
-
T = 6. 9
-
CavitationCloud
~°~N
\
T = 6. 1 ~
Figure 20 Void fraction contours around NACA0015 hydrofoil, a=1 0; a=8. Odeg. ;Re=3.0x 105,
the contour interval is 0. 1 except for the outermost line.
macroscopic view, this new cavity model
treats the inside and outside of a
cavity as one continuum. That is, it
regards the cavitating flow field
phenomenologically as a compressible
viscous fluid whose density varies
greatly. Contour lines of void fraction
can express the cavity shape. In a
microscopic view, a simple LHM (Local
Homogeneous Model) is introduced. This
is a kind of Mean Field Approximation.
This structural microscopic model treats
a cavity as a locally homogeneous bubble
cluster. Assuming bubble density and a
typical bubble radius, a local void
fraction function is given. The BTF
cavity model is significant in the
following points: (l) The BTF cavity
model can investigate the nonlinear
interaction between large-scale vortices
and cavitation bubbles, (2) The BTF
cavity model can consider the effects of
bubble nuclei on cavitation inception,
680
N: -,,,,! _
Hi_
;
Figure 21 Cavity appearance on NACA0015 hydrofoil,
~ =1. 30; a=8. Odeg. ;Re=3. Ox 105.
(3) The BTF cavity model can express
unsteady characteristics of cavitation.
The BTF cavity model, therefore,
includes three essential factors for
cavitation. Those factors are pressure,
nuclei and time. By examining the
- - -
- - -
~:
~ - - - - - -
~ - - -. - - - -
indispensable to improve the computation
scheme of SACT-III to obtain qualitative
agreement with experimental results. The
accuracy improvement of numerical
computations is one of the most
important subjects at present.
This research is partly supported by
the Grant-in-Aid for Co-operative
Research of the Ministry of Education,
Science and Culture. All the
computations in this study were carried
out on the HITAC M-680H computer at the
Computer Center of the University of
Tokyo. The authors express cordial
thanks to Mr. M. Maeda and Mr. M.
Miyanaga for their help. They also
express their appreciation to Prof. H.
Miyata, Dr. Y. Kodama, Dr. N. Baba,
Prof. H. Ohtsubo and Mr. J. Butler for
their encouragement and support.
References
1. Knapp,R.T. et al., Cavitation,
McGraw-Hill (1970).
2. Hutton 9 S.P., Proc.Int.Symp.Cavitation,
Vol.1, Sendai, Japan, pp.21-29 (1986).
3. Kubota,A.et al., J. Fluids Eng.
Trans of ASME, Vol.111(6) (1989).
4. Tulin, M.P., Proc.NPL Symp. Cavil
Hydro., Paper No.16, pp.1-19 (1955).
5. Tulin, M.P., J. Ship Research,
Vol.7, No.3, pp.16-37(1964).
6. Wu, T.Y., J. Fluid Mech., Vol.13,
Part 2, pp.161-181 (1962).
7. Nishiyama, T. et al.,Tech. Rep.Tohoku
,_,_,_~=,=,_, University, Vol.34, pp.123-139 (1969).
_ 8. Furuya, O., J. Fluid Eng., Trans.
ASME, Vol.71, pp.339-359 (1975).
9. Furuya, O., IAHR 10th Symp., Tokyo
Fi gure 22 Veloci ty vectors around NACAOO15 hydrofoi 1, pp.221-241 (1980).
a=l.O;a=8.0deg.;Re=~.OXl05,the bold brokenlo. Yamaguchi, H. et al.,Conf.on Cavil
Iines are void fraction contours of 0.1. Edinburgh, IME, pp.167-174 (1983).
11. Lemonnier, H. et al., J. Fluid
Mech., Vol.195, pp.557-580 (1988).
12. Furness, R.A. et al., J. Fluid Eng.,
Trans. ASME, pp.515-522 (1975).
13. Tulin T.P. et al., 13th Symp.on Naval
Hydr., Tokyo, pp.107-131 (1980).
14. van Houten, R.J., 14th Symp.on Naval
Hydr., Session V, pp.109-158 (1982).
15. van WiJngaaden,L.,Proc.llth Int.Cong
Appl. Mech., Springer, pp.854-861 (1964)
16. van WiJngaaden, L., J. Fluid Mech.,
Vol.33, Part 3, pp.465-474 (1968).
17. M0rch, K.A.,Cavitation and Polyphase
Flow Forum - 1981, pp.1-10 (1981).
18. Chahine, G.L. et al., J. Fluid Mech.
Vol.156, pp.257-279 (1985).
19. d'Agostino, L. et al.,J, Fluid Mech.
Vol.199, pp.155-176 (1989).
computational and experimental results
carefully, the BTF cavity model was
proven to be useful to theoretically
study cavitation characteristics.
Complicated interactions, in particular,
between large-scale vortices and bubble
dynamics were clarified. However,
further effort must be devoted to
improvement of the microscopic cavity
model for development of the BTF cavity
model. For instance, the pressure
gradient of the vortex cavitation cloud
attracts bubbles towards its center. We
must therefore take account of the
convection of bubbles and the slip
between bubbles and liquid to predict
the behavior of cavitation clouds more
accurately. Moreover, it is also
681
~l ~ ~
f —
~'
ol
_l
~.0
\
L 1 ''' ~
1.01} 1 ~ ~ , -
~ /, .
~'
W\, \,,".
~'.,,
, , , X/C
0 0.1 0.2 0.3 0.4
~ / ~ , `\ _
'_~ ~ _
— ~/ - _ \
4-_ ~ ~ \1_
\.
~ _
,. _ _ 74\\ _
/71 ~ _ _, ~Oo
\
\
,`
.0t
,.90
~ ~ I ' ' X/C
0 O.1 0.2 0.3 0.4
Vertical Exaggeration of 4:1
Figure 23 Close-ups of velocity vectors around NACAOO15 hydrofoil, a=1.0; a=8. Odeg. ;Re=3.0X lOs,
the bold broken lines are void fraction contours of 0.1.
20. Kiya, M. et al., J.Fluid Mech.,
Vol.154, pp.463-491 (1985).
21. Roache,P.J., "Computational Fluid
Dynamics", Hermosa Publishers(1976).
22. Kodama, Y. et al., J. Fluids Eng.,
Trans.of ASME, Vol.103(4),pp.557 (1981).
23. Lamb, H., Dover Pub.,
"Hydrodynamics", pp.112 (1932).
24. Kubota, A. et al.,Theoretical and
Applied Mechanics, Vol.36, University
of Tokyo Press, pp.93-100 (1988).
25. Harlow, F.H. et al., Phys. Fluids,
Vol.8, No.12, pp.2182 (1965).
26. Kawamura, T. et al., AIAA paper
84-0341 (1984).
27. Thompson,J.F. et al.,"Numerical Grid
Generation",Elsevier Science Pub.(1985).
28. Abbott, I.H. et al.,"Theory of Wing
Sections", Dover Publications (1958).
29. Kubota, A., Doctral Dissertation,
Department of Naval Architecture,
University of Tokyo (1988).
30. Kodama, Y., J. Soc. Naval Architects
of Japan, Vol.164, pp.1-8 (1988).
31. Hess, J.L. et al., Progress in Aero.
Science, Pergamon Press, Vol.8 (1967).
32. Izumida, Y. et al, Proc. 10th IAHR
Symp, Tokyo, pp.169-181 (1980).
33. Yamaguchi, H. et al., J. Soc. Naval
Architects of Japan, Vol.164, pp.28-42
(1988).
34. Curle, N. et al., Aero. Quart.,
Vol.8, pp.257 (1957).
35. Cebeci, T. et al.,"Momentum Transfer
in Boundary Layer",Hemisphere Co(1977).
36. Fujikawa, S.et al., Proc. Int. Symp.
on Cavitation, Vol.1, Sendai, Japan,
pp.55-60 (1986).
37. Yamaguchi, H. et al., LDV and Hot
Wire/Film Anemometry, pp.29-38 (1985).
38. Kato, H. et al., Proc. Vol.2 18th
ITTC, Kobe, pp.433-437 (1987).
39. Kermeen, R.W., Hydrodynamic Lab.,
CALTEC, Pasadena, California, Report
No. 47-5 (1956).
40. Wade,R.B. et al.,J.Fluid Eng.Trans.
ASME, Vol.88, No.1, pp.273-283 (1966).
41. Alexander, A.J., Conf. on Cavitation,
Edinburgh, IME, pp.27-35 (1974).
42. Shen, Y.T. et al., 12th Symp. on Naval
Hydro., Washington, pp.470-493 (1978).
43. Franc, J.P. et al., J. Fluid Mech.,
Vol.154, pp.63-90 (1985).
682
DISCUSSION
by F. Stern
I would like to congratulate
the authors
on a very ~ nteresting paper which appears to
present a new approach for unsteady
cavitation. The authors discuss the fact that
cloud cavitation is often periodic. In my own
work on uns teady cavi tat ion
(Stern,F.,"Comparison of Computational and
Experimental Unsteady Cavitation on a Pitching
Foil", J. of Fluids Engineering, Vol.111,
No.3, September 1989, pp.290-299), close
correlation was shown between the experimental
cloud-cavitation shedding frequency and the
predicted cavity natural frequency. Do the
authors' results provide an estimate for the
cloud-cavitation shedding frequency and how
does it compare with the experimental value?
Author's Reply
The authors would l ike to thank Prof.
Stern for his valuable discussion.
F.
The computed time period of cavitation
cloud shedding, which is nondimensionalized
based on the uniform f low velocity and the
chord length of the hydrofoi l, is about 1.4 at
a c a v i t a t i o n n u m b e r o f 1 . 0 . T h e
683
dimensionalized shedding period is 1.4x(chord
length)/(uniform flow velocity). It becomes
0.0117sec (85.7Hz) when the chord length and
the uniform flow velocity are 50 mm and
6.0m/s, respectively (the experimental
condition). The authors did not measure the
time period of the cavitation cloud in their
experiments. However, the authors' similar
experiment [ Al ] showed that the
nondimensionalized time period of cloud
cavitation shedding was 2.1. This value is
very close to the present computed result.
The main purpose of this study was to
clarify the generation mechanism of cloud
cavitation. The computed result showed a close
relationship between the behavior of the
separated shear layer and the cavitation
cloud. This is because the cavitation cloud
shedding frequency also depends on the
Reynolds number. Further theoretical
investigation is needed to predict the
cavitation cloud shedding frequency for
arbitrary f low conditions.
[Al ] Kubota,A. et al.: Unsteady Structure
Measurement of Cloud Cavitation on a Foi l
Section Using Conditional Sampling
Technique, ASME J. Fluid Eng., Vol.111,
No.2, 1989, pp.204-210
