Simulation based calculation of ship motions in extreme seas with a body-exact strip theory approach

For all design phases of naval vessels, the fidelity of seakeeping calculations in extreme seas is open to discussion due to the inadequacy of the linear theory of ship motions. Currently the computer-generated time series of ship responses and wave height (the real time computer experiments) are utilized to calculate the distribution of the vertical distortion, shear force and bending moment by means of “ship hydroelasticity theory”. Inspired by these studies a simulation based calculation of symmetric ship motions is performed in long crested irregular head seas, in addition with a body-exact strip theory approach. The scope of this study is limited to the ship motions only. Verification is achieved utilizing the spectral analysis procedure which contains the discrete Fourier transform (DFT) and the smoothing algorithms. The results are compared with the experimental data, and the ANSYS AQWA software results. The simulation results provide adequate data for the extreme responses. This state-of-the-art method in addition with a “body-exact strip theory approach” ensures the consistent assessment of the seakeeping performance in extreme sea condition. As a result, it is evaluated that this calculation method can be used in the design stages of naval platforms.


Introduction
Ship designers aware of that the linear theory is not sufficient for ship motion calculations in extreme seas. Reliability of the ship motion calculations is important for the design and verification of the ship itself and the systems on board. Theoretical non-linear methods for seakeeping calculations mainly include "Non-linear 3-Dimensional (3D) Methods" and "Non-linear Strip Method Extensions" [1]. Non-linear 3D methods are performed in frequency or time domain. Among the 3D methods in the time domain, the "Large Amplitude Motions Program (LAMP)" developed by Lin and Yue [2] is remarkable. In this approach, the perturbation potential and the hydrostatic/Froude-Krylov forces can be found by combinations of solutions "over mean wetted hull surface" and "over the instantaneous wetted hull surface". The solution over mean wetted hull surface is called the "body-linear" solution, while the solution over the instantaneous wetted hull surface is called the "bodynonlinear" approach [3]& [4].
In this paper, it is aimed to use a similar approach to the non-linear strip theory extensions, which are faster than the 3D solution. In this context, as cited in STANAG 4154 [1], studies carried out by Kaplan and Sargent [5], Mansour and d'Oliveira [6] (The computer program ADMASS for the frequency-dependent added mass calculation), Meyerhoff and Schlachter [7], Schlachter [8]& [9] (The computer program DYNBEL) and Guedes Soares [10] on the vertical plane; studies carried out by De Kat [11]& [12] and De Kat and Paulling [13] on the six degrees of freedom are first coming to the forefront. In these approaches, as in the case of the 3D methods, the search for a solution is made "over the instantaneous wetted hull surface". De Kat and Paulling [13] described this solution search as a "geometric" nonlinearity. The solution of the equations of motion in the time domain may require the use of convolution integrals that take into account the memory effects. An example to this approach, in which the nonlinear effects in large amplitude motions are taken into account and the strip theory is used in the time domain with the convolution integrals, is the stud-ies of Fonseca and Guedes Soares [12] [14]& [15]. In these studies the hydrostatic/Froude-Krylov forces are solved over the instantaneous wetted hull surface. This approach, which is also used in studies carried out by Rajendran et al. [16]& [17], which is the continuation of these studies, is described as "body-nonlinear time domain method". However, in the work done by Bandyk [18] and Zhang et al. [19], the solution search over the instantaneous wetted hull surface is called the "body-exact approach". Zhang et al. [19] describes the body-exact approach as "compromise between fully nonlinear computations and linear theory". Furthermore, in the nonlinear approach (which uses conformal mapping) applied by Fan and Wilson [20], the instantaneous wet portion of each section under the wave surface is updated at every time step. In this way, the two dimensional sectional added mass and damping coefficients and wave excitation forces can be calculated depending on the time.
On the other hand, Belik [21], Belik and Gökmen [22] and Başaran [23] investigated the ship responses in ex-treme seas and utilized the computer-generated time series of ship responses and wave height to compute the distribution of the vertical distortion, shear force and bending moment by means of "ship hydroelasticity theory" proposed by Bishop and Price [24]. These simulated time series are termed as "real time computer experiments" [22] & [23]. The significance of these simulations based calculations are their practical and state-of-the-art nature as pointed out by Belik and Gökmen [22]. Without any cost and constraint except the computer power, the desired length of simulation can be performed for all ship responses, as well as the derived responses. Besides this "time domain" approach, Jiao et al. [25] used a linear "frequency-domain" hydroelasticity theory (the Linear Elastic Wave Loads Calculation System, WALCS-LE) to predict the wave induced motions and load responses.
Motivated by the study of Belik [21], Anıl et al. [26] and Anıl [27] performed a simulation based analysis of ship motions (in short crested irregular seas) which focuses   on the clarification of the theoretical background of this process. This paper is the continuation of these previously published researches using the same methodology, but in addition with "a body-exact strip theory approach". The scope of work is limited to the "symmetric" ship motions in long crested irregular head seas.
As clarified by Belik [21] and Newland [28], the discrete Fourier transform (DFT) and smoothing algorithms are utilized for the verification and comparison of the simulated data with the experimental data [29], and the ANSYS AQWA software results of Çekirdekçi [30].

Methodology
In this study, symmetric ship motions in long crested irregular head seas are calculated using the time domain simulation data in addition with a body-exact strip theory approach. The two dimensional sectional added mass and damping coefficients are calculated by the Frank Close-Fit method over the instantaneous wetted hull surface. Hydrostatic/Froude-Krylov forces are also solved over the instantaneous wetted hull surface.
The time series data is generated for both body-linear and body-exact condition. The time series data are trans-ferred to the frequency domain using the spectral analysis method by discrete Fourier transform (DFT) and smoothing algorithms. The DTMB 5415 hull form [31] (in forward motion with a Froude number of 0.41), which is the preliminary form of the DDG-51, is used for the current study. The three Cartesian coordinate systems defined by Beck et al. [32] are used for the simulation. The zaxis is designated to be through the center of gravity, as Salvesen, Tuck and Faltinsen [33] indicated.
The time series data for the body-linear condition is generated from the initially computed frequency domain Response Amplitude Operators (RAOs) and the corresponding motion spectra. The body-linear time series data is used to for the verification of the spectral analysis.
The time series data for the wave surface in moving (inertial) coordinate system is only used for the production of the time series data for the body-exact condition. At each time step, the intersection points of the wave and the ship sections are transformed into the body-fixed coordinate system. The strip theory calculation along with the Frank Close-Fit method is updated for each time step using the instantaneous wet portion of each section as shown in Fig.  1. As seen in figure, the number of the sections is 23 and the distance between the sections are closer at the bow. In each section the number of elements (the discretized two dimensional panels) is 60. Therefore, totally 1380 elements are used for the calculations. The convergence test (which is performed for each frequency) of the added mass and damping coefficients with respect to the total number of elements (maximum 2300) is shown in Fig. 2.
During the simulation, problem is solved by the Runge-Kutta-Fehlberg method (RFK45) [34] in the time domain. The initially calculated linear frequency domain RAOs are used at t = 0 as an "initial condition" but the problem discussed here is not an "initial value problem". Therefore, this approach can be called a "time-varying" strip theory. The flow diagram of the method is presented in Fig. 3 and Fig. 4.

Theoretical Background
The following sub-sections explain the theoretical background of the simulation based calculation of ship motions with a body-exact strip theory approach.

Equations of Motion
The linearized differential equations for the ship motions (in body-fixed coordinate system) can be written in the form as follows [32], [33]: where j and k are the motion index. The components of the generalized mass matrix M jk , added mass A jk , damping B jk and hydrostatic restoring coefficients C jk are given by Salvesen, Tuck and Faltinsen [33] in detail. F j are the complex amplitudes of the exciting forces and moments. Here η k , η . k and η .. k are the amplitude, velocity and acceleration of the k th response respectively.
The amplitude of the k th response can be written as follows [32]: where ω e is the encounter frequency, 1 i = − is the unit imaginary number, and η _ k is the complex amplitude of the k th response: Here |η _ k | and σ k are the magnitude (time independent amplitude) and phase shift angle of the response respectively. On the other hand, the net hydrostatic force is as follows [32]: where F HSj is the hydrostatic and F Gj is the gravitational force. If the net hydrostatic force is taken to the right of the equation, the equations of motion (1) can be written as: Here M is the mass of the ship and I 5 is the moment of inertia for the pitch motion.

Body-Exact Strip Theory Approach
Strip theory was first introduced by Lewis [35], improved by many authors like Gerritsma and Beukelman [36] and Ogilvie and Tuck [37], and get its final form by Salvesen, Tuck and Faltinsen [33]. In this study, with a body-exact approach, the strip theory calculation is updated for each time step using the instantaneous wet portion of each section under the wave surface as shown in Fig. 1. In this way, the two dimensional sectional added mass and damping coefficients and wave excitation forces can be calculated depending on the time. For this reason, the recomputed parts in each time step in the coupled equations of motions (6) According to strip theory [33], the added mass and damping coefficients in the equations of motions (6) are given by the integrals over the length of the ship, which contains two dimensional sectional added mass and damping coefficients for heave. The two dimensional sectional added mass and damping coefficients are calculated by the Frank Close-Fit method over the instantaneous wetted hull surface.
The component of the hydrostatic force acting on the body in the direction j is as follows [32]: Hydrostatic force is also solved over the instantaneous wetted hull surface Here n j is the generalized unit normal outside to the hull surface [32]. The gravitational forces for for heave and pitch motions are as follows [20]: where M is the mass of the ship, g acceleration of gravity and X c is the longitudinal location of the center of gravity.
According to strip theory [33], the exciting force and moment amplitudes in the equations of motions (6) are also given by the integrals over the length of the ship, which contains sectional Froude-Kriloff and diffraction forces. Due to the body-exact strip theory approach sectional Froude-Kriloff and diffraction forces are also solved over the instantaneous wetted hull surface.
The sectional Froude-Kriloff and diffraction forces are given by [33]:  (10) where k is the wave number, ω is the circular wave frequency, χ is the wave encounter angle and C 0 is the immersed contour of the section. Here T * is a mean draft of the section [32].

Frank Close-Fit Method
In order to calculate two dimensional sectional added mass and damping coefficients, Frank [38] solved the boundary-value problem of potential theory, assuming in-compressible and inviscid, and irrotational flow. In this method, namely the Frank Close-Fit method, the surface tension effects and the nonlinear terms of the free surface condition, the kinematic boundary conditions and the Bernoulli equation are neglected. Along the two dimensional section forced into the simple harmonic motion, the hydrodynamic pressures in phase with the displacement and in phase with the velocity, which yields to sectional added mass and damping coefficients, are obtained from the velocity potential using linearized Bernoulli equation [38] & [39].

Generation of the Time Series
Following St. Denis and Pierson's [40] studies, the time series can be generated utilizing the principle of superposition. The details of the time series generation process are explained by Anıl et al. [26] and Anıl [27]. In the earth fixed coordinate system, the time series for the long crested head seas can be produced with the following equation For each component, ζ _ n is the wave amplitude, ε n is the random phase shift angle, and k n is the wave number. If δω is sufficiently small, the component wave amplitude can be written as follows [32]: Here H i is the Response Amplitude Operator (RAO) and σ i,n is the component phase shift angle of the desired response.
The time series of the derived responses can also be produced by the "real time computer experiments". The details of the time series generation process for the derived responses are explained by Anıl et al. [26] and Anıl [27].

Verification
Verification of the time series of all ship responses and the derived responses can be achieved utilizing the spectral analysis procedure which contains the discrete Fourier transform (DFT) and the smoothing algorithms, as clarified by Newland [28] and Belik [21]. The details of the spectral analysis process are explained by Anıl et al. [26] and Anıl [27].
For the verification, the spectral density functions and the corresponding RAOs derived from the "body-linear" and "body-exact" time series data with spectral analysis procedure were compared with the calculated spectra and RAOs by the linear theory in frequency domain. The results are also compared with the experimental data [29], and the ANSYS AQWA software results of Çekirdekçi [30]. The time series were obtained at 18 frequency abscissae. The Nyquist frequency, i.e. the maximum encounter frequency is 2.50 radians per second.
Within the scope of verification studies, the real time computer experiments (simulations) were carried out in two separate environmental conditions (The details of the definitions and the wave spectra for the environmental conditions are explained by Anıl et al. [39] and Anıl [40]).
The simulations were conducted primarily under the environmental condition in which physical experiments [29] were carried out. This environmental condition is represented by the JONSWAP Spectrum [42], [43], with the significant wave height H 1/3 = 2.50 m and the mean wave period T 1 = 8.81 s. This significant wave height value is the beginning of the sea state 5 interval according to the Douglas Scale or World Meteorological Organization standard sea state code (WMO Code Table 3700) [44]. In addition, simulations were performed under the environmental condition represented by the ITTC Two Parameter Spectrum [43] which is a form of the Bretschneider Spectrum [45], with the significant wave height H 1/3 = 5.00 m and the modal wave period T m = 11.74 s. This signifi-cant wave height value corresponds to the sea state 6 [44]. In this paper, only the verification result with the ITTC Two Parameter Spectrum is presented which is consistent with the verification results with the JONSWAP Spectrum (The JONSWAP Spectrum results are given in Anıl [40]).
In the spectral verification graphs: -The black curve represents the motion spectrum calculated by the linear theory in the frequency domain, -The blue curve represents the motion spectrum derived from the real time computer experiment for the "body-linear" condition, -The red curve represents the motion spectrum derived from the real time computer experiment for the "bodyexact" condition.
In the RAO and the non-dimensional acceleration graphs: -The black curve represents the RAO/non-dimensional acceleration calculated by the linear theory in the frequency domain, -The blue curve represents the RAO/non-dimensional acceleration derived from the real time computer experiment for the "body-linear" condition, -The red curve represents the RAO/non-dimensional acceleration derived from the real time computer experiment for the "body-exact" condition, -The green curve represents the RAO/non-dimensional acceleration obtained by ANSYS AQWA software by Çekirdekçi [30], -Dots represent RAO/non-dimensional acceleration obtained from the experiment [29].

Summary and Conclusions
A simulation based calculation of symmetric ship motions is performed in long crested irregular head seas, in addition with a body-exact strip theory approach. At each time step, the strip theory computation along with the Frank Close-Fit method is updated utilizing the instant wet part of each ship section.
In general, all time series obtained for the "body-linear" and "body-exact" condition make it possible to detect the extreme motions and accelerations. Therefore, the simulation results provide satisfactory data about the extreme motions and accelerations which is hard to observe with ship model tests (experiments). This is a huge benefit because real physical experimental studies are already carried out within the limits of linear theory to support numerical studies. Compared to other methods, the presented body-exact method calculates the heave motion higher, but calculates the pitch motion lower. It is evaluated that these differences in heave and pitch movements affect absolute vertical motions and accelerations. However, in general, all body-exact results are consistent with experimental data (By visual estimation, all methods are consistent with each other, or they are equally "off" from the experimental data, especially the acceleration at stern). For this reason, the method can consistently perform seakeeping performance assessments in extreme seas.
The simulation based analysis is remarkably "practical and innovative" since the preferred length of time series record has no cost and constraint except the computer power as indicated by Belik and Gökmen [22]. The simulation results are also suitable for ship-handling, active roll stabilizer fin or helicopter landing on ship simulations, The results are appropriate even for sediment accumulation simulations in double bottom ballast tanks [46]. Moreover, this approach may also be appropriate for ship form optimization studies including especially stochastic optimization.
Since the extreme responses (which is hard to catch by experiments) can be detected and analyzed from the time series record, this state-of-the-art method ensures the assessment of the seakeeping performance in extreme sea condition. As a result, it is found that this calculation method in addition with a body-exact strip theory approach can be used in the design stages of naval platforms in extreme seas.

Funding:
The research presented in the manuscript did not receive any external funding.