**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

**Suggested Citation:**"Numerical Simulation of Three-Dimensional Viscous Flow around a Submersible Body." National Research Council. 1990.

*The Proceedings: Fifth International Conference on Numerical Ship Hydrodynamics*. Washington, DC: The National Academies Press. doi: 10.17226/1604.

Below is the uncorrected machine-read text of this chapter, intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text of each book. Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.

Numerical Simulation of Three-Dimensional Viscous Flow arounct a Submersible Body C-I Yang David Taylor Research Center Bethesda, USA P-M. Hartwich and P. Sundaram NASA Langley Research Center Hampton, USA Abst ract A second-order accurate, implicit, high resolution up- wind scheme has been used to solve the three-di~nensional incompressible Navier-Stokes equations in general curvilin- ear coordinates for the steady-state computation of the flow field around axisym~netric hull geometries at high Reynolds numbers. A hybrid algorithm with relaxation in the stream- wise direction and approximate factorization in the cross- flow plane is used to reduce the temporal splitting error. Three axisymmetric bodies with different stern shapes have been chosen for the present investigation to highlight the boundary layer development and to pinpoint flow separa- tion in the stern regions of the hull geometries. Three- dimensional viscous grids of C-O topology have been gen- erated around each body using a transfinite interpolation technique. Turbulence is simulated using algebraic eddy Hi I tv m Nazi cif Rn1H wi ', T.nm em J ... CAM The computed re- sults are compared with the available experimental data for on-body pressure distribution and radial and axial ve- locity profiles in the afterbody boundary layers with and without propeller in operation. The computed results show close agreement with the measurements in most cases. I. Introduction The flow in the stern region of a submersible is char- acterized by the presence of a thick and possibly sepa- rated turbulent boundary layer. There is an increase in the stern pressure drag and skin-friction drag, called thrust deduction, due to tl~e upstream suction produced by the propeller when operating. Accurate numerical prediction of the flow in the stern region of the hull both with and without the propeller operating is important to evaluate the afterbody thrust deduction, thus limiting the expen- sive resistance and self-propulsion model experiments in a towing tank. Although several numerical algorithms have been so far developed for computing the two- and three- dimensional incompressible Navier-Stokes equations, most of then are prohibitively expensive for very high Reynolds number flows typical in marine hydrodynamics. Several panel methods such as VSAERO (1) ale used In aircraft design and analysis practices. However, their true application in marine flows, particularly for comput- ing the flow field in the stern region, is limited because of the presence of fairly thick boundary layers in the stern . region. A viscous-inviscid interaction procedure Night be once again inadequate because of the possible presence of flow separation in this region. Hence, an accurate compu- tational method for Ellis problem should attempt to solve the inco~npressible Navier-Stokes equations. Aziz and IIellums (2) proposed a vector potential vor- ticity formulation to solve the three-di~nensional inco~n- pressible Navier-Stokes equations. Because of its large stor- age requirements and the necessity to solve three Poisson equations at each time level, the method is not very popu- lar. The direct extension of time-dependellt methods, both explicit and implicit, developed for the compressible Navier- Stokes equations to incompressible flows is not possible be- cause of the 'stiffness' of the physical problem associated with low speed viscous flow. 'ho circumvent this problem, Chorin (3) proposed the use of artificial compressibility when solving the equation of continuity, thus introducing an unsteady term to make the system hyperbolic as in flee case of compressible flows. In this paper, we compute complex flow fields around various axisym~netric bodies by obtaining numerical so- lutions of the incompressible Navier-Stokes equations in primitive variable formulation. A hybrid second-order ac- curate implicit high resolution upwind scheme is used to solve the system of conservation laws in general curvilin- ear coordinates. Tllree axisymmetric bodies with different stern shapes to highlight specific flow details ill that region have been chosen for the present study, primarily because of the availability of well documented experimental results for these bodies (4,5,G,7). Body fitted three-dimensional grid systems of C-O topology have been generated around these bodies using an algebraic grid generator based on transf~lnite interpolation procedure. The computed on-body pressures together with axial and radial velocity profiles are compared with the experimental data. In addition, the flow around the stern region of one of the bodies with a propeller in oper- ation is computed. The propeller is simulated by imbedding body forces in a disk located at propeller plane as suggested by Stern et all. Finally, in order to illustrate the ability of this scheme to simulate a separated flow, the results of a low speed vertical flow about a 3.5 caliber tangent-ogive cylinder at an angle of attack are presented. 59

II. Governing Equations and the Jacobian of the coordinate transformation is given by Using Chorin's (3) artifical compressibility formula- tion, flee incompressible Navier-Stokes equation are written in conservation law form for three-dimensional flow as Qua + (E* - Evict + (F* - Felt + (G*Gall = 0 (1) In Equation (1) the dependent variable vector Q is defined as Q = (p, a, v, W)T and the inviscid flux vectors E*,F*,G* and the viscous shear flux vectors EV,FV,GV are given by E* = Ala, u2 + p, uv, uw)T F = (pv, uv, v2 + p, VW)T G = (pw, uw, vw, W2 + p)T Ev = Re (0, As, sexy, ~z)T Fv = Re~l(O,ty=~7yy'tyz)T Gv = lle (O. The, Tzy'[zz) (2) The coordinates x, y, z are scaled with an appropriate char- acteristic length scale L. In Eq. (2), the Cartesian velocity components a, v, w are nondimensionalized with respect to the free stream velocity TOO, while the normalized pressure is defined as p = (P - poO)/pV2 . The kinematic viscosity z, is assumed to be constant, and the Reynolds number is defined as Re = vow. The artificial compressibility pa- rameter,6 monitors the error associated with the addition of the unsteady pressure terns If in the continuity equa- tion which is needed for coupling the mass and momentum equations in order to make the system hyperbolic. For large values of,0 or when the solution of Eq. (1) leas approached asymptotically a steady state, the continuity equation is accurately satisfied. However, the choice of ,B is also dictated by the important constraint of the stiffness of the partial differential equation as discussed in Ref. 9. Hence, ,5 = 1 has been chosen uniformly in the present computations. To develop the equations in a general curvilinear coor- dinate system, a coordinate transformation of the form ~ = ((x,y,z), ~ = <(x,y,z) and ~ = y~x,y,z) with and where x; y; z: J-i = det x~ ye z: x~ ye zig The Cartesian derivatives of the shear fluxes are obtained by expanding them using chain rule expansions in the (, I, 71 directions. The Jacobians of the inviscid fluxes E, F and G are needed for the flux-difference splitting and for the implicit algorithm. The Jacobian matrices in the different coordi- nate directions are obtained as the linear combination of the Cartesian Jacobian matrices as D = a1A* + a2B* + a3C* where D = A, El, or C with A = All, B = Ale, C = Ail , and (a~,a2,a3) are the row vectors of the T matrix. For the Jacobian matrix A, al = (x/J' a2 = ¢y/J, as = (z/J' and so on. The Jacobians in the Cartesian coordinates themselves are -O BE* 1 A* = = am O O -a B* = i7 = O ~Q 1 O -o * dG* O C = act = O The eigenvalues of D are A = allay (~1,>27~3,>~) 2u v w o v o o o w o 1 0 o a u O u a a n U U 0 2v 0 w u O p~ O u w v 0 2w = diag (U - S,U + S,U,U) where U is the contravariant velocity component in the cor- responding coordinate direction given by U = alu+a2v+a3w has been considered. Eq. (1) is rewritten in strong conser- vation law form as (Q/J)'+ (EEv)< + (FFv)` + (GGv)?7 = 0 (3) (E,F,G)T = [T] (E* F* G*)T (E,F,G)T = [T] (E* F* G*)T and S = [U2 + p~a2 + a2 + at] / = S(U,V,W,ai,a2,a3) The Jacobian matrices in different directions are diag- onalized using a similarity transformation D = RAR-t to obtain the eigenvalues of D. The rows of R-i and the columns of R are computed such that they give an or- thonormal set of left and right eigenvectors. The gener- alized similarity matrices R and R-i are given by (Ref. tT] = [be (y (z ~ R = ~ alp + Ail alp + UA2 -a2 -as BY 71z at+ V>l a2,B + V>2 al + as as ~3,0 + we) as p + W>2 a2 al + a2 60

and A2/p al a2 as -1 1 >1/~ al a2 as R = S2 (~2d1 + Ald2)/p 2$2/kaids 2S2/ka2d5 a3ds L(~2d3 + gilds)/> 2S2/k - ald6 a2d6 - 2S2/ka3d6 where 1; = al + a2 + as al = a' _ etc and d1 = (T31r2l )/k, d2 = (r32 - r22)/1c, d3 = (r41 - r2l )/k, d4 = (r42r22)/k, d5 = d1 + d2, and do = d3 + d4 Quantities such as rid represent the ith row and jth column of the R matrix and S = S(u,v,w,c`~,a2,a3). In hilly Reynolds number flows, it is appropriate to make use of the thin-layer approximation. The justification for such a simplification can be found in Bullions and Steger (10). Consequently, all viscous derivatives in the ~ and direction (along the body) are neglected. The viscous shear flux Gv (normal to the body) and its Jacobian can be derived after Steger (11). That produces the coefficient matrices Gv (Re ]) : ~:lvr' + ¢277 ¢1Wn + ¢29z and with ~O O O O ~ AGE - 1 0 Z22 Z32 Z42 Z=0Q=(Re'J) LOZ39Z33Z43~ O Z42 Z43 Z44 Using (G) we can write (5) as Discretization of the Inviscid Fluxes Consider a system of quasi one-dimensional, hyperbolic partial differential equations (Q/ J)t + He = 0 (4) where ~ = (, (,or y. Defining computational cells with their centroids at I = ,,8~ and their cell interfaces at I ~ 1/2, a discrete approximation to (4) is written as, ( ~ )^Qn+/\lH=0 (5) where At is the time step ,/~Qn = Qn+1 _ Qn and /\I( ) = t( )I+~/2 - ( )I-1/24/~8 · Superscript denotes the time level at which the variables are evaluated. To construct an approximate Riemann solver for the initial value problem in Eq. (5), each variable is regarded as an averaged state in each cell so that the flux difference is preserved in each cell and Eq. (5) can be regarded as an integral rather than a differential law. According to Roe's scheme (12), the flux at interface I it 1/2 call be expressed in terms of the left and right travelling waves, ¢1 = 71~' + 71y + 71~ Z29 = 39= + By + hz Z32 = 39X9y Z42 = 3r1=r/z Z33 = 71x + 371y + hz Z43 = 3 rack Z44 = 71x + By + 3 id ¢2 = 3(~Un + hymn + flown) The inviscid fluxes, the viscous shear fluxes and their Jaco- bians just obtained in Ellis section are ready for discretiza- tion later. III. Numerical Flux Differencing It is well known that upwind schemes possess an in- herent solution-adaptive dissipation that eliminates the ad- dition and fine tuning of artificial dissipation terns for nu- merical stability and accuracy required in schemes based on central diIferencing. In the present approach the invis- cid fluxes are discretized by using Roe's flux-differencing splitting concept (12) and the viscous fluxes are discletized lay the central differencing technique. We shall discuss the numerical disc~etization scheme in one din~e~sio~ arid as- seml~le them together for three-`lin~ensions. Hl+~/2 = Hi :t (OH F I/2) (6) (~`J)~\Q + I\HI+1/2 + /\HI-1/2 = 0 (7) Equation (7) relates the development of Q at the centroid to the waves at the interfaces according to their propagation directions. Defining a mean value matrix according to Roe (12) so that Dl+l/~ = D(QI, Ql+1 ) Dl~l/~Al~l/2Q = AHI+1/2, Using the above we can write (7) in delta form as [(,Nt]) - (Dl+l/2) ~1+1/2 + (Dl+-1/2) ~1-1/2] AQ = (Dl+l/2) /\I+1/2Q (Dl+ l/2) /\I-1/2Q where and with 61 Dl~l/2 = (RAR )I+l/~ = tR(A+A-)R-l~l l/2' (Dl+l/2 ) = (RA-R-1 )1+1/o (Dl+-l/2) = (RA+R )- A =(~A~+A)/2. (8)

The construction of D' /2 is accomplished in two steps. First, the analytic Jacobian D = dH/0Q is formed. D is found to be a function of metric coefficients, of the parameter ill, and of the state vector Q D = D(a, b, c, id, Q ) where a,b,c are the metric quantities <2/J, (y/J, (z/J, etc. Second, the local values aim, balm, cam, and All ~/2 are computed at interfaces and fed into D to form Do ~/2. In addition, in order to ensure that the flux di~er- ences taken over the six bounding surfaces, (i at 1/2,j,~), (i, j ~ 1/2, k) and (i, j, ~ it 1/2) of a three-dimensional cell at (i,j,~-) cancel out completely so that the present finite difference formulation essentially becomes a finite voulme method, a special weighted averaging procedure has been adopted (9), for example, (¢~/J)iJk = t(~k~jY)~j{hZ) (~j{~Y)~{jZ)]i where with ktjY = [({j§~+1 + ([jy)~_1]/2 fj( ~ = [( ~j+1 - ( )j-1]/2 The advantage of the delta form formulation is that the steady-state solution is independent of the time step size At. The justification for the accuracy of Eq. (8) as the approximate Rie~nann solver for the present problem is given in Ref. 9. In order to ensure that the present scheme has a high resolution capability, which is equivalent to preserving To- tal Variation Diminishing (TVD) property for a nonlin- ear equation, it is necessary to construct a quantity sub- ject to TV requirement. For a linear problem such as the one described above, the TV of the characteristic variable W = R-iQ can be forced to diminish in tine. Assuming R,R-i and A to be constant and using Else definition for the characteristic variables, the scheme given in Eq. (8) can be rewritten as [aim ~ Arn^~/2] twin = ~ ~777~/214 In (9' (at= 1,2,3,4) The above equation comprises four scalar linear equations. To enhance the accuracy up to third order, Else characteris- tic variables lF,in at the cell interfaces are reconstructed by using piecewise linear distributions as described in Ref. 9. We then obtain a family of implicit TVD scenes [( ^t Jo r7t iTI/2] Adorn [~\7n ]~ l wj~r7I, + (~l + ~~>r+n) Ghan ~+~ ¢7n ~/4~ /\~+~/2W'n [>m + {(~l In + (<l + W))~777 ) (10) X (urn ~ (~)rn,!~ )/'1] /\!T/2 Mom. where and '¢)m,l = ¢(7 m,|) rob _ ~ (Al_~/2l~m//\l+l/2wm), for /\1~1/21~nt 7t 0 ~ 0' for /\l~l/2wm = 0 Extending this TVD scheme to non-linear hyperbolic con- servation laws, we obtain [( /\~| ) ()m,l+l/2/~1+1/2 Until1/2/~1 - 1/2) ~ am = {'\m,l+l/~2 [(1 W)~,n,l+l/2 + (1 +~)Am,l+l/2] X(¢>m 1+1 am 1)/4} i\l+l/2wm {>m,l+l/2 + |(1 w)~7n,1 - 1/2 + (1 + w))m,l1/2] X(¢m I'¢)n1 1 - 1)/4} i\1 - 1/21/~7n (11) The extension of Eq. (11) to a nonlinear system of conser- vative laws is obtained in two steps. First, the nonlinear equivalent of Eq. (11) is formed. then it is multiplied R throughout from the left to get [( /\t,J ) - (Dl+1/2 i\1+~/2D+ 1/2 /\1-1/2 )] /\Qn ~ where (RK R-l)l+~/2/\l+~/2Qn_ (12 (RK+R-l)i_l/2/\l_l/2Q Kl+~/2 = +fA F+1/2t(1 _ w)~:F+1/~ + (1 + ~)AI~ 1/~] 4 ) t = d~ag(<~ , ~Jr, ,>3, ¢74 )I Discretization of the Viscous Fluxes To discretize the viscous fluxes, the derivatives at flee cell interfaces are approximated lil;e (U()j+~/2 = /\ j+~/2 tl, etc. wherever possible. Otherwise, they are computed from the central differences at the two nearest neighboring points. such as (U`~)i+~/o = 0.25 t(uj+~ - uj_~)i+~ + (uj+i7tj-~)i] The above one di~nensiollal discretization scheme can be used to estimate the interface-flux function in ~nultidi- mensional problem. Tl~e discussions of the extension of the fluxes di~erencing and TVD schemes from one-di~nensiollal to ~nulti-di~nensional application can be found in Roe (13) and Yee (1~1~. Our present strategy is to suln up all inde- pendellt discretizatiolls of the flux derivatives in each coor- dinate direction to forln a three-dimellsiollal forlnuation. 62

Tilde Differencing and the Inexplicit Scheme The backward Euler time di-fferencing of the three- dimensional conserva,tio~ equations with the tl~in-layer a,p- . . . prox~n~at~on Is /\ Q =t/\~(En+i ~ + A`~(Fn+~ ~ + A71l(G7~+iGv+i )] (13) Linearizing it about time level n, we obtain +~9En'~` + (0Fn)~\ + ¢~9Gn _ ~GV )A Amen =f/~(En) + /\~(Fn) + ,/\~(Gn _ Gn)] The left hand side is the implicit part and the piglet hand side is the explicit part of the for~nulatio~. The ex- plicit part is the spatial derivatives in Eq. 3 evaluated at the known time level n; its value diminishes as the steady state solution is approached. Hence, it is also called the residual . The L2 Noreen of the residual is often used as a measure of convergence of a solution. Discretize the in- viscid and viscous fluxes according to upwind difEerencing scheme and central differencing scheme respectively in (, and ~ coordinate direction independently according to Eq. (8) and then assemble then together. Eq. (14) becomes tf /\~J)A:+~/27Ni+~/' + A+ i/~/\i-~/2 BJ+~/2Ai+~/2 + B}_~/~Aj_~/o (C + Z)~+~/2~+~/9 + (C + Z)kI/~kI,/2]n~Qn = - Res(Qn) (15) where i, j, and ~ are spatial indices associate with the (, it and ~ coordinate direction. A+,B+,C+, arid Z are 4 x4 block matrices ~ flux Jacobians) associated with implicit spatial differencing in the coordinate directions by evalu- ating the metric terms at cell interfaces in each direction. Eq.~15) is solved by an implicit hybrid algorithm, where a symmetric planar Gauss-Seidel relaxation is used in the streamwise direction ~ in combination with approximate factorization in the remaining two coordinate directions ~ and y. It is used to avoid the At3 spatial splitting error incurred in fully three-dimensio~a,l approximate fa,ctoriza- tion methods. The hybrid scheme is unconditionally stable for linear systems and offers tl~e advantage of being con~- pletely vectorizable like a conventional three-dimensional approximate factorization algoritl~n. As a result, Eq. (15) becomes [M(B )j+~/2~j+~/2 + (B+ )j-~/2^j-~/24~Q =Res(Qn,Qn+~) iM(C + Z)~+~/~+~/2 + (C + Z)k1/2~kI/2~Q =Mi\Q Qn+1 = Qn + Chin with M = t~J + (A ji+1/~ + (A+)i-l/2] and the residual on the RHS indicates the nonlinear up- dating of Else residual by using Q7~+~ whenever it becomes availabe while sweeping in the ¢ direction back awl forth through the computational do~na,i~. Turbulence Modeling For la~ninar flow computations Else coefficient of molec- ular viscosity ,u = Ill is obtained front Sutherland's law. Turbulence is simulated using the Baldwin-Lon~ax algebraic turbulence model. For turbulent flow computations the lan~inar flow coefficients are replaced by, ~ = pi +pt (14) The turbulent viscosity coefficient lit iS computed by using the isotropic, two-layer Cebeci type algebraic eddy-viscosity model as reported by Baldwin-Lo~nax. In Ellis formulation . . At IS given by ~1 = { (tinner Y < Yc (lit )O?tter lY > Yc where y is the local distance measured normal to the lady surface and Yc is tl~e smallest value of y at Chicle the val- ues front the inner and outer region viscosities are equal. Within the inner region, the Pra,ndtl-Van Driest for~nula- tion is used (let )inner = I ~~ where I = k-yt1- e-(Y /A )] and ~~ is the magnitude of vorticity given by /( Al dv )2 + ( do _ bW )2 + ( Ho _ 0ll )2 N/ By 0,~ fez by ()x 0z and y+ = (rW/,tGw)Y ~ In the outer region, for attached and separated boundary layers, the turbulent viscosity coefficient is given by (lit )outer =It CcpFwa~eFI~ leb(Y) Forage =~7lin(ilmaXFmax, CWk~lmaXUdif /I;37~2ax) FI`leb(Y) =[1 + 5 5(C~leb/Yntax) ] In the above equation k~ and cop are constants and Unlay = 77la~[l~ly(l e (A /A ))] and Y7nax is the value of y at which F''~,aX occurs. Tile quaII- tity Waif is Else difference between ~naxi~nu~n and minilDUTn total velocity in Else profile U4if = \/(U2 + V2 + W2)n7,ax V(tG2 + V2 + W2)77tin The various constants in the model are given in Ref. 15 as ( 16) and 63 A+ = 26, ~ = 0.4, It'= 0.0168, Ccp = 1.6, Cat = 0.25, CI`leb = 0 3

IV. Grid Generation An algebraic grid generation procedure based on trans- finite interpolation technique has been used to generate the viscous grid around the hull geometry y. The three- dimensional grid generated around the body is of C-O type. The C-O grid is particularly needed to adequately resolve the wake region. Details of the theoretical aspects of the transfinite interpolation ~ncthod and the various mapping functions and their behavior in grid control is highlighted in Ref. 10. The method is fundamentally a two body grid generation procedure. The stern shape could be either open or closed. A modified osculating interpolation Injunction has been used in the present programs. V. Results and Discussion The hydrodynamic characteristic of the boundary layer flow around the stern of a ship is quite different with and without a propeller in operation. The action of a propeller produces suction that accelerates the flow upstream. As a result, the pressure and the skill friction drag around the stern increase and the thicl;ness of the boundary layer de- creases. The knowledge of the effective flow profile near the propeller plane and the amount of the added drag is essen- tial for designing an efficient propulsor. In the past, exten- sive efforts have been made to study the interaction between a propeller and a thick boundary layer experimentally and computationally (4,5,C,7,8,17,184. For the reason of sim- plicity, in most cases, axisymmetric bodies were chosen. At present, we perforce the numerical simulations of flow over DTRC Afterbodies 1, 2 and S without a propeller. In addi- tion the flow over afterbody 1 with a propeller in operation is analyzed. The results are then compared with available experimental data. The purpose of the simulations is to develop and validate a numerical scheme for analyzing the complex interaction between boundary and propeller. The details of the afterbodies are shown in Fig. 1, where rmaX is the maximum radius of the body, r is the radial distance measured from the body axis, ~ is the axial distance from the nose and L is the total body length. The afterbody length to ma:;i~num diameter ratios of all three afterbodies are different; they are 4.308, 2.247 and 2.018 for Afterbodies 1 2 and 5, respectively. ITere, the afterbody length is de- fined as the distance between the end of the parallel middle body and the after perpendicular. Furthermore, Afterbody 5 has an inflected stern. The hubs of all three afterbod- ies are identical at the position where the propeller can be mounted. The different stern shapes generate a large data base variation of stern flow suitable for the purpose of val- idation of the computational scheme. Flow over Axisy~n~netric Bodies Flow variables such as pressure and velocity co~npo- nents of an axisymmetric flow are independent of circumfer- ential variation. They can be obtained through a. set of de- generated equations based only on two spatial dimensions. However, in order to validate the three dimensional for~nu- lation we have just presented, full three dimensional com- putations over three different segments of a body ale per- formed. The sizes of the segments are 45°, 90° and 180° . A viscous grid with C - O topology around tlte body is gen- erated by a transfinite interpolation method. Grid domain extends to three body lengths both upstream arid down- stream of the body. Figure 2 shows a partial view of a C - O grid around a generic body. The grid size used for presented simulations is 91 x 25 x 49 in ~ (stream- wise), ~1 (circumferential) and ~ (normal) direction. In the ~ direction, the grid is clustered near the body with the smallest grid spacing 5 x 10-4 of the body length. In the direction, the grid is clustered near the bow and the stern. Grid sizes used for all three segments (45°,90°and 180° ~ are identical. Consequently, we can perform computations based on three different grid densities in the ~ direction. The differences among the solutions are negligibly small, therefore, only the results based on the 90° segment will be presented here. In all our axisymmetric flow computa- tions, the mixing length of the turbulence model has been modified according to IIuang et al (5~. A detailed analysis of Else measurement accuracies is not available. Ilowever, else standard deviations of mea- sured data were estimated from repeat runs. The standard deviation of the measured static wall pressure was less than 5 percent of their mean values and the deviation of the mea- sured velocities was less than 2 percent of the free -stream velocity (19~. The experiments of flow over Afterbodies 1, 2 and 5 were carried out at Reynolds numbers of 6.60 x 106 ,6.80 x 106 and 9.30 x 106, respectively (based on total body length). Figures 3, 4 and 5 show the comparisons of the computed and the measured pressure distributions over the surfaces for Afterbodies 1, 2 and 5, respectively. The body profiles are included in the figures in order to show the relationship between the pressure gradients and the stern shapes. At the region near the end of the body where the stern and the hub meet, the computed pressure shows a sharp decrease followed by a equally sharp increase. Based on some different numerical schemes, this peculiar feature has also been encounted by Chen and Patel (18) and Lee et al (20) in their computations. Chen and Patel attributed such phenomenon to the rapid change of geometry near the hull-hub juncture as well as the upstream influence of the complex pressure interactions in the tail region. A solution we obtained with a panel method exhibited the same fea- ture, as long as enough panel resolution around the hull-hub juncture region was provided. Unfortunately, the experme- ntal data lack the resolution needed to verify this feature. Figures 6, 7 and 8 show the comparisons between the computed and the measured velocity components at sev- eral streamwise locations for afterbodies 1, 2 and 5, respec- tively, where rO denotes the radius of body surface . For Afterbody 1, the agreement is very good. The computa- tion predicted the development of the boundary layer very well. For Afterbodies 2 and 5, the agreement in general is good, except the radial velocity profiles which are less sat- isfactory. The agreement deteriorates as the measurement location moves further downstream. The radial velocity component is relatively small and is more difficult to mea- sure accurately. For both Afterbodies 2 and S. the maxi- mum difference between computed and measured values is about At percent of free-stream velocity. The error bound of 64

the measurement is 2 percent. The measured axial veloc- ities for Afterbodies 2 and 5 are progressively slower than the computed velocities as the end of stern is approached, although the difference is small. Near the end of the stern, the measured axial velocities for both Afterbodies 2 and 5 show an inflection point which is not captured by the com- putations. The pressure gradient of the stern boundary layer is directly affected by the fullness of the stern shape. The stern of Afterbody 1 is the least full among all three bodies studied here. It caused only a mild adverse pres- sure gradient in the boundary layer surrounding the stern region. This may explain the reason why the simulation for Afterbody 1 is the most successful. In each computation presented above, 220 iteration steps were taken. The L2 norm of the residual dropped two orders of magnitude from 10-3 to 10-5 . 220 iteration steps were chosen since further iterations produced only in- significant variations. Each computation requires 17 CPU minutes on a CRAY-YMP machine. Flow over an Axisymmetric Body with a Propeller Flow experiment for an axisyn~netric body with a pro- peller were conducted at DTRC by Huang et al(4,G,7~. A propeller was mounted on Afterbody 1 at x/L = 0.983 . The geometrical and hydrodynamic characteristics of the propeller are given in Huang et al (4) and Huang and Groves (7~. The experiment was performed at a Reynolds number of 6.6 x 106 (based on body length). Numerically, the propeller effect is simulated by imbed- ding body forces in a disl; of finite thickness located in the propeller region. The details of this type of formulation can be found in Stern et al (8~. Distribution of body forces depends on the propeller's characteristics such as, thrust coefficient CT, torque coefficient Icy, advance coefficient J and radial circulation distribution G(r). The axial and circumferential body force per unit volume are obtained from the following equations: CTR9pG(r) MY IR P G(, )? d 4It R3G(r) 1rrJ2 /YX .fR P G(r)' do where f bx and f be are the body forces per unit volume in the axial and circumferential directions, respectively, Rh and RP are the radii of propeller hub and blade, respec- tively, and AX is the thickness of Else disk. The following propeller data were used ill our computation: J = 0.370, It T = 0.227, It q = 0.0453, CT = 0.370 The computed body forces are then incorporated into the right Lund side of Eq. 16 and forms a part of the residual. Based on the identical grid size and grid distribution used previously, a computation was carried out for 220 iteration steps. The L2 Noreen of tl~e residual dropped two orders of magnitude from 10-3 to 10-5 without suffering frown any numerical instability. From the results of the computation, the influence of propeller action can be detected up to about two propeller diameters upstream. The differences between axial veloc- ity components /\u=/VOO Title and without the propeller in operation at two streamwise locations are shown in Fig. 9 . The measurement locations are at X/L = 0.954 and X/L = 0.977 . The agreement between the computed and the measured values is very good. The differences in pres- sure distribution on the afterbody surface upstream of the propeller plane are shown in Fig. 10. The results of the computation are to the same degree of accuracy as those reported by Huang and Groves (7~. The swirl velocity was also computed but due to lack of experimental data to com- pare with, it is not presented here. Flow over Tangent-Ogive Forebody NVith the same numerical scheme discussed above and a modified turbulent model the flow over a 3.5 caliber tangent- ogive forebody was studied (21) at angles of attack of 20° and 30° and at Reynolds Lumbers in the range 0.2 - 3.0 x 106 . The purpose of the study was to investigate the Reynolds number effect on low speed vertical flow and to validate the numerical scheme. A C - O type grid was generated for the purpose of computations. The grid size was 97 x 40 x 91 in ~ (streamwise), ~ (normal) and ~ (circumferential) direction. Fig. 11 shows a comparison of the computed and the measured surface flow patterns at a Reynolds number of 0.S x 106 (based on diameter) and an angle of attack of 20° . The lines indicate that the primary separations are in good agreement. Top views show two distinct regimes in which the surface streamlines appear to collocate. Figure 12 shows a comparison of the computed and the measured circumferential surface pres- sure distributions. The Reynolds numbers are 0.S x 106 and 0.3 x 106 and the angle of attack is 30°. The Reynolds number effect is snore pronounced on the leeward side, near the nose. The agreement is good. VI. Conclusion The 3D incompressible Navier-Stokes equations was discretized by the flux-(lifference splitting and the implicit high resolution schemes. A discretized system of equa- tions was solved by an implicit hybrid algorithm, where a symmetric planar Gauss-Seidel relaxation was used in the streetwise direction in combination with approxi~na- tion factorization in the two remaining directions. The al- gorithn~ is highly vectorizable and suitable for computation on a modern superco~nputer. For the simulation of afterbody boundary layer flows, this method is proven to be effective. Inclusion of the body force propeller model poses no additional problem. To ob- tain a converged solution with the propeller model included, it does not require snore iteration steps in comparison with the case without including a propeller model. The method was also shown to simulate low speed vertical flow with good results. Acknowledgement This study was supported by The Office of Naval Tech- nology. The computing time of CRAY-YMP was provided generously by NASA Antes Numerical Aerodynamic Sin~u- lation (NAS) Program. 65

References 1. Maskew, B.M., "Prediction of Subsonic Aerodynamic Characteristics - A Case for Low Order Panel Meth- ods", AIAA J. of Aircraft, Vol 19 (2), Feb. 1982, ppl57-163 Aziz, K., Hellums, J. D., "Numerical Solution of the Three-Dimensional Equations of Motion for Rectangu- lar Natural Laminar Convection", Physics of Fluids, Vol. 8, Dec. 1965, pp 2182-2189. 3. Chorin, A. J., "A Numerical Metltod for Solving In- compressible Viscous Flow Problems", Journal of Com- putational Physics, Vol. 2, No. 1, Aug. 1967, pp 12-26. 4. Huang, T. T., Wang, H. T., Santelli, N., Groves, N. C., "Propeller/Stern/Boundary Layer interaction on Ax- isymmetric Bodies, Theory and Experiment", DTRC Report 76-0113, 1976. 5. Huang, T. T., Santelli, N.,Belt, G., "Stern Bound- ary layer Flow on Axisy~nn~etric Bodies", Proceedings, 12th Office of Naval Research Symposium on Naval Hy- drodynamics, Washington D.C., 1978, pp 127-157. 6. Huang, T. Y.,Groves, N.C., Belt, G. , "Boundary Layer Flow on an Axisymmetric Body with an Inflected Stern DTRC Report 80/064, Aug. 1980. 7. Huang, T.T., Groves, N.C.,"E~ective Wake: Theory and Experiment", Proceedings, 13th Office of Naval Research Symposium on Naval Hydrodynamics, Tokyo] 1980, pp 651-673. 8. Stern, F., Kim, H.T., Patel, V.C., Chen, H.C.,"A Vis- cous Flow Approach to the Computation of Propeller- Hull Interaction", Journal of Ship Research, Vol.32, No. 4, Dec. 1988, pp 246-262. 9. Hartwich, P.-M., Hsu, C. H., "High Resolution Up- wind Schemes for the Three-Di~nensional, Incompress- ible Navier-Stokes Equations", AIAA Journal, Vol. 26, No. 11, Nov. 1988, pp 1321-1328. 10. Pulliam, T. H., Steger, J. L., "An implicit Finite Dif- ference Solution of Three-Dimensional Flow", AIAA Paper 78-10, Jan. 1978. 11. Stager, J.L.,"Implicit Finite Difference Simulation of Flow about Arbitrary Two-Dimensional Geometries", AIAA Journal, Vol. 16, No. 7, Jul. 1978, pp 679-686. 12. Roe, P.L.,"Approximate Riemann Solvers, Parameter Vectors, and Difference Schemes", Journal of Compu- tational Physics, Vol. 43, No. 2, Oct. 1981, pp 357- 372. 13. Roe, P.L.,"Characteristic-Based Schemes for the Euler Equations", Ann. Rev. Fluid Mech., Vol. 18, 1986, pp 337-365. 14. Yee,H.C.,"Construction of Explicit Symmetric TVD Schemes and Their Applications", Journal of Compu- tational Physics, Vol. 68, 1987, pp 151-179. 15. Baldwin, B.S., Lomax, H.,"Thin Layer Approximation and Algebraic Model for Seperated Turbulent Flows", AIAA Paper 78-257 Jan. 1978. 16. Eriksson, L-E., Rizzi, A.,"Practical Three Dimensional Mesh Generation Using Transfinite Interpolation",VKI Lecture Series 1983-04 Vol. 1, Mar. 1983. 17. Zhou, L-D., Yuan, J-L.,"Calculation of the Turbulent Flow around the Stern and the Wake of a Bocly of Rev- olution", Proceedings, 15th Office of Naval Research Symposium on Naval Hydrodynamics, Hamburg, Ger- many, Sep. 1984 pp 279-290. 66 18 20 BL 8.6 Chen, H.C., Patel, V.C.,"Calculation of Stern Flow by a Time-Marching Solution of the Partially-Parabolic Equations", Proceedings, 15th Office of Naval Research Symposium on Naval Hydrodynamics, Ila~nburg, Ger- many, Sep. 1984 pp 505-522. 19. Huang, T.T., Cox, B.D.,"Interaction of Afterbody Bo- undary Layer and Propeller", Symposium on 'Hydro- dynamics of Ship and Offshore Propulsion System', Hovik Outside Oslo, Mar. 20-25 1977. Lee,Y.T.,Huang, T.T., Sung, C-H., "Flow Prediction about an Appended Body of Revolution", Gth Interna- tional Conference on Numerical Methods in Laminar and Turbulent Flow", Jul. 1989, Swansea, U.K. 21. Hartwich, P-M., Hall R.M.,"Navier-Stokes Solutions for Vortical Flows over a Tangent-Ogive Cylinder", AIAA Paper No. 89-0337, Jan. 1989 2 _ e i_ rrterxrg ~ ~ = - Att~~b~g 2 8.65 8.7 a.7s z/L I 1 1 8.8 e.8s 8.9 Fig. 1 Axisymmetric afterbodies 8.95 1 Fig. 2 Partial View of a CO Grid around a Body 84sr 8.375 _ 8.3 _ 8.225 _ 8.15 _ x x x Experiment - N-S Solution C'' -I x/L Fig. 3 Pressure Distribution on Afterbody 1

B.375 8.3 8.225 8.15 8.875 -. .875 -~.15 -8.225 _a ~ _ ._ -a.375 8.45 9,45 - x x ~c Experiment N-S Solution xT ;;C~8 _ C ~ ~ _ ~1 -8.45 1 1 1 1 ~ 1 1 1 1 1 1 _ -8 ~ 8 8.1 8.2 8.3 8.4 8.5 ~ 8.6 B.7 8.8 8.9 ~ 1 ~ 9.45 x/L Fig. 4 Pressure Distribution on Afterbody 2 O.1 s 0 0 0 Experiment N-S Solution B.375 S _ 0.3 _ B.22s 8.15 e.07s ;338 _ -C.875 _ -~.15 _ -3 .225 _D ~ _._ _ -a.37s _ x x x Experiment N-S Solution ~x. , , , 8 8.t 8.2 8.3 8.4 8.5 8.6 8.7 8.8 1 1 8.9 1 1.1 z/L Fig. 5 Pressure Distribution on Afterbody 5 (6a) Axial Velocity Profiles x/L = 0.846 ~/L = 0.914 ~/L = 0.934 x/L = 0.964 x/L = 0.977 _ ~ , ~,1 0 0.os o.oo 5 , 1 ~'t 1 1 ~7 0 1 0 1 0 1 0 1 0 1 ulV~ 0 0 0 Experiment (6b) Radial Velocity Profiles O. l 5 N-S Solution r/r = n 8~46 z/L = 0.914 =/L = 0.934 0.10 0.05 0.00 x/L = 0.964 z/L = 0.977 ~, ~ ~___ ~ ~ ~ ~ O.O 0.1 O.O 0.1 O.O 0.1 O.O 0.1 O.O . 0.1 ~Ur /V-go Fig. 6 Velocity Profiles on Afterbody 1 67

oO - ~: o :^ . - o - ~ 11 ~ - ~ ~ ~ - - . - o ~ u) P. CQ ~ l co o 11 o ~ - ~2 o l ~ o . . o o ~o o l ~ o . . o o t O U) - C~ o O :- ._ C~ o - ._ ~_ t _' . O ~ U) ~ Z O O O 0 0~2~- - O O o o o 7/(°~~) o o t 00 cg ~ ~ o '=o C 11 ~ ~2 ._ O _ _ O ¢ 11' ~ ~ U- VV--: , _' . CD o C~ 11 0 0 - 0 0 ~ ~ - ~8 co 1°1 - o 0 0 0 0 0 0 - , :: , _ ~ ~ cn _ ~ u, _ ~ Z oo o 1 1l ooooo oooo~ , ~2 0 0 0 0 O~ ' . o o T/(°.1 - .l) co G) o r~ o _~_ 11 ~ [ ° ° ° ° °°od o t G) ,~N /~^ /oOO - o 11 . o n 0 ~: ° _ n ~ 0 U) ~ cn X o o o O o o O U~, 68 :0 O O ' o oo ° o°~ 11 o 1 ~H 1 o ,~~ ~ 0 Ln 0 ~ ~ 0 0 . . . . O 0 0 0 7/(°d.1) 7/(°~~) c~ 0 0 f~ 0 ~ - ° c o 0 cn ~o ~ ~ ~ z 0 ~ 0 ~ ILl ~ P~, 0 .~ 0 0 - 0 ~o · - o o 0 o o c~ o o o o o O ;:s i o o u: o X o u) ~: o :^ ._ c; - 0 tl4 o o o o o

x/L = 0.954 a.06 L 0.06 .05 .04 D...01 _ =/L = 0.977 .05 .04 .03 0.01 _ x x x Experiment N-S Solution 1 x _ _ ~ _~ a l l l l 1 3 0.05 0.1 D.15 ~ 8.05 8.1 0.1S D.2 Au~ /VOO Fig. 9 Measured and Computed Axial Velocity Increase on Afterbody 1 Due to Propeller Suction ~ 4 a) top vlow ~//.~ y compahtlon //~ b) side view O.2 5 _ . _ _ t<)~.1 _ 0.05 x x Experiment N-S Solution . X 1 1 7 D.8 D.S 1 1.1 x/L Fig. 10 Measure and Computed.Surface Pressure Decrease on Afterbody 1 Due to Propeller Suction 1.C x/D = 6.0 .5 .,[ ~ 1 1 1 1 1 1 1 0 20 40 60 80 100 1 = 140 160 180 ~_ 1.C, .s . n 1 x/D = 3.5 _,ia' -.5- ' 1 1 1 1 1 1 1 1 1 0 20 40 60 ~0 1 00 1 20 1 40. 1 60 1 80 '~:x~=2.0 ~ _~ -.751 1 1 1 1 1 1 1 1 1 0 20 40 60 ~ 100 120 1" 160 180 ·7sr .25 y Fig. 11 Surface Flow Paterns - 25 __ _.7~ 69 ~_o ~ a-~ CZ} <experlmcnt computaffon ROD X 1 0 - ~ O 0.8 J~ O 3.0 ~ 1 1 1 1 1 1 1 1 ~ 20 ~o 60 so 1 oo 1 20 1 ~o 1 60 1 80 8, clrcumfcrenffal angle, dc~ Fig. 12. Measured and Computed Circumferential Pressure Distributions