# On the Nonlinear Vibrational Responses of a Large Vessel with a Broad Bow Flare under Wave Excitation: Theory and Experiment

Authors: Haicheng Yu; Huilong Ren; Yingbo Xu; Tiange Li
Journal: Shock and Vibration
 (2018)

---

## Abstract

A fully coupled nonlinear three-dimensional (3D) hydroelastic method is developed to investigate vibrational responses of a large ship with a pronounced bow flare subjected to high seas. This numerical model consists of a 3D boundary element method, 1D Euler-Bernoulli beam model, and a 2D generalized Wagner model. Green water loads were considered. Experimental study was carried out in a towing tank on a self-propelled segmented model with nonuniform steel backbones. The ship model was tested in regular incident waves of large amplitude. Impact pressure and nonlinear vertical bending moments were measured and compared with numerical predictions. The proposed nonlinear model produced similar results to the experimental model. Furthermore, the effects of elastic modes and nonlinearities on the numerical results were analyzed.


---

## Body

## 1. Introduction

Detailed knowledge of wave-induced structural vibration of marine vessels sailing in seaways, more often than not in adverse weather, is of pivotal significance for the evaluation of seakeeping ability. Large vessels, their lowest eigenfrequencies of elastic vibration lying in the vicinity of the harmonics of ocean waves, are highly susceptible to severe transient vibrational responses (whipping) caused by wave impact on either the bow or stern region. This occurs when pitching motion and forward speed lead to emergence of the ship hull above the water surface and subsequent penetration into the wave surface. A pronounced bow flare will further aggravate the severity of whipping due to the large bow-wave contact area. Since the structural vibration of large vessels significantly affects the surrounding fluid flow and therefore the hydrodynamic forces, hydroelastic analysis should be used to analyze ship vibratory behaviors in waves.Beginning with the pioneering work by Bishop and Price [1], significant developments have been achieved in theoretical and experimental methods of predicting the hydroelastic loads acting on ships. Wu et al. [2] extended the 2D linear hydroelastic theory to analyze the vibrational responses of slender ships in waves by coupling Euler-Bernoulli beam theory with linear 2D strip theory. Hermundstad et al. [3] proposed a 2D linear hydroelastic method to study the vibrational responses of high speed ships. Price and Wu [4] presented a 3D linear hydroelastic method to study structural responses of a small waterplane area twin hull (SWATH) in which 3D effects near the bow and stern cannot be omitted.Recently the occurrence of extreme sea conditions is more frequently reported than before, thus highlighting the importance of nonlinear hydroelastic theories. Mikami and Kashiwagi [5] proposed a 2D nonlinear hydroelastic theory to study the characteristics of ship hull vibration. The calculations of incident wave forces (Froude-Krylov forces) and hydrostatic forces are performed for the instantaneous wet hull. In order to further analyze whipping responses, Rajendran et al. [6] presented a comparison of 2D nonlinear computational results and experimental results. The bow flare slamming loads were calculated based on von Kármán model.Y. Kim and J.-H. Kim [7] conducted a benchmark study on motions and loads of a 6750-TEU containership to compare different numerical methods ranging from 2D theories to 3D theories. On average, good agreement was found between computational and experimental results with respect to linear responses. However, for larger wave amplitudes, the computational results were more scattered. Tuitman [8] developed a practical method to calculate the response of the ship structure due to the seakeeping loading including slamming loads. Kim et al. [9] presented a 3D nonlinear hydroelastic method to study whipping responses by coupling a 3D finite element method (FEM) with a 3D Rankine panel based boundary element method (BEM). The nonlinear Froude-Krylov forces and hydrostatic forces were computed from the instantaneous wetted surface. The bow flare impact forces were computed based on the 2D generalized Wagner method. The 2D generalized Wagner method is much more efficient than 3D CFD models and is well tested by drop tests of ship cross sections [10]. However, more experimental validation is needed for the whole elastic vessel model in waves.Experimental studies play a vital role not only in validating applicability of numerical methods but also in discovering new phenomena [11]. Zhu and Moan [12] conducted an experimental study on whipping responses of two large containerships and the experimental results were compared with empirical formulas. The discrepancy found implied that more data from real ships was needed to further improve the accuracy of empirical methods. Jialong et al. [13] presented an experimental and numerical research of hull girder vibration of a large ship. The impact wave forces were calculated based on momentum impact theory. The instantaneous wetted surface was approximately obtained by intersecting the original hydrodynamic mesh with the still water surface. Jialong et al. [14] conducted a large-scale model test in real coastal sea conditions. The bow flare slamming and the subsequent whipping phenomena were captured and analyzed.The gradual march towards a complete understanding of fluid-structure interaction is ongoing. Although various formulations have been developed, some open questions remain on the nonlinear slam-induced responses of large vessels in severe seas, and the available experimental data are not very large. In the present research, a nonlinear time domain hydroelastic method was used to simulate whipping of a large ship in high seas. A 3D BEM, 2D generalized Wagner method, 1D beam model, and 1D dam breaking model were fully coupled in this numerical model. Experimental test of a large ship in severe seas was carried out to validate the numerical method. The numerically and experimentally acquired slamming pressures were compared. A comprehensive comparison and analysis of the measured and simulated structural response of the vessel were provided. In addition, the effects of nonlinearity and elasticity on the response were analyzed.


## 2. Theory

Figure1 demonstrates the general formulation of the problem of wave-induced structural vibration of vessels sailing in regular incident waves. The vessel initially travels at a constant velocity U. The fluid domain is assumed to be infinite. Two inertial frames of reference are introduced to analyze this hydroelastic problem: (1) the earth fixed Cartesian frame Ox∗y∗z∗ with its coordinate plane Ox∗y∗ positioned at the calm waterplane and (2) the translational Cartesian frame oxyz which translates invariantly at the velocity U and remains parallel to Ox∗y∗z∗. Then oxyz is appropriate to describe motion variations of the vessel. At the initial instant, the two frames are identical in space. The Galilean transformation relation between these two frames at time t is shown as follows:(1)xyz=cos⁡βsin⁡β0-sin⁡βcos⁡β0001x∗y∗z∗+-Ut00.Figure 1
Problem description of wave-induced structural vibration.In (1), β is the heading angle and in head seas β is defined to be 0. In Figure 1, ΩF is the infinite fluid domain, R is the infinite radius of ΩF, H is the infinite depth of ΩF, S∞ is the boundary of ΩF, and SF is the free surface. ΩS denotes the structural domain and SB denotes the contact surface between ΩS and ΩF.Figure2 presents three representative states of the vessel: initial state, static hydroelastic state, and the transient hydroelastic state. The initial state corresponds to the initially undeformed vessel in calm water, which serves as the input for the computation of the following static hydroelastic state. The static hydroelastic state represents the deformed configuration of the vessel in calm water due to the uneven distribution of mass. The static hydroelastic state is the ship configuration based on which some hydrodynamic forces are computed. The transient hydroelastic state is the instantaneous deformed configuration of the vessel in waves. Some hydrodynamic forces are calculated based on the transient hydroelastic state, and these are explained later in this section. Figure 2(c) shows the transient hydroelastic state in severe sea conditions, where slamming events and green water phenomena occur.Figure 2
Three representative states of the vessel: initial state (a), static hydroelastic state (b), and transient hydroelastic state (c).
(a) (b) (c)


### 2.1. Structural Domain

The vessel structure as a whole is assumed to be homogeneous, isotropic, and elastic. Instead of adopting complicated and time-consuming 3D FE model, the elastic hull is modeled as a nonuniform Euler-Bernoulli beam. The beam model is more suitable, especially when 3D FE model is not available. In the inertia translational reference frameoxyz, the equation of motion of the elastic vessel can be written as follows [15]:(2)MuN¨+CuN˙+KuN=F,where M, C, and K are mass, damping, and stiffness matrices and F is the external force matrix. uN, uN˙, and uN¨ are the nodal displacement, velocity, and acceleration matrices, respectively.Based on the principle of superposition, nodal variablesuN, uN˙, and uN¨ can be expressed as a sum of the distortion in the principal modes of the vessel:(3)uN=BP,uN˙=BP˙,uN¨=BP¨,where B denotes the principal mode matrix and P denotes the principal coordinate matrix.Substituting (3) into (2), multiplying by BT, the equation of motion expressed in terms of principal coordinates is obtained:(4)mP¨+cP˙+kP=f,where m=BTMB, c=BTCB, k=BTKB, and f=BTF. m, c, k, and f are called generalized mass, damping, stiffness, and external force matrices, respectively.


### 2.2. Fluid Domain

#### 2.2.1. 3D Components of Wave Forces

By assuming that the fluid flow is nonviscid, incompressible, and irrotational, the potential flow theory can be employed to model the fluid in the problem of fluid-structure interaction. In the reference frameoxyz, the fluid velocity potential ∅ of the fluid-structure system in steady state can be written as [16](5)∅x,y,z,t=Re⁡φx,y,zeiωt,where φ is the spatial part of ∅, i is the imaginary unit, and ω is the wave encounter frequency. ω is defined as(6)ω=ω0+Uk0cos⁡β,where ω0 is the natural wave circular frequency and k0 is the wave number.φ can be further decomposed into three parts:(7)φx,y,z=φIx,y,z+φDx,y,z+∑m=1hφRmx,y,z,where φI denotes the incident wave potential, φD denotes the diffraction wave potential, and φRm denotes the radiation wave potential of the mth degree of freedom of the vessel. h is greater than 6 and indices m>6 represent the elastic degrees of freedom. φI, φD, and φR can be solved based on the boundary element method [16].The generalized diffraction forcesFD and radiation forces FR are calculated based on the wet contact surface of the fluid-vessel system in the static hydroelastic state as illustrated in Figure 2:(8)FD=-ρξ∬n⃑·bmiω-U∂∂xφDdS,FR=-ρ∑k=1h∬n⃑·bmiω-U∂∂xφRkpkdS,where ρ is the fluid density, ξ is the incident wave amplitude, n⃑ is the normal vector pointing into the fluid domain, bm is the mth principal mode of the vessel, and dS is the infinitesimal area element of the wet contact surface.The generalized Froude-Krylov forcesFI and hydrostatic restoring forces FHS are computed based on the wet contact surface of the fluid-vessel system in the transient hydroelastic state as demonstrated in Figure 2:(9)FI=-ρξ∬n⃑·bmiω-U∂∂xφIdS,FHS=-ρg∑k=1hpk∬n⃑·bmwkdS-∭ρsgwmdV,where wm is the vertical displacement of the mth degree of freedom and ρs is the mass density of the vessel.


#### 2.2.2. Impact Wave Forces

Slamming loads acting on deformable large vessels with considerable bow flares are calculated based on the 2D generalized Wagner method [17]. Some selected cross sections are chosen to calculate the 2D slamming pressure [10]. Figure 3 displays the description of 2D water entry of cross sections. vt is the velocity of the section, bt and ht are the coordinates of the contact point at time t, and fys is the contour of the cross section. In the present study, vt is the combination of rigid and elastic body motions of the vessel. Note that the effects of pitch on the motion of those cross sections are confined to the vertical components because the ship is slender.Figure 3
Impact problem description of symmetric water entry.The 2D potential flow theory is used to model this fluid-structure problem and the slamming pressure acting on the wet contour of each cross section can be solved by the conformal mapping method. Figure4 demonstrates how conformal mapping is performed. The real free surface contour is linearized as the horizontal line starting from the contact point. The line segments AB and CD from the original real plane ysozs are mapped to intervals -∞,-1 and +1,+∞ in the complex plane ξoη, respectively. The wet contour BC is mapped to the interval -1,+1. The slamming pressure can be obtained as follows:(10)pSLAMξ,t=ρdstdt2pbξ,t+dstdtpaξ,t,where st is the vertical displacement of the cross section. The expressions of paξ,t and pbξ,t can be found in [17]; here we do not explain these two expressions for brevity.Figure 4
Conformal mapping ((a) physical plane and (b) corresponding complex plane).
(a) (b)The solutionpSLAMys,zs,t in the original ysozs plane can be obtained based on the inverse of the conformal mapping. Then the generalized slamming forces can be calculated as(11)FSL=∫n⃑·bmpSLAMdl,where n⃑ is the normal vector pointing out of the wet contour of the cross section and dl is the infinitesimal line segment of the wet contact line of the section.


#### 2.2.3. Green Water Loads

In extreme sea conditions, after the bow impact events, the incident wave may overtop the bow and run up onto the deck. In this paper, green water loads on deck are incorporated and simulated based on the dam-breaking Model. This method models the green water phenomenon as the process of sudden collapse of a dam. Before the sudden collapse, the water stored in the dam is static. Further we assume that during the process the ground exerts no resistant forces to the fluid. We introduce a 2D bow-fixed reference framexGozG in which the origin is located at the tip of the bow. The simplified Saint-Venant equations [18] can be expressed as(12)∂HxG,t∂t+HxG,t∂SxG,t∂xG+SxG,t∂HxG,t∂xG=0,1g∂SxG,t∂t+SxG,tg∂SxG,t∂xG+∂HxG,t∂xG=0,where H and S are the depth and velocity of the green water flow, respectively. g is the gravitational acceleration. Based on the method of characteristics, the depth H and velocity S can be solved as follows:(13)HxG,t=-xG3gt+2H032,(14)SxG,t=23xGt+gH0,where H0 is the initial depth of the water in the intact dam. Velocity S is a linear function of xG, and green water pressure on deck pGW can be easily calculated from S. We assume that pGW does not change in yG direction; then the generalized green water loads can be calculated as(15)FGW=∬n⃑·bmpGWdS.


### 2.3. Hydroelastic Vibration Domain

Based on the analyses of the structural and fluid domain, now the hydroelastic equation of motion can be rewritten as follows:(16)mP¨+cP˙+kP=fI+fD+fR+fHS+fSLAM+fGW,where the right-hand side (RHS) of (4) is decomposed into the six components shown in the RHS of (16). fI, fD, and fR are the generalized incident, diffraction, and radiation wave force matrices, respectively, fHS is the generalized hydrostatic restoring force matrix, and fSLAM and fGW are the generalized slamming and green water force matrices, respectively.For large elastic vessels with pronounced bow flare in severe sea conditions, it is of great importance to appropriately incorporate some nonlinearities in those generalized hydrodynamic force matrices of (16). This may not be purely mathematically consistent because (16) is basically linear. However, this treatment does allow evident improvements to the simulation of fierce fluid-structure interaction. As formulated in Section 2.2, among the six generalized hydrodynamic force matrices, only fD and fR are calculated based on the unchanged static hydroelastic state of the vessel. This is because of the difficulty of recalculating φD and φRn using 3D BEM based on the updated wet configuration of the vessel. In contrast, the other four generalized forces are calculated based on the actual instantaneous wet contact surface because the computational burden is acceptable. In other words, except for fD and fR, all the other generalized external forces are nonlinear.


### 2.4. Necessary Numerical Details

#### 2.4.1. Mesh Generation

Mesh generation is the prerequisite of the computation of the generalized external forces. The hull mesh for the initial state of the vessel is generated by employing the method of cubic spline function [19]. Then the transient elastic hull mesh can be generated according to the superposition of the distortion in the principal modes. Figure 5 displays the initial hull mesh and the transient elastic hull mesh.Figure 5
Mesh generation ((a) initial hull mesh and (b) transient elastic hull mesh).
(a) (b)


#### 2.4.2. NonlinearfI and fHS

The nonlinear generalized force matricesfI and fHS are calculated based on the instantaneous wet contact surface and therefore the acquisition of the exact wet hull plays a significant role. The method used in the present study is called the coarse-fine mesh strategy (see Figure 6), which is illustrated as follows: the initial hull mesh is classified as the coarse mesh. The coarse mesh panels are checked at every time step for their position relative to the incident wave profile. Those panels whose vertices are all below (or above) the incident wave contour are considered as the ordinary fully wet (or dry) panels. Special care is given for those partially wet ones. Each partial wet panel is further subdivided into smaller subpanels which are called the fine mesh. Then the reevaluation of the relative position of the fine mesh with respect to the incident wave profile is carried out. The incident wave profile within every fine mesh panel is approximated as the straight line. Based on this strategy, fI and fHS are calculated by the formulations given in Section 2.2.Figure 6
Coarse-fine mesh strategy for computingfI and fHS.


#### 2.4.3. NonlinearfSLAM

The extraction of the cross sections of the vessel can be carried out based on the information of the hull geometry as shown in Figure7. For bow and stern areas, more cross sections should be extracted not only due to the fact that the slamming phenomena are more notable in these two areas but also for the compensation of the evident 3D effects associated with these regions.Figure 7
Extracted 2D sections along the hull.As mentioned above, the generalized Wagner method assumes that there is no flow separation in the slamming phenomena. Therefore modifications should be made for those sections whose breadth of the outer contour does not increase monotonically with respect to the height of the contour. Figure8 takes one of the cross sections of the bulb bow as an example to demonstrate those modifications.Figure 8
Modification of original outer contour of the bulb bow.For the cross sections where the outer contour consists of two disjoint curves as shown in Figure8, we use the smooth curve called the connector to connect the two separated modified curves. The wet contour of the connector is excluded in the computation of fSLAM.


#### 2.4.4. NonlinearfGW

As formulated in Section2.2.3, the green water phenomenon is assumed to be one-dimensional. Figure 9 shows green water depth distribution at time t along the deck. The mesh for the deck is generated as shown in Figure 5. H0 in (14) is determined by the relative position of the deck and the incident wave profile at time t. Then fGW is computed on the wet deck mesh based on (15).Figure 9
Green water depth distribution along the deck.


#### 2.4.5. Hydroelastic Response

In the present study, we adopt the well tested explicit fourth-order Runge-Kutta algorithm [20] to solve (16). Then the hydroelastic response can be acquired based on the principle of modal superposition:(17)wx,t=∑k=1hpktwkx,Mx,t=∑k=1hpktMkx,Vx,t=∑k=1hpktVkx,where w, M, and V are the displacement, vertical bending moment, and shear force; wk, Mk, and Vk are the kth mode of displacement, bending moment, and shear force, respectively.


## 2.1. Structural Domain

The vessel structure as a whole is assumed to be homogeneous, isotropic, and elastic. Instead of adopting complicated and time-consuming 3D FE model, the elastic hull is modeled as a nonuniform Euler-Bernoulli beam. The beam model is more suitable, especially when 3D FE model is not available. In the inertia translational reference frameoxyz, the equation of motion of the elastic vessel can be written as follows [15]:(2)MuN¨+CuN˙+KuN=F,where M, C, and K are mass, damping, and stiffness matrices and F is the external force matrix. uN, uN˙, and uN¨ are the nodal displacement, velocity, and acceleration matrices, respectively.Based on the principle of superposition, nodal variablesuN, uN˙, and uN¨ can be expressed as a sum of the distortion in the principal modes of the vessel:(3)uN=BP,uN˙=BP˙,uN¨=BP¨,where B denotes the principal mode matrix and P denotes the principal coordinate matrix.Substituting (3) into (2), multiplying by BT, the equation of motion expressed in terms of principal coordinates is obtained:(4)mP¨+cP˙+kP=f,where m=BTMB, c=BTCB, k=BTKB, and f=BTF. m, c, k, and f are called generalized mass, damping, stiffness, and external force matrices, respectively.


## 2.2. Fluid Domain

### 2.2.1. 3D Components of Wave Forces

By assuming that the fluid flow is nonviscid, incompressible, and irrotational, the potential flow theory can be employed to model the fluid in the problem of fluid-structure interaction. In the reference frameoxyz, the fluid velocity potential ∅ of the fluid-structure system in steady state can be written as [16](5)∅x,y,z,t=Re⁡φx,y,zeiωt,where φ is the spatial part of ∅, i is the imaginary unit, and ω is the wave encounter frequency. ω is defined as(6)ω=ω0+Uk0cos⁡β,where ω0 is the natural wave circular frequency and k0 is the wave number.φ can be further decomposed into three parts:(7)φx,y,z=φIx,y,z+φDx,y,z+∑m=1hφRmx,y,z,where φI denotes the incident wave potential, φD denotes the diffraction wave potential, and φRm denotes the radiation wave potential of the mth degree of freedom of the vessel. h is greater than 6 and indices m>6 represent the elastic degrees of freedom. φI, φD, and φR can be solved based on the boundary element method [16].The generalized diffraction forcesFD and radiation forces FR are calculated based on the wet contact surface of the fluid-vessel system in the static hydroelastic state as illustrated in Figure 2:(8)FD=-ρξ∬n⃑·bmiω-U∂∂xφDdS,FR=-ρ∑k=1h∬n⃑·bmiω-U∂∂xφRkpkdS,where ρ is the fluid density, ξ is the incident wave amplitude, n⃑ is the normal vector pointing into the fluid domain, bm is the mth principal mode of the vessel, and dS is the infinitesimal area element of the wet contact surface.The generalized Froude-Krylov forcesFI and hydrostatic restoring forces FHS are computed based on the wet contact surface of the fluid-vessel system in the transient hydroelastic state as demonstrated in Figure 2:(9)FI=-ρξ∬n⃑·bmiω-U∂∂xφIdS,FHS=-ρg∑k=1hpk∬n⃑·bmwkdS-∭ρsgwmdV,where wm is the vertical displacement of the mth degree of freedom and ρs is the mass density of the vessel.


### 2.2.2. Impact Wave Forces

Slamming loads acting on deformable large vessels with considerable bow flares are calculated based on the 2D generalized Wagner method [17]. Some selected cross sections are chosen to calculate the 2D slamming pressure [10]. Figure 3 displays the description of 2D water entry of cross sections. vt is the velocity of the section, bt and ht are the coordinates of the contact point at time t, and fys is the contour of the cross section. In the present study, vt is the combination of rigid and elastic body motions of the vessel. Note that the effects of pitch on the motion of those cross sections are confined to the vertical components because the ship is slender.Figure 3
Impact problem description of symmetric water entry.The 2D potential flow theory is used to model this fluid-structure problem and the slamming pressure acting on the wet contour of each cross section can be solved by the conformal mapping method. Figure4 demonstrates how conformal mapping is performed. The real free surface contour is linearized as the horizontal line starting from the contact point. The line segments AB and CD from the original real plane ysozs are mapped to intervals -∞,-1 and +1,+∞ in the complex plane ξoη, respectively. The wet contour BC is mapped to the interval -1,+1. The slamming pressure can be obtained as follows:(10)pSLAMξ,t=ρdstdt2pbξ,t+dstdtpaξ,t,where st is the vertical displacement of the cross section. The expressions of paξ,t and pbξ,t can be found in [17]; here we do not explain these two expressions for brevity.Figure 4
Conformal mapping ((a) physical plane and (b) corresponding complex plane).
(a) (b)The solutionpSLAMys,zs,t in the original ysozs plane can be obtained based on the inverse of the conformal mapping. Then the generalized slamming forces can be calculated as(11)FSL=∫n⃑·bmpSLAMdl,where n⃑ is the normal vector pointing out of the wet contour of the cross section and dl is the infinitesimal line segment of the wet contact line of the section.


### 2.2.3. Green Water Loads

In extreme sea conditions, after the bow impact events, the incident wave may overtop the bow and run up onto the deck. In this paper, green water loads on deck are incorporated and simulated based on the dam-breaking Model. This method models the green water phenomenon as the process of sudden collapse of a dam. Before the sudden collapse, the water stored in the dam is static. Further we assume that during the process the ground exerts no resistant forces to the fluid. We introduce a 2D bow-fixed reference framexGozG in which the origin is located at the tip of the bow. The simplified Saint-Venant equations [18] can be expressed as(12)∂HxG,t∂t+HxG,t∂SxG,t∂xG+SxG,t∂HxG,t∂xG=0,1g∂SxG,t∂t+SxG,tg∂SxG,t∂xG+∂HxG,t∂xG=0,where H and S are the depth and velocity of the green water flow, respectively. g is the gravitational acceleration. Based on the method of characteristics, the depth H and velocity S can be solved as follows:(13)HxG,t=-xG3gt+2H032,(14)SxG,t=23xGt+gH0,where H0 is the initial depth of the water in the intact dam. Velocity S is a linear function of xG, and green water pressure on deck pGW can be easily calculated from S. We assume that pGW does not change in yG direction; then the generalized green water loads can be calculated as(15)FGW=∬n⃑·bmpGWdS.


## 2.2.1. 3D Components of Wave Forces

By assuming that the fluid flow is nonviscid, incompressible, and irrotational, the potential flow theory can be employed to model the fluid in the problem of fluid-structure interaction. In the reference frameoxyz, the fluid velocity potential ∅ of the fluid-structure system in steady state can be written as [16](5)∅x,y,z,t=Re⁡φx,y,zeiωt,where φ is the spatial part of ∅, i is the imaginary unit, and ω is the wave encounter frequency. ω is defined as(6)ω=ω0+Uk0cos⁡β,where ω0 is the natural wave circular frequency and k0 is the wave number.φ can be further decomposed into three parts:(7)φx,y,z=φIx,y,z+φDx,y,z+∑m=1hφRmx,y,z,where φI denotes the incident wave potential, φD denotes the diffraction wave potential, and φRm denotes the radiation wave potential of the mth degree of freedom of the vessel. h is greater than 6 and indices m>6 represent the elastic degrees of freedom. φI, φD, and φR can be solved based on the boundary element method [16].The generalized diffraction forcesFD and radiation forces FR are calculated based on the wet contact surface of the fluid-vessel system in the static hydroelastic state as illustrated in Figure 2:(8)FD=-ρξ∬n⃑·bmiω-U∂∂xφDdS,FR=-ρ∑k=1h∬n⃑·bmiω-U∂∂xφRkpkdS,where ρ is the fluid density, ξ is the incident wave amplitude, n⃑ is the normal vector pointing into the fluid domain, bm is the mth principal mode of the vessel, and dS is the infinitesimal area element of the wet contact surface.The generalized Froude-Krylov forcesFI and hydrostatic restoring forces FHS are computed based on the wet contact surface of the fluid-vessel system in the transient hydroelastic state as demonstrated in Figure 2:(9)FI=-ρξ∬n⃑·bmiω-U∂∂xφIdS,FHS=-ρg∑k=1hpk∬n⃑·bmwkdS-∭ρsgwmdV,where wm is the vertical displacement of the mth degree of freedom and ρs is the mass density of the vessel.


## 2.2.2. Impact Wave Forces

Slamming loads acting on deformable large vessels with considerable bow flares are calculated based on the 2D generalized Wagner method [17]. Some selected cross sections are chosen to calculate the 2D slamming pressure [10]. Figure 3 displays the description of 2D water entry of cross sections. vt is the velocity of the section, bt and ht are the coordinates of the contact point at time t, and fys is the contour of the cross section. In the present study, vt is the combination of rigid and elastic body motions of the vessel. Note that the effects of pitch on the motion of those cross sections are confined to the vertical components because the ship is slender.Figure 3
Impact problem description of symmetric water entry.The 2D potential flow theory is used to model this fluid-structure problem and the slamming pressure acting on the wet contour of each cross section can be solved by the conformal mapping method. Figure4 demonstrates how conformal mapping is performed. The real free surface contour is linearized as the horizontal line starting from the contact point. The line segments AB and CD from the original real plane ysozs are mapped to intervals -∞,-1 and +1,+∞ in the complex plane ξoη, respectively. The wet contour BC is mapped to the interval -1,+1. The slamming pressure can be obtained as follows:(10)pSLAMξ,t=ρdstdt2pbξ,t+dstdtpaξ,t,where st is the vertical displacement of the cross section. The expressions of paξ,t and pbξ,t can be found in [17]; here we do not explain these two expressions for brevity.Figure 4
Conformal mapping ((a) physical plane and (b) corresponding complex plane).
(a) (b)The solutionpSLAMys,zs,t in the original ysozs plane can be obtained based on the inverse of the conformal mapping. Then the generalized slamming forces can be calculated as(11)FSL=∫n⃑·bmpSLAMdl,where n⃑ is the normal vector pointing out of the wet contour of the cross section and dl is the infinitesimal line segment of the wet contact line of the section.


## 2.2.3. Green Water Loads

In extreme sea conditions, after the bow impact events, the incident wave may overtop the bow and run up onto the deck. In this paper, green water loads on deck are incorporated and simulated based on the dam-breaking Model. This method models the green water phenomenon as the process of sudden collapse of a dam. Before the sudden collapse, the water stored in the dam is static. Further we assume that during the process the ground exerts no resistant forces to the fluid. We introduce a 2D bow-fixed reference framexGozG in which the origin is located at the tip of the bow. The simplified Saint-Venant equations [18] can be expressed as(12)∂HxG,t∂t+HxG,t∂SxG,t∂xG+SxG,t∂HxG,t∂xG=0,1g∂SxG,t∂t+SxG,tg∂SxG,t∂xG+∂HxG,t∂xG=0,where H and S are the depth and velocity of the green water flow, respectively. g is the gravitational acceleration. Based on the method of characteristics, the depth H and velocity S can be solved as follows:(13)HxG,t=-xG3gt+2H032,(14)SxG,t=23xGt+gH0,where H0 is the initial depth of the water in the intact dam. Velocity S is a linear function of xG, and green water pressure on deck pGW can be easily calculated from S. We assume that pGW does not change in yG direction; then the generalized green water loads can be calculated as(15)FGW=∬n⃑·bmpGWdS.


## 2.3. Hydroelastic Vibration Domain

Based on the analyses of the structural and fluid domain, now the hydroelastic equation of motion can be rewritten as follows:(16)mP¨+cP˙+kP=fI+fD+fR+fHS+fSLAM+fGW,where the right-hand side (RHS) of (4) is decomposed into the six components shown in the RHS of (16). fI, fD, and fR are the generalized incident, diffraction, and radiation wave force matrices, respectively, fHS is the generalized hydrostatic restoring force matrix, and fSLAM and fGW are the generalized slamming and green water force matrices, respectively.For large elastic vessels with pronounced bow flare in severe sea conditions, it is of great importance to appropriately incorporate some nonlinearities in those generalized hydrodynamic force matrices of (16). This may not be purely mathematically consistent because (16) is basically linear. However, this treatment does allow evident improvements to the simulation of fierce fluid-structure interaction. As formulated in Section 2.2, among the six generalized hydrodynamic force matrices, only fD and fR are calculated based on the unchanged static hydroelastic state of the vessel. This is because of the difficulty of recalculating φD and φRn using 3D BEM based on the updated wet configuration of the vessel. In contrast, the other four generalized forces are calculated based on the actual instantaneous wet contact surface because the computational burden is acceptable. In other words, except for fD and fR, all the other generalized external forces are nonlinear.


## 2.4. Necessary Numerical Details

### 2.4.1. Mesh Generation

Mesh generation is the prerequisite of the computation of the generalized external forces. The hull mesh for the initial state of the vessel is generated by employing the method of cubic spline function [19]. Then the transient elastic hull mesh can be generated according to the superposition of the distortion in the principal modes. Figure 5 displays the initial hull mesh and the transient elastic hull mesh.Figure 5
Mesh generation ((a) initial hull mesh and (b) transient elastic hull mesh).
(a) (b)


### 2.4.2. NonlinearfI and fHS

The nonlinear generalized force matricesfI and fHS are calculated based on the instantaneous wet contact surface and therefore the acquisition of the exact wet hull plays a significant role. The method used in the present study is called the coarse-fine mesh strategy (see Figure 6), which is illustrated as follows: the initial hull mesh is classified as the coarse mesh. The coarse mesh panels are checked at every time step for their position relative to the incident wave profile. Those panels whose vertices are all below (or above) the incident wave contour are considered as the ordinary fully wet (or dry) panels. Special care is given for those partially wet ones. Each partial wet panel is further subdivided into smaller subpanels which are called the fine mesh. Then the reevaluation of the relative position of the fine mesh with respect to the incident wave profile is carried out. The incident wave profile within every fine mesh panel is approximated as the straight line. Based on this strategy, fI and fHS are calculated by the formulations given in Section 2.2.Figure 6
Coarse-fine mesh strategy for computingfI and fHS.


### 2.4.3. NonlinearfSLAM

The extraction of the cross sections of the vessel can be carried out based on the information of the hull geometry as shown in Figure7. For bow and stern areas, more cross sections should be extracted not only due to the fact that the slamming phenomena are more notable in these two areas but also for the compensation of the evident 3D effects associated with these regions.Figure 7
Extracted 2D sections along the hull.As mentioned above, the generalized Wagner method assumes that there is no flow separation in the slamming phenomena. Therefore modifications should be made for those sections whose breadth of the outer contour does not increase monotonically with respect to the height of the contour. Figure8 takes one of the cross sections of the bulb bow as an example to demonstrate those modifications.Figure 8
Modification of original outer contour of the bulb bow.For the cross sections where the outer contour consists of two disjoint curves as shown in Figure8, we use the smooth curve called the connector to connect the two separated modified curves. The wet contour of the connector is excluded in the computation of fSLAM.


### 2.4.4. NonlinearfGW

As formulated in Section2.2.3, the green water phenomenon is assumed to be one-dimensional. Figure 9 shows green water depth distribution at time t along the deck. The mesh for the deck is generated as shown in Figure 5. H0 in (14) is determined by the relative position of the deck and the incident wave profile at time t. Then fGW is computed on the wet deck mesh based on (15).Figure 9
Green water depth distribution along the deck.


### 2.4.5. Hydroelastic Response

In the present study, we adopt the well tested explicit fourth-order Runge-Kutta algorithm [20] to solve (16). Then the hydroelastic response can be acquired based on the principle of modal superposition:(17)wx,t=∑k=1hpktwkx,Mx,t=∑k=1hpktMkx,Vx,t=∑k=1hpktVkx,where w, M, and V are the displacement, vertical bending moment, and shear force; wk, Mk, and Vk are the kth mode of displacement, bending moment, and shear force, respectively.


## 2.4.1. Mesh Generation

Mesh generation is the prerequisite of the computation of the generalized external forces. The hull mesh for the initial state of the vessel is generated by employing the method of cubic spline function [19]. Then the transient elastic hull mesh can be generated according to the superposition of the distortion in the principal modes. Figure 5 displays the initial hull mesh and the transient elastic hull mesh.Figure 5
Mesh generation ((a) initial hull mesh and (b) transient elastic hull mesh).
(a) (b)


## 2.4.2. NonlinearfI and fHS

The nonlinear generalized force matricesfI and fHS are calculated based on the instantaneous wet contact surface and therefore the acquisition of the exact wet hull plays a significant role. The method used in the present study is called the coarse-fine mesh strategy (see Figure 6), which is illustrated as follows: the initial hull mesh is classified as the coarse mesh. The coarse mesh panels are checked at every time step for their position relative to the incident wave profile. Those panels whose vertices are all below (or above) the incident wave contour are considered as the ordinary fully wet (or dry) panels. Special care is given for those partially wet ones. Each partial wet panel is further subdivided into smaller subpanels which are called the fine mesh. Then the reevaluation of the relative position of the fine mesh with respect to the incident wave profile is carried out. The incident wave profile within every fine mesh panel is approximated as the straight line. Based on this strategy, fI and fHS are calculated by the formulations given in Section 2.2.Figure 6
Coarse-fine mesh strategy for computingfI and fHS.


## 2.4.3. NonlinearfSLAM

The extraction of the cross sections of the vessel can be carried out based on the information of the hull geometry as shown in Figure7. For bow and stern areas, more cross sections should be extracted not only due to the fact that the slamming phenomena are more notable in these two areas but also for the compensation of the evident 3D effects associated with these regions.Figure 7
Extracted 2D sections along the hull.As mentioned above, the generalized Wagner method assumes that there is no flow separation in the slamming phenomena. Therefore modifications should be made for those sections whose breadth of the outer contour does not increase monotonically with respect to the height of the contour. Figure8 takes one of the cross sections of the bulb bow as an example to demonstrate those modifications.Figure 8
Modification of original outer contour of the bulb bow.For the cross sections where the outer contour consists of two disjoint curves as shown in Figure8, we use the smooth curve called the connector to connect the two separated modified curves. The wet contour of the connector is excluded in the computation of fSLAM.


## 2.4.4. NonlinearfGW

As formulated in Section2.2.3, the green water phenomenon is assumed to be one-dimensional. Figure 9 shows green water depth distribution at time t along the deck. The mesh for the deck is generated as shown in Figure 5. H0 in (14) is determined by the relative position of the deck and the incident wave profile at time t. Then fGW is computed on the wet deck mesh based on (15).Figure 9
Green water depth distribution along the deck.


## 2.4.5. Hydroelastic Response

In the present study, we adopt the well tested explicit fourth-order Runge-Kutta algorithm [20] to solve (16). Then the hydroelastic response can be acquired based on the principle of modal superposition:(17)wx,t=∑k=1hpktwkx,Mx,t=∑k=1hpktMkx,Vx,t=∑k=1hpktVkx,where w, M, and V are the displacement, vertical bending moment, and shear force; wk, Mk, and Vk are the kth mode of displacement, bending moment, and shear force, respectively.


## 3. Experimental Setup

To evaluate the capabilities of the methodology presented in this study, a model test of a large vessel with a pronounced bow flare was carried out in the ship model towing tank of Harbin Engineering University. The model was provided by a commercial company (China State Shipbuilding Corporation). The towing tank, as shown in Figure10, is 108 m in length, 7 m in width, and 3.5 m in depth, with the wave absorption beach located at the front end and the electrically driven flap-type wave generator at the back end. The wave generator enables the generation of long-crest regular waves, the periods ranging from 0.4 s to 4 s, and irregular waves with the maximum significant wave height being 0.32 m. A wave probe is positioned close to the wave maker. A motion recorder apparatus, which is both installed in the model and attached to the towing carriage, is designed to obtain the time history of heave and pitch and fix the yaw angle. The maximum velocity of the towing carriage is 6 m/s. The schematic of the experimental setup is shown in Figure 11.Figure 10
Ship model towing tank.
(a) (b)Figure 11
A schematic of the experimental setup.The elastically scaled model consists of seven segmented hulls made from epoxy, three rectangular nonuniform steel beams connecting those segments, a self-propelled system positioned in the stern, and the wired monitoring system used to acquire experimental data. Table1 gives the principal particulars of the elastic model. The propelling system was designed to maintain the expected forward speed of the model in still water. Besides the aforementioned motion recorder, the tethered monitoring units in the model comprise (i) strain gauges glued on the surface of the beams and (ii) impact pressure sensors positioned in the bow area. The main particulars of those units are demonstrated in Table 2. Figure 12 presents necessary details of the arrangement of the beams and pressure sensors. Some complementary and schematic information of experimental setup is displayed in Figure 13.Table 1
Principal particulars of the vessel model.
Item Value Scale 64 Overall length 4.85 m Waterline length 4.30 m Breadth 0.51 m Designed draft 0.12 m Displacement 195 Kg Longitudinal Center of Gravity 1.94 mTable 2
Tethered sensor system instrumentation details.
Measurement Instrumentation Sample rate (Hz) Impact pressure Pressure sensor 5000 Vertical bending moment Strain gauge 200Figure 12
General arrangement of the ship model ((a) beams arrangement and (b) pressure sensors arrangement).
(a) (b)Figure 13
Some details of experimental setup.


## 4. Comparison between Numerical and Experimental Results

The numerical model was established based on the experimental setup. For the nonuniform beam model, 20 elements were utilized and only the first three elastic modes of vertical bending were considered to model the elastic characteristics. The whole ship surface was divided into 2074 panels, 56 of which were the deck mesh for the computation offGW. One coarse panel was subdivided into 32 subpanels when necessary for the computation of fI and fHS. 22 unevenly spaced cross sections were chosen to calculate fSLAM. A series of hydroelastic experiments in regular waves were conducted. The present work focuses on those large wave amplitude conditions. Two representative experimental conditions were picked out for the following comparison and analysis (see Table 3). It should be mentioned here that, for the preservation of commercial benefits of the company that provided the vessel prototype, some parts of the pictorial details in the following experimental photographs were pixelated.Table 3
Two representative experimental conditions.
Vessel speed (kn) Wave height (m) Wave period (s) 5 11.6 13.28 5 16 13.28


### 4.1. Natural Vibrational Characteristics

The elastic modes were calculated based on an improved Euler-Bernoulli beam model. In this treatment, the shell of the vessel was incorporated and modeled as rigid elements assigned with mass. Then the effects of mass distribution and subsequent moment of inertia can be considered in the computation of the natural vibrational characteristics. In order to validate the result, a hammering test, which involves applying an arbitrary impact to the ship model in still water, was conducted before the model tests in waves. The hammering tests were repeated three times at three locations of the model and the averaged wet natural frequencies were compared with those acquired numerically (see Table4). A satisfactory agreement appeared in the results of the first order, with increasingly larger discrepancies in higher orders. This was acceptable because the elastic backbones were designed based on the first-order vibrational characteristics. Moreover, the hydroelastic response components associated with the second and third orders were fairly small compared with those corresponding to the first order. The first three elastic natural modes are shown in Figure 14.Table 4
Comparison of wet natural frequencies at model scale.
Order of mode Numerical results (Hz) Experimental results (Hz) Relative error (%) First 9.709 9.708 −0.01 Second 22.418 24.955 11.31 Third 39.696 48.304 21.68Figure 14
First three elastic beam modes at full scale ((a) displacement and (b) vertical bending moments).
(a) (b)


### 4.2. Wave Calibration

Since the parameters of experimentally generated regular waves may slightly deviate from the designed ones, in the present study, a wave probe was installed near the wave maker to measure the uncertainties associated with the parameters of the generated waves. For the following comparison and analysis, numerical calculation will be based on the measured data from the wave probe to reduce as many preventable errors as possible.Figure15 shows the time history of the experimentally generated waves whose wave height and frequency are 181.25 mm and 1.66 s. The whole process can be divided into two regions: (i) region I describes the duration before the wave front reaches the wave probe and records the first five periods of the waves with ever-increasing wave heights; (ii) region II illustrates the relatively smooth periods of the waves. Among all the periods, 17 periods from 40 s to 70 s as displayed in Figures 15 and 16 are picked out, since those periods correspond to the effective experimental results. The fluctuations of the wave parameters are depicted in Figure 16, thus suggesting strongly the necessity of modifications of the corresponding numerical inputs in the following computation of slamming pressure and bending moments. In this paper, we replace the original wave parameters with the averaged value obtained from the data provided by the wave probe to perform numerical calculation. For brevity, the pictorial details of wave parameters of the second experimental condition are not presented here. The modified wave parameters of the two representative experimental conditions are presented in Table 5.Table 5
Comparison and modification of wave characteristics at model scale.
Wave height in theory (mm) Modified wave height (mm) Wave period in theory (s) Modified wave period (s) Modified encountering wave frequency (Hz) 181.25 189.59 1.66 1.70 0.66 250.00 263.85 1.66 1.70 0.66Figure 15
Time history of the experimentally generated regular waves.Figure 16
Fluctuations of experimental wave parameters at model scale ((a) wave period and (b) wave height).
(a) (b)


### 4.3. Slamming Loads in Bow Area

As shown in Figure12, 12 pressure sensors were positioned on the bow area to measure the impact pressure. Pressure sensor number 2 was picked out as an example to compare experimental and numerical results. Figures 17 and 18, as well as all the other experimental results in this study, reveal the inevitability of the presence of noise in the experimental pressure recordings. This is due to the effects of the surrounding metallic structures and other electromagnetic sources on the pressure sensor recording system. Furthermore, it can be deduced from Figures 17 and 18 that, for those experimental results in which the duration of the impact is comparable to the period of the noise, the common method of filtering the original experimental data to obtain the realistic impact peak value fails to be applicable. Therefore, the exact realistic peak impact value should be replaced by an estimated one, with the error being the peak noise value that is extracted by filtering operation (see the enlarged images in Figures 17 and 18). Table 6 gives the impact peak value intervals. Although the exact peak value cannot be obtained, other information, such as the local structure vibration near the sensors, can be deduced from the pressure oscillations through filtering operation as also shown in those enlarged images.Table 6
Impact peak value ranges.
Operational conditions 1st peak value range (MPa) 2nd peak value range (MPa) 3rd peak value range (MPa) 4th peak value range (MPa) 5th peak value range (MPa) 5 kn, 11.6 m 0.190 ± 0.004 0.244 ± 0.004 0.193 ± 0.005 0.265 ± 0.004 0.221 ± 0.004 5 kn, 16 m 0.646 ± 0.033 0.733 ± 0.034 0.582 ± 0.035 0.877 ± 0.032 0.638 ± 0.033Figure 17
Filtering operation of experimental impact pressure ((a) original and filtered result and (b) enlarged image; corresponding to 5 kn 11.6 m).Figure 18
Filtering operation of experimental impact pressure ((a) original and filtered result and (b) enlarged image; corresponding to 5 kn 16 m).
(a) (b)For both operational conditions (see Figure20), the experimentally obtained impact peak values and periods vary from each other, resulting from the fluctuations of the parameters of experimentally generated waves and the variations of the speed of the test model. On average, the numerically obtained impact peak values are larger than those acquired experimentally. In Figure 20(b), the experimental result shows that, for the severe sea conditions, as exemplified in Figure 19, the real physical characteristics of impact events are actually rather complicated that the impact peak is followed by several smaller peaks and troughs. The current generalized Wagner method cannot fully capture these details. Nevertheless, Figure 20 in general shows reasonable agreement between experimental and numerical results.Figure 19
Emergence of the vessel above the water wave (a) and subsequent severe impact phenomenon (b).
(a) (b)Figure 20
Comparison of impact pressure ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)


### 4.4. Hydroelastic Vibrational Vertical Bending Moments

Experimental and numerical nonlinear structural responses, extracted amidships for the two aforementioned operational conditions, were compared and analyzed in terms of the total vertical bending moments (VBM), low frequency VBM, and high frequency impact VBM (see Figures21–23).Figure 21
Comparison of total VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 22
Comparison of low frequency VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 23
Comparison of impact VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)In general, there are satisfactory agreements between the experimentally and numerically determined VBM for these two operational conditions. In Figure21(a), on average, the experimental results of the peaks slightly outweigh the numerical results, while the opposite phenomenon can be seen for the results of troughs. In contrast both the numerically obtained peaks and troughs outnumber experimental results in Figure 21(b). Figures 21(a) and 21(b) share a common feature that, for both experimental and numerical results, the values of the troughs are bigger than the values of the peaks which can be attributed to the characteristics of impact VBM demonstrated in Figure 23. The low frequency VBM in Figures 22(a) and 22(b) also shows similar asymmetry between peaks and troughs which can be mostly attributed to the nonlinearities associated with the vertical asymmetrical characteristics of the geometry of the hull. In Figure 22(a) the numerical method slightly underestimates peaks of the low frequency VBM, with the observation of relatively larger underestimation of troughs in Figure 22(b). From Figure 23, it can be seen that the numerical method underestimates the peak values while it overestimates the trough values compared with the experimental results. In addition, the impact VBM clearly manifests transient characteristics in these two sea conditions.


### 4.5. Analysis on Effects of Nonlinearity

Linear numerical methods are applicable for the prediction of structural response of vessels in low-to-moderate seas. However, for large vessels in high seas, nonlinearities cannot be neglected. The purpose of this part is to determine the nonlinear effects on the numerical simulation.In the linear case, the 2D impact pressure of (11) is integrated based on the unchanged wet contact line of the cross section of the vessel in the static hydroelastic state. The sectional vertical displacement st of (10) is also calculated linearly. Figure 24 shows that this linear treatment largely underestimates the impact pressure compared with the proposed nonlinear method. In addition, after the impact peak, the following peak caused by the structural vibration feedback is also underestimated by the linear method.Figure 24
Comparison of impact pressure ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Linear and nonlinear results of structural responses are compared (see Figures25–27). In the linear model, all the generalized hydrodynamic forces are calculated based on the unchanged wet contact surface. As shown in these three figures, the linear results are generally smaller than the nonlinear and experimental results. The most evident difference between the linear case and the nonlinear case is that the peaks and troughs determined linearly are nearly the same. This is because the linear model cannot capture the change of the geometry of the wet contact surface caused by the motion of the vessel.Figure 25
Comparison of total VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 26
Comparison of low frequency VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 27
Comparison of impact VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)


### 4.6. Analysis on Effects of Elasticity

A conventional approach takes two steps to calculate whipping responses. The first step assumes that the vessel is a rigid body and the motion and VBM are computed based on this assumption. The second step considers the vessel as an elastic body and the results of motion in the first step serve as the input for the computation of impact VBM. Then the impact VBM is superimposed on the VBM in the first step for the total VBM. This part compares the results determined by this conventional method and the proposed nonlinear method to investigate the effects of elasticity on the impact pressure and structural responses. The comparisons are demonstrated from Figures28–31. The results computed by the conventional model are tagged “rigid result” in Figures 28, 29, 30, and 31.Figure 28
Comparison of impact pressure ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 29
Comparison of total VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 30
Comparison of low frequency VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 31
Comparison of impact VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure28 shows that the impact pressure calculated based on the rigid body outweighs the elastic body based results. Also the rigid body impact pressure decreases smoothly and monotonically after the impact peak. This is because the motion input for the computation of impact pressure is totally determined by the rigid body model and no structural vibration feedback is available. From Figures 29–31, it can be found that, in general, the rigid body method overestimates the structural responses. This is also due to the absence of elastic vibration feedback in the rigid model because the elastic vibration determined by elastic modes will mitigate the uneven external hydrodynamic pressure distribution and the mitigated external loads will reduce the structural responses.


## 4.1. Natural Vibrational Characteristics

The elastic modes were calculated based on an improved Euler-Bernoulli beam model. In this treatment, the shell of the vessel was incorporated and modeled as rigid elements assigned with mass. Then the effects of mass distribution and subsequent moment of inertia can be considered in the computation of the natural vibrational characteristics. In order to validate the result, a hammering test, which involves applying an arbitrary impact to the ship model in still water, was conducted before the model tests in waves. The hammering tests were repeated three times at three locations of the model and the averaged wet natural frequencies were compared with those acquired numerically (see Table4). A satisfactory agreement appeared in the results of the first order, with increasingly larger discrepancies in higher orders. This was acceptable because the elastic backbones were designed based on the first-order vibrational characteristics. Moreover, the hydroelastic response components associated with the second and third orders were fairly small compared with those corresponding to the first order. The first three elastic natural modes are shown in Figure 14.Table 4
Comparison of wet natural frequencies at model scale.
Order of mode Numerical results (Hz) Experimental results (Hz) Relative error (%) First 9.709 9.708 −0.01 Second 22.418 24.955 11.31 Third 39.696 48.304 21.68Figure 14
First three elastic beam modes at full scale ((a) displacement and (b) vertical bending moments).
(a) (b)


## 4.2. Wave Calibration

Since the parameters of experimentally generated regular waves may slightly deviate from the designed ones, in the present study, a wave probe was installed near the wave maker to measure the uncertainties associated with the parameters of the generated waves. For the following comparison and analysis, numerical calculation will be based on the measured data from the wave probe to reduce as many preventable errors as possible.Figure15 shows the time history of the experimentally generated waves whose wave height and frequency are 181.25 mm and 1.66 s. The whole process can be divided into two regions: (i) region I describes the duration before the wave front reaches the wave probe and records the first five periods of the waves with ever-increasing wave heights; (ii) region II illustrates the relatively smooth periods of the waves. Among all the periods, 17 periods from 40 s to 70 s as displayed in Figures 15 and 16 are picked out, since those periods correspond to the effective experimental results. The fluctuations of the wave parameters are depicted in Figure 16, thus suggesting strongly the necessity of modifications of the corresponding numerical inputs in the following computation of slamming pressure and bending moments. In this paper, we replace the original wave parameters with the averaged value obtained from the data provided by the wave probe to perform numerical calculation. For brevity, the pictorial details of wave parameters of the second experimental condition are not presented here. The modified wave parameters of the two representative experimental conditions are presented in Table 5.Table 5
Comparison and modification of wave characteristics at model scale.
Wave height in theory (mm) Modified wave height (mm) Wave period in theory (s) Modified wave period (s) Modified encountering wave frequency (Hz) 181.25 189.59 1.66 1.70 0.66 250.00 263.85 1.66 1.70 0.66Figure 15
Time history of the experimentally generated regular waves.Figure 16
Fluctuations of experimental wave parameters at model scale ((a) wave period and (b) wave height).
(a) (b)


## 4.3. Slamming Loads in Bow Area

As shown in Figure12, 12 pressure sensors were positioned on the bow area to measure the impact pressure. Pressure sensor number 2 was picked out as an example to compare experimental and numerical results. Figures 17 and 18, as well as all the other experimental results in this study, reveal the inevitability of the presence of noise in the experimental pressure recordings. This is due to the effects of the surrounding metallic structures and other electromagnetic sources on the pressure sensor recording system. Furthermore, it can be deduced from Figures 17 and 18 that, for those experimental results in which the duration of the impact is comparable to the period of the noise, the common method of filtering the original experimental data to obtain the realistic impact peak value fails to be applicable. Therefore, the exact realistic peak impact value should be replaced by an estimated one, with the error being the peak noise value that is extracted by filtering operation (see the enlarged images in Figures 17 and 18). Table 6 gives the impact peak value intervals. Although the exact peak value cannot be obtained, other information, such as the local structure vibration near the sensors, can be deduced from the pressure oscillations through filtering operation as also shown in those enlarged images.Table 6
Impact peak value ranges.
Operational conditions 1st peak value range (MPa) 2nd peak value range (MPa) 3rd peak value range (MPa) 4th peak value range (MPa) 5th peak value range (MPa) 5 kn, 11.6 m 0.190 ± 0.004 0.244 ± 0.004 0.193 ± 0.005 0.265 ± 0.004 0.221 ± 0.004 5 kn, 16 m 0.646 ± 0.033 0.733 ± 0.034 0.582 ± 0.035 0.877 ± 0.032 0.638 ± 0.033Figure 17
Filtering operation of experimental impact pressure ((a) original and filtered result and (b) enlarged image; corresponding to 5 kn 11.6 m).Figure 18
Filtering operation of experimental impact pressure ((a) original and filtered result and (b) enlarged image; corresponding to 5 kn 16 m).
(a) (b)For both operational conditions (see Figure20), the experimentally obtained impact peak values and periods vary from each other, resulting from the fluctuations of the parameters of experimentally generated waves and the variations of the speed of the test model. On average, the numerically obtained impact peak values are larger than those acquired experimentally. In Figure 20(b), the experimental result shows that, for the severe sea conditions, as exemplified in Figure 19, the real physical characteristics of impact events are actually rather complicated that the impact peak is followed by several smaller peaks and troughs. The current generalized Wagner method cannot fully capture these details. Nevertheless, Figure 20 in general shows reasonable agreement between experimental and numerical results.Figure 19
Emergence of the vessel above the water wave (a) and subsequent severe impact phenomenon (b).
(a) (b)Figure 20
Comparison of impact pressure ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)


## 4.4. Hydroelastic Vibrational Vertical Bending Moments

Experimental and numerical nonlinear structural responses, extracted amidships for the two aforementioned operational conditions, were compared and analyzed in terms of the total vertical bending moments (VBM), low frequency VBM, and high frequency impact VBM (see Figures21–23).Figure 21
Comparison of total VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 22
Comparison of low frequency VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 23
Comparison of impact VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)In general, there are satisfactory agreements between the experimentally and numerically determined VBM for these two operational conditions. In Figure21(a), on average, the experimental results of the peaks slightly outweigh the numerical results, while the opposite phenomenon can be seen for the results of troughs. In contrast both the numerically obtained peaks and troughs outnumber experimental results in Figure 21(b). Figures 21(a) and 21(b) share a common feature that, for both experimental and numerical results, the values of the troughs are bigger than the values of the peaks which can be attributed to the characteristics of impact VBM demonstrated in Figure 23. The low frequency VBM in Figures 22(a) and 22(b) also shows similar asymmetry between peaks and troughs which can be mostly attributed to the nonlinearities associated with the vertical asymmetrical characteristics of the geometry of the hull. In Figure 22(a) the numerical method slightly underestimates peaks of the low frequency VBM, with the observation of relatively larger underestimation of troughs in Figure 22(b). From Figure 23, it can be seen that the numerical method underestimates the peak values while it overestimates the trough values compared with the experimental results. In addition, the impact VBM clearly manifests transient characteristics in these two sea conditions.


## 4.5. Analysis on Effects of Nonlinearity

Linear numerical methods are applicable for the prediction of structural response of vessels in low-to-moderate seas. However, for large vessels in high seas, nonlinearities cannot be neglected. The purpose of this part is to determine the nonlinear effects on the numerical simulation.In the linear case, the 2D impact pressure of (11) is integrated based on the unchanged wet contact line of the cross section of the vessel in the static hydroelastic state. The sectional vertical displacement st of (10) is also calculated linearly. Figure 24 shows that this linear treatment largely underestimates the impact pressure compared with the proposed nonlinear method. In addition, after the impact peak, the following peak caused by the structural vibration feedback is also underestimated by the linear method.Figure 24
Comparison of impact pressure ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Linear and nonlinear results of structural responses are compared (see Figures25–27). In the linear model, all the generalized hydrodynamic forces are calculated based on the unchanged wet contact surface. As shown in these three figures, the linear results are generally smaller than the nonlinear and experimental results. The most evident difference between the linear case and the nonlinear case is that the peaks and troughs determined linearly are nearly the same. This is because the linear model cannot capture the change of the geometry of the wet contact surface caused by the motion of the vessel.Figure 25
Comparison of total VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 26
Comparison of low frequency VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 27
Comparison of impact VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)


## 4.6. Analysis on Effects of Elasticity

A conventional approach takes two steps to calculate whipping responses. The first step assumes that the vessel is a rigid body and the motion and VBM are computed based on this assumption. The second step considers the vessel as an elastic body and the results of motion in the first step serve as the input for the computation of impact VBM. Then the impact VBM is superimposed on the VBM in the first step for the total VBM. This part compares the results determined by this conventional method and the proposed nonlinear method to investigate the effects of elasticity on the impact pressure and structural responses. The comparisons are demonstrated from Figures28–31. The results computed by the conventional model are tagged “rigid result” in Figures 28, 29, 30, and 31.Figure 28
Comparison of impact pressure ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 29
Comparison of total VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 30
Comparison of low frequency VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure 31
Comparison of impact VBM ((a) 5 kn 11.6 m and (b) 5 kn 16 m).
(a) (b)Figure28 shows that the impact pressure calculated based on the rigid body outweighs the elastic body based results. Also the rigid body impact pressure decreases smoothly and monotonically after the impact peak. This is because the motion input for the computation of impact pressure is totally determined by the rigid body model and no structural vibration feedback is available. From Figures 29–31, it can be found that, in general, the rigid body method overestimates the structural responses. This is also due to the absence of elastic vibration feedback in the rigid model because the elastic vibration determined by elastic modes will mitigate the uneven external hydrodynamic pressure distribution and the mitigated external loads will reduce the structural responses.


## 5. Conclusions

In the present research, impact pressure and vibrational responses of a large vessel with a pronounced bow flare were studied by a nonlinear hydroelastic method, and experimental tests were conducted to validate the proposed numerical model. Main findings are concluded as follows:(1)
The fluctuations in the parameters of the experimentally generated regular waves seemed inevitable. Therefore, for the following comparison between experimental and numerical results, the corresponding numerical wave parameters needed to be modified based on the measuring data of the wave probe.(2)
The existence of noise in the experimental impact pressure makes the filtering operation necessary. In some cases where the period of noise was very close to the duration of the impact, the exact peak impact pressures were difficult to obtain and should be treated as estimated values.(3)
The impact pressures in the two representative operational conditions were overestimated by the fully coupled nonlinear method. As for the vibrational responses, the numerical and experimental results were in reasonable agreement.(4)
The linear method underestimated the impact pressure and structural responses compared with the proposed nonlinear method and experimental model and produced similar peaks and troughs in VBM. This indicates that nonlinear effects should be considered in severe sea conditions.(5)
The conventional rigid body method overestimated the impact pressure and the structural responses in comparison with the results determined by the proposed hydroelastic model and experiments. This is because of the inability of the rigid model to consider the elastic vibration feedback of the vessel. Therefore the hydroelastic model is much more suitable than the rigid body model in terms of computing structural response for large vessels sailing in high seas.



---

*Source: 7054768-2018-02-27.xml*