| 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 882
Combining A ccuracy and E ffciency with R obustness
in S hip S tern Flow Computation
A. van der Ploeg (Maritrme Research nstitute, T h e N etherlands)
L. Eca (Instituto Superior Tecnico, Portugal)
M. Hoekstra (Maritime Research nstitute, T h e Netherlands)
Abstract
Usmg c matme RANS-solver for ship tem flows
as c starting pomt, we consider some possibilities to
erJvmce eflflciency md robustness wifhout sarff cing
acuray Among fhem me alterm~tive domcm de-
compositions md vector izab le preconditi onmg tech-
niques It is show fnat cgoodbahmcebetw enfl:x-
~bility, robust ess md efflciency c m be obtained md
fnat it is ju tffed to negdect dfff sion in mcin sheam
di ection md some of fhe cor~ctive terms
1 Introducidon
For the use of iscous flow computations in fhe
desigm process of ships, robust ess of the m merical
solution pmcedure is c pre~equisite Praticclly this
me ms that mdemecsoruible dema ds on g id qua lity,
convergence of the iterctive sohtion is warrmted in
fl w simoktions mder t pical ship opercti g condi-
tions How wv, robust ess is m sewval methods ob-
tained by sacrffcing a uray (for example by using
fl st -or der di scret isati on) A d mger ous route , be mse
fhe ~esults of c robust code have no pratical value ff
not c certam level of acuray c m be ahieved it is
one recson why considerale mphcsis is p t nowc-
days on veriflcati on md val idat ion Mm e over, in view
of usual time conshamts in ship desigm, it is desirale
fnat c good quality sohtion ca be obtamed et low
cost m c small am o mt of time (eflflciency) his pcper
cim s to sh w how c good bahmce betwen robustness,
acuray md eflflciency c m be obtamed in fhe compu-
tation of the flow aro md c ships stern by solving fhe
RAN S < qua t ions m c non-convent iom~l way
Most methods for the c mputation of the vis-
cous flow aro md ships employ eifher fhe pressure-
correction or the artiflcicl comp~ bility method
Recogmi mg fnat particularly for flows et ve y high
Rey olds n mbers fhese itemtive medhods have dif-
flculties in, respectively, restormg fhe couplmg be-
ad reahing c stedy-state sohtion, m alterm~tive
aproah hcs beff~ chosen his hcs g own out to
fhe code PARNAS505, which is MARINs propri-
etary RANS-solver for ships Hoekshc & Egc 1998,
Hock trc 1999) his code is cunently in use for
quality csses ment of hull desigms on ~equest of ship
yards, r~vies md of her customers ~ our cpproah
fhe coupling between the equations is mcmtcined in
fhe iterctive solution, m impo t mt fator for robust-
ness Reduction of the si:D: of fhe equation ys-
tem is ahieved by dividmg the computation domcin
mto mb-domcins Eah mb-domcin is c g id pkme
roughly perpff~dicohr to the mcin st~eam clong the
hull ~ fhe iterctive process these sub-domcim are
visited in d wn tream di ection, which is ob iously
fhe be t choice By c specicl update of fhe pres-
sure fleld cfter such c weep f ough the domcin, the
convergence rcte is erJvmced Between 40 md 100
weeps me usually sufflcient to ~educe the maxim m
ncrm of the chmges of the non-dimensiom~l varictles
betw en successive iterctions to below 10 4 Mme-
over ehbomte veriflcation st dies have ben done et
model scale R by Hoekst~a & Bgc 1999 md et full
sccleR byEgc&Hoekshal999
After cbrief outlme of fhe origincl fectmes of PAR-
NASSOS in Section 2, this pcper focuses on c n mber
of new elements which bring us closer to m optim m
m combini g robust ess, acmay md eflflciency
I Fir t, the poss~bility is inchded to choose k ger
sub-domcim he single-g id-phne sub-domcin
hcs originclly ben chosen to minimize mem-
ory requi ements, but on present-day mahmes
memory is cmply cvaibble So it is now pos-
sible to choose fhe size of the mb-domcins in
the r mge fi om one g id pla to even the com-
plete domcin ~ Section 3 I w will describe
this techmique m mme detail
2 he systems of Imear equations are solved with
preconditioned GMPES he p~econditioning
OCR for page 883
technique, which performs well on c vector pro-
cessor md is suitable for implementation on par-
cllel mahmes, is partly re ponsible for eff-
ciency md flexibility: it allows to choose She size
of the sub-domcms between only one g id plume
md She entire compnhtioncl domain it will be
described m some d tail m Section 3 3
3 Originally, we discarded She dfff sion m She
mcmtrresm direction md also c part of She con-
vection terms m reversed-flow zones With
larger mb-domcins as described alder 1 These
terms c m readily be included implicitly, while in
the pkme-by-pla version they c m be inco po-
rcted explicitly (Section 3 4) it will be show
that these temms are msig dic mt indeed inrep~e-
sentative circ m st Aces, but that nevertheless The
behaviour of the iterative solution is improvedby
including Them
Section 4 w illustrate The effects of the clove cd-
ditions to PARNASSOS by showing re mlts of cppli-
cation to the flow aro md The stern of The HSVA md
Mystery takers Section 5 d scribes some current
developments Croat how c flee no Ace m be mcor-
pomted, md in The flncl Section 6 The main ahieve-
ments are summarized
2 Original solution strategy
As mentioned before, PARNASSOS solves The
f fly coupled steady moment m md mass conserve-
tion equations in primitive variable fomm, without ~e-
sortmg to pressure conection methods or The artff-
cicl compressibility cpproah Additional t mspo t
equations associated with the turbulence model me
h ected es uncoupled from mass md moment m equa
tions he g ids are less med to be body-fltted, but
may be generally non-or6hogomr1 hey are sketched
towards the hull m order to resolve The g cdients in
The bo mdary bye' he gowxning equations are inte-
g cted dow to the wall (no wall functions used) even
for full-sccle Rey olds n mbers he i herently ve y
high c pect ratio of the cells near The hull puts high
d m mds on the solver for The linear systems, which
is one of The reasons to mcmtam the cm phngt~w en
The equations m The iterative solution
For turbulence modelling in PARNASSOS, The
concept of the isotropic eddy viscosity is used it
is possible to choose from c wide ridge of tmbu-
lence models, varymg fr m algebraic models to sev-
e~al two-equation turbulence models Ego & Hoeksh c
2000)
A detailed description of the math maticcl model,
The computaticrurl g id md the PDEs in curvilinear
coordinates is given m R f Hoekshc 1999 Here w
reshict ourselves to c sho t description of The dis-
cretization ~ PARNASSOS The foll wing Smite dif-
ference approximations are used:
e For The m mericcl evaluation of th g id met-
ric terms w use second-order, central difference
schemes
e in the continuity equation w use c second-order,
th Copout backward scheme for the mcinsh eam
derivative, md c thi d-order four-pomt scheme
with c -i ed bits in nommcl a d gir6h-wise di~ec-
tion
e For The derivatives of The velocities that occur
in The convective terms w use c second-ord r
upwind scheme m sheamwise dinection, ad c
thi d-order upwind scheme for The nommcl md
gir6h-wise direction
e For The g adients of the pressure m the momen-
t m conservation equations w use c third-order
scheme For stability reasons, the bias ht. to be
opposite to that of The 'conespondmg' derivative
in the continuity equation Hoek Ha 1999)
e All second derivatives in The diffusive temms are
dismetized by second-order central-dffference
schemes
vation equations are dismetized et let t second-ord r
accurately To avoid negative turbulence q entities,
only m The h mspo t equations in the turbulence mod-
els we use c -i t-order upwind sch me for convection,
while for diff sion The same dismetization is applied
as in The m oment m conservation equations
Application of The cfore mentioned dismetization
leads to c huge set of non-linear tlgebr tic equations,
to which quasi-Newton linearization is applied h or-
der to reduce the size of The die rete equation yst m,
PARNASSOS ut s a marching solution t h me The
velocities md presume of one g id plume across the
mamsOetm direction (say, a Bacon t mt pane), are
solved simult meously The g id plumes are visited
m d wn tr mm order, while the elliptic character of
The RANS-equations is m merically nYmered by it-
eration Each step of this iteration t h me inchdes
not onh 6X dowmOetm sw ep th ough The computa-
tion d ma:, m which The eddy viscosity, the veloci-
ties md the pressme are updated, but also a upwind
sweep in which only the pressme is updated in the
t quel of 6 is pap r, this scheme will be denoted as
The 'global iteration'
The plare-ts-plane marching yields of comse
enormous savings m memory requirements But the
OCR for page 884
price to be pcid is tbat the upstr cm commumication of
d wn tream occunen es may proceed mther sl wly
vie the gdobal iterction process Neverdbeless, c robust
solution procedu e bas be n ahieved u ing c proper
pr conditioning tecbmiqu
3 New elements in solution proce-
dure
3.1 Larger sub-domains
If memory is not c serious limitation, fLe si:D: of
fLe sub-domcms c m be in ecsed So, in the n w cp-
proah, eah sub-domcin is defimed by c number (say,
gg~l)ofconsecutiv gid-pimes Incsub-domcin
fLe variables ce solv d simulta ously We con truct
c ystem of lin ar equations, fr m which w obtam
c corr ction for the pressu e md v locities et cll g id
nodes m fLe cunent sub-domcm be r suiting coef-
ficient matri is g times bigger tb m in the tmditiorul
pimeby-pla cpproah m PA tNASSOS, ad con
tains terms representmg the coupling betwe n the g
pimescov redbytb sub-d mcm
~ the planeby-pi me solution trctegy we are used
to cpply severcl steps of the quasi-Newton iterction
m c sub-domcin m order to let the global iterction
safely cor rge it tmn d out tbat wibb g > I it is
possible to mcke only on tep of fLe qursi-N wton
itemtion befor going to the next sub-domcm Hen e
fLe global itercti m acoumts completely for tb non
lin arity, while c sig ific mt saving in computation
time is obtamed
Let us denote fLe m mber of g id nodes m sheam-
wise, nommcl md gi thwise di ecti m by NX, NY md
NZ,rspectiv Iy~fLeextrmecaseofg NX(mb-
domcm = complete domcm), the global itemtion be-
comes the quasi-Newton iterction for the whole prob-
lem ~ fnat particular case, there is no n ed for
fLe pressu ecorr ctmg UPSh eam sw ep, which is nor
mclly p ut of c global iterction cycle
We fu fLer note fnat, even if g NX, fLe tmbu-
len e m odel equation(s) is (are) still h ected es u m-
pled fr m the of ber equations Suppose tbat fLe tb ee
v locity components, fLe pr ssme md fLe eddy vis
cosity wibb g id indices i, j, md k ce denoted, r -
pectiv Iy, by u~jt, v~jt, w~ jt. P~jt md (Vt)~/t ben
fLe complete non-linear pro'blem tbat is obtam d cfter
discr tization c m be w itten es
F{u~jt.V~it Wlit Ptit {V[i,jt) O
mwhichF isav ctor-valu dfmctionwibb 5 xNXx
NY u NZ compon nts if n coumts fLe number of
global iterction steps, on cycle of the global iterction
process consi ts of the following teps:
I compute the eddy viscosity {v~),jt es c fu tion
of u,jtl,v,jtl, mdw,jtl Inobberwords: solv
fv,),jtaproximatelyfiom
F{u,jtl v,jt 'W,jt -P,it .{V,i,jti O'
If c one- or two-equation tmbulen e model is
used. the h m port equrtion(s) are till solv d us-
ing c pi meby-plane stmtegy
2 Imearize tb tANScquations, using {v~),jt
computed ctov, md c mpute u,jt, v,jt, w,/t,
md P,jt by solvmg the remltmg yst m of Im-
ear equrtions to c prescribed acu ay
Hen e w c umot expect c quadrctic speed of cor r-
gence, not even when cll cor ctiv terms bav been
linearized second-order acu ctely
3.2 The linear system solver
As bas be n mention d ctov, fLe coeflflcient me-
h i is g times larger tb m in the phneby-phne v r-
sion of PA tNASSOS, md it contains coeflflcients
which acoumt for the couplmg between the separcte
6=const mt phnes bis matrix will be denoted by
A Let fLe entries of A be g ouped in blocks so fnat
all elements multiplyi g fLe varictles m c ~ =constat
pkme form c block Let su h c block, which itself is
c square mch ix of size 4 u NY u NZ, be represented
byon e tryA~:,inwhichi md/arerow mdcolumn
mdices, with i for convenien e taken equal to fLe g id-
mdex of the relevat ~ =con t mt pkme wibbin the cu -
r nt mb-domcin As c result of the chosen discretiza
tion, fLe block A,l ca onlybe non-7mo if 1/—il ~ 2
Hen e A bas the following pentc-diagorul sh u tme:
dl cl fi o o o .
c: d: e: fi o o .
b3 c3 d3 e3 f3 0
O b4 C4 14 e4 f4
, . 0 0 be Ce de
· (1)
Usually, fLe coeflflcients e~ md f only contam ele-
ments coming from the dicretization of the pressure
derivativ in mcin str cm dir ction O Iy in case of
fl w separction, these matrix elements c m clso con-
tain mcll conh ~butions fi om the derivativ of the v -
locities m mcin sheam di ection We ce faed now
wibb the problem to solv the yst m of Imear equ~-
tions Ax b Dir et medbods cc out of the qu s-
tion sin etheyrequi e too much storage adfloctmg-
pomt operctions ber for, w prefer to use a iter-
ctiv medbod, suited for non-symmetric matrices, I kc
OCR for page 885
Bi-CGSTABO ( m der Vcrst 1992) or G~PS~
(Sasd & Schultz 1986) Usuclly, the most expensive
operctions in such medhods are th computation of
i merpro ducts , vector up date s, md 6he matrixwect or
multiplicetion in our particohr situation, 6he mch i -
vector multiplication is much more expensive th m m
i merpod ct or vector update Hen e w should mini-
mi:m the m mber of matrixwector multiplications, md
6herefore we choose G~PS~ es the linear system
solver (where ~ indicates 6~t G~PS is restarted
eve y Mit step)
3.3 Preconditioning
The sped of cor rgence of G~PS trongdy de-
p nds on 6he eiger~lue spech m of A Ideclly, 6he
matrix should have cll eigemvlues clu ter d; 6hen A
would be close to the id ntity metri md G~PS
will rapidly cor rge in order to brmg c given spec-
h m closer to thet idecl situation, on c m use c pr -
conditionmg techmique Instecd of solving Az b,
on solvesM~IAz M-lb,inwhichMissomeap-
proximation of A The mch ix M should be ec y to
constn t md for c given vector z the comp htion of
M-i z should be cheap There me several possibili-
ties: on c m hy to cppro imete the mverse directly,
for exemple, by using c poly omicl of 6he origirul
matrix A How ver, the totcl n mber of matrixwector
multiplicetiom will then not be less th m with full
G~PS Therefore, w have chos n c dfffere t kind
of preconditioning techmique: we exploit 6he fat that
6he equations me clmost parabolic in 6he str emwise
di ection, byn glectmgthetemms e~ adf m (1) dur-
mg the comtn tion of the pr condition r We note
6~t G~PS still hcs to solve 6he ystem of linear
equations wi6h the coefEcient matrix (1) 6~t in hdes
6hese block This is the essenticl dffferen e wi6h 6he
aproah m th phneby-pkme version of PARNAS-
SOS, m which 6he systems of lin ar equations do not
have my entries to aco mt for the coupling betw en
6he separcte 6=constat pkmes, sin e 6hey have cll
moved to 6he right-h md side vector b
The mchi M is conshucted m mch c way that it
has the block shuctme
ml O O
c: m O
b3 c3 m3
O b4 c4
~ O
O O O
O O O .
O O O .
m4 0 0 .
o
O
.. O
O 2)
0 0 be Cc me
ad it is mapproximation of
dl O O O
c: d: O O
b3 c3 d3 0
O b4 c4 d4
~ O
O O O
O O ...
O O ...
O O ... O . (~)
O ~
0 be Cc de
Hen e m~ should be m cpproximation of 6he block d~
To that end we con tn t c coupled mcomplete LU-
decomposition of d~ Then m~ i~U~ m which the
fators i~ md U~ are sparse md low r- md uppe tri-
ag lar, respectively The in omplete decomposition
is constn ted on c 4 x 4 block-level Eah 'enDy' in
i~ md U~ consi ts of c 4 >~ 4 block
Bodh 6he con truction md the tri mg kr solves us-
mg 6he fators i~ md U~ should not co t too mmy
floeting-point operctions A other requi ement is 6~t
bodh the conshuction of i~ md U~ ad 6he trimgmbr
solves ca be vectorized The last requi ement im-
plies 6~t th sparsity pattern of i~ md U~ should be
r gmbr, so 6~t indi ect cdd essing c m be avoided
Let 6he enhies of d~ be g oup d m square blocks
of si:m 4 so that all elements multiplying the velocity
com p on nts md pre s sure in c g i d p o mt fomm c b lock
L t mch c block be repr sented by on entry The
mch ix d 6hen hcs c b lock par s ity par em 6~t c orr -
ponds to the followmg 9-point discretization ten il
(4)
The coefificients m C,~v, Ch3. CSTV md C5z are
r htively smell, ad 6herefore, 6hose block me n -
glected durmg the mcomplete LU-decomposition of
d~ Fur6hermore, cll fill-m blocks me neglected
Hen e the mch i i~ U~ hcs c block sparsity pattern
6~t corre ponds to the following 5-point sten il
. U
is
(5)
Bodh in 6he con tructi m of i~ md U~ md m the tri-
agmlar solves there are no r curr n e rehtions when
looping clong c diagonal of 6he g id This me ms 6~t
by diagorul order mg, we c m w~ct ori:m the lo op s Ti i s
m be regarded es c kmd of hype phne ordering
Wheng 1,w usecspeciclprconditioningtech-
nique, m which only the coefificients m C,~,, ad Cs~v,
OCR for page 886
togedber with cll fill-m blocks me n glected Tbis pr -
conditi onmg tecbmi qu bas not been vect ori:D:d, s mce
computations in which g I are more likely to rum on
scclar mahin s, instecd of on c vector processor or c
parcllel computer
As mention d ctove, it is not n essuy to solve
6be yst m of n m-lmear equations for eah domcm:
it s ffices to cpply only on step of the qu~si-N wton
itemtion ff g ~ 2 As c remit, 6be r quir d CPU-
time bas be n re du cd by cpproximate ly c fat or of
10, compmed to the pkmeby-plane implementation
of PARNASSOS if m clgebrcic tmb ien e model
is used 6be redu tion m CPU-time c m be even more
sin e 6be number of global iterctions d crecses by m-
crecsmg g ff c on -equation cr two-equation tmbu-
len e model is used, this favou ctle effect bas not yet
occu red, which is probably du to 6be sh ongly non
lin ar sou e temms
~crecsing the size of the sub-d mcms (ie m-
crecsmg the valu of g) leads to higher memory r -
qui ements Bel w we estimate the requi cd mem
o y by listmg 6be most memory-comummg fields of
varictles Herem N denotes the total number of g id
pomts
e Tbe velocities ad pr ssme, coordinates of 6be
g id pomts, pr ssme of the pr vious sw ep,
qumtities for 6be tmbulen e model: I IN vari-
ctles
Tbe mch ix: 4 x 26 x 84x
~ Tbe precondition r: 4 x 20 x S4x
e Right-b md side md solution vector of linear sy -
tem: 4 x 2 x S4x
e For GMRPS, ff r started ewxy 25tt tep: 4 x
27X~z
Togetb r wibb sewxcl fields of varictles tbat are not
so large, we arrive et 6be followmg erLmote of 6be
number of REAL*8 varictles to be stored:
{~] 350g)N
(6)
order to in ecse 6be Adfiop-mte on on processor
of the Cmy C90, memory bak co fiicts should be
avoided Tber for, ff the lecding dimension of krge
arrcys is c pow r of two, on the Cmy C90 it is m-
crecsed by on Tbis gives m cdditional mcrecse of
6be memory usage by cpproximately 15%
A very robu t shategy is to avoid the use of c
marching scheme by choosmg g NX How ver,
fi om m mericcl experiments it bas cppeared tbat the
number of global itemtion steps does not mcrecse
markedly ff g equals NX,, 2 mstecd of NX, while the
memory requi ements c m ckmo t be balved
3.4 Including 'negligiblet terms
Tbe fiow aroumd mor or less slender bodies I kc
ships md cir rafts et high R y olds m mber is char-
aterised by c domi mt fi w dir ction D ffusion in
6~t dir ction is pratically n glig~ble, md we bave
discarded 6be cssocicted temms m the moment m con-
servation equations in seri g up PARNASSOS Also
6be cor tion temms for 6be domi mt fiow dir c-
tion bave clways be n tr cted somewbat differ ntly
fiom 6be cor ction terms m the Ictercl dir ctions
Although PARNASSOS bas been decling wibb fiows
with b oumd uy Icyer separct ion in f Le sense tbat there
iscrgionwher thefi wisagainstthedominmtfi w
di ection, w bave cssumed fnat fLe reversed-fiow r -
gion will be of relatively small extent md fLerefore
with low cor ctive h m port cgamst fLe mcin fiow
As c consequ n e we bave neglected the corre pond-
mg cor ction term
Wbile in c planeby-pkme marchmg scheme the
'n gligible' terms bave to be in luded explicitly, i e
m the right-bmd side w~ctor, m c sub-domcm of
g ecter extension fLe neglected terms c m recdily be
mcluded implicitly Accordingly, w bave ir ti-
gated fLe role of these terms to confirm fnat fLey are
negligible mdeed md to verffy their effect on fLe con-
vergen bebaviou
4 Numerictil results
Tbe present comparisons are pe fommed for fLe ccl-
culation of fLe fiow aroumd the HSVA-t mker et model
sccle R y olds number, R 5 x I o6, md aroumd the
Mystery Tmker both et model sccle R y olds m m-
ber, R S x I o6, md et full sccle R y olds m m-
ber, Rn 2 x 109 Except for fLe HSVA-taker,
wher C beci-Smith's algebmic tmbulen e model is
used (Cebeci & Smith 1984), the r suits showninthis
section w re g n rcted usmg Menter's on -equation
model ~nter 1997) Tbe g id aroumd fLe HSVA-
tmker consists of 137 x 81 x 37 nodes in longitudi-
rul, wall-nommcl md gi thwise di ection, r spectively
For the cclcoktion of the fiow aroumd tb Mystery
tmker et model sccle the g id consi ts of 161 x 81 x
41 nodes For fLe f 11-sccle computati m m in r csed
number of 121 g id nodes is used in wall-normcl di-
r ction, brmging the total number of nodes for 6~t
case to ckmo t I milli m Although the g ids used to
genercte the results shown in 6bis section cover only
OCR for page 887
6he ~egion aro md 6he afterbody of the ship, the cur-
~ent impleme htion of PA tNASSOS zll ws to :x-
tend the domam to mclude the b w part zlso
he inlet ad outlet p me are z=con t mt pla s
he inlet phne is located zt z 0.5l (midship) md
6he outlet p me zt z I.25l he exterm~l b mdary
is m elliptical cylmder, giv by:
y ~ (z—0.056l\2 l:
~0.149) ~ 0.140 /
he ~emainmg bo mdaries are 6he still water surfae,
p me z 0.0 5 6l, the symm etry p me of the ship, y
0, ad the hull surfae
0 6he hull smfae w specify Di ichlet conditions
for 6he th e v locity components, zs dictated by 6he
no-slip ad impemmeability conditions On the sym-
met y phne of th ship, z Di ichlet condition for 6he
v locity component normal to 6he bo mdary is used,
ad Ne marm condition for bodh t mgential v locity
components in 6his secti m, fiee surfae effects me
not t kff~ into aco mt, h nce on 6he still water sur-
fae bo mdary we use similar symmeby conditions
(doublebody c mputation) Section 5 presents some
current dev lopments aout how to incorporate flee
surfae effects Bodh zt 6he outlet phne md 6he exter-
rul bo mdary, a e at condition is not a~ilable for
z computation domam of finite extent Hence 6he lo-
cation of 6hese bo mdaries determines for z g et deal
6he quality of the bo mdary cond6tions At 6he exter-
rul bo mdary we impose Dirichlet bo mdary condi-
tions for 6he pressure ad for the v locity components
t mgential to th bo mdary, deriv d from m i viscid-
fl w computation At the outlet p me, w take as 6he
bo mdary cond6tion that 6he g zdie t of the pressure
m longit dim~l di ection v mishes At 6he mlet p me,
w use Di ichlet bo mdary conditions for the v loc-
ity components f diff sion m mam sheam di~ection
is not t ken into aco mt, this completes 6he set of
bo mdary conditions if zll diff siv terms are t ken
mto aco mt, we use zs ext~a bo mdary condition at
6he outlet phne that 6he normal derivativ s of th v -
locities v mish
f z one-equation turbulence model is used, suit-
zble values for 6he eddy viscosity (deriv d wi6h 6he
zid of m zlgebraic tmbulence model) are imposed on
6he inlet bo mdary, while v~ is set to :D:ro on z no-slip
bo mdary as well zs on 6he exterm~l bo mdary At 6he
odher bo mdaries, w use sim ilar bo mdary conditions
m 6he tmbulence mod is zs for mass md moment m
equations
A illustmtion of z poss~ble effect of choosing large
mbdomams is show m Fig I ad 2, which sh w 6he
convergence behavior for 6he computation of 6he flow
aro md th HSVA tarker on model-scale Rey olds
n mber f g NX,,2 the whole separation region
is contained in one mb-domain, in tead of bei g dis-
h~buted ov r more, ad z considemble improv ment
m the convergence behavior c m be obtzmed We
lik to point out that on the v rtical axis the maxi-
m m L~)-nomm md not the more common L -nomm
(which would hav mggested z much smoother con-
v rgence behavior) is measmed he CPU-time on
one processor of z C zy C90 for the computation cor-
~espondmg with Fig 2 is zpproximately I O minutes
100 i
ID .
10
10 '
10'
i,\
T~_,,
`` ;
0~2 ~ ~~
0 ' ~ , ~\
dp
- - - dw
dev
10-8 ic 20 i~ lo so
Numbar of 9 oba t rat ons
Figme 1: Convergence behavior for HSVA-taker,
C beci-Smi6h's turbulence model g 4 Rn 5 x
Ne t, we st dy the effect of t kmg into a o mt zll
convectiv terms on 6he convergence behavior of the
global itemtion Fig 3 sh ws 6he convergence be-
haviour for 6he computation of 6he fl w aro md the
Mystery tar ker on m odel-scale Rey olds n mber 6~t
m main tream di ection is d opped, whe~eas Fig 4
shows the cor~rgence behavior obtained when zil
convectiv temms me t ken mto aco mt It zppears
6~t t king all convectiv terms mto aco mt has z
tzbilising i fluence on the global convergence be-
havim his is c msed by 6he second-order upwind
scheme for th g zdients m 6he t~eamwise convec-
tion, which giv s z positiv conh ibution to the main
diagorurl of the coeffcient matri For comparison,
Fig 5 shows the convergence behavior m case g I,
also for model-scale Rey olds n mber It zppears 6~t
6he n mber of global itemtion steps slightly mcreases
ff g is inmeased from I to 4 How v r, 6he requi cd
CPU-time is reduced considerably smce zt eah do-
OCR for page 888
100 _
10-'
10-
10
10 '
10-'
10'
10 ~ _
10'
"~:.'.~\
, _ \~
~ F~ I I
10 zo 30 do so
Numbol otg obal wratons
dp
dv
dw
dw
Figme 2: Convergence behmior for HSVA-tarDcer,
C beci-Smif 's tmbulence model g NX, 2 Rn
sxlo6
~n ~
M nter's tmbulence model No terms oglected g
4 Rn s,~lo6
Figme 3: Co orgence behwior Mystery tarD~r,
M nter's turbulence model St oomwise co oction
neglectedinr or od-flowzo o g 4 Rn s,~lo6
10' ,
1
n
10
10' ~
Figme 5: Co orgence behmior Mystery tar~cer,
M nter's tmbulence model No terms oglected g
I Rn s~lo6
OCR for page 889
OCR for page 891
OCR for page 892
OCR for page 893
OCR for page 894
OCR for page 895
OCR for page 896
Representative terms from entire chapter:
convergence behavior
main only one step of the quasi-Newton method is
applied per global itercti m tep Fig 6 sh ws She
converge Be behavior obtained by i xeasmg the size
of She mbdomams g 16) it appears mat, m con-
ha t wish the results show m Fig 2, She speed of
converge Be is not improved by choosing She value
of g larger th m 4 This c m be explained by She
fact At in Menter's model She converge Be behav-
ior c m be d teriomted by trongly non-linear so t Be
terms, whereas in Cebeci-Smith's tubule Be model
She eddy viscosity is c mpletely delerromed by She
velocity field
The CPU-times on one processor of a C ay I: 90 for
She computations cone pondmg with Fig 4 md 6 me
approximately one hour We note that in practice it
is s fflcient to red Be the maxim m of the ch mge in
She pressure to below 10~5 He Be in that case She
amo mt of req ired CPU-time is sig die mtly smaller
Figure 6: Converge Be behavior Mystery tarDer,
M nter's turbule Be model No terms neglected g
16 f s5 1lo6
Next, w study She i flue Be of the Rey olds m m-
ber on the global convergence behavior Comparing
Fig 7 with Fig 4 shows At She speed of converge Be
d creases somewhat, but the mcrease m the m ml er of
iterations is not ve y d cmatic This is in accord mce
wish previous experie Be that She corergence peed
of PARNASSOS is essentially mdependent of f in
6 is case, the requited CPU-time is approximately i
hour, the mcrease being mamly due to the mcreased
n mber of g id nodes
~ order to st dy She elf at of the terms that have
l en neglected so far m PA tNASSOS, some integ al
ID
10
ID
0 _
0- ~
1
is
0
,~,
Number of 9 oba t rat ons
Figure 7: Converge Be behavior Mystery tacker,
M nter's tubule Be model No terms negl Ted g
4 f ': 1~:~
q mtities md maxim m,mmim m mdmemvalues
of the urJcoow s me given for the My te y tarDcer in
Tables I md 2 The different versions of PA tNAS-
SOS have been mdicated as follows
e PNS: sheamwise dfff sion md convection in
streamwise flow separctim legions, fUI < 0),
are neglected
e PNS: Onlyneglects sheamwise dfff sion
e PNS : No temms in the RANS
The integrals in the definitions above are calculated
with Gaussian quadrature rules assuming a bi-linear
variation of the unknowns between the grid nodes.
The area Q for the calculation of the wake fraction
is the propeller disc which has been defined by:
· x = 0.989L
· Axis of the propeller at z = 0.0166L
· Hub radius equal to zero.
· Propeller radius R = 0.015L.
The maximum, (max.)' minimum, (min)' and mean,
(med)' values are obtained for the following flow vari-
ables:
· Us, axial velocity.
· Vw, cross-stream velocity, ~/(U2~2 + (USA.
· UOT, axial velocity at the outlet plane.
· Cp, pressure coefficient, Cp = p/ 2 p (U ~ )2.
vt, eddy-viscosity scaled by UiL.
Version
CD X 104
CD X 104
Umed
VWmed X 10
CPmed x 10
(Vt~med X 104
Wf
Umin
Umax
(Umin)O
V Wmax
CPmin
Cpmax
(Vt~maX X 104
Table 1: Model Scale Rn = 5 x 106
PNS++
1.606
0.741
0.652
0.632
-0.361
0.262
0.619
-0.068
1.064
0.686
0.311
-0.111
0.148
1.261
. PNS
1.605
0.747
0.653
0.631
-0.361
0.262
0.624
-0.121
1.064
0.687
0.310
-0.111
0.149
_ 1.267
PNS+
1.605
0.747
0.653
0.634
-0.364
0.262
0.621
-0.071
1.064
0.687
0.311
-0.111
0.149
1.265
The results in these tables show that the effect of
the "missing" terms is very small, especially for full
scale Reynolds number. However, in case of flow re-
versal, the missing convective term has a larger effect
on the maximum reversed flow velocity Umin.
Fig. 8 shows the limiting streamlines computed by
PNS++ both for model scale and full scale. A 0 in the
figure indicates where flow separation was detected.
The results of PNS and PNS+ are not shown in these
Version
CD X 104
CD X 104
Umed
VWmed X 10
CPmed x 10
(Vt)med X 104
Wf
Umin
Umax
(Umin)O
V Wmax
CPmin
Cpmax
(Vt)maX X 104
f ? ( ~V) ~
Table 2: Full Scale, Rn = 2 x 109
PNS
0.791
0.537
0.657
0.633
-0.304
0.109
0.353
0.000
1.073
0.804
0.410
-0.155
0.203
0.626
PNS+
0.791
0.537
0.657
0.633
-0.304
0.109
0.353
0.000
1.073
0.804
0.410
-0.155
0.203
0.626
T PNS++
0.791
0.535
0.657
0.632
-0.303
0.109
0.353
0.000
1.073
0.803
0.411
-0.155
0.203
0.626
Figure 8: Limiting streamlines Mystery Tanker,
above: Rn = 5 x 106, below: Rn = 2 x 109
nn4
nn~
0.02
\ 0.0 1
].QO~
0.01~
0.02-
0.oo 0.01 0.02 0.03 0.04
y/L
+ ozo
030
040
K O 50
090
070
9 0 RO
Fig;ure 9: Isolines Ul for ~ 0.989l, solid 1ine: PNS,
dotted: PNS , d4shed dotted: PNS Rn 5 >~ I O9
O.E
0.02
\ 0.01
. .
1.00
.04 ' / ~ ~ + 029
~ / / / ~ 030
0.03 1 / ~ o 040
1 1 1 / ~ oso
11 / ~ 070
~// ~ 050
// ~ 9090
N .
0.01
0.02
ooo 001 002 003 004
y/L
Fig;ure 10: Isoli es Ul for ~ 0.989l, solid ii e:
PNS, dotted: PNS, dashed dotted: PNS Rn
2xlO9
003
oo
\ oo~
nno
A n~
. :~/
~ : ~
O.OZ
ooo 001 002 003 004
Ott
ot4
ot6
017
ot9
02t
024
025
027
y/L
Fig;ure 11: Isoli es Cp for ~ 0 989l, solid ii e:
PNS, dotted: PNS , d4shed dotted: PNS Rn
S~lo6
0.04 `~ 4
003 // / /
002
M ///~ ,/'//~ ~
000 ~ /?
001 ~ ~
0.02
ooo 001 002 003 004
Ott
014
ot6
ot7
ot2
02t
024
025
027
y/L
Fig;urel2: Isoli es Cp for ~ 0 989l, solid ii e:
PNS, dotted: PNS , d4shed dotted: PNS Rn
2~109
on
nn
on
\ 0.0
o.o
0.02
000 001 002 003 004
y/L
405
415
425
435
445
455
465
475
485
495
505
F'gme 13 1 olin s g~ofv ~ fm 89 ohd
nn!
0.00 t
].0
0.02
000 001 002 003 004
y/L
Figme 14: Isolines LoglOfv~) for ~ 0.989l solid
line: PNS dotted: PNS dashed dotted: PNS
R 2>~109
004
0.03
00
001
onc
nn~
0.00 0.01 0.02 0.03
y/L
+ so
45
0 40
~5
30
zs
zo
Figme 15: Logl O of th di fe~ences between PNS md
PNS in the mcg itude of the tr msv rse v bcitv
vctorforphnex 0.989l Rn Sx106
0 ~4~5 003 ~ 0 ~40
~ 445 1 ~ 30
~4~55 0.02 )~ 6 ~25
485 ~ \1~/ /
485 \ 0.01 |j~,// /
000 :~
00 ) \\\
002 ~ ~ \ +
0.00 0.01 O.OZ 0.03 0.04
y/L
Figme 16: Logl O of th di fe~ences between PNS md
PNS m the mcg itude of the tr msv rse v bcitv
vctorforphnex 0.989l Rn 2>~109
figures, because these limiting rresmlines me vi t -
ally identical O Iy m the neighbourhood of The
prop I ler planes the computed velocities md pressure
differ slightly Therefore, Fig 9-14 sh w the com-
putedresultsindhepropellerpla x 0.989l These
results show that the effect of The neglected terms is
not large, md for model scale R y olds m mber, The
d fferences are located mainly in the region that con-
tains -I w separation From These figures it appears
That for f 11-scale R y olds m mber the i fiuence of
The rrsdionslhy neglected terms in PA BASSOS is
even mcllerthmfor mod 1-~t R y olds n mber
D spite These small d he~encei, The velocity fields
computed by PNS md PNS sh w hardly my dif-
ference, es is clearly demon trctedby Fig 15 md 16
They show The isolines of The Logl 0 of The differences
m the mcgmitude of the transverse velocity vector
Fig 15 shows the results for model-sccle R y olds
n mber in The propeller plume, whereas Fig 16 shows
similar re mlts for full-sccle R y olds m mber
5 Current developments: Incor-
porating the free surface
By basing The m odellmg a d The solution cpproah
on The physics of the problem considered, m accurate,
efficient md robust method has thus been obtained A
similar philosophy md riles our crorenr research on
mcorpomtion of The tree smfae
Viscous effects on the wave pattern of c ship me
Restively limited, md are essentially cord ad to The
tem are if we tee the tem waves md wake
out of taco mt for c moment, experience from vali-
datiom shows That the nonlmear fiee-surfae potff~-
ticl -I w model cm produce a mate predictions of
The wave pattern A tANS-sohtion subject to The
kinematic md dynamic fiee-surfae bo mdary condi-
tions FSBC) should reproduce the same _9e pat-
tern everywhere except m The tem melt Conversely,
ff w -i st use c potential flow code to compute The
wave pattern, md Then solve for The viscous flow m-
der That given wave surface (imposing just The kme-
matic FSBC), we should get The same solution, i e
The dynamic condition should be automatically satis-
fied es w 11 Fig 17 ilhstmtes The result of 6 is com-
posite cpproah Wmdt & Hess en 2000) based on The
codes tAP D Oven 1996) md PA BASSOS, for
6heSeries60modelatFn 0.316. Rn 34x106.
The bottom huff shows isolines for The wave elevation
aro md the cfterbody, The top half sh ws the cone-
ponding isolines for The en or m the dynamic FSB C,
expressed es c wave elevati m difference The error
is sigma mt only m The tem area md is negligi-
ble eve wh re else Therefore, for mo t of the do-
mcm w have Weedy fo md c solution equivalent to c
f 11 RANS FS solution, with a effort ju t marginally
Urger thm that for c doublebody RANS solution
This easy md most efficient procedure will soon ~e-
pkce The doublebody Afro toh for m ost of our p~a-
tical computations
R m tins, however, the stern are a where viscous ef-
fects on The wave pattern me present, es show by the
Urger pressme differences et the wave surface Lo-
cally, c full RANS FS method is still needed While
most such methods use c time-depffndent solution cp-
proah, msteady waves may deny the cpproah to
teddy state considerably We prefer to solve the
steady form of The fiee-surfae problem, by ite~a-
tion rather 6 m by time-stepping Refs Raven & m
Br mmelen 1999ad mBrummelff~&Rcven2000
derive c new fomm of The FSBC Citable for iterative
solution, md show The superior performance for c 2D
test case The iterative process converges ve y la t,
ad The total effort is just 2 or 3 times That for cor-
respondmg problem without tree smfae The same
aproah is n w bemg implemented m PARNASSOS
to be tested in 3D, ad will hop huh suppl merit the
composite procedure jn t described
6 Conclusions
The re mlts for The My te y t mker have show
That the effect of the temms that have been negdected
so far m PARNASSOS is indeed very small B pe-
cially at full scale R y olds m mber, it is justified to
omit diffusion in main tream di ection altogether md
to discard part of the convection terms m regions with
-I w separation H wever, especially The Utter temms
m have a positive i fiuence on the omen 3mx:e be-
hsviour 11 ~relore, it is better to take all convective
terms mto taco mt It The cunent production version
of PA BASSOS, we also take The sheamwise dfffih
non into a o mt Although we have d monshated
That The effect of this temm on both the convergence
behavior md the computed results is ve y small for
The Mystery t mker, taking all temms mto taco mt re-
moves my do ht in this re pect
The new feature in PARNASSOS, a spa -
marching solution with sub-domains of h get size
6 m a single g id plume has improved The efficiency
ad The fi:xibility: it is n w possible to choose the
sin: of The via domains m the rmge from one-g id-
phme to even the complete domain The phmeby-
phme solution sh ateg -has The advmtage that it is ve y
efficient with respect to memo y usage A fill scale
compohtion using I million g id pomts c m be done
onaPCwi6h only 64Mhyte of mainmemo y
n1L
02 04
06 08 1
xlL
PP
Figme 17: Cpmposite sohtion for viscous flow mder wave suriape Bottom hcff: wave elevation aro md stern
(mm p~d) Top half: error m dynamm FSBC (vi pous), expressed as wave elevation differe pe Contour interval
O OOOSl.
The multi-pkme tmtegy, on the odher h md, hcs 6he
cdvmtage that it is faster, mcmly bec mse the qua si-
N wton Imearization process pm be shffted to 6he
global iterction in praptice, the best tmtegy cm
be chosen depff~dent on 6he mcchine a~ihtle for
6he computations Treating 4 phnes simultmeously
seems c good choice in mmy pmpticcl situations in
6~t case, 6he memo y requir ments mcrease only by
cpproximately c fcotor of 4 compared to the phne-
by-phne solution shategy Even c flow compubtion
at f 11-socle Rey olds m mber c m 6hen till be done
on c workstatipm or c PC 6~t contsins et lecst 256
M ytes of m mo y
ff enough computer memory is a~ihtle md ff
m algeb~aic tmbule pe model I kc the Cebeci-Smith
model is used, mcrecsing 6he size of 6he sub-domcms
is c good option smce it red pes the m mber of global
itemtion teps ffMenter's turbule pe model is used,
6he n mber of global itemtion steps does not decrecse
sigmffcmtlybymcrecsingg h 6henearflure,w
hope to remedy 6his phenomffmoby hectingthe non-
linear sompe temms diffe~ently
More 6 m 80 % of the CPU-time is pent on 6he
solution of the systems of linear equations, usmg p~e-
conditioned GMPES The Mflop-mte that is obtamed
on c ptor computer like 6he Cmy C90 is quite high:
on one processor of 6his mcchme, cpproximately 300
Mflop/sec is cthined This is cchieved by choosing
c suitdole ordering m 6he computation which cl c m-
vents recmre pes in 6he consh ption md cpplication
of the p~econditioner
We are co fldent that 6he erJumcements i po po-
mted in PA tNASSOS es descobed m this pep r have
brought us closer to the optim m combim~tion of fl:x-
~bility, cocurapy, effciff~py md robustness in ship
tem flow compubtions
Acknowledgements
The mthorswishto6 mkH Rcvenformmystim-
ubting discussions md m merous mggestions for im-
proving 6he presentation of 6he pcper Th y clso th mk
6he miversity of Utr pht for the fort m routine 6~t
implements GMPES~
References
[1] TCebedandA.M.O.Smdth,Ano/ysuofTurhr
lent Boupdoy yens AcademicP'ess,November
1984
[2] L. E;a and M. Hoekstra, On the m mericcl ver-
ffcation of ship stern flow calcoktions h poper
pres pted ot the MA NET w rkshop, Borcelono,
Spoin, 1999
[3] L. Epa and M. Hoekstra, N merical Predic-
tion of Sccle Effects in Ship Stern Flows with
Eddy-Viscosity Turbule pe Models h Twmy
Th ind Symp s ium on Now I Hyd mdyn mics Appl i
cohons t Ship Flow o d Hull Fo m Derign, Os
oko, Jop p, September 2000
[4] M. Hoekstn' and L. Epa, PA tNASSOS: An ef-
flcient medhod for ship stern flpw cclculation in
Thind Osok CollogulumonAdw ced CFDAppli
cohons to Ship Flow and Hul l Fo m Design, Os
ok, Jopan,1998
[5] M. Hoekstra and L. Epa, An example of error
q mtification of ship-'ehted CFD results in Sew
m6h Confe~mce on Numericol Ship Hydnodyn m
ics, Nantes, Fnance, 1999
[6] M. Hoekstra, Numencol Simulation of Ship
St rn Flow w ith o Spoce Morch ing Novier Sto kes
Method Ph thesis,MARIN,Wzgenmgen,1999
[7] FR Menter, Eddy Viscosity Trasport Equa
tions md Thel Rehtion to 6he k—~ Model Joum
nol of Fluid Engineering, 119:876 884, L cem-
ber 1997
[8] H.C.Raven, ASoluLonMethodfw6heNonlin
eor Ship Wove Resut nce Pnoblem Ph 6hesis,
MARIN, Wzgenmgen, 1996
[9] H.C. Raven and H. van Brummelen, A
new aproah to computmg steady flee-surfae
viscous fl w problems h Ist MARNET
CFD w rk hop, Boncelono To be downlooded
f om http:/wwwmorin nVpublicotions/pgmetu
t nce html, l 999
[10] Y. Saad and M.H. Seh~dtz, Agenerali:mdmm-
lmal residual zlgori6 m for solvmg nonsymmet-
ric linear yst ms SIAM J. S i Stotut Comput,
7:856 869,1986
[11] H. van Brummelen rmd H.C. Raven,
N merical sohtion of steady flee-surfae
Navier-Stokes flow h 156h Int Work
shop on Woter Woves and Flooting Bodies,
Coesoneo, Isnoel To be downlooded f~m
h ttp://ww~ morin n Vpublicotions/pg~esu
t nce html, 2000
[12] H.A. van der Vorst, Bi-CGSTAB:A fast md
smoothly cor~rgmg varlmt of Bi-CG for the so-
hti m of nonsymmehic linear systems SIAM J.
Sci Stotut Comput,13(2):631 644,1992
[13] J. Wlndt and H.C. Ravea, A composite proce-
dme for ship viscous fl w with fiee surfae 2000
DISCUSSION
A. Sanchez-Caja
VTT Manufacturing Technology, Finland
In your very interesting paper you mention that
the outlet plane for your computations is
located at x= 1.25L, midship being located at
x=0.5L. This means that the position of the
outlet plane relative to the aft-perpendicular is
0.25L, which seems to be too small. Have you
checked the influence of the outlet plane
location on the results, for example on the
value of the drag coefficient?
AUTHOR'S REPLY
The choice of the location of the outlet plane is
based on our experience and is a compromise
between maximizing grid density and
minimizing the influence of approximate
boundary conditions. The inexactness of the
boundary conditions on the outlet plane is in
the first place related to the condition imposed
for the pressure. Where in most methods a
Dirichlet condition (undisturbed pressure) is
applied, we use a less restrictive Neumann
condition that allows the selection of a
smallest domain size, without significant loss
of accuracy.
As an illustration, we give some results for the
Mystery tanker at model scale Reynolds
number. The table below summarizes a few
meaningful results for three positions of the
outlet plane. It appears that the computed data
hardly changes when the position of the outlet
plane relative to the aft-perpendicular is moved
from 0.25L to 0.365L. Furthermore, the plots
of the isolines of the pressure and velocity
distributions obtained with the outlet plane at
x/L=1.365 are virtually identical to the plots
shown in Section 3.4, where the position of the
outlet plane is atx/L=1.25. On the other hand,
by moving the outlet boundary to, say, x/L=2.0
and keeping the same number of grid nodes we
would observe a much more drastic effect on
the solution. So positioning the outlet plane at
xlL= 1.25 is a deliberate choice.
Table 3: Computed values for several locations
of the outletplane.Rn=5xl06
xlL 1 1.135 1 1.250 1 1.365
CD X 10 1 1.605 1 1.606 1 1.606
CD X 10 1 0.688 1 0.741 1 0.744
Wf 1 0.624 1 0.619 1 0.619
U min 1 -0.069 1 -0.068 1 -0.068
U max 1 1.064 1 1.064 1 1.064
Vwmax 1 0.310 1 0.311 1 0.311
Cpmin 1 I. 1 10 1 -O. 1 1 1 1 -O. 1 1 1 l
Cpmax 1 0.148 1 0.148 1 0.148 1