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 439
Nonlinear Ship Waves
Y.H. Kim (David Taylor Research Center, USA)
T. Lucas (University of North CarolinaCharlotte, USA)
ABsTRAcr
A boundary element method is presented for solving
a nonlinear free surface flow problem for a ship moving a
constant speed on a calm sea. The method of solution is
based on distribution of simple Rankine type~llr) sin
gularities on the body and on the true location of the free
surface, which must be obtained by iteration. The key is
the iterative algorithm for determining the free surface and
wave resistance using a new one parameter family of up
stream finite difference methods A verification of numeri
cal modeling is made using Wigley hull and validity of
computer program is examined carefully by comparing the
details of wave profiles and wavemaking resistance with
Series 60, CB =0.60 model.
=RoDucrIoN
Many of problems facing engineers involve such
difficulties as a nonlinear governing equation and nonlinear
boundary conditions at known or unknown bound
aries.The determination of surface ship waves is a non
linear problem wherein the governing equation is linear but
there is a nonlinear free surface condition imposed at an
unknown boundary. In practice the solution is approxi
mated using numerical techniques, analytical techniques
and combinations of both. Foremost among the analytical
techniques is the systematic method of perturbations in
terms of a small parameter.
In the past, problems of this type were generally
treated after the boundary condition on the free surface was
linearized. Limited studies on nonlinear boundary condi
tions have been conducted for simplified twodimensional
cases. With the recent advent of supercomputers and
numerous new techniques in computational fluid dynamics
there is growing interest in solving the threedimensional
nonlinear free surface problems by various schemes.
Dommermuth and YueL11 reported their study on the non
linear threedimensional gravity waves created by a
moving surface pressure distribution using a higherorder
spectral method. During the 5th International Conference
on Numerical Ship Hydrodynamic a significant amount of
work was been reported on this subject including
Jensen[2], Musker[3l, and Bai[41. SSPA's series of
reports by XiaL5l, Nit6], and KimE71 also deal with this
topic.
439
Many different numerical schemes have been pro
posed. The method proposed here is a boundary element
method with Rankine type source as the kernel function.
The drawback to this approach is, as well known, the
necessity of proper numerical implementation of the radi
ation condition required since the potential flow model
contains no viscosity. The upstream finite difference oper
ator first introduced by DawsonE8] successfully satisfies
the radiation condition by producing an inherent local
damping of the free surface waves by unwinding. Since
then numerous researchers have tried to improve
Dawson's method and to extend it to the nonlinear free
surface problem. SWIFT (Ship Wave Inviscid Flow
Theory) is a computer code developed by the authors and
S.H. KimE9] based on the linearized free surface con
dition. SWIFT uses a higher order panel method with
either constant or linearly varying singularity strength
distribution across a curved panel . The stability and accu
racy of the computational results have been examined and
found to be quite robust.
In the present paper the SWIFT approach is extended
to the nonlinear free surface problem and, though devel
oped independently, is very similar to Nit6] and Kiml7] of
SSPA. The convergence of the method has been major
focal point among researchers who have been working on
the nonlinear problems. XiaL5] experienced severe con
vergence problems during iteration. NiE6] used a a higher
order method and test results were converged if an under
relaxation coefficient was used to modify the wave eleva
tion and velocity distribution for the next iteration. KimE7]
obtained convergent solutions successfully by satisfying
the kinematic and dynamic boundary conditions simulta
neously without a relaxation factor. They all used the three
point finite difference operator instead of Dawson's four
point operator. Musker[3] extensively studied stability and
accuracy of a nonlinear code, and experienced the usual
difficulties in convergence with the four point scheme. In
this paper we introduce a new one parameter family of four
point schemes which can give an arbitrarily large
damping. The four point method of Dawsont8] and three
point method of Xia[5] are special cases of this approach.
We have experienced greatly enhanced stability for all
tested ship hullforms over a ~broad range of speeds.
Wigley mathematical hull is used to verify the numerical
modeling of the problem. The measured wave profiles and
wavemaking resistance of Series60, CB=0.60 are used to
validate the developed computer program.
OCR for page 439
MATHEMATICAL FORMULATION
We consider a ship which moves in the positive x
direction with constant forward speed UO in Figure 1. Let
Oxyz be a Cartesian coordinate system with Oz opposing
the direction of gravity and z=0 coincides the undisturbed
free surface, the positive xaxis in the direction of the
ship's forward velocity. The fluid is assumed to be
inviscid and incompressible and it's motion is irrotational
such that the velocity field of the fluid Vcan be defined as
v~x,y,z) = V¢(x,y,z) `1'
where ¢(x,y,z) is the velocity potential and satisfies the
Laplace equation
v2¢ = o (2)
in the fluid domain D and the boundary condition
An =Uon (3)
on the body surface SO where n=(nX,ny,nz) denotes the
outward unit normal vector on the boundary.
On the free surface z = ((x,y), where ~ is the free
surface elevation, the kinematic boundary condition and
dynamic boundary conditions can be given respectively as
~x~x+¢y~y¢z = 0 ~4'
go + ~ (V) ~ V0 Uo2) = 0
(5)
where g denotes the gravitational acceleration. Combining
both dynamic and kinematic conditions Eqs. (4) and (5)
becomes
V¢eV:2(V¢~23 + gig = 0
(6)
Finally, energy consideration requires that velocity poten
tial approaches the uniform onset flow potential and that
there be no waves far upstream of the ship, and that waves
always travel downstream.
The problem described in equations (1) thru (5) is
nonlinear, since the free surface boundary conditions(4)
and (5) themselves are nonlinear and should be satisfied on
the true free surface which is unknown and should be
obtained as a part of solution. In their exact form the
problem is difficult to solve. In order to be able to treat the
problem, the equations are approximated by ones which
are more tractable.
Following linearizations are similar to NiL6] and
Kim(~71. First we define
D~o,() = Add + tY(Y ~ (z
D2~,,() = go + 2(Vt V)  Uo2)
We denote ¢° and (° the velocity potential and the
wave profiles, respectively, obtained from the previous
step and introduce two small parameters 6o and b: with
respect to the previous solution. Take Taylor series
expansion of D~(c,,() and discard higher order terms
greater than Act and S(, then the kinematic free surface
condition becomes
D~o,() ~ D~°,(°) + ha D,(o,(°~6c; + a3` D~(c,°,(~
=40~0 +~0~0_40
+ 6¢xCx + bPy~y ~ Adz
+ ¢°~( + ¢°~` +~°xz(°x + (yz(Y ~ (zz)~( (7)
We redefine a new velocity potential ~ = 0° + ~ and the
linearized kinematic boundary condition can be expressed
as follows:
D~o,() ~ (x~x + (Y(Y  (z
+ {x6(x + (y6(y + (¢xz (x + tyz By ~ (zz )~( = 0 ~8'
Similarly the dynamic boundary condition can be linearized
as follows;
D2(c',() ~ 6` (1 + g(¢xtxz + (Y¢YZ + (zfzz)) + ~
 ~ {Uo2 + ¢° + ¢° + ¢°  2~°¢X + ¢°¢y + LIZ LIZ)} = 0
(9)
Now the new linearized boundary value problem is for
mulated. The singularity distribution and wave profile are
continuously updated until the solution converges and
satisfies original Eqs. (4) and (5~. Once the singularity
strength distributions are determined, the wave resistance
can be computed by
Rw = J~J. pnx4S
SO (10)
where SO is the wetted ship hull surface and the fluid
pressure p is given by the Bernoulli equation
p = PI ~v¢~2 u2~ pgz
NUMERICAL SCHEME
(1 1)
The method used here is a boundary element method
in which simple Rankine sources are distributed across
each panel. The boundaryvalue problem formulated above
then reduces to a determination of an unknown singularity
distribution over the boundary surface of the fluid domain.
Once the singularity distribution is determined, the hydro
dynamic quantities of interest, the velocity and the pres
sure, are also determined.
440
OCR for page 439
The present method uses a panel scheme to achieve a
numerical solution to the problem of Eqs. (1) thru (9~.
Panel schemes proceed by first dividing the boundary
surface into panels. Source distributions are assigned to
each panel. These distributions are expressed in terms of
unknown singularity parameters Hi associated with the
panel and neighboring panels. A finite set of control
points (equal in number to the number of singularity
parameters) is selected at which the boundary conditions
are imposed.
The construction of each network requires numerical
development in these areas: A. Surface geometry defi
nition; B. Singularity strength definition; and C. Control
point selection and boundary condition specification.
Essential features of the computational scheme in each of
these areas are as follows
A. Geometry input for a network is assumed to be
a grid of corner point coordinates partitioning the network
surface into panels. Panel surface is obtained by fitting a
paraboloid to corner points in an immediate neighborhood
by the method of least squares.
B. Discrete values of singularity strength at certain
points on each network are assigned as singularity pa
rameters. Singularity splines are constructed for each
network by fitting a linear distribution on each panel of the
network by the method of least squares.
C. Certain standard points on each network are
assigned as control points. These points include panel
center points as well as network edge points. Boundary
conditions involving the specification of potential or veloc
ity are applied at panel center points for the purpose of
controlling local properties of the flow.
Curved Panel Approximation
We use curved panels to patch the ship and free
surface. Let S be a true surface element bounded by four
corner points four corner points P~,P2,P3, and P4 as
shown in Fig. 2. We assume that S may be represented in
the approximate form
All, rig = (0 + (~ t + (~ rat + 2 (~ ~2 + (~ ~q + 2 (~2
where ~ (,q,() are orthogonal coordinates local to S. The
six coefficients ~ (0, Hi, All, Aid, (an, (~ ~ are obtained
by requiring that the approximate surface given by Eq.
(12) pass through its four corner points exactly and
through its neighboring points approximately in a least
square sense. That is, for each panel the expression
R = ~W,`(~,31~)  (,`)
2 ,` (13)
is minimized with respect to the six unknown coefficients
appearing in Eq. (12~. The summation in Eq. (13) ranges
over the N neighboring grid points. The choice of the N
points and the weight Wk has been made in an attempt to
.. . . ..
minimize any 1rregu arltles t net may appear in t ne
paraboloid approximation on the true surface. In the
present formulation N=16 and Wk=1 or 108 are used
depending on the location of the panel in the surface panel
network. Here we used the term network to represent the
collection of quadrilaterals.
The choice of these features is required by the
argument given by Hess[10] that certain computing
methods which use flat panels do not obtain increased
accuracy over the basic zero order method~constant source
on flat panels) even though they use higher order source
distributions. For purposes of computational efficiency
Eq. (12) can be reduced by an appropriate coordinate
transformation to the canonical form
~=at2+bll2 (t,11)~I; (14)
where ~ is the quadrilateral formed by the projection of
the corner points on the (mu) plane. In the case of flat
panel, a and b are both zero.
Singularity Strength Definition
The true singularity distribution will be approximated
by a truncated Taylor's series on each panel. Such a rep
resentation is valid on any interior part of the network
providing the paneling there is sufficiently fine. We will
consider a linearly varying singularity distribution on
source panels. There may be an advantages in using even
higher order distribution, but as pointed out by Hess[10],
it would be necessary to consider a higher order panel
geometry definition for the sake of consistency.
Specifically, we assume that the singularity strength ~ at
a point ~ A, 11, ~ ~ on a panel S is given by
o(t,11) = co + Cl`t + still (15)
The coefficients in Eq. (15) are determined by using a
method of weighted least squares over the panel and up to
eight of its neighboring panels. This method requires that
the form of Eq. (15) gives exactly the singularity value at
its centroid and approximate values in the least square
sense at centroids of the neighboring panels.
Induced Velocity Potential
The velocity potential at P=(x,y,z) induced by a
singularity distribution on S is given by
~ =  4 J ds (16)
where r = ~ (~x)2+~11y)2+~z)2 ~ I/2 and ~ is given in
Eq. (15~. We are now in a position to evaluate the integral
Eq. (16) using the relations of Eqs. (13), (14), and (15~.
The ~ can be expanded in closed form by using a com
bination of closed form calculations and recursive rela
tionships. Details are fully included in Johnsonfll] and
Kim et al.~91. Derivatives of velocity component along z
direction in Eqs. (8) and (9) are new and their final closed
forms are summarized in the Appendix.
METHOD OF SOLUTION
As usual we discretize the fluid boundary into a finite
number of panels. Because of symmetry, half of the ship
as well as the free surface domain is used as an input. A
typical free surface domain is truncated at a half ship
length(L=2.0) ahead of the ship and is extended more than
one ship length aft of the ship and a half ship length or
more to the side. Since the convective term in the free
surface condition involves x and y directional derivatives,
OCR for page 439
we use an algebraically generated waterlinefitted coordi
nate system, independent of the double model streamline
coordinate, to create the panel network. As shown in Fig.
3 the transverse lines are straight and vertical to the center
line of the ship. Our numerical experiments indicate that
the solution becomes more stable and accurate when the
aspect ratio of the free surface panel away from the ship
closes to unity. The same finding was discussed by
Musker(1989~.
During iteration, the panels are rearranged such that
the projection on the horizontal plane z=0 is unchanged,
i.e. x and y coordinates are fixed and the zcoordinate is
allowed to move up and down. The method of solution
requires that the linearized free surface conditions (8) and
(9) are satisfied for a given 0° and (° which were
obtained from the previous step. The iteration starts from
the solution which satisfies the wellknown linearized free
surface condition. In each iteration the free surface panels
as well as body surface panels are adjusted according to
the newly computed wavy surface and boundary con
ditions are applied to the updated boundaries. The itera
tion continues until the solution satisfies the exact free
surface conditions. In practice, we provide a certain
convergence criteria in the computation that requires the
iteration to continue until the residual error satisfies the
provided tolerance . This part will be discussed more
extensively with numerical results later.
One Parameter Family of Advection Methods to Enhance
Convergence
In this section we derive a one parameter family of
advection methods to enhance convergence for the non
linear free surface problem. A number of researchers have
been concerned with the question of convergence of
various forms of the nonlinear method. Xia had severe
convergence problems with his approach and recom
mended the three point finite difference operator instead of
Dawson's more accurate four point operator. Ni and Kim
have developed more powerful methods. Recently Musker
has completed a detailed stability and accuracy test of a
nonlinear code, and experienced the usual difficulties in
convergence with the four point method. Musker used two
approaches to enhance convergence: he used a spline
scheme involving four upwind points and he elevated the
free surface panels. His nonlinear iterations, unlike those
of Ni and Kim, are all on the flat surface.
In this paper we introduce a new one parameter
family of four point advection schemes which can give an
arbitrarily large damping. The four point method of
Dawson and the three point method are special cases of
this approach, and without the use of elevation of the
panels (which we have not examined) we have experienced
gready enhanced stability. We recommend the inclusion of
this technique into existing codes. For most programs this
would require only a minor modification of a few lines of
code, as will be seen below.
The four point method is:
Slopei = ( A4fi + B4fi ~ + C4fi2 + D4fi3 ~ (17)
where for an uniform mesh spacing of size 1, A4 = 1.6671
1, B4 = 2.511, C4 = 11, D4 = 0.16711. The related
three point scheme is of the form
Slopei = ( A3fi + B3fil + C3fi2 ~ (18)
where for an uniform mesh spacing of size 1, A3 = 1.5ll,
B3 = 2.011, C3 = 0.511. In considering all possible four
point schemes for approximating f ' one might like to
require that they be exact for quadratics. These three
restrictions on a formula of the type
Slopei = ( ANt + BNfil + CNfi2 + DN03 ~ (19)
leave only one degree of freedom and since both the three
point and four point methods satisfy a quadratic exactly,
the general derivation of this general fam~ly is very simple:
AN = ~ 1  Qmul ~ A4 + Qmul A3
BN = ~ 1  Qmul ~ B4 + QmulB3
CN = ~ 1  Qmul ~ C4 + Qmul C3
DN= (1Qmul ~D4
(20)
where Qmu~ is a parameter to be chosen greater than or
equal to zero. The case Qmu~ = 0 gives the familiar four
point formula while Qmu~ = 1 returns the more stable
three point formula. We have found that Qmul = 1 is often
fine for lower speeds but for higher speeds Qmul = 3 (or
even sometimes Qmul= 5 ~ eliminates most of the
problems with convergence for reasonably regular com
putational grids and with little loss of accuracy.
To analyze the effect of the one parameter method
and to get some feel for the meaning of the parameter
Qmul, it is instructive to examine the coefficient of fi"'
which is the leading term in the Taylor expansion of the
error for both the three point and four point formula.
(Recall that the basis of the four point formula was not that
the fi"' term was eliminated but rather that the fi""
term was required to have a zero coefficient.) Here we
will take an uniform spacing of size 1. The coefficient of
fi"' in the three point method is of the magnitude 12/3
while in the four point method it is 12/6, leading again to
the well known conclusion that the superior damping of
the three point method will frequently enhance conver
gence over that of the four point method, an important
finding of Xia for example. Thus the one parameter Qmu~
method has a coefficient for fi"' of the magnitude
12( (1  QN,UI)/6 + Q~ /3 ~ = 12( 1 + Ql~lUI)/6
giving again 12/6 and 12/3 for Qmu~ = 0 or 1 as before.
When denser meshes are used, or for large speeds where
f"' may be smaller, the us~e of larger values of Qmul such
as Q~ul = 3 or 5 gives the (larger) coefficients of fi"' of
magnitude 2/312 or 12 giving in these cases twice or triple
the damping of the 3 point method.
We now examine briefly the spline scheme consid
ered by Musker, which for an uniform mesh has the formf
refer to (201) with A4 = 1.5551l, B4 = 2.17711, C4=
0.6891l, D4 = 0.06671~. This method does not quite
make the quadratic term zero, so is not included in the
above family, but it is very close with a value of Qmul of
about 0.6. Thus it has much of the improved stability
associated with the three point scheme, giving a coeffi
cient of f" of the magnitude 0.2712, closer to the 0.3312
of the 3 point method than the 0.16612 of the 4 point
method. Thus we would be led to conjecture that the 3
point scheme would be more stable than the spline scheme,
at negligible loss of accuracy.
RESULT AND DISCUSSION
To facilitate comparisons, the following nondimen
sionalization has been made:
442
OCR for page 439
Cw~wavemaking resistance)  RWI(1I2 pU2S)
((wave profile) = ~ /(L/2)
= ~(U2/~2g)) (used in Figures)
Although the higher order panel method gives better
results, we used linearly varying source strength distri
bution across curved panel for the ship hull and constant
source strength across flat panel for free surface during
computation because of the following reasons:
(1) The higher order panel method requires a much
longer computing time, normally 3 or 4 times longer with
the same number of panels.
(2) The finite difference scheme used for the
convective terms in the free surface condition requires such
small panel sizes that the advantages of higher order panel
scheme cannot be fully utilized over the free surface
domain.
Numerical experiments have indicated that the differences
in results are negligible between using higher order panel
and flat panel on the free surface domain. On the other
hand, a Neumann condition imposed on the hull surface
does not require any numerical differentiation and the
higher order panel method gives better results than flat
panel result with even fewer panels. Furthermore, the
wavemaking resistance is small quantity and to get this
quantity accurately the integration of the pressure distri
bution over the entire wetted hull surface has to be carried
out precisely. It is therefore required that the pressure
distribution as well as panel surface be precisely approx
~mated.
Qmul Test with Wigley Hull
Wigley parabolic hull was used to verify the pro
posed one parameter family of advection methods. The
following computational conditions were used:
Wigley model L=2.0,L/B=10.0,D/B=0.625
Bow = +1.0, Stern =1.0
Computation domain x= 2.0, +2.0
(S SPA recommended) y= 0.75
Number of panels 132 on ship (23 Stax7VVls)
400 on free surface (41xl 1)
Froude numbers 0.316, 0.408
Qmul's 0, 1,3,5
Computer Cray XMP/24
We used the computational domain that XiaL5] and other
SSPA reports recommended and also matched the panels
on the ship hull and free surface. We selected higher
Froude numbers because several researchers have experi
enced convergence problems particularly at higher speeds.
Table 1 shows the computed wave making resistance using
four different Qmul's: 0, 1, 3, and 5. The solutions
diverged when Qmul = 0 and 1, and converged when
Qmul=3 and 5. The experimental values of wavemaking
resistance coefficients were reported at the 17th I1lCE12]
and Cw ranges between 1.7x103 and 2.3x103 at
Fn=0.408. Qmul=0 is the four point scheme and the solu
tion diverged after 2 iterations. At Qmul=l, the three point
scheme, the solution diverged after 5 iterations. Qmul=3
required only 4 iterations and Qmul=5 required just 2
iterations to get a converged solution. The magnitude of
the wave making resistance is slightly smaller than that of
Qmul=3. It is likely that the larger damping by using
Qmul=5 caused the difference.
Table 1 Qmul Test
(Wigley Hull at EN = 0.408)
Cw(Experiment) = 2.0
linear
2.248
2.304
2.317
2.231
. 1 1st
 Iter
1.562
1.890
2.009
_ 1.965
2nd
Iter
1.863
2.067
1.976
_ 1.939
3rd
Iter
Divrg
1.683
1.982
1.939
.
4th
Iter
1.726
1.978
1.940
5th
Iter
1.899
1.978
1.940
Cwx103
Non
linear
_
Divrg
Convg
Convg
The two small parameters 8o and 6( used when
linearizing the free surface conditions were continuously
monitored during the computation. In Table 2 EITor, the
summation of the square root of bt squared terms, is tabu
lated at each iteration for different Qmul's and in Table 3,
that of 8~. At Qmull, Error oscillates during the itera
tions until the solution diverges. With Qmul=3 or 5, the
magnitude of Error decreases steadily as the iteration
continues. In practice, the convergence rate after several
iterations is rather slow so we assumed the solution to be
converged when the ratio of the current Error to the first
iteration Error is smaller than 1% for both Act and by.
Table 2 ~ Convergence Test
(Wigley Hull at EN = 0.408)
N
errorl= ~)2,
i
Qmul
o
3
5
1st
Iter
8.248
4.199
.
2.955
l 2.859
2nd
Ite r
4.550
.
4.493
.
1.857
0.695
3rd
Iter
152.5
5.488
0.507
0.298
l
N = 400
4th
Iter
Divrg
2.683
0.333
0.033
error1 x1 n2
5th
Iter
4.975
0.036
._ 1
Non I
linear I
1
1
Divrg 
Convg 
3
Table 3 ~ Convergence Test
(Wigley Hull at EN  0.408)
N
error2= ~)2,
J
~ 1
l
..4
443
1st 1
Iter 
1.10 1
1
0.62
0.43
3
2nd
Iter
1.30
1
1.1 0
o~o9 1
~3
3rd
Iter
Diverg
1 .00
0.062
0.012
N = 400
1
1 4th
I Iter
1
1 0.58
I 0.017
1
l
1 5th
I Iter
1
1 0.84
0.009
0.001
1
 Non
 linear
 Diverg
 Convrg
 Convrg
OCR for page 439
Fig. 4 depicts the wave profiles at the hull surface
during iterations including the wave profile obtained from
the linearized free surface conditions. As the iteration
progresses the bow wave grows and the first trough
becomes deeper with the phase shifting slightly towards
the bow. Only insignificant differences are observed aft
during the iterations.Figure 5 shows the wave profile
comparisons between the experimental, linear and non
linear computations. The nonlinear solution is now much
closer to the experimental.
Series60, CB=0.60
The Series60 with block coefficient 0.60 has ex
tensive model experimental results, including wave
profiles and wavepattern resistance. Experimental values
used here were determined by Kim and Jenkins[131.
Experimental data only for model fixed were considered in
comparisons. Wave profiles were marked at every station
along hull with a grease pencil during the run and were
read after each run. Wave resistances were obtained from
longitudinal wave profile measurements. The following
computational conditions are used:
Series60 model L=2.0, L/B=7.5, D/B=0.4
Model is fixed
Bow = +1.O, Stern= 1.0
Computation domain x = (3.5, 2.0)
y = (1.5, 0.0)
Number of panels192 on ship(25 Sta x 9WLs)
696 on free surface(59 x 13)
Froude numbers0.22,0.25,0.28,0.30,0.32,0.35
ComputerCray XMP/24
Figure 6 shows the comparison of wave profiles along the
hull between computed and measured. The wave profiles
obtained from linearized free surface condition are also
plotted together. The nonlinear solution significantly
improves the wave profiles, in particular, at the bow and
the first trough and after that the difference between linear
and nonlinear results seems insignificant. Qmul = 1 is used
for Fn = 0.22, 0.25, and 0.28, and Qmul = 3 for Fn =
0.30, 0.32, 0.35.
Table 4. Wave Resistance Coefficient Comparison
Series60, CB = 0.60 (Model Fixed)
Fn~
0.22
0.25
0.28
0.30
0.32
0.35
linear
0.369
0.664
1.581
2.144
1.815
2.148
1 st
Iter
0.278
.
. 0.398
1.123
1.548
1.401
. ~
1.886
2nd
Iter
0.271
0.357
1.049
1.474
1.333
1.800
3rd
Iter
0.268
0.354
1.043
1.467
1.350
1.816
4th
Iter
0.266
0.355
1.043
1.460
1.354
1.828
(CWX1 03)
.
Nonlin
ear
0.255
0.354
1.045
1.460
1.357
1.836
Table 4. summarizes the wavemaking resistance
computation results at each iteration and in Figure 7 the
linear and nonlinear results are compared with the
experimental data. It is interesting to see that linear solution
consistently predicts higher wavemaking resistance
throughout the Froude number range. Nonlinear results
show very close agreement with the measured data. Even
though all the runs converged successfully, at the lower
Froude numbers 0.22 and 0.25 agreement on wave
profiles as well as wavemaking resistance were rather
poor. It is likely that the panel size may be too large at
lower Froude number where wave length is shorter and is
inadequate to resolve the flow phenomena accurately. At
FN=0.22 the estimated characteristic wavelength is 0.608
which is less than 1/3 of the ship length. In usual practice
between 10 and 12 points per wave length are considered
to get acceptable numerical accuracy for wave making
resistance. This means at least 33 evenly distributed panels
along the ship hull have to be provided. No further attempt
was made to increase the panel density at the lower Froude
numbers, because it exceeds Cray XMP memory
restriction 3 million words. This particular run requires 2.6
million words computer memory~maximum usable
capacity is 3 million+) and approximately 70 CPU seconds
were used at each iteration.
As seen in Table 4. the most sizable improvement
was made at the first iteration and in subsequent iterations
Improvement was small. In practice we may need one
more iteration over the linear solution to have reasonable
results. Convergence of the solution is always monitored
carefully during computation. Table 5 shows one exam
ple.
Table 5. Convergent Test
series60, CB=0.60 at Fn = 0.32
error1
0.9147
0.1783
0.0513
0.0219
.
0.0153
error2
0.670
0.170
0.01 7
0.037
0.013
SOUrCe
2.9383
3.0646
3.1 1 57
3.1411
3.1 52 1
3.1521
_ ~
1 .8 1 5
1 .401
1 .333
1 .350
1 .354
1 .357
Here errorl and error2 are defined in Tables 2 and 3 and
NF
NF is the total number of free surface control points. For
example, if we average errorl by NF=696, average wave
height increment ~hl at each free surface panel after the
first iteration is 0.001314 and after the fifth iteration 6~5
becomes 0.000022 which are nondimensionalized by half
of the ship length. If Series60 hull length is 200m,
6~1=131.4mm becomes 6~5=2.2mm after five iterations.
Fig. 8 shows two pairs of the contour plots using the
locally developed visualization program WAVEt143. The
Froude number is 0.32. The nonlinear solution was
obtained after 6 iterations. The first pair shows contours of
the wave height and source strength distribution on the free
444
OCR for page 439
surface for the linear step, and the second pa* for the
nonlinear solution. The wave height is plotted for the
positive yplane and the source strength for the negative y
plane for convenience. 0.03 uniform spacing is used
between contour lines for both plots. Solid lines indicate
positive wave height and positive source strength, dotted
lines negative wave height and negative source strength,
and dashed lines zero. We observe the following: (1) the
source strengths on the body toward the bow of the
ship~not plotted here) are large and positive, but at the free
surface near the bow a locally concentrated sink dism
bution is noticed, with a strong and wide source distribu
tion following immediately downstream,, (2) at the stern
there is only sink distribution on the ship surface and the
free surface, (3) the humps and hollows of wave height are
out of phase to those of source strength. Zeroes of wave
height occurs near the mean of either source or sink distri
bution along xaxis.
CONCLUSIONS
1. The introduction of a one parameter family of
advection methods has significantly enhanced the con
vergence. This new scheme can control the magnitude of
damping in the numerical solution, and proper choice of
~ ] may eliminate the numerical instabilities that several
researchers have experienced when using either the three
point or the four point method.
2. Contour plots similar to Fig.8 for example, are
strongly recommended when studying such a complex
numerical model. Using the contour plots we can more
easily visualize the quality of computational results and
choices of truncation regions and panel densities which
often may be misjudged just looking the local flow charac
teristics such as the wave profile or streamline tracing
along the ship hull or especially a scaler quantity such as
wave resistance alone. We have found that numerical
instability sometimes starts at the edge of a boundary with
negligibly small magnitude, and it progressively grows
and spreads as iteration continues.
3. The most significant improvement was made
between the linear solution and the first iteration. It is
suggested that at least for the first iteration an under
relaxation coefficient either obtained from Eq. (22) or
provided with small value such as 0.5 used by Nit6] be
used. For practical purposes, 2 to 3 iterations may be
sufficient.
4. It now apples that the long search for programs
to accurately estimate wave resistance is now drawing to a
successful conclusion, and the next emphasis should be on
the more difficult problems of simulating the wave height
on the hull and more generally the wave heights on the free
surface. Here the comparison with experimental results
combined with flow visualization methods is essential.
ACKNOWLEDGMENT
The first author was funded by the Applied
Hydromechanics Research(AHR) Program supported by
the Office of Naval Research and administered by David
Taylor Research Center(DTRC).The second author was
supported in part by an IPA grant from DTRC and super
computer grants from the National Science Foundation
(Grant Number ECS 8515174) and the North Carolina
Supercomputer Center. The authors express their gratitude
to Dr. WenChin Lin and Dr. Michael Wilson of DTRC for
their valuable suggestions and encouragement during the
course of this work. The authors also thank Mr. Peter
Chang and Mr. Steve Fisher for their kind support.
REFERENCES
1. Dommermuth, D.G. and D.K.P. Yue, "The Nonlinear
ThreeDimensional Waves Generated by a Moving
Surface Disturbance," Proceedings of the 1 7th
Symposium on Naval Hydrodynamics, the Hague,
Netherlands, 1988
2. Jensen, G., V. Bertram, and H. Soding, "Ship Wave
Resistance Computation," Proceedings of Fifth
International Conference on Numerical Ship
Hydrodynamics. Hiroshima, Japan, Sept., 1989.
3. Musker, A.J., "Stability and Accuracy of a Non
Linear Model for the Wave Resistance Problem,"
Proceedings of Fifth International Conference on
Numerical Ship Hydrodynamics. Hiroshima, Japan,
Sept., 1989.
4. Bai, K.J., J.W. Kim, and Y.H. Kim, "Numerical
Computations for a Nonlinear Free Surface Flow
Problem," Proceedings of Fifth International Conference
on Numerical Ship Hydrodynamics. Hiroshima, Japan,
Sept., 1989.
5. Xia, F., "Numerical Calculations of Ship Flows, with
Special Emphasis on the Free surface Potential Flow,"
Chalmers University of Technology, Goteborg, Sweden,
1986.
6. Ni, S.Y., "Higher Order Panel Methods for Potential
Flows with Linear or Nonlinear Free Surface Boundary
Condition," Chalmers University of Technology,
Goteborg, 1987.
7. Kim, K.J., "Ship Flow Calculations and Resistance
Minimization," Chalmers University of Technology,
Goteborg, Seden, 1989.
8. Dawson, C.W., "A Practical Computer Method for
Solving Ship
Wave Problems," Proceedings of the 2nd International
Conference on Numerical Ship Hydrodvnamics,
Berkeley, Calif., 1977.
9. Kim, Y.H., S.H. Kim, and T.R. Lucas, "Advanced
Panel Method for Ship Wave Inviscid Flow Theory,"
DTRC89/029, 1989.
10. Hess, J.L., "Consistent Velocity and Potential
Expansions for Higher Order Surface Singularity
Method," MDCJ 6911, pp.l30, June, 1975.
11. Johnson, P.T., "A General Panel Method for the
Analysis and Design of Arbitrary Configuration in
Incompressible Flow," NASA Contract Rpt 3079,
Boeing Com. Airplane Co., Seattle, 1980.
12. Norrbin, N.H.(editor), "The Proceedings of the 17th
International Towing Tank Conference," SSPA,
Goteborg, 1984
445
OCR for page 439
13. Kim, Y.H. and D. Jenkins, "Trim and Sinkage Here
Effects on Wave Resistance with Series60, CB=0.60,"
DTNSRDC/SPD101301, 1981.
14. Robinson, C.R., "Ship Wave Color Graphics Using
DI3000." Computer Science Senior Project, U. of North
Carolina at Charlotte, 1986.
APPENDIX
Denvatives of Velocity Components in Zdirection
¢,~z =  (t )~( ~ sec(( · n~d~d
(yz = IJ (ll Y)~( ~ sec(( · n~d~d
(21)
(22)
(yZ= (~;s ) sec((·n~d~dq~sec((·n~d~dll
where cs((,ll) is given Eq.~15), r= (~E,x)2+(lly)2+~
zy2~ 1/2, and ~ shown in Fig. 2.
We assume that surface panel is sufficiently fine that the
term sec(( · n) becomes zero. After algebraic manipu
lation, the following expressions are obtained:
0~z = ooK,~,~(l,l) + o~,~xK~,,~(l,l) + K
os
o~
c~
o.o
o~
0.1
n.e
t
0.2
~0.8 0.6 0.4 0.2 0.0 0~ 0.4 0.6 0.8
z/(L/2)
Fig. 4 . Wave Prof iles at Hull Surface
Wigley Hull at Fn=O. 316
1.0  .8 0 6 0.4 o~ o.o o~
z/(L/2)
Fig. 5. Wave Profile Comparisons
Wigley Hull at Fn=O. 32
co
o
X3
_
~ Linear Theory
O Nonlinear Theory _
1=~
T
.~
\ Lo
\ o
\ ~
_ _
~6 
at Hull Surf ace
a5
c~
2.5'
0.5
Fig. 6. Wave Making Resistance Comparisons
Series60, CB=0 . 60
.
Linear
· Theory /
Nonlinear /
· Theo ~/
~ Experiment 77
~ / ~
1 ~ /
. ~
_ ~
1
L\ '/
C  Rw
_ w 2pU2S
n~
O.l
nn
n l
0.2
_n ~
O1 0.1
_ . . , . . . . . . ·
0.2 0.25 0.3 0.35 ,,~
Fn= U/(gL)1/2
0.3
447
. ~
0.4 0.2 0.0 0.2 0.4
x/(L/2)
Fn = 0.22
1 1 1 1
~r I I T
~:
. = 1 1 1
 7` Linear Theory L
I O Nonlinear lheory I
I Experiment 1rith Nodel Fixed I
L I
0.2 0.0 0.2
x/(L/2)
Fn = 0.28
I ~linear lheory L
I O Nonlinear 1heory 1
I Experiment 1rith Yodel Fixed
 1.0  .8 0.6 0.4 0.2 0.0 0.2
x/(L/2)
Fn = 0.30
1.0
0.(
O.
1.0
F7FI
~1 ~'o'll 1
1 1 11
Ll 11
0.6 0.8 1.0
OCR for page 439
0.3
0.2
0.1
Cal 0.0
0.1
0.2
0.3
e
_
 1.0  .8 0.6 0.4 0.2 0.0 0.
x/(L/2)
Fn = 0.32
 ~Linear Theory _
O Nonlinear Theory
E3rper~ment Tenth Yodel Fixed
~ V ~
0.6 0.8 1.0  1.0 0.8
Figure 7. Wave Profile Comparisons at Hull Furface
Series 60, Cg = 0.60
Aim,\ ( \~ \ \ five U\ight
_ _
~linear Theory
O Nonlinear Theory _ _
I Experiment with Nodel Flax _
/.
0.6  .4 0~ 0.0 0.2 0.4 0.6
x/(L/2)
Fn = 0~35
my',
\~
. _ .
'> 1 1 1 1 1 1 T 7 We___ 7 1
Aft//; //~ /Z~/L~Jsource Strength /
LINEAR SOLUTION
`~W \N \\ \ C\ ~\ Wa~eight In,
FIR
~ 111
_~ 29
. __ _
__
0.8 1.0
~ 1
3 . 0
art /; Yew' ale / source Strength '/~
NONLINEAR SOLUTION
Figure 8. Contour Plots of Wave Height and Source Strength
Series60, CB=0.60 at Fn=0.32
448
. ~
OCR for page 439
DISCUSSION (SESSION V, PAPER 5)
John V. Wehausen
University of California at Berkeley, USA
The agreement between the authors' ~nonlinear. computations and
measurement is impressive, but I wonder is one should not be
suspicious of such good agreement. Since the computations are based
on the irrotational flow of an inviscid fluid, shouldn't a measured
wave resistance based upon residuary resistance, or even wavepattern
resistance, show at least some discrepancy with the computed
resistance? The wave resistance is a rather sensitive test, for it is the
difference between two large numbers and one expects to see some
effect of viscosity in the stern region. I shall be curious to see the
results for similar computations for Series 60, CB = 0.80. I also
have a question concerning the numerical procedure. ~Convergence,.
or lack of it, appears to depend upon the choice of Qua, but Q,~,, =
3 or 5 both lead to convergence, but to different values. Which
should I believe? and why? Even if they were identical, can one
give an estimate of the difference between the obtained value and the
solution to the exact problem (when 1 _ 0 in some fashion)?
DISCUSSION
Ronald W. Young
University of California at Berkeley, USA
Your results show remarkably good agreement in terms of wave
resistance and wave/profile when compared with experimental data.
Jensen and Soding (1988, Jahrbuch der Schiffbautechnische
Gesellschaft) have also obtained similarly good agreement for the
same formulation but using discrete sources. They also have very
good conveyance in their iterative scheme. I was wondering if it will
be possible for you to compare your results with theirs. To my
recollection, their stern wave profiles do not have as good an
agreement with experiments as yours. Of course, for an identical
model, and identical mathematical problem, there should be only one
set of correct results. For a blunt bow ship, Jensen & Soding had
difficulty obtaining convergent and accurate results in the bow area.
You may like to try a similar test with your code. I want to
congratulate you and your colleagues, including Dr. S. H. Kim, for
the successful development of this code.
DISCUSSION
Hauime Marno
University of California, Santa Barbara, USA
The examples, which are employed in this paper, are Wigley hull and
Series 60 hull. In these cases the effect of nonlinearity is
comparatively small, and the advantage of the nonlinear computation
may not be well demonstrated. The nonlinear effect is remarkable in
the case of hull forms as was employed in the paper by Maruo &
Oginara (1985), for which the nonlinear computation is indispensable.
DISCUSSION
Hoyle Raven
Maritime Research Institute Netherlands, The Netherlands
Your method shows a very good convergence for the cases studied.
However, this convergence depends on the numerical (artificial)
damping caused by the difference scheme. This damping does not
vanish upon convergence but remains present in the result. For the
rather mild cases shown, its influence does not seem to be large, but
for full hull forms the required value of AN could be much higher
and the damping could be more detrimental. Similarly, a reduction
of the panel size is likely to ask for higher Q i. Therefore, I believe
that other inherently more robust formulations of the iterative
procedure remain something to be searched for. Could you comment
on this?
DISCUSSION
Kuniharu Nakatake
Kyushu University, Japan
I congratulate you on your fine results. You tested a new family of
upstream finite difference schemes in order to satisfy the radiation
condition. As well known, the calculated results vary with choice of
different scheme. I think it is not reasonable. We are using a phase
shift method which we call KU (Kyushu Univ) method.' The feature
of this method is nonuse of difference scheme. According to our
experience, it is applicable to very narrow computation domain. And
a theoretical proof for it is to be presented by Dr. H. Seto
(Mitsubishi Heavy Ind.) in Nov. 1990. I recommend you to test KU
method too.
[1] J. Ando and K. Nakatake, Tran. of WestJapan Society of Naval
Architects, Vol. 75, 1987.
DISCUSSION
Dimitris Nakos
Massachusetts Institute of Technology, USA
The design of iteration schemes, based on linearization about the
previous solution step, is the most critical issue in numerical solutions
of nonlinear steady ship waves. Such schemes are desired to be
convergent within a few iterations, due to the large computational
effort involved at each solution step. It is fair to state that the part
of the wave flow which is expected to resist convergence is associated
mostly with the diverging wave system. The authors have
demonstrated effectively that excessive numerical damping is able to
Filter out. all short diverging waves and consequently accelerate the
convergence (see e.g. Figure 8). The question that arises is the
following: Is it proper to disregard essential features of the flow in
the name of numerical convergence, in particular, in light of the fact
that the relatively long transverse waves are mostly unaffected by
nonlinearities?
AUTHORS' REPLY
We would like to express our thanks to all of the participants who
have shown their interest in our paper, made comments, and
expressed their critiques. Several questions have a common interest
and we have grouped them without advance consent.
Reply to Professor Wehausen:
The experimental values of the wavemaking resistance we used to
make comparisons were obtained by the longitudinal wave profile
measurement suggested by Eggers[l]. In his derivation the effect of
the viscosity is neglected and the corresponding computer code has
been reported by Reed[l]. Table A summarizes the linear, nonlinear
and measured wavemaking resistance for Series 60, CB = 0.60. At
the lower Froude number, the discrepancy between the measurement
and the nonlinear solution is significant and even at the higher Froude
number the difference is still noticeable.
449
OCR for page 439
Reply to Professors Wehausen, Marno, Yeung, and Dr. Raven:
They have shown common interest in the nonlinear effects on a ship
with blunt bow and large block coefficient. We have selected Series
60, CB = 0.80 which has a half entrance angle of 43 degrees as an
additional test case. Figures A and B show the linear and nonlinear
wave profiles together with a series of wave profiles obtained during
each iteration for Fn = 0.20 and 0.25 with Q,l,'., = 3 and 5
respectively. No comparisons were made with measurements because
of lack of availability. As Prof. Marno pointed out, the nonlinear
effects on Series 60, CB = 0.80 is more significant than those on
Wigley hull or Series 60, CB = 0.60. One noticeable difference
observed between ships during the iteration is that for Series 60, CB
= 0. 80 the wave profile after the first iteration has the largest hump
and hollow near bow region and from the second iteration the wave
profiles oscillate within the linear and the first wave profile band
until the solution converged. On the other hand, for Wigley hull or
Series 60, CB =0.60, the wave profiles always progressively grow as
the iteration number increases until the solution converged. For Fn
= 0.30, we increased Qua from 5 to 7 to 9, but failed to get a
converged solution. As the iterations continued, the wave amplitude
grew continuously and after the third iteration the wave amplitude
near the bow region exceeded the ship draft and the solution
diverged. Without having measurements, it is difficult to evaluate the
quality of our converged solutions for a fullform ship. Whether the
Fn = 0.30 case is beyond the capability of our code would require
further study, but as the A,,,, method is very new, and has now been
shown to converge nicely for lower Froude numbers even for this
fullform ship, with further research we would anticipate convergence
even at Fn = 0.30.
Reply to Prof. Young:
Dr. Jensen compares his wavemaking resistance computations with
Ogiwara's experimental results which were conducted at free sinkage
and trim conditions. Our code can be extended to handle the free
sinkage and trim conditions. There are no plots on the wave profile
in Jensen's paper.
Reply to Prof. Nakatake:
The phase shift and the shift of collocation point may be different
terminology for the same approach. A shifting also introduces
numerical damping and is algorithmically similar to the use of
upstream finite difference formulae. We had earlier implemented the
shift of the collocation method in our computer code and tested it.
In our numerical experiments, we shifted the collocation points
forward between 10% and 50% of the characteristic length of each
panel. Our numerical results showed that the shifting apparently
removed the oscillation of wave profiles, in particular near the bow
region, but decreased the wave amplitude proportionally as we
increased the shifting. According to our numerical experience, we
haven't found significant advantages to applying the shifting logic to
our code.
Reply to Prof. Wehausen, Dr. Nakos, and Dr. Raven:
It is well known that the panel method using a Rankine type
singularity with an upstream finite difference scheme to satisfy the
free surface and radiation conditions provides good results near the
ship [4]. But because of an inherent local damping introduced by
unwinding, the wave amplitudes attenuate significantly as they
propagate and the results get poorer far downstream. Our aim was
to develop a quality computer program which is practical for ship
designers in the shipbuilding business. The linear codes of the past
made many simplifications. We have provided an approach to solve
the original nonlinear problem that eliminates these, and gives results
that compare well with those of experiment, with a modest use of
computer time. We haven't examined the relations between panel
size, Proude numbers, and A s carefully at this point, however we
will be considering these questions in the future. Your comments and
critiques are very much appreciated.
References:
[1] Eggers, K.W.H., cuber die Ermittlung des Wellenwilderstandes
eines Schiffsmodells durch Analyse seines Wellensystem,. I.II
Schiffstechnik 9, pp. 79 84, 1962; Disk 85; 10, pp. 93106, 1963.
[2] Reed, A.M., "Documentation for a Series of Computer Programs
for Analyzing Longitudinal Wave Cuts and Designing Bow Bulbs,.
DTNSRDC/SPD082001, June, 1979.
[3] Sclavanous, P.D. and D.E. Nakos, Stability Analysis of Panel
Methods for FreeSurface Flows with Forward Speed, Proc. 17th
Symp. on Naval Hydrodynamics, The Hague, The Netherlands, 1988.
[4] Lindenmuth, W.T., T.J. Ratcliffe, and A.M. Reed, Comparative
Accuracy of Numerical Kelvin Wake Code Prediction'WakeOff',"
DTRC/SHD1260101, May, 1988.
450
OCR for page 439
Table A WaveMaking Resistance Comparisons
Series60, Cg=0.60(Model Fixed)
En i Linear I Nonlinear
0.22 0.369 0.255
.
0.25 0.664 0.354
0.28 1.581 1.045
0.30 2.144 1.460
0.32 1.815 ~1.357
0.35 2.148 1.836
(CWX103)
Experiment
0.176
0.230
1.011
1.375
1.316
1.780
Fig.A. Series60, CB=O.BO at Fn=0.200
0.5
1.0
__ _
Symbol SWIFT
~ Linear Have
_ I E}  Iteration No.1
O Iteration No.2
O Iteration No.3
I ~ Iteration No.4
_   ~   Nonlinear
l
~ ~ = 11 .
_N = 1 1 ~
_= === =.
.
Q ~
===~N
~ =~
_~ fit _
==~= _
=__ _
_ _
, ~
_
 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4 0.6 0.8 1.0
x/(L/2)
451
OCR for page 439
Fig.B. Series60, CB=O.8O at Fn=0.250
.o
0.5

cot
~ o.o
_'
._
o
at
0.5
1.0
_ I S~nbolSVIFr
l ~Linear Have
_ l ~Iteration No. 1
O Iteration No.2
O Iteration No.3
l ~Iteration No.4
_ l    ~   Nonlinear
. .
14~:
~ W~+ _0
_
_
=== ==
I I I AL
===~
r
~
=== 
 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4
x/L/2)
452
0.6 0.8 1.0