Computation of viscous and inviscid flow around ship hulls

Introduction

The development of new ships, especially in the commercial field, still is a slow evolution-ary process. Series are generally small and time for development of a new ship is generally too short. There is no opportunity for modifications after intensive testing of a prototype, siiice the prototype is immediately delivered to the customer and has to meet specifications. Modifications to further ships in a series are often not possible, because the keel of the last vessel ha.s been laid before service experience of the first ship is available. These are the rea-sons for development in small steps to minimize technical risks. The implications for the field of ship hydrodynamics are that the next generation of ships will he modified only to a very limited extent. In addition, modifications will often not be investigated systematically due to limitations in time and budget. Fore- and aftbody, as well as propeller and overall dimensions are modified at the same time so that only the overall effect, but not the causes in detail, are known. CFD (computational fluid dynamics) methods, which are often cheaper and quicker than experiments and deliver more detailed information, will help the naval industry to over-come this problem in the future.

Inviscid free-surface flow

2.1. Introduction

Panel codes are commonly used in the aerospace and automotive industry to preoptimize aerodynamic characteristics before wind tunnel tests are performed. More recently there have been also reports of similar commercial applications in the field of naval architecture. The main difficulty in adapting panel methods developed in the aerospace industry to ship flow problems lay in the free water surface. Here a nonlinear condition lias to be fulfilled at an a priori unknown location - at the wave system elevation created by the ship. A first breakthrough was Daw.son's (1977) method even though only a linear approximation of the free-surface condition and other numerical shortcomings caused limited accuracy. It took another 10 years of research before the first solutions for the nonlinear problem overcame these shortcomings, Ni (1987,,), Jensen (1988), Kirn and Lucas (1990,). Research still progresses in this field as problems persist for strong flare, high-speed applications and sufficient accuracy for quantative p redic t io ii.

2.2. The wave resistance problem

For the wave resistance problem, water is considered to be incompressible, irrotational and inviscid. Surface tension is neglected. The ship's hull is assumed smooth. Appendages and propeller are neglected. Furthermore, we exclude breaking waves. These assumptions limit us in essence to displacement ships of Froude number F, < 0.4. Conventional cargo ships are not affected by this restriction.

Incompressible potential flow is governed by Laplace's ecluation for the velocity potential which holds everywhere in the fluid domain. Furthermore, for a unique solution of the potential and the a priori unknown position of the free water surface, we state boundary conditions on all boundaries:

- Water does not penetrate thewetted hull surface (Hull condition)

- Water does not penetrate the water surface (Kinematic condition) 'Bramfelder Str 164, 22305 Hamburg, Germany

Schiffstechnik Bd. 40 - 1993 / Ship Technology Research Vol. 40 - 1993 115

Volker Bertram, Jochen Laudan, HSVA'

TECHNISCHE UNIVERSIiEFI

Computation of Viscous and Inviscid Flow around SIij

mecharca

Archiof

Mekefweg 2, 228 CD De!t Tel: Ol575873-F _O15781CC3

Water does not penetrate the sea bottom (Bottom condition) Water does not penetrate the side walls (Side wall condition)

At the water surface there is atmospheric pressure (Dynamic condition) Waves created by the ship do not propagate ahead (Radiation condition)

Far away from the ship the disturbance caused by the ship has vanished (Decay condition)

- Waves pass unreflected through the artificial boundary of the computational domain (Open-boundary condition)

The ship is in equilibrium (Equilibrium of forces)

If the ship moves in unrestricted waters, the bottom and side wall conditions are covered by the decay condition.

The wave resistance problem has been tackled in the last loo years by both mathematicians and engineers. On top of the physical simplifications stated so far, the classical methods require further mathematical simplifications. The boundary condition at tlìe free surface is linearized

by assuming that the wave elevation is small compared to the ship's length and that the devia-tion of the local flow velocity from the parallel flow is small. This crude approximadevia-tion implies certain restrictions for the ship's geometry. The ship is assumed to be thin (width/length small) or slender (width/length, draft/length small). Unfortunately real ships usually do not meet these requirements. As a consequence, the agreement between experiments and computations is unacceptably poor, hindering commercial applications of these methods. Better computer equipment allowed the development of panel methods. Panel methods discretize surfaces where a boundary condition has to be enforced numerically into little segments (panels). This allows the treatment of arbitrary shapes and complicated boundary conditions. Appendix A gives some mathematical background on free-surface l)aIel methods and panel types.

2.3. Results

Most panel methods give as a direct result only the source strength of the discrete panels and in a next step the velocities at the collocation point. Bernoulli's equation then gives pressures and surface elevation again a.t discrete points. Integrating either pressure or wave elevations gives forces and moments which in turn might be used to determine dynamical trim and sinka.ge (squat). Selecting the maximal value of the normal velocity at the free surface or the residual of the nonlinear free-surface condition gives one scalar criterion of the quality of the iterative nonlinear solution.

Generally speaking, the very abstract information of the potential a.t a large miumber of collocation points (typically 1000) has to be aggregated to allow a civality check of the coin-putation and to derive information concerning possible hull improvements. An inadequate use of colour plots meant only for quality control in customer reports and presentations has led to nick-naming CFD colourful displays. Lack of experience in using the tool CFD lets engineers who are not CFD experts be happily snowballed with plots where in fact only 2 or 3 suffice for the expert to derive his recommendation.

Pressure distributions at the bow are usually quite accurate. A colour plot of interpolated isobars serves to identify critical areas of high velocities. These areas are usually at the forward shoulder resulting in deep wave troughs and at the bottom of the forebody where a reduction

of the critical area can reduce the viscous resistance. Interpolation gives a more realistic impression of the pressure distribution. However, time grid resolution determines the accuracy of the solution. To allow judgement. a plot of the hull grid should therefore always be included.

8

-8

Fig.1: Velocity plots give an indication of streamlines (Example: container vessel)

Fig.2: Wave patterns are mainly suited to check plausibility (Example: SWATH ship)

CW,6R + o 00 0. OOcg 0 ₀

0

Fig.3: Wave and residual resistance for container vessel 'Ville de Mercure'

0 experiment, +CFD old grid, CFD new grid

The wave profile along the hull shows how various wave systems interact. _{It can thus be} used to recommend hull changes for an improved interaction to reduce wave resistance. For clarity, the scale in the vertical direction is usually stretched. Again, interpolation often gives the illusion of high accuracy. However, wave profiles on the hull can be calculated at least by three methods: (1) for each panel on the hull the corresponding wave height can be calculated; if this height cuts the panel, a resulting line is drawn or its center marked; (2) the innermost line of collocation points at the water surface is plotted; (3) the water surface is extrapolated in

the y-directidn. A comparison of options (1) and (2) showed noticable differences only directly at the bow for Froude numbers higher than 0.2.

Velocity plots. Fig.1, give an indication of the streamlines and represent a numerical equi-valent to flow visualisation with tufts. Unproper grids tend to produce unreasonably large velocities on single panels which indicate unreliable results. Pressure plots, however, prove to

be an even more sensitive tool to discover such short-comings.

Other popular plots include isolines of the surface elevation showing the typical wave pattern. Fig.2. as well as pressure and velocity plots at the afthody or oblique views of the water surface. Although generally fasciiiating or pleasing to the eye all these plots are only suited to clìeck general plausibility of the computations and have little practica.1 value otherwise.

Nonlinear computations for real ships typically give differences of 25% in the wave resistance compared to linear computations, improving accuracy considerably. Although the pressure distribution at the bow is believed to be quite accurate, the wave resistance might still show errors in the order of 50% or more for common discretizations of 400 to 500 panels on the ship hull of a container vessel, unlike test computations for simple geometries such as the

parabolic Wigley hull which show excellent agreement with experiments. The reason for the still unsatisfactory accuracy for real ships is found in the high relative error in the integration, the viscous interaction in the aftbody and other residual resistance components. If the hydrostatic pressure is integrated to a force in the x-direction, it is found that iii many cases the resulting force has an order of magnitude of the wave resistance rather than being zero. By subtracting this force from the computed wave resistance the discretization error can be reduced. However.

it is still felt that Rankine panel methods with usual discretizations should not be used for quantitative resistance prediction. On the other hand, if forces on the forebody or pressure distribution for different bow shapes with comparable discretization are compared, valuable qualitative information can be obtained.

Generally, better grids improve accuracy. Trial computations for a modern container vessel ("Ville de Mercure") with about 400 elements on the hull gave errors in the wave resistance of about 100% for the lowest investigated Froude number which was the most realistic by today's standards, Fig.3. A new grid used about 25% more elements giving a finer resolution mainly on the bow. Furthermore, the new grid was generated automatically from a CAD surface giving a more regular (smoother) distribution of elements. Results were drastically improved.

The force in the vertical direction is not strongly affected by numerical error propagation as no difference of numbers of same magnitude is involved. Subsequently dynamical sinka.ge is predicted with high accuracy. Trim again is tricky as numerical errors and the influence of viscous effects lead to significant differences at higher speeds. Restricted water might amplify this effect to the point where the sign of the trim angle is predicted differently by computations and experiments. Finally, the maximum error in fulfilling the free-surface condition expressed either by the left-hand side of (6) or by the normal velocity at the free surface should be reduced by one order of magnitude in one or two iterative steps. The reduction of the error indicates if iterations succeeded well and can be regarded as a general quality control.

In conclusion, commercial CFD reports for inviscid flows should contain: - Information for quality control:

Plots of hull grid and water surface grid to help judge the accuracy of the interpolated results

Number of iterative steps and the error criteria at the free surface to judge how well the nonlinear solution was approximated

Plot of isolines of the free-surface elevation to see if implausible results appear at the outer boundary of the computational domain (e.g. reflections) or directly at the ship (usually near bow or stern)

Information for hull improvement:

Pressure isobars (preferably in colour) in views showing all relevant areas. Oblique views from below and above have been found to be valuable.

Wave elevation at the hull especially at the forward half of the ship showing to see the interaction of the wave systems. Information on how the wave elevation was obtained should be given.

Velocity distribution on the hull showing flow directions.

In case of comparisons of different hull versions. relative changes of the wave resis-tance showing if wave resisresis-tance has been improved

In special cases: dynamical trim and sinkage for example prediction of grounding in shallow water

All plots should include information about the ship geometry like sections and undisturbed waterline to allow comparisons to the lines plan.

3. Viscous Flow ('Navier-Stokes' Solver)

The inviscid part of the resistance accounts only for a fraction ( < 30%) of the total resistance. Knowledge of the viscous part is therefore of high interest for the evaluation of a ship hull geometry. Unfortunately all existing methods are unable to predict it with adequate accuracy.

However, such unwanted flow phenomena as separation, vortex generation and nonuniformity of the wake field are dominated by viscous effects. Therefore application of viscous flow codes makes sense, as qualitative insight of the flow is possible already today. Viscous investigations of the flow in the aftbody region serve in judging the propulsive properties and are used as input for propeller design.

The Navier-Stokes and the continuity equations are generally considered to be sufficient to describe in principle all real fluid physics for ships. Analytic solutions for the resulting system of nonlinear partial differential equations can only be achieved for some simple cases. For real ship geometries analytic solutions are impossible. This leaves only numerical methods of integration. Even if the influence of the free surface is neglected, full numerical solutions are still not possible even on the most powerful computers. For a ship speed of 20 knots the smallest eddies have a length scale of approximately 1irn and a fluctuation period of i0s. The computational domain covers approximately 106m3. To perform a meaningful time-average also over the largest eddies, the integration time has to be approximately lOs. This discretization of time and space leads to an extremely large nuniber of cells (10's to 1020, which cannot be handled in a reasonable time, Larssorz et ai. (1990).

Therefore 'Navier-Stokes' solver split velocities and pressures into a temporal mean and fluc-timation part to allow simpler numerical treatment of the equations. The resulting Reynolds-averaged Navier-Stokes equations (RANSE) need for a determinant solution additional equa-tions to describe the turbulence. A simple turbulence model uses algebraic relaequa-tions between turbulent R.eynolds stresses and mean velocities, Baldwin and Lomax (1978). This model has two advantages: No additional differential equation lias to be solved, and the model allows

direct enforcement of the no-slip condition on the hull. The (lisadvantage is that the displace-ment thickness is used as a length scale. This length is difficult o define in the aftbody region and wake. Nevertheless Hoekstm (1989) and Kodaina (1989) report successfull computations based on this turbulence model.

The commonly used k-r model is another turbulence model, Patankar and Spalding (197). The eddy vicosity is modelled by the turbulent kinetic energy k and its rate of dissipation r. k and r are described by partial differential equations with empirical constants. The method has two advantages: No length scale is needed and k and E are scala.r quantities, thus independent from the coordinate system. This allows the use of non-orthogonal coordinate systems. The disadvantage is that the model is only valid iii regions of fully turbulent flow which excludes the laminar sublayer along the body surface. Therefore the introduction of a wall function is necessary to substitute for the no-slip condition. Generally the logarithmic law for thin two-dimensional boundary layers is used as wall function coupling the wall-shear velocity with a nondimensional distance y from the wall. The introduction of the wall function allows a reduction of numbers of cells in the vicinity of the wall. A severe disadvantage of the combina-tion of k-r model and logarithmic law of the wall is the inahilitity to predict flow separacombina-tion at curved surfaces. ('hen and Pate! (1987) use the k-r model in the total computational domain with the exception of the cells lying at the wall. Here the no-slip condition is enforced by the introduction of an algebraic expression for the equation describing s. The results obtained with this model can be considered as good.

We use the commercial code STAR. CD choosing the k-r model for the eddy viscosity. Ap-pendix B gives more details of the equations and numerical techniques involved in solving the Reynolds- averaged Navier- Stokes equations.

4. Applications

Inviscid CFD computations were successfully used to determine the flow about a 280,000 tdw VLCC. For the low Froude number considered, wave resistance is predicted very inaccurately and has no informative value. However, flow details such as local pressures give valuable insight to the designer. The results of the first version were analysed to draw conclusions for hydrodynamic hull improvement. Results for the initial design and a modification were compared, Fig.4.

The following main changes were observed:

The low pressure (indicating high flow velocity) between bottom and sides near stations 18 and 19 is significantly reduced. The hull modifications in this area had the desired effect and should reduce viscous resistance.

The wave trough near the forward shoulder is reduced as is the low pressure area there. For the same VLCC project, viscous flow computations were performed for two hull ver-sions. Distributions of velocities and pressure in the aftbody were calculated. In this region potential flow calculation lose their validity due to the strong viscous effects. The cylindrical computational domain covered the region from 0.5L forward of the ship to 1.OL behind the ship. The cylindrical radius was 0.5L. The structured grid used 147 x 30 x 30 cells. Fig.5 shows part of the grid. A control computation with a finer grid using 147 x 30 x 50 cells showed no significant differences in the results.

The computations were performed for model conditions (model scale À = 32.2) using the k-s-model and the logarithmic law of the wall. The ship speed was V0 = 15.5 kn corresponding to a model Reynolds number Re = 1.17 IO. A computation for the real ship would require grids with considerably more cells if the nondimensional distance y of the innermost grid

Fig.4: Pressure distribution for two versions of a VLCC tanker. The low-pressure were reduced.

Fig.5: Grid for VLCC tanker (140000 cells)

Fig.6: Viscous pressure distribution in the aftbody region of original (left) and modified (right)

hull. lsoha.rs spaced by Cp = 0.1. The low-pressure areas Iìave been reduced by the hull

modification.

\

\\

_\\

\\

_\\

_\

\

Fig.7: Velocity distribution in plane 0.2 diameter forward of propeller: Lines of constant axia.l velocity and velocity vectors projected onto the plane .A vortex appears in the upper region of the l)ropel1eI circle.

point at the wall and the aspect ratio of the cells is kept constant. The value of y was chosen as approximately 1000. Ju and Patel (1991) show that even for y = 1400 the logarithmic law of the wall is valid.

Fig.6 shows the pressure distributions in the aftbody region for the original and the mod-ified hull forni. The 111111 modification reduced the low-pressure region at the aft shoulder considerably. This leads to smaller waves and reduced danger of flow separation. Fig. 7 shows the wake at a section 0.2 diameters forward of the propeller. A vortex appears in the upper region of the propeller circle. The full hull shape will certainly have separation which was not predicted by the computations. We see the reason in the combination of the k-E model aiid the wall function. A different turbulence niodel or wall function could improve results but would require grid refinements.

The investigation of varions bow forms of a container vessel were among the first practical applications of panel methods for real ships in Germany. Recently HSVA performed a CFD preoptimization of the foresliip of a slender container vessel (CB = 0.61). The wave elevation at the hull, Fig.8, implied that the concavity of waterlines at station 16 should be nioved forward.

-2 3 5 0 7 8 9 0 1 2 3 1 9 8 7 8 9 20

Fig.8: Wave elevation at container ship design (solid line) with CB = 0.61 and alternative design (dashed line) with CB = 0.58

Fig.9: Inviscid pressure distribution for container ship design (heft CB = 0.58, right C'B = 0.61) The uneven pressure distribution, Fig.9, implied a removal of the concavity of waterlines at the bulbous bow. The ship owner asked an external consultant for another design. He was wiffing to sacrifice some TEU for a lower power requirement. The subsequent design of block coefficient CB = 0.58 was investigated using a very comparable grid. Surprisingly, wave resistance increased and pressure distribution was clearly less favourable, Fig.9. Wave systems of bow and forward shoulder interacted rather infavourably heading to the significant trough at station 12. Fig.8. After reviewing the CFD results, the recommendations for the minor hull

modifications of the better design were accepted and incorporated in the line drawings before model tests commenced.

References

AANESLAND, V. (1989), A hybrid model for calculating wave-making resistance, 5th lut. Conf. Nurn. Ship Hydrodyn., Hiroshima

ANDO, J.; NAKATAKE, K. (1988), A method to calculate wave flow by Rankine source, Trans. West Jßpan Soc. Naval Architects

BALDWIN, B.: LOMAX, H. (1978), Thin layer approximation and algebraic model for separated tur-bulent flows prediction code, AIAA paper 78-257

BERTRAM, V.; JENSEN, G. (1987), A new approach to non-linear waves generated by a body moving steadily at a free surface, IUTAM-Symp. Non-linear Water Waves, Tokyo

BUSCH, S. (1990), Numerische Berechnung des Wellenwiderstands eines schnellen Sportbootes, ilS-Diplom thesis, Univ. Hamburg

CAO, Y.-S.; SCHULTZ, W.W.; BECK, R.F. (1990), Three-climensioiìal insteady computations of tian-linear waves caused by underwater disturbances, 18th Symp. Naval Hydrodyn., Ann Arbor

CHEN, H.C.; PATEL, V. (1987), Practical near-wall turbulence models for complex flows includiug separation, AIAA paper 87-1300 DAWSON, W.C. (1977), A practical computer method for solving ship wave problems, 2nd Tnt. Conf. Nurn. Hydrodyn., Berkeley

DIMIRDZIC, I; GOSMAN, A.D., ISSA, RI. (1980), A finite-volume method for the prediction of turbulent flow in arbitrary geometries, 7th Tnt. Conf. Num. Methods in Fluid Dyn., Stanford

HESS, J.; SMITH, A.M.O. (1962), Calculation of non-lifting potential flow about arbitrary three-dimensional bodies, Douglas Aircraft Division Report No. E.S.40622

HOEKSTRA, M. (1989), Recent developments in a ship stern flow prediction code, 5th Tnt. Conf. Nurn. Hydrodyn., Hiroshima

ISSA, RI. (1986), Solution of the implicitly discretized fluid flow equations by operator-splitting, J. Comp. Physics 62

JENSEN, G. (1988), Berechnung der stationären Potentia!strörnung um ein Schiff unter

Berücksichtigung der nichtlinearen Randbedingung an der freien Wasseroberfläche, IfS-Report 484, Univ. Hamburg

JENSEN, G. (1991), A panel method using numerical integration, 7th GAMM Seminar Num. Techn. for BEMs, Kiel

JENSEN, G.; MI, Z.-X.; SÖDING, H. (1986), Rankine source methods for numerical solutions of the steady wave resistance problem, 16th Symp. Nava! Hydrodyn., Berkeley

JENSEN, G.; BERTRAM, V.; SÒDING, H. (1989), Ship wave-resistance computations, 5th lut. Conf. Num. Ship Hydrodyn., Hiroshima

JU, S.; PATEL, V.V. (1991), Stern flows a.t full-scale Reynolds numbers, J. Ship Research 35/2 KIM, Y.-H.; LUCAS, T. (1990), Nonlinear ship waves. 18th Symp. Naval Hydrodyn., Ann Arbor KODAMA, Y. (1989), Grid generation and flow computation for practical ship hull forms and propellers using the geometrical method and IAF scheme, 5th Imt. Conf. Num. Hydrodyn., Hiroshima

LARSSON, L.; BROBERG, L.; KIM, K.-J.; ZHANG. D-H. (1990), A method for resistance and flow

prediction in ship design, SNAME

MUSKER, A.J. (1988), A panel method for predicting ship wave resistance, 17th Symp. Naval Hydro-dyn., Den Haag

NAKOS, D. (1990), Ship wave patterns and motions by a three-dimensional Rankine panel method, Ph.D. thesis, MIT, Cambridge

NAKOS, D. (1991), Transverse wave cut analysis by a Rankine panel method, 6th Tnt. Workshop Water Waves and Floating Bodies, Woods Hole

NI, S.Y. (1987), Higher order panel methods for potential flow with linear or nonlinear free surface boundary conditions, Chalmers Univ., Goteborg

OGIWARA, S. (1988), Numerical solution for steady ship waves by Rankine source method - A review of current studies, 19th ITTC R&FC., ist Meeting, A.R.E. 88. 5.24-25

PATANKAR, S.V.; SPALDING, D.B. (1972), A calculation procedure for heat, ma.ss and momentum transport in three-dimensional parabolic flows, hit. J. Heat Mass Transfer 15

PIERS, W.J. (1983), Discretization schemes for the modelling of free water surface effects in first-order paiiel methods for hydrodynarnic applications, National Lucht- en Ruimtevaartlah. NLR TR 83093 RAVEN, H.C. (1988), Variations on a theme by Dawson; recent improvements of potential flow calcu-lation method for ships, 17th Symp. Naval Hydrodyn., The Hague

RAVEN, H.C. (1990), Adequacy of free surface conditions for the wave resistance problem, 18th Symp. Naval Hydrodyn., Ann Arbor

RAVEN, H.C. (1992), A practical nonlinear method for calculating ship wavemaking and wave resis-tance, 19th Symp. Naval Hydrodyn., Seoul

SPALDING. D.B. (1972), A novel finite difference formulation for differential equations involving both first and second derivations, mt .i. Num. Methods Eng. 4

THOMPSON, J.F., WARSI, Z.U.A., MASTIN, C.W. (1985), Numerical Grid Generation, North-Holland

WEBSTER, W.C. (1975), The flow about arbitrary three-dimensional smooth bodies, J. Ship Res. 19/4

Appendix A: Potential free-surface flow

The equations are formulated in a Cartesian coordinate system with z pointing downwards. The velocity potential meets Laplace's equation in the whole fluid domain:

The no-penetration condition ou the hull, bottom and side walls can be expressed as

where is the unit normal vector on the respective surface. The no-penetration condition on the sea bottom is easily enforced for a flat bottom by using mirror images of all panels with the sea bottom as plane of reflection.

In the special case of a ship with transom stern we assume the transom to he dry. At the edge zt of the transom stern we require the pressure to be equal to the pressure at the free surface. For flow dominant in x-direction this leads to:

= -LJ1/i -

The dynamic condition at the free water surface z = ç is derived from Bernoulli's equation:

- gz = U2 ₍₄₎

Time kinematic free-surface condition gives at z =

Vç' = (5)

where the index indicates a partial derivative (For simplification we write (x, y, z) with &(/&z = 0.) The unknown surface height ç' is eliminated from these two boundary conditions, resulting at z_{= ( in:}

VV(V)2

This free-surface condition is nonlinear. An important quality criterion of a method for computing the wave resistance is the approximation of this condition. Cla.ssical methods linearize the difference between the real flow and the uniform flow (Kelvin condition). Panel methods adopt this very crude approximation only for reasons of numerical comparison, e.g. Raven (1990), or as an initial value for further iterative solutions. Dawson (1977) introduced a linearization about the double-body flow. How-ever, his often copied formulation is faulty, Jensen (1988), Raven (1988). Double-body linearization is used e.g. by Aanesland (1989), Nakos(1990) and Raven (1990). Musker (1988) and Ogiwara (1988) present nonlinear approximations for (6). However, these methods should not he confused with methods that iteratively solve the exact condition (6), like Ni(1987), Jensen (1988,), Kim and Lucas (1990,), and Raven. (1992). A correct solution procedure fulfills repeatedly a linearization about an arbitrary ap-proximation of the potential at an arbitrary approximationZ of the free surface elevation. Consistent linearization leads to the linearized free-surface condition at z= Z, Jensen eI al. (1985,:

VV [(V)2 + V

V} + VV(V)2

-ri

4[(VF)2+2V.vU2]gZ

+ [VV(V

- g - .

=0

(7)

The denominator in the last term is zero mf the vertical particle acceleration is equal to gravity acceleration g; this is a stability limit of the approximate flow . A higher acceleration would violate our initial assumption of no breaking waves and cause a break-down in the iteration.

The decay condition states that far from the body the flow field tends towards uniform flow: hm VØ = (LT,0,0)

x2y + z2 -co

The radiation condition must be enforced by a numerical specialty. Most panel methods use special upstream finite-difference (FD) operators for second derivatives of the potential in (6). Dawson's (1977,) 4-point FD-operator introduced unwanted dissipation and for usual discretization underpredicted wave lengths by about 5%. Pzers (1983,), Nm (1987) and K2?n and Lucas (1990,) give improved FD-operators using 3 to 5 points. Another improvement is achieved by the use of splines as by Musker (1988) and

Nakos (1990). We use the completely clifferemit approach of staggered grids, Jensen ei al. (1986,.). An extra row of sources downstream and an extra row of collocation points upstream effectively enforce the radiation condition. For equidistant grids this can also he interpreted as shifting the source grid versus the collocation point grid. Amido and Nakalake (1988) also show that this technique works well using Hess&Smith type panels. Its advantage over FD-operators lies in simpler programs and higher accuracy as demonstrated in 2-dimensional test cases.

The open-boundary condition must also be specially enforced. Methods based on FD-operators use 2-point operators at the downstream end of the grid which strongly damp the waves. At the upstream end, no waves should appear. Here often a simple condition like zero deviation from uniform flow speed in the x-direction is applied. Nakos (1990) sphine-based method requires that waves do not reach the outer boundary of the grid. This leads to rather large grids at the free surface. Aanesland (1989,) tries to reduce the size of the grid by coupling the near-field panel solution to a far-field thin-ship theory solution. However, computation time was not significantly reduced whereas complexity of programming increased considerably. The staggered grid technique of Jensen et al. (1986) is the most elegant solution. Waves leave the computational domain at the downstream end and the outer houiidary without reflection. A only disadvantage is the need to analytically compute derivatives of the potential of one order higher than for FD or spline methods. However, the source distribution is smooth which is not the case for techniques using 2-point FD-operators, Ando and Nakatake (1988).

The forces and moment acting ou the body are determined by integrating pressure over the actual wetted hull surface. An alternative to pressure integration is the analysis of the wave energy behind the ship. Busch (1990) and Nakos (1991) used transverse cut wave analysis successfully hut the computation time is large due to the large grid necessary behind the ship. Busch (1990,) shows that wave cut analysis gives wave resistance coefficients that correctly tend towards zero for small Froude numbers. Integration of pressure gives a rather constant finite value for small F, which is anindication of a numerical discretization error.

The equilibrium floating position of the body is determined from the condition that the forces and moments on the body are zero. Forces and moments considered are those due t.o the potential fluid pressure, the weight, and the towing force.

(8)

Panel methods discretize surfaces where a boundary condition has to he enforced numerically into little segments (panels). Each panel acts as a distributed Rankine source. Panel methods therefore fulfill automatically Laplace's equation and the decay condition. Form and source strength distribution differ for various panel codes:

- point source

It is used by almost every panel code as an approximation if the distance between panel and point on the boundary is sufficiently large.

- point source cluster

A point cluster consists of sveral point sources of equal source strength. This approximation by Jensen et at. (1989,) gives good results at the water surface in short CPU time.

- plane quadrilateral, constant source strength

This panel type was developed by Hess and Smith. (1962) and is commonly used, e.g. Aan.es-land(1989), Musker (1988,), Nakos (1990) and Raven (1990).

- plane triangle, bilinear varying source strength

This panel type was developed by Webster (197.5,,. Personal experience proved it to he unsuitable for ship geometries, Bertram and Jensen (1987,).

- arbitrary polyeder, constant source strength

This type is used by us since Jensen (1988,). It is described in detail in Jensen (1991). - curved panel with bilinear varying source strength

This higher order panel gives higher accuracy for the same number of unknowns and should be suited to approximate also second derivatives of the potential at the panel center, Ni (1987,), Kim

and Lucas (1990).

We use the arbitrary polyeder as panel type for the hull with a collocation point at its center. This panel type has been found to be very flexible in discretizing ship geometries. The integrals for the induced velocity have a singular integrand that precludes straightforward numerical integration. Here we apply a special mathematical trick: we deduct a simple function with the same singular behaviour, namely the effect of a source distribution of constant strength on a sphere that touches the hull in the panel center/collocation point. The radius of the sphere should be as large as possible but with the center of the sphere still lying within the hull. The residual function is simply evaluated by numerical integration. Desingularization, i.e. a slight removal of panels from a slightly curved surface improves the numerical accuracy, Cao, Schultz and Beck (1990). For hull surfaces this technique can not be recommended due to the usually sharp bow form. However, for the free water surface it has been proved to he advantageous,

Jensen et al. (1986), Musker (1988), Raven (1992).

Appendix B: Navier-Stokes solver

The continuity and Na.vier-Stokes equations for incompressible flows are written in a Cartesian coordi-nate system using tensor notation:

1. continuity equation

2. Navier-St.okes momentum equation

3u ôuztj i dP

+fi

dt ôx p òx d.cjò.v3

where u is the velocity, p the pressure, u the kinematic viscosity and f a body force.

Splitting the velocity into a temporal mean and fluctuation part allows a simpler numerical treatment. of the Navier-Stokes equation:

u =)i.+u

p=p+p'

(9)

(10)

where the bar denotes the mean value and the prime the fluctuation. Equations (8) and (9) then become:

axt (12)

&u auu i a ô2u

+V,

+f

ut ox1 p ox ox uxux1

The difference between these Reynolds-averaged Navier-Stokes equations (RANSE) and the equations

for laminar flow lies in the Reynolds stresses puu. For a determinant solution, further equations

for the Reynolds stresses have to he introduced. A simple turbulence model uses an algebraic relation between Reynolds stresses and mean velocities:

5u Du1 2

pu =

+ ) +

with the turbulent kinetic energy

U' ti'

the eddy viscosity Vt, and the Kronecker symbol 6jj. The dependence of the eddy viscosity on the location is described by an algebraic empirical function, Baldwin and Lomax (1978).

In the k-e model the eddy viscosity is modelled by the turbulent kinetic energy k and its rate of

dissipation e:

0. 09k2

Vt

with 0.09 as empirical constant. k and e are described by two partial differential equations:

dk

= () + ¿'t( + )-- C

dt Ôxj 1.ODxj dxj Ôxj 3x1

de Ô u De 1.44ve Ôu Du Dv.t 1.92e2

=(---)+

(_±L)__

(1)

dt Ôx 1.2 3x k àxj Ôxj Òx k

These equations contain four further empirical constants (1.0, 1.20, 1.44 and 1.92).

The k- model is coupled to a wall function. Generally the logarithmic law for two-dimensional turbulent boundary layers is used as wall function:

1

ln(9.0y)

u1. 0.42

with 0.42 and 9.0 as empirical constants, y a nondimensional distance from the wall and u the

wall-shear velocity.

The differential equations describing conservation of mass and momentum are discretised in finite vo-lumes, Dimirdzic et al. (1980) and Spalding (1972). Finite volume methods integrate over a single cell before dependent variables are approximated by the cell-centered nodal values. The RANSE can l)e described by the general transport equation:

+ V (p

fV) = s

where th is one of the dependent. variables u, k, and e. F are the corresponding diffusion terms and

s the source terms. The transport equation for stationary flows consists of three parts: convective

flow, diffusive flow and source. Gauss' law transforms the integration of the convective and tue diffusive terms into integrations over the surface of the cells. For the general transport equation the convective

flow is approximated over a cell volume by:

b j=i

(p)j

Here j represents the average value of the variable at the cell face aiea S. The mass flows F are assumed to he known (i.e. the values of the previous iterative step are used). The values for j

have to be expressed by the unknown values at the cell centers (cell nodes). This can be done by various interpolation schemes. The simpLest case uses the value of the next cell node upstream (upwind

differencing). This interpolation scheme ensures that physical limits for are not trespassed, but leads very often to numerical complications. The linear upwind differenciug method uses the nodal values of the next two upwind cells. The central differencing interpolation method uses the two closest nodal values regardless of the flow direction. The two latter methods are numerically more stable hut might give values with no physical meaning. We used in our computations linear upwind differencing interpolation.

For the computation of the diffusive terms the gradient at the cell side has to be expressed in the derivatives of in all three coordinates. For non-orthogonal grids this adds linear terms for S which are treated explicitly in solving the system of equations. This leaves only the values of the cell center and its 6 immediate neighbors as unknowns in the discretized equations. The right-hand side of the equation is integrated over the volume. lt is assumed that the source term s at the cell center represents a mean of the total cell. The additional stress ternis in the momentuni equations represent. an exception. They are a divergence of a vector and are integrated over the cell sides just as the diffusive terms. The fluai form of the integrated and discretized transport equation is:

App +

The coefficients A cover the contributions of the implicitly discretized parts of the convective and diffusive terms and Q covers the source terms and the explicitly treated parts of convection and diffusion. The system of equations can he solved efficiently in an iterative way.

The vector equation for the conservation of momentum gives three scalar equations. These iteratively determine the velocities for a given pressure which is estimated in the beginning. Usually these velocities do not fulfill the continuity equation. The introduction of a coupling term for the pressure allows a correction of pressure and subsequently velocities. Commonly used for coupling pressure and velocities is the SIMPLE method, Palankav and Spa/ding (1972,) and the PISO method, Issa (l986. The SIMPLE method is fast hut has numerical instabilities for suboptimal grids. The PISO method is numerically stable but increased CPU time by a factor of 5 in one case. As CPU time is neither a critical time or cost factor in practical CFD computations, we recommend the use of the PISO method.

In addition to the governing equations of the fluid, boundary conditions have to be formulated. The upstream boundary plane of the computational domain is defined as the inlet. The x-componerìt of the velocity is set to negative ship speed. The other components are zero. The turbulent energy k is estimated to k 0.001 as better experimental data are not available. The dissipation rate then

follows from the semi-empirical relation

Schiffstechnik Bd. 40 - 1993 I Ship Technology Research Vol. 40 - 1993 129

A plausible value for the turbulence length I is the mixing length which depends ou the distance y from the wall:

1

(24)

At the inlet I is set to approximately 0.01 times the radius of the cylindrical computational domain as no distance from the wall can he defined here.

The downstream boundary plane of the computational domain is defined as the outlet. The gradients in the x-direction of all variables are set t.o be zero at the outlet. The downstream bouiidary plane, the outer boundary, the center plane and the water surface are treated as planes of symmetry, i.e. the gradients in the normal direction of all variables are zero. On the ship hull the wall function substitutes for the no-slip condition.

Time and costs involved in grid generation for RANSE calculations is often underestimated. In a first step grids are generated on the boundaries of the computational domain. The computational domain is then divided into planes normal to the x-axis of the ship. Grid generation within the planes follows

Thompson eÍ aI. (1985,). The elliptic grid generation method solves the Poisson equation V2() = P. P are control functions, influencing distance and orientation of coordiiiate lines. The generated grid lines are interpolated in a second step using a spline function. This ensures that for highly convex sections the distance of the innermost cells from the wall (loes not become too large.

0.09314k312

KAVITAT1ON

Ein wichtiges Standardwerk

zur Schiffbautechnik

In 3., erweiterter Auflage

KAVITATION

von Dr-Ing. W-H. Isay,

Professor für Strömungsmechanik

am Institut für Schiffbau der Universität Hamburg

441 Seiten mit 209 Abbildungen, Format DIN A5

DM 64,80, inkl. MwSt. zuzüglich Versandkosten.

Best.-Nr. 726

Die nun vorliegende dritte Auflage ist durch Einarbeitung der wichtigsten in den vergangenen fünf Jahren erreichten Forschungsergebnisse nochmals erweitert und ergänzt worden.

Zum Inhalt: Zunächst werden die strömungsmechanischen und thermodynamischen Grundlagen der

Kavitation in einer Keime enthaltenden Zwei-Phasen- und Zwei-Komponenten-Strömung sowie Fragen der Keimbildung und der Stabilität, des Wachstums und Zusammenfalls von Keimen und Blasen verschiedener

Formen (auch Porenkeime und Blasensysteme) behandelt. Weitere Kapitel sind den Problemen des Kavitationseinsatzes und dem Einfluß der Strömungsgrenzschicht, der Oberflächenrauhigkeit und der Turbulenz gewidmet sowie der Vorausberechnung der Kavitationserscheinungen an Strömungskörpern,

insbesondere PropeIlerflügen. Die Bedeutung der Kavitation (u.a. auch in Spitzenwirbeln) für die von

Propellern an der Schïffsaußenhaut induzierte Vibrationsbelastung und die verschiedenartigen

Maßstabsef-fekte zwischen Modeilversuchen mit Propellern und den Criginalbedingungen werden abschließend erörtert.

Dieses Buch ist aus Manuskripten der Vorlesungen über Kavitation entstanden, die der Verfasser seit 1976 regelmäßig im Fachbereich Physik der Universität Hamburg abhält. Auch die Vorlesun-gen über Propellertheorie im Institut für Schiffbau enthalten wesentliche Teile des in dem Buch verarbeiteten Stoffes.

>

Ich/Wir bestelle(n) E per Rechnung E per Nachnahme E Scheck anbei

über die Buchhandlung Exemplare

KAVITATION, DM 64,80

zuzüglich Versandkosten, inkl. Mehrwertsteuer Name/Firma

Anschrift

Datum/U ntersch rift

X

Schiffahrts-Verlag ,,Hansa"

C. Schroedter & Co. Postfach 520365 22593 Hamburg