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.
Towards the Simulation of Seakeeping and Manoeuvring based on the Computation of the Free Surface Viscous Ship Flow A. Cura Hochbaum and M. Vogt (Hamburg Ship Model Basin, Germany) ABSTRACT The simulation of ship motions due to rud- der manoeuvres and/or to incoming waves based on free surface viscous flow computations will be one of the most fascinating challenges in ship hydrody- namics for the coming years. The present paper deals with the steps we are taking to achieve this long-term objective. The numerical method used is described, predicted forces and moments on the hull of a ship in a variety of cases involving ba- sic manoeuvring flow situations are compared with measurements and results of simulations for ships moving against harmonic waves coming from ahead are presented. INTRODUCTION There have been research activities in this field since the 21th ONR Symposium in Trond- heim 1996 and the 7th Int. Conf. on Numerical Ship Hydrodynamics in Nantes 1999. McDonald et al. [6] and Takada et al. [12] simulated manoeuvring based on viscous flow computations. Wilson et al. t16], Kinoshita et al. t5] and Sato et al. [9] simulated ship motions in waves coming from ahead. Miyake et al. t7] also treated oblique incoming waves. They all generated regular waves prescribing velocities at the inlet according to the potential wave theory. Because the techniques used in the publica- tions mentioned above have shown to be eEective and promising, we extended our numerical method in a similar way in order to simulate the motions of a ship model moving straight ahead against incom- ing waves. Due to the assumed symmetry, the ship Ship motions and loads are traditionally cat- only performs surge, heave and pitch motions. An culated using methods based on the strip theory or accurate computation of forces and moments on the with panel methods. These methods rely on po- hull during steady oblique and turning motions, is tential flow theory and take viscous effects into ac- certainly a prior requirement for a reliable simula- count in a very coarse manner. For this reason, motions with a significant viscosity influence, e.g. surge and roll, can hardly be predicted well. More- over, most of the methods are not able to deal with steep waves nor with breaking waves, especially be- cause they linearise the motions and/or the water free surface, impeding the treatment of cases where slamming and green water occur. The long-term objective of this work is to achieve numerical simulations of ship motions based on the calculation of the free surface viscous flow around a manoeuvring ship in calm water or in ocean waves. For most tasks concerning seakeeping and manoeuvring, conventional methods will still be used in future. However, the present method may become useful for the special cases mentioned above, because it takes all relevant aspects involved into account. For the present first step however, we have focussed on much simpler cases than the really interesting problems in practice. tion of motions of a manoeuvring ship. In order to carry on the validation of our code for these cases, extensive work is being performed parallel to the simulation of motions in waves. This includes the computation and comparison of results with model test measurements for a large combination of yaw rate, drift angle, rudder angle and heel angle. Re- sults for zero rudder and heel angle are presented at the end of this paper. METHOD Coordinatesystem In order to describe the motion of the ship we define the earth fixed coordinates X, Y. Z and the ship fixed coordinates x, y, z as well as the Eu- ler angles ,o, I, ¢. The origin of the ship fixed coordinate sys- tem is positioned mid ship, on the symmetry plane
X,I x,~ f Y ~ Projection of x ~ ~ ~ b_\ -Y,J Z k ~ , _ v Z K m[(u+w2wTV) XG (W2 +W3) YG (W3~1~2) + ZG (~2 + ~1 W3)] = Fit m [(v~1 W + W3 U) + XG (W3 + w1 ¢~2) YG (~1 + ~3) ZG (~1 ~2 W3)] = FY m [(w + ~1 v~2 U) XG (~2~1 W3) + YG (W1 + ~2 ~3) ZG (W1 + ~2)] = FZ IX ~1 - (IYIN ) W2 ~3IXY (W2~1 W3) IT (~2 _ ~3 )IXZ (~3 + ~1 W2) +m[yG(W+~tlVMU) Figure 1: Coordinate systems ZG (V~1 W + (1)3 U) ] = K at the undisturbed free surface. The hybrid base i, r', u and the projection of the x-axis on the X' O Y' plane serve to define the rotations and Euler angles cp (roll angle), ~ (pitch angle) and ~6 (yaw angle). The unit vector ~ is defined as the product of K and ~ for K x ~ 76 0, see Fig. 1. With this def- inition of the rotations the position of the ship is unique, independent of the order of rotations. The components of the angular velocity in the ship fixed coordinate system are: we = ~~ sine . . ~2 = %6cos~sin~+§cos~ (1) . . w3 = ~ Cost COsy'~ sing The ship fixed and the earth fixed coordi- nates are related to each other as follows: (, Y \) = (\ YO \\) + T ~ y ~ (2) where XO, YO, ZO are the coordinates of O and T is the transformation matrix. The elements of the rows of T are the components of the earth fixed unit vectors expressed in the ship fixed system. Motion equations of the ship The motion of the assumed rigid ship is de- scribed by the momentum and the moment of mo- mentum equations. These equations are written in the ship fixed coordinate system: (3) IY ~2 - (\IZIX ~ Wit ~3 IXY (<W1 + ~2 W3 Iy: (w3~1 W2)Ixz (~3W +m f ZG (A +~2 wW3 v) He (w +~ vw2u) ~ = M Izw3 - (`IxIy~wlw2 IXY(~2_~)2) Iyz (W2 + ~J1~3)I=z (~1~2 ~3) + m [L XG (V~1 W + W3 U) YG(U+~2W~3V)] = N (4) m is the mass of the ship and XG' YG, ZG are the (constant) Cartesian components of the center of gravity. a, v, w and At, w2, W3 are the Cartesian velocity components of O and the angular velocity components respectively and u, v, w and At, w2, we their time derivatives. On the RHS of Eq. (3) are Tic (longitudinal force), Fy (side force) and Fz ('vertical' force) the components of the resulting external force in the ship fixed system. Id, IY, IZ are the moments of inertia and Icy, In, Iyz the moments of deviation for the axis through O. On the RHS of Eq. (4) are K (roll moment), M (pitch moment) and N (yaw moment) the components of the moments of the external forces. The state of the moving ship is defined by the position of O (earth fixed coordinates XO, YO, ZO), its velocity (u, v, w), the Euler angles be, it, ¢) and the angular velocity All, ~2, ~3~. To solve Eqn. (3) and (4) we need the forces and moments on their RHS. They are written here:
Fat = For + mg cosoe1 + Fop Fy = Fyh + mg costar + Fyp Fz = Fzh + mg COS(X3 + Fzp = Kh + mg (COSC\3 YGcostly ZG) + KP (5) M = Mh+ mg (cosa1 ZGCOSCX3 XG) + Mp (6) N = Nh + mg (cosot2 XGCOSOY1 YG) + Np The directional cosine of the gravity force, in the ship fixed coordinate system, costar, i = 1, 2, 3 coincide with the components of the earth fixed ba- sis vector K and as such, with the elements of the third row of the transformation matrix T. By inte- grating the pressure and shear stresses on the hull we get the hydrodynamic forces Fish, Fyh' Fzh and moments Kh, Mh, Nh in Eqn. (5) and (6) where 9 is the acceleration due to gravity. Fop, Fyp' Fzp and Up, Mp, Np are the propeller contributions to the forces and moments and are here approximated with a simple model. Equations of motion of the fluid The fluid flow around a moving ship is con- sidered to be a flow of one fluid with two immisci- ble phases (water and air). The interface between both phases represents the free surface. Because each phase is incompressible, the global flow can as such be handled as incompressible. The governing equations are the RANSE and continuity equation. Reynolds stresses are approximated with an eddy viscosity modell. In the Cartesian ship fixed coor- dinate system xi=x, x2=y, X3=Z the conservative form of the governing equations can be written: dvi + 0(vivj) _ 1 0(p + 2/3 r: k) (7' lit {Jxi rat Ski + F 2 + Ski [( Rn ( fJ /; ~ ) (~i + (ijk (l)j Ok) (ijk DJj~km7' Am X 2 (ijk Wj ok f ijk JO X ~ i = 0 (8) All variables have been nondimensionalised by the ship speed Uo, the ship length L and the wa- ter density Pw vi is the relative velocity component in the Cartesian direction i, xi the corresponding Cartesian coordinate and t the time. Fin = To// is the Ffoude number and Rn = Uo L/Z/W is the Reynolds number. p is the mean pressure and ~ the turbulent kinetic energy. The eddy viscosity zig = k/w where w is the specific dissipation rate of ~ is determined by using the kw turbulence model from Wilcox t154. \Vall functions were used in all computations presented in this paper. r, = P/Pw and r2 = ',/Z'w are the ratios between the local den- sity and kinematic viscosity and the corresponding values for water. These ratios have constant values P'/Pw and i1/uw in air, pi and i'` being the density and viscosity of air, and both equal one in water. f ijk iS the permutation symbol. Surface tension has not been taken into account in Eq. (7~. The terms in the third and forth row of Eq. (7) are the inertial forces stemming from the ac- celeration and rotation of the coordinate system used. The first term represents the acceleration of the origin, the second and third the centrifugal and the Coriolis forces respectively. The last term is identical to zero for steady motions e.g. steady oblique motions. Different from our earlier work, where the motion of the ship was given, the inertial forces now have to be updated in the course of the simulation. The transport equations for the turbulent quantities k and ~ for a incompressible fluid with two immiscible phases in Cartesian coordinates are: 0t + fJ2;i foci [(Rn 0~£ ] + (0Vi + i9Vj) Levi _ id* k 0t + deli = hi [(Rn +a Gil (10) ~ (49~i + Ski) If - jB,-4,2 a* = ~ = 0.5 it* = 0.09 ,3 = 0.075 By = 0.555 The transport equations for k and ~ Eqn. (9) and (10) and the continuity equation Eq. (8) remain without changes valid in the ship fixed co- ordinate system. These equations and the momentum equa- tions Eq. (7) are solved in the whole computa- tional domain. Special treatment of the two phases
and boundary conditions on the free surface is not needed. The inertial forces in Eq. (7) are calcu- lated with the velocities pi and wi and their time derivatives obtained from Eqn. (3) and (4~. Level Set Method The free surface is taken as the interface be- tween the two phases water and air. This is cap- tured with the Level Set method see e.g. Osher and Sethian t8] and Sussman et al. t114. Among the first ship hydrodynamics related applications of the Level Set method was reported in t133. The interface can be expressed by the equa- tion Phi, t) = 0. Since the function o is always zero in all points on the interface its material deriva- tive is zero there. This also applies for compli- cated topologies e.g. breaking waves. In the Level Set technique, a scalar field ~ having a zero level set ~ = 0 coinciding with the interface is defined throughout the computational domain. We have air where ~ is negative and water where o is posi- tive. Extending the condition ~ = 0 to all points in the domain leads to the transport equation for the level set function: do + 0(Vi ¢) o (11) In every point (cell center) of the compu- tational domain, the initial value of o is chosen to be the distance (with sign) to the initial posi- tion of the interface. According to Eq. (11), o is exclusively transported by convection and always remains positive in one phase and negative in the other. For the case of water and air the density is P = Pw (water) for ~ > 0 and p = pi (air) for 0 < 0. Also the molecular viscosity of the fluid changes its value suddenly when ~ changes sign. In order to avoid numerical problems, these discontinuities are somewhat smoothed and the density and the viscosity are determined as follows: D Pw V LAW = (1c) Pi + c Pw 1 o 2 [1 + sir (26)~ > or (water) <or (air) of _ 0 < ~ The thickness of the transition region can so, be prescribed with the parameter or. As long as ~ is the nondimensional distance from the considered point to the interface, this thickness is 2 or, as can be seen in Eq. (12~. We choose or to be roughly the vertical grid spacing near the waterline of the ship. During the time marching procedure, first the momentum equations, then the pressure-cor- rection equation (mass conservation) followed by the transport equations for the turbulence param- eters are solved assuming a known o distribution. After that a new approximation of ~ is computed. COMPUTATIONAL GRID The quality of the numerical grid used is cru- cial for the accuracy of the results and for the con- vergence behaviour of the method. While the loss of accuracy due to a poor grid can (in principle) be reduced increasing the grid resolution, conver- gence problems mostly remain. Usual grid require- ments are smoothness and the proper resolution to capture all interesting aspects of the flow. The latter usually only depends on the Reynolds num- ber (boundary layer, wake) and the Froude num- ber (free surface). However, in order to calculate ship motions using a ship fixed grid, the expected motions themselves have to be taken into account when generating the grid. The free surface (includ- ing incoming waves) changes the position in the course of the simulation. Even small pitch motions for instance, lead to large vertical shifts of the free surface in the grid somewhat away from the ship. To avoid related numerical problems, without dras- tically reducing the time step of the time integra- tion, we simply generate the grid in such a way that the waves are always inside regions foreseen for this purpose in front of and behind the ship. Grids with satisfactory quality have in the past been generated with our own elliptical grid generator. The grid for the KRISO container ship of the CFD Workshop Gothenburg 2000, [4], and for the Series 60 ship in this paper, Figs. 10 and 11, were constructed with this tool. However, due to the poor user interface, the grid generation takes several days using our program. Therefore, we now
have begun to use a commercial grid generator, re- ducing this time considerably, at acceptable costs regarding grid quality. The grid for the container ship C-Box in this paper was generated in this way. ;: ~ ~ 'in Figure 2: Partial view of the coarse grid at the cen- ter plane of the C-Box container ship The grid has roughly 1 million cells and 34 blocks. Because of the symmetry of the case con- sidered, only the starboard side of the ship was meshed. Most of the simulation results presented in this paper were obtained on a coarse version of this grid, in which every other grid line was skiped. Fig.2 shows a cut at the center plane of the ship. Because of the very intensive computational time, fine grid computations have been restricted to the model fixed condition till now. NUMERICAL METHOD The mathematical model described in the previous section was implemented by the authors in the RANSE code Neptun in the course of several research projects. The implementation details are certainly crucial for the accuracy of the resulting code. Most of these details are common to many other codes and/or have been described in previous publications e.g. A. Therefore, the overall numer- ical method is briefly outlined here and we focus only on a few particular aspects of our method. The conservation equations of momentum (7), mass (8) and turbulence parameters (9), (10) are discretised with a Finite Volume Method. Con- vective terms are approximated with the Linear Upwind Differencing Scheme (LUDS) and diffusive terms with the Central Differencing Scheme (CDS). The SIMPLE method is used to iteratively solve the resulting set of equations in each step of the time marching procedure used. The convection equation of the Level Set Function Eq. (11) is also discretised with LUDS in a finite volume manner and explicitely integrated in order to get the new position of the free surface. Then the derivatives of the ship velocities (u, v, w) and of the angular ve- locities Owl, w2, ~3) are eliminated in the equations of motion of the ship Eq. (3) and Eq. (4) and ex- plicitely integrated in order to get the position and velocities for the next time step. Hereby, known velocities and angular velocities from the previous time step (or initial condition) are used together with the forces and moments on the hull, deter- mined when solving the fluid equations. Note that a more sophisticated time integra- tion, e.g. as recommended by Soeding [10], was not necessary in the cases considered, because the aris- ing accelerations were relatively small. Because the integration of the ship motion equations does not take much computational time, significant time is not saved performing the integration with a larger time step than for the fluid equations. Thus, as long as the time step needed when solving the flow is smaller than the time step needed for stability reasons during the explicit integration of the ship motion equations, it is not worthwhile changing to an implicit integration, which demands more com- putational effort. Unbalance of Discretised Terms Because all variables are only stored for the cell centers, pressures and velocities at cell sides are interpolated linearly along grid lines to determine fluxes at the cell faces. If the mass term coso~i/F2 in Eq. (7) is treated in non-divergence form, while the pressure term is treated in divergence form, an unbalance between the discretised expressions of these terms can arise on curvilinear grids, leading to large numerical disturbances. L==~ Figure 3: Velocity field at a cross section of the KRISO container ship at rest. Usual algorithm (left) and improved algorithm (right)
In case of an uniform flow, with hydrostatic pressure distribution away from the ship, both terms should cancel. Fig.3 illustrates this effect at the cross section of the KRISO container ship at rest. The exact velocity distribution vi = 0 is strongly disturbed after some hundred iterations on the left hand side of the figure. A remedy of this problem consists in transforming the mass term in diver- gence form before integrating this term over the considered finite volume, determining the needed coordinates at the centers of the cell faces through linear interpolation like the pressures. The result- ing expression is then identical to the pressure term in case of a hydrostatic pressure distribution, avoid- ing the flow disturbancies almost completely, as shown on the right hand side of Fig.3. A similar effect could arise from an unbal- ance between the discretised forms of the convec- tive terms and the inertial terms in Eq. (7~. In a steady turning motion for instance, these terms should cancel each other, as well as the mass and pressure terms, in the undisturbed regions away from the ship. Despite no attempt having been performed to ensure the balance in this case, no disturbances have been detected in all cases anal- ysed untill now, because the angular velocities dur- ing ship motions are quite small. W~ ~ uO /~` // ~~\ N\ // ~~N N\ /~> << Figure 4: Comparison of the computed (black) and the undisturbed velocity field (red) in the water- plane of a twin rudder RoRo ship in steady turning The computed vector field at the waterplane of the twin rudder FSG RoRo ship in steady turn- ing (black) is compared with the exact velocity in the far field (red) in Fig.4. The non-dimensional Eq. (13) is explicitely integrated in the pseudo yaw rate, based on the ship speed Uo and ship time A. S(o) = I/ is a smoothed sign length L, is r' = 0.4. The agreement between both function. We choose ~ to be roughly as large as or. vector fields is very satisfactory in those regions When the steady state is achieved, the first term in where they are expected to agree. It is interesting to note that departures are restricted to the near neighbourhood of the ship and the wake. This can be seen more clearly in Fig.5, where the magnitude of the difference of both vector fields are depicted. o.' 1 QO9 ~ 0.08' 0.07 ~ 0.06' 0.05! 0.04 ~ 0.03' 0.02' 0.01 ~ o 1 Figure 5: Magnitude of the difference between the computed and the undisturbed velocity fields in the waterplane of the same RoRo ship in steady turning Reinitialization of ~ The convection equation (11) moves the free surface in the correct way but 0 does not remain a distance function. Keeping o a distance function is crucial to assure a constant thickness of the tran- sition zone between water and air, which is needed for stability and for limiting mass losses. Thus, the Level Set function ~ is 'reinitialized' before starting a new time step, i.e. replaced by a new distribution ¢, which in each point again represents the distance to the free surface. Because in all points on the in- terface the distance is zero, the interface remains unchanged while doing so, i.e. the isosurfaces ~ = 0 and o=0 coincide. The reinitialization is here performed deter- mining ~ as the steady state solution of the follow- ing equation, with the initial condition ¢(xt, 0) = ¢(xi, t): '9r (<b) (~- 1 ) = 0 (13)
Eq. (13) vanishes. Thus, it must be ~~ = 1. This is exactly the constraint which must be fulfilled by o to be a distance function. Usual convection schemes are not suitable for the discretisation of Eq. (13~. Here, the 1st order ENO scheme described in [11] for Cartesian grids was extended for 3D curvilinear grids. After integrating Eq. (13) in time, we get: +i = onfir S(O) (:D2 + D2 + D2 _ 1) (14) D denotes a discretisation operator for the magnitude of the x, y and z derivatives, see below. The non-dimensional pseudo time step b~ is set to 0.01 - 0.001 times or (smaller values for finer grids). Using the ENO scheme on a Cartesian grid, the x derivative writes: Do = max (max(D2, 0), min(D+, 0~) 0 > 0 Do = max (max(D+, 0),min(D2 ~ 0~) 55 < 0 D2 = (°i,j,k(iI,j,k) /^X D+ = (°i+l,j,k°i,j,k) /~X (15) The indices i, j, k denote the position of the grid point considered. DO and D+ are backward and forward finite differences. The role of the min/ max expressions is mainly to build the derivative with those points, which are closer to the interface o = 0. For y and z, analog expressions are valid. In order to extend the scheme to curvilin- ear grids, with coordinates (i, we denote the com- ponents of the contravariant base vectors, e.g. if], with As. and the Jacobian with J. The backward x-derivative can then be written: _ ~ (D-(2)Af +D-(~)A~ +D~ (~)A<) (16) The upper(x) in the RHS means, that the corresponding backward finite difference has to be built regarding to x. In order to determine how to proceed, we consider the components al of the first covariant base vector: D; ( ) = (¢i,j,kOiI,j,k)/~( for al > 0 Do ( ) = (¢i+l,j,k(i,j,k)/~( for al < 0 (17) All other derivatives are built in the same way. When solving Eq. (14), needed values of v at boundaries are set equal to the direct neighbour inside the domain. Figure 6: Isolines of the Level Set function at a cross section of the KRISO container ship at the beginning (left) and end (right) of the computation Fig.6 shows the result of the reinitialization at a cross section of the KRISO container ship men- tioned before, moving straight ahead in calm water. On the left hand side of the figure, the initial o dis- tribution is shown. The distribution at the end of the computation on the right hand side of the figure agrees satisfactorily with the distance to the free surface represented by the first isoline from above. Because o has to be reinitialized in the neighbour- hood of the free surface only, just 3 to 5 pseudo time steps solaces when solving Eq. (14~. BOUNDARY CONDITIONS The boundary conditions on the boundaries of the computational domain have remained almost unchanged compared to previous work, see e.g. [1], the only difference being that they now have to be expressed in a moving coordinate system. On the hull the no slip condition is enforced. Boundary conditions on the free suface are not necessary in conection with the Level Set method for two-phase flows. The level set function and the velocities are given at the inlet. At all other boundaries the <b values are set equal to the value in the center of the cell. On the lateral boundaries of the domain a mixed boundary condition is used for unsymmet- rical flow situations, e.g. steady turning motions, and the free slip condition otherwise. In the case of ship motions in waves we use the mixed condi- tion also on the top and bottom boundaries. This mixed boundary condition consists in setting the free stream values only for the convecting veloci- ties (mass flux) in the convection fluxes, where the fluid leaves the domain. The convected velocities
are still extrapolated by LUDS. Where the fluid enters the domain, all velocities are given as usual. On the outlet the pressure is given. For free surface computations it is set to the hydrostatic pressure relative to the undisturbed free surface: .o.o~s .o.o, -0.005 0.005 0.01 0.015 Poutiet = Fn2 (ZO + tat x + t32 y + t33 z) (18) 0.02 t3i, i = 1, 2, 3 are the components of the third row of the transformation matrix T. Due to ship motions and incoming waves, the boundary conditions in a point (x, y, z) at the inlet are in our coordinate system: (>irliet = Zo + tat x + t32 y + t33 z + (w (19) = = V3 inlet = (A~3y+W2Z) + (ill US + t21 VW + t31 WW) (v + w3 xw1 z) (20) + (tl2Uw+t22Vw +t32WW) (w w2x+wly) + (tl3Uw +t23vw +t33WW) (w is the wave elevation relative to the earth fixed coordinate system. Uw, Vw and Ww are the earth fixed velocity components in the generated incoming waves. The values at the inlet boundary consist of two parts, the contribution from the rigid motion of the ship fixed coordinate system and the absolut values in the far field. According to the potential flow theory we have in deep water: N (w = Mew sin(kjX+=jt+olj) (21) j=1 N Uw = ~Wj~we-kiz sin(kjX+wjt+o~j) j= VW = 0 (22) N Ww = ~ wj (w e kj Z cost; X + Wj t + A) j=1 N ist the number of superimposed harmonic waves. w; is the frequency, (w the amplitude, kj the wave number and off the phase of the wave component j. Each individual wave is assumed to run in negative earth fixed direction X. In order to completely express the inlet values in terms of t\,'>~: 1 0.75 0.5 0.25 ~ All- r -f TV,- - -I- - -it -A - - ~ - - - - ~ \ ~ / 1~ ~ g: ;~ ~ _ _ jet _ I_ ~ ~ V~ ~ '~ ~ 2 -0.25 -0.5 -0.75 - -1.25 -1.5 1.75 x/L Figure 7: Profile and velocities of a harmonic wave with amplitude (w/L = 0.01 and length A/L = 1. Computed profile (red) and potential theory (blue) o 0.05 0.1 n 1.~ 0.2 _ x/L Figure 8: Profile and streamlines of a harmonic wave with amplitude (w/L= 0.01 and length A/L= 1 ship fixed variables, the earth fixed coordinates X and Z of the boundary point considered have to be transformed according to Eq. 2. Using this technique, the waves are consid- ered to be fully developed already at the inlet. They can be seen as having been generated far in front of the grid and pass without disturbance the inlet into the computational domain. The main problem when generating waves with a finite vol- ume method is to avoid an excessive damping due to numerical dissipation. The implemented Level Set technique seems to do a good job in this regard. In Fig. 7 a snapshot of a computed harmonic wave (in red j superimposed to an uniform flow is compared with the corresponding profit and veloc- ity distribution of the potential theory (in blue). The given non-dimensional wave lenght and ampli- tude were (w/L = 0.01 and A/L = 1 in this case. The computational domain extends from x/L = 1 to x/L = -5, but a strong cell stretching was used starting from x/L = - 2 in order to damp out the waves toward the outlet. Only 100 cells were lo- cated in the shown region of the Cartesian grid used, in order to simulate a situation, which can be af- forded when performing ship motion simulations. Nevertheless, the computed wave elevation and ve- locities agree well with the theory. The streamlines
-0.01 5 -0.01 -0 ~5 n 0.005 0.01 0.015 , _ , I , , it, I , I , ~ tar - - - r - -~/~A,~ a- _ _ ~ _ _ _ ~ ~ ~ ~ ~ 5, \ ~ I'm ~t i-~ ~ r ~ ~ . >~ /~ ' \ ~ ~ \ ' ~12~ ~ ' ~ z~ IN ~ \ 0.02 -~,wk~l~l/,<~~,~]l&~-,~l/~',-,l~-<l<,/~<,~N 2 1 0.75 0.5 0.25 0.25 -0.5 -0.75 -1 -1.25 1.5 -1.75 x/L Figure 9: Profile and velocities of two superim- posed harmonic waves with amplitudes (W/L = 0.005 and lengths AI /L = 1 and A2/L = 2. Com- puted profile (red) and potential theory (blue) of the computed velocity field (after subtracting the superimposed uniform flow) are shown in Fig. 8. As expected no mass transport in x direction can be noticed. Superimposed harmonic waves with ampli- tude (W/L=0.005 and lenghts )~/L=1 and A2/L=2 are shown in Fig.9. The same grid resolution and computational domain than before were used. The agreement is very satisfactory also in this case. The time step and accuracy when solving the momentum and continuity equations are decisive in order to correctly predict the wave amplitude and lenght. The time step having more influence on the amplitude and the accuracy, on the wave lenght. Non-dimensional time steps of the order 0.001 and 3 to 5 SIMPLE iterations were sufficient for achiev- ing the very satisfactory results obtained. APPLICATIONS The method described in the previous sec- tions has first been used to calculate the flow around ship models moving steadily straight ahead against incoming waves. All degrees of freedom were sup- pressed in these cases. After having got some in- sight in the resulting flow and having compared the computed forces on the hull with measurements, we felt confident enough to go over to the simulation of the ship motions in head waves. Additionally, we have intensified our efforts in predicting forces and moments on the hull fore more realistic ma- noeuvring situations than in the past. We regard this as a requirement for reliable manoeuvring sim- ulations and we now focus on modern ship forms. Results for the Series 60 ship Simulations were done for the Series 60 ship CB = 0.6 at Froude number Fir = 0.316 moving steadily straight ahead against incoming harmonic waves of amplitude (w/L = 0.01 and the wave length A/L = 1, yielding a non-dimensional en- countering period of 0.442 time units. The compu- tations were performed at model-fixed condition at Reynolds number Rn = 5.27 106. A relatively small computational domain, ex- tending from X/L = 1.5 at the inlet to X/L = - 1.5 at the outlet was discretised with a grid having roughly 400000 cells. The time step was set to 0.001 and 3 SIMPLE iterations were done per time step. After having achieved the steady state in calm water, the wave generation was started. The computation was continued until the wave passed the outlet of the domain, demanding roughly 40 hours on a small workstation. Calm water Waws I I Figure 10: Computed wave patterns in calm water and in head waves. Series 60 ship at Fir = 0.316 In Fig. 10 the wave pattern in calm water and in waves are compared. Due to the high Froude number, the pronounced wave system of the ship is moderately distorted by the incoming waves. Two snapshots from the simulation are shown in Fig. 11. The first snapshot shows the instant when the wave has almost reached the (not depicted) ship bow after 0.66 time units, according to the result- ing relative phase velocity. In the second snapshot the wave front has already passed the outlet. Results for the C-Box container ship Computations in model-fixed condition were also performed for the C-Box container ship of the
Figure 11: Two snapshots of the simulation of the flow around the Series 60 ship in head waves german ship yard FSG, in order to exploit available force measurements for validation purpose. Simu- lations were done at Fr' = 0.23 and Rr' = 8.2106, A again with (WIL = 0.01 and A/L = 1. In Fig. 12 the wave patterns computed on the coarse grid and on the fine grid are compared for a given instant. Due to the relatively low Troupe number, the wave system of the ship can be seen less clearly now. Nevertheless, it is interesting that on the coarse grid having roughly 125 000 cells, in- stead of 1 million as for the fine grid, almost all wave features are captured. In Fig. 13 the cor- responding wall shear stresses obtained on both grids are compared. Now, the poorer quality of the coarse grid results is more obvious. Fig. 14 compares the computed and mea- sured longitudinal force acting on the hull of the model. Because of the very high computational time, simulations on the fine grid could not be suf- ficiently continued to compare the obtained time history with experiments. Therefore, the predicted force on the coarse grid is shown instead. The agreement is surprinsingly good, taking into ac- count the grid resolution. The encountering pe- riod has been captured almost exactly and even the amplitude is well predicted. The computed mean value of the longitudinal force departes roughly 10% from the measured value, but confirmes the very low added resistance compared to the calm water value when all degrees of freedom are suppressed. A snapshot from the simulation on the fine grid at a moment when a wave crest is at the bow is compared with a photo from experiment in Fig.17. Due to non- hydrodynamical considerations, this fine gnd Hi_ coarse grid Figure 12: Wave patterns in head waves computed in the coarse grid (bottom) and in the fine grid (top). C-Box container ship at En = 0.23 _ ~ ~ ; , ~ ~ I, , ' ~ ~ ' ~ ~ ~ ', ',' ~ ~ ~ ~ ~ ~ ~ ~ ~ , 'I , = 'a ; ~ _<,3,~ ~~ . it. . .~ ... 5 ~ 5 ~ ~-3__= ~ + ~,;~ ~-.3-~ A= -= ~ I =, ;, ..~ ~ I.= ~ ~~ ::: .:::,~ .: _ ~=.~ :-~ >3 += zest (~ ~~:-~,-:~ .~P~:~ Figure 13: Wall shear stress field at the hull of the C-Box container ship in head waves computed on the coarse and on the fine grid container ship has got a relatively large block co- eflicient CB = 0.74 and was not provided with a bulbous bow. As a result, pronounced wave break- ing occur at the bow. This is captured, at least qualitatively, by the simulation. Unfortunately, the side walls of the ship were simply vertically pro- longed when generating the grid, regardless of the real form of the back. This seems to have an im- portant influence on the wave breaking process. The above computations were then repeated for a condition where the ship model is free to pitch, heave and surge. During the simulation the pro- peller thrust was kept constant at the value corre- sponding to the calm water situation. Again the simulation started from the calm water solution without taking dynamic sinkage and trim into ac- count. Therefore, when the wave generation be- gins, the initial conditions does not represent an equilibrium situation. This explains why the ship carries out small heave and pitch motions already
.i,IN] -1 60 -120 -80 .dn 40 80 on Figure 14: Comparison of the computed (solid red line) and the measured (dashed blue line) longitudi- nal force at the C-Box container ship in model-fixed condition in head waves before the waves reach the hull. Moreover, the gen- erated waves are not harmonic at the beginning, neither in the towing tank nor in the simulation. For these reasons, a comparison of computed and measured complete time histories is not possible. A selected interval of the computed pitch angle time history is compared with experiment in Fig. 16. Again, the period is well predicted, but the am- plitude is strongly underestimated. The reason for this discrepancy most probably relies on a stronger influence of the coarse grid in this case then for the model-fixed condition, leading to an increased wave damping. In Fig. 15 the computed axial velocity in the propeller plane is compared with the results in calm water. The axial velocity field shows much higher values now. The chosen instant corresponds to the situation where a wave trough is at the stern, lead- ing to a strong acceleration of the flow there. Finally, we show some snapshots of the sim- ulation performed for the C-Box container ship in harmonic head waves. The sequence depicted in Fig. 18 covers a periode of the surge, heave and pitch motions, starting on the top of the figure from a slightly positiv trimed position. The wave break- ing mentioned above and observed in the model tests in both, the model-fixed and the model-free condition can certainly not be predicted very real- istically on the coarse grid used, but are at least present. ~ o.s 3 0.8 0-7 0.6 0.5 0.4 0-3 0.2 0.1 1 o - :~ Figure 15: Computed axial velocity in the propeller plane of the C-Box container ship in calm water (top) and in head waves (bottom) Results for the FSG RoRo ship Computations for many combinations of the non-dimensional yaw rate r, = rL/Uo and the drift angle ~ at the main section were performed for a modern RoRo ship of the FSG. The Froude number En = 0.20 and the Reynolds number Rn = 1.53 107 were determined with a suitable speed for steady turnings. The predictions of forces and moments on the hull of this twin screw, twin rudder ship model are compared with experimental data ob- tained with the Computarized Planar Motion Car- riage (CPMC) in the towing tank of the Hamburg Ship Model Basin. The ship motions were prescribed and the free surface was not taken into account. Because of the low Fioude number, we do not expect this to have a significant influence on the results. The computations were performed on a non-matching block-structured grid with roughly 1 million cells. Between 1000 and 2000 SIMPLE iterations were necessary to achieve a good convergence, demand- ing a CPU time of 7 to 14 hours on a PC.
16 1R Figure 16: Comparison of computed (solid red line) and Reassured pitch motion (dashed blue line) in head waves. C-Box container ship at Fit = 0.23 __ ~ = 1·1 1 _ - _ - ~ 51 1 _ __ ~ , ~ ~ _ - _e ~ 1 - - i - Figure 17: C-Box container ship in model fixed con- dition in head waves. Comparison of a snapshot of the simulation on the fine grid with a photo of the model test The flow features at the fore body of the RoRo ship turning with constant yaw rate r'=0.4 and drift angle p=10°, which represents a realis- tic situation, look very similar to those of a ship in straight motion. Because the bow shows towards the center of the circle, the local drift angle at the bow is practically zero. This yields a very sym- metric pressure distribution at the bulbous bow as shown in Fig.19. Contrary, on the rear part of the ship both the contributions of r' and ,B to the lo- cal drift angle yield a large cross flow which can be discern on the pressure distribution on the skeg. The computed non-dimensional side force Y' Ye/ Uo2 L2) is compared with measurements in Fig.20. Although the overall trends are captured well, discrepancies are somewhat larger than in ear- - - Figure 18: Sequence of the simulation of the flow around the C-Box container ship in head waves, coarse grid tier cases. The latter especially applies to the com- parison of the non-dimensional yaw moment N' = N/~2 Uo2 L3) with experiments shown in Fig. 21. These larger discrepancies are probably dun to the propellers, shafts and V-brackets present in the ship model, Fig. 22, which were not taken into account in our computations. To check this and to improve our predictions we now are going to in- clude shafts and brackets in the computational grid and repeat the computations. Next, computations will include rudder deflections and heel angles. ACKNOWLEDGEMENT This work was supported by the German Ministry of Education and Research.
4 2 4 , -6 Figure 19: Computed pressure coefficient on the hull of the FSG RoRo ship at steady turning with r' = 0.4 and ~ = 10° - -10 . , ,_____,_ ...... , i, ,~ : _1 ___ _' - ~ - ~ .......... ... ~r'=0.4~ ~ .................................................. I I ...................... , ~ '^? ., . ~, /~ 1 2. .~ 1 ' ............................................ 1 -. ~ : ~ - 1 3 1.5 1 5 o -A t Figure 21: Comparison of measured and computed non-dimensional yaw moment on the FSG RoRo ship for different yaw/drift combinations Figure 20: Comparison of measured and computed non-dimensional side force on the FSG RoRo ship for different yaw/drift combinations Figure 22: Stern of the model of the FSG RoRo ship with rudders and propellers, shafts and V- brackets
References [1] Cura Hochbaum A. (1994), " A Finite-Volume Method for Turbulent Ship Flows", Ship Technology Research Schiffstechnik, Ham- burg, Germany t2] Cura Hochbaum A. (1998), " Computation of the Turbulent Flow Around a Ship Model in Steady Turn and in Steady Oblique Motion", 22nd ONR Symp. on Naval Hydrodyn., Wash- ington D.C., USA t3] Cura Hochbaum A. and Schumann C. (1999), " Free Surface Viscous Flow Around Ship Mod- els", 7th Int. Conf. on Numerical Ship Hydro- dyn., Nantes, France t4] Cura Hochbaum A. and Vogt M. (2000), " Flow and Resistance Prediction for a Container Ship", A Workshop on Numerical Ship Hydro- dynamics, Chalmers University of Technology, Gothenburg, Sweden [5] Kinoshita T., Kagemoto H. and Fujino M. (1999), "A CFD Application to Wave Induced Floating Body Dynamics", 7th Int. Conf. on Numerical Ship Hydrodyn., Nantes, France t6] McDonald H. and \Yhitfield D. (1996), "Self- Propelled Maneuvering Underwater Vehicles", 21st ONR Symp. on Naval Hydrodyn., liond- heim, Norway t7] Miyake R., Kinoshita T. and Kagemoto H. (2000), " Ship Motions and Loads in Large Waves", 23rd ONR Symp. on Naval Hydro- dyn., Val de Reuil, France [8] Osher S. and Sethian J.A. (1988), " Fronts Propagating with Curvature-Dependent Speed: Algorithms Based on Hamilton-Jacobi Formulations", Journal of Computational Physics Vol.79 [9] Sato Y., Miyata H. and Sato T. (1999), "CFD Simulation of 3-Dimensional Motion of a Ship in Waves: Application to an Advancing Ship in Regular Heading Waves", Journal of Marine Science and Technology Vol.4, pp. 108-116 [10] Soeding H. (2001), " How to Integrate Free MA tions of Solids in Fluids", 4th Numerical Tow- ing Tank Symposium, Hamburg, Germany [11] Sussman M., Smereka P. and Osher S. (1994), "A Level Set Approach for Computing So- lutions to Incompressible Two-Phase Flow", Journal of Computational Physics Vol.114 t12] Takada N., Miyata H. and Sato T. (1999), "CFD Simulation of 3-Dimensional Motion of a Vehicle with Movable Wings", 7th Int. Conf. on Numerical Ship Hydrodyn., Nantes, France t13] Vogt M. (1998), "A Numerical Investigation of the Level Set Method for Computing Free Sur- face Waves", Licensiate Thesis at Chalmers University of Technology, Gothenburg, Swe- den [14] Vogt M. and Cura Hochbaum A. (2002), "Nu- merical Simulation of Ship Motions in Head Waves with a RANSE Method", HSVA Re- port 1649 (to appear), Hamburg, Germany t15] Wilcox D.C. (1993), "Turbulence Modeling for CFD", DCW Industries, La Canada, Califor- nia, USA t16] Wilson R., Paterson E. and Stern F. (1998), " Unsteady RANS CFD Method for Naval Combatants in \Vaves", 22nd ONR Symp. on Naval Hydrodyn., \Vashington D.C., USA
