|
Items in cart [0] |
Questions? Call 888-624-8373 |
| Copyright © 2009. National Academy of Sciences. All rights reserved. Terms of Use and Privacy Statement |
Below are the first 10 and last 10 pages of uncorrected machine-read text (when available) of this chapter, followed by the top 30 algorithmically extracted key phrases from the chapter as a whole.
Intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text on the opening pages of each chapter.
Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.
Do not use for reproduction, copying, pasting, or reading; exclusively for search engines.
OCR for page 185
New Viscous and Inviscid CFD Techniques for Ship Flows
L. Larsson~ and L. Broberg
FLOWTECH International AB4, Sweden
K. J. Kim
Daewoo Shipbuilding & Heavy Machinery LTD, Korea
D. H. Zhang
Wuhan University of Water Transportation, China
ABSTRACT
The research on CFD for ship flows
carried out at SSPA and Chalmers Uni-
versity of Technology in recent years
is summarized. Methods for potential
flow calculations including first and
higher order theories with linear or
non-linear free surface boundary con-
ditions are presented. The importance
of higher order effects and non-
linearity is discussed. A viscous flow
method based on the Reynolds-averaged
Navier-Stokes equations is also intro-
duced. The method may be run either in
the partially parabolic mode or fully
elliptically. A comparison is made
between results of the two modes. Grid
dependence studies and validation
against experiments are presented for
all methods.
1. INTRODUCTION
Computational Fluid Dynamics (CFD)
has been a major research area at SSPA
for a number of years. During the
seventies several computer programs for
calculating laminar and turbulent
boundary layers were developed. These
methods turned out to be very useful
for thin boundary layers but failed
completely near the ship stern [1].
During the eighties a strong effort has
therefore been made to develop more
- Also Chalmers University of Techno-
logy
Presently at FLOWTECH International
AB
Presently at SSPA Maritime Consul-
ting AB
FLOWTECH International AB is a SSPA
subsidiary specialized in CFD for
external flows.
185
accurate methods, such as the higher
order boundary layer integral method
[2], [3], the ADI method in body-fitted
coordinates [4] and the streamline
curvature method [5]. All of these
produced improved results near the
stern, but none of them could accura-
tely predict the near wake and the
viscous/inviscid interaction. A more
accurate approach to these problems was
obviously required.
The situation was anticipated around
1980, when the development of a Navier-
Stokes method was started. After
testing and analysing various routes to
find the best possible coordinate sys-
tem and solution procedure a method was
developed in 1985-1988. Several diffe-
rent programs were written based on
somewhat different approximations of
the Navier-Stokes equations [6], [7]
[8] and [9]. Introduction of the pro-
peller effect has been made recently
[10].
In the early eighties interest was
also directed towards the free surface
potential flow, and the first methods,
based on Dawson's theory, were deve-
loped during 1983-1986 [11], [12],
[13]. The call for more accurate solu-
tions in certain cases prompted the
development of methods with a more
exact free surface boundary condition
and better numerical accuracy. These
methods are presented in detail in
[14], [15], [16], [17], [18] and [19].
The purpose of the present paper is to
introduce the new development, to give
some examples of the application of the
programs and to compare the results
from different levels of approximation,
in the potential flow as well as in the
viscous flow.
OCR for page 186
2 . POTENT IAL FLOW METHODS
2.1 Governing equations
In this class of methods, the flow is
considered steady, inviscid, irrota-
tional and incompressible. A right-
handed coordinate system Oxyz is
employed with the origin on the mean
free surface, x positive in the
direction of the uniform flow, and z
positive upwards. A ship, piercing the
free surface, is assumed to be in a
uniform onset flow of velocity U.-. (See
Fig 1.) The flow field around the ship
may be described by a velocity
potential ¢, which is generated by a
distribution of sources on a surface S
and by the uniform onset flow in the x-
direction.
z fx
u=~ ~
p{fj,llj,Cj)
Fig 1 Coordinate systems in poten-
tial flow methods
¢(x,y,z) = ~ o(q)/r(p,q) dS + U-x
S
where o(q) is the source density on the
surface element dS and r(p,q) is the
distance from the point q to the field
point p(x,y,z) where the potential is
being evaluated.
The potential ~ given in Eq (1) is
governed by the Laplace equation.
V2¢ = 0
in the fluid domain and satisfies the
regularity condition at infinity
(2)
Vat = > (U-,O,O) as r -- (3)
The source density ~ is to be
determined from the boundary condi-
tions on the hull and free surface. On
the wetted hull surface the solid
boundary condition is
On = 0
where n denotes the outward normal.
(4)
At the free surface, two boundary
conditions must be imposed, i e the
flow must be tangent to the surface
xhx + Achy _ ~ = 0
and the pressure constant
gh + 2 ( Vat Vat - UP ) = 0 (6)
Further, no upstream waves must be
generated.
The exact problem formulated above is
nonlinear, since the equations (5) and
(6) are nonlinear and are to be satis-
fied on the wavy surface z=h(x,y),
which is unknown. In the present work
an iteration procedure is applied and
the free surface boundary condition in
each iteration is linearized about the
previous solution. When the process
converges the exact solution is
approached.
Unknown sources ~ on the hull and
wavy surface z=h(x,y) will induce a
potential ~ and a wave elevation h,
which satisfy the boundary conditions
(5) and (6).
Dl(~'h) = ~xhx + Why - Liz - O
D2(o,h) = h ~ 2g [ UP - (¢x+¢y
+¢z ) ] = 0
Small perturbations 5o and ah of the
previous solution are introduced, and
the equations expanded to first order
in a Taylor series.
Dl(a,h) ~ Dl(~°,h°) + ADl(~'h°)
+ ~Dl(~°,h)
~ Dl(~°,h° ) + aa Dl(Cl,h° )~o
ah D1(~ ,h)6h ~ 0
D2(a,h) ~ D2(o ,h ) + ~D2(~'h )
+ AD2(~°,h)
~ D2(o ,h ) + at D2(~,h° )
+ a D2(~°,h)ah - O
186
(8)
OCR for page 187
Here the superscript, °, corresponds
to the previous solution, which is
assumed to be known. Thus,
D1(~° h°) = ~xhx + Why - fez
D2(~° h°) = h° - 2g [Up - (~x
+ By + Liz)]
where ~ is the potential induced by
a° on the free surface z = h°(x y).
Introducing the perturbation 5¢, the
partial increments of D1 and D2 may be
written
AD1(o h°) = 6¢xhx + 5¢yhy Fez (10)
AD1(~° h) = ~X5hx + Why + (~xzhx
+ Achy - ~zz)5h
AD2(o h ) = g(~x5¢x + byway Z z
(11)
(~° h) = ah + g(~x~xz + ~Y.YZ
+ ~z~zz)6h
Therefore the linearized free surface
boundary conditions are
xhx + Why - fez + ~X5hx + Achy +
(~Xzhx+~yzhy ~zz)5h = 0 (12)
{ 1 + g(.x~xz+.y~yz+~z~zz) }ah =
1 { U 2 ~ 2 ~ 2 ~ 2
(~x6¢x+~y5¢y+~z6¢Z) } -h°
(12) and (13) are to be applied on
z=h°(x,y)
(13)
To start the iteration process the
double model case is first computed.
This yields ~ and ho, and a linear
solution ¢, h may be obtained from
(12) and (13). However, for this
particular case Liz, as well as z-
derivatives of fix and by are equal to
zero, for symmetry reasons, so (12) and
(13) are simplified to
h ° + ~ h ° - ~ + ~ ah
+ Why = 0
(14)
Sh - {u2+~x+~y - 2(~x~x+~y~y)} / 2g
These are the equations solved in
the linear methods. They may be applied
either on z=0 [11], [19] or on
z=h°(x,y) [15], the so called Bernoulli
wave.
To obtain a solution satisfying the
exact boundary conditions (non-linear
in ~ and h) the iterations continue
from the linear solution using (12) and
(13). After 5-10 iterations ah and 5¢
are usually sufficiently small, and the
process may be considered converged.
2.2 Numerical Solution
2.2.1 Free Surface Grid
In the numerical implementation the
first step is to discretize the
integral domain S in (1) into a large
number of small panels. The domain
consists of the hull surface and part
of the free surface. In principle the
free surface panel distribution should
be extended over the entire region
where there are significant waves, but
for practical reasons the region of
free surface panels has to be limited
to the area near the hull. The
influence of the truncation will be
discussed later.
To generate the free surface grid,
two different principles have been
used. Following Dawson [26], Xia Ill]
generated a grid based on the stream-
lines of the double model flow. This
approach has the advantage that a
somewhat simpler form of the boundary
conditions (14) and (15) may be used in
the linear case. Dawson showed that the
equations can be rewritten as ordinary
differential equations in l, when l is
the arc length along a streamline.
Using a streamline grid thus enables
simple numerical differentiation
(Dawsons derivation was not quite
correct, but the error is small).There
is, however, one serious disadvantage
of the streamline grid, as noted by Xia
& Larsson [12]: the grid gets very
coarse at the ends of the hull, where a
fine grid is required. Therefore
another approach has been chosen in the
more recent methods t12]-[19]. The
longitudinal lines are now generated
numerically by interpolation between
the edge of the panelled region and the
hull. The transverse lines can be
187
OCR for page 188
either hyperbolas [11]-[17] or straight
lines [18], [19]. A typical grid is
shown in Fig.2. For non-linear calcu-
lations the free surface panels are
moved to the wavy surface in each
iteration.
tax
Fig 2 Panel distribution. SSPA
Ro-Ro ship, medium bulb
To obtain derivatives of a function
in the x and y directions, derivatives
along the L (longitudinal) and T
(transverse) directions are obtained
from
gx = aXlgT - aX2gL
gy = gT
where
all = fL
r ~
ax2 = - i1 + fL2
and y = fL(x) is the equation for a
longitudinal line.
(16)
Derivatives with respect to L and T
are obtained either from three point or
four point finite difference operators
gLi = CAigi + CBigi_NL + CCigi-2NL
gTi = GAigi + GBigi+l + GCigi+2
or
gLi = CAigi + CBigi_NL + CCigi-2NL
+ CDigi_3NL
gTi = GAigi + GBigi+1 + GCigi+2
+ GDigi+3
(17)
where NL is the number of longitudinal
strips on the free surface.
The coefficients of the upstream
operator (CA, CB, CC and CD) and the
downstream operator (GA, GB, GC, GD)
are calculated from the distance bet-
ween successive control points on the
free surface panels, see [11].
Dawson conjectured that the use of
upstream differences along the stream-
lines should prevent upstream waves.
This has been true of all cases com-
puted by the present methods even
though the streamlines have been re-
placed by the body-fitted longitudinal
lines.
During the iteration process in the
nonlinear procedure the calculation
domain is changed in each iteration.
The free surface sources are moved to
the previous wavy surface in every step
and the boundary condition is applied
there. This means that the hull surface
has to change as well, so a special
panelization procedure has been intro-
duced in the program. A typical panel
distribution is shown in Fig. 3.
Fig 3
Hull panel distribution in
the non-linear method. Series
60, CB = 0.60. Fn = 0.32
2.2.2. Higher Order Effects
Xia's methods tll]-[14] were of first
order i.e. the panels on the hull and
free surface were flat and with a con-
stant source strength. To obtain a more
accurate solution without increasing
the number of panels Ni tl5]-[18]
introduced a higher order technique in
which the panels are parabolic and the
source distribution linear. The tech-
nique, which is similar to the one pro-
posed by Hess t20], has been further
improved by Kim [19].
A parabolic panel is defined in
following form
(' = ZO + Ao: + BOn
+ po:~2 + 2Qo:'nt + ROn'2
(18)
where a panel element co-ordinate
system (:', n ' , ~ ~ ) is constructed using
the four corner points. Thus, the
188
OCR for page 189
origin is defined as the average of the
corner points, while the ;'-direction
is normal to the two diagonal vectors.
The six parameters (Zo,Ao,Bo,Po,Qo and
Ro) are determined by: (a) requiring
the panel to pass through the corner
points of the panel (four conditions),
and (b) requiring the panel to pass as
closely as possible, in a least squares
sense, to the eight additional input
points for the adjacent panels (two
additional conditions).
In the present method the control
point is the point on the curved panel
closest to the average point of the
four input points. This control point
is characterized by the condition that
the vector from the origin is parallel
to the local normal vector
[: ei+n ej+( (:',p')ek] x t-(' t
-A 'n .ej+ek] = 0
Eq (19) is equivalent to the two
scalar equations
G(:' ,n ' ) = ~ '+; '(a ' ,n ~ )~: ~ (;',n ' ) = 0
(20)
H(~t,n ~ ) = n t+` ~ (a ~ ,n ~ )¢n ~ (a ~ ,n ' ) = 0
These nonlinear equations are solved
by Newton-Raphson iteration.
Once the control point is determined,
the panel co-ordinate system (I', n ' , ~ '
is transformed to a new projected flat
panel which is tangent to the parabolic
panel and the tangent point is both the
control point and the origin of new
panel co-ordinate system (gin, ~ ) see
Fig 1. The equation of the panel may be
written
F(~'n'() = (_ [p:2 + 2Q:n + Rn2]
= 0
The source density distribution is
assumed to be linear on each panel.
() = ~° + ~~: + Ann (22)
This is fitted in the least squares
sense to the values of source density
at the control points of the M adjacent
panels. Thus the source derivatives on
the panel in question may be expressed
in terms of the unknown values of the
source density on the adjacent panels
in the form
M (I)
=ok k
M (a )
= ~ Ck Ok
n k=o
(23)
The desired source density coefficients
Cal) and Cin) are obtained by
minimizing the error function
Min ~
M
(l9) ct?)Ckn ) k-o[ k (~°+~:k+~n~k)]2 (24)
According to Eq (1) the perturbation
potential Bid at the i-th field point
( ~ i' n i';i) induced by the j-th panel on
which the source density ~ is
distributed is
ij = |Aj(~/r) dS (25)
Since in the higher order method the
integration is to be carried out over a
curved surface, expressions for 1/r and
dS are different from the corresponding
first order ones. The expression for
Aid is further complicated by the
variation of ~ on each panel. In [18]
the following expression for Ail is
derived
hit = .(j°) j + { P.ijP) + 2Q~.jQ)
+ Ouija) }aj + ~(jl6)~ + ~(1~)~
= A* ~ + ~ { C( ~ )~( 14 )
(21) J J k=1 k lo
+ Ckn)~ij1n)}~k (26)
Here M is the number of adjacent
panels.
189
OCR for page 190
All the ~ terms in (26) can be
calculated numerically and they may be
interpreted as follows: ()
corresponds to a flat panel with a
constant source density, .(~), ~(Q)
and () are caused by the parabolic
panel shape, () and ~(\n] come from
the linear variation of the source
density. All these terms depend only
on the geometry of the j-th panel, but
the curvature terms P. Q. R and the
source derivative coefficients Cog) and
con) depend on the surrounding panels.
as
NE
· = ~
X1 j=1
The individual It's in Eq (26) are ~yi
Fiji) = iAj r dads
Fiji ) = jAj (i;~/rf dads
Fiji) = iAj (i~j~j/rf dada (27)
(R) =
iJ
|Aj (ing/r3f dads
(14)= |Ajfj/rf duds
~ij(ln) = lAj nj/rf dads
where rf is the distance between
(ti' ni' (i) and the j-th control point
( Id, nj, O) on the flat panel (see Fig
1). Then the velocity induced by the
panel is obtained by taking the
gradient of the corresponding
potential.
At= ~~aj + E {ck(~)~(l4)+ck(n)~(ln)}~k
n nisi +k[-l{ck(~)~(l()+ck(n)~(ln)}~
(28)
E1{Ck(~)~(l;)+Ck(n)~(lM)
The subscript ij is omitted in the
equation for simplicity, *and new
velocity terms (~*, .*, ~ ) are intro-
duced. They are velocity Components
induced at the i-th control point by a
unit source density on the j-th panel.
All the induced velocities in (28) can
be obtained analytically.
The velocities the reference
coordinate can be calculated from the
coordinate transformation
:~x) all a21 all ~ 4]
Bye = al2 a22 a32 ~ ~~ t (29)
zJ al3 a23 a33 l rJ
These velocities may also be written
am a21 a31
a19 all aR2
Xi j (Sj ~ Uco
NE
~ Yi-~-
j=1 J J
NE
~ Zing.
j=1 J J
(30)
Here the velocity influence coeffi-
cients Xij, Yij and Zij are the
velocity components in the reference
coordinate system (x, y, z) at the i-th
control point, induced by a source
density which is unity at the control
point of the j-th panel. In the first
order method the velocity induced by a
panel depends only on the panel itself.
The essentially new feature of the
higher order method is that the
velocity induced by a panel depends on
the values of source density also at
the control points of adjacent
elements. Thus the influence
coefficients for a panel depend not
only on the geometry of that panel but
also on the geometry of adjacent
panels.
The second order derivatives of Fig
required in (12) and (13) are more
complicated than the first order ones.
In the present work, the magnitude of
the second derivative terms related to
the curvature and linear variation of
source density are assumed to be small
and vanish rapidly during the
iteration. Therefore the second order
derivative terms are calculated from a
source velocity which corresponds to
the flat panel with constant source
density.
xz NE also
~yz j-1 a23 J
fizz i a33 j
| all a21 a3l
~13 a23 an
190
fall a22 a32 ~ nt Ann Ant (31)
; din ~~
ij
OCR for page 191
NE is here again the total number of
panels.
2.2.3 Solution Procedure
The boundary conditions (4) and (12),
(13) or (14), (15) may be transformed
into a set of linear equations in ~
using the relations above for the first
and second derivations of A.
[ A ] { ~ } = { B } (32)
where tA] is a NE x NE matrix and {B}
is the right-hand side array.
The upper part of the A-matrix and B-
vector, corresponding to the hull
boundary condition, is obtained in a
straightforward way, while the lower
part, corresponding to the free
surface, is considerably more complex.
The equations are given in [18] and
[19]. Unfortunately, the lower part of
the matrix is not diagonally dominant,
so iterative methods cannot be used in
the solution. This is therefore found
by a Gaussian elimination procedure.
2.2.4 Wave Resistance
Having obtained the source strengths,
velocities may be computed from (30)
and the pressures from the Bernoulli
equation. To obtain the resistance, two
principally different methods have been
employed. Either the force is found by
pressure integration over the hull or
by integration of the sources over the
free surface. In the former case,
different formulas are used depending
on whether or not the higher order
technique is employed, i.e.
NB NB
W i CPi NXi Nisi / ~ Nisi (33)
where NB is the number of panels and Nx
and ASi are the x-component of the
normal and the area respectively of a
panel, or
NB
JC i E ( P° PC ~ Purl)
(a11N(; + a21Nrl + a31N; ) (34)
|NB
( 1+2~2 )d~dq| ~ |Ei ( 1+2¢ 2 )d~dn
The latter formula is derived
assuming that the pressure varies
linearly over the panel in a similar
way as the source strength. Cp: and Cp~
are the slopes of the pressure
variation, all, a21 and a31 are the
first column of the transformation
matrix between the panel and reference
coordinate systems, N; E, N. E and N: E
are the components of the Formal vector
in the panel coordinate system and ~ 2
is a term which takes into account the
change in panel area due to curvature.
The formula for obtaining the
resistance from the free surface
sources is derived from the momentum
theorem applied to a control volume
outside the hull, under the undistrubed
free surface and extending to infinity
sidewards, downwards and longitudi-
nally. The detailed derivation is given
in [18]. It yields
NF NB
Cw = (4n ~ UBici~Si) / (USA ASi) (35)
1 i
It should be noted that this formula
is valid only for the linear case.
2.3 Validation
The methods presented have been
validated in a number of ways. A grid
dependence study was first carried out
for the linear first order method to
investigate the sensitivity of the
solution to the panelization on the
free surface. Most likely the results
from this study are valid also for the
other methods. The improvement in
accuracy when introducing the higher
order method was then tested against
two analytical results without a free
surface.
Comparison with measurements have
been made for a number of cases such as
the Wigley hull, the Series 60, CB=0.60
hull, the HSVA-tanker, the high speed
ship ATHENA, a Ro-Ro ship designed at
SSPA and two sailing yachts. All these
comparisons are described in detail in
[11]-[19]. As an example only the
calculations for the Ro-Ro ship will be
presented here. This case is parti-
cularly interesting since a study of
several bulb alternatives was made. It
is also the only case for which calcu-
lations have been carried out using the
first and higher order linear, as well
as the higher order non-linear method.
Interesting comparisons between all
three methods may thus be made.
2.3.1 Grid Dependence Study
Since one of the most important
questions in connection with the
present methods is the dependence of
the solution on the panel distribution,
particularly on the free surface, a
systematic variation was carried out.
191
OCR for page 192
The test case was the Wigley hull.
Employing the coordinate system defined
in Fig 1, the equation for the hull
surface reads
y = 10 ( 1 - x2) ( 1 ~ Z64 ) (36)
In table 1 (last page) the number of
panels on the hull and on the surface,
the extent of the panelled part of the
free surface in the x- and y-
directions, the Froude number and the
wave resistance coefficient are given.
In run No 1, only approximately 1/4
of the half ship length, 1, was used
for the free surface extent in the y-
direction. The calculation gave an
unresonably high wave resistance at Fn
= 0.266 and broke down at Fn = 0~45. In
runs of No 2, 3 and 4, the free surface
region was enlarged and divided into 5,
10, and 15 strips respectively. At Fn =
0.266, the resistance changed signi-
ficantly from run to run, but at Fn =
0.45 there was no change between the
last two runs.
The free surface extent in the x-
direction was tested in runs No 3, 5,
6, 7 and 8. From No 3 to No 5 the free
surface region was extended by 1/4 1.
The wave resistance did not change very
much, but as the free surface region
was enlarged by 1/2 from run No 6 to
No 7, the wave resistance changed,
particularly for the highest Froude
numbers. From run No 7 to run No 8,
part of the free surface region was
moved from the bow to the stern. The
wave resistance did not change at Fn =
0.266, whereas it was increased by 10% viva
at Fn = 0 45
There were no systematic changes in
the panel distribution on the hull, but
two slightly different hull panel
distributions were used in the
successive tests. This should not cause
large changes in the final results. The
considerable difference in the wave
resistance between runs No 3 and No 6
is likely to have been caused by the
change of the panels on the free
surface. More even panels were placed
on the free surface for run No 6 than
those for case No 3.
The results obtained imply that the
choice of the panels depends on the
Froude number to be evaluated. Smaller
panels should be used for low Froude
numbers, while a larger free surface
portion is required for high Froude
numbers. It is wise to vary the panel
arrangement for different Froude
numbers.
A panel distribution similar to case
No 9 of Table 1 has been used in most
work afterwards.
2.3.2 Analytical Test Cases
To investigate the improvement in
accuracy when introducing the higher
order technique the flow around a
sphere was calculated. In Fig 4 the
velocity distribution along a generator
through the stagnation point is shown.
Two grids were tested, one with a
distribution similar to a typical ship
case, i.e. with 22 longitudinal and 40
transverse strips on one quarter of the
sphere and one with 5xlO strips. It is
seen that both higher order solutions
are extremely accurate, while there is
some error also for the 880 panel first
order case. Small as this error may
seem, it may have a significant in-
fluence on the drag. Integrating the
pressure on one half of the sphere
yields the following results for the
drag coefficient CX
First order Higher order Exact
Cx 0.0550 0.0622 0.0625
Error 12% 0.5% -
~ Anolytic solution
+ First order solution
0 Higher order solution
First order solution
O Higher order solution
.0
0.5
r'
B80 effective
panels
50 effective
panels
,.
V:,': ~
v<
, , ~
30 60 (3 90
Fig 4 Velocity distribution on a
sphere
192
OCR for page 193
The higher order calculation is thus
considerably more accurate. It should
be noted that the computer time for a
given number of panels increases by
about 10% for a case of this size.
z.3.4 The SSPA Ro-Ro Ship
The body plan of the SSPA Ro-Ro ship
is shown in Fig 5 with a small, a
medium and a large bulb. A comparison
between the measured residuary
resistance and the computed wave
resistance (pressure integration) using
two linear (first and higher order)
versions and one non-linear (higher
order) version of the same program [19]
is made for the medium bulb in Fig 6.
The panel distribution is shown in Fig
2.
Fig 5 Body plan, SSPA Ro-Ro ship
with different bulbs
C W X 1 0 3
3.0,
2.0
1.0
O _
Fig 6
0 0.1
-—— Measurements
Linear solution (first order)
· Linear solution (higher order)
· Non - I inear solution
_.
.1
,~^
0.2 En
i
l
,, 0.3
Predicted wave resistance
compared with measured
residual resistance. SSPA
Ro-Ro ship, medium bulb
A common problem in most evaluations
of wave resistance calculations is that
the wave pattern resistance is seldom
measured. In the present case the form
factor has been determined using a
Prohaska plot, and the viscous
resistance based on the ITTC-57 corre-
lation and the form factor has been
subtracted from the total resistance to
get the residuary resistance.
Cwx10-3
0.6 r
0.5
~ J
0.4
03
02
01
F.,
0.15
021
025
1 ~
Medi um
bulb
Otto
0~0
O 1 2 3 4 5 6 7
I terati on
Fig 7 Convergence history. Wave
resistance of SSPA Ro-Ro
. ship, medium bulb
As appears from Fig 6 the difference
between the linear and non-linear
higher order methods is quite small for
this case. This is also apparent from
Fig 7 where the convergence history of
the wave resistance is shown. The final
converged value after 7 iterations is
not very much different from the first
value, corresponding to the linear
solution. Larger differences have been
noted for bluffer ships [16] and non-
linearity should be extremely important
for hulls with large overhangs such as
sailing yachts.
193
OCR for page 194
2h Calculated
0.02 ..
Bulb
/;
/>\\ ~
~ _,
-0.01
No bulb
L Measured
0.02
Hi
o.olL
I
I
I
,0 2X/L
0;0 2 X/L
Fig 8 Wave profiles at En = 0.21.
SSPA Ro-Ro ship
CWX 10-3
an,
.
2.0
1.0
Type Measurement Calculation
No bulb
Small bulb ___
Medium bulb
Large bu l b
1.
11'
hi:
/i,1
:/'1
~/~///
0.1
an
Fig 9 Wave resistance. SSPA Ro-Ro
ship
The first order solution of Fig 6 is
not as accurate as the higher order
solutions, but even the first order
results are not grossly in error.
Interesting comparisons between
results for the three bulbs and the
case without a bulb can be made in Fig
8 and 9 based oh higher order linear
solutions. The measured wave profiles
of the four cases are very much
different as appears from Fig 8, and
the differences are quite well
predicted. The absolute values of the
wave resistance, Fig 9, are not as
accurate in all cases, but the method
is able to rank the cases in the right
order.
3 VISCOUS FLOW METHODS
Although r as indicated in the
Introduction' a number of different
methods for computing the viscous flow
(boundary layer/wake) around ships
have been developed at SSPA/CTH, only
the Navier-Stokes methods will be
presented in the present review.
3.1 Coordinate systems and grid
In the SSPA Navier-Stokes methods
the coordinate systems of Fig 10 are
employed. These systems, which differ
from the ones of the potential flow
methods, are more convenient for tensor
notation which is required for the
fully transformed Navier-Stokes
equations in curvilinear nonortho-
gonal coordinates. ;~ are the global
Cartesian reference coordinates, while
xi are the grid oriented ones, as
appears from Fig 10. The relation
between the two systems is of the form
V2xi = fi (i = 1, 2, 3) (37)
where v2 is the Laplacian.operator in
orthogonal coordinates (ha) and fir are
control functions which may be assigned
appropriate values to yield the desired
stretching of coordinate surfaces. The
equations are inverted by exchanging
the dependent and independent
variables
siiaiaj:1 + fkak~l = V251 (38)
The Laplacian operator may be written
2 1 3 hlh2h3 ~
V hlh2h3 j z- 1 a a; j (
2 - ~
194
(39)
OCR for page 195
: ~ ~
. 1 1 1 1 1 1 1 1
Fig 10 Coordinate systems in Navier-
Stokes methods
here hi, h2, h3 are metric coefficients
in the frame id. After selecting
41(xl) and assigning proper values for
the control functions fin, equations are
solved numerically. A more detailed
description of the grid generation is
given in [9]. Fig 11 shows the
coordinate system for the SSPA 720
Model.
2 Governing equations
The governing equations solved are
the Reynolds equations obtained from
the Navier-Stokes equations by
averaging over the largest, turbulent
time scale. In this way the
unsteadiness due to turbulence is
removed at the expense of having to
introduce a turbulence model for the
Reynolds stesses. In the coordinate
system above the incompressible
momentum equations may be written (see
[6] for a derivation).
Otui + h(k)ukuik = - h(i)gikakP +
(gjkulj + h(i)h(k)g u,j)SkvE
VE{gjkakajui + gikOk(ulQjil) + gjkQ
gj kr1 ui }
J , (40)
q ~
F
l
1 l
1 2
1 0 l
n l
0 ~
0 4
0.2
View of the water plane
Fig 11 Longitudinal and transverse
views of the grid. SSPA
720 Model
The incompressible continuity
equation may be written
J 1ak(Jh(k)uk) = 0 (41)
where J is the Jacobian
us are the physical components of the
velocity vector, i.e. the components
based on unit vectors parallel to the
covariant base vectors go, whose
magnitude igi~ is denoted h(i), gij is
tlhe contravariant metric tensor and
r i j is a Christoffel symbol of the
second kind. a i is to be interpreted as
a partial derivative with respect to xi
and ul,j is defined from the relation
195
OCR for page 200
- The pressure is computed from (60)
using the same algorithm.
- The transport equations for k and
are solved using the line by line tech-
nique.
- The same steps are repeated at the
next time level.
In both methods the solution is
assumed converged when the global
divergence is below a critical value.
3.4 Boundary Conditions
The calculation domain is shown in
Fig 11.
(1) Upstream and initial condition
The profiles of the velocity
components us and the turbulent
quantities k and ~ are required at the
inflow section. Only standard boundary
layer profiles for ul, k and ~ (U2 = us
= 0) have been used in the calcualtions
so far. No condition for the pressure
is required.
Lagging has been used during the
first time step for the velocity and
the turbulent quantities and the
pressure field is assumed to be con-
stant and equal to zero.
(2) Outflow condition
In the present calculations the
outflow section is located in the far
wake and a zero pressure gradient
condition alp = 0 is used.
(3) Outer edge
The outer plane may be located far
from the body, where the flow is nearly
undisturbed and thus
p = 0, ul = 1,u3 = 0 , a2k = a2 ~ = 0
and the normal velocity component u2 is
calculated from the local continuity.
Alternatively, p, u1 and us may be
taken from a potential flow solution.
(4) Symmetry planes
The undisturbed water surface and the
centerline plane are treated as symme-
try planes, i e
us = 0, 33u1 = 33u2 = ask = a3E = 0
(5) Wake center plane
The following conditions have been
applied at wake center planes that are
degenerated to a single line
200
u2 = us = 0, a2u1 = a2k = a26 = 0
otherwise
U2 = 0, 32U1 = a2u3 = a2k = a2c = 0
For the body boundary condition there
are two alternatives: either the
computational domain may be extended
into the viscous sublayer and the no-
slip condition applied (this requires a
two-layer turbulence model) or the wall
law may be used. In the former case the
conditions are simply
ul = u2 = us = k = 0
If wall functions are used a
procedure similar to the one of Chen &
Patel [21] has been employed. The
procedure requires that the two points
closest to the wall are located in the
logarithmic layer. It is assumed that
the velocity at the innermost point is
parallel to the wall which yields
u2 = 0 at x2 = 2 and that the velocity
vector does not rotate between the
first and second point in the plane
parallel to the wall. The last assump-
tion yields
( q ~ x2=2 =(q ~ X2=3
(q X2=2 =(q- x2=3 (63)
where q is the magnitude of the velo-
city.
The procedure is iterative and starts
by calculating the magnitude of the
velocity q2 at x2 = 2 with an assumed
skin friction velocity us from the wall
law
q2 = ~ 1 lny+2 + B
where K= 0.42,B = 5.45 and y+ = Re use.
The velocity components at x2 = 2 can
then be calculated from equation(63).
The components serve as boundary condi-
tions for the momentum equations. q3
can be calculated from the solution and
an updated u can be calculated from
the wall law applied at point x2 = 3.
The procedure is repeated until conver-
gence.
The boundary conditions at x2 = 2 for
the k- and E-equations are
uT2
k2 = ~
up
62 =
~Y2
OCR for page 201
3.5 Validation
3.5.1 Grid dependence study
A study was carried out to investi-
gate how the solutions depend on the
grid and what the proper number of
clusters is in each coordinate
direction for ship stern flow calcu-
lations. Due to computer limitations
the study had to be rather restricted,
however.
The systematic variations were made
using four different grids to investi-
gate the grid dependence in the three
coordinate directions. For all four
grids the calculation domain was the
same, that is the inlet plane was
placed at the midship, 2x/L = 0.0; the
outlet plane at 2x/L = 10 in the far
wake; the outer edge boundary on a cir-
cular cylinder located one ship length
from the centre line; and the two
innermost grid points were located
within the logarithmic layer.
In case 1 the grid has a number of
clusters equal to LLxMMxNN = 60 x 19
x 9, where LL, MM and NN2are the number
of clusters in the x1, x and x3
directions respectively.
The number of clusters in the
transverse direction x3 in case 1 might
be too small to capture the variation
of the cross flow and the rapid pres-
sure change near the region of the
keel. Thus in case 2, NN was increased
to 15. The grid in case 2 has the
number of clusters equal to LLxMMxNN
= 60 x 21 x 15.
In case 3, the grid was designed to
test the grid dependence in the pre-
dominant flow direction x1, LL was
decreased to 30. The grid had the
number of clusters equal to LLxMMxNN =
30 x 21 x 15.
The velocity profiles and the pro-
files of other physical quantities in
the normal direction might be quite
influenced by the number of clusters in
the normal direction x2. Based on the
grid of case 3, we increased the number
of clusters in the normal direction to
MM = 41, getting the grid case 4: LLx
MMxNN = 30 x 41 x 15.
A summary of the four cases is given
in Table 2.
Case LL
1 60
2 60
3 30
4 30
MM NN
19 9
21 15
21 15
41 15
Total number of
clusters
10 260
18 900
9 450
18 450
Table 2. The four grid variation cases
The calculations for the four cases
were carried out using standard boun-
da2ry later profiles for ul, k and ~
(u = u = 0) as the inlet profiles.
\ Of
2X/L = 0.8
o
N
l
n
__
i
a
80X19X9
- - 80X21X15
30XZiXi5
...... 30X41X15
, l l l
.00 0 . 05 0. 10 0. 15 0 . 20
GIRTH LENGTH
2X/L= 0.8
— 60XI9X9
— — BOX21X15
· 30X21X15
..... - poX41X15
o . Go 0.05 0. 10 0. 15 0.20
GIRTH LENGTH
Fig 13 Grid dependence. Distribution
of Cp and up at 2x/L = 0.8
The girthwise pressure and friction
velocity variations at 2x/L = 0.8 for
the four grids are shown in Fig 13. The
influence of the cluster numbers in the
transverse direction can clearly be
seen from the figure. Comparing the
grid case 1, 60 x 19 x 9, with other
grids it is obvious that nine clusters
are too few to resolve the problem in
the transverse direction. The Cp
variation is almost the same for the
other three grids, since the same
number of clusters was used in the
transverse direction. Owing to the
201
OCR for page 202
limitation of the computer size we did
not increase the number of clusters in
the transverse direction to more than
15. So it is difficult to say whether
or not the solution would be improved
if more than 15 clusters were used.
This number seems to be the minimum
required, however.
O _ . ~X4iX15
- - 30~1Xi5
o 0° Waterl i ne
N _ ~
O DO ~ 1.00
U
0 - - X41X15
- - ~~lX15
45° /
0 00 0.~ t.00
U
-
~X41Xi5
_ ~~iXi5
90° Keel
o.oo 0.50 t.oo
u
it, ~
Grid 60 x l9 x 9
I'
/////
/
-
/
Grid 60 x 21 x 15
////// /
Grid 31 x 41 x 15
/~/~
/////
/
11//
A//
Grid 31 x 21 x 15
Fig 14 Grid dependence. Velocity
profiles at 2x/L = 0.8 Fig 15 Grid dependence. Iso-wakes at
Waterline, mid-girth and keel 2x/L = 0.9
202
OCR for page 203
The influence of the number of
clusters in the normal direction can be
tested by comparing the results of grid
3 and grid 4. The total velocity pro-
files at 2x/L = 0.8 along the water
line, keel and at the bilge are shown
in Fig 14. There are small discrepan-
cies between the results for the two
grids. Also in Fig 13 the girthwise up
variation of the grid 30 x 41 x 15 is
only slightly different than the
others. This means that there is a very
small grid sensitivity in the normal
direction. If we look at the isowakes
from the results of the four grids at
2x/L ~ 0.91 shown in Fig 15, they
correspond quite well with each other.
The pressure and friction velocity °
variations along the water line and the
keel from the results of grid 2 and 3
are shown in Fig 16. Grids 1 and 4 are
not shown in this figure, since their
(xl, x2)-planes are slightly different.
60X2iXt5
30X2iXt5
C1
CO
o
n
CO
O _
o
O _
tu
O _
,,\
Watertine
1 1 1 1 1 1 1 1 1
0.0 O.Z 0.4 0.6 O.B 1.0 1.2 I.4 i.B 1.8 2.0
2X/L
60X2iX15
—· 30X21X15
Keel
0 0 0.2 0.4 0.6 O.B 1.0 t.2 1.4 1.6 1.B 2.0
2X/L
~ 1 1
0.0 O.Z 0.4 0.6 O.B
ZX/L
o
o
*
,_, o
~ A
o
CU
o
80XZlXIS
—30XZlX15
Waterline
1 1 1
1.0
60XZiX15
· 30XZlX15
Keel
1 1 1
0.4 0.6
ZX/L
O.B 1.0
Fig 16 Grid dependence. Cp and up
along the waterline and keel
The influence of the number of
clusters in the xl-direction is tested
from this comparison. The results are
almost the same in the whole region.
Increasing the number of clusters from
30 to 60 gives only very small changes.
It is therefore concluded that 30
clusters is enough in this direction.
These limited systematic studies of
grid-dependence indicate that there is
some grid sensitivity in the cross
section, but not in the longitudinal
direction. Special attention must be
paid to regions of abrupt change in
geometry, such as the stern region and
near the keel. Concentration of
clusters to these regions may be an
adequate alternative to fine grids.
Considering this possibility a grid of
30 x 20 x 15 would be proper for
practical use.
203
OCR for page 204
3.5.2 Location of the outer edge of the
computational domain
As indicated above, the outer edge of
the computational domain can be located
so far from the body that free stream
conditions may be specified. Another
alternative is to obtain velocities
(except for the normal component) and
pressure from a potential flow solu-
tion. However, outside the boundary-
layer there is a region in which the
viscous and inviscid flow interact. To
apply the boundary condition here an
iterative procedure between the poten-
tial and viscous solutions must be in-
troduced. To investigate the location
of the region of viscous-inviscid
interaction for the SSPA Model a series
of calculations was carried out using
different locations of the outer edge.
The largest distance to the edge mea-
sured from the keel is 1.08 (dimension-
less by L/2) for a grid with 60 x 19 x
15 clusters. The distance was varied by
dropping clusters in the normal
direction so that the distance was re-
duced to 0.51, 0.24, 0.11 and 0.05
respectively. The last grid consists of
60 x 11 x 15 clusters. For each calcu-
lation the boundary conditions at the
outer edge and at the inflow (outside
the boundary layer) were calculated
from the potential flow solution. The
resulting pressure-fields are
normalized at a point by subtracting a
constant in the following comparisons.
Fig 17 shows pressure profiles at four
different x-stations, distinguished by
symbols. The symbols are placed at the
corresponding outer edge. The line
types distinguish different calcula-
tions, i e different locations of the
outer edge. The profiles shown in Fig
17 (N = 14) are from the region close
to the keel where the boundary layer is
very thin, and the profiles shown in
Fig 17 (N = 7) are from bilge region,
where the boundary layer is thick. Both
figures show that the profiles from the
calculation with the outer edge located
at 0.05 are different from the others,
indicating that the boundary conditions
have been applied in the interaction
region. Thus, the interaction region
extends out to a surface which is
located between the outer edges of the
calculation domain 0.11 and 0.05.
A calculation of the potential flow
and the laminar and thin turbulent
boundary-layer normally precedes the
stern-flow calculation. The cost for
using the potentialflow condition is
therefore low. The gain would be a
reduction of clusters in the normal
direction and a convenient way to
handle blockage-effects in tunnel
simulations.
Variation of outer surf ace (N-14)
2X/L-0.7
ZX/L-O .8
2)(/L-0 .8
2X/L-0.95
· 2Are/L-t . 03
— 2Are/L-0.51
2Are/L-0.24
- 2Are/L-O.~t
~ 2Are/L-0.05
0.0 O. 1 O.Z
Cp
N
Var I at ion of outer sun f ace tN-7)
-
x.
Or
ID
0 _
~ _
x 2X/L-0 .7
~ 2X/L-D . ~
+ 2x/L-o . g
0 2X/L-0 .85
· 2Are/L-I. OB
2Are/L-0 . 51
— — 2Are/L-0.24
2Are/L-0 . ~ 1
- 2Are/L-0 . 05
Be.! 0.0 0.t 0.2
Cp
-
Fig 17 Variation in outer edge
location. Cp distribution in
the x2 direction
3.5.3 Comparison between oarabolized
and elliptic solution
The comparion between the
parabolized and elliptic methods was
carried out for the SSPA 720 Model with
a grid consisting of 60 x 19 x 15
clusters. The inlet profiles were
standard boundary-layer profiles,
normal and transverse velocity
204
OCR for page 205
components were zero and a constant
velocity u1 = was outside the boundary-
layer. The solution after one sweep
from the partially-parabolic method was
used as input for the fully-elliptic
algorithm.
The pressure and friction velocity
variation in the girthwise direction at
2x/L = 0.8 are shown in Fig 18. The
partially-parabolic and fully-elliptic
algorithms give almost the same
results. No noticeable difference
between the solutions can be seen.
This is not surprising since the only
simplification made in the parabolized
algorithm, as mentioned before, is that
the second derivative ( a 1 a 1 ¢) in the
predominant direction is neglected.
This higher order term is of signi-
ficant importance only for bluff stern
flow with separation, but not for
slender sterns, like the one of the
SSPA model. Thus the two algorithms
should give the same results.
The CPU time for the partially-
parabolic algorithm is about 70% of
that of the fully elliptic one.
O SSPA720 Cp AT 2X/L - O.B
0 0 0 . 1 0 . ~ 0 . 2 0 . 2 N
GIRTH LENGTH T
0.0 0.2 0.4 0.6 O.B 1.0 1.2 t.4 1.6 I.8 2.0
2X/L
cat
o
l
0
lo
lo
0
lo
lo
-
O
0 _
ELLIPTIC
— — Pi"B~IC
SSPA720 Ut AT ZX/L - O.B
~
ELLIPTIC
— — PARABOLIC
cam
L,
' l l 1 lo
0 0 0.! 0.1 0.2 0.2 ,o
GIRTH LENGTH r
Fig 18 Comparison between results of
elliptical parabolized
methods
3.5.4 Comparison with experimental
data
The grid 60 x 21 x 15 was used to
make the final calculations. The inlet
profile at midship (2x/L = 0.0) for the
longitudinal velocity component U
within the boundary layer was deter-
mined from the wall law using the
measured boundary layer thickness and
friction velocity from Larsson. The
normal component V and the transverse
component W were assumed to be zero.
Outside the boundary layer the
distribution of U was calculated using
Hess & Smith's potential flow calcula-
tion method. The k and ~ within the
boundary layer were calculated
according to Klebanoff & Bradshaw as
presented in [24] and set to zero out-
side.
The calculated results for the Reynolds
number Re = UoL/ = 5.0 x 106, were
compared with the experimental data of
Larsson [25]. Owing to the different
coordinate systems used in the
calculations and the experiments it is
difficult to make direct comparisons of
all quantities of interest. Therefore
only limited comparisons are presented
here.
fir
o SSPA720 Cp AT O DE6.
mu
lo
o
0 _
~ I~
tic \
- 60X21Xt5
0 SSPA720 Cp AT 90 DE6.
ru
0 _
o
0 _
—60X2IXI5
l
0.0 0.2 0.4 0.6 O.B I.0 1.2 I.4 1.6 t.8 2.0
2XIL
Fig 19 Comparison with experiments.
Cp distribution along water-
1lne and keel
205
OCR for page 206
The pressure distribution on the hull
surface along the water line and the
wake centerline as well as along the
keel are prsented in Fig 19 along with
the measurements. The original
experimental data of Larsson[25] are
subject to windtunnel blockage effects.
Corrections have been introduced in the
figures using the formula given by
Larsson[25], which was formulated after
carrying out potential-flow calcu-
lations with and without tunnel
constraints. Fig 20 shows the girthwise
pressure distributions at three
stations. It is clear from these
figures that the agreement between
calculations and measurements is quite
good.
o
° r SSPA720 Cp AT 2XJL - O .70 l
._ ~/?
c MEASUPEH£NT
60X2IXt5
O 00 0.05 0.10 0.15 O. 20
GIRTH LEN6TH
SSPA720 Cp AT 2X/L - O.BO
- ...
1° _
to
N
O' I _
0.00 0.05
~ MEASUREMENT
- ~nx2lxts
_ . _
0.10 0.15 0.20
GIRTH LENGTH
SSPA720 Cp AT ZX/L - O.S
o
to
o
can
to
o
0 _,
o
to
up SSPA720 Ut AT O DEG.
BOX21X15
0.2 0.4 0.6 0.8 1.0
Z)(/L
SSPA720 Ut AT 90 DEB.
· Rnx~~Xls
010 O'Z 0'4 0~6 ·11 110
2X/L
Fig 21 up distribution along water-
line and keel
Fig 21 shows the calculated wall
shear velocity distributions along
water line and keel and in Fig 22 the
girthwise distributions at three
transverse cross sections are given
along with the measurements. It can be
seen from the figures that in the stern
region the calculations overpredict the
wall shear velocity along the water
line, but underpredict it along the
keel. This indicates that the
calculated velocity profiles in these
regions might differ somewhat from the
measured profiles. It is difficult to
make comparisons of velocity profiles,
owing to the different coordinate
systems used in calculations and
experiments, as mentioned above.
SSPA720 Ut AT 2X/L - O . 70
g
o
to
_~
\ lo/ ~ htE AS UREME NT
\_/ 60X2IXt5
n, 00 O .05 0. 10 0 . 15 O. ZO
BIRTH LEN6TH
Fig 20 C distribution at three
stations
. 00 0 .05 0 . 10 0. 15
GIRTH LEN6TH
2~
0.20
OCR for page 207
SSPA7ZO Ut AT ZX/L - O.BO
to
o
o
\
O MEASUREMENT
60XZiXt5
:~ ~ ~
O MEASUREMENT
0 60XZiXt5
to
0 I , I
0.00 0.05 0.10 0. t5 0.20
6IR1H LENGTH
172~ ~ AT 2Xh _ n P' Fig 23
s~
BOX2IX15
0.00 0.05 0.10 0.15 O.
GIRTH LENGTH
Fig 22 up distribution at three
stations
Wake contours are given for two
stations in Fig 23. A reasonable
correspondence with measurements i.
noted at 2x/L = 0.8, while at 0.9 the
outermost iso-wake is displaced
outwards. Close to the keel the
boundary layer thickness is over-
predicted at both stations. The results
were obtained with the wall law as the
inner boundary condition. Preliminary
results from the latest version of
NASPAR, where the wall law is removed,
indicate that the prediciton of the
flow near the keel can be much im-
proved in this way.
Measurements Colculotion
Iso-wakes measurements x=O.R u/U=0.6, 0.7, 0.8, 0.9, 0.95
calculation x=0.81 q/U=0.6, 0.7, 0.8, 0.9, 0.95
Measurements Colculation
:~
/
Iso-wakes measurements x= 0.9, u/U=0.6, 0.7, 0.8, 0.9, 0,95, 1.0
calculation x=0.&9, q/U=0.6, 0.7, 0.8, 0.9, 0.95
Iso-wakes at two stations
4 FINAL COMMENTS
Methods for computing free surface
potential flows as well as complex
viscous flows have been presented. The
development has taken place over a long
period of time and includes five PhD
projects. Several computer programs
have been written, but they are now
combined into one program for the
potential flow and one for the Navier-
Stokes flow. These flow programs,
together with a boundary layer code
(for laminar and turbulent flows as
well as transition) represent a
powerful tool, named SHIPFLOW, for
investigating the flow and resistance
properties of ships. The SHIPFLOW
system is also equipped with state-of-
the-art pre- and postprocessors, which
makes it easy to use. For preliminary
project investigations the system
represents an alternative to tank
testing, i.e. it could be used as a
numerical towing tank. It should be
noted that there are several features
of the programs which have not been
presented in this survey.
Such features are the transom stern
option and the sinkage and trim calcu-
lation in the potential flow as well as
the propeller effect in the Navier-
Stokes code. Quite interesting is also
the most recent development, in which a
mathematical optimization scheme is
linked to the flow programs, enabling
optimization of hull forms from a
resistance point of view to be carried
out automatically. Owing to space limi-
tations these features have had to be
left out, but they are explained in the
references [10], Ill] and Elf].
207
OCR for page 208
run no
panel arrangement
(half range)
hull free surface
Table 1. Grid dependence study.
Wigley hull.
extent of the free surface F
(in half ship length) n
x-direction y-direction
1 24 x 6 34 x 5 -1.50 to 1.50 0.28 0.266
0.450
2 24 x 6
3 24 x 6
34 x 5 -1.50 to 1.50 0.775
34 x 10 -1.50 to 1.50 0.775
4 24 x 6 34 x 15 -1.50 to 1.50 0.775
5 24 x 6 36 x 10 -1.75 to 1.75 0.775
6 22 x 6 36 x 10 -1.50 to 1.50 0.775
7 22 x 6 32 x 10 -2.00 to 2.00 0.775
8 22 x 6 40 x 10 -1.50 to 2.50 0.775
0.266
0.313
0.266
0.313
0.350
0.400
0.450
0.266
0.450
0.220
0.310
0.350
0.450
0.266
0.450
0.266
0.450
0.266
0.450
9 22 x 6 40 x 10 -2.00 to 2.00 0.775 0.266
References
1.
Larsson, L. (editor), "The SSPA- 4.
ITTC Workshop on Ship Boundary
Layers 1980", SSPA Publication No.
90 (1981).
2. Larsson, L., & Chang, M.-S.,
"Numerical Viscous and Wave
Resistance Calculations Including 5.
Interaction", 13th Symposium on
Naval Hydrodynamics, Tokyo (1980).
Kim, K.-J. & Larsson, L.,
"Comparison Between First and
Higher Order Methods for Computing
the Boundary Layer and Viscous
Resistance of Arbitrary Ship
Hulls", International Symposium on
Resistance and Powering Perform-
ance, Shanghai (1989).
v
w
x 103
1.204
0.598
1.437
0.785
1.495
1.212
1.849
3.054
0.847
3.054
0.603
1.551
1.199
3.046
0.859
2.788
0.931
3.317
0.931
3.662
0.913
(0.966)
3.168
(3.235)
Broberg, L. & Larsson, L., "A
Calculation Method for Ship Stern
Flows Using and Analytic Body-
Fitted Coordinate System", 15th
Symposium on Naval Hydrodynamics,
Hamburg (1984).
Larsson, L. & Johans son, L.-E., "A
Streamline Curvature Method for
Computing the Flow near Ship
Sterns", 14th Symposium on Naval
Hydrodynamics, Ann Arbor (1982).
Broberg, L., "Calculations of
Partially Parabolic Stern Flows",
SSPA Report 2801-2 (1988).
7.
208
Broberg, L., "Numerical Solution
of the Time Averaged Navier-
Stokes Equations for Ship Stern
Flow", SSPA Report No. 2803-1
(1988).
OCR for page 209
8. Broberg, L., "Numerical Calcula- 18.
tion of Ship Stern Flow", Ph.D.
Thesis, Chalmers University of
Technology (1988).
9. Zhang, D.H., Broberg, L. &
Larsson, L., "A Numerical Method
for Generating Smooth Body-Fitted
Coordinate Systems Suitable for
Ship Stern Flow Calculations",
International Symposium on
Resistance and Powering Perform-
ance, Shanghai (1989).
10. Zhang, D.H., "Numerical Prediction
of Propeller-Hull Interaction in
Viscous Flow", SSPA Report No.
2963-1 (1989) .
11. Xia, F., "Calculation of Potential
Flow with a Free Surface", SSPA
Report No. 2912-1 (1984) . See also
Report No. 65, Dept of Marine
Hydrodynamics, Chalmers University
of Technology, Gothenburg, Sweden
(1984) .
12. Xia, F. & Larsson, L., "A
Calculation Method for the Lifting
Potential Flow around Yawed,
Surface-Piercing 3-D Bodies", 16th
Symposium on Naval Hydrodynamics,
Berkeley ( 1986) .
13. Xia, F., "Numerical Calculation of
Ship Flows with Special Emphasis
on the Free Surface Potential
Flow", Ph.D. Thesis, Chalmers
University of Technology ( 1986) .
14. Xia, F., "A Study on the Numerical
Solution of Fully Non-Linear Ship
Wave Problems", SSPA Report No.
2912-3 (1986) .
15. Ni, S.Y., "A Higher Order Panel
Method for Double Model Linearized
Free Surface Potential Flows",
SSPA Report No. 2912-5 (1987) .
16. Ni, S.Y., "A Method for Calculat-
ing Non-Linear Free Surface
Potential Flows Using Higher Order
Panels", SSPA Report No. 2912-6
(1987).
17. Ni, S.Y., "Higher Order Panel
Methods for Potential Flows with
Linear or Non-Linear Free Surface
Boundary Conditions", Ph.D.
Thesis, Chalmers University of
Technology (19 87) .
2
Ni, S.Y., Kim, K.-J., Xia, F. &
Larsson, L., "A Higher Order Panel
Method for Calculating Free
Surface Potential Flows with
Linear Surface Boundary Condi-
tions", International Symposium on
Resistance and Powering Perform-
ance, Shanghai (1989).
19. Kim, K.-J., "Ship Flow Calcula-
tions and Resistance Minimiza-
tion", Ph.D. Thesis, Chalmers
University of Technology (1989).
20. Hess, J.L., "A Higher Order Panel
Method for Three-Dimensional
Potential Flow", Douglas Report N
62269-77-C-0437 (1979).
21. Chen, H.C. & Patel, V.C.,
"Calculation of Trailing Edge,
Stern and Wake Flows by a Time-
Marching Solution of the Partially
Parabolic Equations", Iowa Inst.
of Hydraulic Research, Report No.
285 (1985).
22. Patankar, S.V., "Numerical Heat
Transfer and Fluid Flows", McGraw
Hill (1980).
23. Chen, H.C. & Patel, V.C.,
"Practical Turbulence Models for
Complex Flows Including Separa-
tion", AIAA 19th Fluid Dynamics,
Plasma Dynamics and Lasers
Conference, Honolulu (1987).
24. Markatos, N.C., "The computation
of Thick Axisymmetric Boundary
Layers and Wakes around Bodies of
Revolution", Proc. Inst. Mech.
Engineers, Vol 198C (1984).
25. Larsson, L., "Boundary Layer of
Ships. Part III: An Experimental
Investigation of the Turbulent
Boundary Layer on a Ship Model",
SSPA, Goteborg, Report No. 46
(1974).
26. Dawson, C.W.: "A Practical
Computer Method for Solving Ship-
Wave Problems", Proceedings of the
Second International Conference on
Numerical Ship Hydrodynamics,
Berkeley (1977).
OCR for page 210
DISCUSSION
by C.M. Lee
This is one of the rare papers which cover
both inviscid and viscous flows about ships. I
believe that a rational approach to ship flow
problems should be based on an intelligent
combination of potential and viscous flow
solutions. This paper apparently shows the
efforts in that direction. The authors deserve
congratulation for their efforts. Bulbous bows
certainly can be cited as one of the useful
inventions achieved by the ship hydrodynamic
research community. However, we are still
uncertain of its effects on flow around ship
bows, particularly on the viscous flows. Does
a bulbous bow really contributes to reducing
the total resistance of a fat ship like a
super tanker at low speeds? I think it is
about the right time to examine the bow flow
for ships with bulbous bow in detail to find
the stagnation point on the bow as a function
of ship speed in order to develop a correct
model for computing ship flows more
scientifically.
Author's Reply
We agree with Prof. Lee that a thorough
investigation of the flow around the bulb of a
fat hull would be a most interesting exercise.
Ideally, such an investigation should include
an experimental and a numerical part. The
experiments should include detailed velocity
measurement, such as the ones by Fry and Kim,
presented at the 15th Symposium on Naval
Hydrodynamics using an LDV, and flow
visualization to find the stagnation point.
The numerical part should include free surface
as well as viscous flow calculations. For this
detailed investigation the best approach would
be to use a free surface Navier-Stokes method,
like the one developed by Prof. Miyata and his
co-workers. In such a case the wave breaking,
which is likely to be important in this case,
could also be considered.
DISCUSSION
by S. Ogiwara
There are several ways in prediction of
wave resistance by the panel method, namely,
the pressure integral over a hull surface,
computation of momentum change and so on. If
you have experiences to compare the results
between them, especially for the case of the
series models with protruding bulb, would you
comment about that?
Author's Reply
We have compared the two ways of computing
the resistance for several hulls. Results are
presented in reference [19] of the paper. For
the thin Wigley hull the differences between
the two methods are very small, while for the
Series 60, CB=0.60 hull, differences of about
10% were noted. The momentum approach yielded
smaller values than the higher order pressure
integrations. The Series 60 results are in our
experience typical for ships of moderate
thickness. However, for the very full HSVA
tanker we experienced the opposite trend with
the momentum results higher than the ones from
the pressure integration. The differences were
again around 10% .
DISCUSSION
by F. Stern
Fig.17 of your paper is very surprising
and contradictory to our work on viscous-
inviscid interaction, ie, our work indicates a
much greater extent of the interaction. Is
the inviscid-flow calculation for the bare or
displacement body? If for the former, please
comment on the lack of dependency in your
solution to the location of the outer boundary
even for outer-boundary locations near the
boundary-layer edge.
Author's Reply
It should be noted that the distance to
the outer edge of the grid is measured from
the keel, where the boundary layer is very
thin. Fig.A1 may be helpful. We claim in the
paper that there is an interaction at 2x/L =
0.05 but not at 2x/L = 0.11. As appears from
Fig.A1, the latter edge is about three times
the maximum boundary layer thickness away from
the hull. When reading the paper one might get
the impression that the 2x/L = 0.11 grid
extends only to the boundary layer edge, in
which case the results would have been most
questionable.
\
Fig.A1
210
11\ N:
\ | Boundory \//
\\\\ ~:
\ No os -
\~1 ----~_
\, ~
Representative terms from entire chapter:
boundary layer
