| 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 192
24th Symposium on Naval Hydrodynamics
Fukuoka, Japan, 8-13 July 2002
Unstructured Nonlinear Free Surface Simulations for the
Fully-Appended DTMB Model 5415 Series Hull Including
Rotating Propulsors
Clarence 0. E. Burg, Kidambi Sreenivas, Daniel G. Hyams, Brent Mitchell
(Computational Simulation and Design Center, Mississippi State, MS)
ABSTRACT
Nonlinear free surface simulations around realistic ge-
ometries, such as the DTMB Model 5415 Series hull,
are a necessary step to achieve the goal of simulation of
maneuvering surface vessels. As part of this effort, re-
searchers in the Computational Simulation and Design
Center at Mississippi State University have developed
and parallelized a three-dimensional unstructured code,
U2NCLE, which solves the incompressible Reynolds-
averaged Navier-Stokes equations using a node-based
flux differencing finite volume method with Roe-averaged
variables. The simulations produced by this code real-
istically capture turbulent flow and vortices arising from
bulbous bows and tips of propulsors and rudders and can
model the rotation of the propulsors accurately. Because
the code uses unstructured grids, the code is applicable to
a wide range of geometries, including those with shafts
and struts, and can solve the flow through small gaps such
as those between the rudder and the boot. Current re-
search efforts include the addition of nonlinear free sur-
face capabilities to this code. A nonlinear free surface is
obtained by solving the kinematic free surface boundary
condition, imposing the resulting hydrostatic pressure dis-
tribution onto the boundary for the Navier-Stokes solver
and moving the grid to match the free surface while con-
forming to any geometry that intersects the free surface.
The current free surface implementation is loosely cou-
pled with the Navier-Stokes solver and makes the assump-
tion that the flow is steady-state.
INTRODUCTION
Modeling the free surface generated by surface ships and
by submerged vessels near the water surface is important
for proper understanding of the flow around these ves-
sels. In particular, the free surface changes the resistance
of the vessels through the water by changing the pres-
sure distribution and the wetted surface area. The free
surface affects the location and magnitude of the vortices
that originate from various locations, including the bow,
appendages and propulsors, which can greatly affect the
performance of the propulsors. The wave signatures of
ships and submarines is also of importance. Surface ves-
sels under investigation by ship designers attempt to miti-
gate the negative influences of the free surface interaction
on cavitation, power requirements and wave signature and
attempt to use the free surface to improve the performance
in other ways. Because of the novelty of many of these
designs, which often are not direct extensions of current
designs, extensive sets of experimental data may not yet
exist, and designers are more apt to rely on data from nu-
merical simulations. Furthermore, the effects on the ship's
performance caused by maneuvering through a change in
direction, or through various sea states, especially in re-
gards to sea-keeping and cavitation in littoral water are of
importance in the designs and can not be readily tested by
experiment.
A long-term goal of the research into numerical
simulations is to develop the ability to study the perfor-
mance of a full ship design, including the interactions of
the various appendages, sonar domes, rudders, shafts and
propulsors, as the ship maneuvers in response to changes
in the settings of the rudders and propulsors. To accom-
plish this goal, efficient and accurate free surface simula-
tions are needed. Other necessary components to achieve
this goal include capabilities to monitor the onset of incip-
ient cavitation, to perform simulation of several vessels in
motion relative to one another and to adapt the grid to cap-
ture vortices in the water.
Currently, the implicit unstructured code being
developed by researchers in the Computational Simu-
lation and Design Center at Mississippi State Univer-
sity solves the three-dimensional, incompressible Navier-
Stokes equations, via an edge-based, flux-differencing fi-
nite volume method with Roe-averaged variables. Hyams
(2000a, 2000b) successfully parallelized the code and has
ported the code to a wide variety of high-performance ma-
chines. The unstructured grid generator (Marcum, 1998)
has the capability of building highly stretched high-aspect
ratio grids near viscous surfaces, which include higher or-
der elements such as pyramids and prisms. Other areas of
code development include realistic rudder and propulsor-
induced motion and full-scale Reynolds number simula-
tion.
OCR for page 193
At each time step, the nonlinear free surface al- modeled under these assumptions.
gorithm solves the kinematic free surface equation
BY BY
ot + (u—V=) ~ —(v - Vy) + (w—Vz) ~ = 0 (1)
where Y = Y(x,z,t) is the free surface defined as a
single-valued function over the xz plane, (u,v,w) are
the velocity components in the coordinate directions and
(V~, Vy7 Vz) are the grid velocities. A flow-throughbound-
ary condition for the free surface based on characteristic
variable boundary conditions is used, with the pressure on
the free surface set to P = ~ where P is the pressure
and Fr is the Froude number. After several time steps, the
grid is moved to match the free surface while conforming
to the geometry, using a three dimensional extension of
Farhat's torsional spring method (Farhat, 19981. This grid
movement algorithm is quite robust, allowing for moder-
ate to extreme distortions as required by the free surface
simulations. As the solution converges, the flow along the
free surface becomes tangent to the free surface, and the
nonlinear free surface solution is obtained.
Much of this work is an extrapolation of the free
surface algorithm used within the three-dimensional struc-
tured code UNCLE also developed at Mississippi State
University. Beddhu (1994)developed the structured free
surface solver, using a modified artificial compressibility
formulation. This algorithm was applied to steady and
unsteady flow around the Wigley hull (Beddhu, 1998a),
to the barehull model 5415 series hull (Beddhu, 1998b),
to the Series 60 CB = 0.6 ship (Beddhu, 1998c) and to
a more detailed stern analysis for the model 5415 series
hull (Beddhu, 1999), and these results were presented at
the Gothenburg conference (Beddhu, 2000~. Initial veri-
fication and validation exercises of the unstructured non-
linear free surface algorithm has been presented by Burg
(2002), which includes a grid refinement study for flow
over a submerged NACA0012 hydrofoil, flow around in-
viscid and viscous Wigley hulls, and flow around the bare-
hull model 5415 series hull.
In this paper, results for the DTMB Model 5415
Series hull, which has a transom stern, are presented. The
numerical results are compared against experimental re-
sults which include the hull profile for the barehull 5415
and the stern region topologies for the unpowered and
powered fully-appended 5415. Because of the use of un-
structured grids, the actual geometry with the shafts, struts,
rudder, the gap between the rudder and the boot and the
DTMB propellers 4876 and 4877 is used to generate the
grid. For the unpowered case, the propulsors are removed,
so a steady-state simulation is obtained; for the powered
case, the propulsor is rotating providing an unsteady simu-
lation. Even though the nonlinear free surface implemen-
tation assumes a steady-state flow, it is assumed that the
unsteadiness produced by the propulsor can be adequately
In the following section, the free surface algo-
nthm is presented, which includes the method for solving
the kinematic free surface equation, the imposition of the
hydrostatic pressure on the Navier-Stokes equations and
the grid movement algorithm. Then, the results for the
DTMB Model 5415 series simulations are presented.
NAVIER-STOKES SOLUTION ALGORITHM
Ihree different sets of equations are solved in the process
of simulating the flow around surface vessels and sub-
merged vessels near the water surface. The unsteady in-
compressible Reynolds-averaged Navier-Stokes equations,
which are presented here in Cartesian coordinates and in
conservative form, are solved to determine the velocity
and pressure within the computational domain. Either the
Spalart-Allmaras or the q—w turbulence model is used
to simulate the turbulent viscosity primarily within the
boundary layer, and the kinematic free surface equation
is solved to advance the free surface in time.
GOVERNING EQUATIONS
Assuming that gravity acts in the y-direction (i.e., the ver-
tical direction), the incompressible Navier-Stokes equa-
tions can be expressed in dimensional form, denoted with
superscript *, as
0~* ~v* dw*
dx* by* ~z*
~t* + `'x* (u*2 + p ) + ,'0* (u*v*)
+ ~ (u*w*~=~*V2u*
dv* + ~ ~u*v*y + 00* (V*2 + p + 9Y
+~*(v*w*~=~*V2v*
oot* + oo* (u*w*) + o0* (v*w*)
+ Ad (W*2 + P ) = p*V2w*
oz* ~ POO,
(2)
The variables in the preceding equation are normalized
with respect to a characteristic length scale (L) and free
stream values of velocity (UOO), density (pOO), and vis-
cosity (,uOO). Thus, the Reynolds number is defined as
Re = UooLI z/OO. Pressure is normalized with P = (P* +
Poo9Y*—poO)/pOOU2, where P* is the local dimensional
static pressure. Following (Chorin, 1967), an artificial
time derivative term (63pa/5t, where Pa = P//~3) has been
added to the continuity equation to cast the complete set
of governing equations into a time-marching form. The
OCR for page 194
nondimensionalized equations can be written in integral
form as
The solution algorithm consists of the following
basic steps: reconstruction of the solution states at the
`, _ 1 control volume faces, evaluation of the flux integrals for
/0 Qd;V+ ~ F ndA = /^ G. n~dA (3) each control volume, and the evolution of the solution in
(fit JO JAQ Re JOQ each control volume in time.
where n is the outward pointing unit normal to the con-
trol volume V. The vector of dependent variables and the
components of the inviscid and viscous flux vectors are
given as
Q=~v I
IW
~ 3(~-at)
F n = up + nap
wE) + n:P
o
G · n = need + ray + nz~Z
n~z + nary + nz~z
news + news + nzrzz
where ~ is the artificial compressibility parameter (typi-
cally 15 in this work), it, v, and w are the Cartesian ve-
locity components in the x, y, and z directions, and no,
no, and n: are the components of the normalized control
volume face vector. (3 is the velocity normal to a control
volume face:
Ed = not + nyv + nzw + al (7)
where the grid speed al = —(Vent + Vyny + V:n:).
Note that Vs = Vat + Van + V:k is the control volume
face velocity. The viscous stresses given in Equation 6
are
Tij = (A + fit) (0Xj 0~z )
where ,u and ,~b~ are the molecular and eddy viscosities,
respectively.
NUMERICAL APPROACH FOR NAVIER-STOKES
EQUATIONS
The baseline flow solver is a node-centered, finite volume,
implicit scheme applied to general unstructured grids with
nonsimplical elements. The flow variables are stored at
the vertices and surface integrals are evaluated on the me-
dian dual surrounding each of these vertices. The nonover-
lapping control volumes formed by the median dual com-
pletely cover the domain, and form a mesh that is dual to
the elemental grid. Thus, a one-to-one mapping exists be-
tween the edges of the original grid and the faces of the
control volumes.
RECONSTRUCTION
A higher order spatial method is constructed by extrapo-
`4~ rating the solution at the vertices to the faces of the sur-
rounding control volume. The unweighted least squares
method (solved via QR factorization (Anderson, 1994~) is
used to compute the gradients at the vertices for the ex-
trapolation. With these gradients known, the variables at
(5) the interface are computed with a first order Taylor series
as:
Qf = Qo + VQo r `9y
where r is a vector that extends from node 0 to the mid-
`6' point of the edge associated with the control volume face
In question.
RESIDUAL EVALUATION
The governing equations are discretized using a finite vol-
ume technique; thus, the surface integrals in Equation 3
are approximated by a quadrature over the surface of the
control volume of interest. So, the numerical discretiza-
tion of the spatial terms associated with the control vol-
ume surrounding vertex 0 results in
Sqo + ~0 = 0 (10)
where the spatial residual ~ contains all contributions from
the discrete approximation to the inviscid and viscous terms
(~ = Kiev + Evict. Also, the quantity q is defined as
q = By Qd;V.
SPATIAL RESIDUAL
The evaluation of the discrete residual is performed sepa-
rately for the inviscid and viscous terms given in Equation
3. The Roe scheme (Roe, 1981) is used to evaluate the in-
viscid fluxes at each face of the control volume. The alge-
braic flux vector is replaced by a numerical flux function,
which depends on the reconstructed data on each side of
the control volume face:
~ = 2 (IF f QI ~ + FTQR))—2 A (QR—QI ~ (1 1)
~ ~ ~ ~
where A = CAR-. The matrix R is a matrix constructed
from the right eigenvectors of the flux Jacobian, and A is
a diagonal matrix whose entries contain the absolute val-
ues of the eigenvalues of the flux Jacobian. The eigensys-
tem used in this work is based on that reported in (Tay-
lor, 1991~. Note that A is evaluated with Roe-averaged
OCR for page 195
variables, which is simply the arithmetic average between
left and right solution states in the case of incompressible
flows.
For general element grids, it is expedient to use
only edge-local information to compute the viscous fluxes.
This allows the evaluation of viscous fluxes on each face
of the control volume without regard to the varying ele-
ment types of the mesh. An algorithm in which no ele-
ment information is used outside of metric computations
is termed a "grid transparent" algorithm (Haselbacher, 1999~.
To this end, the viscous fluxes are evaluated directly at
each edge midpoint using separate approximations for the
normal and tangential components of the gradient vector
to construct the velocity derivatives (Marcum, 1997~. Us-
ing a directional derivative along the edge to approximate
the normal component of the gradient and the average of
the nodal gradients to approximate the tangential com-
ponent of the gradient leads to the following expression
(Hyams, 2000b):
VQij ~ VQ + (Qj - Qi - VQ As] -2 (12)
The weighted least squares method is used to evaluate the
nodal gradients in the preceding formula.
TEMPORAL RESIDUAL
After the spatial terms have been suitably discretized, the
time derivative term appearing in Equation 10 must be ap-
proximated. A general difference expression is available
for this purpose (Beam, 1978) (Taylor, 1993), and is given
as follows:
/`qn = B1/Yt 0(Aqn)
1+62 fit (13)
~t ~ (qn) + 82 Con—1
where /`qn = qn+l _ qn. A first order accurate in time
Euler implicit scheme is given by the choices b1 = 1,
62 = 0. Correspondingly, a second order time accurate
Euler implicit scheme is given by B1 = 1, 82 = 1/2.
Since 81 = 1 for both time discretizations used in this
work, Equation 13 can be further simplified:
In = /\t 63 (qn+l) + 82 Akin—1 (14)
Using Equation 10 to replace the time derivative,
Ax — 1+92 my 1 ~n+1 <15y
At 1 + 82 where
By the definition of q, one can write q = QV, where Q is ~0 + (Q + )
the volume averaged solution variable vector 1/V iso QdV.
Then, the following two identities can be formed:
/`qn—1 = An—1~\Qn—1 + Qn~`vn—1 (17)
Inserting the above two identities into Equation 15, one
arrives at the following expression:
Vn+1 /\ Qn _ 82 in—1 /\ Qn—1
/\t
Q [ At ] (18)
+ 1 ~pn+1 = 0
1+82
Now, one must consider the Geometric Conservation Law
(GCL). This statement relates the rate of change of a phys-
ical volume to the motion of the volume faces:
Alit = ~ V · V~dV = ~ Vs ·ridA (19)
According to Thomas and Lombard (1978) and later Janus
(1989), the solution of the volume conservation equation
must be performed in exactly the same manner as the flow
equations to ensure that GCL is satisfied. This proce-
dure ensures that spurious source terms caused by volume
changes are eliminated. Using the same time differencing
expression (Equation 14) to approximate Equation 19,
vn_ 82 In—1 1
1+82 = 1 + ~ GIL (20)
where GIL = I(0) VnO+i1 · ~oi+l' and \(0) is a list
of nodes surrounding node 0. Note that the left hand side
of the preceding equation is exactly the bracketed term in
Equation 18. Replacing the bracketed term and rearrang-
ing slightly gives the final form of the discretization of the
time derivative:
(1 + 6~2 )vn+1 is Qn—82Vn—1 ~ Qn—1
2\t
QngRn+1 + jRn+1 = 0 (21)
For incompressible flows, a divergence-free velocity field
at the end of each Newton iteration is desirable. To this
end, the contribution of the time derivative and GCL terms
to the residual are removed for the continuity equation.
TIME EVOLUTION
A Newton iterative time evolution scheme is applied, which
requires the solution of a sparse linear system at each non-
linear subiteration:
duct n+1 m
_~On+l~m = U0 ~~Qn+l,m (22)
( 1 + 62 ) VOn+l i\ QOn—§2 VOn—1 i\ QOn—1 ~23'
At
/`qn = Vn+l /`Qn + Qn~vn (16) + Q ~O,GCL + ~0
OCR for page 196
Qn+1 m Qn+l,m+1 _ Qn+l,m. Now, ex-
panding the terms and performing the required differen-
tiations of ~ results in the following expression for New-
ton's method:
~ (1 + 82 judo + (Qon+ltm _ QOn)—82vn—1~Qn—1
_
At
+ QO.~O,GCL + ~ Hoi
ick(O)
(1 + 82)Vo +1I
L At +
~ `~j~t am . rlOni+1 1 ~ +1
-
The initial guess to beg~n the iteration is /\Qn+l,° = 0.
In the above formulas, Hoi represents the combined invis-
cid and viscous flux vector from node O to node i. \~0)
represents the set of neighbors for node 0, N~.(O) is the
list of neighbors such that the node label earl < Eros, and
~u(O) is the list of neighbors such that the node label
earl ~ egos
BOUNDARY CONDITIONS
. rl,r~+1 I =
Viscous conditions are enforced by modifying the linear
system such that no change is allowed in the velocity,
and the pressure is driven according to the imbalance in
the continuity equation in the boundary control volume
(Anderson, 1994~. Farfield conditions are handled via a
characteristic variable reconstruction; all boundary condi-
tions are handled in an implicit fashion. The free surface
boundary conditions will be presented below - in essence,
~~Hn+l,m -n+1 ~ at steady-state, the velocity is tangent to the free surface
iEN~(O) ~ Oi no, ~Qu+ltm] (24' diver, ~ nfs = 0) Ad Ale pressure is P = An, where Y
where /`QOn-1 = QOn _ QOn-l. For notational conve-
nience, both the inviscid and viscous terms are collapsed
into a single flux function H. Note that the iteration can
be started by using an initial guess of Qn+1,0 = Qn.
Also, performing only one iteration of Newton's method
per time step (with 1st order time discretization and no
GCL terms) is equivalent to a time linearization of the
spatial terms only. However, writing the method in this
framework is more general than a straightforward time
linearization of the nonlinear terms.
To solve the resulting linear system, a bidirec-
tional Gauss-Seidel solution algorithm is adopted. Split-
ting the matrix into diagonal, upper triangular, and lower
triangular parts as [~] = [L] + [~] + tU], the Gauss-Seidel
sweeps may be written as the following two-step process
per subiteration (k is the linear subiterative index):
tr + v~ ~Qk+2 + tu~ ~Qk = ~n+ltm (25)
tv + Up ~Qk+1 + tr~ Auk+- = ~n+~,m (26)
where it is understood that /`Q is evaluated at the rt +
1, m + 1 time level and Newton subiteration. The diago-
nal, lower, and upper operators are defined as
D ~1 + 82)VoI +
L At
dHo7li+l~m . rl,On+l ~
ic~O) Who ~ ()
`,Hn+l,m -n+i
= ~ Oi · Hi (.)
L= ~
iced (o)
is the elevation of the free surface above the undisturbed
waterline and the Froude number is Fr = Add.
A symmetry plane boundary condition imposes
V; · b = 0 for any arbitrary variable ¢. In addition, no
flow is allowed through a symmetry surface; so, like a
solid wall condition, ~ = 0. To simulate symmetry condi-
tions, a layer of phantom entities that is a mirror image of
the entities inside and connected to the symmetry plane is
created; one element layer is mirrored. Thus, control vol-
umes on the symmetry plane are closed and behave just
as interior control volumes. Phantom nodes created by
the mirroring process are updated in the same hierarchy as
described for phantom nodes generated by the paralleliza-
tion. Note that care must be taken to copy scalars, mirror
vectors, and mirror tensors appropriately to the symmetry
phantom nodes such that no fluxes are allowed through
the plane.
TURBULENCE MODEL
Both the one-equation turbulence model of Spalart and
Allmaras (1992~(Sheng, 1999) and the q—w turbulence
model (Coakley, 1985) are available within the solver. Con-
stants are taken from the version II q—w model given in
(Coakley, 1985~. The interested reader is referred to these
references.
The diffusive terms in both turbulence models
are discretized in the same manner as the viscous terms
for the mean flow, and the convective terms are computed
<27y via pure unwinding. Appropriate consideration is given to
maintain positive operators in the formation of the Jaco-
(28) bian matrix for the implicit solution of the transport equa-
tionfs). The respective turbulence models are incoroo-
. I . ] . 1 ~ ~
n+l~m nr,.+1 rated with the mean flow solution in a "loosely-coupled"
Oi Em oi `.y (29) procedure;thatis,the core governing equations are solved
OCR for page 197
I ~ ~ ~ ~ ~ ~ in, rQ ~q At
I ;::: lj';^IA ~ ~ ~ :: ~ ~ L ~ ~ J
s~ ma : ~ : :: : ::~:__
r ~ t
~ em l ~ jut
T;m" atnr. ~ ~ i:: :.: . :: : : : ~ : :~: ~ l
~ ~ ~~ v—_
' ~ ~ ~ ~ Ha [Q. VQ, <*
I Cycle over local subdomains~:~:: i: ~ I r- *
~ Lit, lit]
~~: : A? ~ T:~ ~ ?~: : ~
_ NewtonIteration : ~ :~:~ : ~~
Compute VQ em VQt
Compute ~
Compute [A]
~ I
Linear ~ subiterations: ~: : ~ ~ ~ :1
_-
Matrix-Vector Muliiply - _ A <)t
~ ~ ~ :~: ~: : : ~ ~ ::~: ~ :: ~ : ~ ~ ::: :: :~: ~ 1
~~SolYe~I\frbulence~:~Model : ::~: ,,,,B,,~
: :: : ~ [u pt]t
~ ~~ ~ ~~ ~ ~ ~ 1 ~
~~ol~re:~::~nemat~c free Menace - ~ t
_ _—[Yf~~Pfs]
* sequential updating hierarchy
t concurrent updating hierarchy
Figure 1: Iteration hierarchy used for the parallel unstruc-
tured solution algorithm
first, then the turbulence model is solved independently.
This procedure allows for easy interchange of the turbu-
lence models.
PARALLELIZATION
For quick turnaround time in a design environment, it is
essential to parallelize the flow solution algorithm. The
present parallel unstructured viscous flow solver is based
on a coarse-grained domain decomposition for concurrent
solution within subdomains assigned to multiple proces-
sors. The solution algorithm employs iterative solution of
the implicit approximation, with the concurrent iteration
hierarchy as shown in Figure 1. Also, domain decompo-
sition takes place with each node in the domain uniquely
mapped to a given task. The code employs MPI message
passing for interprocessor communication.
In general, the parallelization of an existing val-
idated flow solver should satisfy several constraints. First
and most important, the accuracy of the overall numerical
scheme must not be compromised; i.e., the solution com-
puted in parallel must have a one-to-one correspondence
with the solution computed in serial mode. For the cur-
rent numerical algorithm, this ability has been shown in
(Hyams, 2000~. Also, the code must be efficient in its use
of computational resources. This characteristic is mea-
sured in terms of memory usage and scalability, as well
as the fact that the parallel code should degenerate to the
serial version if only one processor is available. Finally,
the consequences of the inevitable domain decomposition
should not seriously compromise the convergence rate of
the iterative algorithm.
NONLINEAR FREE SURFACE ALGORITHM
A nonlinear free surface is obtained for a steady-state sim-
ulation by solving the kinematic free surface boundary
condition at each time level and imposing the hydrostatic
pressure distribution based on the new free surface eleva-
tion onto the free surface boundary within the mean flow.
After several time steps, typically on the order of 200, the
grid is moved to match the free surface elevations while
conforming to the surfaces that intersect the free surface,
such as the hull of a ship or the sail of a submarine; and as
the loosely coupled interaction between the free surface,
the Navier-Stokes equations and the grid movement al-
gorithm converges, the solution approaches the nonlinear
free surface solution.
GOVERNING EQUATIONS
In deriving the kinematic free surface boundary condition,
the free surface ~(x, y, z, t) = 0 is considered as a mate-
rial surface (i.e., no flow travels through the surface), so
D71 = 0 (30)
or
lilt + (a- Vx) ,~,71 + (v- Vy) '971 + (w-Vz) ,0~71 = 0 (31)
By making the assumption that the free surface 77 can be
expressed as 71(X, y, z, t) = y—Y(x, z, t) = 0, this equa-
tion reduces to
lit + (a - V=) `~ - (v - Vy) + (w - Vz) ,, = 0 (32)
At steady-state, the kinematic boundary condition becomes
(u+at) (~,-1' Liz) 0 (33)
(U+at)ootrifs =0
where or = )/1 + ( OoY ) + ~ BY ~ , so the flow is tangent
to the free surface at steady-state.
The pressure in the numerical simulation at the
free surface is set based on the free surface elevation. For
the original, dimensional equations, the pressure at the
free surface is the atmospheric pressure POO- Using the
nondimensionalization, p = P +Pp°°9u2-P°° with P* =
POO and ~ = Y. the pressure is clearly set to P = An,
where the Froude number is up .
NUMERICAL APPROACH FOR FREE SURFACE
The kinematic equation for the linear free surface is solved
via a Galerkin approach. In this approach, an algebraic
equation is obtained for each node j by multiplying the
OCR for page 198
governing equation by a weight function ~j(x,z), ap-
proximating the terms in the governing equation by using
linear and bilinear interpolating functions Minx, z) and in-
tegrating over the computational domain, or
~Q~i(X~Z) (at + (U-V~) o~ - (V_Vy)
~ \
+(W-Vz),g )=0 (34)
where the interpolated free surface elevation is Y(x, z) =
Pi Yi (i (x, z) and the velocities are of the form USA, z) =
Ui~i(X~Z).
For the triangular elements, the resulting equa-
tions are integrated either exactly or using 7 point Gauss
quadrature, and for the quadrilateral elements, the inte-
grals are solved either using 3 or 6 point Gauss quadra-
ture in both directions or by decomposing the quadrilat-
erals into 4 overlapping triangles and evaluating them as
triangles. The examples given herein use 7 point Gauss
quadrature for the triangular elements and 6 point Gauss
quadrature for the quadrilateral elements. A first-order
backward time discretization of Equation 34 is used for
the temporal derivative, where CITY = Y lit Y and the
spatial terms evaluated at time level r' + 1. This dis-
cretization results in a linear algebraic system of the form
~fs(Yn+t, Y~) = 0. Since these equations are linear in
the unknown variable Yin+i at each time level (i.e., the ve-
locities are frozen), the Jacobian matrix is calculated only
once per time level. A Gauss-Seidel iterative solver simi-
lar to the one used to solve the discretized Navier-Stokes
equations is used to identify the solution of the kinematic
free surface equation by solving
of the viscous surface is on the order 10-6 for model scale
simulations and 10-8 for full scale simulations.
FREE SURFACE BOUNDARY CONDITION
At the free surface boundary, the pressure is set to the hy-
drostatic pressure of P = ma. Determining the velocity
to impose on the boundary is more difficult. If the ve-
locity were forced to be tangent to the free surface (i.e.,
(n + at) · Ifs = 0), then the free surface equation would
reduce to MY = 0, and the free surface could not evolve.
Thus, this tangency condition is relaxed to allow the free
surface to evolve.
This boundary condition is implemented using
a characteristic variable boundary condition, similar to
the derivation used for the farfield boundary condition.
For hyperbolic systems, flow information travels along
characteristics; and for three-dimensional incompressible
flow, three characteristics originate upstream, and one orig-
inates downstream. If the flow is traveling out of the do-
main (i.e. (u + at) · n > 0), then only one characteristic
needs to be specified from the outside, which is derived
from the hydrostatic pressure. If the flow is traveling into
the domain, three characteristics need to be specified from
the outside. The characteristic variables are derived from
considering the inviscid fluxes across the the boundary
face, via
~Q + VF rl = 0 (36)
where r' is the normal to the grid for each boundary face.
The inviscid flux term can be rewritten as
VF = L3Q VQ = RAR-I VQ (37)
where A is a diagonal matrix with the eigenvalues 0, 0,
`~ (3 + c and (3—c, where (3 = uric + vr~y + wry + at and
n+t,m _ ~ ~yn+i,m yin ~35> ~ 2
~yn+i Y —— fey ~ J ~ J c=~(~)_ I) +~. PremultiplyingbyR-t,
where yn+i,m+i = yn+i,m + /`yn+i,m and yn+i = -~0Q AR-IV —0 38
yn+itm+i when ~~Ayn+~,m~ < tolerance. R ~ + Q ~—
Within the viscous boundary layer, special care
is needed in solving the kinematic free surface equation.
On viscous surfaces, the flow velocity plus the grid ve-
locity (i.e., u + all is set to zero, which prevents the free
surface from moving at the viscous surface. Physically,
however, surface tension effects force the free surface to
rise and fall along the viscous surfaces. Since the solver
does not simulate surface tension, another method must
be employed to move the free surface near the viscous sur-
faces. Following the method used in the structured solver,
the free surface at a certain distance from the viscous sur-
face is extended at a constant height to the viscous surface.
This distance is typically on the order of 2 x 10-4 times
the characteristic length L, whereas the point spacing off
Freezing the matrix Ro ~ = R-t f~Qin`),
ARE Q +~V (Ro 1Q) ri = 0 (39'
Using the definition W(`Q`) = Ro-iQ' Equation 36 decou-
ples into four equations for the four characteristic vari-
ables W as
0~9~ +AVW ri = 0 (40)
Three of the characteristics are obtained from upstream
values, and one of the characteristics is obtained from
downstream values. In the case of a boundary face, there
are a set of characteristics associated with the inside flow
OCR for page 199
I.
Wln
W = win
w4ut
characteristic variables are chosen as
we) and a set of characteristics associated with the into the volume grid, the use of a linear spring analogy
outside flow W°Ut (Q°Ut). The flow variables on the bound- where each edge is a spring whose stiffness is determined
ary face are determined from the appropriate characteris- by the length of the edge, a torsional/linear spring anal-
tic variables via Q = RoW. For flow traveling out of the ogy where the angle between the edges affects the stiff-
domain, the characteristic variables are chosen to be ness of the springs in the mesh, or solving the linear elas-
i~ ticity equations to propagate the perturbations within the
wag grid. The spring analogy is computationally efficient but
W = w2n <41y is not robust. Both the spring analogy and the Laplacian
wout solver generate elements with negative volumes for mod-
4 crate amounts of movement, on the order of the size of 3
because the fourth eigenvalue (3—c is negative, while to 4 elements. Solving the linear elasticity equations to
the other eigenvalues are positive. For flow traveling into move the grid is robust, but the computational cost asso-
the domain, only the third eigenvalue is positive, so the elated with this method is the primary drawback.
Thus, following the work of (Farhat, 1998) and
extending it to three dimensional grids, the torsional/linear
spring analogy has been developed and applied to the prob-
`42' lem of moving the grid to match the linear free surface el-
evations. For tetrahedral grids, this method is quite robust,
allowing severe distortions on the surface while maintain-
ing positive volume elements. In practice this algorithm
For farfield boundary conditions the characteristic vari- . .'
' allows up to approximately 80% compression of elements.
ables associated with the outside flow are taken as the . .
Within the boundary layer, (led, near the viscous surfaces),
freestream variables. For free surface boundary cond~-
.. . . the grid moves the same as the points on the boundary.
tons the characteristic variables for the outside flow are
~ This grid movement method Is computat~onally costly but
determined from the hydrostatic pressure and the only avail-
prov~des excellent robustness for the free surface solver.
able velocity information which is the velocity on the ~n-
~ The linear so arms method for moving an unstruc-
side. However this velocity is modified by removing a . . . ~ ~
' tured arid Is presented In (Bating, 1989) where each edge
portion of the velocity component that is not tangent to . ~ . . . . . .
In the cad Is replaced by a linear spring whose stiffness Is
the grid, via
I-WO?lt
wout
wont
~ = (A + at)—by ((n + at) n) (n + at) (43)
where
:_ ~—(Nat) ~ if —1 < (~+at) ~ < 0
t1 if (Nat) ~ < -1
(44)
When y is 1, the velocity imposed from the outside is tan-
gent to the grid (i.e., clout ~ = 0~. When the velocity
was not modified or constrained in some fashion, the free
surface algorithm became unstable for flow into the do-
main due to the inconsistency of using downstream infor-
mation. This modification provides enough control over
the velocity to maintain stability of the algorithm, and if
the grid is allowed to move to match the free surface, then
this modification ensures that the flow will be tangent to
the free surface at convergence.
GRID MOVEMENT ALGORITHM
After several time steps, the grid is moved to match the
free surface while conforming to any solid surfaces inter-
secting the free surface, with displacements on the surface
being propagated into the volume grid. Several methods
are available for this grid movement, including the use of
a Laplacian solver to propagate the surface perturbations
inversely proportional to the length of the edge. Thus, for
the edge connecting nodes i and j, the stiffness kid of the
. .
spring IS
( (Xi—X j ) 2 + (hi - y j ) 2 ) p/2 = ],[ p (45)
where the coordinates of nodes i and j are (xi, pi) and
(Xj ~ yj) and p is a predetermined coefficient, usually 1 or
2, and lid = \/((xi—Xj)2 + (Pi—yj)2). Given a set
of nodal displacements or forces acting on the boundary
of the computational grid, the following equations for the
interior displacements are solved iteratively until all the
~ . . . .
forces are In equip. ~ Pram
i ~ k
i ti
~yn+1 = E; kid EVE
i 3}
1
(46)
where i is summed over all edges connected to node j.
This method is easily extended to three dimensions and
is computationally efficient requiring only a few Jacobi
iterations to achieve an acceptable level of accuracy. Sev-
eral researchers have used this method to move nodes in
unstructured grids for simulating flows around objects in
relative motion (Singh, 1995) and for design optimization
(Anderson, 1997).
OCR for page 200
./~a
(xi,yi) Xj-Xi
~ (Xj ,yj )
Yj-Yi 2
Figure 2: A Displacement in y-direction results in dis-
placements in x and y directions
Technically speaking, these equations do not sim-
ulate the behavior of a network of springs because there
is no interaction between the x and y coordinates. A dis-
placement in one coordinate will not influence the loca-
tion in the other coordinate, as would be the case for a
network of springs. An example of this coupled inter-
play is shown in Figure 2, where a displacement in the
y-direction is applied to node j which results in a displace-
ment in both directions for node i. To simulate the behav-
ior of a network of springs, the following set of equations
for each spring must be summed over all springs:
forces in the x-direction:
kij [~Xi—Axj) cos2 ~
+ (~\yi—/\yj) COS Ol sin Ol]
forces in the y-direction:
kij [taxi—A~j) cos ~ sin
+ (phi—Ayj) sin2 a]
For either representation, however, this method
often fails for complicated geometries and for large changes
in the boundary, because negative areas are generated when
a node crosses over an edge in the grid. The creation
of negative areas is illustrated in Figure 3, where node
1 is pushed downwards and node 4 is forced to cross the
edge between nodes 2 and 3. The reason for this failure
is that the stiffness in the linear spring method prevents
two nodes from colliding but does not prevent a node
from crossing over an edge. As two nodes get closer to-
gether, the stiffness increases without bound preventing
the collision; but there is no mechanism to prevent a node
from crossing an edge, because these crossovers can occur
without the stiffness increasing without bound.
To provide a more robust movement algorithm
for unstructured grid, Farhat (1998) developed an algo-
rithm to prevent a node from crossing an edge by using
torsional springs around each node, as shown in Figure 4.
~ ~L~ ~\~
3 2 3 2 4 3
1
Figure 3: Negative Areas Produced by Linear Spring
Method
k
Linear Spring :~
~ \
~ \
~)X/\e ilk \
ire ~ .
~_J
Torsional Spring Linear Spring
Figure 4: Placement of Linear and Torsion Springs
The stiffness of the torsion spring Cijk is inversely pro-
portional to the sine of the angle, so that the stiffness
grows without bound as the angle decreases towards zero
or increases towards 180°, or
Cijk = j 2 Jo,, (48)
where 83jk iS the angle centered at node i formed from
the edges ij and ik. Thus, as a node moves towards an
edge, the angle goes towards zero, and the stiffness of the
torsion spring grows. The sine of the angle is squared to
prevent a negative stiffness.
In (Farhat, 1998), the equations for the torsional
spring methods are derived for two-dimensional grids, and
an iterative solution method is presented. In the present
work, these methods have been extended into three di-
mensions, which provides an adequate level of robustness
for the grid movement algorithm. The derivations of the
three-dimensional torsional spring equations are beyond
the scope of this paper and will be presented elsewhere.
EXAMPLES
The two examples given in this paper deal with the fully-
appended DTMB Model 5415 series hull, one without
propulsors and the other with propulsors rotating. Ex-
perimental data for comparison was provided by Ratcliffe
OCR for page 201
at the Naval Surface Warfare Center Carderock Division
(NSWCCD) (Ratcliffe, 1990, 2000), which consisted pri-
marily of the profile of the free surface intersection with
the hull and free surface elevations in the stern of the hull.
The Model 5415 series hull has a transom stern, which
produces a sizeable "rooster tail" in the water surface be-
hind the stern. The parameters for the various runs are
given in Table 1.
r Parameter | Value l ~
I 1 1
Speed TOO (m/sec) 2.06
Ship Length L (m) 5.72
Propulsor RPM 436
Advance Ratio (Rev. per L) 20.177
Froude Number (Fr) 0.28
Reynolds Number (Re) 11.655 Million
Omega 126.59 red/second
At for 1.5 degrees 2.068089 x 10-4
Table 1. Flow parameters for Model 5415 Series
simulations.
The grids used within the following simulations were gen-
erated by AFLR3 (Marcum, 1998) and consist of prisms
in the boundary layer generated from the triangles on the
viscous surfaces, tetrahedra and pyramids in the transition
region at the edge of the boundary layer grid and tetrahe-
dra outside of the boundary layer. On the free surface
grid, quadrilaterals exist within the boundary layer region
and triangles exist outside of the boundary layer. The use
of tightly packed prisms within the boundary layer allows
for good resolution of the boundary layer effects and typ-
ically generates y+ values on the order of 1.0.
Information about the size of the grids used in
these simulations is provided in Table 2. For the powered
and unpowered simulations, the same grid is used except
for the regions around the two propulsors. For the un-
powered grid, the propulsors were replaced by hubs and
were not rotated. For these simulations, the grid included
both sides of the hulls, so that maneuvers using these grids
could be studied in the future. The hull, struts, rudders
and propulsors were symmetric about the centerline, but
because of the use of unstructured grids, the grid was not
symmetric about the centerline. These simulations were
run locally on a LINUX-based cluster running MPI-Pro.
I l Unpowered l Powered l
Hull Profiles for Model 5415 Hull
0.018
0.015
0.012
0.009
0.006
0.003
o
-0.003
-0.006
0 0.2
Fr=0.28, Re=11.355 M
' I ' I ' 1 ' 1 '
Barehull (Structured)
- Unpowered
- Powered
0 DTMB Experimental Data
0.4 0.6 0.8 1
Distance from Bow of Ship
Figure 5: Free Surface / Hull Intersection
which the simulation was run at a constant time step, so
that a time accurate simulation could be achieved for the
powered case. To ensure a valid comparison, both simu-
lations were run at a constant time step. For the powered
case, the images and plots use instantaneous data rather
than time-averaged data.
UNPOWERED FULLY-APPENDED 5415
The unpowered fully-appended 5415 simulation used the
same grid as the powered fully-appended 5415 simulation
with the exception that the grid around the propulsor/hub
was changed and did not rotate. Hence, the free surface
grid was identical and the forces acting on the free surface
near the bow should be almost identical between the sim-
ulations. The profile of the intersection between the free
surface and hull is given in Figure 5. This profile com-
pares the results for the powered and unpowered simula-
tions with the experimental results from the Naval Surface
Warfare Center Carderock Division runs (Ratcliffe, 1990,
2000) and with the computational results of the structured
flow solver UNCLE (Beddhu, 1999, 2000~. The profile
for the powered and unpowered cases agree with each
other nicely until approximately y/L = 0.6 at which point
the unsteady influence of the propulsors are felt. The
bow wave in the unstructured grids are higher than for
the structured grid and seem to agree more with the bow
wave elevation of the experiment. The drop behind the
bow wave is more severe for the unstructured simulations
Nodes 4251444 5437074 - this difference may be caused by a lack of refinement in
Tetrahedra 8195672 9912263 this region.
Pyramids 23592 35101 Figure 6 shows the bow wave for the unpowered
Prisms 5446868 7160111 simulation head on, as well as the bulbous bow underneath
Processors 32 42 the free surface. In Figure 7, the surface grid near the bow
Table 2: It formation abo at Grids. is shown. This region undergoes the largest displacement
due to the changes in the free surface, and yet the grid
In each case, the simulation was run with a constant CFL movement algorithm is able to produce nice quality in the
until the gross features in the flow field developed, after
OCR for page 202
Figure 6: Bow Wave for Model 5415
Figure 7: Grid Moved to Match Free Surface Elevation
near Bow
distorted grid.
The hullform was symmetric about the center-
line, but due to the use of unstructured grids, neither the
volume nor the surface grids were symmetric. This asym-
metry in the grid produces an asymmetry in the discretiza-
tion error within the code, but the free surface, even in the
stern region, regains symmetry as the solution approaches
steady-state. In Figures 8 and 9, contour plots of the free
surface elevations around the unpowered 5415 are pre-
sented showing the Kelvin wake pattern originating from
the bow and the "rooster tail" in the stern region. For the
unpowered case, the solution is nearly symmetric. Due to
the coarseness of the grid away from the body, the free
surface dissipates approximately one body length behind
the hull. By extending the resolution on the free surface,
the free surface waves could be carried further. The free
..f ~ ,
Figure 8: Contours of Free Surface for Unpowered 5415
Figure 9: Contours of Free Surface for Stern Region of
Unpowered 5415
OCR for page 203
Unpowered Model 5415 Stern Topology
O.05
-D.05
—1.05 1.1 X,L
y
0.~e
0.~07
0.006
0.005
O.~04
0.003
0.002
0.00 1
.~D 1
.002
.003
Figure 10: Comparison with Experimental Data for Un-
powered 5415
surface in the stern region was also compared Wltn tne ex-
perimental data from Ratcliffe at NSWCCD and is shown
in 10. The agreement for the unpowered case is remark-
able, with similar features and elevations. There is a no-
ticeable shift forward for the computational results.
Finally, the free surface-stern intersection for the
powered and unpowered simulations are given in Figure
11, as well as the profile for the structured barehull simu-
lation. The structured barehull simulation was run using a
symmetry plane, so its profile was mirrored about z = 0.
The unpowered simulation is almost symmetric about the
stern, whereas the unsteadiness in the powered simulation
has created some asymmetry (this asymmetry is probably
a result of not being as highly converged as the unpowered
case). For the powered simulation, this stern profile is an
instantaneous cut, and is not as converged as the unpow-
ered simulation. The profiles for the simulations differ
greatly. The reason for this difference lies primarily with
the difference in the flow field behind the rudders. For the
barehull, the struts, shafts and rudder have been removed,
so the how has not been slowed nor straightened, thus
the flow passing the stern is faster for the barehull case
than for the unpowered case. For the unpowered case, the
flow has been slowed by these appendages and straight-
ened by the rudder, so that the flow is slower, and the free
surface is not sucked down. And for the powered case,
the propulsors have accelerated the flow past the rudders,
which sucks the free surface down. In Figure 12, the ver-
tical velocity on the free surface and for a cutting plane
just aft of the stern are presented from a view beneath the
hull, showing that the velocity is relatively uniform be-
tween the rudders and the vertical velocity is nearly uni-
. .. . . . . · ..
Free Surface Stern Intersection
0.002
,, /
,,
-0.002
I A
v
-0.004
-0.04 -0.02
i.\ ,j
~ ~ , ~ —,j ,~
~ if ,' -it—·_,,,
., "
ll
Powered
Unpowered
---- Barehull
I'
.
0 0.02
z/L
Figure 11: Free Surface / Stern Intersection
0.04
Figure 12: Vertical Velocity for the Unpowered 5415.
form along the free surface behind the stern. In Figure 13,
the instantaneous vertical velocity for the powered case is
shown from the same location and on the same scale. The
swirl in the velocity field generated by the propulsors are
clearly seen and indicates that the propulsors are rotating
outboard over the top. The velocities are much more ex-
treme at the free surface behind the stern, which results in
the drop in the free surface elevation along the stern.
POWERED FULLY-APPENDED 5415
The powered fully-appended 5415 simulation was per-
formed using the same parameters as the unpowered case,
with the exception that the powered case rotated the ac-
tual propulsors used in the experiment. As a result, a pe-
riodic solution is expected for the powered case. Figure
14 shows the drag on the port propulsor as well as the to-
tal drag on the hull, struts, shafts and rudders and on the
OCR for page 204
0.0008
~ 0.0004
c,
n
I I l ~ / / /
Figure 13: Vertical Velocity for the Powered 5415 (Instan- Figure 15: Contours of Free Surface for Stern Region of
taneous). Powered 5415
0.0012 —
Drag (Thrust) for Powered 5415
Nondimensional
-
Powered Fully-Appended 5415 Stern Topology
51 52 53 54
Rotation Count
Figure 14: Drag (Thrust) on Powered 5415.
Overall 5415. The drag on the propulsors is negative be-
cause the propulsors are introducing thrust into the simu-
lation. The drag is clearly periodic, approximately repeat-
ing every 1/5 of the propulsor rotation, since the propul-
sors had 5 blades. The overall drag is reducing slightly in-
dicating that a highly converged solution has not yet been
achieved.
A contour plot of the free surface in the stern re-
gion for the powered case is shown in Figure 15 and shows
that the flow field is nearly symmetric. Finally, the stern
region for the powered case is compared against the exper-
imental results from NSWCCD in Figure 16. The exper-
imental results are time-averaged, whereas the computa-
tional results are instantaneous. The height of the "rooster
tail" for the computational results is higher than for the
experimental results. This difference may be caused by
V3
· 0.0080
0.0072
0.0064
0.0056
0.0049
0.004 1
0.0033
0.0025
0.0017
0.0009
0.000 1
-0.0006
-0.00 1 4
-0.0022
_ -0.0030
55 56 1.D5 1.1 x`,Ll .15 1.2
Figure 16: Comparison with Experimental Data for Pow-
ered 5415
the assumption in the nonlinear free surface algorithm that
the free surface is steady-state when the flow is clearly un-
steady.
CONCLUSIONS
A robust, nonlinear free surface algorithm coupled with an
incompressible Navier-Stokes solver on three-dimensional
unstructured grids has been presented. The underlying un-
structured solver has been developed to simulate turbulent
boundary layer effects and is capable of performing rud-
der and propulsor induced maneuvering. The addition of
nonlinear free surface capabilities is a necessary step to-
wards the simulation of maneuvering surface vessels.
The nonlinear free surface algorithm is an ex-
tension of the algorithm used in the structured code UN-
OCR for page 205
CLE developed by researchers at Mississippi State Uni-
versity. The kinematic free surface equation is solved at
each time level, and the hydrostatic pressure is imposed
on the Navier-Stokes solver via a characteristic variable
boundary condition. After several iterations, the grid is
moved to match the free surface while conforming to the
underlying geometry, using a three dimensional extension
of the torsional spring method (Farhat, 1998~. By allow-
ing the grid to move, a nonlinear free surface solution is
obtained.
In this paper, the free surface around the fully-
appended DTMB Model 5415 series hull with and with-
out rotating propulsors has been simulated. The hull pro-
files for both simulations upstream of the propulsors is in
close agreement with each other and adequate agreement
with the experimental data and with computational results
using the structured solver UNCLE. The contours of the
free surface elevations for both simulations show that the
free surface simulation is symmetric even though the grid
was not symmetric. And the flow in the stern region for
both geometries is symmetric and agrees well with the ex-
perimental data generated by Ratcliffe at NSWCCD. Dif-
ferences in the free surface at the stern can be attributed
to differences in the velocity field between the propul-
sors and between the rudders. More work is necessary
in demonstrating that the powered simulation does have a
periodic free surface, while the current work has shown
that the drag on the geometry is nearly periodic.
Future work includes performing these same sim-
ulations at full scale Reynolds number, analyzing the ef-
fects of stern flaps on the Model 5415 series hull forms
and Series 60 CB = 0.6 ships, simulating flow around
unconventional hull forms and performing prescribed and
propulsor/rudder induced maneuvers.
ACKNOWLEDGMENTS
This research was sponsored by the Office of Naval Re-
search with Dr. L. Patrick Purtell grant monitor. This sup-
port is gratefully acknowledged. Special thanks to Toby
Ratcliffe at Naval Surface Warfare Center Carderock Di-
vision for providing the data files for the experimental
data and to Stephen Nichols at the Computational Sim-
ulation and Design Center for many hours of conversa-
tions about free surface issues, especially pertaining to the
structured code UNCLE.
References
[1] Anderson, W. K., and Bonhaus, D. L., "An implicit
upwind algorithm for computing turbulent flows on
unstructured grids," Computers in Fluids, Vol. 23,
No. 1, Jan. 1994, pp. 1-21.
[2] Anderson, W. K., and Venkatakrishnan, V., "Aero-
dynamic Design Optimization on Unstructured
Grids with a Continuous Adjoint Formulation,"
Proceedings of the 35th AIAA Aerospace Sciences
Meeting end exhibit, AIAA Paper 97-0643, Jan.
1997.
[3] Batina, J. T., "Unsteady Euler Airfoil So-
lutions Using Unstructured Dynamic Meshes,"
27th AIAA Aerospace Sciences Meeting AIAA Pa-
per 89-0115, Jan. 1989.
[4] Beam, R. M., and Warming, R. F., "An implicit
factored scheme for the compressible Navier-Stokes
equations," AIAA Journal, Vol. 16. No. 4, April
1978, pp. 393-402.
[5] Beddhu, M., Taylor, L. K., and Whitfield,
D. L., "A Time Accurate Calculation Pro-
cedure for Flows with a Free Surface Using
a Modified Artificial Compressibility Formula-
tion," Applied Mathematics and Computation, Vol.
65, 1994, pp. 33-48.
[6] Beddhu, M., Jiang, M. Y., Taylor, L. K., and Whit-
field, D. L., "Computation of Steady and Un-
steady Flows with a Free Surface Around the Wigley
Hull," Applied Mathematics and Computation, Vol.
89, 1998,pp.67-84.
[7] Beddhue, M., Jiang, M. Y., Whitfield, D. L., Taylor,
L., K., and Arabshahi, A., "CFD Validation of the
Free Surface Flow Around DTMB Model 5415 Us-
ing Reynolds Averaged Navier-Stokes Equations,"
Proceedings of the Third Osaka Colloquium on Ad-
vanced CFD applications to Ship Flow and Hull
Form Design, 1998.
[8] Beddhu, M., Nichols, S., Jiang, M. Y., Sheng, C.,
Whitfield, D. L., and Taylor, L. K., "Comparison
of EFD and CFD Results of the Free Surface
Flow Field about the Series 60 CB = 0.6 Ship,"
Proceedings of the Twenty-Fifth American Towing
Tank Conference, 1998.
[9] Beddhu, M., Jiang, J. Y., Whitfield, D. L.,
and Taylor, L. K., "Computation of the Wet-
ted Transom Stern Flow over Model 5415,"
Proceedings of the Seventh International Conference
on Numerical Ship Hydrodynamics 1999.
[10] Beddhu, M., Pankajakshan, R., Jiang, M. Y.,
Taylor, L., Remotigue, M. G., Briley, W. R., and
Whitfield, D. L., "Computation and Evaluation
of CFD Results for Practical Ship Hull Forms,"
Proceedings of Gothenburg 2000: A Workshop on
Numerical Ship Hydrodynamics, 2000.
OCR for page 206
[11] Burg, C. O. E., Sreenivas, K., Hyams, D., and [20]
Mitchell, B., "Unstructured Nonlinear Free
Surface Solutions: Validation and Verification"
Proceedings of the 32nd AIAA Fluid Dynamic
Conference and Exhibit AIAA Paper 02-2977 June
, ,
2002.
[12] Chorin, A. J., "A numerical method for
solving incompressible viscous flow problems,"
Journal of Computational Physics, Vol. 2, 1967, pp.
12-26.
[13] Coakley, T. J., Hsieh, T., "A comparison be-
tween implicit and hybrid methods for the
calculation of steady and unsteady inlet flows."
Proceedings of the 21st AIAA/SAE/ASME/ASEE
Joint Propulsion Conference, AIAA Paper 85- 1125,
July 1985.
[14] Farhat, C., Degand, C., Koobus, B., and Le-
sionne, M., "Torsional Springs for Two-
Dimensional Unstructured Fluid Meshes,"
Computer methods in Applied Mechanics and Engin- [24]
sing, Vol. 163, 1998, pp. 231 -245.
[15] Haselbacher, A., McGuirk, J. J., and Page, G. J., "Fi-
nite volume discretization aspects for viscous flows
on mixed unstructured grids," AIAA Journal, Vol. [25]
37, No. 2, Feb. 1999, pp. 177-184.
[16] Hyams, D. G., Sreenivas, K., Sheng, C., Bri-
ley, W. R., Marcum, D. L., and Whitfield,
D. L. "An investigation of parallel implicit
solution algorithms for incompressible flows
on multielement unstructured topologies,"
Proceedings of the 38th AIAA Aerospace Sciences
Meeting and Exhibit, AIAA Paper 2000-0271,
January 2000.
[17] Hyams, D. G., "An Investigation of Parallel Implicit
Solution Algorithms for Incompressible Flows on
Unstructured Topologies," PhD thesis, Mississippi [27]
State University, May 2000.
[18] Hyams, D. G., Sreenivas, K., Sheng, C.,
Nichols, S., Taylor, L. K., Briely, W. R.,
and Whitfield, D. L., "An Unstructured
Multielement Solution Algorithm for Com-
plex Geometry Hydrodynamic Simulations,"
Proceedings of the Twenty-Third Symposium on
Naval Hydrodynamics, September, 2000.
[19] Janus, J. M. " Advanced 3-D CFD Algorithm for
Turbomachinery," PhD thesis, Mississippi State
University, May 1989.
Marcum., D. L. Private communications, Engineer-
ing Research Center for Computational Field Simu-
lation, Mississippi State, MS, December 1997.
[21] Marcum, D. L. "Unstructured Grid
Generation Using Automatic Point
Insertion and Local Reconnection,"
The Handbook of Grid Generation, edited by J.F.
Thompson, B. Soni and N.P. Weatherill, CRC Press,
Boca Raton, 1998, Section 18.
t22] Ratcliffe, T. J., and Lindenmuth, W. T., "Kelvin
Wake Measurements Obtained for Five Surface Ship
Models," DTRC-89/038, January, 1990, David Tay-
lor Model Basin.
[23] Ratcliffe, T. J., "An Experimental and Compu-
tational Study of the Effects of Propulsion on
the Free-Surface Flow Astern of Model 5415,"
Proceedings of the Twenty-third Symposium on
Naval Hydrodynamics, September 2000.
Roe., P. L., "Approximate Riemann solvers,
parameter vectors, and difference schemes,"
Journal of Computational Physics, Vol. 43, pp.
357-372, 1981.
Sheng, C. Hyams, D., Sreenivas, K., Gaither, A.,
Marcu, D., Whitfield, D., and Anderson, W. K.,
"Three-dimensional incompressible Navier-Stokes
how computations about complete configurations
using a multiblock unstructured grid approach,"
Proceedings of the 37th AIAA Aerospace Sciences
Meeting and Exhibit, AIAA Paper 99-0778, January
1999.
[26] Singh, K. P., Newman, J. C., and Baysal, O., "Dy-
namic Unstructured Method for Flows Past Multiple
Objects in Relative Motion," AIAA Journal, Vol.33,
No. 4, April 1995, pp. 641-649.
Spalart, P. R., and Allmaras, S. R., "A one-equation
turbulence model for aerodynamic flows," AIAA
Paper 92-0439, 1992.
[28] Taylor, L., K., "Unsteady Three-Dimensional In-
compressible Algorithm based on Artificial Com-
pressibility," PhD thesis, Mississippi State Univer-
sity, 1991.
[291 Taylor, L. K., Busby, J. A., Jiang, M. Y., Arab-
shahi, A., Sreenivas, K., and Whitfield, D. L.,
"Time accurate incompressible Navier-Stokes
simulation of the fiapping foil experiment."
Proceedings of the Sixth International Conference
on Numerical Ship Hydrodynamics, August 1993.
OCR for page 207
[30] Thomas, P. D., and Lombard, C. K., "Geometric
conservation law and its application to flow compu-
tations on moving grids" AIAA Journal, Vol. 17,
No. 10, 1978, pp. 1030-1037.
OCR for page 208
DISCUSSION
L. J. Doctors
The University of New South Wales, Australia
Thank you for a fascinating exposition of a very
challenging and practical problem. I would like
to ask you how you handle the low-speed flow
calculations in the neighborhood of the transom
stern.
In this region, the water level on the surface of
the transom must drop in a complicated fashion
as the speed of the vessels increases, until some
point is reached where the flow will separate
clearly along the girth of the transom.
AUTHORS' REPLY
We choose a Froude number where the transom
remained wetted. Thus, our initial "ridding
efforts started with a planar free surface at the
original undisturbed water line. Then, as the
flow develops, the free surface sets up and the
grid is moved to match the free surface. Around
the corner of the transom, the free surface drops
quite rapidly. If the grid were to match the free
surface, the points on the transom would need to
be compressed substantially. In this region, we
relax the condition that the grid match the free
surface, to enable us to have a robust grid
movement/ nonlinear free surface algorithm.
If we chose a Froude number which results in a
dry transom, we would need to use a different
initial grid. If we started with a planar free
surface, as above, then we would need to
compress the points along the transom down
towards the bottom of the transom. Since we
would not be able to compress them completely
to an edge, the grid movement algorithm would
fail and we would not be able to simulate this
flow. However, by starting with a nonplanar free
surface grid, we hope to be able to simulate the
flow around a dry transom, as well as for other
similar geometries.
DISCUSSION
Jacquin Erwan
Bassin d'Essais des Carenes, France
I congratulate the authors for the quality of their
work, and for their presentation. Calculations
using unstructured mesh, with nonlinear free
surface conditions with rotating propulsors is a
hard challenge, and the qualitative results you
have shown seems very realistic. Nevertheless,
no results are published on quantitative
quantities that are available on this hull form in
the literature (pressure and frictional resistance,
boundary layer velocities, ...~. Did you make
comparisions between your results and
experimental data? Thank you for your answer.
AUTHORS' REPLY
You have clearly pointed out a weakness in our
paper, in that we do not compare with
experimental quantities such as drag or velocity
profiles. The goal of the work that was
presented was to develop efficient, robust and
stable algorithms for the nonlinear free surface
capabilities of our unstructured flow solver. In
the future, we plan on analyzing the ability of
our code to predict the forces on the hull form
and the propulsors.
Currently, our unstructured flow solver does not
capture the hook feature in the wake region
directly upstream of the propulsors. Our
structured flow solver is able to capture this
feature reasonably well. This failure is of great
concern to us and we are seeking to remedy this.
We also have concerns about the number of
nodes required to achieve a grid independent
solution. We have run the model 5415 series
hull with several different grids and have noticed
some significant differences in the flow field,
along the hull behind the bow wave and in the
stern region. Thus, we expect that we will need
to run with significantly more points that we are
currently using.
Concerning the velocity distributions, we would
not expect to match the free surface along the
hull, if the velocity distribution near the hull was
not accurate. And as we move away from the
hull, which means that the grid is becoming less
refined, we would expect that the velocity
distribution would be smeared as is the free
surface. By adding more points in these regions,
we would probably maintain the velocity
distribution better.
Finally, the free surface in the transom region is
highly dependent on Reynolds number which
suggests that the free surface is dependent on the
turbulent features in the flow. To adequately
capture the flow features in this region, we will
probably need to convect the turbulent features
OCR for page 209
upstream of the transom better, which includes
simulating the hook feature in the nominal wake.
As a predictive tool for forces and moments, we
have much work ahead of us, before we and
others will be confident of the results.
DISCUSSION
Toby J. Ratcliffe
Naval Surface Warfare Center, Carderock, USA
The authors have presented an impressive set of
results from the Unstructured Grid code
developed at MSU. Clearly, these results show
much improvement in predicting the surface
wave elevations behind the transom of Model
5415, in both the unpowered and powered
conditions, over the structured grid solutions. At
the Froude number of 0.28 for which these
comparisons are made, the transom of the model
is partially wetted. At higher speeds, the stern
becomes completely dry. This dry transom
condition has presented problems for modelling
in a "traditional" structured grid code. Would the
unstuctured RAN S code, detailed in this paper,
with the non-linear free surface capabilty, be
better suited toward solving the dry transom
condition?
AUTHORS' REPLY
We have not run cases with a dry transom stern
using either our structured or unstructured RANS
solvers. The main issue for both codes is how to
build the grid to satisfy the dry transom
condition, and for the unstructured case, this
issue is complicated with the issue of how to
terminate the viscous packing off the viscous
surfaces.
For the wetted case, the viscous packing can
terminate into the water surface without any
special considerations, as shown in the diagram
labeled Wet Transom Case. However, for the
dry transom case, if the viscous packing were to
terminate into the water surface, the grid quality
at the transom corner would be quite poor.
Other options for terminating the viscous
packing have disadvantages as well. Hence, we
would suggest wrapping the grid around the
transom for one element, so as to allow the
viscous packing to terminate cleanly into the
water surface, as shown in the diagram labeled
Dry Transom Case. By doing this, no special
considerations for either the free surface solver
or the grid movement algorithm would be
required, and the points in the free surface along
the transom stern would still be allowed to move
up and down to a small extent and would then be
a better indication of whether the numerics
suggest that the transom is truly dry.
Water Surface
,, ......... at
Transom
Hull Boy
Wetted Transom Case
r
Transom |
Hull Bo¢t~
DISCUSSION
Dry Transom C:ase
Michel Visonneau
Ecole Centrale de Nantes, France
/
J
water Surface
Concerning the 3D regrinding tool, do you use
the three dimensional generalization purposed by
Farhat recently or did you build your own
regridding tool?
AUTHORS' REPLY
In the fall of 1999, Dr. Farhat visited our
university and shared with us his work in 2D
torsional springs. We implemented his 2D
algorithm and achieved similar excellent results.
After several unsuccessful attempts to extend the
method to 3D, we developed an algorithm which
performs well in 3D. This algorithm was
developed independently from Dr. Farhat during
the winter of 2000. We are currently working on
a paper to describe this algorithm as applied to
viscous mixed element type grids but have not
OCR for page 210
published it because improvements are under
development in regards to computational
efficiency.
Representative terms from entire chapter:
nonlinear free