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 805
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 805
Computation of Nonlinear Turbulent Free Surface Flows Using the
Parallel Uncle Code
M.Beddhu, R.Pankajakshan, M.-Y.Jiang, M.Remotigue, C.Sheng, L.Taylor, W.Briley, D.Whitfield (Mississippi State
University, USA)
ABSTRACT
A numerical approach is presented in this work suitable for the computation of nonlinear free surface flows over
complex geometries such as ship hulls in a fast, reliable and robust manner. The governing equations solved are the
incompressible Reynolds Averaged Navier-Stokes (RANS) equations coupled with the free surface kinematic condition
and a two equation turbulence model. Simple no normal-gradient dynamic boundary conditions are used at the free
surface. The governing equations are cast with respect to an unsteady (non-inertial) general curvilinear coordinate system.
The numerical approach uses the modified artificial compressibility formulation. The governing equations are discretized
using a finite volume approach where the numerical fluxes at cell interfaces are obtained using Roe's inviscid flux
averages coupled with van Leer's MUSCL formulation for higher order flux extrapolation. Viscous fluxes are averaged
using central differencing. Time is discretized implicitly using the first order Euler backward differencing. The resulting
non-linear algebraic equations are solved using the discretized Newton-relaxation (DNR) approach with symmetrical
Gauss-Seidel sweeps. To speed up the solution process a parallel implementation of the numerical algorithm that uses
MPI for message passing is used. In order to accelerate the solution convergence process a multilevel approach coupled
with the traditional multigrid approach is taken. The resulting algorithm has been applied to various ship geometries and
comparisons with the sequential code solutions and experimental results are presented. The results show that the parallel
version of the free surface UNCLE code accurately reproduces these earlier results.
INTRODUCTION
Nonlinear, turbulent free surface flows represent an important class of problems with immediate naval applications
especially when such flows occur in the vicinity of a body. These problems are very challenging from a Computational
Fluid Dynamical point of view in that they demand a very robust numerical algorithm, very large computational resources
and very large amounts of computing time. More importantly, in addition to robustness, the CFD algorithm must model
the correct physics. The original UNCLE (for UNsteady Computation of fieLd Equations) code (Taylor (1991), Whitfield
and Taylor (1991)) was developed based on first principles to solve the unsteady Reynolds averaged Navier-Stokes
(unRANS) equations without any further simplifying assumptions. This sequential version was further extended in
(Beddhu et al (1994), (1999)) to include the effects of a free surface. In Beddhu et al (1999) the free surface governing
equation was formulated in terms of surface curvilinear coordinates introduced on the actual evolving free surface for the
first time. Previous efforts have introduced the surface curvilinear coordinates on a flat surface. The formulation
introduced in Beddhu et al (1999) allows for the computation of steep and breaking waves. The sequential version of the
free surface UNCLE code has been applied to various geometries with a good degree of success (see Beddhu et al (1999)
and the references therein). However, these computations took enormous computing time.
With the motivation of reducing the run-time of the sequential UNCLE code when applied for complex
configurations, a parallel version of the UNCLE code without the free surface capability was developed by Pankajakshan
and Briley (1996). It is designed to operate at coarse/medium grain parallelism levels for optimum performance and is
based on Single Program Multiple Data (SPMD) model. It uses MPI for message passing. This version represents one of
the first parallel CFD codes that has been used for solving practical problems in a routine manner. Numerous test cases
were used to test the robustness and accuracy of the code (Pankajakshan (1997)). The effort presented in the present paper
incorporates the free surface capability into the parallel version of the UNCLE code. This involved a complete re-write of
the free surface code using FORTRAN-90 which is now included as a separate module in the parallel code.
The approach adopted to compute free surface flows is to cast the governing equations with respect to a non-inertial
frame so that the free surface can always be made to coincide with a coordinate surface. Thus, Navier-Stokes equations
are cast with respect to a set of general unsteady curvilinear coordinates and are solved
the authoritative version for attribution.
OCR for page 806
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 806
along with the non-linear free surface kinematic boundary condition. The modified artificial compressibility method
Beddhu et al (1994) is used to march the solution in time. Approximate inviscid dynamic boundary conditions are used at
the free surface since the wind stresses and surface tension are neglected. A two equation turbulence model is used for
closing the momentum equations. Characteristic Variable Boundary Conditions (CVBC) are used on the inflow and
outflow boundaries.
The numerical scheme in Taylor (1991) and Whitfield and Taylor (1991) uses a finite volume formulation and uses
the Roe scheme for obtaining the first order fluxes and the MUSCL scheme of van Leer to obtain higher order fluxes at
the cell faces. The flux Jacobians which are required in an implicit scheme are obtained numerically. The viscous flux
Jacobians are obtained using the thin layer approximation. The resulting system of algebraic equations are solved using
the Discretized-Newton Relaxation (DNR) scheme which is comprised of Newton iterations with symmetric Gauss-Seidel
sub-iterations in each Newton iteration. A Full Approximation Scheme (FAS) multigrid algorithm is used to accelerate
the convergence. The free surface kinematic boundary condition is solved using the same approach.
Since the present approach falls under the category of front tracking methods, it is necessary to regenerate the free
surface and the underlying grid at each time step. The free surface is assumed to be of the form y=Y(ξ, ζ, τ) where ξ and ζ
are surface curvilinear coordinates and τ denotes time. Thus, in solving the kinematic boundary condition one obtains an
increment ∆y at every point on the free surface for an increment ∆τ in time. Using the ∆y's the new free surface position
is obtained and the flow grid is updated using a background grid.
The parallel code was verified by running the same cases as the sequential code and the results were compared with
experiments. The majority of the results shown are for Series-60 CB=0.6 and Model 5415. The flow parameters for
Series-60 CB=0.6 were specified by the experiments to be as follows: a Froude number of 0.316 and a Reynolds number
of 4.2×106. For Model 5415, the Froude number was given as 0.2756 and the Reynolds number as 12.02×106.
Some of the precursory works in free surface flow computations about surface ships were performed by Miyata and
Nishimura (1985), Kodama (1989), Hino (1989) and Farmer et al (1993). Recent works are that of Tahara and Stern
(1996) and Beddhu et al (1994), (1998), (1999). Tahara and Stern (1996) use a Poisson equation for pressure and use the
pressure-implicit split-operator (PISO) algorithm to solve the mean flow equations and the Beam and Warming approach
with artificial dissipation for solving the free surface equation. They use the so-called finite analytic method to discretize
the governing equations. In addition, they use separate grids for solving the mean flow and the free surface and use
interpolation for transferring data between the grids. Beddhu et al (1994) introduced the modified artificial
compressibility method for solving the mean flow equations and use an explicit method for calculating the free surface
motion. For the computation of mean flow they use a finite volume, implicit scheme patterning their numerical scheme
after compressible flow solvers. Thus, Roe scheme is used for obtaining first order numerical fluxes and van Leer's
MUSCL approach is used to obtain third order corrections. The viscous terms are approximated using second order
central differences and the time derivative is approximated using either first order or second order backward differences.
In contrast to Beddhu et al (1994), Beddhu et al (1999) use an implicit method for calculating the free surface. The
method was chosen to be the same as the one used for the mean flow. In addition, a novel way of tracking the free surface
in time was introduced which preserves the shape of the free surface at each time step during various grid operations.
Since the geometries of actual ships are very complex such an approach is necessary for advancing the free surface along
curvilinear coordinates. In essence, a background grid was used which is fixed in time. Grid points on the free surface are
allowed to move along a particular family of coordinate lines, chosen a priori, of the background grid. The portion of the
coordinate lines below the free surface is then used to rebuild the actual grid at the next time level using the arclength
distribution of the actual grid line at the current time level.
GOVERNING EQUATIONS
The governing equations in this work are cast with respect to a non-inertial frame. In tensor invariant form the
continuity equation in the modified artificial compressibility method (Beddhu et al (1994), (1999)) is given by
(1)
where β is the artificial compressibility parameter. The momentum equation for viscous, incompressible flows in a
non-inertial frame of reference, in a gravitational field, in non-dimensional, vector invariant form is given by (Beddhu et
al (1994), (1999))
(2)
the authoritative version for attribution.
OCR for page 807
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 807
In Eqs. (1) and (2), u=u*/U0, is the non-dimensional absolute velocity vector, v=u+w is the non-dimensional relative
velocity vector, w is the non-dimensional grid velocity vector, τ=tU0/L, is the non-dimensional time, is the Stokes
tensor and p′=p+χ/Fr2 where is the non-dimensional pressure and χ is the body force potential due to
gravity. Re0 is the Reynolds number, Re0=ρ0U0L/µ0, where, ρ0 is a reference density, U0 is a reference velocity, L is a
reference length, and, µ0 is a reference coefficient of viscosity. Fr is the Froude number given by where, a is
the acceleration due to gravity. Note that when the positive y-direction is aligned in the direction opposite to the gravity
vector, one obtains χ=y. A tilde over a quantity denotes that it is a tensor and an underscore denotes that it is a vector. The
Stokes tensor is given by
Fig. 1. Schematic illustrating the Cartesian and curvilinear coordinates chosen.
(3)
where, µ=µ*/µ0, is the non-dimensional coefficient of viscosity. The superscript ‘T' in Eq. (3) denotes the transpose
operation.
Casting the governing equations (1) and (2) with respect to a curvilinear coordinate system (ξ, η, ς, τ) and using the
so-called partial transformation in which the vectors and tensors that appear within the divergence operator are expressed
with respect to the underlying Cartesian coordinates whereas the divergence operator itself is expressed in terms of the
curvilinear coordinates, one can obtain the so-called numerical vector form which is as follows:
(4)
where,
the authoritative version for attribution.
u, v, and w, are the components of the absolute velocity vector with respect to a Cartesian coordinate system, σxx,
etc., are the Cartesian components of the Stokes tensor, ξx, ξy and ξz, are the Cartesian components of the contravariant
base vector grad ξ. Expressions for G and H are similar to F and can be obtained from F by replacing ξ by η and ζ
respectively. Similarly GV and HV can be obtained from FV. Figure 1 indicates the coordinate systems chosen.
Equation (4) is the system of non-linear equations to be solved numerically, using a set of physical and numerical
boundary conditions. The physical boundary conditions include the noslip condition
OCR for page 808
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 808
(v=0 ⇒ u=−w) on the body and the kinematic and dynamic free surface boundary conditions.
The kinematic condition used in this work can be derived as follows. Assume that the free surface is represented by
the function Then the kinematic condition implies (Warsi (1993))
(5)
In terms of a set of inertial curvilinear coordinates Eq. (5) can be written as
(6)
u1 u2 u3
where the contravariant velocity component is defined below Eq. (4) and and can be similarly defined. On
the other hand, in terms of a set of non-inertial curvilinear coordinates Eq. (5) can be
written as
(7)
where w1, w2 and w3 are the contravariant components of the grid velocity vector. Choosing η=y, Eq. (7) becomes
(8)
Up to this point the notion of surface curvilinear coordinates is not needed. Note that there is no restriction placed on
the coordinates ξ and ς. In general, ξ=ξ(x, y, z, t) and ς=ς(x, y, z, t). Now, let the surface be represented by
(9)
Substituting, Eq. (9) in Eq. (8), one obtains
(10)
In Eqs. (9) and (10) the curvilinear coordinates ξ and ς need to be interpreted as surface curvilinear coordinates.
However, ξ and ς are still of the form ξ=ξ(x, y, z, t) and ς=ς(x, y, z, t). The velocity components u1 and u3, and, the grid
velocity components w1 and w3 are now surface contravariant components.
Note that Eq. (10) is identical inform to the one obtained by Farmer et al (1993). However, they have explicitly
assumed that ξ=ξ(x, z) and ς=ς(x, z). The formulation introduced here is completely general and is valid as long as Y
remains a single valued function of ξ and ς. Thus, the formulation introduced here is valid for tracking breaking waves up
to the point of reentry (see Fig. 2). The conceptual difference between the present formulation and that of Farmer et al
(1993) is that Farmer et al (1993) introduce the curvilinear coordinates on a flat surface (i.e. xz-plane) whereas the present
formulation introduces them on the actual free surface.
Fig. 2. Breaking wave representation.
Y(x) is multi—valued for x0≤x≤x1.
Y(ξ) is however single—valued until reentry.
Any book on differential geometry, Warsi (1998) for example, can be consulted for obtaining the metrics of the
surface curvilinear coordinates. From these metrics the surface contravariant components of the flow velocity and the grid
velocity can easily be obtained as outlined in Beddhu et al (1999).
Since, Eq. (10) is cast in terms of curvilinear coordinates, the numerical scheme for solving it can be patterned after
that of Eq. (4).
Two two-equation turbulence models are available in the UNCLE code. They are the modified Shih and Lumley
model (Shih and Lumley (1993), Yang et al (1995), Liou and Shih (1996)) and the q-ω model (Coakley, 1983). The
governing equations of the model are:
(11)
the authoritative version for attribution.
(12)
OCR for page 809
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 809
where C1=1.44, C2=1.92, In this formulation the eddy viscosity is
defined as usual as
(13)
However, fµ and Cµ are defined as
The governing equations of the q-ω model are:
(14)
(15)
In this formulation, and the eddy viscosity is defined as before as
(16)
where Cµ=0.9 and fµ=1−exp(−0.0065qy/v). The other constants and functions are defined as C2= 0.92, C1=0.045
+0.405fµ, and σω=1.3.
In the absence of surface tension, continuity of the stress vector across the interface is the exact dynamic free surface
boundary condition that one needs to impose. This condition was originally obtained by Hirt and Shannon (1968). An
efficient way of implementing these exact conditions is outlined in Beddhu et al (1997). However, in the present study
exact dynamic boundary conditions are not used. Since the external tangential stresses are neglected and only the
atmospheric pressure is considered in the normal stress component and also since the grid near the free surface is not fine
enough to resolve the weak surface layer, the dynamic boundary condition has been approximately implemented by many
authors for computing ship related flows. For example, commonly used approximate conditions are
where, z is the direction opposite to the direction of the gravity vector. In the present work, a characteristic variable
based approach is used which can be outlined as follows. First, every term in Eq. (4) is neglected except the unsteady
term and the inviscid term in the η-direction. The resultant equation is linearized and cast in a diagonal form as outlined in
Taylor (1991). Thus, one obtains
(17)
is the left eigenvector of the flux Jacobian ∂G/∂Q and is
where, is the characteristic variable vector,
the diagonal matrix containing the eigenvalues of the flux Jacobian ∂G/∂Q. The subscript ‘0' in denotes that it is
treated as a constant matrix. The eigenvalues are given by
(18)
where,
On a η=ηmax boundary which is typically chosen as the free surface, λ1 and λ2 are zero, λ3 is positive and λ4 is
the authoritative version for attribution.
negative. Thus the characteristic variable W4 needs to be prescribed and W3 needs to be extrapolated from within the
computational domain. However, since λ1 and λ2 are zero one can either extrapolate or specify W1 and W2. In this work,
W2 and W3 are extrapolated from within the computational domain, and instead of prescribing W1 and W4 one uses the
conditions ηt+uηx+vηy+wηz=0 and p′=y/Fr2. A 3x3 matrix is solved for the velocity components at each grid point using
W2, W3 and the kinematic condition. Even though, this is an inviscid approximation to the exact viscous dynamic
boundary condition, it works quite well for the cases considered and very good agreement has been obtained with
measured wave profiles.
The numerical boundary conditions are imposed on artificial (outer) boundaries which are
OCR for page 810
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 810
introduced to truncate the computational domain to a finite size so that the resulting problem can be solved using a
computer. The assumption is that the artificial boundaries are far removed from the physical body that they won't affect
the accuracy of the solution in the vicinity of the body. In this work the fluid flow is assumed to be along the positive x
direction with the body fixed at the origin. Thus, far upstream of the body one uses a characteristic variable based inflow
boundary condition and a characteristic variable based outflow boundary condition far downstream of the body. The up
stream and downstream boundaries are located at least at 5 body lengths each from the origin. At the far away side
boundary and the bottom boundary either characteristic variable based inflow or outflow boundary condition is used
depending upon the local velocity vector. The side boundary and the bottom boundaries are located at 5 body lengths
from the origin. At the z=0 boundary, flow symmetry is imposed. In addition, the computational domain is subdivided
into arbitrary sub-domains and at the boundaries of these domains a two point symmetry condition is used.
The following procedure is used to advance the solution from time step n to time step n+1:
1) Solve for the interior of the flow field using the Discretized Newton Relaxation method.
2) Update the boundary conditions on the flow variables on all surfaces.
3) Use the kinematic condition to find the new position of the free surface.
4) Find the intersection of the background η lines with the updated free surface using a nearest point
approximation.
5) Recreate the volume grid from the free surface obtained in step 4.
6) Obtain the new metric coefficients and grid speeds.
7) Update the free surface dynamic boundary condition once again. This is done since the shape of the free
surface has changed (in step 2 the free surface at the previous time level n was used) and is found to enhance
the stability of the scheme.
8) go to step 1.
An integral quantity of immense practical interest is the total resistance of a ship. The total resistance is calculated
using the following expression:
where
In the actual computation, the integral is replaced by summation over all the surface grid cells. nx is the component
of the unit normal in the x direction at any point on the body. The quantity |∇ς|dξdη is the surface elemental area where
Note that the body is represented by a ς=constant surface everywhere except the transom surface
which is a ξ=constant surface. This fact is properly accounted for in the results presented.
NUMERICAL PROCEDURE
The numerical scheme used in this study is similar to that proposed by Pan and Chakravarthy (1989) and is discussed
in detail by Taylor (1991), and, Whitfield and Taylor (1991). An extensive discussion of the methodology has been
presented by Whitfield and Taylor (1994) applicable to two dimensional flows. The approach taken in this work is to
solve Eq (4) implicitly using the Discretized Newton-Relaxation (DNR) scheme (Ortega and Rheinboldt (1970),
Whitfield and Taylor (1991)), where the fluxes at the cell faces are obtained using the Roe scheme (1981) with higher
order accuracy achieved using the MUSCL approach (van Leer (1979); Whitfield and Taylor (1991)). Writing Eq. (4) in
discrete form,
(19)
where and so on. Note that for a higher order flux representation depends on
the authoritative version for attribution.
as well. If Eq. (19) is expanded for each grid cell, a system of algebraic equations are obtained in terms of qn
and
+1 at each grid cell where Strictly speaking Fn+1 is a function of both qn+1 and the metrics at n+1. Since
the metrics at n+1 are known, no linearization needs to be done with respect to the metrics. Hence
OCR for page 811
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 811
Eq. (19) is regarded as a function of qn+1 alone. In functional form, Eq. (19) can be represented as
(20)
Solving Eq. (20) involves finding the roots of a system of non-linear algebraic equations. Using Newton's method,
the solutions of Eq. (20) are obtained from the following linear equations:
Fig. 3 Schematic of the Parallel Iteration Hierarchy
(21)
In order to limit the band width of the matrix, the operator (∂X/∂q) is obtained using higher order fluxes in a special
manner and is rearranged along the lines of Whitfield and Taylor (1994), into a strong diagonal form. The viscous flux
Jacobians are obtained using the thin layer approximation whereas the residue X(qn+1, m) contains all the viscous terms.
Within each Newton iteration symmetric Gauss-Seidel passes are used. The resulting algorithm is termed the Discretized
Newton-Relaxation (DNR) procedure. When the iteration in m converges, qn+1 is obtained and the calculation procedure
is extended to the next time level. As the iteration in m converges, the LHS of Eq (21) goes to zero. Hence time accuracy
is introduced into the scheme by multiplying the local time derivative term in X(qn+1, m) by a conditioning matrix Ia where
Ia= diag (0, 1, 1, 1). The inviscid fluxes on the RHS of Eq (21) are obtained using a third order MUSCL-type flux and the
viscous fluxes using a second order central differencing.
PARALLEL ITERATION HIERARCHY
The parallel implicit algorithm used here was developed by Pankajakshan and Briley (1996) and Pankajakshan
the authoritative version for attribution.
(1997) in terms of an iteration hierarchy in which the Newton, multigrid and relaxation schemes are re-examined in the
context of a domain decomposition that exploits coarse/medium grain parallelism. The basic iteration hierarchy was
extended with the addition of multilevel acceleration (Whitfield and Taylor (1998)) which is used with an embedded
multigrid scheme at the finer grid levels. A schematic of the parallel iteration hierarchy is shown in Figure 3. The
linearized approximations at each time step are solved using a block Jacobi symmetric Gauss-Seidel (BJ-SGS) relaxation
procedure that involves communication at each sub-iteration level. The sequencing of forward and backward sweeps
interspersed with data exchange is shown in Figure 4.
The code uses MPI for all communication primarily through point-to-point messages in asynchronous non-blocking
mode. There are separate message-passing routines for the solver, two-equation model and the free surface module.
Parallel runs typically have communications overheads of 10–15%. Scalability studies using heuristic performance
estimates that account for algorithm characteristics and system resources indicate that parallel efficiencies of 80% can be
OCR for page 812
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 812
achieved on up to 500 processors using appropriately sized grids of up to 50 million points.
MULTIGRID STRATEGY
The multigrid approach (Sheng et al (1995)) used to accelerate the iterative convergence at each multilevel level is
described schematically in Figure 5. It is seen that a time-dependent solution is obtained through the use of Newton sub-
iterations, which are the most costly part of the overall work. However, this effort can be reduced with the aid of an
efficient multigrid method (FAS). The unsteady multigrid method is basically the same as the steady, except that both the
fine grid and coarse grid equations must be solved at the same time level to ensure temporal consistency, while in the
steady multigrid approach, time is advanced in the fine grid as well as the coarse grid to achieve full efficiency. The two-
level multigrid method for a general unsteady equation can be briefly described as follows:
1. Iterate Nh(Qh)=0 times on the fine grid h by Newton's method.
2. Restrict the residual and solution to the coarser grid 2h, and iterate N2h(Q2h)=τ2h times, where τ2h=N2h(I2hhQh)
2h (NhQh) is the relative truncation error between the grids h and 2h.
−R h
3. Interpolate the correction from the coarser grid to the fine grid and update the solution
Fig. 5. Schematic of the Multigrid Cycle
Fig. 4. Sequence of Operations for One BJ-SGS/1 Sweep
in Asynchronous Mode
times at the same time level, using Qh as the new approximation to Qn+1.
4. Repeat steps 1~3 for
In the above procedure, is the number of Newton sub-iterations for the fine grid and coarser grids, and is
the number of multigrid cycles implemented at each time step. Choosing different values of and may form
different multigrid strategies and result in different effects. However, an important fact is that the cost of CPU time is
proportional to the multigrid cycles ( ) at each time step.
GRID GENERATION
The grids presented in this work are generated using Graphical Unstructured Multi-Block (GUMB) structured grid
system (Jiang (2000), Remotigue (1999), Jiang and Remotigue (1998)) which is being developed in house. GUMB
retained and further enhanced both the initial General Topology Model (GTM) data structure and geometry engine based
on Non-Uniform Rational B-Splines (NURBS) which were originally developed for NGP (Thompson (1992), Remotigue
(1994)). They are both coupled with a structured grid
the authoritative version for attribution.
OCR for page 813
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 813
generation library. A schematic of the organization of the libraries of GUMB is shown in Figure 6.
Fig 6. GUM-B Schematic of Library Layout
GUMBO, Graphical Unstructured Multi-Block Omnitool, is a graphical user interface, developed at Mississippi
State University, that has flexible and general repartition and manipulation algorithms that automatically account for the
boundary conditions and connectivity. Furthermore, the application of boundary conditions and connectivity is simple and
straight forward. In addition, grid assurance and quality measures are included to validate the grid, boundary conditions,
and connectivity. The schematic of GUMBO can be seen in Figure 7.
Fig. 7. GUMBO Schematic of Library Layout
For all the cases presented, load balancing through domain decomposition, specification of boundary conditions and
automatic detection of inter-block boundaries were done using GUMBO.
For nonlinear free surface calculations it becomes necessary in the present method that the flow grid be adapted to
the free surface such that the top surface of the flow grid coincides with the free surface. Beddhu et al (1998b) took an
intersection approach where the every η-line of the background grid is intersected with the actual free surface to
redistribute points on the free surface in order to generate the flow grid for the next time step. Thus, this approach
preserves the shape of the free surface during the grid regeneration process and is suitable for unsteady free surface flows.
Since the present work is focussed on steady free surface flows, the following approach is taken to redistribute points on
the free surface: For each point (i, k) on the free surface find the closest point on the corresponding η-line of the
background grid. Note that this approach distorts the free surface. Since only steady state is of interest, one expects the
free surface to move little as the iteration count increases and thus one expects the free surface to reach the same shape as
it would have reached via the intersection approach. However, in practice, though very little free surface movement is
noticed at large iteration counts, the final shape is slightly inferior to the one obtained through the intersection approach.
Further testing and an intersection approach based on NURBS are underway.
Fig. 8 Hull Profile Comparisons for Wigley Hull
the authoritative version for attribution.
Fr=0.289; Re=3300000
● Expt ——— Sequential — - — Parallel
RESULTS
The free surface version of the parallel UNCLE code described in the earlier sections was used to compute the free
surface flow fields around three popular hulls: (1) Wigley, (2) Series 60 CB=0.6 and (3) Model 5415. These results are
presented below.
Wigley hull is a fairly simple geometry and has an analytical description to it. Referring to the coordinate system in
Fig. 1, the equation for the hull is given
OCR for page 814
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 814
by: where B is the breadth, L is the length and D is the draft of the ship. The following
ratios were taken L/B=10 and L/D=16. The flow conditions are Fr=0.289 and Re=3.3×106. Computed results using the
sequential UNCLE code for the Wigley hull were reported in Beddhu et al (1998a). Hull profile for Wigley hull is shown
in Fig. 8. Note that the parallel code produces a profile that is in very good agreement with the sequential code and the
experiment.
Fig. 9 Comparison of Hull Wave Profiles for Series 60 CB=0.6
Fr=0.316; Re=4.02×106
● Expt ——— Sequential — - — Parallel
Fig. 10 Transverse Wavecut Comparison in the Bow
Region for Series 60 CB=0.6
Fr=0.316; Re=4.02×106
● Expt ——— Sequential — - — Parallel Fig. 11 Transverse Wavecut Comparison in the Stern
Region for Series 60 CB=0.6
Fr=0.316; Re=4.02×106
● Expt ——— Sequential — - — Parallel
Fig. 12 Comparison of Experimental and Computed (Parallel UNCLE) Wave Contours for Series 60CB=0.6
the authoritative version for attribution.
Fr=0.316; Re=4.02×106
Series 60 CB=0.6 is a geometry for which a lot of experimental data is readily available (Longo et al (1993)). The
flow conditions chosen for running the parallel UNCLE code were as follows: Fr=0.316 and Re=4.02×106. Results using
the sequential version
OCR for page 815
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 815
of the UNCLE code were reported in Nichols (1998) and Beddhu et al (1998c). Comparison of the hull profile between
the sequential and parallel versions with the experiment is shown in Fig 9. It can be seen that the prediction of the parallel
version is similar to the sequential version. This similarity can be further seen in the transverse wavecuts presented in Figs
10 and 11. In Fig. 12, comparison of the overall wave contours between experiment and computation using the parallel
UNCLE code is shown. This figure shows that there is reasonable overall agreement. In Fig. 13 the u-velocity contours on
the propeller plane of the parallel code is compared with the experiment. Here, it can be seen that the agreement is quite
good. A similar agreement for the sequential code is shown in Beddhu et al (1998c).
Fig. 13 Comparison of Experimental and Computed (Parallel UNCLE) u-velocity contours on the Propeller Plane for Series
60 CB=0.6
Fr=0.316; Re=4.02×106
Fig. 14 Comparison of Hull Wave Profiles for Model 5415
Fr=0.2756; Re=12.02×106
Fig. 15 Comparison of Longitudinal Wavecuts at z=0.0965 for Model 5415
Fr=0.2756; Re=12.02×106
● Expt ——— Sequential — - — Parallel
the authoritative version for attribution.
Fig. 16 Comparison of Wave Contours for Model 5415
Fr=0.2756; Re=12.02×106
The next test case to be presented is the Model 5415. Extensive experimental results for this model is reported in
Ratcliffe and Lindenmuth (1990) and Oliv
OCR for page 816
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 816
eri et al. The flow conditions are Fr=0.2756 and Re=12.02×106. The computed wave profile has been plotted against the
available experimental wave profiles in Fig. 14. Note that the experimental wave profiles themselves do not fall in each
other's uncertainty limits. Thus, the implication of such disagreements between experimental results to the efforts for
verification and validation of numerical results is not clear. However, the computed profiles agrees reasonably well with
either measurement. In Fig. 15, the comparison of the longitudinal wavecut at z=0.0965 between experiment and the two
versions of the UNCLE code is shown. It can be seen that the agreement is quite reasonable. The overall comparison of
the computed and experimental wave contours shown in Fig. 16 also shows good agreement. Fig 17 compares the
contours of u-component of velocity between the parallel version and the experiment. It can be seen that the agreement is
excellent. Computed and measured total resistance comparisons are shown in Fig. 18. It can be seen that though the
sequential and parallel versions follow different paths of convergence they both converge to the experimental value in
about 4000 cycles. The difference in the paths of convergence is due to algorithmic differences between the two versions.
For example the sequential version has no multilevel capability. Finally, the computed and measured stern wave patterns
are compared in Fig. 19. Since the transom is wet at this Froude number, computing and matching the stern wave system
is supposed to be a tough challenge. It can be seen that the parallel version of the UNCLE code does a very good job of
predicting this tough flow feature. A similar comparison with the sequential code is shown in Beddhu et al (1998b, 1999).
Fig. 17 Comparison of Contours of U-Velocity Component at the Propeller Plane for Model 5415
Fr=0.2756; Re=12.02×106
Fig. 18 Resistance Comparisons for Model 5415 Fr=0.2756; Re=12.02×106
Fig. 19 Comparison of Stern Wave Contours for Model 5415
Fr=0.2756; Re=12.02×106
• Expt ——— Sequential — - — Parallel
CONCLUSION
The free surface version of the parallel UNCLE code was designed to be a production mode code that is fast, robust
the authoritative version for attribution.
and reliable. It has been tested against various geometries and has produced similar agreements with experiments as the
free surface version of the sequential UNCLE code. A NURBS based intersection algorithm is being developed that will
improve the accuracy of the computation of free surfaces. This code has been successfully transitioned to DTMB and has
undergone extensive user trials. Currently, efforts are being focussed in solving unsteady free surface flows.
ACKNOWLEDGEMENTS
This work was supported by grant N00014–97–1–0959 from the Office of Naval Research. The grant monitor is Dr.
Edwin Rood. This support is greatly appreciated.
REFERENCES
Beddhu, M., Taylor, L.K. and Whitfield, D.L., “ATime Accurate Calculation Procedure for Flows with a Free Surface Using a Modified Artificial
Compressibility Formulation,” Applied Mathematics and Computation, Vol. 65, 1994, pp. 33–48.
OCR for page 817
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 817
Beddhu, M., Jiang, M.Y., Whitfield, D.L, Taylor, L.K., and Arabshahi, A. “Computational Physical Oceanography—A Comprehensive Approach based
on Generalized CFD/Grid Techniques for Planetary Scale Simulations of Oceanic Flows,” MSSU-EIRS-ERC-97–5, Mississippi State, MS
39762, Feb., 1997.
Beddhu, M., Jiang, M.Y., Taylor, L.K. and Whitfield, D.L., “Computation of Steady and Unsteady Flows with a Free Surface Around the Wigley Hull,”
Applied Mathematics and Computation, Vol. 89, 1998a, pp. 67–84.
Beddhu, M., Jiang, M.Y., Whitfield, D.L.Taylor, L.K. and Arabshahi, A., “CFD Validation of the Free Surface Flow Around DTMB Model 5415 Using
Reynolds Averaged Navier-Stokes Equations,” Proceedings of the Third Osaka Colloquium on Advanced CFD applications to Ship Flow and
Hull Form Design, 1998b.
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, 1998c.
Beddhu, M., Jiang, M.Y., Whitfield, D.L., and Taylor, L.K., “Computation of the Wetted Transom Stern Flow over Model 5415,” Proceedings of the
Seventh International Conference on Numerical Ship Hydrodynamics, Nantes, France, July 19–22, 1999.
Coakley, T.J., “Turbulence Modelling Methods for the Compressible Navier-Stokes Equations,” AIAA 16th Fluid and Plasma Dynamics Conference,
July 12–14, 1993, Danvers, Massachusetts.
Farmer, J., Martinelli L., and James on A., “A Fast Multigrid Method for Solving Incompressible Hydrodynamic Problems With Free Surfaces,” AIAA
Journal, Vol. 32, No. 6, 1993, pp. 1175–1182.
Hino, T., “Computation of Free Surface Flow Around an Advancing Ship by the Navier-Stokes Equations,” In Proceedings, Fifth International
Conference on Numerical Ship Hydrodynamics, 1989, pp. 103–117.
Hirt, C.W. and Shannon, J.P., “Free-Surface Stress Conditions for Incompressible-Flow Calculations,”, Journal of Computational Physics, Vol. 2, 1968,
pp. 403–411.
Jiang, M. and Remotigue, M., “GUM-B Grid Generation Code and Applications,” Proceedings of the 6th International Conference in Numerical Grid
Generation in Computational Field Simulations, London, England, July 1998, pp. 823–832.
Jiang, M., “GUM-B Applications to Submarine and Surface Vessels,” Proceedings of the 7th International Conference in Numerical Grid Generation in
Computational Field Simulations, Whistler, British Columbia, Canada, September 2000.
Liou, W.W. and Shih, T.H., “Transonic Turbulent Flow Predictions with New Two-Equation Turbulence Models,” NASA CR 198444, Jan., 1996.
Longo, J., Stern, F., and Toda, Y., “Mean flow Measurements in the Boundary Layer and Wake and Wave Field of a Series 60 CB=0.6 Ship Model—
Part 2: Scale Effects on the Near Field Wave Patterns and Comparisons with Inviscid Theory,” Journal of Ship Research, Vol. 37. No 1., 1993,
pp. 16–24.
Kodama, Y., “Grid Generation and Flow Computation for Practical Ship Hull forms and Propellers Using the Geometrical Method and the IAF
Scheme,” Proceedings, Fifth International Conference on Numerical Ship Hydrodynamics, 1989, pp. 71–85.
Miyata, H. and Nishimura, S., “Finite Difference Simulation of Nonlinear Ship Waves,” Journal of Fluid Mechanics, Vol. 157, 1985, pp. 327–357.
Nichols, III, D.S., “Calculation of Free Surface Wave Forms and Flow Field about the Series 60 CB=0.6 Ship,” Masters Thesis, Department of
Aerospace Engineering, Mississippi State University, Mississippi State, May 1998.
Ortega, J.M. and Rheinboldt, W.C., “Iterative Solution of Nonlinear Equations in Several Variables,” Academic Press Inc., New York, 1970.
Pan, D. and Chakravarthy, S., “Unified Formulation for Incompressible Flows,” AIAA-89–0122, Jan., 1989.
Pankajakshan, R. and Briley, W.R., “Parallel Solution of Viscous Incompressible Flow on Multi-Block Structured Grids Using MPI”, Parallel
Computational Fluid Dynamics—Implementations and Results using Parallel Computers, Eds: S.Taylor, A.Ecer, J.Periaux and N. Satafuca,
Elsvier Science, B.V., Amsterdam, 1996, pp. 601–608.
Pankajakshan, R., “Parallel Solution of Unsteady Incompressible Flow Using Multi-Block Structured Grids,” Ph. D. Dissertation, Mississippi State
University, Dec. 1997.
Pankajakshan, R., Taylor, L.K., Jiang, M.Y., Remotigue, M.G., Briley, W.R. and Whitfield, D.L., “Parallel Simulations for Control-Surface Induced
Submarine Maneuvers,” AIAA Paper 2000–0962, 2000.
Olivieri, A., Palini, M., and Penna, R., “Studio sperimentale del campo fluidodinamico intorno ad una care
the authoritative version for attribution.
OCR for page 818
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE 818
na dislocante veloce. Parte II: rilievo della formazione ondosa a Fr=0.28 e Fr=0.41,” Rapporto Tecnico INSEAN (in Italian).
Remotigue, M. and the NGP Team, “The National Grid Project: Making Dreams into Reality,” Proceedings of the 4th International Conference on
Numerical Grid Generation in Computational Fluid Dynamics and Related Fields, Swansea, Wales, 1994, pp. 429–440.
Ratcliffe, Toby J., and Lindenmuth, William T., “Kelvin Wake Measurements Obtained on Five Surface Ship Models,” DTRC-89/038, January 1990.
Remotigue, M.G., “Structured Grid Technology to Enable Flow Simulation in and Integrated System Environment,” Ph.D. Dissertation, Mississippi
State University, December 1999.
Roe, P.L., “Approximate Riemann Solvers, Parameter Vectors, and Difference Schemes,” Journal of Computational Physics, Vol. 43, 1981, pp. 357–372.
Sheng, C., Taylor, L.K. and Whitfield, D.L., “Multigrid Algorithm for Three-Dimensional Incompressible High Reynolds Number Turbulent Flows,”
AIAA Journal Vol. 33, No. 11, November, 1995, pp. 2073–2079.
Shih, T.H. and Lumley, J.L., “Kolmogorov Behavior of Near-Wall Turbulence and its Application in Turbulence Modelling,” Comp. Fluid Dyn.. Vol 1,
1993, pp. 43–56.
Tahara, Y., and Stern, F. “A Large-Domain Approach for Calculating Ship Boundary Layers and Wakes and Wave Fields for Nonzero Froude Number,”
Journal of Computational Physics, Vol. 127, 1996, pp. 398–411.
Taylor, L.K. “Unsteady Three—Dimensional Incompressible Algorithm Based on Artificial Compressibility,” Ph. D. Dissertation, Department of
Aerospace Engineering, Mississippi State University, Mississippi State, May 1991.
Thompson, J.F., “The National Grid Project,” Computing Systems in Engineering, Vol. 3, Nos 1–4, pp. 393–399, 1992.
van Leer, B., “Towards the Ultimate Conservation Difference Scheme V, A Second Order Sequel to Gudunov's Method,” Journal of Computational
Physics, Vol. 32, 1979, pp. 101–136.
Whitfield, D.L. and Taylor, L.K., “Discretized Newton-Relaxation Solution of High Resolution Flux-Difference Split Schemes,” AIAA-91–1539, June
1991.
Whitfield, D.L. and Taylor, L.K., “Numerical Solution of the Two-Dimensional Time-Dependent Incompressible Euler Equations,” MSSU-EIRS-
ERC-93–14, Mississippi State University, Mississippi MS 39762 April, 1994.
Whitfield, D.L. and Taylor, L.K., “Variants of a Two-Level Method for the Approximate Numerical Solution of Field Simulation Equations,” MSSU-
EIRS-ERC-98–09, Mississippi State University, Mississippi MS 39762 July, 1998.
Warsi, Z.U.A., “Fluid Dynamics: Theoretical and Computational Approaches,” 1st. ed., CRC Press, Boca Raton, p. 10, 1993.
Warsi, Z.U.A., “Ruminations on Applied Tensor Analysis,” MSSU-EIRS-ASE-98–1, Mississippi State, MS 39762, Oct. 1998, p. 78.
Yang, Z., Georgiadis, N., Zhu, J. and Shih, T.H., “Calculations of Inlet/Nozzle Flows Using a New k- model,” AIAA Paper 95–2761, 1995.
the authoritative version for attribution.
OCR for page 819
About this PDF file: This new digital representation of the original work has been recomposed from XML files created from the original paper book, not from the original typesetting files. Page breaks are true to the original; line
lengths, word breaks, heading styles, and other typesetting-specific formatting, however, cannot be retained, and some typographic errors may have been accidentally inserted. Please use the print version of this publication as
the authoritative version for attribution.
S.Cordier
None Received.
Bassin d'Essais des Carenes, France
DISCUSSION
AUTHOR'S REPLY
for example). Can you tell us what your plans are to address this limitation?
COMPUTATION OF NONLINEAR TURBULENT FREE SURFACE FLOWS USING THE PARALLEL UNCLE CODE
819
The use of grids fitted to the free surface has limits when simulating non-linear free surface flows (braking and jets
Representative terms from entire chapter:
free surfae
