Abstract

Wave loads estimation and structural strength evaluation are the fundamental work at the ship design stage. The hydroelastic responses and slamming strength issues are also concerned especially for large-scale high-speed ships sailing in harsh waves. To accurately predict the wave-induced motions and loads acting on the ship sailing in regular waves, a fully coupled 3D time-domain nonlinear hydroelasticity theory is developed in this paper. The vibration modal characteristics of the flexible hull structure derived by the 3D finite element method (FEM) and simplified 1D nonuniform Timoshenko beam theory are firstly described. The hydrostatic restoring force and hydrodynamic wave force are calculated on the real-time wetted surface of hull to address geometric nonlinearity due to the steep wave and large amplitude motions. The bow slamming and green water loads acting on the ship in severe regular waves are estimated by the momentum impact method and dam-breaking method, respectively. Moreover, a small-scaled segmented ship model is designed, constructed, and tested in a laboratory wave basin to validate the hydroelasticity algorithm. The results predicted by theoretical and experimental approaches are systemically compared and analyzed. Finally, future work for predictions of ship hydroelasticity and slamming loads in irregular waves is prospected.

1. Introduction

The evaluation of the ship seakeeping performance and structural strength is a key element in the design, optimization, and operation guidance of seagoing ships. Fundamental to this is the predictions of the environment loads acting on ships when sailing in seaways [1]. To date, considerable efforts have been devoted to predicting ship motion and load responses in waves. The hydrodynamic responses of the ship at zero or low speed in small amplitude wave conditions are well addressed by the classical 2D and 3D linear frequency-domain potential flow theory [2]. In fact, with the development of ships towards large-scale, high-speed, and light-weighted, hydroelasticity theory should be employed to predict the load responses of large hull structures operating in waves [3]. In addition, as illustrated in Figure 1 (downloaded from the Internet [4, 5]), high-speed seagoing vessels with broad flare bow operating in harsh conditions are more susceptible to slamming impact. The slamming impact has pronounced influence on hull girder hydroelastic vibrations. Therefore, the slamming loads should be evaluated associated with the hull girder hydroelastic responses so as to accurately predict the wave-induced global motions, external loads, and structural responses of the large ship sailing in severe waves [6].

The hydroelasticity is a phenomenon concerned with the mutual interactions among inertial, hydrodynamic, and elastic forces. The concept of hydroelasticity was first proposed in the field of aerodynamics by Heller and Abramson [7]. Ship hydroelasticity theory was then widely developed since the 1970s [8]. The 2D ship hydroelasticity theory was first established on the basis of strip theory and linear beam structural dynamics to predict the symmetry (i.e., vertical bending) and antisymmetry (i.e., coupled horizontal bending and twisting) vibration responses of beamlike ship in waves [9, 10]. Wu and Price [11, 12] extended the hydroelasticity theory to 3D by establishing the boundary condition of flexible body with respect to fluid domain (the so-called Price-Wu boundary condition). To date, the 3D frequency-domain and time-domain linear or nonlinear hydroelasticity theories are widely adopted in the field of naval architecture and ocean engineering [13].

Transient impulsive loads, induced by bottom slamming, flare slamming, and green water on deck, should be considered during hydroelastic analysis especially for ships with pronounced flare bow sailing in harsh waves or at high speed. Severe slamming loads can result in not only local structure damages but also transient global whipping responses. Up to now, considerable efforts have been made on predicting the slamming loads acting on ships [14, 15]. A good review of slamming and whipping loads has been undertaken in the ISSC load committee [16]. Generally, the Wagner model [17] and momentum impact theory [18] are the two classical methods in use for the estimation of ship slamming loads. Moreover, green water loads can be predicted by many approaches, for example, simplified empirical approach, theoretical calculation, and Reynolds Averaged Navier–Stokes (RANS) simulation [19]. Generally, the dam-breaking model and flood wave model are the two classical methods in use for the estimation of green water loads together with the potential flow theory [20].

Although remarkable achievements have been made on developing the ship hydrodynamics and hydroelasticity theory, physical experiments constitute an indispensable tool in the field of naval architecture and ocean engineering. Experiments do not only reflect what really happens on the ship, but are also used to validate the numerical algorithm [21]. Small-scaled segmented model tank test is an effective tool in the measurement of wave-induced ship global motions, sectional loads and structural responses [22, 23]. Jiao et al. [24] developed a seven-segment model and a set of tubular structure steel backbones with varying cross section are used to connect the segmented hulls. The local slamming pressure can be also measured by mounting pressure sensors at bow or aft areas of the scaled segmented model [25]. Hong et al. [26] investigated the bow-flare slamming load behavior of a 10,000 TEU containership in regular waves. Lavroff et al. [27] studied the wave slamming loads on wave-piercing catamarans operating at a high speed by tank model test.

This paper aims at presenting a fully coupled 3D time-domain nonlinear hydroelasticity theory to estimate the motions and wave loads acting on ships sailing in severe regular waves. Small-scaled segmented model tests are also conducted to validate the numerical algorithm. The nonlinear behavior of ship motion and wave load responses under different conditions is comprehensively analyzed based on the numerical and experimental results. Moreover, our recent work regarding theoretical and experimental investigation on ship hydroelasticity and slamming loads in both long-crested and short-crested irregular waves is prospected.

2. Hull Structure Vibration Modal Analysis

The dry mode analysis of the flexible hull structure is the fundamental work prior to ship hydroelastic loads computation. Real hull girder is a highly nonlinear vibration system due to the complexity of the hull structure and external environment loads. The 3D finite element method (FEM) has good accuracy in the simulation of the vibration response of complex hull structure; it is, however, complex, expensive, and time-consuming [28]. On the other hand, for simplification, the 1D nonuniform Timoshenko beam system is usually used to simulate the vibration mode characteristic of hull girder. The mode analysis methods, including the 3D FEM and 1D beam theory, are concluded as follows.

2.1. Structure Vibration Response by 3D FEM

The hull structure is assumed to be homogeneous, isotropic, and linear elastic. According to the knowledge of 3D structural finite element (FE) theory, the dynamic governing equation of flexible hull structure is expressed as follows:where , , and are, respectively, structural mass, damping, and stiffness matrices; denotes nodal displacement matrix; and denotes external force matrix.

The vibration characteristic equation of hull structure is expressed as follows:

The rth (r = 1, 2, …, m) intrinsic frequency ωr and modal shape {Dr} of elastic structural system can be calculated by solving Equation (2). The rth order intrinsic modal shape is expressed as follows:where denotes the generalized displacement vector of the node with ID j in rth order vibration and is expressed as .

The overall modal shape of elastic structural system with first m order mode considered is expressed as

Based on the principle of mode superposition, the nodal displacement matrix is expressed aswhere denotes generalized modal principal coordinate matrix, denotes the principal coordinate vector corresponding to the rth order mode and can be decomposed into time-dependent term and space-dependent term using the following equation:

Substituting Equation (5) into Equation (1) and multiplying [D]T at left side of each term, the dynamic governing equation of flexible hull structure is rewritten as follows:where the generalized structural mass, damping, and stiffness matrices are expressed as [m] = [D]T[M][D], [c] = [D]T[C][D], and [k] = [D]T[K][D], respectively. The generalized external force matrix is expressed as [f] = [D]T[F].

2.2. Hull Girder Response by Timoshenko Beam Theory

The 1D nonuniform hull girder is dispersed into a number of beam elements, and each element can be regarded as a uniform Timoshenko beam. The vibration mode characteristics of the overall 1D hull girder can be solved by transfer matrix method (TMM). The acquiring of vertical vibration mode of hull girder by TMM theory is described below, which is also applicable for hull horizontal vibration mode analysis. The vertical vibration equation of a simplified uniform beam element, with the influence of shearing deformation and moment of inertia considered, is expressed as follows:where E denotes elastic modulus, G denotes shear modulus, I denotes sectional vertical bending moment of inertia, Iy denotes the moment of inertia of section with respect to y axis, AS denotes sectional shear area, μ denotes mass of beam element, ky denotes the sectional shape-related coefficient, denotes the vertical vibrational displacement, and the subscript i denotes the dispersed element number.

The vertical displacement of the element can be further decomposed into space-dependent term and time-dependent term:where denotes the amplitude of vertical vibration displacement of the element, ω denotes the natural frequency of the hull girder system, and ε denotes the phase angle.

By substituting Equation (9) into Equation (8), the following equation can be obtained:

The solution of Equation (10) can be generalized as follows:where A, B, C, and D are the undetermined constant coefficients, the expressions of λ1 and λ2 are given in Appendix.

According to beam theory, the sectional rotation angle θ(x), vertical bending moment (VBM) M(x), and vertical shearing force (VSF) V(x) of beam element can be obtained by taking the derivative of vertical displacement in Equation (11). Here, the state vectors at the left and right sides of the ith element are expressed by and , respectively. And the relationship between them can be expressed as

As observed from (12), it is clear that the matrix transfers the motion and force state information from the left side to the right side of the ith beam element. The detailed expression of is presented in the Appendix. By employing the continuity condition, the transfer functions can be induced as follows:where the 4 × 4 dimensional matrix Π denotes the overall transfer matrix of the hull girder system, and it is determined by multiplying of all the Fi (i = 1, 2, …, n).

Therefore, the relationship between the left side boundary condition (at stern side) and right side boundary condition (at bow side) of hull girder is written as follows:

Since the hull girder can be regarded as a beam with two free ends, the sectional VBM and VSF acting on the two ends are zero. Therefore, we obtain M0 = V0 = 0 and Mn = Vn = 0. Then the Equation (14) can be simplified as

Due to the fact that and θ0 are not always zero for a free vibration beam, the determinant of coefficients in Equation (15) should be zero. Therefore, the natural frequency of hull girder in different orders can be obtained by solving the following equation:

Once the natural frequency ω(k) (k = 1, 2, …) is obtained, the principal mode shape for vertical displacement , rotation angle θ(k)(x), VBM M(k)(x), and VSF V(k)(x) at different vibration order k can be obtained by substituting the natural frequency into the matrix Fi (i = 1, 2, …, n).

When the coupled effect between horizontal bending and twisting is not considered, the horizontal vibration mode of hull girder can be obtained in a same manner as vertical bending vibration. While for torsional deformation of hull girder, the relationship between the left side boundary condition (at stern side) and right side boundary condition (at bow side) of hull girder is written as follows:where φ denotes sectional twist angle and T denotes torsional moment.

According to the boundary condition of two free ends (T0 = T0 = 0), we know that Π21 = 0. Then the natural frequency and modal shape for torsional deformation can be solved. For the coupled horizontal and torsional vibration case, the mode analysis method by TMM can be found in Reference [29].

3. Time-Domain Nonlinear Hydroelasticity Theory

A 3D time-domain nonlinear hydroelasticity algorithm for ship motion and load predictions is developed in this section. In the hydroelasticity algorithm, the potential flow theory-based 3D boundary element method (BEM) and structural response-based 3D FEM (or simplified 1D Timoshenko beam theory) are combined to estimate the external loads acting on the elastic hull surface. The calculation of hydrostatic restoring force, wave excitation force, wave diffraction force, and radiation force are performed on the instantaneous wetted surface of elastic hull. The impact loads caused by bow slamming and green water on deck are included in the time-domain hydroelastic governing equation. The fourth order Runge–Kutta method is used to solve the time-domain nonlinear differential equation and the hull hydroelastic responses are obtained by the modal superposition principle.

3.1. Basic Assumption

In the potential flow theory, the fluid is assumed to be ideal, that is, incompressible, inviscid, continuous, and irrotational. The boundary conditions of fluid domain for the fluid-flexible structure interaction (FFSI) issue are illustrated in Figure 2, where Ω, SB, SF, SH, and S denote the fluid domain, body surface condition, free surface condition, bottom condition, and infinity condition, respectively.

The ship is advancing in waves at a constant speed of U, and the influence of propulsion system on the flow field is not considered. The incident waves are assumed to be linear regular waves. In this study, the head wave condition is defined as 0° heading angle, while 90°, 180°, and 270° are for starboard beam wave, following wave, and port beam wave, respectively. In order to describe ship motion in waves, a total of three right-hand coordinate systems (Figure 3) are adopted, and they are described as follows.

(1)Coordinate System 1: the space-fixed system O-XYZ is an inertia system with the plane OXY lying on the undisturbed free surface, OX directing with positive in the oncoming direction of incident waves and the OZ pointing vertically upward.(2)Coordinate System 2: the plane movement system o-xyz is a moving coordinate system with the same moving speed and direction as the ship forward speed U and wave heading β. The origin o is coincident with the origin point O initially. The plane oxy is coincident with the undisturbed free surface all the time.(3)Coordinate System 3: the body-fixed system G-xbybzb is fully fixed on the moving ship. G is coincident with the center of gravity (COG) of the ship. The xb, yb, and zb axes are directed toward the bow, the port side, and the sky, respectively.

For consistency, the relationship between values in these three coordinate systems is expressed as follows:where t denotes the time and zG denotes the height of COG of the ship with respect to the undisturbed free surface.

3.2. Solution of Fluid Velocity Potential

Since the fluid is ideal, the fluid velocity vector can be expressed by the gradient of the velocity potential:

The steady velocity potential corresponding to the steady wave-making flow component is ignored in this study. The disturbance velocity potential around the ship can be further decomposed into incident wave potential, diffraction wave potential, and radiation wave potential in the plane movement system o-xyz:where ζa denotes wave amplitude, ωe denotes encounter frequency, ϕ0(x, y, z) denotes incident wave potential under unit wave amplitude, ϕd(x, y, z) denotes diffraction wave potential under unit wave amplitude, and ϕr(x, y, z) denotes radiation wave potential under unit principal coordinate amplitude for rth order motion mode. It is noted that for r = 1∼6, the motions are rigid body six degree-of-freedom (DOF) motion, while for r = 7∼m, the motions are structural elastic deformation.

The incident wave potential ϕ0 is fully determined by the known incident wave field. For finite water depth h and infinite water depth conditions, the incident wave potentials are, respectively, expressed as follows:where k denotes wave number.

For the convenience of formulism, the diffraction potential ϕd(x, y, z) is replaced by ϕm+1(x, y, z) in the following statements. The radiation potential ϕr(x, y, z) (r = 1∼m) and diffraction potential ϕm+1(x, y, z) should satisfy the following fluid domain and boundary conditions. The fluid domain continuous condition, free surface condition, wetted body surface condition, bottom boundary condition, and infinite boundary condition are summarized as follows:where n denotes inward-directed unit normal vector on the wetted hull surface and ur denotes displacement vector of hull surface caused by the rth order motion or deformation.

In this study, the source distribution method is used to solve the flow field velocity potential. The velocity potential can be obtained by integrating the distributed source over the body surface:where σ(j)(q) denotes distributed source strength; j denotes the boundary element panel number on the body surface; and G(p,q) denotes frequency-domain Green function which satisfies all the boundary conditions above-mentioned except for the body surface boundary condition, the field point p is on the hull surface, the source point q is in the fluid field.

According to the hull body surface boundary condition, the distributed source strength should satisfy the following boundary integral equation:

The Green functions for both finite water deep and infinite water deep conditions are employed for calculation. The Hess-Smith panel element method is used to solve the distributed source strength in Equation (24). Once the distributed source strength is determined by solving the Green function, the velocity potential components can be obtained. The added mass and damping coefficient can be obtained by the known radiation potential. The pressure components for incident, radiation, and diffraction wave acting on elastic wetted hull surface can be obtained by the Bernoulli equation once the velocity potentials are obtained.

3.3. Nonlinear Hydrostatic and Hydrodynamic Forces

The nonlinear effects of hydrodynamic forces for the ship sailing in severe waves at high speed are obvious. Therefore, the hydrostatic and hydrodynamic forces, including the hydrostatic restoring force, incident wave force, diffraction wave force, and radiation force, are calculated in real-time to take into consideration the nonlinear effects attributed to body geometry.

The real-time hydrostatic restoring force acting on the elastic hull surface is calculated as follows:where S(t) denotes the instantaneous wetted body surface of the elastic hull, pka is the amplitude of mode principal coordinate, denotes the vertical displacement component of uk (i.e., (, , )), the hydrostatic restoring force coefficient element Crk(t) is further expressed as follows:

The incident wave force acting on the instantaneous wetted elastic hull surface is expressed as

The diffraction wave force acting on the instantaneous wetted elastic hull surface is expressed aswhere the diffraction potential ϕd(x, y, z, t) is for the value of steady state solution of the instantaneous wetted surface of hull.

The radiation force acting on the instantaneous wetted elastic hull surface is expressed as

The hydrodynamic coefficients in real-time are computed based on the instantaneous static wetted body surface condition by using frequency-domain Green function method. The real-time hydrodynamic coefficients are expressed as follows:

3.4. Roll Damping Correction

As aforementioned, the potential flow theory ignores the fluid viscosity. In fact, the roll damping of the ship is associated with strong nonlinear behavior and should be considered when the ship sailing in beam or oblique waves. The roll damping coefficient B44 calculated by Equation (30), which is only responsible for wave-making radiation force component, only accounts for small proportion of the total roll damping. The following linearized roll damping correction coefficient is adopted to reflect the nonlinear behavior of roll motion. The nonlinear roll damping moment is expressed as follows:where p4 denotes angular displacement of roll motion, a and b are specified coefficients which can be determined by experiment or experiential method, and kφ denotes linearized equivalent roll damping coefficient.

In this study, Muller’s method is used to determine the roll damping coefficient [30]. The linearized equivalent roll damping coefficient is expressed aswhere p4a denotes amplitude of roll motion in specified regular waves and k1 and k2 are expressed as follows:where L denotes ship waterline length, B denotes moulded breadth, T denotes draft, Cb denotes block coefficient, Fr denotes Froude number, hx denotes transverse metacentric height, and Rb, lb, and Ab are bilge keel related coefficients, Cv = 4.85 − hx1/2 and C44 = ρghx∇.

3.5. Impact Loads for Slamming and Green Water on Deck

Generally, Wagner theory and momentum impact theory are the classical methods in use to acquire the slamming pressure and slamming force, respectively. In this study, the momentum impact theory is adopted to calculate the sectional impact force. In the theory, the slamming force is determined by the momentum changing rate of the fluid around the bow. Thus, the sectional slamming force is obtained as follows:where, m(x) is the sectional added mass in the heave mode at infinite frequency and is the vertical displacement of the hull relative to the wave surface elevation.

The selected hull sections for 2D slamming estimation are shown in Figure 4. The overall flare slamming loads acting on the whole elastic hull are obtained by integrating the 2D sectional slamming force longitudinally over the ship length:where is vertical normal vector component of n.

In harsh wave conditions, the incident wave may overtop the bow and run up onto the deck after the occurrence of bow slamming event. In this study, the green water load on deck is estimated by 1D dam-breaking model. This method treats the green water flow as the fluid subject to a sudden collapse of a dam, the phenomenon of which is illustrated in Figure 5. The water is considered to be static before the collapse of dam. The dam-breaking flow is governed by the 1D Saint-Venant’s equations, which are written as follows:where u denotes the horizontal velocity component of the green water flow and h denotes the water height on the deck, and they are, respectively, determined by equations as follows:where H0 is the equivalent dam height before breaking.

The green water pressure on deck pgw can be obtained based on the Bernoulli equation. Then the overall green water loads can be calculated by integrating the pressure over the whole deck surface:

3.6. Solution of Time-Domain Motion Governing Equation

The time-domain motion governing equation of the elastic hull in regular waves, including both rigid motion and elastic distortion, can be expressed aswhere [a], [b], and [c] are generalized structural mass, damping, and stiffness matrices, respectively; [A], [B], and [C] are generalized fluid added mass, damping coefficient, and hydrostatic restoring matrices, respectively; {pr(t)} is the rth mode principal coordinate column matrix; {FI(t)}, {FD(t)}, {Fslam(t)}, and {Fgw(t)} are incident wave force, diffraction wave force, slamming loads, and green water loads column matrices.

The structural damping has significant influence on high-frequency hydroelastic vibrations of the ship [31]. The structural damping coefficient not only affects the decaying time of whipping loads, but also influences the amplitude of springing loads. The structural damping matrix [b] (assumed diagonal) can be determined experimentally or empirically. In this study, the logarithmic decrement of peak value of high-frequency vibrations concluded by Hirowatari is used:where ωr denotes rth order natural frequency.

The diagonal element in the structural damping matrix is obtained by:where arr denotes rth diagonal element in the structural mass matrix.

The fourth order Runge–Kutta method is employed to solve the time-domain nonlinear differential equation of Equation (39). Once the principal coordinates for different modes are obtained, the displacement, moment, and shearing force of a given cross section can be obtained by the following modal superposition method:where ωr(x), Mr(x), and Vr(x) are the rth natural mode of the model’s displacement, moment, and shearing force, respectively.

4. Tank Model Experimental Setup

In order to validate the hydroelasticity theory algorithm, a segmented model was designed and manufactured according to a bow-flare ship to allow the conduct of the corresponding physical experiment. The experiment was conducted in a laboratory wave basin for regular head wave tests. The experimental details including model design, tank facility, and testing procedure are reported as follows.

4.1. Scaled Model Design

The design of scaled model hydrodynamic test is guided by Froude’s similitude law. A 1:50 scaled model is designed according to geometric similarity, kinematic similarity, dynamic similarity, and vertical bending vibration modal similarity principles. The longitudinal distribution of weight and stiffness of model is also designed to be similar with the prototype. Main dimensions of the prototype and model ship are listed in Table 1.

The model has 20 stations and is cut into seven parts of segments by divisions of 2nd, 4th, 6th, 8th, 10th, and 12th. A flexible steel backbone is used to connect the discontinuous segments. The backbone is elaborately designed with varying cross section so as to match the natural frequency and longitudinal stiffness distribution of the prototype in vertical bending vibration mode. The fluid forces subjected by the segmented hulls are fully transferred to the continuous backbone beam. Gaps of 15 mm wide are provided between adjacent segments to prevent their contact due to elastic deformation of backbone in waves. The gaps are sealed by latex rubber for waterproofing. The sectional wave and slamming loads at cut divisions are measured by strain gauges mounted on the backbone surface. An array of pressure sensors is mounted on the hull surface at the bow-flare area to measure the bow slamming pressures. Three accelerometers are mounted on bow, middle, and aft area of deck to measure the vertical acceleration. A set of propulsion system is installed at stern area from divisions 13th to 20th. The overall arrangement of model setup is illustrated in Figure 6.

The model hull shell is made by Fiberglass-Reinforced Plastics (FRP). The backbone beams are made by hollow rectangular and cylindrical steel tubes, and they are fixed on the segment’s bases rigidly by aluminum fixing plate. The backbone beam is calibrated by applying known static loads prior to the assembly of segments. The ballast iron is installed at the desired positions to adjust the model’s weight distribution and moment of inertia. View of the segmented model before tank measurement is shown in Figure 7.

4.2. Laboratory Tank Facility

The regular wave experiments were conducted in the towing tank of Harbin Engineering University (HEU). The tank has a dimension of 108 m long, 7 m wide, and 3.5 m deep. Both long-crested regular and irregular waves can be generated by a single flap type hydraulically driven wave-maker. There is a wave absorbing beach at the opposite side of the wave-maker. A resistive wave probe is installed in situ near the wave-maker to measure the surface elevation of the actual generated waves. Another wave probe is installed onboard the carriage to measure the encountered waves suffered by the model. View of the tank facilities is shown in Figure 8.

The model’s motions are measured by an in-housed developed 5-DOF (i.e., heave, pitch, roll, sway, and surge) motion measurement device, which is installed on the measurement bridge of towing carriage. The model is attached to the 5-DOF motion measurement device by two heave sticks. The motion measurement device guides the model sailing heading and acts as speed reference during model running. The self-propelled model’s forward speed is achieved by its four propellers. The scheme of model setup on the 5-DOF measurement device is shown in Figure 9. The measured data, including ship motions, vertical accelerations, sectional loads, and impact pressures, are recorded and stored by a 32-channel data collector.

4.3. Testing Scheme and Procedure

The regular wave testing scheme is determined by ship speed, wave states (wave height, wavelength, or frequency), and wave heading angle. The regular head wave tests are conducted for both zero-speed and forward speed conditions in this study. The testing schemes involved are summarized in Table 2. As is seen, the tests over a wide range of wave frequency (λ/L = 0.4, 0.6, 0.8, 0.9, 1.0, 1.1, 1.2, 1.5, 2.0, and 2.5) are conducted in low wave height conditions in order to obtain the linear Response Amplitude Operator (RAO) of ship, whereas the tests are conducted around the resonant frequency range (λ/L = 0.8, 0.9, 1.0, 1.1, and 1.2) in high wave height conditions to obtain the ship responses in harsh and extreme testing conditions.

The tests were conducted by a team of experienced researchers and staff in the towing tank laboratory of HEU. The testing procedure for each set of running measurement is described as follows:(a)The carriage and model are situated at the start side (opposite the wave-maker) of the towing tank at first. The generated waves will propagate towards the model once the wave-maker started working.(b)The carriage will start running when the waves are about to reach the model. The start time of carriage is comprehensively determined by the model testing speed and wave frequency so as to reduce the influence of wave reflection and also achieve maximum effective testing measurement distance. The data collector will start to record the model’s responses prior to model running.(c)During the acceleration phase of carriage, the model is towed by ropes tied on its centerline at bow and stern. When the carriage reaches the desired speed, the ropes will be released. Then the model forward speed is achieved and maintained by its propellers.(d)When the carriage is about to approach the end of the running, the ropes are recovered and tightened, after which the carriage starts deceleration. Finally, the carriage and model will return to the start point and wait for calm water and prepare for the next running measurement.

5. Overview of Obtained Numerical and Experimental Results

The necessary details during the conduct of numerical simulation and experimental measurement are reported in this section. The vertical vibration mode characteristics are reported at first. Then overview of the obtained typical numerical and experimental results in regular waves is presented.

5.1. Vibration Mode Analysis

For simplification, the flexible hull girder is described by a 1D nonuniform Timoshenko beam. The combination of 1D beam theory with 3D BEM turned out to be applicable and reliable for ship hydroelastic analysis as a compromise between the calculation accuracy and efficiency [32]. In the following illustrative theoretical calculation, only the first three order vertical bending modes are considered for the ship in head waves. The dry mode characteristics of full-scale ship in vertical bending mode are calculated by TMM to obtain the natural frequency and corresponding principal mode shape. The obtained 1st∼3rd order principal mode shapes for sectional distortion and force of prototype with 20 sections of equal length in dry condition are shown in Figure 10. The modal shapes for displacement and rotation angle are normalized by setting the value at 0th station (forward perpendicular) as unit; the 1st and 3rd order VBM and 2nd order VSF modal shapes are normalized by setting the value at 10th station (amidships) as unit; and the 2nd order VBM and 1st and 3rd order VSF modal shapes are normalized by setting the value at 6th station (quarter ship length) as unit. The vibration natural frequencies of scaled model are obtained based on those of full-scale according to the similitude law. The model’s backbone system is designed to match the target frequency of two-node (1st order) natural vibration.

Impact hammer test was performed on the model to identify the actual natural frequency of model vibration in wetted condition. The test is performed by hitting the bow of model in calm water and recording the subsequent stress decaying curve. The corresponding frequency-domain information can be obtained by applying the fast Fourier transformation (FFT) on the recorded stress decaying curve. The measured vertical bending stress decay time series at the six different measurement stations (i.e., 2nd, 4th, 6th, 8th, 10th, and 12th) and the corresponding frequency-domain results are presented in Figure 11. The vertical vibration natural frequency of the model in the wetted condition can be clearly read from the peak frequency: the two node (1st order) natural frequency is 3.87 Hz, and the three node (2nd order) natural frequency is 9.38 Hz. They are close to the designed values of wetted natural frequency.

5.2. Numerical Simulation Results

A 3D hydroelastic analysis code was in-house developed using FORTRAN language based on the nonlinear time-domain hydroelasticity theory established in Section 3. The 3D hull hydrodynamic grid is generated by means of the cubic spline curve method [33]. Figure 12 shows the generated hull surface meshes, and there are 1802 hydrodynamic panels in total. The hydrodynamic hull model was divided into 21 blocks by 20 divisions, and the sectional loads at these divisions are calculated. The mass of each block is input so as to calculate and verify the COG and inertia of the whole ship. The sectional area, position of shear center, and centroid are also input for sectional load calculation.

Typical numerical simulation results for the full-scale ship sailing at a forward speed of 5 knots in the harsh head wave condition (wave height H = 12 m and wavelength/ship length λ/L = 0.8) are illustrated as follows. Figure 13 shows the motion results at COG of the full-scale ship in regular head waves. The rigid body global motions (heave and pitch) are shown in Figures 13(a) and 13(b), and the elastic deformations (1st and 2nd order vibrations) are shown in Figures 13(c) and 13(d). As is seen, the motion responses trend to be stable approximately since the 500 s. Figures 13(e) and 13(f) show the comparison of vertical displacement between rigid body motion (heave) and elastic deformation (1st∼3rd order vertical vibrations). The comparison indicates that the amplitude of vibration is much smaller compared with the rigid body motion, and the amplitude of vibration decreased rapidly from the 1st order vibration to the 3rd order vibration. Therefore, the 1st order vibration should be considered in the simulation of large flexible ship’s motion and load responses, while the high-order vibrations can be ignored in ordinary engineering application.

The calculation results for VSF and VBM at amidships (the station 10) of the ship in regular head waves are shown in Figures 14(a) and 14(b), respectively. The corresponding short-term zoom view of the time series is shown in Figures 14(c) and 14(d), respectively. As seen from the results, the responses trend to be stable approximately since the 500 s. To be different with the rigid body global motion, there are obvious high-frequency load components associated with the wave frequency loads due to the bow slamming effects in severe waves, and the slamming loads have more influence on VBM than VSF. Moreover, the slamming loads largely increase the sagging load of VBM.

5.3. Experimental Measurement Results

The experimental results for the same condition scheme as that presented in Section 5.2 (i.e., full-scale ship speed V = 5 knots, wave height H = 12 m, and wave length/ship length λ/L = 0.8) are illustrated as follows. It is noted that the results are presented at the model scale without extrapolation. The incident waves recorded by the in situ and onboard wave probes are shown in Figures 15(a) and 15(b), respectively. As seen from the waves measured by the in situ wave probe, the wave field is very stable and the wave parameters (i.e., wave height and period) coincide with the expected values. However, the wave peaks measured by the onboard wave probe fluctuate in a certain range due to the nonlinear second-order wave drift effects, wave reflection effects, and ship wave-making interference.

The heave and pitch motions at COG of the model in the above wave condition measured by the 5-DOF measurement device are shown in Figure 16. As is seen, the heave and pitch motions of the model show the same tendency as the wave time series measured by the onboard wave probe. The peak fluctuates over a certain range, and an envelope can be clearly observed. The reasons are as follows: (i) the model forward speed is low at 0.364 m/s, (ii) the wave length is long (λ = 4.672 m), and the wave is steep (H/λ = 0.051), (iii) the resonant frequencies of heave and pitch are close to the wave encounter frequency, and (iv) the surge motion of the model is released. Thus, the model needs much time to reach a stable motion state. Moreover, nonlinear second-order wave force, wave reflection effects, and experimental uncertainty may also contribute the irregularity of measured time series. The motions reveal relatively steady tendency since the time of the 40 s. Therefore, the average amplitude values (crest and trough values) after 40 s are calculated for further frequency-domain response analysis.

The measured time series of VBM loads at different stations (i.e., 2nd, 4th, 6th, 8th, 10th, and 12th) are summarized in Figure 17. As is seen, the load signal reveals obvious nonlinearity characteristics due to the high-frequency slamming load components. The high-frequency loads, mainly contribute the sagging loads of VBM, are very pronounced especially at stations 2 and 4 due to the occurrence of bow slamming. Moreover, the slamming loads at stations 2 and 4 have a strong randomicity, and thus, the sagging load for each slamming event fluctuates significantly even for a regular wave-testing scheme.

The extreme values (includes peak and valley values) of total VBM at different stations are extracted based on the time series during a steady run range of 20∼45 s. The amplitude values (includes hogging and sagging values) of linear wave-frequency VBM at different stations are obtained by Fourier’s filter and read from the corresponding filtered time series. Moreover, the ratio coefficient of the extreme value of total load to the amplitude value of wave-frequency load is calculated. The distribution of extreme values of total load and amplitude values of wave-frequency load for VBM along ship length is shown in Figure 18. As is seen, the sectional VBM increases almost linearly from station 2 at bow to station 10 at amidships; and the largest sectional VBM occurred at station 10 for both peak value and amplitude value. The ratio coefficient of hogging VBM fluctuates slightly around 1 over the whole ship length; while the ratio coefficient of sagging VBM decreased rapidly from station 2 to station 4 due to the large slamming load components at the bow area.

The bow slamming and green water phenomena recorded by video cameras are shown in Figure 19. As seen from the carriage fixed video camera’s recordings, during one slamming period, the bow of the model emerged from the water due to large amplitude motion of hull in pitch and heave modes (Figure 19(a)) and then reentered the water at a relative high velocity (Figure 19(b)). The bow slamming and green water on deck phenomena recorded by the deck fixed camera are shown in Figures 19(c) and 19(d), respectively.

6. Comparison of Numerical and Experimental Results

The results for ship motions and loads in different conditions obtained by theoretical and experimental methods are compared in this section. The numerical algorithm is also well validated by comparing with the model testing results.

6.1. Global Motion Responses

The RAOs for ship sailing in head wave conditions at a forward speed of 5 knots are compared between calculation and measurement. The wave frequency in numerical simulation ranges from 0.02 rad/s to 2 rad/s with a step of 0.02 rad/s. The testing wave frequency is selected at points λ/L = 0.4, 0.6, 0.8, 0.9, 1.0, 1.1, 1.2, 1.5, 2.0, and 2.5. Figure 20 shows the comparison of pitch and heave RAOs between calculation and measurement. As is seen, generally, the results show good agreement between calculation and experiment. The heave experimental results are generally a little larger than the calculated over the whole λ/L range. The pitch experimental results are in very good accordance with calculation apart from that they are slightly higher at point λ/L = 2.5.

The comparative full-scale pitch and heave time series between calculation and measurement in different wave height conditions (V = 18 knots; H = 4 m and 12 m; λ/L = 0.8) are shown in Figures 21 and 22. As is seen, the time series show good agreement between calculation and measurement especially for the H = 4 m condition. While for the H = 12 m condition, the peaks of the motion signal by experiment fluctuate over a small range due to the experimental uncertainty caused by wave nonlinearity and ship high forward speed.

6.2. Sectional Load Responses

Figure 23 shows the comparison of VBM responses between calculation and measurement for the ship in head wave conditions at a forward speed of 5 knots. Figure 23(a) shows the RAO for VBM amidships, and Figure 23(b) shows the longitudinal distribution of VBM response amplitude along the ship for the specified λ/L = 1.0 condition scheme. As seen from Figure 23(a), the VBM increases shapely from λ/L = 0.4 to 0.8, and the largest VBM occurred in the range of λ/L = 0.8∼1.0. The measured sagging values are generally a little larger than the hogging values especially around the peak frequency region due to nonlinearity effects of waves and hull structure. As seen from Figure 23(b), the measured values are slightly lower than the calculated values ahead of the 8th station while becomes higher than calculation for the values at stations 10 and 12. On all account, the VBM results show good agreement between calculation and measurement.

The comparison of time series of sectional VBM amidships (at station 10) between calculation and measurement for different wave height conditions (V = 18 knots; H = 4 m and 12 m; λ/L = 0.8) is shown in Figure 24. As seen from the results, the sectional loads show good agreement between calculation and measurement. To be different from the global motion signals, the sectional loads comprise high-frequency slamming components besides the low frequency wave loads. And the high-frequency slamming loads are more obvious for the H = 12 m condition due to the severe slamming events occurred in high waves.

The comparison of different VBM load components (low-frequency wave loads and high-frequency slamming loads) at amidships for schemes (V = 18 knots; H = 4 m and 12 m; λ/L = 0.8) are plotted in Figures 25 and 26. It is noted that the low-frequency wave loads and high-frequency whipping loads are separated from the total loads by employing Fourier filtering. As is seen, the wave-frequency loads are almost linear to the incident wave height, since the double amplitudes (peak-to-peak value) of wave-frequency VBM are about 2.07 GN·m and 6.51 GN·m for 4 m and 12 m wave heights, respectively. The slamming loads are minor in the H = 4 m condition when compared with the wave loads. While in the H = 12 m condition, severe slamming phenomenon occurred and the whipping loads are almost in the same magnitude as the wave loads. The whipping loads decayed rapidly after slamming impact due to the structure damping effects. To summarize, the total loads, wave loads, and whipping loads show satisfactory agreement between calculation and measurement.

6.3. Slamming Pressure

The comparative bow-flare slamming pressures at positions of sensor nos. 1 and 2 (Figure 6) in the severe wave condition (V = 18 knots; H = 12 m; λ/L = 0.8) by the generalized Wagner theory calculation [32] and measurement are shown in Figure 27. In the numerical simulation, the hull girder is, respectively, assumed to be rigid and flexible in order to investigate the influence of hull elasticity on slamming pressure responses. The results indicate that the slamming pressure peak by elastic hull theory is slightly lower than that by rigid body theory. This is due to the fact that the elastic effect of hull girder will buffer the transient impact pressure. Therefore, the results by the hydroelasticity method are more reasonable and reliable, and they are much closer to the measurement results.

7. Perspective of Ship Hydroelasticity in Irregular Waves

Although this paper presents detailed theoretical and experimental approaches for the prediction of ship hydroelasticity in regular waves, the prediction of ship hydroelasticity in irregular waves will be more insightful and practical. Accurately prediction of motion and loads performance of the ship in irregular waves is of great importance for full-scale ship design and evaluation. Our research project is still in progress. We are devoted to investigating ship hydroelasticity and slamming loads in both regular and irregular waves theoretically and experimentally. In our recent work, the theoretical and experimental investigations of ship hydroelasticity in both long-crested and short-crested irregular waves are undertaken, which is briefly prospected as follows.

7.1. Theoretical Approach

To date, ship motion and load responses in random waves are usually predicted on the assumption that the incident waves are long-crested irregular waves [3436]. The realistic sea waves are however short-crested irregular waves. Real practice reveals that the ship motion and load responses induced by short-crested irregular waves are different from those induced by long-crested irregular waves; and the wave directional spreading has significant influence on ship motion and load responses. With this in mind, in our recent work, a scheme of 3D time-domain nonlinear hydroelasticity theoretical methodology is proposed for the prediction of ship motions and loads in both long-crested and short-crested waves.

According to the linear superposition theory, the 3D irregular wave field is composed of a number of regular waves with different frequencies and different spreading directions:where ωi denotes the frequency of the ith (I in total) component regular wave, ki denotes the wave number of the ith component regular wave, βj denotes the spreading direction of the jth (J in total) component regular wave, aij denotes the amplitude of the ijth component regular wave, and εij denotes the phase of the ijth component regular wave.

A nonlinear hydroelasticity theory is established to predict ship motions and loads in irregular waves with instantaneous wetted surface condition, nonlinear wave forces, slamming, and green water loads considered. The flow field velocity potential around the elastic hull is solved by the Rankine panel method to consider ship forward speed effects. The fourth order Runge–Kutta algorithm is used to solve the nonlinear governing equation in time domain. The nonlinear hydroelastic governing differential equation of the ship sailing in irregular waves is expressed as follows:where [A] denotes the infinite frequency-added mass matrix, [B] denotes the infinite frequency damping coefficient matrix, K(τ) denotes the time retardation function matrix, τ denotes the time interval, and other symbols are the same meanings as those defined in (39).

The real-time hydrostatic restoring force acting on the elastic hull surface in 3D irregular waves is calculated by integrating the pressure over the instantaneous wetted hull surface:

Considering the 3D effects of the flexible hull structure and irregular wave field, the incident wave force is obtained by integrating the contribution of the component wave at different frequencies and different spreading directions. Thus, the time-domain convolution integral method is used to obtain the incident wave force induced by short-crested irregular waves:where the impulse response function of the incident wave for the rth motion mode can be obtained by applying the Fourier transformation using the frequency-domain response function (i.e., RAO) in regular waves:

The diffraction wave force of the ship in short-crested irregular waves can be obtained in the same manner:

The radiation wave force of the ship in short-crested irregular waves is obtained by the indirect frequency-domain method. To take the wave memory effects into account, the radiation wave force is expressed aswhere the element Krk(τ,β) in the time retardation function matrix can be obtained by applying the Fourier transformation using the frequency-domain damping coefficient Brk(ω,β):

The slamming loads Fslam(t) and green water loads Fgw(t) of ship sailing in short-crested irregular waves is a reason of concern. The slamming loads can be also obtained by the momentum impact theory given in (34). The vertical relative velocity of the ship with respect to the 3D wave surface can be obtained by derivating the relative displacement:where denotes the rth mode coordinate of the ship surface under the plane movement system and the wave elevation ζa(x, y, t) can be obtained by (43).

7.2. Experimental Measurement

The tank model experiments were conducted in three different tanks using the small-scale segmented model, described in Section 4, to fulfill all the testing conditions, which include long-crested regular and irregular wave testing conditions. These three tanks are (i) the towing tank (108 m long by 7 m wide by 3.5 m deep) of HEU, (ii) the comprehensive deep ocean basin (50 m long by 30 m wide by 10 m deep) of HEU, and (iii) the high-speed hydrodynamic tank (510 m long by 6.5 m wide by 4 m deep) of Aviation Industry Institute No. 605. Views of model setup in the three tanks are shown in Figure 28. The testing conditions for head and following regular waves were conducted in the first tank (reported in Section 4.2). For beam and oblique wave testing conditions, the tests were conducted in the second tank with the help of the two orthogonal carriages so as to provide arbitrary sailing heading. While for long-crested irregular wave testing conditions, the tests were conducted in the third ultralong tank to achieve enough signal samples for statistical analysis [37].

To allow the experimental investigation of ship hydroelasticity and slamming loads in short-crested irregular waves, a large-scale segmented model was constructed at a scale of 1/25 according to the same prototype [38]. For the comparative investigations, model setup and sensor arrangement of the large-scale model are designed to be the same as the small-scale model (at a scale of 1/50), as described in Section 4. A directional wave buoy is used to measure the short-crested sea waves. The sea waves, ship navigational state, ship global motions, sectional loads, bow slamming pressure, and green water phenomenon are measured and recorded during the sea trials for postvoyage analysis. The large-scale testing model and sea trial scene are shown in Figure 29.

8. Conclusions

This paper focuses on theoretical and experimental investigations of hydroelastic responses and slamming loads behavior of large ships sailing in regular waves. This study also lays a foundation and provides fundamentals for the ship structure design and strength evaluation. The following conclusions can be made from this study:(a)The developed fully coupled 3D time-domain nonlinear hydroelasticity theory, which combines 3D BEM, 3D FEM(or 1D nonuniform Timoshenko beam model), 2D momentum impact model, and 1D dam-breaking model, is an effective tool in the estimation of the nonlinear wave and slamming loads of the bow-flare ship sailing in regular waves.(b)Although the 3D FEM has good accuracy in the simulation of the vibration response of the complex hull structure, it is complex, expensive, and time-consuming. The 1D nonuniform Timoshenko beam turned out to be applicable and reliable for ship modal analysis as a compromise between the calculation accuracy and efficiency.(c)The established segmented model and tank experimental system are reliable and capable for the measurement of global motions, wave loads, slamming loads, and impact pressure acting on the ship sailing in waves. The numerical results are successfully validated by the tank testing results.(d)The time series of motions are smooth and stationary even for the ship sailing in high wave states. However, the sectional loads comprise considerable high-frequency components besides the wave-frequency component especially for high speed or high wave state conditions. Therefore, the slamming loads should be concerned during the ship structure design and operational guidance.

Appendix

The expressions of λ1 and λ2 in the vertical displacement equation in (11) are expressed as follows:where the coefficients in the equation are expressed as follows:

The transfer matrix Fi in (12) is expressed as follows:where the coefficients in the matrix elements are expressed as follows:

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that there are no conflicts of interests regarding the publication of this paper.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant no. 51709118), the China Postdoctoral Science Foundation (Grant nos. 2017M622696 and 2017M612669), the Foundation for Distinguished Young Talents in Higher Education of Guangdong, China (Grant no. 2017KQNCX004), the Funds for Economic Development of Guangdong, China (Grant no. GDME-2018B003), and the Science and Technology Program of Guangzhou, China (Grant no. 201804010482).