Marine Structures 33 (2013) 71–99
Contents lists available at SciVerse ScienceDirect
Marine Structures journal homepage: www.elsevier.com/locate/ marstruc
Fully coupled BEM-FEM analysis for ship hydroelasticity in waves Kyong-Hwan Kim a, e, Je-Sung Bang b, Jung-Hyun Kim a, Yonghwan Kim a, *, Seung-Jo Kim c, Yooil Kim d a
Department of Naval Architecture & Ocean Engineering, Seoul National Univ., 1 Gwanak-ro, Gwanak-gu, Seoul 151-744, Republic of Korea Systems Engineering Research Division, Korea Institute of Machinery and Materials, 156 Gajeongbuk-Ro, Yuseong-Gu, Daejeon 305-343, Republic of Korea c School of Mechanical and Aerospace Engineering, Seoul National Univ., 1 Gwanak-ro, Gwanak-gu, Seoul 151-744, Republic of Korea d Department of Naval Architecture & Ocean Engineering, Inha University, 253 Yonghyun-dong, Nam-gu, Incheon 402-751, Republic of Korea e Korea Institute of Ocean Science & Technology, 1312-32, Yuseong-daero, Yuseong-gu, Daejeon 305-343, Republic of Korea b
a r t i c l e i n f o
a b s t r a c t
Article history: Received 16 July 2012 Received in revised form 17 March 2013 Accepted 13 April 2013
This paper considers the problem of ship hydroelasticity, which is an important technical issue in the design of ultra-large vessels. For the analysis of fluid-structure interaction problems, a partitioned method is applied. The fluid domain surrounding a flexible body is solved using a B-spline Rankine panel method, and the structural domain is handled with a three-dimensional finite element method. The two distinct methods are fully coupled in the time domain by using an implicit iterative scheme. The numerical results of natural frequency and the motion responses of simple and segmented barges are computed to validate the present method through comparisons with experimental and numerical results. This study extends to the application to two real ships, 6500 TEU and 10,000 TEU containerships, for more validation and
Keywords: Ship hydroelasticity Fully coupled analysis Rankine panel method Finite element method Direct time integration
* Corresponding author. Tel.: þ82 2 880 1543; fax: þ82 2 876 9226. E-mail addresses:
[email protected],
[email protected] (K.-H. Kim),
[email protected] (J.-S. Bang),
[email protected] (J.-H. Kim),
[email protected] (Y. Kim),
[email protected] (S.-J. Kim),
[email protected] (Y. Kim). 0951-8339/$ – see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.marstruc.2013.04.004
72
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
also observation on the practicality of the present method. Based on this study, it is found that the present method provides reliable solutions to linear ship hydroelasticity problems. Ó 2013 Elsevier Ltd. All rights reserved.
1. Introduction Recently, commercial ships have been becoming larger and faster due to the increase of global trade. As a result, the natural frequencies of ship structures become lower, approaching the frequency range of waves in open sea. As ships become faster, the encounter frequency of the ships tends to increase, moving closer to the frequency range of hull girder vibration. These trends contribute to the increased occurrence of ship springing phenomena, and the risk of fatigue failure is also increased. Ship hydroelasticity should be taken into account in the design of large vessels such as ultra-large containerships. In previous studies, the hydroelastic behaviors of ships have been successfully analyzed based on the modal superposition method, which uses the Eigen modes of an elastic body to represent the elastic deformation of ships [1–3]. Ship structure is idealized to several beam elements, and the hydroelastic response of the ship is analyzed based on a modal approach. Strip theory in the frequency domain has been used to handle the fluid domain surrounding a flexible body in most previous studies. The coupling approach of both the one-dimensional beam theory for structural analysis and the strip theory for hydrodynamic analysis is very useful in the initial stage of design due to its simplicity and low computational burden. Recently, a time domain method is applied [4,5], and a nonlinear analysis for hydrodynamic problem is carried out [6]. As computational capacity has increased, three-dimensional FEM (Finite Element Method) as a tool for structural analysis has been acclaimed rather than one-dimensional beam analysis [7,8]. It is regarded that one-dimensional beam analysis is vulnerable in the prediction of localized and transversal stress. Three-dimensional FEM needs much more computational time compared with beambased analysis, but localized and transversal stress on hull girders can be accurately evaluated. Moreover, it has no limitation in terms of the complexity of hull geometry, such as with containerships that have large open mid-ship sections. It is not easy to analyze complex hull geometry using onedimensional analysis. In particular, a proper modeling of offshore structures is very hard using beam elements. Accordingly, it is necessary to apply three-dimensional FEM to the hydroelasticity problem. Recently, three-dimensional Rankine panel method has widely been applied to ship motion problems, as it is able to reflect three-dimensional effects which cannot be included in strip theory. It is easy to import nonlinear effects induced by free surface and body shape in this method, which are hard to consider in the application of wave Green function method. Furthermore, this method needs much less computational time than methods based on CFD (Computational Fluid Dynamics). Despite such merits, it is hard to find research that adopts Rankine panel method in ship hydroelasticity problems, especially coupling with three-dimensional FEM. Malenica et al. [7] applied a wave Green function method by coupling with three-dimensional FEM. In their study, only the first few lowest modes computed by three-dimensional FEM were considered to represent the elastic behavior of ship structure. It is a clever approach to tackle the low modes of hull girder vibration induced by hydroelasticity. However, it is hard to obtain the subset of modes necessary for accurate representation of localized stress on ships in their approach. In spite of a huge effort based on CFD to analyze the ship motion problem, it is regarded that CFD-based analysis is not yet mature enough for the prediction of ship motion. Rankine panel methods can be an alternative for accurate ship motion prediction and could reduce computational burden in its current status. As mentioned, a coupling approach using three-dimensional FEM and Rankine panel method with a direct integration scheme can be an applicable technique for the analysis of global elastic behavior and localized stress of ships, because modal superposition has limitations in the prediction of stress with high accuracy since it uses a finite number of mode shapes. This limitation of modal approaches motivates the present study. For the analysis of fluid-structure interaction problems, the fluid and structural domains are decomposed and analyzed using three-dimensional BEM (Boundary Element Method) and FEM, respectively. WISH, a computer program for linear and nonlinear Wave-Induced SHip motion and
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
73
structural loads, is applied for the hydrodynamic portion. The program is based on a time domain Rankine panel method [9]. For the structural part, IPSAP (Internet Parallel Structural Analysis Program) is used, which is a three-dimensional structural analysis program based on FEM [10]. Both methods are fully coupled in the time domain by using an implicit iterative scheme. The developed partitioned approach is verified through comparisons of the wet mode natural frequency, motion responses, and stress histories of a simple barge, with the other computational results obtained by a beam-based solver, WISH-FLEX, which was developed by Kim et al. [11,12]. For validation purposes, the natural frequency and motion responses of a segmented barge are compared with experimental data [13,14]. In order to investigate reliability on a real ship application, the motion responses and stress histories of a 6500 TEU containership are observed through simulation in regular waves. Moreover, the motion responses and structural load of a 10,000 TEU containership are observed. Based on these computational results, the applicability and strengths of the present method for ship hydroelasticity are discussed. 2. Mathematical background 2.1. Fluid-structure interaction problem Let us consider a freely floating body with moving speed (U) in incident waves. A mean-body fixed coordinate system (o xyz), which is advancing with the same speed of the body, can be employed to define the body motion, as shown in Fig. 1. This body has a finite value of hull rigidity, so that the body structure is flexible. SB, SF and SN are the wetted-body surface, free surface, and boundary surface at infinite distance from the body, and, UF and US are fluid and structural domains, respectively. The wave amplitude and frequency of an incident wave are represented by A and u, and the wave heading angle is b. To solve fluid-structure interaction problem, the two domains should be considered at the same time. Both domains can be treated by applying the domain decomposition approach, also known as the “partitioned approach,” where the two field equations are handled separately, allowing only variable exchange between the two fields at every time step. The governing equations of fluid and structural domains can be represented f and s, as follows:
_ u; 4f ¼ 0 f ¼ p F u;
(1)
s ¼ u_ SðpÞ ¼ 0
(2)
Fig. 1. Definition of coordinate system.
74
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
where p is the hydrodynamic pressure acting on the body surface and 4f is the velocity potential on free surface. u and u_ are the deformation displacement and velocity at the considered location. Eqs. (1) and (2) represent that the displacement and velocity are a set of inputs for the fluid field equation, while structural responses are dependent on dynamic loads due to fluid pressure. The dynamic pressure is used as input for structural analysis which provides the kinematic response as a part of the solution. This coupling is essential to the fluid-structure interaction problem.
2.2. Structural domain The equation of structural motion can be derived by the principle of virtual work. Once the virtual displacement, (du), is assumed, the external virtual work should be equal to the internal virtual work.
Z
Z
T
fdug fFV gdV þ V
T
fdug fFS gdS ¼
Z
T T € g þ fdugT kd fug _ dV fdεg fsg þ fdug rB fu
(3)
V
S
where ε and s are strain and stress, respectively. FV is the volume force acting on the volume of the structure, and FS is the surface force acting on the body surface. rB and kd represent the density of the body and the damping coefficient, respectively. The displacements are expressed by a nodal displacement vector and shape function as follows:
fug ¼ ½Nfud g _ ¼ ½Nfu_ d g fug € g ¼ ½Nfu €d g fu
(4)
€ d g are displacement, where [N] is a function of space, the so called shape function, and fud g, fu_ d g and fu velocity, and acceleration vectors at nodes, respectively. By substituting Eq. (4) into Eq. (3), the following equation is derived:
2 T4
Z
½BT fsgdV þ
fdud g
V
Z
rTB ½NT ½NdVfu€d g þ
V
Z
3
Z V
kd ½NT ½NdVfu_ d g
Z
½NT fFV gdV
V
(5)
½N fFS gdS5 ¼ 0 T
S
[B] is the strain-displacement relationship multiplied by the shape function Equation (5) is satisfied for any fdud gT , and this leads to a system of equations as
o n € d g þ ½Cfu_ d g þ F int ¼ F ext ½Mfu
(6)
where [M] and [C] are mass and damping matrices. Those are defined as follows:
Z
rB ½NT ½NdV
(7)
kd ½NT ½NdV
(8)
½M ¼ V
Z ½C ¼ V
The vectors of internal force, fF int g, and external force, fF ext g, are defined as follows:
n
F int
o
Z ¼ V
½BT fsgdV
(9)
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
ext F ¼
Z
½NT fFV gdV
V
Z
½NT fFS gdS
75
(10)
S
The external volume and surface forces are changed to the nodal force by Eq. (10). If a small displacement is assumed, the internal force can be represented by the linear stress–strain relationship as follows:
n
F int
o
¼ ½Kfud g Z
½K ¼
(11)
½BT ½E½BdV
(12)
V
where [K] is the stiffness matrix, and [E] is the stress–strain relationship. Eventually, the equation of motion can be derived as follows:
€ d g þ ½Cfu_ d g þ ½Kfud g ¼ ½Mfu
ext F
(13)
The above equations are derived for an element, so the global equation is formulated by assembling the element equations. 2.3. Fluid domain Under the assumption of incompressible and inviscid flow with irrotational motion, the velocity potential can be introduced, which satisfies the Laplace equation in the fluid domain.
V2 4 ¼ 0
in UF
(14)
The velocity potential satisfies the boundary conditions as follows:
! ! v! v4 u ! ¼ U$n þ $ n on SB vn vt
d þ V4$V ½z zðx; y; tÞ ¼ 0 on SF dt
d4 1 ¼ g z V4$V4 dt 2
on SF
(15)
(16)
(17)
! ! where U is the velocity vector of a body and n is and a normal vector on the body surface. z and g are ! the wave elevation and gravity, respectively. In addition, d=dt ¼ v=vt U $V. The boundary value problem can be linearized by decomposing the total velocity potential and wave elevation as follows:
! ! ! ! 4ð x ; tÞ ¼ Fð x ; tÞ þ 4I ð x ; tÞ þ 4d ð x ; tÞ
(18)
zðx; y; tÞ ¼ zI ðx; y; tÞ þ zd ðx; y; tÞ
(19)
F is the basis velocity potential which satisfies the double-body linearization proposed by Dawson [15]. It can also be Ux in the Neumann-Kelvin formulation. 4I and 4d are incident wave and disturbed wave potentials, respectively. The order of the basis potential is O(1), and the incident and disturbed wave potentials are O(ε).
76
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
The linearized body boundary condition takes the following form [16]:
!
! !
! ! v4d v u ! v4I ! ! ¼ þ ð u $VÞ U VF $ n þ $n U VF $V u $ n on SB vn vn vt
(20)
The linearized free surface boundary conditions are as follows:
!
vzd ! v2 F v4 z þ d þ U VF $VzI U VF $Vzd ¼ vt vz vz2 d
on z ¼ 0
!
! v4d ! vF 1 U VF $V4d ¼ g zd þ U $VF VF$VF þ U VF $V4I 2 vt vt
(21)
on z ¼ 0
(22)
where SB represents the wetted-mean-body surface.
3. Numerical implementation A partitioned method for the fluid-structure interaction problem is widely applied to many fields in engineering and applied science. The partitioned method allows the use of well-established discretization and solution methods for each sub problem. The coupling can be treated fully explicitly or implicitly, or in a mixed explicit-implicit manner. The implicit treatment is applied due to its stability [17,18]. In the present study, the fluid and structural domains are decomposed and handled by threedimensional BEM and FEM, respectively. Those methods are fully (strongly) coupled by using an implicit iterative scheme in the time domain. 3.1. Structural domain The hull motion including elastic deformation is analyzed by using a three-dimensional FE analysis program, IPSAP. It is a large scale, high performance structural analysis program developed by the Aerospace Structures Laboratory in Seoul National University. For the time integration of the equation of motion, the Newmark-b method [19] which takes the following form is applied:
tþDt d 1 C þ K M þ fud gtþDt ¼ F ext aDt aDt 2
1 1 1 t t t _ € u þM þ þ 1 u f f fu g g g d aDt d 2a aDt 2 d
d d Dt d € d gt 1 fu_ d gt þ 2 fu þC fud gt þ aDt a 2 a
(23)
here, t and t þ Dt indicate the current and next time steps, respectively. The velocity and displacement of structure can be approximated at the next time step by using the solutions of the previous time step as follows:
n
n
€ dtþDt u
D u_ dtþ t
o
o
¼
n o 1 n tþDt o t 1 n to 1 _ € td u u 1 u u d d d aDt 2a aDt 2
¼
n o n oi h n o € dtþDt € td þ d u u_ td þ Dt ð1 dÞ u
(24)
(25)
where a ¼ 0:25; d ¼ 0:5. In the present direct integration method, Eq. (23) is directly handled to obtain the equation of motion. The external force is obtained at the time step t þ Dt, therefore this scheme should be considered in an implicit manner. In the present study, the structural damping term is not modeled.
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
77
3.2. Fluid domain To solve the prescribed boundary value problem in the fluid domain, the WISH program for rigid body dynamics is revised. The WISH program solves the motion responses and wave load on a rigid hull by using a higher-order Rankine panel method. By using Green’s 2nd identity, Eq. (14) is converted to the integral equation as follows:
Z 4d þ
4d SB
vG dS vn
Z SF
v4d GdS ¼ vn
Z
v4d GdS vn
SB
Z 4d SF
vG dS vn
(26)
The body surface and free surface near the body are discretized and the Rankine source is distributed on the boundary surface. The B-spline source potential is employed as follows:
P ! ! 4d ð x ; tÞ ¼ 4d;i ðtÞBi ð x Þ i
P zd ð! zd;i ðtÞBi ð! x ; tÞ ¼ xÞ
(27)
i
X v4 v4d ! ! d ð x ; tÞ ¼ ðtÞBi ð x Þ vn vn i i
! where Bi ð x Þ is a B-spline basis function, and 4d, zd and v4d/vn are the potential coefficient, wave elevation coefficient and the normal flux of the potential coefficient, respectively. By substituting Eq. (27) into Green’s second identity and the boundary value problems described above, a matrix equation with unknown coefficients can be assembled. By solving the matrix equation, the normal flux on the free surface and the velocity potential on the body surface can be obtained. In the present study, the kinematic and dynamic free surface boundary conditions are solved using the mixed explicit-implicit Euler scheme proposed by Nakos et al. [20]. In order to satisfy the radiation condition, the concept of numerical damping zone is applied, which applies artificial wave damping near the truncated free surface and damps outgoing waves numerically. Then the free surface boundary conditions on z ¼ 0 are written as follows:
n n znþ1 zd n v4 d ¼ P Fn ; z d ; d Dt vn
n
2nzd þ
n2 g
4nd
4nþ1 4nd nþ1 d ¼ Q zd ; Fnþ1 ; 4nþ1 d Dt
(28)
(29)
where n indicates the numerical damping strength which has non-zero values inside the artificial damping zone. In this study, the hydrodynamic force including the Froude-Krylov force is computed by using the following equation:
v ! ! U VFj $V 4I;j þ 4d;j $ n j Dsj Fj ¼ r vt
v Fj ! 1 ! U $VFj þ VFj $VFj $ n j Dsj r 2 vt
(30)
the subscript j represents the j-th panel, and Dsj is the panel area. In fact, there is extra term which is proportional to the spatial derivative of the second term in Eq. (30). Such extra term vanishes when the Neumann-Kelvin linearization is applied, but it is not zero when the double-body linearization is used. However, according to our experience, the extra component is ignorable compared with Eq. (30) despite the difficulty in numerical evaluation. Therefore, the extra term is not included in the present numerical computation.
78
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
The restoring force, so called hydrostatic stiffness, can be applied as follows:
! ! Fres;j ¼ rgdz$ n j Dsj rgz$d n j Dsj
(31)
! z and n are water head and normal vector at its mean position. d denotes its leading order variation. The first and second terms in Eq. (31) are called pressure variation and normal vector variation, respectively. The evaluation of hydrostatic stiffness acting on a flexible body is not as straightforward as for a rigid body [21–23]. In the present study, only two terms of Eq. (31) are considered by following Kim et al. [11]. The linearized body boundary condition can be derived as follows [24]:
6 x X v j v4d v4 ¼ nj þ xj mj I vn vt vn j¼1
on SB
(32)
equation (32) is derived under the assumption that the ship is a rigid body. For shell element model, it can be assumed that each shell element of body works as rigid one. Based on this, Eq. (32) can be applied, here, x1 ; x2 ; x3 mean the translational motions and x4 ; x5 ; x6 mean the rotational motions of each shell element. This assumption is reasonable because, first of all, the size of shell element is normally very small. Moreover, the deformation gradient of each shell element would be very small thanks to the relatively high stiffness of general ship structure. 3.3. Fluid-structure coupling method The initial boundary value problems of the two different domains can be written as Eqs. (1) and (2). To solve the coupled equations in a fully coupled manner in the time domain, an implicit iterative scheme should be involved. To this end, a fixed-point iteration scheme is introduced. This scheme is simple but powerful to handle the coupling problem in an implicit manner. When an iterative scheme based on fixed-point iteration is considered, iteration starts with the assumed initial pressure field or the converged solution in the previous time step. First, the structural problem is solved with this pressure field and the structural deformation and deformation velocity are obtained, which are subsequently passed to the fluid field equation. Then the new pressure field can be obtained from this equation. This whole procedure is repeated until a solution converges. This procedure is summarized by the following equations:
k u_ kþ1 tþDt ¼ S ptþDt kþ1 kþ1 _ kþ1 ¼ F ; u ; 4 pkþ1 u D t tþ tþDt tþDt f ;tþ Dt iterate until u_ ktþDt u_ kþ1 tþDt < ε
(33)
where k means the number of iterations, and ε indicates an error tolerance. In order to increase the convergence speed and obtain numerical stability, a relaxation scheme presented by Kim [25] is used in this study. During the fixed-point iteration, the external force and displacement of the ship surface should be transferred to the structural grid and panel model, respectively. In this kind of coupling analysis, the accuracy of computational results can be dependent on the accuracy of data transfer between two domains. Moreover, unbalanced force, which can be generated by inaccurate mapping of the wave load to the FE model, can disturb the numerical solution. In the present study, hydrodynamic pressure on the panel model is transferred to the skin mesh of FE model by matching the closest node between the skin mesh and the panel model. If the resolution of the structural mesh is fine enough, this kind of approach can be applicable. Usually, the structural mesh of the ship is made sufficiently fine by comparing with the hydrodynamic mesh. In Fig. 2, the examples of pressure distribution on hydrodynamic panels and FE skin meshes are shown, which look very consistent. The time histories of vertical force obtained from the panel and FE models are shown in Fig. 3, showing good agreement
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
79
Fig. 2. Pressure distribution on panel model and skin mesh of FE model of 6500 TEU containership.
between the forces on the panel and the FE models. There is insignificant unbalanced force in between the vertical forces of the panel and FE models. The transfer of displacement from the FE skin mesh to the hydrodynamic panels is carried out by following the same method of transferring pressure. In Fig. 4, the displacement distributions on the skin mesh and the panel model are shown, and they look reasonably transferred. However, for better accuracy, the interpolation scheme for pressure and displacement mapping should be more accurate.
4. Comparison with beam model for a simple barge 4.1. Test model A simple barge model is applied to verify the developed program through comparison with other numerical results obtained by WISH-FLEX. This program is a computer program for the analysis of ship hydroelasticity problems using a B-spline Rankine panel method and Vlasov beam theory. These are fully coupled in the time domain by using an implicit iterative scheme. The model was systematically validated by Kim et al. [11,12]. The panel and FE models of the barge are shown in Fig. 5, and the principal dimensions of the barge are presented in Table 1. 3700 panels for half domain are distributed for hydrodynamic analysis, while
20000 Panel model FE model
F3 [kN]
10000
0
-10000
-20000 0
20
40
60
80
100
Time [sec] Fig. 3. Comparison of vertical forces on panel and FE models of 6500 TEU containership, Fn ¼ 0.0485, b ¼ 180deg, A/L ¼ 1/286.3, l/L ¼ 1.0.
80
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
Fig. 4. Displacement distribution on panel model and skin mesh of FE model of simple barge, Fn ¼ 0.0, b ¼ 180deg, A/L ¼ 1/60, l/L ¼ 1.0.
the structural meshes have 22,900 nodes and 26,100 elements. The number of nodes is determined through Eigenvalue analysis that converges up to the 20th mode. Shell elements are distributed on the whole outer surface of the FE model. The characteristic of this model can be said to be globally soft and locally stiff like a beam. This characteristic is achieved by distributing very stiff bulkheads inside the barge. The bulkheads have no mass, so these do not affect the mass distribution.
Fig. 5. Panel and FE models of the barge.
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
81
Table 1 Principle dimensions of the barge. Length (L) Breadth (B) Depth (D) Displacement Draft at AP/FP KG
60.0 m 10.0 m 4.0 m 1200 ton 2.0 m/2.0 m 2.0 m
A beam model that has structural responses equivalent to the three-dimensional shell model is generated for comparison. In Table 2, the natural frequencies of the beam and shell models for several mode shapes are presented, and the corresponding mode shapes are shown in Fig. 6. The components from the first to the sixth modes indicate the vibrations of a rigid body, and the higher mode shapes are flexible vibrations. It is clear that the beam and shell FE models have almost equivalent structural responses. 4.2. Natural frequency on wet mode The natural frequency of the wet mode can be evaluated through a hammering test, which involves an arbitrary impact applied to the body in still water. The same magnitude of impact is given to both the beam and shell element models of the barge, and the impact is applied to AP (after perpendicular). The time histories of vertical displacement at FP(forward perpendicular), mid-ship, and AP of the barge are shown in Fig. 7. Some transient signal is found after the impact, and eventually the signal shows a regular vibration. Basically, the oscillating period of vibration is the fundamental natural frequency of the barge in wet condition. The natural frequencies computed by the present method and WISH-FLEX are 0.6084 Hz and 0.6091 Hz, respectively, and those are quite similar. 4.3. Motion response The motion responses of the barge with zero speed computed by the present method are compared with those of WISH-FLEX. The vertical displacements of the barge at FP, mid-ship, and AP are shown in Fig. 8, for when the wave length and ship length are the same. Here, un represents the natural frequency of the body computed in the hammering test, and ue indicates the encounter frequency. With this wave length, the rigid body motion is dominant. The present computational result shows a tendency similar to the result of WISH-FLEX. An asymptotic rigid body oscillation can be computed by increasing the stiffness of elements. By removing the rigid body motion from Fig. 8, pure flexible deformation can be obtained as shown in Fig. 9. The flexible vibration oscillates with the frequency which is higher than the incident wave frequency. The results of the present method show very acceptable agreement with those of WISHFLEX. As time goes by, time lag between the results of both methods is increased. This is quite natural because a small difference can exist between the two beam and shell models. Moreover, difference in numerical approaches between the two methods can be the reason of this time lag. Table 2 Natural frequencies for several modes of beam and shell models (dry mode). Mode no.
7 8 9 10 11 12
Natural frequency (Hz) Beam model (WISH-FLEX)
Shell model (IPSAP)
Shell model (NASTRAN)
0.92 1.72 2.10 2.22 3.76 3.88
0.92 1.72 2.13 2.23 3.77 3.89
0.92 1.72 2.13 2.22 3.77 3.89
82
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
Fig. 6. Natural mode shapes of the barge (computed by IPSAP).
The contours of disturbance waves near the barge are shown in Fig. 10. The wave patterns computed by the present method and by WISH-FLEX look similar, even though the structural analyses were carried out by different manners. If the wave frequency is close to the natural frequency of the ship structure, the body can experience a resonance vibration, or “springing.” The vertical displacements at FP, mid-ship and AP of the barge are shown in Fig. 11. In this case, the incident wave frequency is close to the natural frequency which was obtained in the hammering test. From the motion histories, it is clear that the body experiences resonance vibration of two-node vertical bending. The amplitudes of resonance vibration computed by the two methods are different in this case. In the region near the resonance frequency, the amplitude of motion can be dramatically changed, and therefore the amplitudes obtained by the two methods can be different.
0.02
FP Mid-ship AP
w [m]
0.01 0
-0.01 -0.02 0
5
10
15
20
25
30
35
Time [sec] Fig. 7. Time history of vertical displacement of barge with hammering test.
40
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
83
2 WISH-FLEX Present method
w [m]
1
0
-1
-2 0
10
20
30
40
50
Time [sec]
(a) FP 2 WISH-FLEX Present method
w [m]
1
0
-1
-2 0
10
20
30
40
50
Time [sec]
(b) Mid-ship 2 WISH-FLEX Present method
w [m]
1
0
-1
-2 0
10
20
30
40
50
Time [sec]
(c) AP Fig. 8. Time histories of vertical displacement of barge, Fn ¼ 0.0, b ¼ 180deg, A/L ¼ 1/60, l/L ¼ 1.0 (ue/un ¼ 0.264).
4.4. Stress result An accurate prediction of the stress response of hull girders is essential in the prediction of the fatigue life of a ship. In many previous studies, the stress RAO of the ship was predicted by using beam theory, which involved idealizing three-dimensional structure into one-dimensional beam elements.
84
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
WISH-FLEX Present method
w [m]
0.2
0
-0.2
0
10
20
30
40
50
Time [sec]
(a) FP WISH-FLEX Present method
w [m]
0.2
0
-0.2
0
10
20
30
40
50
Time [sec]
(b) Mid-ship WISH-FLEX Present method
w [m]
0.2
0
-0.2
0
10
20
30
40
50
Time [sec]
(c) AP Fig. 9. Time histories of vertical displacement of barge (flexible component), Fn ¼ 0.0, b ¼ 180deg, A/L ¼ 1/60, l/L ¼ 1.0 (ue/ un ¼ 0.264).
As mentioned, the one-dimensional beam theory is very useful in the initial stages of ship structural design, but it has limitations in the handling of complex hull structure. Three-dimensional FEM and a direct integration scheme are therefore very desirable for an accurate prediction of the stress response of a ship structure.
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
85
Fig. 10. Instantaneous counters of disturbance waves near the barge, Fn ¼ 0.0, b ¼ 180deg, A/L ¼ 1/60, l/L ¼ 1.0 (ue/un ¼ 0.264).
0.1 AP Midship FP
w
0.05
0
-0.05
-0.1 0
10
20
30
40
30
40
Time [sec]
(a) Present method 0.1 AP Midship FP
w [m]
0.05
0
-0.05
-0.1 0
10
20
Time [sec]
(b) WISH-FLEX Fig. 11. Time histories of vertical displacement of barge, Fn ¼ 0.0, b ¼ 180deg, A/L ¼ 1/60, (ue/un ¼ 0.993).
86
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
In order to verify the stress results computed using the present method, the stress of the barge structure in regular waves is compared with the results of WISH-FLEX. Two spots on the barge are considered to compute stress: at 2/3L from AP (Pt.1) and mid-ship (Pt.2) on deck of the barge. The time histories of axial stress on Pt.1 and Pt.2 are shown in Fig. 12. The axial stress of WISH-FLEX is computed by integrating the vertical bending moment and axial force along the x-direction from the ends of beam elements. The stress results computed by both methods show reasonable correspondence. Particularly, oscillations with high frequency components at Pt.1 are quite similar in both the results of the two methods. As can be seen in Fig. 12, a time lag appeared in the time histories of stress. The transient oscillation remains in the stress histories because structural damping is not modeled in the present study. 5. Comparison with experiment for a segmented barge 5.1. Test model
σ
σ
The validation of the developed program can be carried out by comparison with experimental results. In this paper, the comparison with experimental data for a segmented barge model is introduced. The elastic barge has been experimented with by Malenica et al. [13] and Remy et al. [14]. Their experimental model is shown in Fig. l3. This barge was made of 12 segmented pontoons, and each pontoon is 0.19 m long, 0.6 m wide and 0.25 m deep. The overall length is 2.445 m and the draft is
Fig. 12. Comparison of axial stress of barge, Fn ¼ 0.0, b ¼ 180deg, A/L ¼ 1/60, l/L ¼ 1.0(ue/un ¼ 0.264).
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
87
0.12 m. The foremost pontoon is slightly modified as shown in Fig. 13. Twelve segmented pontoons were tied together with two steel plates placed side by side on the deck level in Malenica’s model (Model A). In Remy’s model (Model B), steel plates were replaced by steel rods with the same pontoon arrangement. The barge model is quite flexible in order to produce high levels of hydroelastic phenomena. The panel model and three-dimensional FE model of the segmented barge are shown in Fig. 14. 2600 panels for half domain are distributed for the panel model. The structural mesh has 9500 nodes and 9700 elements. For the computation of hydrodynamics, the gaps between the pontoons are not considered. Therefore, hydrodynamic properties such as added mass, damping and fluid restoring can differ slightly from those of the real situation. 5.2. Asymptotic solution with very large rigidity Fig. 15 shows the comparison of the present results for very high rigidity with those of WISH (i.e. rigid body solution) and experimental data for a rigid body. Reasonable agreement is observed among the three results, and therefore it is clear that the asymptotic solution recovers the rigid body case as the rigidity of the body structure approaches infinity. Some discrepancy with experimental data at low frequency is found, though. 5.3. Natural frequency on wet mode The wet mode natural frequency of the segmented barge is computed through hammering test in still water. An arbitrary vertical impact is given to Pt.12, and the time history of vertical displacement at the mid-ship of the segmented barge is shown in Fig. 16. The oscillating frequency is the natural frequency of two-node vertical bending. The natural frequencies obtained by experiment and computation are 1.063 Hz and 1.010 Hz, respectively. The natural frequency computed by the present method is close to the value measured by experiment, but there is some difference. As mentioned above, the gaps that exist between the pontoons are not modeled in the computation, which can result in some differences in the added mass and restoring. This difference may be a primary source of the difference in natural frequencies.
Pt.11
Pt.9
Pt.7
Pt.5
Pt.3
Pt.1
25 cm
Pt.12
2445 cm
(a) Overall geometry
10 cm
60 cm
19 cm
5 cm 5 cm
5 cm
(b) Foremost pontoon Fig. 13. Experimental model [10,11].
88
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
Fig. 14. Panel and FE models of segmented barge.
5.4. Flexible motion response in waves Senjanovic et al. [26] computed the motion responses and wave loads on the segmented barge by using their beam-based approach. They presented that there are discrepancies between experimental and numerical results near resonance frequency, and they adjusted the numerical results by applying an arbitrary damping. Remy et al. [14] also applied a damping mechanism in their computation in order to consider structural and viscous damping. In the present study, the damping effect is imposed as an external force on the outer surface of the pontoons, which is proportional to the body velocity. The damping force should be imposed, because the viscous damping effect due to fluid motion can be significant in the segmented pontoon model which has sharp edges. Moreover, the gaps between the pontoons can be the source of fluid damping.
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
89
2 6DOF Experiment Present method
w/A [m/m]
1.5
1
0.5
0
2
3
4
5
ω [rad/sec]
6
7
8
(a) At Pt.1 (bow) 2 6DOF Experiment Present method
w/A [m/m]
1.5
1
0.5
0
2
3
4
5
ω [rad/sec]
6
7
8
(b) At Pt.7 (mid-ship) Fig. 15. Vertical response of rigid body, Fn ¼ 0.0, b ¼ 180deg.
0.002
w [m]
0.001
0
-0.001
-0.002 0
5
10
15
Time [sec] Fig. 16. Time history of vertical displacement at mid-ship of segmented barge with hammering test.
90
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
3
FP Midship AP
w/A [m/m]
2 1 0 -1 -2 -3 0
5
10
15
20
Time [sec]
(a) Without damping 3
FP Midship AP
w/A [m/m]
2 1 0 -1 -2 -3 0
5
10
15
20
Time [sec]
(b) With damping Fig. 17. Time histories of vertical motion, Fn ¼ 0.0, b ¼ 180deg, u ¼ 6.0 rad/s (ue/un ¼ 0.9).
Fig. 17 shows the vertical displacement of the segmented barge in regular wave conditions. The wave frequency is close to the natural frequency of the barge structure. The motion response shows two-node vertical bending vibration as naturally expected. If the damping force is applied to computation, the motion amplitude is decreased.
Pressure variation Normal vector variation
40
F3 [N]
20 0
-20 -40
0
5
10
15
20
Time [sec] Fig. 18. Time histories of restoring force on segmented barge, Fn ¼ 0.0, b ¼ 180deg, A/L ¼ 0.0041, u ¼ 6.0 rad/s (ue/un ¼ 0.9).
Experiment (Remy et al., 2006) Computation (restoring1+damping) Computation (restoring2+damping) Computation (restoring2)
2.5
2.5
2
w/A [m/m]
w/A [m/m]
2
1.5
1.5
1
0.5
0.5
4
6
8
0
10
4
ω [rad/sec]
(a) At Pt.1 (bow) 2.5
w/A [m/m]
w/A [m/m]
2
1.5
1
1
0.5
0.5
4
6
ω [rad/sec]
8
10
0
4
6
8
10
ω [rad/sec]
(d) At Pt.12 (stern)
Fig. 19. Vertical response of segmented barge, Fn ¼ 0.0, b ¼ 180deg.
91
(c) At Pt.9
10
Experiment (Remy et al., 2006) Computation (restoring1+damping) Computation (restoring2+damping) Computation (restoring2)
2.5
1.5
0
8
(b) At Pt.3
Experiment (Remy et al., 2006) Computation (restoring1+damping) Computation (restoring2+damping) Computation (restoring2)
2
6
ω [rad/sec]
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
1
0
Experiment (Remy et al., 2006) Computation (restoring1+damping) Computation (restoring2+damping) Computation (restoring2)
92
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
Fig. 20. Panel and FE models of a 6500 TEU containership, WILS JIP Phase-I model
Fig. 21. Segmented model of 10,000 TEU containership in WILS JIP Phase-II and comparison of size with 6500 TEU containership.
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
93
As mentioned above, the application of hydrostatic stiffness for hydroelasticity problems is not straightforward. In the present study, the pressure variation and normal vector variation is only considered as shown in Eq. (31). In Fig. 18, the magnitude of each variation term is presented. The magnitude of normal vector variation is quite smaller than the one of pressure variation. It can be said that the normal vector variation does not play an important role in the prediction of vertical motion response. The vertical motion RAOs at each measuring point are shown in Fig. 19. Here, ‘restoring1’ means that only the pressure variation is considered, while ‘restoring2’ means that the normal vector variation is also considered, and “damping” means that artificial damping is applied. The damping force is tuned to achieve minimal discrepancies between the numerical and experimental results. Without damping, the computational results show large magnitudes of vertical displacement at Pt.1 and Pt.3 in a frequency range close to the natural frequency (un ¼ 6.66 rad/s). By applying the damping force, better agreements are observed between the experimental and numerical results. Particularly, the vertical displacement at u ¼ 6.0 rad/s is dramatically decreased. In the real situation, the motion of the segmented barge can be affected by structural damping as well as fluid damping. Based on the present results, it can be said that the damping plays an important role for this segmented barge. At Pt.12, the discrepancy between the experimental and numerical results looks larger than in other locations. It can be understood that the aft pontoon is affected by more severe damping during heave and pitch motions in head sea condition. The effect of normal vector variation is not severe like in Fig. 18.
4 WISH-FLEX Present method
w [m]
2
0
-2 0
20
40
60
80
100
Time [sec]
(a) FP 4 WISH-FLEX Present method
w [m]
2
0
-2 0
20
40
60
80
100
Time [sec]
(b) AP Fig. 22. Vertical displacement of 6500 TEU containership, Fn ¼ 0.0485, b ¼ 180deg, A/L ¼ 1/286.3, l/L ¼ 1.0.
94
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
6. Real ship application The reliability of the developed program for application to a real ship is investigated by observing the motion and stress responses of a 6500 TEU containership and the motion and structural load of a 10,000 TEU containership. This model was tested in WILS JIP Phase I and II [27,28] for the evaluation of nonlinear motions and wave loads. In WILS JIP Phase I, the 6500 TEU containership was tested in rigid body condition, while the 10,000 TEU containership was tested in WILS JIP Phase II, which is a segmented model with elastic beam. The solution panels and FE models of the 6500 TEU containership
(a) Location of stress computation 6E+07 WISH-FLEX Present method
Normal stress
4E+07 2E+07 0
-2E+07 -4E+07 -6E+07 0
10
20
30
40
50
Time [sec]
(b) Normal stress in longitudinal direction
von Mises stress
6E+07 WISH-FLEX Present method 4E+07
2E+07
0 0
10
20
30
40
50
Time [sec]
(c) Von Mises stress Fig. 23. Time histories of stress of 6500 TEU containership, Fn ¼ 0.0485, b ¼ 180deg, A/L ¼ 1/286.3, l/L ¼ 1.0.
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
95
applied to the present computation are shown in Fig. 20. The number of panels is 7200 for a full domain. The structural meshes have 13,700 nodes and 39,500 elements. The segmented model of the 10,000 TEU containership tested in WILS JIP Phase-II is shown in Fig. 21. The advantage of the present method is the direct computation of stress at any location of ship structure, regardless the geometric complexity. However, the available experimental stress data of real ship structure, particularly which can be compared with the results of the present study, is very hard to find. In many experimental studies which backbone models are used, sectional structural loads such as vertical bending and torsion moment have been of primary interests. Therefore, the direct comparison of local stress between computation and experiment is not an easy task. In this study, the stress signal of a 6500 TEU containership is compared with the solution of 1D FEM, and the structural load on the 10,000 TEU containership is compared with experimental data of WILS JIP Phase-II. In Fig. 22, the vertical displacements of the 6500 TEU containership at FP and AP are presented. The normalized wave length is 1.0, and therefore the rigid motion is dominant in this wave range. The results of the present method and WISH-FLEX show a similar tendency. The time histories of normal stress in the longitudinal direction and von Mises stress are presented in Fig. 23. The wave condition is the same as in Fig. 22. The stresses are computed at a hull point near mid-ship, as shown in Fig. 23(a). The normal stress and von Mises stress computed by the present method and WISH-FLEX show favorable agreement. On a hot spot where stress is concentrated, the stress computed by the beam-based approach and this three-dimensional FEM can be more different. The wet mode natural frequency of two-node vertical bending can be roughly estimated by using infinite-frequency added mass, and which was found to be 3.77 rad/s in the present ship. In the case that the wave frequency is close to the natural frequency, resonant vibration should appear, as well shown in Fig. 24. In Fig. 25, the deformation and distribution of von Mises stress are presented for the wave condition of Fig. 24. Large stress is found near mid-ship on the side hull. By considering the deformation of the ship, the stress is reasonably distributed. As shown in Fig. 25, the three-dimensional FEM allows direct access to the structural responses (stress and strains) at any required location. Moreover, structural discontinuity such as in a large opened deck can be directly considered. By increasing the resolution of the structural mesh, a localized stress on an arbitrary hot spot can be accurately evaluated using this method. The wet mode natural frequency of 10,000 TEU containership was obtained in the WILS JIP Phase-II. Fig. 26 shows the vertical displacement of the 10,000 TEU containership in a hammering test. Two-node vertical vibration is observed which is more complicated than the barge models. The natural frequency evaluated from the numerical simulation is about 0.441 Hz. The targeted wet natural frequency of WILS JIP Phase-II was 0.45 Hz. Therefore, the present method provides reliable solution.
0.1
AP Midship FP
w [m]
0.05 0
-0.05 -0.1 0
5
10
15
20
25
30
Time [sec] Fig. 24. Vertical displacement of a 6500 TEU containership, Fn ¼ 0.0485, b ¼ 180deg, A/L ¼ 1/286.3, ue/un ¼ 0.998.
96
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
Fig. 25. Deformation and von Mises stress on a 6500 TEU containership, Fn ¼ 0.0485, b ¼ 180deg, A/L ¼ 1/286.3, ue/un ¼ 0.998.
Figs. 27 and 28 show the motion response and structural load of the 10,000 TEU containership in waves. The forward speed of ship is 5knots. The computation results of the present method show reasonable correspondences with the experimental data and the other numerical solution. The present partitioned method for the fluid-structure problem of ships is quite time-consuming, and therefore it is not yet preferable for practical purposes. The CPU time of WISH-FLEX takes several hours, usually. However, the CPU time of the present method takes several days or one week, depending on the number of solution panels and FE meshes and test condition, by using an Intel Xeon
Fig. 26. Time history of vertical displacement of 10,000 TEU containership with hammering test.
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
97
Fig. 27. Motion responses of 10,000 TEU containership at mid-ship, ship speed ¼ 5knots, b ¼ 180deg.
CPU which is one of the fastest CPU in the world. Moreover, the parallel computation is embedded in the FEM solver of present method. Nevertheless, this method can be used to investigate the complicated hydroelastic behaviors of ships and offshore structures, which cannot be approximated by beam elements. Furthermore, in the future when more powerful computational capacity is available, this method will surely be applied for practical ship design.
98
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
Fig. 28. Vertical bending moment of 10,000 TEU containership at mid-ship, ship speed ¼ 5knots, b ¼ 180deg.
7. Conclusions A fully coupled three-dimensional BEM-FEM approach for hull girder hydroelasticity was introduced. This method is based on time domain formulation, adopting a higher-order Rankine panel method for hydrodynamic problems and three-dimensional FEM for structural problems. The solvers of fluid and structure domains are coupled by using an implicit iterative scheme. From numerical tests, validation, and applications for real merchant ships, the following findings and conclusions are made: The present time domain method can be a good analysis tool to understand global and local hydroelastic behaviors of ship structure. Particularly, the present method provides direct access to structural responses at any required position. Moreover, a complicated hull structure can be directly handled. It is observed that the differences in global motion and stress responses in head sea conditions are not very significant between the beam element based model and three-dimensional shell element model applied in this study. This implies that the beam-based approach is valid for global ship responses, but the prediction of local structural responses requires the present three-dimensional approach. In the analysis of the segmented barge, an arbitrary damping should be involved for an accurate prediction of the motion response and hydroelastic behavior of flexible bodies. The pressure variation term of hydrostatic stiffness plays an important role rather than the normal vector variation in the prediction of vertical motion response. Only head sea condition was considered. Further research on oblique sea cases is required in order to investigate more complicated hydroelastic behavior. Second-order springing resonance is mainly caused by nonlinear force acting on the hull structure, and therefore the present method must be extended to nonlinear analysis solvers. Acknowledgments Authors appreciate Prof.Molin and Dr.Malenica for their kindness to provide their experimental data. Their data were very helpful for the present study. This study has been carried out as a part of a
K.-H. Kim et al. / Marine Structures 33 (2013) 71–99
99
project jointly funded by The LRET*-Funded Research Center at SNU for Fluid-Structure Interaction and Office of Naval Research (Grant No. N00014-08-1-1085), and also as a part of WISH-FLEX JIP funded by Daewoo Shipbuilding and Marine Engineering, Hyundai Heavy Industry, Korean Register, Samsung Heavy Industry, STX offshore and shipbuilding, Korea. Their support is acknowledged. (*LRET: The Lloyd’s Register Education Trust) Also the administrative support of the Engineering Research Institute (ERI) and Research Institute of Marine System Engineering (RIMSE) should be acknowledged. References [1] Bishop RED, Price WG. Hydroelasticity of ships. London: Cambridge University Press; 1979. [2] Jensen JJ, Pedersen PT. Bending moments and shear forces in ships sailing in irregular waves. Journal of Ship Research 1981;24(4):243–51. [3] Wu MK, Moan T. Linear and nonlinear hydroelastic analysis of high-speed vessel. Journal of Ship Research 1996;40(2): 149–63. [4] Taghipour R, Perez T, Moan T. Hybrid frequency-time domain models for dynamic response analysis of marine structures. Ocean Engineering 2008;35:685–705. [5] Taghipour R, Perez T, Moan T. Time-domain hydroelastic analysis of a flexible marine structure using state-space models. Journal of Offshore Mechanics and Artic Engineering 2009;131. [6] Shao YL, Faltinsen OM. Towards development of a nonlinear perturbation method for analysis of springing of ships. Korea. In: Proceedings of the 23rd international workshop on water waves and floating bodies 2008. [7] Malenica S, Tuitman JT, Bigot F, Sireta FX. Some aspects of 3D linear hydroelastic models of springing. France. In: Proceedings of the 8th international conference on hydrodynamics 2008. [8] Oberhagemann J, Moctar O. Numerical and experimental investigations of whipping and springing of ship structures. US. In: Proceedings of the 21th international offshore and polar engineering conference 2011. [9] Kim KH, Kim Y, Kim Y. WISH JIP project report. Seoul, Korea: Seoul National University; 2008 [written in Korean]. [10] http://ipsap.snu.ac.kr, July 4, 2012. [11] Kim Y, Kim KH, Kim Y. Analysis of hydroelasticity of floating ship-like structure in time domain using a fully coupled hybrid BEM-FEM. Journal of Ship Research 2009;53(1):31–47. [12] Kim Y, Kim KH, Kim Y. Springing analysis of seagoing vessel using fully coupled BEM-FEM in the time domain. Ocean Engineering 2009;36:785–96. [13] Malenica S, Molin B, Senjanovic I. Hydroelastic response of a barge to impulsive and non-impulsive wave loads. UK. In: Proceedings of the 3rd international conference on hydroelasticity in marine technology 2003. [14] Remy F, Molin B, Ledoux A. Experimental and numerical study of the wave response of a flexible barge. China. In: Proceedings of the 4th hydroelasticity in marine technology 2006. [15] Dawson CW. A practical computer method for solving ship-wave problem. US. In: Proceedings of the 2nd international conference on numerical ship hydrodynamics 1977. [16] Timman R, Newman JN. The coupled damping coefficients of symmetric ships. Journal of Ship Research 1962;5(4):34–55. [17] Matthies HG, Steindorf J. Partitioned but strongly coupled iteration schemes for nonlinear fluid-structure interaction. Computers & Structures 2002;80:1991–9. [18] Longatte E, Verreman V, Souli M. Time marching for simulation of fluid-structure interaction problems. Journal of Fluids and Structures 2009;25:95–111. [19] Bathe KJ. Finite element procedures in engineering analysis. New Jersey: Prentice-Hall; 1982. [20] Nakos DE, Kring D, Sclavounos PD. Rankine panel methods for transient free-surface flows. Iowa, US. In: Proceedings of the 6th International conference on numerical ship hydrodynamics 1994. [21] Riggs HR. Comparison of formulations for the hydrostatic stiffness of flexible structures. Journal of Offshore Mechanics and Arctic Engineering 2009;131(2). [22] Senjanovic I. Reply to Prof. Riggs’ discussion on paper “Investigation of ship hydroelasticity” by I. Senjanovic, S. Malenica, S. Tomasevic. Ocean Engineering 2008;35:1287–8. [23] Senjanovic I, Tomic M, Tomasevic S. An explicit formulation for restoring stiffness and its performance in ship hydroelasticity. Ocean Engineering 2008;35:1322–38. [24] Ogilvie TF, Tuck EO. A rational strip theory of ship motion: part I. Report No 013. Michigan, US: University of Michigan; 1969. [25] Kim Y. Time-domain analysis on hull-girder hydroelasticity by fully coupled BEM-FEM approach. PhD dissertation. Seoul, Korea: Seoul National University; 2009. [26] Senjanovic I, Malenica S, Tomasevic S. Investigation of ship hydroelasticity. Ocean Engineering 2008;35:523–35. [27] MOERI/KORDI. Report of wave induced loads on ships joint industry project. Daejeon, Korea: Maritime and Ocean Engineering Research Institute; 2007 [limited distribution]. [28] MOERI/KORDI. Report of wave induced loads on ships joint industry project-II. Daejeon, Korea: Maritime and Ocean Engineering Research Institute; 2010 [limited distribution].