Below are the first 10 and last 10 pages of uncorrected machineread 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, chapterrepresentative 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 103
Computation of a Free Surface Flow around
an Advancing Ship by the NavierStokes Equations
T. Hino
Ship Research Institute
Tokyo, Japan
Abstract
The finitedifference solution method for the
NavierStokes equations with nonlinear free surface con
dition is applied to the simulation of flow field with a
free surface around an advancing ship. The bodyfitted
coordinates system is used in order to cope with a ship
of an arbitrary hull form. The coordinates system does
not fit to the free surface geometry which must be de
termined as the part of solution in the time marching
procedure. The nonlinear free surface condition is im
plemented in the numerical scheme. The algebraic tur
bulence model is used together with the wall function on
the body boundary condition to simulate high Reynolds
number flow. The numerical results are compared with
the experimental data.
1. Introduction
A viscous flow field around a ship is strongly non
linear even when the free surface effects are neglected.
When a free surface deformation is taken into account,
the geometry of free surface boundary should be deter
mined as a part of solution by the nonlinear free sur
face condition and flow field becomes more complicated
in both physical and numerical aspects. A number of
efforts have been made to solve this nonlinear prob
lem. Among them, the finitedifference solution meth
ods for the NavierStokes equations with free surface
effects t1,2] seem to be most promising because of their
generality. However, even the NavierStokes solvers for
doublemodel flow around a ship [3,4] have not been
established well in the case of turbulent flow simula
tion because of various problems. There are much more
difficulties to be overcome in the development of free
surface flow solvers, such as the treatment of free sur
face boundary condition, grid generation strategy and
turbulence model.
In the present paper, the finitediRerence method
for the NavierStokes equations with nonlinear free sur
face condition developed in Reference t13 is extended to
high Reynolds number flow simulation around an ad
vancing ship. Nonlinear free surface condition is imple
mented in the scheme. The algebraic turbulence model
103
is used together with the wall function method for the
body boundary condition.
The outline of numerical scheme is described in
Section 2. The numerical results for Wigley's parabolic
hull and Todd's Series 60, Cb=0.6 are shown and com
pared with measured data in Section 3 and 4, respec
tively. The concluding remarks is given in Section 5.
2. Numerical Scheme
2.1 Governing Equations
The governing equations are the Reynolds aver
aged NavierStokes equations and the continuity equa
tion for mean velocity of unsteady threedimensional
incompressible fluid. They are written in dimensionless
form as follows;
at ~ mu + vuy + wuz
=—p+(l/Re+z~t)(u2+uy,~uzz) (la)
Vt + us + vvy + wvz
=—PI + (1/Re + lJt`)(~V~ + VYY ~ VZZ) (
OCR for page 103
In Eqs. (1)(4), subscripts, x, y, z and t mean the
partial differential.
The bodyfitted curvilinear coordinates system
(if, 71, ¢) is introduced to cope with the body boundary
of an arbitrary form, where ~ is the direction from fore
to aft, ~ the direction from a ship or a center plane to
the side outer boundary and ¢, the girth direction from
keel to deck. As same as the previous method [1], the
computational coordinates do not fit to the free surface
shape, so they are not timedependent. The coordinates
transformation is given as follows;
= ((a, y, a), 7~ = 71(z, y, Z), ¢ = ¢(z, y, Z), t = t (3)
The momentum equations (1) and the continuity
equation (2) are transformed through Eqs. (3) as
up + Us + Vu + We
=  (~` + rears + ¢~) + (l/Re + z~t)(V2u) (4a)
vet + Us + Vv,7'+ We
—((y¢~+0y¢~7+¢y¢~)+(l/Re+L/t)(V2v) (4b)
we+ Us +Vw77+ Ww:
=  ((z¢< + r~z¢,7 + Czar) + (1/Re + z~t)(V2w) (4c)
BUD + n,,~u,7 + ¢,ru~ + (yin + ~yV77 + ¢yVc
+(zw~ + ~zw,7+¢zw: = 0 (5) where
where (U. V, W) are the unscaled contravariant velocity
components and defined
U = flu + (yv + Saw (6a)
V = YOU + 71yV + 77zW (fib)
W = Mu + ¢yV + ¢zw (6c)
~ is pressure from which hydrostatic component is ex
tracted;
~ = p + z/Fn2 (7)
V2 is the transformed Laplacian operator and defined
as
V2q = (~2 + f2 + f2)q66 + (~2 + 02 + 02)q
+(¢+ + ¢y + ¢Z )q(<
+2( + By fly + (zNz)q6n
+2(r~¢+ + py¢y + Hz¢z)quc
+2(¢2~+¢y~y+¢zfz)q(E
+(~m + (§9 + (zz)q~ + (em + H99 + ~zz)
+(~= + ¢99 + ¢zz)q~
where q is arbitrary scalar quantity. A, by and so on
appeared in Eqs. (4)(8) are the metrics of the grid.
2.2 Basic Algorithm
The basic algorithm is same as that of the MAC
method [5~. The discretization is made in the non
staggered grid, that is, all variables are defined in the
intersections of grid lines. The present method is based
on the time marching procedure and is divided into two
stages.
On the first stage, velocity is updated by the mo
mentum equations (6). The forward difference is used
in time. The spatial differences are the thirdorder up
stream difference by Kawamura and Kuwahara t6] for
the convection terms, the secondorder central differ
ence for the pressure gradient terms and for the diffu
sion terms and the fourthorder central difference for
the grid metrics terms.
On the second stage, pressure on the next time step
is computed so that the velocity field on the next time
step may satisfy the continuity condition. By taking di
vergence of the momentum equations (6), the following
Poisson equation for pressure is derived.
v2¢ = —IKE—~Ku—DIE
—Lyle—Kyle—¢yLc t9
—(zM~—~zMn—¢zM~
—Dt
K = US + Van + Wuc—(l/Re + ut)(V2u)
L = US + Vv77 + Wvc—(l/Re + L/t)(V2v)
M = US + Vw,'+ Ww:—(l/Re + z~t)(V2w)
and
D = Cup + mud + ¢u; + love + ~yv,7 + ¢yvc
+(zw~ + Hzwn + ¢zw:
The righthandside of Eq.~9) is evaluated by the values
at the present time step. The spatial differences for
K, L and M are same as that for Eqs.~6). The time
differential appeared in the last term is expressed by
the forward difference. Then D, divergence of velocity,
on the next time step is set zero from the continuity
condition, while D on the present time step which is not
necessarily zero is evaluated by the secondorder central
difference to avoid accumulation of numerical errort53.
The lefthandside of Eq.~9) is evaluated by the second
order central difference and solved iteratively by the
Jacobi method.
The initial condition is still state, that is, velocity
and wave elevation are zero and pressure is hydrostatic
(8) in the whole domain of computation. The constant ac
celeration is made by adding the inertia force to the
momentum equation in xdirection, Eq.~4a), until the
inflow velocity becomes unity.
04
OCR for page 103
2.3 Free Surface Conditions
When the effects of viscosity and free surface ten
sion are neglected, the free surface conditions consist of
the following two conditions. One is the pressure con
dition that means that pressure on the free surface is
equal to atmospheric one. The other is the kinematic
condition that tells the fluid particles on the free surface
keep staying on it.
Because the grid points are not on the free surface
in the present grid system, it is not easy to satisfy the
free surface conditions on the exact location of the free
surface. The pressure condition is implemented in the
solution process of the Poisson equation for pressure. To
give the boundary condition at the intermediate point
between grid points, the 'irregular stars method' used in
the SUMMAC method [73 is extended to the curvilinear
coordinates system.
The kinematic condition is used to determine the
free surface shape in the time marching process. The
wave elevation is defined in the computational coordi
nates as
~ = h((, 0, t) (10)
The kinematic condition is written as 2.6 Turbulence Model
ht + US + Vhn—W = 0 on ¢ = h (11)
Eq.~11) is transformed into the finitedifference form in
the same manner as that for the momentum equations
(6). Velocity (U. V, W) on the free surface is extrapo
lated equally from the value at the lower grid points.
2.4 Body Surface Conditions
For the body surface condition, the wall function
approach is used to reduce the computation time. With
the noslip condition, the minimum grid spacing in the
direction normal to the body surface should be small
enough to resolve the viscous sublayer of the boundary
layer on the body. Because the explicit scheme in time
is used, the time increment is limited by the CFL condi
tion and should be also small in proportion to the grid
spacing. The total computational time to convergence
would be too large for the present computer power.
The wall function used here is the general logarith
mic law, that is,
q/ur = 1/,c In y+ ~ B (12)
where q is velocity magnitude, u, the friction velocity
and ye the normal distance from the wall normalized
by I//UT. The constants tic and B are set 0.4 and 5.5,
respectively.
Following Chen and Patel [8i, the two velocity
points ~ j = 2 and 3, where j = 1 is the wall ~ are
assumed to be located in the logarithmic region. From
q and the normal distance from the wall at j = 3, UT
is determined from Eq. (12). Then q at j = 2 can
be calculated from UT and the normal distance of the
point. Velocity at j = 2 is treated as the boundary
value in the velocity updating process. The direction of
velocity at j = 2 is assumed parallel to the wall. In the
accelerating period, the noslip condition is imposed for
velocity on the body. Pressure and wave elevation on
the wall (j = 1) is set equal to those at j = 2.
By the use of the wall function, the minimum grid
spacing can be more than ten times as large as that in
the case of the noslip condition. In the present compu
tations, typical value of y+ at j = 3 is taken about from
20 to 100. The computing time is reduced drastically
by this procedure.
2.5 Other Boundary Conditions
On the inflow boundary, velocity is uniform flow
in xdirection after the acceleration and pressure is hy
drostatic with zero wave elevation. On the outflow and
side boundaries, pressure, velocity and wave elevation
are extrapolated with zero gradient from the inside.
Turbulence model used is the twolayer algebraic
model by Baldwin and Lomax t93. It is widely used in
the aerodynamic computations and also in the incom
pressible flow computation around a ship by Kodama
t33. In the present study, flow is enforced to be turbulent
from the fore end of a ship. The free surface effect on
turbulence is not included in the model. There has not
been any turbulence model that can be applied to the
boundary layer and wake of a surfacepiercing body like
a ship. Further investigation in both computation and
experiment is required to establish a turbulence model
under the free surface effects.
3. Computation for Wigley's Hull
3.1 Computational Condition
The first computational results are for Wigley's
parabolic hull. The waterlines and the frame lines of
the hull geometry are defined by the parabolic lines.
The computations are made with two, coarse and fine,
grids. The grid generation scheme based on the geo
metrical method is used.
The coarse grid is shown in Fig.la. The number of
grid points is 51, 20 and 18 in the ((, A, A) directions,
respectively. The HO grid topology is adopted. The
grid points are clustered near the body and near the
still water surface. The computational domain in the
dimensionless coordinates ( where x = 0 is the midship,
y = 0 the center plane and z = 0 the still water level)
IS
105
OCR for page 103
Fig. la Coarse computational grid for Wigley's hull.
Fig. lb Fine computational grid for Wigley's hull.
1 < x < 1, 0 < y < 0.5, 0.5 < z < 0.0555
It should be noted that the domain includes the region
above the undisturbed free surface, that is, z ~ 0. The
grid points below the still water surface is 51 x 20 x 10.
The number of grid points inside the fluid varies as the
wave field develops. The minimum grid spacing in
direction is 0.001.
The fine grid shown in Fig.lb has 100 x 20 x 38 grid
points in ((, ~1, ¢) directions, respectively. The computa
tional domain is same as that for the coarse grid, except
that—0.5 < z < 0.036. The grid points under the free
surface is 100 x 20 x 30 and the minimum grid spacing
in redirection is 0.0008 in this case.
The Froude number is 0.25 and the Reynolds num
106
her is 106 in both computations. The acceleration is
made in the first 500 time steps. The dimensionless
time increment is 0.0005 for the coarse grid computa
tion. In the fine grid case, the dimensionless time incre
ment is 0.0005 from 1 to 2000th time step and 0.0003
from 2001th time step for stabilization of computation.
3.2 Accuracy Analysis

The time derivatives in the momentum equations
(6) and in the free surface kinematic condition (11)
are replaced by the forward onesided difference, that
means accuracy in time of the present method is the
firstorder. For the spatial differences, the pressure gra
dient terms and the diffusion terms have the second
order accuracy. The thirdorder difference is used for
the convection terms, in which the leading error is the
fourthderivative term and does not affect the diffusion
OCR for page 103
/: '\\ ~~ : ~ ~ ~ ~ ~ :~ ~
~ ... i,:.:::''' i,\~/'~ 1
Fig. 2a Time evolution of computed wave pattern around Wigley's hull
with the coarse grid. Contour interval is 0.02 x 2gh/UO. Dotted lines show
negative values. Top; 6000th step ~ t = 3.0 ), middle; 8000th step ~ t = 4.0
and bottom; 10000th step ~ t = 5.0 ).

'' """" .'~).J'. J
Fig. 2b Time evolution of computed wave pattern around Wigley's hull
with the fine grid. Contour interval is 0.02 x 2gh/UO. Dotted lines show neg
ative values. Top; 6000th step ~ t = 2.2 ), middle; 8000th step ~ t = 2.8
and bottom; 11000th step ~ t = 3.7 ).
107
OCR for page 103
1 an' W
Fig. 3a Pressure distribution on hull surface, center plane and free sur
face around Wigley's hull computed with the coarse grid. Contour interval is
0.02Cp. Dotted lines show negative values.
J ~2'  ~
' —l'' W<.~=''2''''\
Fig. 3b Pressure distribution on hull surface, center plane and free surface
around Wigley's hull computed with the fine grid. Contour interval is 0.02Cp.
Dotted lines show negative values.
terms of the momentum equations. The grid metrics
are evaluated by the fourthorder difference. Numerical
errors due to these finite differencing are the function
of the time increment and the grid spacing.
Other factor that determines accuracy is conver
gence. As for convergence of the Poisson solution for
pressure, the residual is typically 0~104) after 20 iter
ations. In the time integration process of the present
method, the quantity that converges most slowly is the
wave elevation. Therefore, convergence of the solution
is examined by steadiness of the wave patterns.
The comparison of the numerical results with the
fine and coarse grids provides information concerning
grid density effect. Fig. 2a shows the time sequence
of the wave patterns around Wigley's hull computed
with the coarse grid. The wave pattern has not reached
the steady state at 10000th time step when the dimen
108
sionless time is 5.0. The grid resolution seem to be not
sufficient to get convergence. In the case of the fine grid
shown in Fig. 2b, on the other hand, the wave pattern
has become almost steady at 1100~th time step (the
dimensionless time is 3.49~. The waves far from the
body in the coarse grid case are less steep than those in
the fine grid case. The numerical dissipation due to the
finite differencing error decreases the wave amplitude
when the grid spacing is large.
Fig.3a and 3b show the pressure distribution on
body surface, center plane and free surface in the coarse
and fine grid cases, respectively. Pressure value on the
free surface is identical to the wave elevation, because
hydrostatic component is extracted from static pres
sure. Strong wiggles of pressure can be seen on the body
surface in the coarse grid case. That may be one reason
why the solution has not converged. In the fine grid
OCR for page 103
Fig. 4 Perspective view of computed waves around Wigley's hull. Wave
height is three times magnified.
2*9*h/(U*U)
0.40
0.20
0~~
0.20
Computed
Measu red
J
J
~.5 x/L
Fig. 5 Comparison of computed and measured wave profiles on ship sur
face of Wigley's hull. En = 0.25, Re = 106 in computation and Re = 5 X 106
in measurement.
case, however, the wiggles are limited in the restricted
regions such as near the fore end and the aft end or near
the free surface and their magnitude is small.
3.3 Results
Hereafter, only the fine grid results are shown.
Fig.4 shows the perspective view of the free surface
configuration, where the wave amplitude is magnified
by three times. Waves far from the ship hull cannot be
seen clearly. The grid spacing is still too large, partic
ularly far from the ship hull, to diminish the numerical
dissipation effects. Fig.S is the comparison of the com
puted and the measured [10] wave profiles on the body
surface. Agreement is very well in wave amplitude and
in wave length. The orgin of wave making, apart from
the propagation of waves which is affected by the nu
merical dissipation, is simulated properly.
109
In Figs.3b and 5, stern wave generation which is not
clear in the previous computation for the low Reynolds
number laminar flow [1] is simulated in the present
results. Stern waves are related to pressure recovery
at the stern region and this becomes higher as the
Reynolds number increases.
Computed pressure distribution on body surface is
shown in Fig.6 together with the measured one [103.
Pressure patterns are in good accordance with each
other, except for the region of wiggles described above.
These wiggles seem to come from inappropriate treat
ment of the boundary condition. Pressure recovery at
the stern in the high Reynolds number flow described
above is simulated well.
Figs.7 show the wake (u velocity) contours and the
cross flow vectors (v and w velocity) at various stations.
The vertical velocity component appears beneath the
OCR for page 103
Computed
0.1 0 0.1 0
_., . .
o
AP
Fig. 6 Comparison of computed and measured pressure distribution on
ship surface of Wigley's hull. Contour interval is 0.02Cp. Dotted lines show
negative values. En = 0.25, Re = 106 in computation and Re = 3.4 x 106 in
measurement.
:,\~\
x=0 .5 (AP)
Fig. 7 Computed wake contours and cross flow vectors at various stations
of Wigley's hull. Contour interval is 0.1 x u. Top left; x =  0.5 (F.P.), top
right; x = 0 (midship), middle left; x = 0.3, middle right; x = 0.4, bottom
left; x = 0.5 (A.P.), bottom right x = 0.6.
110
OCR for page 103
x=0.4
x~0.5 (AP) x=0.6
Fig. 8 Computed eddy viscosity distribution at various stations of Wigley's
hull. Contour interval is 4.0 x 106ut. Top left; x = 0.3, Top right; ~ = 0.4,
bottom left; x = 0.5 (A.P.), bottom right x = 0.6.
free surface. In particular, large upward velocity is seen
at F.P. station, which corresponds with the generation
of bow waves. The wake contours at A.P. station seems
to be too thick though the corresponding measured data
is not presented. The improvement of the turbulence
model and/or the wall function approach is required.
In Figs.8, the eddy viscosity distributions at var
ious stations are presented. The discontinuity in the
distribution comes from the fact that the eddy viscos
ity is determined line by line. The strong eddy viscosity
regions spread widely near the free surface at x = 0.4
and 0.5 stations. This may be because the turbulence
model is affected by the free surface motion and it does
not have physical meaning.
4. Computation for Series 60, Cb=0.6
4.1 Computational Condition
Fig. 9 Computational Grid for Series 60, Cb=0.60
The second result is for the practical ship hull form,
Series 60, Cb=0.6. The computational domain is
1 < x < 1, 0 < y < 0.5, 0.5 < z < 0.0384
and the number of grid points is 100 x 20 x 38 which is
same as the fine grid case for Wigley's hull. The com
putational grid is shown in Fig. 9. The minimum grid
spacing in direction is 0.0008. The Froude number is
0.22 and the Reynolds number is 106 in this case. The
acceleration is made in the first 500 time steps. The di
mensionless time increment is made smaller gradually
as the time step increases, otherwise the computation
cannot be stable. From 1st to 1200th time step, the
time increment is set 0.0005, then 0.0003 from 1201th
to 2000th time step, 0.00025 from 2001th to 3000th
time step and 0.0002 after that.
111
4.2 Results
OCR for page 103
t=7 .09
t=7.49
I
. . .
. .
/
~~ :~ G
,.~ :~.~1 .''~,
t=789 ~ )J{~
Fig. 10 Time evolution of computed wave pattern around Series 60 Con
tour interval is 0.02 x 2gh/UO. Dotted lines show negative values. Top; 31000
th step ~ t = 7.09 ), middle; 3300~th step ~ t = 7.49 ~ and bottom; 35000th
step ~ t = 7.89 ).
The time evolution of wave patterns are shown in
Figs.10. The wave field has not reached the steady
state at 35000th time step (the dimensionless time is
7.89~. It takes very large time to get converged solution
for free surface problems by the time marching proce
dure, though the use of the wall function approach con
tributes to time saving to some extent. Improvement
of the numerical scheme is required to get faster con
vergence. One reason why it takes longer to get con
vergence in the case of Series 60 than in the case of
Wigley's hull may be the difference of flow complexity.
The hull geometry of a practical ship, such as Series
60, is more complicated than that of Wigley's hull and
flow around the complicated geometry does not become
steady in a short time.
Hereafter, the numerical results at 3500~th time
step are shown, because the flow field near the ship can
be considered as almost steady.
Fig. 11 shows the pressure distribution on body sur
face, center plane and free surface. Wiggles of pressure
can be seen on the body surface near the fore and aft
ends and near the free surface. Wiggles beneath the
free surface is partly due to the discontinuity of grid
spacings near the free surface which comes from the
constraint of the grid generation scheme.
Fig.12 shows the perspective view of the free sur
face configuration, where wave amplitude is magnified
by three times. Fig.13 is the comparison of the com
puted and measured wave profiles t11] on the body sur
face. The slight unphysical oscillation of wave elevation
is found in the aft part. The hull geometry of Series
60 is more complicated, particularly in the stern part,
than that of Wigley's hull and the grid lines around the
stern are more distorted. Numerical error due to the
distorted grid causes the oscillation of waves. Except
that, the computed result agrees well with measured
data in wave amplitude and in wave length.
In Fig.14, the computed pressure distribution on
body surface is compared with the measured one t124.
Pressure patterns are in good accordance with each
other except for the slight wiggles of computed results.
Free surface effects on pressure distribution beneath the
free surface, that is, the low pressure zone below the
wave trough and the high pressure zone below the wave
crest, can be seen in both computation and measure
ments.
Figs.15 show the wake (u velocity) contours and
the cross flow vectors (v and w velocity) at various sta
tions together with the measured data t123. At the F.P.
station, large upward velocity appears as well as in the
case of Wigley's hull. The wake contours at the midship
station becomes thin around the bilge circle in the mea
surement and this is well simulated in the computation.
At the A.P. station, longitudinal vortices can bee seen
OCR for page 103
Fig. 11 Pressure distribution on hull surface, center plane and free sur
face around Series 60. Contour interval is 0.02Cp. Dotted lines show negative
values.
Fig. 12 Perspective view of computed waves around Series 60. Wave
height is three times magnified.
f: 0.20
nit
~  ~ .vv
—0.5  0.4 ~—O. 1 ~
0.20
2*9*h/(U*U)
0.4OI
Computed
M ensured
04 0.5 x/L
Fig. 13 Comparison of computed and measured wave profiles on ship sur
face of Series 60. En = 0.22, Re = 106 in computation and Re = 1.39 x 107
in measurement.
113
OCR for page 103
Computed
~ ;g' it' ~
. .
0.1 ~ ~
Fig. 14 Comparison of computed and measured pressure distribution on
ship surface of Series 60. Contour interval is shown in figure. Dotted lines show
negative values. Fn = 0.22, Re = 106 in computation and Re = 7.7 X 106 in
measurement.
Computed
rid rid
/
Measured Computed
,Ii,,~m!, j 1 ~ 1, '\
x=0.25
Measured / Computed
0.9 / 0.9
I (~'\'''"N'\~.~:,_=
\
a\\\'\\\\
1111111\ \ \
x=0.5 (AP)
Measured Computed
x=0 (midship) 
Measured Computed /
0.9 /0.9
(go\ 1 ~
·,.lU
Measured
Computed
O .9
1
/: / —f '
Fig. 15 Computed and measured wake contours and cross flow vectors at
various stations of Series 60. Contour interval is 0.1 x u. En = 0.22, Re = 106
in computation and Re = 7.7 X 106 in measurement. Top left; x =  0.5 (F.P.),
top right; x = 0 (midship), middle left; x = 0.25, middle right; x = 0.4, bot
tom left; x = 0.5 (A.P.), bottom right x = 0.6.
114
OCR for page 103
/
/
1 west 7
1_~
Fig. 16 Computed eddy viscosity distribution at various stations of Series
60. Contour interval is 4.0 x 106L,t. Top left; ~ = 0.25, top right; x = 0.4,
bottom left; x = 0.5 (A.P.), bottom right x = 0.6.
in both measurement and computation at the almost
same position. However, the computed wake contour is
thicker than measured one. At the section of x=0.6, the
positions of the longitudinal vortices are different from
each other.
Figs.16 show the eddy viscosity contours at various
stations. Nonphysical eddy viscosity beneath the free
surface appeared in the case of Wigley's hull cannot be
seen in this case.
5. Concluding Remarks
The finitediderence solution method for the
Reynolds averaged NavierStokes equations is applied
to the simulation of high Reynolds number flow with a
free surface around an advancing ship. The numerical
results for Wigley's hull and Series 60, Cb=0.6 show
good agreement with the experimental data in wave
profiles and pressure distributions on the hull surface.
The further efforts are required in the turbulence mod
eling under the free surface effect to obtain quantitative
agreement for the wake distribution. The improvement
of the numerical scheme is also needed to get converged
solution faster.
Acknowledgment
The author would like to express his sincere grati
tude to the members of the CFD group at Ship Research
Institute for their valuable discussions. In particular, he
is grateful to Dr. Yoshiaki Kodama, Ship Research In
stitute, for his suggestions and for providing the author
with the graphic output code.
The computations were carried out by ACOS 910
with Scientific Attached Processor at the Computer
Center, Ship Research Institute. The CPU time was
115
about two hours per 1000 time steps in the fine grid
case.
References
[1] T. Hino, "Numerical Simulation of a Viscous Flow
with a Free Surface around a Ship Model", Jour
nal of The Society of Naval Architects of Japan,
Vol.161, ~ 1987 ).
t2] T. Sato, H. Miyata, N. Baba and H. Kajitani,
"FiniteDifference Simulation Method for Waves
and Viscous Flows about a Ship", Journal of
The Society of Naval Architects of Japan, Vol.160,
(1986~.
t3] Y. Kodama, "Computation of High Reynolds Num
ber Flows Past a Ship Hull Using the IAF Scheme",
Journal of The Society of Naval Architects of
Japan, Vol. 161, ( 1987~.
[4] A. Masuko, Y. Shirose, Y. Ando and M. Kawai,
"Numerical Simulation of Viscous Flow around a
Series of Mathematical Ship Models", Journal of
The Society of Naval Architects of Japan, Vol.162,
(1987~.
F.H. Harlow and J.E. Welch, "Numerical Calcu
lation of TimeDependent Viscous Flow of Fluid
with Free Surface", The Physics of Fluids, Vol.8,
(1965~.
t6] T. Kawamura and K. Kuwahara, "Computation
of High Reynolds Number Flow around a Circular
Cylinder with Surface Roughness", AIAA paper,
840340, (1984~.
OCR for page 103
[7] R.K.C. Chan and R.L. Street, "A Computer Study
of FiniteAmplitude Water Waves", Journal of
Computational Physics, Vol.6, ( 1970~.
[8] H.C.Chen and V.C. Patel, "Calculation of Trailing
Edge, Stern and Walce Flows by a TimeMarching
Solution of the PartiallyParabolic Equations",
lIHR Report, No.285, (1985~.
[9] B. Baldwin and H. Lomax, "ThinLayer Approxi
mation and Algebraic Model for Separated Turbu
lent Flows", AIAA paper, 78257, (1978~.
[10] "Cooperative Experiments on Wigley Parabolic
Models in Japan", 1 7th ITTC Resistance Commit
tee Report,(1983~.
[11] H. Takeshi, T.Hino, M. Hinatsu, Y. Tsukada and
J. Fujisawa, "ITTC Cooperative Experiments on a
Series 60 Model at Ship Research Institute", Pa
peTs of Ship Research Institute, Supplement No.9,
(1987~.
[12] F. Mewis and H. J. Heinke, " Untersuchungen
der Umstromung eines Modells der Serie 60 mit
Cb=0,60", Sciffbaujorschung, Vol.23, No.3, (1984~.
6
OCR for page 103
DISCUSSION
by S. Ju
1. In the wall function approach, y+
values of third grid points from the body
between 20 and 100 look too small, since y+
values of 2nd grid points from the body will
be less than 10 and in the region of laminar
sublayer.
2. How the transition from the laminar to
turbulent flow is treated?
Author's Reply
1. Yes. The grid spacing near the wall in
this calculation is a little too small. This
should be improved in the future calculation.
2.Though the original BaldwinLomax
turbulence model can cope with the transition,
flow is assumed to be turbulent from the fore
end of a ship in this computation.
DISCUSSION
by R.C. Ertekin
I would like to ask two questions.
1) Because you are solving a symmetric
problem, and therefore, using the symmetry
condition on the center plane, the normal
vector at the point where the center plane
meets the bow or stern is double valued or it
has discontinuous first derivatives. As a
result, the Jacobian of the transformation
matrix vanishes there. How did you cope with
this problem? Where were the stagnation points
on the bow and stern? Have you had any
difficulty with the noslip boundary condition
at these points?
2) I am surprised not to see any
comparisons with the care in which the linear
free surface boundary condition is supposed.
Why can you solve the nonlinear problem but
not the linear one? I don't also understand
why you can not calculate the resistance
experienced by the hull. After all, shouldn't
the objective of such calculations be the
determination of power requirements? By the
way when I say the linear problem I mean the
linear (viscous) governing equations and
boundary conditions!
Author's Reply
1. Along the line of mapping singularity,
flow quantities, are not computed but set the
average of values of the adjacent grid points.
The present scheme does not use the noslip
boundary condition but use the wallfunction
approach. Anyway, I don't have any difficulty
on the bow or stern.
2. The implementation of the linearized
free surface condition in the present scheme
is not impossible, but apparently the solution
becomes less accurate. We cannot compute the
wavemaking resistance separately from the
numerical results but, of course, we can
compute total resistance. Please see the reply
to Dr. Musker.
DISCUSSION

by A.J. Musker
I should be interested to know whether you
have made any efforts to calculate total
resistance. The pressure contours appear to be
well predicted and you presumably also have
calculated data on wall shear stress. Perhaps
you could comment.
Author's Reply
Pressure resistance and friction
resistance are computed by the integration of
numerical results as follows.
. Ship Fn Re rp rF rT
_
Wigley O .25 101 ~ 1. 44xlOd 3.73x103 5.17xlO:
Series O .22 106 0.41x103 4.13x103 4.54x103
where rp = (Pressure Resist. ) / 2 ,oU2S
rF = (Friction Resist. ) / 12 ,oU2S
rT = rp + rF
OCR for page 103