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 699
An Interactive Approach for Calculating Ship Boundary
Layers and Wakes for Nonzero Froude Number
Y. Tahara, F. Stern
(Iowa Institute of Hydraulic Research, The University of Iowa, USA)
B. Rosen (South Bay Simulations Inc., USA)
ABSTRACT
An interactive approach is set form for calculating
ship boundary layers and wakes for nonzero Froude
number. The Reynoldsaveraged NavierStokes equa
tions are solved using a small domain with edge condi
tions matched with those from a sourcedoubletDawson
method solved using the displacement body. An
overview is given of both the viscous and inviscidflow
methods, including their treatments of the freesurface
boundary conditions and interaction procedures. Results
are presented for the Wigley hull, including comparisons
for zero and nonzero Froude number and with available
experimental data and the inviscidflow results, which
validate the overall approach and enable an evaluation of
the waveboundary layer and wake interaction.
NOMENCLATURE
A¢,B¢,etc.
Aij,Bij,aj,bjk
i
bj
CD,CP,CU,Cnb finiteanalytic coefficients (nb = NE,
NW, SE, etc.) 2
friction coefficient (= 2~W/pUo)
pressure coefficient
residuaryresistance coefficient (=
2R/pSU2o)
Froude number (= U
Historically, inviscidflow methods have been
used to calculate wavemaking and viscousflow methods
the boundary layer and wake, in both cases, without
accounting for the interaction. Recent work on wave
making has focused on the solution of the socalled
NeumannKelvin problem using both Rankine and
Havelocksource approaches. Methods implementing
these approaches were recently competitively evaluated
and ranked by comparing their results with towingtank
experimental data [11. In general, the methods under
predicted the amplitude of the divergent bow waves,
were lacking in high wavenumber detail in the vicinity
of the bowwave cusp line, and overpredicted the ampli
tudes of the waves close to the stern. These difficulties
were primarily attributed to nonlinear and viscous
effects. The methods using the Havelocksource
approach generally outperformed those using the
Rankinesource approach, except with regard to the near
field results (i.e., within one beam length of the model)
for which one of the the lager methods [2] was found to
be far superior.
Considerable effort has been put forth in the
development of viscousflow methods for ship boundary
layers and wakes. Initially, threedimensional integral
and differential boundarylayer equation methods were
developed; however, these were found to be inapplicable
near the stern and in the wake. More recently, efforts
have been directed towards the development of Navier
Stokes (NS) and Reynoldsaveraged NavierStokes
(RANS) equation methods; hereafter both of these will
simply be referred to as RANS equation methods. At
present, the status of these methods is such that practical
ship geometries can be considered, including complexi
ties such as appendages and propellers. Comparisons
with experimental data indicate that many features of the
flow are adequately simulated; however, turbulence
modeling and grid generation appear to be pacesetting
Issues with regard to future developments (see, e.g., the
review by Patel [3] and the Proceedings of the 5th
International Conference on Numerical Ship
Hydrodynamics [43~.
Relatively little work has been done on the inter
action between wavemaking and boundary layer and
wake. Most studies have focused separately on either
the effects of viscosity on wavemaking or the effects of
wavemaking (i.e., waves) on the boundary layer and
wake. Professor Landweber and his students have both
demonstrated experimentally the dependence of wave
resistance on viscosity and shown computationally that
by including the effects of viscosity in inviscidflow cal
culations of wave resistance beKer agreement with exper
imental data is obtained (most recently, [5~. Such
effects have been confirmed by others, including other
more detailed aspects of the flow field such as surface
pressure distributions and wave profiles and patterns [61.
Most studies concerning the effects of waves on
boundary layer and wake have been of an approximate
nature utilizing integral methods and assuming small
crossflow conditions (see Stern [7] for a more complete
review, including references). In [7,8], experiment and
theory are combined to study the fundamental aspects of
the problem utilizing a unique, simple model and
computational geometry, which enabled the isolation and
identification of certain important features of the wave
induced effects. In particular, the variations of the wave
induced piezometricpressure gradients are shown to
cause acceleration and deceleration phases of the
streamwise velocity component and alternating direction
of the crossflow, which results in large oscillations of
the displacement thickness and wallshear stress as com
pared to the nowave condition. For the relatively simple
geometry studied, f~rstorder boundarylayer calculations
with a symmetrycondition approximation for the free
surface boundary conditions were shown to be satisfac
tory; however, extensions of the computational approach
for practical geometries were not successful [91.
Miyata et al. [101 and Hino [11] have pursued a
comprehensive approach to the present problem in which
the NS equations (subgrid scale and Reynolds aver
aged, respectively) are solved using a large domain with
approximate freesurface boundary conditions. In both
cases, the basic algorithms closely follow those of MAC
[12] and SUMMAC [131. However, [10] uses a time
dependent freesurface conforming grid, whereas [111
uses a fixed grid which does not conform to the free sur
face. The results from both approaches are promising,
but, thus far, have had difficulties in accurately resolving
the boundarylayer and wake regions and, in the case of
[10], have been limited to low Re.
The present interactive approach is also compre
hensive. Two of the leading inviscid [2] and viscous
flow [14] methods are modified and extended for inter
active calculations for ship boundary layers and wakes
for nonzero Fr. The interaction procedures are based on
extensions of those developed by one of the authors for
zero Fr [151. The work of [7,8,15] is precursory to the
present study. Also, it should be mentioned that the
present study is part of a large project concerning free
surface effects on boundary layers and wakes. Some of
the related studies under this project will be referenced
later.
In the following, an overview is given of both
the viscous and inviscidflow methods, with particular
emphasis on their treatments of the freesurface
boundary conditions and the interaction procedures.
Results are presented for the Wigley hull, including
comparisons for zero and nonzero Fr and with available
experimental data and inviscidflow results, which
validate the overall approach and enable an evaluation of
the waveboundary layer and wake interaction. In the
presentation of the computational methods and results
and discussions to follow, variables are either defined in
the text or in the NOMENCLATURE and are
nondimensionalized using the ship length L, freestream
velocity UO, and fluid density p.
COMPUTATIONAL METHODS
Consider the flow past a shiplike body, moving
steadily at velocity UO, and intersecting the free surface
of an incompressible viscous fluid. As depicted in figure
1, the flow field can be divided into four regions in each
of which different or no approximations can be made to
the governing RANS equations: region 1 is the inviscid
flow; region 2 is the bow flow; region 3 is the thin
boundary layer; and region 4 is the thick boundary layer
and wake. The resulting equations for regions 1 and 3
and their interaction (or lack of one) are well known.
Relatively little is known about region 2. Recent exper
iments concerning scale effects on nearfield wave
patterns have indicated a Re dependency for the bow
wave both in amplitude and divergence angle [161;
however, this aspect of the problem is deferred for later
study. Herein, we are primarily concerned with the flow
700
OCR for page 699
in region 4 and its interaction with that in region 1. As
discussed earlier, the description of the flow in region 4
requires the solution of the complete RANS equations
(or, in the absence of flow reversal, the socalled
partiallyparabolic RANS equations, however, this
simplification will not be considered here).
There are two possible approaches to the solution
of the RANS equations: a global approach, in which one
set of governing equations appropriate for both the
inviscid and viscousflow regions are solved using a
large solution domain so as to capture the viscousinvis
cid interaction; and an interactive approach, in which dif
ferent sets of governing equations are used for each
region and the complete solution obtained through the
use of an interaction law, i.e., patching or matching
conditions. Both approaches are depicted in figure 1.
The former approach is somewhat more rigorous
because it does not rely on the patching conditions that
usually involve approximations. Nonetheless, for a
variety of reasons, both types of approaches are of inter
est. In [151, both approaches were evaluated for zero Fr
by comparing interactive and largedomain solutions for
axisymmetric and simple threedimensional bodies using
the same numerical techniques and algorithms and turbu
lence model. It is shown that both approaches yield sat
isfactory results, although the interaction solutions
appear to be computationally more efficient. As men
tioned earlier, the present study utilizes the interactive
approach. This takes advantage of the latest develop
ments in both the inviscid and viscousflow technolo
gies; however, a largedomain solution for the present
problem is also of interest and a comparative evaluation
as was done previously for zero Fr is planned for study
under the present project for nonzero Fr.
ViscousInviscid Interaction
Referring to figure 1, there are two primary dif
ferences between the interactive and largedomain
approaches with regard to the solution of the RANS
equations: (1) the size of the solution domain, i.e., the
placement of the outer boundary SO; and (2) the bound
ary (i.e., edge) conditions specified thereon. For the
largedomain solution, uniformflow and waveradiation
conditions are appropriate, whereas the interaction solu
tion requires the specification of the match boundary
(i.e., SO) as well as an interaction law, and also a method
for calculating the inviscid flow.
In the present study, solutions were obtained
with the match boundary at about 2b, where ~ is the
boundarylayer and wake thickness. The interaction law
is based on the concept of displacement thickness 8*. A
threedimensional 5* for a thick boundary layer and
wake can be defined unambiguously by the two require
ments that it be a stream surface of the inviscid flow
continued from outside the boundary layer and wake and
that the inviscidflow discharge between this surface and
any stream surface exterior to the boundary layer and
wake be equal to the actual discharge between the body
and wake centerplane and the latter stream surface. A
method for implementing this definition for practical
geometries is presently under development [171; how
ever, in lieu of this, an approximate definition is used in
which twodimensional definitions for a*, i.e.
~ = `: (MU) dr (1)
for the keelplane and waterplane at each station are con
nected by a second order polynomial.
In summary, the inviscidflow solution is
obtained for the displacement body 5*. This solution
then provides the boundary conditions for the viscous
flow solution, i.e.
U(So) = Up(So) = Ue W(SO) = Wp(SO) = We
p(So) = pp(So) = pe (2)
Because S* and Vp (SO) are not known a priori, an initial
guess must be provided and the complete solution
obtained by iteratively updating the viscous and
inviscidflow solutions until the patching conditions (1)
and (2) are satisfied.
Viscous Flow
The viscous flow is calculated using the large
domain method of Patel et al. [14], modified and
extended for interactive calculations and to include free
surface boundary conditions. The details of the basic
method are provided by [141. Herein, an overview is
given as an aid in understanding the present
modifications and extensions.
Eguations and Coordinate System
The RANS equations are written in the physical
domain using cylindrical coordinates (x,r,8) as follows:
au + ~ a ~rV' + ~ aW = 0 (3)
Dt = ~ ~ (¢ + Hi) ~ ear Ad) ~ r am (uw)
UV+ReV2U (4)
DDV W2 =  ad `~y ~9 +~  r are (vw)
 r (A  ww) + Re (V2Vr~ at rim) (5)
DW VW ~ a 
D +=  ~ (uw) ~ ar (vw)
r at (6 + ww)  r (vw) + Re (V2W
+ ~ DO ~ r=) (6)
with
701
OCR for page 699
DD` =~+Uaa +vaa +_ a
and v2=~+~+ ~ a + ~ a2
Closure of the RANS equations is attained
through the use of the standard k£ turbulence model
without modifications for freesurface effects. The lim
ited experimental data available for surfacepiercing
bodies [18] indicate that, near a free surface, the normal
component of turbulence is damped and the longitudinal
and transverse components are increased. This effect
has also been observed in openchannel flow [19] and in
recent measurements for freesurface effects on the wake
of a submerged flat plate [20] and a plane jet [211. Such
a turbulence structure cannot, in fact, be simulated with
an isotropic eddy viscosity turbulence model like the
present one; however, this aspect of the problem is
deferred for later study.
In the standard k£ turbulence model, each
Reynolds stress is related to the corresponding mean rate
of strain by the isotropic eddy viscosity vt as follows:
au av 1 au aw
uv=vt(~r + ~X) Uw=vt(r as + ax)
1av aw w
Vw=vt(r all + fir ~ r )
 UU =Vt (24)  3k  VV =Vt (2 ear) 3 k
 ww = vt (r aa~ 2 Vr ) 3 k (7)
Vt is defined in terms of the turbulent kinetic energy k
and its rate of dissipation £ by
k2
Vt = Cot
£
where Cp is a model constant and k and £ are governed
by the modeled transport equations
Dk_ a (1 ak
fit ax Rkax)
+ r Or(Rkraark) +~aa~ (Rika~ )+ G  £ (9)
Do a 1 as 1 a 1 as
Dt = aX (R aX' + r ar (R£ r ar)
+~a~ (Ready + C£1 k G  C£2 k (10)
G is the turbulence generation term
G=Vtt2~(aaXU)+(aavr)+(aa~w3+vr)
+ (aU+av'2+ ~ aU+aaxw' +(rDa~v+aDw r~ }(11)
The effective Reynolds number Rib is defined as
Rip Re ~, (12)
in which ~ = k for the kequation (9) and ~ = £ for the £
equation (10~. The model constants are: Cal = .09, C£1
= 1.44, C£2 = 1.92, t7U = <'V = oW = ok = 1, C7£ =
1.3.
The governing equations (3) through (12) are
transformed into nonorthogonal curvilinear coordinates
such that the computational domain forms a simple rect
angular parallelepiped with equal grid spacing. The
transformation is a partial one since it involves the coor
dinates only and not the velocity components (U,V,W).
The transformation is accomplished through the use of
the expression for the divergence and "chainrule" defi
nitions of the gradient and Laplacian operators, which
relate the orthogonal curvilinear coordinates x1 = (x,r.,8)
to the nonorthogonal curvilinear coordinates 41=
(t,q,(~. In this manner, the governing equations (3)
through (12) can be rewritten in the following form of
the continuity and convec~avetransport equations
am (bin + b2V + b3W) + as (b2U + b2V + bow
+ Ha: (b3U + b3V + b3W) = 0 (13)
all a2¢ + g22 a2¢ + g33 02~¢ = 2At a:
(8, + 2B~ ad) + 2C~ am + If am + so (14)
Discretization and VelocityPressure Coupling
The convectivetransport equations (14) are
reduced to algebraic form through the use of a revised
and simplified version of the finiteanalytic method. In
this method, equations (14) are linearized in each local
rectangular numerical element, AK, = Aq = A: = 1, by
evaluating the coefficients and source functions at the
interior node P and transformed again into a normalized
form by a simple coordinate stretching. An analytic
solution is derived by decomposing the normalized
equation into one and twodimensional partialdifferen
tial equations. The solution to the former is readily
obtained. The solution to the latter is obtained by the
method of separation of variables with specified bound
ary functions. As a result, a twelvepoint finiteanalytic
formula for unsteady, threedimensional, elliptic equa
tions is obtained in the form
702
OCR for page 699
1 8
(P= R {) Cnb~nb
l+CpECU+CD+  ~ 1
+ Cp(CU¢U + CD¢D +~  S)1
(15)
It is seen that up depends on all eight neighboring nodal
values in the crossplane as well as the values at the
upstream and downstream nodes MU and ~D, and the
values at the previous time step ~ I. For large values of
the cell Re, equation (15) reduces to the partially
parabolic formulation which was used previously in
other applications. Since equations (15) are implicit,
both in space and time, at the current crossplane of
calculation, their assembly for all elements results in a set
of simultaneous algebraic equations. If the pressure field
is known, these equations can be solved by the method
of lines. However, since the pressure field is unknown,
it must be determined such that the continuity equation is
also satisfied.
The coupling of the velocity and pressure fields
is accomplished through the use of a twostep iterative
procedure involving the continuity equation based on the
SIMPLER algorithm. In the first step, the solution to the
momentum equations for a guessed pressure field is cor
rected at each crossplane such that continuity is satisfied.
However, in general, the corrected velocities are no
longer a consistent solution to the momentum equations
for the guessed p. Thus, the pressure field must also be
corrected. In the second step, the pressure field is
updated again through the use of the continuity equation.
This is done after a complete solution to the velocity field
has been obtained for all crossplanes. Repeated global
iterations are thus required in order to obtain a converged
solution. The procedure is facilitated through the use of
a staggered grid. Both the pressurecorrection and pres
sure equations are derived in a similar manner by substi
tuting equation (15) for (U,V,W) into the the discretized
form of the continuity equation (13) and representing the
pressuregradient terms by finite differences.
Solution Domain and Boundary Conditions
The solution domain is shown in figure 1. In
terms of the notation of figure 1, the boundary condi
tions on each of the boundaries are as follows. On the
inlet plane Si, the initial conditions for ~ are specified
from simplellatplate and the inviscidflow solutions.
On the body surface Sb, a twopoint wallfunction
approach is used. On the symmetry plane Sk, the
conditions imposed are 8(U,V,p,k,£~/a~ = W = 0. On
the exit plane Se, axial diffusion is negligible so that the
exit conditions used are a2~/aX2 = 0, and a zerogradient
condition is used for 19. On the outer boundary SO, the
edge conditions are specified according to (2), i.e.,
(U,W,6) = (Ue,We, be) and a (k,£~/ar = 0, where
(Ue,We,pe) are obtained from the inviscidflow solution
evaluated at the match boundary SO.
On the freesurface S,~(or simply 11), there are
two boundary conditions, i.e.
and
V n = 0
*
tijnj = tiinj
(16)
(17)
where n is the unit normal vector to the free surface and
tij and All are the fluid and externalstress tensors
respectively, the latter, for convenience, including
surface tension. The kinematic boundary condition
expresses the requirement that ~ is a stream surface and
the dynamic boundary condition that the normal and
tangential stresses are continuous across it. Note that ~
itself is unknown and must be determined as part of the
solution. In addition, boundary conditions are required
for the turbulence parameters, k and £; however, at
present, these are not well established.
In the present study, the following approxima
tions were made in employing (16) and (17~: (a) the
external stress and surface tension were neglected; (b)
the normal viscous stress and both the normal and tan
gential Reynolds stresses were neglected; (c) the curva
ture of the free surface was assumed small and the tan
gential gradients of the normal velocity components were
neglected in the tangential stresses; and (d) the wave ele
vation was assumed small such that both (16) and (17)
were represented by firstorder Taylor series expansions
about the mean waveelevation surface (i.e., the water
plane Sw). Subject to these approximations, (16) and
(17) reduce to the following:
(UX11 X + Vylly  WZ) ~ = 0 (18)
Sw
19(Sw)=ll/Fr2~] afi I (19)
az ~
low
Ia~v,k,£'
=0
US ~
NEW
(20)
where Cartesian coordinates (x,y9z) have been used in
(18) and (19~. Conditions (18) through (20) were
implemented numerically as follows. The kinematic
condition (18) was used to solve for the unknown free
surface elevation ~ by expressing the derivatives in
finitedifference form and 11 in terms of its difference
from an assumed (or previous) value. A backward dif
ference was used for the xderivative, a central difference
for the yderivative, and the inviscidflow lip was used
as an initial guess. The dynamic conditions, (19) and
(20), were used in conjunction with the solution for ~ in
solving the pressure and momentum and turbulence
model equations, respectively. Backward differences
were used for the z and ~derivatives.
Inviscid Flow
The inviscid flow is calculated using the method
of Rosen [2], i.e., the SPLASH computer code. The
method is an extended version of the basic panel method
of Maskew [22,231 originally developed for the predic
tion of subsonic aerodynamic flows about arbitrary con
figurations, modified to include the presence of a free
surface and gravity waves both for submerged and
surfacepiercing bodies. As is the case with the basic
703
OCR for page 699
method, lifting surfaces and their associated wake
treatments as well as wall boundaries are included;
however, the present overview and calculations are for
nonlifting unbounded flow (see [24] for SPLASH
results for lifting flow). The details of the basic method
are provided by [22,231. Herein, an overview is given
as an aid in understanding the extensions for the
inclusion of the free surface and gravity waves and the
present interaction calculations.
The flow is assumed irrotational such that the
governing differential equation is the Laplace equation
V20 = 0 (21)
where 0 is the external perturbation velocity potential,
i.e.
Vp = UOx + Vo (22)
A solution for ~ may be obtained by defining also an
internal perturbation potential ¢~ and applying Green's
theorem to both the inner and outer regions and
combining the resulting expressions to obtain
¢=  J FOUR~) + R~} dS (23)
sb
where RpQ is the distance from the surface point Q to the
field point P and ~ = ¢~  O and c, = D(¢  o~/anQ are the
dipole and source strengths, respectively. In [22], the
nature of solutions to (23) is investigated for two dif
ferent specifications for hi, i.e., ¢~ = 0 and UOx. In both
cases, (23) is solved for the surface potential (i.e.,
¢(Sb)) by representing the body by flat quadrilateral
panels over which ~ and ~ are assumed constant and uti
lizing the farfield a ~ O and body D¢/an =  UOnx
boundary conditions. The zero internal perturbation
potential formulation Cot = 0) is shown to produce
"results of comparable accuracy to those from higher
order methods for the same density of control points."
In this case, the velocity normal to the external surface
vn is
Vn = UOnx + 80/0n = UOnx + ~ (24)
and, the velocity tangent to the external surface V' is
V~=Uotx+a~l~t=uotx~/8t (25)
where tx is the xcomponent of a tangent vector and t is
arclength in a tangential direction. For solid surfaces, Vn
is usually zero, but it may be a specified nonzero value to
simulate body motion, boundarylayer growth, inflow
and outflow, controlsurface deflection, etc. Hence, in
the basic method, (24) is used to evaluate the source
strengths directly. The corresponding doublet strengths
are then given by solution of the discretized form of
(23~. Values of V' are subsequently computed using
(25) with a central difference for the tderivative. It
should be recognized that the socalled zero internal
perturbation formulation is, in fact, equivalent to
methods based on Green's third formula applied directly
to the external perturbation potential (e.g., [251~.
In the SPLASH code, the internal zeroperturba
tion boundary condition is satisfied not only inside the
submerged portion of the configuration, but also on the
"other side" of a finite portion of the free surface. Both
are represented by sourcedoublet singularity panels and
flow leakage from one side of the freesurface to the
other, at the freesurface outer boundary is assumed to
be negligible. This assumption is valid if the outer
boundary of the free surface is sufficiently far from the
configuration, and if wave disturbances are eliminated
before reaching the freesurface outer boundary. In this
case, the discretized form of (23) is
(i = Hi, Aij ale + Hi, Bij Hi = 0 (26)
Sb + Sw Sb + Sw
The freesurface shape is determined by
representing the undisturbed free surface by panels,
whereupon freesurface boundary conditions linearized
with respect to zero Fr are imposed [261. The zero Fr
velocities, UO, VO, and WO, are obtained by first
considering all freesurface panels as solid and fixed (in
contrast to a traditional approach which employs the
double panel or image model). The nonzero Fr velocities
are then expressed as small increments to those for zero
Fr. The velocities tangent and normal to a freesurface
panel are, respectively
Ux ~ UO + AU Vy ~ VO + TV (27)
and
since WO = 0 for a freesurface panel. Through
Bernoulli's equation, the pressure on freesurface panels
is a function of local velocity, and is approximated by
retaining only firstorder incremental velocity terms
Vn = Wz ~ WO + AW ~ AW (28)
pa = ~ ~ 1 (U2 + Vy + W2) ~
~2 1(U20+V2)t{Uo~U+Vo~V }
~ 2 ~ ~  (U20 + V20 ~ ~  { Uo (Ux  Uo)
+Vo (Vy~ Vo) } (29)
Freesurface boundary conditions are linearized
in a similar manner, retaining only h~rstorder incremental
velocity and surfaceelevation terms. The kinematic free
surface boundary condition (18) is approximated by
Wz = Vn ~ UO 11X + VO lly ~ (UO + VO)1/2 also (30)
where the subscript so denotes differentiation along a
zero Fr streamline. The dynamic freesurface boundary
condition (19), after differentiation along so, and substi
tuting for 1lsO from (30), becomes
aSO Fr2 BUS + V2~1/2 (31)
704
OCR for page 699
A fivepoint backward difference is used in the ~ and ~
directions, and the freesurface grid metrics, are used to
compute the pressure gradient
COOP + vO8P
an ax ay
asO (U2+v2~ll2
=
uO (aP a; + ap all) + V (aP at + ap all)
at ax ~ ax at by Do, by
 (32)
(U2 + V2)1/2
The pressuregradient algorithm is structured to permit
the use of any blocked freesurface grid arrangement.
Also, using less than a f~vepoint backward difference
tends to dampen wave amplitudes. This wavedamping
mechanism is employed on panels near the outer
boundary of the finite freesurface model, so that wave
disturbances are eliminated before reaching the free
surface outer boundary.
At this point, a sufficient number of linear
dependencies have been established to permit the elimi
nation of the unknown freesurface source strengths in
(26), i.e., (24) relates source strength to panel normal
velocity, (31) relates freesurface panel normal velocity
to streamwise pressure gradient, (32) with backward
differences relates streamwise pressure gradient to free
surface pressures, (29) relates freesurface pressure to
freesurface panel tangential velocities, (25) relates panel
tangential velocities to the local surface gradient of dou
blet strength, and central differences relate the local
surface gradient of doublet strength to doublet strengths.
Hence, freesurface source strengths can be expressed as
a linear combination of freesurface doublet strengths,
A.
<,jaj +2 . bjk ink
w
Substituting for Hi from (33) into (26) yields
(i  ~ Aij ~j + ~ Bij ~j .
Sb + Sw sb
+v
sw
Bij (aj +2 bjk ink)
w
With freesurface source strengths eliminated,
and source strengths on the solid body evaluated
directly, solution of (34) yields the corresponding dou
blet strengths. The freesurface source strengths are then
given by (33), and (24) and (25) are used to compute the
resulting velocities on both body and freesurface panels.
Pressures on freesurface panels are given by (29~. A
similar linearized formula is used for pressures acting on
body panels, and configuration forces and moments are
obtained by panel pressure integration.
For interactive calculations, the SPLASH code
calculates the inviscid freesurface flow about the
equivalent displacement body resulting from the previous
viscous calculation. For this purpose, the equivalent
displacement body is treated as a solid fixed surface.
The inviscid flow velocities required for the next viscous
flow calculation, at offbody points on the viscous grid
outer boundary SO, are obtained using the computed
sourcedoublet solution and velocity influence
coefficients. A subpanel velocity influencecoefficient
algorithm was developed which utilizes a bilinear varia
tion of source and doublet strength across each panel.
The continuous variation of source and doublet strength
on each panel, and across panel edges, enhances the
accuracy of offbody velocity calculations at points close
to any body andlor freesurface panels.
WIGLEY HULL GEOMETRY AND
EXPERIMENTAL INFORMATION
The Wigley parabolic hull was selected for the
initial calculations since the geometry is relatively simple
and it has been used in many previous computational and
experimental studies. In particular, it is one of the two
hulls, the other being the Series 60 CB = .6 ship model,
selected by the Cooperative Experimental Program
(CEP) of the Resistance and Flow Committee of the
International Towing Tank Conference [27] for which
extensive global (total, wave pattern, and viscous resis
tance, mean sinkage and trim, and wave profiles on the
hull) and local (hull pressure and wallshearstress distri
butions and velocity and turbulence fields) measurements
were reported. It was for these same reasons that the
Wigley hull was selected as the first test case of the basic
viscousflow method [141, including comparisons with
some of the zero Fr data of the CEP. Herein, compar
isons are made for zero Fr with this same data and for
nonzero Fr with the appropriate data of the CEP. As will
be shown later, the nonzero Fr data is not as complete or
of the same quality as that for zero Fr, which was the
motivation for a related experimental study for the Series
60 CB  .6 ship model [28] for which calculations and
comparisons are in progress. However, the
comparisons are still useful in order to validate the
present interactive approach and display the
shortcomings of both the computations and experiments.
The coordinates of the Wigley hull are given by
y = 2~4x(1  x)~1  (z/d)2} (35)
where B = .1 and d = .0625. Waterplane and typical
crossplane views are shown in figure 2.
RESULTS
In the following, first, the computational grids
(figures 2 and 3) and conditions are described. Then,
some example results are presented and discussed for
zero Fr, followed by those for nonzero Fr, including
wherever possible comparisons with available experi
mental data, and, in the latter case, with inviscidflow
results. The convergence history of the pressure is
shown in figure 4. Figure 5 provides a comparison of
the largedomain and interactive solutions. The freesur
face perspective view and contours, wave profile, and
surfacepressure profiles and contours are shown in fig
ures 6 through 10, respectively. The axialvelocity con
tours, crossplanevelocity vectors, and pressure, axial
vorticity, and turbulent kinetic energy contours for sev
eral representative stations are shown in figures 11
through 13. Lastly, the velocity, pressure, and turbulent
705
OCR for page 699
kinetic energy profiles for similar stations are shown in
figures 14 through 16. On the figures and in the discus
sions, the terminology "interactive" refers to results from
both the interactive viscous and displacementbody
inviscid solutions. When the distinction is not obvious it
will be made. The terminology "inviscid" or "bare
body" refers to the noninteractive inviscid solution.
Computational Grids and Conditions
The viscousflow computational grid was
obtained using the technique of generating bodyfitted
coordinates through the solution of elliptic partial differ
ential equations. Because of the simplicity of the present
geometry, it is possible to specify the axial fit and cir
cumferential f3 control functions as, respectively, only
functions of I, and (; however, in order to accurately
satisfy the bodysurface boundary condition and resolve
the viscous flow, f2 = f2~t,ll,(~. Partial views of the
grids used in the calculations are shown in figures 2a,b
for a longitudinal plane and typical body and wake cross
planes, respectively. Initially, a largedomain grid was
generated. Subsequently, a smalldomain grid was
obtained by simply deleting that portion of the large
domain grid that lay beyond about r > .2. The outer
boundary for the smalldomain grid is shown by the
dashed line in figure 2. For the largedomain grid, the
inlet, exit, and outer boundaries are located at x =
(.296,4.524) and r = 1, respectively. The first grid point
off the body surface is located in the range 90 < y+ <
250. 50 axial, 30 radial, and 15 circumferential grid
points were used. As already indicated, the small
domain grid was similar, except 21 radial grid points
were used. In summary, the total number of grid points
for the large and smalldomain calculations are 22,500
and 15,750, respectively.
The inviscidflow displacementbody and free
surface panelization is shown in figure 3. 423 panels are
distributed over the displacement body and 546 over the
free surface for a total number of 969 panels. The panel
ization covers an area corresponding to 1 ship length
upstream of the bow, 1.5 ship lengths in the transverse
direction, and 3 ship lengths downstream of the stern.
This panel arrangement was judged optimum based on
panelization dependency tests E161.
The conditions for the calculations are as follows:
L= l;Uo= l;Re=4.5x 106;Fr=0and.316; end on
the inlet plane the average values for ~ and Us are .0033
and .0455, respectively. These conditions were selected
to correspond as closely as possible to those of the
experiments of the CEP with which comparisons will be
made [5,29,301.
Initially, largedomain calculations were per
formed for zero Fr. A zeropressure initial condition
was used and the values for the time a`, pressure ap, and
transport quantity Ad (where ~ = k and £)
underrelaxation factors and total number of global
iterations were .05 and 200, respectively. Next, small
domain calculations were performed, first for zero Fr,
and then for nonzero Fr. For zero Fr, the interaction
calculations were started with a zeropressure initial
condition and freestream edge conditions (Ue = 1,We=
Pe = 01. After 200 global iterations, the edge conditions
were updated using the latest values of displacement
thickness. Subsequently, the edge conditions were
updated every 200 global iterations until convergence
was achieved, which took three updates. For nonzero
Fr, the calculations were started with the zero Fr solution
as the initial condition and with nonzero Fr edge condi
tions obtained utilizing the zero Fr displacement body.
This solution converged in 200 global iterations. Most
of the results to be presented are for this case; however,
some limited results will be shown in which the nonzero
Fr edge conditions were obtained using an updated
nonzero Fr displacement body. The values for at, up,
and At (where ~ = k and £) used for the smalldoma~n
calculations were the same as those for the largedomain
calculations; however, for nonzero Fr, in addition, a
value of .01 was used for Al (where ~ = U) for grid
nodes near the outer boundary. The a~9/3z term in (19)
was found to have a small influence and was neglected in
many of the calculations; however, this may be due, in
part, to the present grid resolution. The calculations
were performed on the Naval Research Laboratory
CRAY XMP24 supercomputer. The CPU time required
for the calculations was about 17 minutes for 200 global
iterations for the viscousflow code and 1 minute for the
inviscidflow code.
Extensive grid dependency and convergence
checks were not earned out since these had been done
previously both for the basic viscousflow method [14]
and for other applications. However, some calculations
were performed using both coarser and finer grids.
These converged, respectively, more rapidly and slower
than the present solution. Qualitatively the solutions
were very similar to the present one, but with reduced
and somewhat increased resolution, respectively. The
convergence criterion was that the change in solution be
less than about .05% for all variables. Usually the solu
tions were carried out at least 50 global iterations beyond
meeting this criterion. Figure 4 provides the conver
gence history for the pressure and is typical of the results
for all the variables. In figure 4, the abscissa is the
global iteration number it and the ordinate is the residual
R(it), which is defined as follows:
imax imax
R(it)= Ad, ~ptitl)  pkit) i/ At, Pith ~ (36)
i=1 i=1
where ill and imax are the total number of iterations and
grid points, respectively. Referring to figure 4, global
iterations 1  200 correspond to the final iterations of the
zero Fr solution and global iterations 200  400 to those
for the nonzero Fr solution.
Zero Fr
Figure 5 provides a comparison of the zero Fr
largedomain and interactive solutions and experimental
data. The two solutions are nearly identical and show
good agreement with the data, which validates the pre
sent interactive approach. The agreement with the data
for the largedomain case is, of course, not surprising
since this was already established in [14] for a similar
grid and conditions, i.e., the present zero Fr solution is
essentially the same as that of [141. Some additional
aspects of the zero Fr solution are displayed in figures 11
through 15 for later comparison with the nonzero Fr
solution. Reference [14] provides detailed discussion of
the zero Fr solution, including comparisons with the
available experimental data. In summary, there is a
downward flow on the forebody and an upward flow on
the afterbody in response to the externalflow pressure
gradients. The boundary layer and wake remain thin and
attached and the viscousinviscid interaction is weak;
706
OCR for page 699
however, on the forebody, the boundary layer is rela
tively thicker near the keel than the waterplane, whereas
the reverse holds true on the afterbody and in the near
wake. The stern vortex is very weak. In the intermedi
ate and far wake, the flow becomes axisymmetric. As
indicated in figures 5 and 14 through 16, the agreement
between the calculations and data is quite good;
however, there are some important differences, which
are primarily attributed to the deficiencies of the standard
ke turbulence model with wall functions. In particular,
the axial velocity and turbulent kinetic energy are
overpredicted near the stern and there is a more rapid
recovery in the wake.
Nonzero Fr
Figure 5 also includes nonzero Fr results for
comparison. On the waterplane, the surface and wake
centerplane pressure displays very dramatic differences,
the wallshear velocity shows similar trends, but with
reduced magnitude, and the wake centerplane velocity
indicates faster recovery in the intermediate and far
wake. As will be shown later, the first closely follows
the wave profile, the second is due to an increase in
boundarylayer thickness near the waterplane for the
nonzero Fr case, and the third can be explained by the
waveinduced pressure gradients. On the keel, all three
of these quantities are nearly the same as for zero Fr.
The freesurface perspective views (figure 6) and
contours (figure 7) vividly display the complex wave
pattern consisting of both diverging and transverse wave
systems. The bow and stern wave systems are seen to
initiate with crests and the shoulder systems initiate with
troughs, which conforms to the usual pattern described
for this type of hull form. Very apparent is the reduced
amplitude of the stern waves for the interactive as com
pared to the inviscid solution. Also, the diverging wave
system is more pronounced and at a smaller angle with
respect to the centerplane. Note that the axial and trans
verse waveinduced pressure gradients can be discerned
from these figures, but with an appropriate phase shift,
i.e., increasing and decreasing wave elevations imply,
respectively, adverse and favorable gradients. The wave
profile along the hull is shown in figure 8, which, in this
case, includes experimental data for comparison. On the
forebody, the two solutions are nearly identical and
underpredict the amplitude of the bowwave crest and the
first trough. On the afterbody, the interactive solution
indicates larger values than the inviscid solution, with the
data in between the two. The wave profile for the
nonzero Fr displacement body (figure 3b) is also shown
in figure 8. The differences are minimal on the
forebody, whereas, they are significant on the afterbody
and depart from the data. It appears that the present
simple definition ( 1 ) is insufficient for "wavy"
displacement bodies.
The surfacepressure profiles (figure 9) show
similar tendencies as just discussed with regard to the
wave profile. On the forebody, the two solutions are
nearly identical, but, in this case, in very close agreement
with the data. The pressure on the forebody shown by
the dashed line is that obtained from the inviscid dis
placementbody solution. On the afterbody, here again,
the interactive solution indicates larger values than the
inviscid solution, with the data in between the two. The
waveinduced effects are seen to diminish with
increasing depth and the agreement between the two
solutions and the data on the afterbody shows
improvement. The surfacepressure contours (figure 10)
graphically display the differences between the two
solutions and the data. Note that the axial and vertical
surfacepressure gradients can be discerned from these
figures, i.e., increasing and decreasing pressure imply,
respectively, adverse and favorable gradients. The larger
wave elevation and pressure on the afterbody for the
interactive solution results in the closed contours near the
stern displayed in figure 10b. As already mentioned, the
viscousinviscid interaction is weak for the Wigley hull,
which is the reason that the inviscid and viscous pressure
distributions are quite similar. However, it appears that
the interaction is greater for nonzero as compared to zero
Fr.
Figures 11 through 13 show the detailed results
for several representative stations, i.e., x = .506, .904,
and 1.112, although the discussion to follow is based on
the complete results at all stations. Note that for zero Fr
the upper boundary shown is the waterplane, whereas
for nonzero Fr, it is the predicted free surface. Also, the
axialvelocity, vorticity, and turbulent kinetic energy
contours are not shown for the inviscid solution since, in
the former case, their values are all very close to 1 and,
in the latter two cases, they are, of course, zero. Solid
curves indicate clockwise vorticity.
On the forebody (figure 11), the boundary layer
is thin such that many aspects of the solutions are simi
lar; however, there are some important differences. The
nonzero Fr pressure fields show local and global effects
of the free surface, i.e., near the free surface, regions of
high and low pressure coincide with wave crests and
troughs, respectively, and at larger depths, the contours
are parallel to the free surface. Also, for nonzero Fr, the
crossplanevelocity vectors are considerably larger,
especially for the interactive solution. The inviscid solu
tion clearly lacks detail near the hull surface. The extent
of the axial vorticity is increased for nonzero Fr and is
locally influenced by the free surface. In both cases, as
expected, the direction of rotation is mostly anticlock
w~se.
On the afterbody (figure 12), almost all aspects
of the solutions show significant differences. The
boundary layer is thicker near the waterplane for nonzero
as compared to zero Fr. This behavior begins at x ~
.825, which coincides with a region of adverse axial
waveinduced pressure gradient (see figure 7~. The
differences for the pressure field and axialvorticity
contours are similar as described for the forebody;
however, in the case of the crossplanevelocity vectors,
there is an additional difference that near the free surface
the interactive solution displays downward flow. This is
consistent with the fact that the freesurface elevation is
above the waterplane and the pressure is generally higher
near the free surface than it is a larger depths, i.e., 11 ~ O
and ap/az ~ 0. Note that, as expected, in both cases, the
direction of rotation for the axialvorticity is mostly
clockwise. The turbulent kinetic energy contours are
nearly the same for both Fr.
In the wake (figure 13), the solutions continue to
show significant differences. Initially, the lowvelocity
region diffuses somewhat and covers a larger depthwise
region; then, for x > 1.2, recovers quite rapidly. A
similar behavior was noted earlier for the wake centerline
velocity for x ~ 1.2, both of which, as already
mentioned, are consistent with the wave pattern. The
zero Fr pressure field is nearly axisymmetric and fully
707
OCR for page 699
recovered by the exit plane. The nonzero Fr pressure
field continues to show freesurface effects, i.e., the
contours are parallel to the free surface, but also fully
recovered by the exit plane. Note the considerably larger
wave elevation near the wake centerplane for the inviscid
as compared to the interactive solution, which was
pointed out earlier with regard to figures 6 and 7. Here
again, the crossplanevelocity vectors are larger for
nonzero as compared to zero Fr, especially near the wake
centerplane for the interactive solution. The interactive
and inviscid solutions display differences near the free
surface, which appear to be consistent with the
differences in their predicted wave patterns. The zero Fr
axial vorticity decays fairly rapidly, whereas, for
nonzero Fr, the decay is slow with a layer of nonzero
vorticity persisting near the free surface all the way to the
exit plane. The turbulent kinetic energy contours are
similar for both Fr, but recover faster for the nonzero
case.
Figures 14 through 16 show the velocity,
pressure, and turbulent kinetic energy profiles for similar
stations as for figures 11 through 13, i.e., x = .5, .9,
and 1.1. Also, included are both zero and non zero Fr
experimental data. At the largest two depths, z = .05 and
.0625, data for both Fr are available, whereas, at the
waterplane, z = 0, only zero Fr data are available. At the
intermediate depths, data are available for both Fr, but
for different z values. Since the interest here is primarily
nonzero Fr and the zero Fr data and comparisons were
already displayed in [14], only nonzero Fr data are
shown for z = .0125, .025, and .0375. For zero Fr, a
corrected pressure is also shown which includes a
constant (= .03) referencepressure correction as
described in [141. Turbulent kinetic energy data are only
available for zero Fr.
At x = .5, consistent with previous discussions
the differences between the two solutions are quite small
and the agreement with the zero Fr data is good.
However, the nonzero Fr data show some unexpected
differences. In particular, the axialvelocity profile has a
laminar appearance and the boundarylayer thickness is
relatively large, the vertical velocity is upward, and the
pressure shows considerable scatter. It is pointed out in
[5] that the pressuremeasurement error was appreciable.
At x = .9 and 1.1, here again, consistent with
previous discussions the differences between the two
solutions are significant and the agreement between the
zero Fr solution and data is good, except for the
aforementioned discrepancies. The nonzero Fr solution
shows larger axial velocities than the measurements for
the inner part of the profiles. Here again, the measured
profiles have a laminar appearance and the boundary
layer is thick. However, no doubt, a part of the
difference is due to the calculations, i.e., as is the case
for zero Fr, due to deficiencies of the ke turbulence
model an overprediction of the velocity near the wall and
wake centerplane is expected. The transverse velocity is
small and with similar trends for both calculations and
measurements. The calculations indicate downward
vertical velocities near the free surface and upward
values for the midgirth region and near the keel. The
agreement with the data near the keel is satisfactory, but
in the midgirth region and near the free surface the data
display greater upward flow than the calculations. In the
wake, the nonzero Fr data show surprisingly small
vertical velocities near the wake centerplane. Here again,
the nonzero Fr pressure data shows considerable scatter
and is difficult to compare with the calculations.
Consistent with earlier discussions the turbulent kinetic
energy profiles are nearly the same for both Fr.
Lastly, Table 1 provides a comparison of the cal
culated pressureresistance coefficient and experimental
values of the residuaryresistance (i.e., total  frictional)
coefficient. The experimental values cover a range of
Re, including the present value, and clearly show a
dependency on Re. Interestingly, the inviscid result
compares well with the data at the highest Re, whereas
the interactive result is close to that that the data implies
at the present Re.
WAVEBOUNDARY LAYER AND WAKE
INTERACTION
The comparisons of the zero and nonzero Fr
interactive and inviscidflow results with experimental
data enables an evaluation of the waveboundary layer
and wake interaction. Very significant differences are
observed between the zero and nonzero Fr interactive
results due to the presence of the free surface and gravity
waves. In fact, the flow field is completely altered.
Most of the differences were explicable in terms of the
differences between the zero and nonzero Fr surface
pressure distributions and, in the latter case, the addi
tional pressure gradients at the free surface associated
with the wave pattern. The viscousinviscid interaction
appears to be greater for nonzero as compared to zero Fr.
It should be mentioned that other factors undoubtedly
have important influences, e.g., waveinduced
separation, which are not included in the present theory.
The interactive and inviscid nonzero Fr solutions
also indicate very significant differences. The inviscid
solution clearly lacks "realfluid effects." The viscous
flow close to the hull and wake centerplane is clearly not
accurately resolved. The interactive solution shows an
increased response to pressure gradients as compared to
the ~nv~sc~d solution, especially in regions of low
velocity. Also, the inviscid solution overpredicts the
pressure recovery at the stern and the sternwave
amplitudes.
CONCLUDING REMARKS
The present worlc demonstrates for the first time
the feasibility of an interactive approach for calculating
ship boundary layers and wakes for nonzero Fr. The
results presented for the Wigley hull are very
encouraging. In fact, in many respects, the present
results appear to be superior to the only other solutions
of this type available, i.e., [10,111. This is true both
with regard to the resolution of the boundarylayer and
wake regions and the wave field. Furthermore, it
appears that the present interactive approach is
considerably more computationally efficient than the
largedomain approaches of [10,111. This is consistent
with the previous finding for zero Fr [151. However, a
complete evaluation of the present method was not
possible. In the former case, due to the limited available
experimental data. As mentioned earlier, a related
experimental study for the Series 60 CB = .6 ship model
[28] was recently completed for which extensive
measurements were made at both low and high Fr for
which calculations and comparisons are in progress. In
the latter case, due to the considerable differences in
numerical techniques and algorithms and turbulence
models between the present methods and those of
708
OCR for page 699
[10,11]. As mentioned earlier, the pursuit of a large 10.
domain approach to the present problem is also of
Interest and will enable such an evaluation.
Finally, some of the issues that need to be
addressed while further developing and validating the
present approach are as follows: further assessment of
the most appropriate freesurface boundary conditions;
improved definition and construction of displacement
bodies; the inclusion and resolution of the bowflow
region; extensions for lifting flow; and the ever present
problem of grid generation and turbulence modeling.
Also, of interest is the inclusion of nonlinear effects in
the inviscidflow code.
ACKNOWLEDGEMENTS
This research was sponsored by the Office of
Naval Research under Contract N0001488K0113
under the administration of Dr. E.P. Rood whose sup
port and helpful technical discussions are greatly appre
c~ated.
REFERENCES
1. Lindenmuth, W., Ratcliffe, T.J., and Reed, A.M.,
"Comparative Accuracy of Numerical Kelvin Wake
Code Predicitions  "WakeOff"," DTRC/SHD
12601, 1988.
2. Rosen, B.,"SPLASH FreeSurface Code:
Theoretical/Numerical Formulation," South Bay
Simulations Inc., Babylon, NY, 1989 (proprietary
report).
3. Patel, V.C., "Ship Stern and Wake Flows: Status
of Experiment and Theory," Proc. 17th Office of
Naval Research S ymposium on Naval
Hydrodynamics, The Hague, 1988, pp. 217240.
Proc. 5th International Conference on Numerical
Ship Hydrodynamics, Hiroshima, 1989.
5. Shahshahan, A., "Effects of Viscosity on
Wavemaking Resistance of a Ship Model," Ph.D.
Thesis, The University of Iowa, Iowa City, IA,
1985.
6. Ikehata, M. and Tahara, Y., "Influence of
Boundary Layer and Wake on Free Surface Flow
around a Ship Model," J. Society of Naval
Architects of Japan, Vol. 161, 1987, pp. 4957 (in
Japanese).
Stern, F., "Effects of Waves on the Boundary
Layer of a SurfacePiercing Body," J. of Ship
Research, Vol. 30, No. 4, 1986, pp. 256274.
8. Stern, F., Hwang, W.S., and Jaw, S.Y., "Effects
of Waves on the Boundary Layer of a Surface
Piercing Flat Plate: Experiment and Theory," J. of
Ship Research, Vol. 33, No. 1, 1989, pp. 6380.
9. Stern, F., "Influence of Waves on the Boundary
Layer of a SurfacePiercing Body," Proc.4th
International Conference on Numerical Ship
Hydrodynamics, Washington, D.C., 1985, pp.
383406.
709
Miyata, H., Sato, T., and Baba, N., "Difference
Solution of a Viscous Flow with FreeSurface
Wave about an Advancing Ship," J. of
Computational Physics, Vol. 72, No. 2, 1987, pp.
393421.
11. Hino, T., "Computation of a Free Surface Flow
around an Advancing Ship by the NavierStokes
Equations," Proc. 5th International Conference on
Numerical Ship Hydrodynamics, Hiroshima,
1989.
12. Harlow, F.H. and Welch, J.E., "Numerical
Calculation of TimeDependent Viscous Flow of a
Fluid with Free Surface," The Physics of Fluids,
Vol. 8, 1965, pp. 21822189.
13. Chan, R.K.C. and Street, R.L., "A Computer
Study of FiniteAmplitude Water Waves," J. of
Computational Physics, Vol. 6, 1970, pp. 6894.
14. Patel, V.C., Chen, H.C. and Ju, S., "Ship Stern
and Wake Flows: Solutions of the FullyElliptic
ReynoldsAveraged NavierStokes Equations and
Comparisons with Experiments," Iowa Institute of
Hydraulic Research, The University of Iowa, IIHR
Report No. 323, 1988; also J. of Computational
Physics, Vol. 88, No. 2, June 1990, pp. 305336.
15. Stern, F., Yoo, S.Y. and Patel, V.C., "Interactive
and Large Domain Solutions of HigherOrder
ViscousFlow Equations," AIAA Journal, Vol. 26,
No. 9, 1988, pp. 1052 1060.
16. Longo, J., "Scale Effects on NearField Wave
Patterns," M.S. Thesis, The University of Iowa,
Iowa City, IA, 1990.
17. Black, R., "Definition of ThreeDimensional
Displacement Thickness Appropriate for Ship
Boundary Layers and Wakes," M.S. Thesis, The
University of Iowa, Iowa City, IA, expected 1991.
18. Hotta, T. and Hatano, S., "Turbulence
Measurements in the Wake of a Tanker Model on
and under the Free Surface," Fall Meeting of the
Society of Naval Architects of Japan, 1983.
19. Rodi, W., "Turbulence Model and Their
Application in Hydraulics," Presented at the IAHR
Section on Fundamentals of Division II:
Experimental and Mathematical Fluid Dynamics,
1980.
20. Swean, T.F. and Peltzer, R.D., "Free Surface
Effects on the Wake of a Flat Plate," NRL Memo
Report 5426, Naval Research Laboratory,
Washington D.C., 1984.
21. Ramberg, S.E., Swean, T.F., and Plesniak,
M.W., "Turbulence Near a Free Surface in a Plane
Jet," NRL Memo Report 6367, Naval Research
Laboratory, Washington D.C., 1989.
22. Maskew, B., "Prediction of Subsonic
Aerodynamic Characteristics: A Case for Low
Order Panel Methods," Journal of Aircraft, Vol.
19, No. 2, 1982, pp. 157163.
OCR for page 699
23. Maskew, B., "A Computer Program for
Calculating the NonLinear Aerodynamic
Characteristics of Arbitrary Configurations,"
NASA CR 166476, 1982.
24. Boppe, C.W., Rosen, B.S., Laiosa, J.P., and
Chance, B., Jr., ''Stars & Stripes '87:
Computational Flow Simulations for
Hydrodynamic Design," The Eighth Chesapeake
Sailing Yacht Symposium, Annapolis, MD., 1987.
25. Stern, F., "Comparison of Computational and
Experimental Unsteady Cavitation on a Pitching
Foil, ASME J. Fluids Eng., Vol. 111, 1989, pp.
290299.
26. Dawson, C.W., "A Practical Computer Method for
Solving ShipWave Problems," Proc. 2nd
International Conference on Numerical Ship
Hvdrodvnamics, Berkeley, CA., 1977, pp. 3038.
27. "Report of the Resistance and Flow Committee,"
Proc. 18th Int. Towing Tank Conf., Kobe, Japan,
1987, pp. 4792.
28. Toda, Y., Stern, F., and Longo, J., "MeanFlow
Measurements in the Boundary Layer and Walce
and Wave Field of a Series 60 CB = .6 Ship Model
for Froude Numbers .16 and .316," Iowa Institute
of Hydraulic Research, The University of Iowa,
IIHR Report No. xxx, 1990 (in preparation).
29. Sarda, O.P., "Turbulent Flow Past Ship Hulls 
An Experimental and Computational Study," Ph.D.
Thesis, The University of Iowa, Iowa City, IA.,
1986.
30. Kajatani, H., Miyata, H., Ikehata, M., Tanaka,
H., Adachi, H., Namimatsu, M., and Ogiwara, S.
"The Summary of the Cooperative Experiment on
Wigley Parabolic Model in Japan," Proc. 2nd
DTNSRDC Workshop on Ship WaveResistance
Computations, 1983, pp. 535.
Table 1. ResiduaryResistance Coefficients
L(m) T(°C) UO(m/s)
Fr Re CR
ExperimentIHI 6 12.8 2.423 0.316 11.9 x 106 1.803 x 103
Expenment SRI 4 10.6 1.978 0.316 6.14 1.998
Experiment UT 2.5 17.3 1.564 0.316 3.6 1.866
Inviscid ~0.316  1.79
Interactive ~0.316 4.5 x 106 1.92
710
OCR for page 699
DIVERGING WAVE
80W WAVE '~ST E R N WAVE
REGION 2. BODY SURfACE, Sb / / TRANSVERSE / /
Bow FLOW yt ~ \ / / WAVES ~ ~
Uo ~ /4~MMETRY PLANE, Sk x
S > ~\l~ I
(INTERACTION} REGION 3 ~ ~ ~ 1 ~ S. \
THIN BOUNDARY LAYER \ / o
REGION 4. ~tINTERACTtoN,
THICK BOUNDARY LAYER ~ WAKE
REGION 1. INVISCID FLOW
INLET BOUNDARY, S.
(LARGE DOMAIN: EXIT BOUNDARY, Se
OUTER BOUNDARY, S
(LARGE DOMAIN)
f
(a) (x,y)plane
(a) Fr = 0
O.
(b) Fr=.316
Figure 3. Displacement bodies.
a' REGION 1. INVISCiD FLOW
ARC
I' \ So ,/
_ {lNTERACTlON) /
OUTER BOUNDARY, So  ,/
(LARGE DOMAIN)
SYMMETRY PLANE, Si
REWATER PLANE, W
(b) (y,z) plane
Figure 1. Definition sketch of flowfield regions and
solution domains.
° VATERPLANE
o _ .
~ .

o
to
~ Fr=o
 F`r=0.316
%~% 
%% _
)
50 100 150 200 250 300 350 400 44 to
Global Iteration Number
Figure 4. Convergence history.
(a)
 Fr=0.318
 i'r=O:luLeracLive
O ['r=O:LargeDoulain
O t:xp.l'stmuff & J.
_9~
_ _ _ _
~_ _ _ _
0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0
X
(a) longitudinalplane
1=7 X= 0.~ff
1=~D X= 1.002
l
0.6 0.8 1.0 .2
 F.r=0,316
 F`r=O:lnteractive
O Fr=O:LareeDomain
O Esp.Ylatmuff dc J.
_~
, , ,
D.4 0.5 0.6 0.7 0.8 0.9 1.0
(C)
 FY=O.SIG
 F`r=O:lDteractive .
O Fr=O:LargeDomain
O Exp.Sarda
 I , T ,
0.00 0.05 0.10 0.15 0.20 0.00 0.05 0.10 0.15 0.20
Y Y
(b) body and wake crossplanes
Figure 2. Computational grid.
1.50 1.75 2.00
X
Figure 5. Comparison of interactive and largedomain
solutions for the watelplane: (a) surface and
wake centerplane pressure, (b) wallshear
velocity; and (c) wake centeIplane velocity.
711
OCR for page 699
(a) interactive ~111111111111111111111
(b) inviscid
Figure 6. Freesurface perspective view.
..~
0.2 0.4 O.B 0.8 1.0
X
(a) z/d=.04
 Present cal. Inviscid
 Present cal.VjBCOUB
O Bare body  Inviscid
O Exp.IHI
"o~ODOOODOaC
 Present cal. Invincid
 Present cal. Viscous
O Bare body  loviscid
O Esp.  IHI
[ ~&
l l
).0 0.2 0.4 0.6 0.8
X
(b) zld = .92
Figure 9. Surfacepressure profiles
0.2 0.0 O.Z 0.4 0.6 0.8 1.0 1.2 1.4 1.6
X
(a) interactive

or~
~,
o
= 0.001 "
X
(b) inviscid
Figure 7. Freesurface contours.
l
 Disp.body  O_Fr
  Disp. body  Non_O_Fr
O Bare body
O Esp.
~ =~__= =~;
I ~ 1 ~
0.0 0.2 0.4 0.6 O.B 1.0
X
Figure 8. Wave profile.
·c,
~0
1 ~
0.0 0.2 0.4 0.6
X
(a) expe~iment
..
0.0 o.z
~ . ~ .
: :
. .
. .
0.8 1.0
0.4 0.6 0.8
X
(c) inviscid
Figure 10. Surfacepressure contours.
.015
1.0
OCR for page 699
Interactive
(a)
O_
o
U~
o
o

(C) O
o
(d
(e)
o
.
o
o
C~
O
o
o

U]
o
o
...
,.
o

o
O ~
Inviscid
Fr = 0 Fr = .316 Fr = .316
F.:
Do.soo
(0.900
,. l l
0.00 0.05 0.10
(b o
,f ' ..
o_o:oi5
.. oo,olo
, I
~ U.05 0.10 0.15 0.20
. , ,"
. ,' ~ , ~
,.' .' ' ,
~' 0.1
. r
, ,
O.00 0.05 0.10 0.15 0.20
\~J/~
0.00 0.05 0.10
. .
0.00 0.05 0.10
`~de
/' , .,.' ' b,_'o,,05o.,
/ . b, 0 02S.
/'''' '
·  o~o;OZO
0.00 0.05 0.10 0 15 0;20
  .n n2n
~0,030
. . .
0_o.076
. . ,
1
0.1
r
, ', , ,
0.00 0.050.10 0.15 0.20 0.000.05 0.10 O.15 0.20
y
0.00 0.05
y
Figure 11. Comparison of solutions at x = .506: (a) axialvelocity contours; (b) pressure contours; (c)
crossplanevelocity vectors; (d) axialvorticity contours; and (e) turbulent kinetic energy contours.
713
OCR for page 699
(o)
1
all
(b) ~
o
~ .
 ,/
~ O
a
a
O .  ~. . O . . . .
0.000.05 0.10 0.15 O.ZO o.eo 0.05 0.10 0.15 O.ZO
(C) an_ _
(d 
o
0
. .
n an n no n In
. . .
O.UO 0.05 0.10
Y
Disco
Or = 0 ~ = .316 ~ = .316
.....
. . . .
O.00 0.()5 0.10 0.15 O.3O
/
o
O_
o
6 ~
O_ I 1
0.00 0.05 0.10
0.00 0.05 0.10 0.15 O.ZO
, .. . . . . . .
0.00 0.05 0.10 0.15 0.20 0.00 0.05 0.7 O.l5 O.ZO
~0
Y
Fugue 12. Carson of spurns ~ x = .9~: (a) ~i~vel=1~ cone; (b) Is cares; ~)
s~l=~ve~1~ vow; (d) =1~v~c1~ congas; ad (e) opulent anemic entry congas.
714
OCR for page 699
(a)
o
o
o
~
L~
o
o
o

O_ _
(b) o
o
o
o
o

o
in

o
Interactive
_. Inviscid
Fr= 0 Fr = .316 Fr = .316
7J
Al=
J
o
1 1 0_ ~
0 t1() 0.05 0.10 0.00 0.05 0.10
_._
_~005
it//
_ / /
i 15 /
.010
I I ~I
0.00 0.05 0.10 0.15 0.20 0.00 0.05 0.10 0.15 0.20
(C) o
o .
U:
_ .
. _
(d)
(e)
0.00 0.05 0.1 0 (). 15 C).''0
 1
0.00 0.05 0.10
y
. ~1
. ..
0.00 0.05 0.10 0.15 0.20 0.00 0.05 0.10 ().15 0. 'U
y
1
0.00 0.05 0.10
y
Figure 13. Comparison of solutions for x = .1.112: (a) axialvelocity contours; (b) pressure contours; (c)
crossplanevelocity vectors; (d) axialvorticity contours; and (e) turbulent kinetic energy contours.
OCR for page 699
l

 CAL. [Fr=0.316]
 CAL. [Fr=O]
EXP.Sarda [Fr=O]
EXP.Ali S. [Fr=0.313]
o
o
_ ~Z= 0.0000
o
_ ~~~~~ ~Z= 0.0125
0 ~
o
o
o
o
c
,~Z= 0.0250
Q':;~ z=0.0375
i~£~ z= 0.0500
_ ~ Z=0.625
0.00 0.01 O.OZ 0.03 0.04 0.05
y
o
u~
0
CAL. [Fr=0.31~;]
 CAL. [Fr=O]
O EXP.Sarda [Fr=O]
EXP.Ali S. [Fr=0.313]
~_ ~ Z= 0.0000
0.0250
Z= 0.03~5
Z= 0.0500
~X:
Z= 0.0625
CAL. [Fr= 0. 316]
 CAL. [Fr=O]
O EXP.Sarda [Fr=O]
~ EXP.Ali S. [Fr=0.313]
0
o
_ ~Z= 0.0000o
~c ~z= 0.01250
. ~o
~~°~~~~~~ ~Z= 0.02500
e ~Z = 0 . 0 3 7 50
__ .
~:A Z= 0~0500 0
Z=0.0625 0
cz
1 1 1 1 O
0.00 0.01 0.02 *0.03 0.04 0.05
y
CAL. [Fr= 0 . 316]
  CAL. [Fr=O]
O EXP.Sarda [Fr=O]
X EXP.Sarda:CORRECTED
o ~ EXP.Ali S. [Fr=0.313] o
o  o
0 ~0
0 .~_ Z= 0.0000 o
.~ ~ ,,< Z= 0.0125 o
o _~__ ~Z= 0.0250 o
0 0
0 ~ ~_0_ ~Z~ 0.0375 0
0 ~Z= 0.0500 0
>~ 0~__,
~w.~
0 ~Z= 0.0625 o
~ ~ ~ ~ ,~, .
O ~ ~=';7941~W' ~
~0
O 1 1 1 1 0_
0.00 0.01 O.OZ .0.03 0.04 0.05
y
CAL. [Fr=0.316]
  CAL. [Fr=O]
O EXP.Sarda [Fr=O]
~' Z= 0.0000
Z= 0.0125
~~ Z= 0.0250
Z= 0.0375
~ Z= 0.0500
~ Z= 0.0625
1 1 1 1
0.00 0.01 O.OZ *0.03 0.04 C) 05
Figure 14. Velocity, pressure, and turbulent kinetic energy proD~les at x = .5.
716
0.00 0.01 0 02 0.03 0.04 0.05
y
OCR for page 699
 
In

 CAL.tFr0.316]
 CAL. [FrO]
O FXP.Sarda [l`'r=O]
EXP.Ali S. [Fr=0.313]
~, V9~Z=O.OOOO
 1, _ Z=0.0125
~Z=0.0250
:;~ Z= 1).0375
~, :~Z= 0.0500
_ ~ Z= 0.0625
O 1 1 1 1
0.00 0.01 O OZ *0.1)3 0.04 0.05
oF
1
o
o
0
of
O F
o
o
CAL. [Fr=0.316]
 CAL. [Fi=O]
O EXP.Sarda [Fr=O]
X EXP.Sarda:CORRECTED
EXP.Ali S. [Fr=0.313]
~Z= 0.0000
~;
o

o
o
o
o
c~
~~~~~ Z= 0.0125
~Q^ ^~ Z= 0.0250
~Z=0.0375
~Z= 0.0500
 z= 0.0625
(}~
1 1 1 1
0.00 0.01 0.02 *0.03 0.04 0.05
y
 CAL. [Fr=0.316]
 CAL. [Fr=O]
O EXP.Sarda [F'r=O]
EXP.Ali S. [Fr=0.3 13]
1\ _
~a~
Z= 0.0000
Z= 0.0125
Z= 0.0250
~9~
Z= 0.0375
Z= 0.0500
9~
~ Z= 0.0625
,~
~ 0.00 0.01 0.02 *0.03 0.04 0.05
y
CAL. [F'r=0.316]
 CAL. [Fr=O]
O EXP.Sarda [Fr=O]
Z= 0.0000
 CAL. [Fr=0.316]
 CAL. [Fr=O]
O EXP.Sal da LFI =()]
EXP.Ali S. [Fr=n.313]
. ~Z= 0.0000
__ _ = Z= 0.0125
__~ Z= 0.0250
Z= 0.03~5
 ,~4Z= 0.0500
Z= 0.0625
. ~5L=1 1
O.00 0.01 0 02 0.03 0.04 0 05
*
y
~~ Z= 0.0125
Z= 0.0250
~Z=0.0375
= 0.0500
Z= 0.0625
1 1 1
0.00 0.01 0.02 .0.03 0.04 0.05
y
Figure 15. Velocity, pressure, and turbulent kinetic energy proD~les at x = .9.
OCR for page 699
o
.
o
o
 ar
o
o(

o
o
o 
o
o
_  N
o
o
C'AL. [Fr= 0 . 3 1 6]  CAL. [Fr= 0. 31 6]
 CAL. tFr=O]  CAL. [Fr=O]
': EXP.Sarda [Fr=O] O EXP.Sarda [Fr=O]
 i\ EXP.Ali S. [Fr=0.313]  ~A EXP.Ali S. [Fr= .313] 
~)Z= o.oooo I r =c Z= .0000 ~O.
Z=O.OIZ5 1 °~ Z= .0125 1 °
=oo250 1 0~:z= .0250 1 °~
:Q ~, z= 0.0375 1 °~=~\ z= 0375 1 °
~ ~:~ ~bo Z = 0 . 0 5 0 0 ~° ~ = / ~ ~ Z = . 0 5 0 0 ~
_ ~ ~ Z= 0.0625 o ~=~ ~ Z= 0.0625 o
0.00 0.01 0.02 0.03 0.04 0.()5 0.00 0.01 0.02 0.03 0.04 0.05
* *
Y Y
o

0 ~
o
o
o
o
o
o
P~
o
o
o
1
~4
o
o
o
o
o
o
o
o
o
o
o
o
1 1 1'1 1 °
0.00 0.01 0.02 *0.03 0.04 0.05
y
 CAL. [F'r=0.316]
 CAL [Fr=O]
O EXP.Sarda [Fr=O]
X EXP.Sarda:CORRECTED
~ EXP.Ali S. [Fr=0.313]
~ Z= 0.0000
Z= 0.0125
,:^,~ ~ Z= 0.0250
_ _ e _ _ _ _ _ _ _ ~
Z= 0~0375
~Xz= 0 0500
A Z= 0.0625
 CAL. LFr=0.3 1~3]
  CA L. lFr= 0]
O EXP.Sarda tPr=()]
EXP.Ali S. [Fr=0.313]
Z= 0.0000
~= Z= 0.0125
~= Z= 0.0250
Z= 0.03~5
~3~Z= 0.0500
~_ _Q ~ ~  A Z= 0.0625
. . . .
0.00 0.01 0.02 *0.03 0.04 0.05
y
CAL. tFr=0.316]
 CAL [Fi=O]
O EXP.Sarda [Fr=O]
~Z= 0.0000
Z= 0.0125
Z= 0.0250
Z 0.0375
~~
Z= 0.0500
 11~  O 0625
Z .
0.00 0.01 O.OZ *0.03 O.04 0.05
y
Figure 16. Velocity, pressure, and turbulent kinetic energy profiles at x = 1.1.
718
OCR for page 699
DISCUSSION
Kuniharu Nakatake
Kyushu University, Japan
I appreciate your interesting paper. In your calculation, the flow
field and pressure distribution around wave surface are not taken into
account. If you distribute source panel on the calculated wave surface
and determine its strength from the kinematic condition on the wave
surface, you may obtain more plausible results. This is possible in
the linearized framework.
AUTHORS' REPLY
We agree that a higherorder treatment of the freesurface boundary
conditions would be preferable. We hope to make progress in this
area in the future.
DISCUSSION
Hoyle Raven
Maritime Research Institute Netherlands, The Netherlands
This valuable paper addresses the difficult problem of prescribing
freesurface boundary conditions inside the viscous domain. The
authors solve the wave elevation from integration of the kinematic
condition, and thus from the velocities at the undisturbed free surface.
Prescribing boundary conditions for these velocities here is, therefore,
critical. If I understand it correctly, in (20) not only the stresses but
also the associated strain rates have been neglected, which are of
leading order in the wake elevation. This is not a definite solution
of course. What alternatives for the freesurface condition do you
consider?
AUTHORS' REPLY
The approximations (a) through (d) utilized in deriving equations (18)
through (20) have been clearly stated in the text, and therefore do not
require further explanation. An alternative treatment of the free
surface boundary conditions within the present overall framework is
to retain more terms in the approximations (a) through (c) and higher
order terms in approximation (d).
DISCUSSION
William B. Morgan
David Taylor Research Center, USA
In the authors' presentation of their problem they stated that they
were using an "Unsteady RANS Code using a Kit turbulence model."
I believe this statement is inconsistent in that a RANS Code with a
kin model cannot be unsteady. I can understand attempting to use
this code in a Quasisteady" sense, but I believe it is not applicable
to the unsteady problem. Would the authors please comment?
AUTHORS' REPLY
Although not discussed in the present paper, we recognize the
limitations of the kit turbulence model with regard to simulating
unsteady flow. Issues concerning this point were discussed in one of
our earlier papers on another topic [31]. In our presentation, we only
wished to emphasize the general capability of the IIHR basic viscous
flow method for unsteady flow notwithstanding the limitations of
current turbulence models for such applications.
[31] Stern, F., Kim, H.T., Patel, V.C., and Chen, H.C., "A
ViscousFlow Approach to the Computation of PropellerHull
Interaction,. Journal of Ship Research, Vol. 32, No. 4, December
1988, pp. 246262.
DISCUSSION
Kazuhiro Mori
Hiroshima University, Japan
1. It must be important to be of the same order in approximations
when the viscous effects are taken into account in the freesurface
computation. From this standpoint of view, the use of the
displacement thickness method is not consistent where the thickness
is calculated exactly. This may be crucial when the 3D separation
is dominant. My suggestion is that the viscous flow should be taken
into account as the double model flow directly in the inviscid
computation.
2. According to our computation and experiment, the separation of
the flow at stern is much affected by the bow wave system. This
means that we should not expect precise discussions on the interaction
between the viscosity and the free surface by the iterative procedure
as done in the present study.
AUTHORS' REPLY
We are not clear as to the precise meaning of your questions;
however, we would like to point out that viscous effects have been
included directly in the inviscidflow computation through the use of
the displacement body. As discussed in the text, the limitations of
this approach are not yet known. The present results are
encouraging, but a complete evaluation requires further validation
through comparisons with experimental data and a large~omain
approach. Work along these lines is in progress.
We thank both the oral and written discussers of our paper for their
pertinent remarks.
719
OCR for page 699