The motion of a ship/offshore platform at sea is governed by a coupled set of nonlinear differential equations. In general, analytical solutions for such systems do not exist and recourse is taken to time-domain simulations to obtain numerical solutions. Each simulation is not only time consuming but also captures only a single realization of the many possible responses. In a design spiral when the concept design of a ship/platform is being iteratively changed, simulating multiple realizations for each interim design is impractical. An analytical approach is preferable as it provides the answer almost instantaneously and does not suffer from the drawback of requiring multiple realizations for statistical confidence. Analytical solutions only exist for simple systems, and hence, there is a need to simplify the nonlinear coupled differential equations into a simplified one degree-of-freedom (DOF) system. While simplified methods make the problem tenable, it is important to check that the system still reflects the dynamics of the complicated system. This paper systematically describes two of the popular simplified parametric roll models in the literature: Volterra GM and improved Grim effective wave (IGEW) roll models. A correction to the existing Volterra GM model described in current literature is proposed to more accurately capture the restoring forces. The simulated roll motion from each model is compared against a corresponding simulation from a nonlinear coupled time-domain simulation tool to check its veracity. Finally, the extent to which each of the models captures the nonlinear phenomenon accurately is discussed in detail.

## Introduction

The dynamics of any rigid body are governed by the Newton's second law of motion. For a ship or an offshore platform involving multi-DOFs, the application of Newton's law results in a coupled set of nonlinear differential equations. In general, analytical solutions for such systems do not exist and hence require either linearizing to solve it analytically or performing nonlinear time-domain simulations to get a numerical solution [1].

The numerical solution of ship motions in irregular seas is time consuming and each simulation captures only a single realization of the many possible responses. In order to get a statistical confidence, multiple realizations are needed which require further time and computational resources.

In a basic design phase when the concept design of a ship/platform is being iteratively changed, simulating multiple realizations for each interim design is not only impractical but also extremely time consuming. On the contrary, an analytical approach which provides the answer almost instantaneously is more preferable. Another advantage with the analytical models is that they do not require analysis of multiple realizations for statistical confidence.

However, the analytical solutions are only possible for simple models. Thus, for a possibility of analytical methods, the equations of motion must be simplified to a model involving a few manageable parameters. The simplifying assumptions often come at the cost of rendering model incapable to capture the complete dynamics of the system. It is imperative that the designer, during the design phase, make the right compromise between the computational time involved in a simulations and simplification of the model.

The parametric roll of ships in head seas is a complex nonlinear phenomenon and many researchers have approached this problem with a combination of simplified models and numerical simulations. Initial investigators preferred to study this phenomenon only in regular waves. Paulling [2] was one of the first investigators of this phenomenon who suggested a few simplified approaches, but these methods were insufficient to study large amplitudes of motion. Neves and Rodriguez [3] proposed a coupled model with third-order restoring coefficients and showed that the second-order model was insufficient in capturing the dynamics in regular waves. The analytical expressions for the coefficients in regular waves were obtained based on the offset data of the ship. Neves and Rodríguez [4] later extended this method to calculate the new stability boundaries for the Hill's equation. However, as the model is improved from second- to third-order, both the complexity and the number of coefficients to evaluate increase tremendously. There also is a speculation [5] that such detailed geometrical data might not always be handy. Spyrou et al. [6] suggested analytical and probabilistic techniques to predict the susceptibility of a hull form to parametric roll in regular waves. However, these methods were based on the assumption of a Mathieu-type instability and were not applicable as the mathematical models become more realistic [7]. Moideen suggested a Hill's equation approach instead of a Mathieu equation for modeling the parametric roll in regular waves and developed the three-dimensional stability charts to predict the roll amplitude in regular waves. The advantage of using Hill's equation was to capture nonharmonic but periodic variation of hydrostatic stiffness in regular waves [8,9].

While the investigations into regular wave parametric roll continued, many researchers also started investigating the problem in irregular seas. One of the first ideas was provided by Grim [10] to extend the regular wave models to the irregular waves using an equivalent regular wave representation of the irregular wave profile. Umeda et al. [11] developed a simplified model for regular waves in which the roll restoring stiffness was calculated using a Grim's effective wave concept. The Grim's effective wave reduced the problem of calculating roll stiffness over an irregular wave profile to an equivalent regular wave whose crest/trough was fixed at midship. The model was quasi-statically balanced on the regular wave to match the calm water displacement. However, Hashimoto et al. [12] compared this model with the experiments to find that the two did not agree well suggesting that an improvement was needed. Bulian [5] improved upon Umeda's model by incorporating the improved Grim effective wave (IGEW) into the roll stiffness. The improved Grim effective wave allowed the equivalent regular wave to have crest/trough at any point along the length of the hull [13] and resulted in a better fit of the irregular wave profile. This model predicted a better quasi-static equilibrium trim than the original Grim effective wave model. However, even the improved model showed only a some what reasonable agreement with the experiments [5]. The SCAPE committee as a part of its investigation also documented the comparison of the Grim effective wave and IGEW models against the experimental results and suggested the use of a vulnerable criterion and direct assessment for the design against parametric roll [14]. It concluded that although the simplified 1DOF model could be used as an assessment, a complete nonlinear time-domain simulation tool would still be required to gain a thorough understanding of the phenomenon.

While the Grim's effective wave was becoming popular, an alternate simplified model was proposed by Hua et al. [15], which modeled the *GM* variation in waves using a system of Volterra transfer functions. This method was further investigated by Moideen et al. [16] and Somayajula et al. [17] in detail. The Volterra series method allows the *GM* variation to be obtained through frequency-domain transfer functions while including the effect of dynamic heave and pitch of the vessel in waves. Its advantage includes using the exact irregular wave profile instead of using an equivalent regular wave. In addition, the Volterra model also incorporates the dynamic effect of heave and pitch motion on the roll restoring moment. However, its drawback is the exclusion of the time-varying cubic restoring stiffness. A further improvement of the model to include constant cubic restoring stiffness was also investigated by Somayajula and Falzarano [18].

Continuing to investigate the problem of parametric roll, Somayajula and Falzarano [1] also developed a 6DOF nonlinear time-domain simulation tool called SIMDYN which included the nonlinear Froude–Krylov forces and nonlinear hydrostatic forces and solved for the large amplitudes of motion. SIMDYN, however, still uses linear radiation and diffraction forces calculated by three-dimensional radiation and diffraction program using the Green function approach [19–21]. Somayajula and Falzarano [22] compared the simplified model described in Ref. [18] against the nonlinear time-domain simulation tool and found that there was still a missing piece in the single-DOF models which rendered them incapable of capturing the complete dynamics of the phenomenon of parametric roll.

This paper continues upon the above investigations to systematically compare the various simplified models for simulating parametric roll of ships in head seas. The simplified models will also be compared against the coupled nonlinear time-domain simulations to ascertain the limitations of each model. Specifically, we will be looking at two different simplifying approaches which have gained much popularity in the recent times to reduce the problem into a 1DOF uncoupled roll model—Volterra series method of calculation of *GM* variation and the improved Grim's effective wave approximation. The Volterra series method of *GM* calculation used by both Hua et al. [15] and Moideen [23] included a time-varying *KG* in the calculation of *GM*. In this paper, we discuss this assumption and suggest an improvement to consider *KG* to be time invariant.

For the purpose of simulation, a standard hull form, i.e., APL China with minor modifications (C11 hull form), is chosen. The C11 hull form is known to exhibit parametric roll instability [24] and is chosen specifically due to its susceptibility. The particulars of the ship used are shown in Table 1, and its body plan is shown in Fig. 1.

Particulars | Variable | Value |
---|---|---|

Length between perpendiculars (m) | L_{pp} | 262.00 |

Breadth (m) | B | 40.00 |

Depth (m) | D | 24.45 |

Mean draft (m) | T | 12.30 |

Displacement (ton) | Δ | 74,520 |

Vertical center of gravity (m) | KG | 18.111 |

Transverse metacentric height (m) | GM_{t} | 2.075 |

Roll natural period (s) | $T\varphi $ | 25.27 |

Critical damping ratio (%) | ζ | 0.79 |

Particulars | Variable | Value |
---|---|---|

Length between perpendiculars (m) | L_{pp} | 262.00 |

Breadth (m) | B | 40.00 |

Depth (m) | D | 24.45 |

Mean draft (m) | T | 12.30 |

Displacement (ton) | Δ | 74,520 |

Vertical center of gravity (m) | KG | 18.111 |

Transverse metacentric height (m) | GM_{t} | 2.075 |

Roll natural period (s) | $T\varphi $ | 25.27 |

Critical damping ratio (%) | ζ | 0.79 |

## Roll Equation of Motion

*GM*(

*t*) can be expressed as a sum of the calm water metacentric height

*GM*

_{0}and the time-varying component

*δGM*(

*t*) as shown in the following equation:

*K*

_{3}(

*t*) is independent of time to make the problem tenable. It is common to estimate the value of

*K*

_{3}by curve fitting the calm water

*GZ*curve. However, it will be later shown that this assumption leads to significant deviation from the actual physics and needs to be used with extreme caution. Dividing Eq. (3) throughout by [

*I*

_{44}+

*A*

_{44}(

*ω*

_{0})] leads to the following equation:

The linear and quadratic damping coefficients in Eq. (5) are generally calculated from a free decay experiment. In case a free decay experiment is unavailable, they may be estimated from experimentally measured roll motion time series using a system identification approach [26,27]. However, when the experimental data are unavailable, the linear and quadratic damping coefficients can be obtained using the empirical Himeno approach [28,29]. The critical damping ratio in Table 1 is computed as the ratio of effective linearized damping to the critical damping. The linearized effective damping is calculated using the linear and quadratic coefficients for a chosen roll amplitude of 10 deg.

It is clear from Eq. (5) that in order to solve for roll motion time series, we need to estimate the *GM* variation in waves. In this paper, we will investigate two different methods of calculating the *GM* variation.

*GM* Variation by Volterra Series Approach

This section describes the method of estimating the *GM* variation from a Volterra series approach. Although this method is based on the works of Hua et al. [15] and Moideen et al. [16,23], there are significant differences which are detailed below.

The ship is separated into a number of strips along the length as shown in Fig. 2. When the ship is encountering irregular waves, at every time instant, the local draft at each of these sections is different due the effects of the wave, instantaneous heave, and instantaneous pitch of the ship. The local breadth *B*(*x*, *T* + *z*) and moment of sectional underwater area about the keel *M*(*x*, *T* + *z*) at any section *x* and local draft *T* + *z* can be expanded into a Taylor series as shown in Eqs. (10) and (11), respectively

*GM*for any free floating structure is given by the following equation:

*BM*and

*KB*can be expressed as integrals along length of the ship in terms of the relative waterline

*r*(

*t*,

*x*) as shown in the following equations:

*GM*in waves. Hua et al. [15] also suggested a time-varying

*KG*across the length of the ship. However,

*KG*of the vessel only depends on the mass distribution of the vessel and is immune to the motion of the ship. Hence in this paper, we shall consider

*KG*to be time invariant and accordingly modify our Volterra series model

*GM*variation can be expressed as a sum of various orders of

*GM*variation as shown in the following equation:

The plot of *G*_{1}(*x*) and *G*_{2}(*x*) along the length for the C11 hull form is shown in Fig. 3.

*a*,

_{m}*k*,

_{m}*ω*, and

_{m}*ε*represent the amplitude, wave number, encounter frequency, and phase of the

_{m}*m*th regular wave component:

*ξ*

_{3}(

*ω*) and

_{m}*ξ*

_{5}(

*ω*) are the heave and pitch response amplitude operators (RAOs) at encounter frequency

_{m}*ω*, the corresponding relative wave elevation incorporating the dynamic heave and pitch of the vessel can be expressed as follows:

_{m}*f*(

*ω*) is the first-order transfer function given by

*GM*variation is given by the following equation:

*g*

_{1}and

*g*

_{2}are the second-order transfer functions given by

Figure 4 shows a plot of the first-order transfer function *f*(*ω*) as defined in Eq. (30). Similarly, the second-order transfer functions *g*_{1}(*ω*, *ω*) and *g*_{2}(*ω*, *ω*) are shown in Figs. 5 and 6, respectively.

In order to verify the veracity of the Volterra series formulation, the first-order *GM* variation is compared against the actual *GM* variation in waves. The actual *GM* variation is obtained by calculating the instantaneous relative waterline which is then used to evaluate the instantaneous *GM* at every time step. Figure 7 shows the comparison between the first-order *GM* variation and the actual simulated *GM* in the same input wave realization from a Bretschneider spectrum with a significant wave height *H _{s}* = 3 m and modal period

*T*= 14.1 s. The comparison shows a decent agreement which suggests that the

_{p}*GM*variation obtained by Eq. (29) can be used as a reasonable approximation for

*δGM*(

*t*) in Eq. (5) to simulate the parametric roll motion of the ship.

## Grim's Effective Wave Approach

This section describes a second approach to obtain the hydrostatic variation in waves. In this approach, instead of using the irregular wave profile across the ship, the problem is reduced to a ship in a regular wave. The original idea of reducing the irregular wave profile across a ship into an equivalent regular wave was proposed by Grim in 1961 [10].

*η*is shown in Eq. (37). Note that the Grim's effective wave has a restriction that the crest/trough of the wave be located at the origin. Assuming that the origin is at midship, the equivalent regular wave is always assumed to have the crest/trough at midship

_{g}*a*(

*t*),

*η*(

_{c}*t*), and

*η*(

_{s}*t*) can be found by solving the least-square fit problem described in Eq. (40) and are shown in Eqs. (41), (42), and (43), respectively

The square of the transfer functions *f _{a}*(

*k*),

_{m}*f*(

_{c}*k*), and

_{m}*f*(

_{s}*k*) are shown in Fig. 8. Figure 9 shows the quality of fit of IGEW and Grim effective wave to the actual wave elevation in space at a particular time instant. It can be seen that the flexibility of IGEW to move the crest along the length clearly results in a better fit than the original Grim effective wave.

_{m}With the improved Grim effective wave, the problem of calculating the hydrostatic stiffness of a ship on an irregular wave profile is reduced to that on an equivalent regular wave. Within this approach of using a regular wave approximation, there are again two different techniques to model the hydrostatic stiffness.

*GM* Variation Using IGEW.

The first method is to evaluate the instantaneous *GM* by imposing the regular wave profile over the ship. Bulian [31] suggests allowing the ship to trim freely to a quasi-static equilibrium before calculating the *GM*. However, in this paper, we will consider a model which uses IGEW and the dynamic heave and pitch of the vessel to evaluate the instantaneous waterline and *GM*. The inclusion of dynamic heave and pitch is an improvement over quasi-static trim approach and is implemented to make a fair comparison with the Volterra model.

Figure 10 shows the comparison of the *GM* calculated by using IGEW with the actual simulated *GM* for the same input wave realized from a Bretschneider spectrum with significant wave height *H _{s}* = 3 m and modal period

*T*= 14.1 s. In contrast to the Volterra series method, the power spectral density of

_{p}*GM*variation obtained using IGEW has a shifted peak as compared to the simulated

*GM*variation.

*GZ* Variation Using IGEW.

*GM*varies with time, the hydrostatic stiffness is modeled using a time-varying

*GZ*in regular wave. The roll equation of motion corresponding to this model is given by the following equation:

At every time step, the instantaneous wave profile across the ship is approximated by *η _{g}*(

*t*,

*x*) as given by Eq. (38). The instantaneous

*GZ*at each time step is interpolated from a lookup table. The lookup table for

*GZ*values is calculated prior to the simulation for various regular wave amplitudes, crest location along the length of the ship, and different roll angles [5,32]. These

*GZ*values for various cases are usually obtained using any standard hydrostatic calculation software. A sample

*GZ*curve for a 4 m high regular wave for various crest positions is shown in Fig. 11(a).

It is important to note that this model does not include the effect of dynamic heave and pitch on the change in instantaneous waterline. While calculating the *GZ* lookup table, in each of the cases, the model is allowed to trim freely. Thus, a lookup table similar to a *GZ* can also be calculated for quasi-static trim angle, which can be used to compute the equilibrium trim angle time series. The quasi-static trim values for various heel angles and crest positions for a 4 m high regular wave are shown in Fig. 11(b). Figure 11(c) shows the comparison of dynamic pitch with the quasi-statically calculated trim time series.

## Comparison of Simulated Roll Time Histories

So far, we have only discussed the differences in all three roll models in terms of how they model the roll restoring moment. However, in order to understand whether a model captures the physics accurately, it is important to compare the simulated roll time series with either experiments or coupled nonlinear time-domain simulations. Due to the lack of extensive experimental data in irregular waves, in this paper, we will compare the simplified models against a nonlinear time-domain simulation program, SIMDYN, developed in house at Texas A&M University [1].

This time-domain simulation tool has been validated against the regular wave experiments performed by Silva et al. [33]. Silva et al. documented the observed parametric roll in regular wave tests performed at Canal de Experiências Hidrodinâmicas de El Pardo (CEHIPAR), Spain as a part of the HYDROLAB III project. Twelve regular wave tests were performed which spanned three wave heights of 6, 8, and 10 m and four wavelengths corresponding to *λ*/*L* of 0.8, 1.0, 1.2, and 1.4. As the roll decay tests were unavailable, the simulations used a standard linear and quadratic damping model proposed by Himeno [28]. Each simulation was started from a zero initial roll condition, and the steady-state solution was measured after visually examining that the transient had decayed. The comparisons between simulations and experiments for all 12 cases are shown in Fig. 12.

It can be seen that the SIMDYN simulations agree quite well with the experiments. Both the roll amplitudes and the disappearance of parametric solution with higher *λ*/*L* are captured well. Based on these comparisons, SIMDYN is chosen as the benchmark to compare the simplified models.

Figure 13 shows the comparison of the roll response from the three models described against the results from the nonlinear time-domain simulation. The same input wave realized from a Bretschneider spectrum with significant wave height *H _{s}* = 3 m and modal period

*T*= 14.1 s has been used for all models to be consistent (Fig. 13(a)). Unlike SIMDYN, the 1DOF models require an initial perturbation for the instability to manifest. Thus, each simulation (SIMDYN and 1DOF models) was started from an initial condition of 20 deg roll angle for consistency. While many realizations were investigated as a part of this study, it was found that all of the realizations demonstrated the same qualitative results as shown in Fig. 13. For brevity, only one of these realizations shown in Fig. 13 will be discussed here. It can be seen that all the three models are far off from the nonlinear time-domain simulations. The

_{p}*GZ*variation model performs a bit better than the other two models in terms of the maximum roll angle. This clearly shows that the simplified models are unable to capture certain key aspects of the phenomenon of parametric rolling.

Although the heave and pitch of the vessel significantly affect the roll motion through the variation of instantaneous hydrostatic restoring moment, the effect of roll on pitch and heave is negligible. Even in cases of severe parametric roll, the pitch and heave motions of the ship are nearly linear. This has been demonstrated from both the experiments [12] and simulations [1,34]. This indicates that the major difference between the coupled nonlinear simulations and the simplified models is in the handling of nonlinear stiffness and Froude–Krylov forcing.

SIMDYN uses a nonlinear Froude–Krylov force vector which is obtained by integrating the Wheeler stretched pressure [35] around the hull. It also uses a nonlinear hydrostatic force vector accounting for the instantaneous position of the hull with respect to the incident waterline. These are the two main factors which are not completely accounted for in the simplified models.

To ascertain the dominant factor between them, the SIMDYN simulations in regular waves were performed with and without the nonlinear Froude–Krylov force options. In both cases, the nonlinear hydrostatics, large angles of rotation (Euler angles), and convolution integral were included. Two cases with different wave heights of 2 m and 3 m but the same period of 14.1 s were investigated. The results are shown in Figs. 14(a) and 14(b), respectively.

It is interesting to see that including the nonlinear Froude–Krylov force option in the simulation does not change the final steady-state amplitude of roll angle for both wave heights. Including the nonlinear Froude–Krylov moment only changes the transient behavior and does not affect the steady-state roll amplitude reached by the system. These results suggest that the nonlinear hydrostatics (including the dynamic effect of heave and pitch) is the key to enhancing the simplified models to capture the physics of the problem accurately.

Comparison of the three simplified models with the nonlinear time-domain simulations shows that the simplified models overpredict the roll motion. Both the Volterra and IGEW *GM* models account for only the time-varying linear stiffness while assuming a constant cubic restoring coefficient. It can be seen that the *GM* variation predicted by the Volterra model compares better with the actual *GM* variation than the IGEW *GM* model. However, its inability to capture the final roll angle accurately as compared with the IGEW *GZ* model suggests that the nonlinear restoring component in waves plays a significant role. This suggests that if a Volterra type of model can be implemented for the higher-order (polynomial) coefficients, then it might lead to a more accurate modeling and predictions closer to the time-domain simulations.

## Conclusions and Future Works

Although parametric roll of ships has been known to researchers for a long time, it is still an active area of investigation. Various simplified models have been proposed by researchers which include varying degrees of complexity. In this paper, an effort has been made to study in detail three of the more popular approaches and compare them against the coupled nonlinear time-domain simulations to assess the strengths of the different models. The first method is based on modeling the *GM* variation in waves using a Volterra series approach. Since the *GM* variation is reduced to a set of transfer functions, this method is quite suitable for analytical and probabilistic techniques. The second and third models use the Grim's effective wave as a basis for reducing the hydrostatic restoring moment in irregular head seas to that in an equivalent regular wave problem.

All three models were compared against the coupled nonlinear time-domain simulations. It was found that between the two *GM* variation models, the Volterra model has a better accuracy in predicting the actual *GM* variation. However, as the current Volterra model does not include time-varying higher-order stiffness coefficients, it compares worse than the IGEW *GZ* model.

All the three models investigated in this work are found to overpredict the roll response when compared with the nonlinear time-domain simulations. While the overprediction of roll response of IGEW *GZ* model has been previously reported [13], the corresponding comparisons of other models have not been investigated to the authors' knowledge. However, the comparative analysis shows that the *GM* variation is better predicted by the Volterra series method as compared to the IGEW *GM* model. This hints at the possibility of achieving a better model if the Volterra model can be extended to higher-order stiffness coefficients. This approach will be investigated further in the future.

## Acknowledgment

The authors would like to thank Dr. Frans van Walree of MARIN for making available to us the hull form description for our analysis. The authors thank Dr. Sergio Ribeiro e Silva for making available the parametric roll experimental data from HYDROLAB III project for comparisons. This work has been funded by the Office of Naval Research (ONR)—ONR Grant N000-14-16-1-2281. The authors thank Dr. Paul Hess for facilitating the funding for this work.

## Nomenclature

*A*_{44}=frequency-dependent roll added moment of inertia

*B*_{1}=linear viscous roll damping coefficient

*B*_{2}=quadratic viscous roll damping coefficient

*B*_{44}=frequency-dependent roll radiation damping

*BM*=metacentric radius—distance between center of buoyancy and metacenter

*g*=acceleration due to gravity

*GM*=metacentric height—distance between center of gravity and metacenter

*GZ*=roll restoring arm

*GM*_{0}=calm water metacentric height

*I*_{44}=roll mass moment of inertia

*K*_{3}=cubic coefficient of

*GZ*approximation*KB*=distance of vertical center of buoyancy from keel

*r*(*t*,*x*) =relative wave elevation at station

*x*at time*t**T*=draft at calm water

- ∇ =
displacement of the vessel

*ρ*=density of water

- $\varphi $ =
roll angle

- $\varphi \u02d9$ =
roll velocity

- $\varphi \xa8$ =
roll acceleration

*ω*_{0}=roll natural frequency of the ship