A Simple and Efficient Method for the Coupled Motions of Multiple Bodies

HUANG Ming-han,ZOU Zhi-li

(State Key Laboratory of Coastal and Offshore Engineering,Dalian University of Technology,Dalian 116024,China)

1 Introduction

The analysis of hydrodynamic problems is an important subject of ship security and design.The hydrodynamics of a single structure in waves,such as a ship or a production platform,is usually calculated by the three-dimensional source-sink distribution method[1-2].But difficulties in actual computations may arise for calculations of the interaction of multiple bodies[3-6],since a lot of source panels are needed for this kind of numerical calculations,especially for a high wave frequency case.The computation work increases with the increase of size and number of structures.So for multiple-bodies case,a simple and efficient method is called for.The method proposed by Lean et al[7]is a possible choice for such a purpose.This method approximates a ship cross section by a rectangle with the same area of the original ship section,and divides the flow field into inner domain below the ship bottom and the outer domain beside ship hull.The inner-domain velocity potential is given by a simple analytic solution,and the outer-domain velocity potential,given by a source-sink distribution along the contour of ship water area,is expanded over water depth as a Fourier series.

The present study extends Lean’s method to solve the interaction of multiple ships and piers.The calculation for a vertical pier is treated as the case of the ship with a zero gap between ship bottom and seabed.By this method,the present method couples the calculation for multiple ships and the calculation for multiple piers and gives a simple and efficient method for the interaction of multiple ships and piers in waves.

2 Mathematical models

In order to describe the mathematic model for wave motions induced by multiple bodies,define a global coordinate system oxyz and a set of local coordinate systems oMxMyMzM（M=1,2,…,NS+N C）,in which NS and NC denote the numbers of ships and piers,respectively.The former is fixed on earth and the latter is fixed on each body,as shown in Fig.1.The global coordinate system has the z-axis pointing vertically upward and the oxy plane rests on still water surface.The local coordinate system has its origin locating at the ship’s center of gravity or on still water surface for the pier,xM-axis parallel to the body’s longitudinal direction(to the ship’s bow)and zM-axis pointing vertically upward,β the incident wave angle(β=0°for head sea).

Fig.1 Description of coordinates

Fig.2 Ship section and equivalent rectangle in calculations

Suppose the fluid is ideal,incompressible,inviscid and irrotational.Consider a ship of the length L,beam B and draught D in constant water depth h.The incident wave amplitude and the ship motion are assumed to be small to meet the linearity assumption.The flow field can be described by the velocity potential

in which‘Re’means taking the real part.The complex,time-independent velocity potential φ in the fluid field is written for the inner and outer domains in the form

The first expression in Eq.(2)is the velocity potential in the outer domain.φ0and φ70are the incident wave and diffraction wave potentials,respectively.φjNis the radiation potential of the jth mode of motion of the Nth ship（N=1,2,…,NS）.The second expression in Eq.(2)is the velocity potential in the inner domain for a ship.For piers,there is not the inner potential.φujNis the radiation potential of the jth mode of motion of the Nth ship（N=1,2,…,NS）,φu70is the inner-domain potential induced by diffraction waves.（x′,y′,z′）are the inner domain coordinates.ζjNis the amplitude of the jth mode of motion of the Nth ship.j=1,2,…,6 denotes surge,sway,heave,roll,pitch and yaw,respectively.The incident wave potential is given by

where ζ0is the incident wave amplitude,ω the incident wave frequency,g the acceleration of gravity,and k0wave number which satisfies the usual dispersion relation:ω2=gk0tanhk0h.

Based on linearized potential theory in the frequency domain,the governing equation and linear boundary conditions to be satisfied by the diffraction and radiation potentials of inner and outer domains,φjNand φujN(for convenience of writing we take diffraction potential and radiation potential coming together to increase a group number j=7 and N=0),are summarized as follows:here υ=ω2/g,r2=x2+y2,R2=x2+y2+z2,SMis the surface of the Mth body,is the jth component of the normal vector:where the normal vector⇀nis defined as positive when directing into the fluid from the boundary surface.The condition[M]denotes matching condition on the interface S′belonging to both the inner domain and the outer domain.

3 Velocity potential in the outer domain

For the outer domain,the velocity potential can be obtained by distributing source-sink only along the contour of ship water area.For the calculation of ships,we approximate the ship cross section with a rectangle with the same area in Fig.2,and the flow field is divided by the inner domain below ship bottom and the outer domain outside ship hull.Here we assume that the gap δ between ship bottom and seabed is small quantity compared to water depth,which is the case in a harbor or coastal region.So the velocity potential in the inner domain can be expressed in a simple analytic form.For the outer domain,the velocity potential is expressed by a source-sink distribution method along the contour of ship water area.Due to the rectangle form of the simplified ship cross section,the lateral sides of the hull are vertical.Therefore the Green function of the flow can be expressed by a Fourier series expansion over water depth.The inner and outer domains potentials are matched at the interface of two domains beneath the ship hull.To satisfy the governing equation and boundary conditions(4)-(9),the resulting potential has the following form for ships:

where p=（xp,yp）and q=（xq,yq）are the field point and source point respectively,andare the source intensities on the water surface of Mth body due to the jth mode of motion of the Nth ship.There can be determined by the boundary conditions(6)and(10);anddenote the first kind of Hankel function and the second kind of modified Bessel function;k0and kiare the roots of the dispersion equations k0tanhk0h=ω2/g and kitankih=-ω2/g.

Solving the case of piers is equivalent to the case of ships with a zero gap between ship bottom and seabed,and there is no inner domain.The flow round the sides of the pier can be represented by distributing sources along the contour of pier water area.The Green function of the source-sink distribution can be expressed by a Fourier series expansion along water depth.The resulting potential has the following form for piers:

where p=（xp,yp）and q=（xq,yq）are the field point and source point respectively,is the source intensity on the water surface of Mth body due to the jth mode of motion of the Nth ship.There can be determined by the boundary conditions(6)and(11);denotes the first kind of Hankel function;k0is the roots of the dispersion equations k0tanhk0=ω2/g.

4 Velocity potential in the inner domain for ships

For the inner domain,some simplifications are introduced to the description of the hull shape to make solving the water flow easier:we assume the hull is everywhere rectangular in cross-section with the sides vertical and the bottom flat.We can assume that the water flow is essentially horizontal and lateral,but lateral flow velocities are much greater than vertical and longitudinal velocities.So the three-dimensional flow under ship bottom is simplified by twodimensional cross-flow along the ship width(see Fig.3),neglecting longitudinal-flow along the ship length.And we can have a simple method to match potentials between the inner domain and the outer domain on the lateral sides of the hull.Now we give the concrete method.

Fig.3 The flow along ship width direction in inner domain

Fig.4 Grid distribution along the contour of body water surface

Integrating the two-dimensional continuity equation from water bottom to ship bottom,this leads to the following equation:

substituting

into Eq.(13),and integrating gives

where f1（x′,z′）and f2（y′,z′）can be solved by substituting Eq.(15)into the following vertical two-dimensional Laplace equation

This process also relies on the expression of,which is a constant(-iω,for heave motion;0,for surge,sway,yaw and diffraction motions)or linear function for roll and pitch motions.According to above results we can get the expression φujN:

where subscript u means under the ship bottom;p′is the field point;FujN（p′,z′）is velocity potential concerned about two-dimensional flow due to the jth mode of motion of the Nth ship;cross-flow under-keel velocity VujN（p′）and mean under-keel potential φujN（p′）can be solved by the matching condition.

5 Solving unknown source intensities and velocity potentials

Then we get the integral equations for the ship below:

Integrating matching condition(8)along the gap between seabed and ship bottom,we obtain matching equation

For the case of the pier,we take potential(11)into body surface condition(6),multiply vertical function coshk0（z+h）and have integral equations for the piers below:

where the subscripts 0 and i of(20a,b)correspond to the subscripts 0 and i of(10),N2=N1+1,substituting N1(N1equals to 5 in this paper)for∞in(10),the others are as follows:

For j=1,2,…,6;N=1,2,…,NS:

For j=7;N=0:

In order to solve above integral equations,the contour of Mth body on free water surface is discretized into 2n（M）grids.The distance between two adjacent grids is dx（M）.n（M）is the function of number for the Mth body.Therefore Eqs.(20a,b),(21)and(22)can be expressed as

here K1 is the unit serial number of field point peK1=（xeK1,yeK1）which denotes p′eK1=（x′eK1,y′eK1）in the local coordinate system;K2 is the unit serial number of source point qeK2=（xeK2,yeK2）which denotes q′eK2=（x′eK2,y′eK2）in the local coordinate system;other parameters can be written as

when the field point and the source point belong to the common section,A0and Aican be expressed in the form:

6 Linear pressure forces and ship motions

Once the velocity potentials on the body surface are obtained,it is straightforward to compute the first-order hydrodynamic forces.The linearized unsteady pressure on the body surface can be computed directly from the obtained velocity potentials.Integrating this pressure multiplied by the ith direction of the Mth ship can be computed.Likewise the restoring forces can be computed from the change of hydrostatic pressure due to the displacement of a ship from its equilibrium position.These results may be expressed in the form:

in which the horizontal（i=1,2,6）direction

the vertical（i=3）direction

the horizontal and vertical（i=4,5）directions

here ρ is density of water,FWiMis the wave exciting force in the ith direction of the Mth ship,and μijNMand λijNMare the added mass and damping coefficients,respectively,in the ith direction of the Mth ship due to the jth mode of motion of the Nth ship.CijMin Eq.(26)denotes the restoring-force coefficients,which are free from hydrodynamic interaction among NS ships and thus almost zero except for CjjM（j=3,4,5）and C35M=-C53M.

Denoting the generalized mass matrix of the Mth ship with mijMand using Eq.(26)for the hydrodynamic and hydrostatic forces on the Mth ship,the coupled motion equations of NS ships are written in the form

The complex motion amplitude ηjMcan be determined by solving these coupled equations;thereby the first-order solution will be completed.

Fig.5 Small scale body plan

7 Validation of the present method

In order to validate the numerical results,the results for a ship are compared with experimental result of Oortmerssen[1].The experimental result includes two parts:one giving the responses of a single ship,and the other giving the hydrodynamic coefficients for a ship against a quay.For these two parts,the ship hull is the same,which is a 200 000 TDW tanker.Particulars of the ship are listed in Tab.1.Fig.5 shows the cross sections of the ship at 21 equally spaced stations along its length starting at the stern:stern sections appear on the left and bow sections appear on the right.

Added mass μiiand damping coefficients λii(ith mode of motion),wave exciting forces Fiand responses ζi(ith direction),period T and frequency ω are presented in non-dimensional form as shown below:

here V is the volume of ship displacement,S is the area of free surface water of the ship.The comparison of the present numerical results with the experimental data is given in the following two parts.

Tab.1 Tanker details

7.1 Comparison for a single ship

Fig.6 shows ship responses for a single ship case obtained from numerical calculation and experiments.The results are given for three wave directions 180°,225°and 270°and sixdegrees of freedom of the ship motions.The results show that numerical results are in better agreement with experimental results.

Fig.6 Ship responses

7.2 Comparison for a ship against the quay

In the experiment,the quay wall is parallel to the ship for the case of a ship against the quay.The distance between the ship hull and the quay wall is b=16.5 m,as shown in Fig.7.We calculate this case by two different ways in order to validate the present method.One is to apply the image method and the other is to consider two ships symmetric about the quay wall and remove the quay wall,which is shown in Fig.8.For the image method,the ship is reflected with respect to the quay wall,and this can be treated in numerical calculation by changing the Green function,i.e.doubling the source strength of the Green’s function for the case of a single ship without the quay.So for this method we need only to deal with one ship.But for the second method,we apply the method mentioned in previous sections for the calculation of multiple ship case.In doing this,we need to deal with the two ships separately,with each ship being represented by source-sink distribution.So the Green function does not need to be changed for this method.Therefore the element number of the second method is doubled that of the first method.The comparison of results of two methods will give the validation of the present method,if the two different methods give the same results.

Fig.7 A ship against the quay

Fig.8 Two ships symmetric about the quay wall

Figs.9 and 10 show added mass and damping coefficients in sway and heave for the wave direction perpendicular to the ship.It can seen from the figures that there is a peak in the curves,and this is because the water trapped in the clearances between the ship and the quay can go into resonance[8].The agreements among the results given by the two-ship calculation,the image method and the experimental data are good.This gives a validation of the present method for the coupled motions of multiple ships.

Fig.9 Added mass in sway and heave

Fig.10 Damping coefficients in sway and heave

8 A calculation for the interaction between ships and a harbor

In the previous section,we have discussed the case of a ship against the quay,but the actual harbor is usually semi-closed.So here we study a more complex case:a ship in the U-shaped harbor,as shown in Fig.11.The U-shaped harbor consists of three quays:the middle quay of 1 040 m long and the two lateral quays of 700 m long.All the quays are 20 m in thickness.The ship considered here is the same as that mentioned in section 7.The distance between the quay wall and ship hull is also the same as before(16.5 m).Since the quays have large size,the calculation of interaction between the ship and the quays will involve large computation works if the conventional three-dimensional source-sink distribution method is used.As shown below,the computation can be largely simplified when the present method is applied for this case.It is known from the previous sections that the simplification introduced by the present method comes from the simple velocity potentials given by Eqs.(10)and(11)for the calculation of wave motions induced by ships and piers,which will be used here to calculate the wave motions induced by the ships and quays.

Fig.11 A ship in the U-shaped harbor

Fig.12 Added mass and damping coefficients in sway and heave

From the following results,it can be seen that a ship in the U-shaped harbor is influenced on the standing waves between the two lateral quays,so the responses of a ship in this kind of semi-closed harbor are different from the case of a ship against the quay shown in Fig.7.Figs.12-14 show hydrodynamic coefficients,wave forces and ship responses for wave direction perpendicular to the ship.In order to make comparison,the results for the case of a single ship and a ship against the quay are also given in the figures.It can be seen in the figures that there are many peaks in the results for a ship in the U-shaped harbor,but this is not the case for the cases of a single ship and a ship against the quay.The peaks are caused by the occurrence of standing waves excited by the ship between the two lateral quays.The surface elevation of the standing wave between two quay walls has the following form:

Fig.13 Wave forces in sway and heave

Fig.14 Ship responses

here Anis the wave amplitude,wave number kn=nπ/d which has the relation with wave frequency ωn:ω2n=gkntanhknh.d is the length between the two lateral quays.Since the presence of the ship will influence the standing waves,standing waves across the ship-occupied region will deviate from that given by Eq.(39).In order to show the relation between the peaks in the figures and the standing waves,we calculate the corresponding wave number knand n which correspond to the frequency of the peaks in the figures.The values of n calculated by this way have been shown in Figs.12-14.From these values of n,we can judge the form of standing waves corresponding to the peaks in the figures:If the value of n is an integer or close to an integer,the peak in the figures will correspond to the standing wave given by Eq.(39).For example,the values 3.1,4.9,5.8 of n in Fig.12 are not integers,but close to integer n=3,5,6,so the peaks corresponding to these values of n are produced by the standing waves given by Eq.(39)with n=3,5,6.Otherwise,the value of n is not an integer or far from an integer,this means that the peak in the figures is produced by the standing wave that deviates from the standing waves given by Eq.(39),that is,it is produced by the standing wave appears over the ship-occupied region.

From above analytic method by using the values of n shown in the figures,we can know that the multiple peaks in the figures are reflection of the effect of the standing waves due to the presence of the two lateral quays.This shows that there is quite strong interaction between the ship and the U-shaped harbor.This interaction is presented in the hydrodynamic coefficients,wave forces and the ship motion responses.

9 Conclusions

This paper proposes a simple and efficient method for the interaction of multiple bodies in waves,such as the interaction among multiple ships,piers,and quays.This method combines the method of ships and the method of piers to solve the coupling problems of multiple bodies.The simplification introduced by the present method comes from the simple velocity potential given by Eq.(11)for piers and the simple velocity potential given by Eq.(10)for ships.The form is used in the present method to calculate the wave motions induced by the breakwater or quays and can reduce the computation work greatly.The later can also greatly simplify the calculation for multiple ships,since it only needs to distribute the sources along the contour of ships water area.The comparisons of experimental data show that the present method has good accuracy.By applying the method to calculate the interaction between a ship and the U-shaped harbor,it is found that there are many peaks in the curves of results.This reflects that there is strong interaction between a ship and a semi-closed harbor.The simplicity and the efficiency shown by this application of the present method demonstrate that the present method is a very useful method for the calculation involving large size harbor and multiple bodies(ships and piers).

[1]Gerard van Oortmerssen.The motions of a ship in shallow water[J].Ocean Engineering,1976,3(4):221-255.

[2]Gao Zhiliang,Zou Zaojian,Wang Huaming.Numerical solution of the ship wave-making problem in shallow water by using a high order panel method[J].Journal of Ship Mechanics,2008,12(6):858-869.

[3]Chen Gungrong,Fang Mingchung.Hydrodynamic interactions between two ships advancing in waves[J].Ocean Engineering.2001,28:1053-1078.

[4]Choi Y R,Hong S Y.An analysis of hydrodynamic interaction of floating multi-body using higher-order boundary element method[C]//Proceedings of the 12th International Offshore and Polar Engineering Conference.ISOPE,Kitakyushu,Japan,2002,3:303-308.

[5]Masashi Kashiwagi,Kazuaki Endo,Hiroshi Yamaguchi.Wave drift forces and moments on two ships arranged side by side in waves[J].Ocean Engineering,2005,32:529-555.

[6]Sun Zhaochen,Wang Huanyun.Coupling calculation of hinged multi-body floating structure by analytical method[J].Journal of Ship Mechanics,2009,13(1):34-40.

[7]Lean G H,Bowers E C,Spencer J M A.Underkeel:A computer model for vertical ship movement in waves when the underkeel clearance is small[R].HR Report SR 160,1988.

[8]Huang Minghan,Zou Zhili.A simple and efficient method for multiple ships and piers in waves[J].Journal of Ship Mechanics,2010,14(11):1227-1240.