11institutetext: Keiji Katsumura 22institutetext: Taiga Kono 33institutetext: Yuuki Maruyama 44institutetext: Masahiro Sakai 55institutetext: Atsuo Maki66institutetext: Department of Naval Architecture and Ocean Engineering,
Graduate School of Engineering, Osaka University,
2-1 Yamadaoka, Suita, Osaka, Japan
Leo Dostal
77institutetext: Institute of Mechanics and Ocean Engineering, Hamburg University of Technology, 21043 Hamburg, Germany

Validation of Theoretical Estimation Methods and Maximum Value Distribution Calculation for Parametric Roll Amplitude in Long-Crested Irregular Waves

Keiji Katsumura    Leo Dostal    Taiga Kono
Yuuki Maruyama
   Masahiro Sakai    Atsuo Maki
(Received: date / Accepted: date)
Abstract

Parametric rolling is a parametric excitation phenomenon caused by GM variation in waves. There are a lot of studies of the estimation the conditions, the occurrence, and the amplitude of parametric rolling. On the other hand, there are relatively few cases in which theoretical methods for estimating parametric roll amplitudes in irregular waves have been validated in tank tests. The primary objective of this study is to validate theoretical estimation methods for the parametric roll amplitude in irregular waves and improve their accuracy. First, the probability density functions (PDF) of the parametric roll amplitude obtained from the model ship motion experiment in irregular waves are compared with that obtained from theoretical estimation methods. Second, the method to improve the accuracy of estimation of the roll restoring variation in irregular waves is suggested. Third, the method to estimate the distribution of the maximum amplitude of parametric rolling in irregular waves. As a result, the PDFs of the roll amplitude obtained from the experiments differ from the results of theoretical estimation. After that, by correcting GM variation, the results of theoretical estimation are closer to the experimental results. Moreover, by the theoretical estimation method using the moment equation, the qualitative estimation for the PDF of the maximum roll amplitude is succeeded.

Keywords:
Model test Probability density function Grim’s effective wave concept Moment equation Extreme value GM Variation
journal: Journal of Marine Science and Technology

1 Introduction

Parametric rolling is a parametric excitation phenomenon caused by GM variation in waves. It is particularly likely to occur on container ships, whose transom stern and bow flare cause significant changes in the secondary moment of the waterline near the stern of the ship. For example, the accident of a C11 class container ship caused by parametric rolling in 1998 is well known France2003 . Moreover, the accident of a pure car and truck carrier (PCTC) was also occurred in recent years Rosen2012experience . In order to prevent such dangerous phenomenon, parametric rolling, it is necessary to estimate its conditions, occurrence, and amplitude.

1.1 Related Research

There is a long research history on theoretical estimation methods for the conditions, the occurrence, and the amplitude of parametric rolling. The parametric rolling in regular seas has been theoretically explored by Kerwin Kerwin1955 , Zavodney et al. Zavodney1989 , Francescutto Francescutto2001 , Bulian Bulian2004approximate , Spyrou Spyrou2005paramet , Umeda et al. Umeda2004nonlinear , Maki et al. Maki2011parametric , and Sakai et al Sakai2018 . Furthermore, in particular, since the 1980s, there have been a lot of studies on theoretical estimation methods for parametric rolling in irregular waves. In order to study parametric rolling in irregular waves, systems excited by colored noise must be treated. One of methods to treat such systems theoretically is through a probabilistic approach. Spyrou Spyrou2005paramet , Dostal Dostal2012non and Maki Maki2023 discussed the occurrence of parametric rolling in terms of the destabilization of the origin due to the roll restoring variation in irregular waves. For the estimation of parametric roll amplitudes, studies using stochastic averaging methods have been conducted. Using the stochastic averaging method, Roberts Roberts1982 derived the stochastic differential equation for a phase and an amplitude. In Roberts’ results, the probability density function (hereinafter referred to as PDF) of the roll amplitude reflected the characteristics of the damping component and those of the cubic restoring component. Subsequently, Roberts et al. Roberts1986 ; Roberts2000 proposed an energy based methodology and attempted to reflect the restoring component. Dostal et al. Dostal2011 proposed an energy-based stochastic averaging method using the Hamiltonian. Furthermore, Maruyama et al. Maruyama2022 developed Dostal’s method. They compared the results of Roberts Roberts1982 and Dostal Dostal2011 ’s methods with those of Monte Carlo Simulation (hereinafter referred to as MCS). As a result, they reported a discrepancy in the tail of the PDF and proposed a method called the Simulation-Based Stochastic Averaging Method to solve the discrepancy. In addition, Maruyama et al. Maruyama2022_moment proposed a method for estimating parametric rolling using the moment equation Bover1978 (hereinafter referred to as the moment method). This method also allows to obtain some quantitative agreement with the result of MCS for the PDF of the parametric roll amplitude. On the other hand, there are relatively few cases in which these analytical methods for estimating roll motion have been validated in tank tests.

Moreover, estimation of the roll restoring variation is important in the theoretical estimation of parametric rolling. In the studies of the theoretical estimation of parametric rollong in irregular waves by Maruyama et al. Maruyama2022 Maruyama2022_moment Maruyama2023_momentamp mentioned above, the roll restoring variation was estimated by considering only the wave component based on the Froude-Krylov assumption under quasi-statically balancing heave and pitch in waves and introducing Grim’s effective wave concept Grim1961 . However, Hashimoto et al. Hashimoto2004 suggest that the estimation method of the roll restoring variation considering only the wave component based on the Froude-Krylov assumption under quasi-statically balancing heave and pitch in waves might overestimate the risk of parametric rolling due to the large roll restoring variation during wave passage through the captive model test. Yu et al. Yu2019 compare the degree of this overestimation using five models that estimate nonlinear restoring forces and Froude-Krylov forces. Furthermore, Hashimoto et al. Hashimoto2006 measured the roll restoring variation in a captive model test in irregular waves and compared it with the results of the calculation of the roll restoring variation using Grim’s effective wave. The results suggested that the accuracy of the roll restoring variation using Grim’s effective wave remained a problem.

1.2 Object and Scope

The primary objective of this study is to validate theoretical estimation methods for the parametric roll amplitude in irregular waves and improve their accuracy. The contributions of this study are as follows:

  1. 1.

    Validation of theoretical estimation methods for parametric roll amplitude by comparison with corresponding experimental result,

  2. 2.

    Suggestion for the method to improve the accuracy of estimation of the roll restoring variation in irregular waves, and

  3. 3.

    Suggestion for the method to estimate the distribution of the maximum amplitude of parametric rolling in irregular waves.

In this study, we first conducted the model ship motion experiment in order to obtain the PDF of parametric roll amplitude in long-crested irregular waves in the towing tank of Osaka University. After that, the accuracy of the theoretical methods is validated by comparing the experimental results with the PDFs of the parametric roll amplitude calculated using three theoretical calculation methods; Roberts’ stochastic averaging method Roberts1982 , the energy-based stochastic averaging method Dostal2011 (hereinafter referred to as ESAM), and the moment method Maruyama2022_moment ; Maruyama2023_momentamp . Furthermore, we attempted to improve the accuracy of the theoretical estimation method by using the equation of motion reflecting the correction of the roll restoring variation based on the captive model test by Kono et al.

Extreme value theory is also widely used in risk assessment in various engineering fields. It would be significant to introduce the concept of extreme value theory to the risk assessment of parametric rolling motion. On the other hand, accurate estimation of the distribution of maximum values requires accurate estimation of the PDF of the parametric rolling amplitude down to the tail. Therefore, we pay particular attention to the behavior of the tail section of the PDF of the parametric roll amplitude in order to estimate extreme values. Finally, we propose a method for calculating the distribution of the maximum roll amplitude based on the theoretical method.

2 Subject ship

The subject ship is a C11 class post-Panamax container ship. In this study, a geometric similarity model with a scale of 1/100 is used for model tests at the towing tank of Osaka University. The main features of the model are shown in  Table 1.

Table 1: Principal Particulars of C11 at model and full scale.
Items Value
Scale Model Full
Lpp(m)subscript𝐿ppmL_{\mathrm{pp}}\,\left({\mathrm{m}}\right)italic_L start_POSTSUBSCRIPT roman_pp end_POSTSUBSCRIPT ( roman_m ) 2.62 262
B(m)𝐵mB\,\left({\mathrm{m}}\right)italic_B ( roman_m ) 0.4 40
D(m)𝐷mD\,\left({\mathrm{m}}\right)italic_D ( roman_m ) 0.2445 24.45
d(m)𝑑md\,\left({\mathrm{m}}\right)italic_d ( roman_m ) 0.115 11.5
W(kg)𝑊kgW\,\left({\mathrm{kg}}\right)italic_W ( roman_kg ) 67.247 6.7247×1076.7247superscript1076.7247\times 10^{7}6.7247 × 10 start_POSTSUPERSCRIPT 7 end_POSTSUPERSCRIPT
Cbsubscript𝐶bC_{\mathrm{b}}italic_C start_POSTSUBSCRIPT roman_b end_POSTSUBSCRIPT 0.56 0.56
GM(m)GMm\mathrm{GM}\,\left({\mathrm{m}}\right)roman_GM ( roman_m ) 0.019299 1.9299
Tϕ(s)subscript𝑇italic-ϕsT_{\phi}\,\left({\mathrm{s}}\right)italic_T start_POSTSUBSCRIPT italic_ϕ end_POSTSUBSCRIPT ( roman_s ) 2.44 24.4

The extinction coefficients used in the calculations were calculated from the results of free-rolling tests using the model in a water tank.

3 Towing Tank Experiment

By conducting model ship motion experiments in long-crested irregular waves and obtaining the PDF of the roll amplitude, we validate the theoretical calculation method as described in Section 1. This section describes the experimental method used in this study.

The bow and stern of the model ship is connected to the towing dolly with a rubber strap to gently restrain the model ship. The tension of the rubber strap should be adjusted appropriately so as not to constrain the surge or roll motion too much. The stern rubber strap was not used in most of the experiments in head-wave conditions. The state of the model ship is shown in Fig. 1. The long-crested irregular head waves are generated using ITTC spectrum. The angle of occasional rolling of the model ship is measured by a fiber-optic gyro sensor mounted at the center of gravity of the model ship. Wave heights of irregular waves are measured by a wave height meter installed at the front of the model ship. , The experimental conditions for irregular waves are shown in Table 2.

Table 2: Wave condition at full scale
Items Value
H1/3(m)subscript𝐻13mH_{1/3}~{}\left({\mathrm{m}}\right)italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT ( roman_m ) 5.0 7.0
T01(s)subscript𝑇01sT_{01}~{}\left({\mathrm{s}}\right)italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT ( roman_s ) 10.0 10.0
Num. of realization 24 24

The producing time for each irregular wave is 300s300s300~{}\mathrm{s}300 roman_s. Here, the data between 270s270s270~{}\mathrm{s}270 roman_s from 70s70s70~{}\mathrm{s}70 roman_s to 340s340s340~{}\mathrm{s}340 roman_s after the start of producing waves were considered for the analysis.

Refer to caption
Figure 1: Schematic view of the experiment
Refer to caption
Figure 2: The situation of the experiment

4 Experimental Results and Discussion

4.1 Analysis methods and theoretical estimation methods for roll amplitudes

An example of a time series of the roll angle obtained in the experiment is shown in Fig. 3.

Refer to caption
Figure 3: Time history of roll angle, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s

In this study, the following two methods of analyzing the roll amplitude are used.

4.1.1 Experimental Zero Crossing method

In the time series data of the roll angle obtained in the experiment, The time series data of roll angle is obtained by the experiment. The maximum value between every zero up crossing and zero down crossing and the minimum value between every zero down crossing and zero up crossing are recorded. The absolutes of them are the roll amplitudes. The PDF is calculated from the distribution of the obtained roll amplitudes.

4.1.2 Experimental Envelope method

The Hilbert transform of the time series of roll angles obtained in the experiment is used to create an envelope. Then, the absolute value of the envelope of the roll angle at each sampling point is recorded as the amplitude. The PDF is calculated from the distribution of the values of roll angles.

4.2 Theoretical estimation methods

Theoretical estimation of the amplitude of parametric rolling in irregular waves requires the estimation of the corresponding the roll restoring variation. Therefore, in this study, we introduce Grim’s effective wave theory Grim1961 for irregular waves. Moreover, considering only the wave component based on the Froude-Krylov assumption under quasi-statically balancing heave and pitch in waves, we obtain an equation relating the displacement of the wave at the center of the hull to the amount of GM variation Maruyama2022 . In this way, we estimate GM variation in irregular waves. The specific procedure of this transformation method is described in detail in Appendix Appendix A. Based on the roll restoring variation estimated in this way, we estimate the amplitude of parametric rolling in irregular waves using three theoretical methods, Roberts’ stochastic averaging method Roberts1982 , ESAM Dostal2011 , and the moment method Maruyama2023_momentamp .

4.3 Experimental results

The results of comparing the PDFs of the roll amplitudes obtained from experiments and various theoretical calculation methods are shown in Figs. 4 and 5, and Figs. 6 and 7. Here, the black dots in the figure are the PDF calculated by the experimental zero crossing method described in  Section 4.1.1 and the red dots are the PDF calculated by the experimental envelope method described in  Section 4.1.2. The blue line is the PDF obtained by Roberts’ stochastic averaging method  Roberts1982 , the yellow line is the PDF obtained by ESAM Dostal2011 , and the purple line is the PDF obtained by the moment method Maruyama2023_momentamp . The green dots are the PDF obtained by MCS. For MCS, the initial value of the roll angle was set to 5deg5deg5~{}\mathrm{deg}5 roman_deg, and the initial value of the roll velocity was set to 0deg/s0degs0~{}\mathrm{deg/s}0 roman_deg / roman_s. The time for one trial of simulation was 3600 s. The number of trials was 104superscript10410^{4}10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT times.

Refer to caption
Figure 4: PDF of roll amplitude, H1/3=5.0msubscript𝐻135.0mH_{1/3}=5.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 5.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, linear scale
Refer to caption
Figure 5: PDF of roll amplitude, H1/3=5.0msubscript𝐻135.0mH_{1/3}=5.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 5.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, logarithmic scale
Refer to caption
Figure 6: PDF of roll amplitude, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, linear scale
Refer to caption
Figure 7: PDF of roll amplitude, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, logarithmic scale

4.4 Consideration

The PDFs of the roll amplitudes obtained from the experiments differ quantitatively and qualitatively from the PDFs obtained from the simulation and from each theoretical prediction method. Around 0deg0deg0~{}\mathrm{deg}0 roman_deg, the simulation and theoretical results show that the PDF converges to 00, while the experimental results are non-zero. The experimental values have a local maximum from 10deg10deg10~{}\mathrm{deg}10 roman_deg to 20deg20deg20~{}\mathrm{deg}20 roman_deg, and then they are smaller than the theoretical values. In the experiments, an intermittent rolling behaviour was observed, whereas pure parametric rolling behavior resulted from the theoretical methods and the MCS of the equations of motion. In Figs. 6 and 7, the form of the PDF from the experiments obtained here, which takes non-zero values around 0deg0deg0~{}\mathrm{deg}0 roman_deg and has one more hump, is not consistent with the two forms of the theoretical PDFs presented in previous studies Maruyama2023_momentamp ; Maki2019 . We currently believe that the reason for such PDFs is the change of yaw angle during the experiment. Further investigation is necessary in the future.

Next, we compare the three PDFs obtained by the theoretical method with that obtained by MCS. In Fig. 6, the PDFs obtained by Roberts’ stochastic averaging method and ESAM agree well with the result obtained by MCS in the amplitude range from 0deg0deg0~{}\mathrm{deg}0 roman_deg to 10deg10deg10~{}\mathrm{deg}10 roman_deg. However, in the log scale in Fig. 7, the PDF obtained by the moment method agrees better with the result of MCS in the tail (large amplitude part) of the PDF. As mentioned in Section 1, the estimation of the maximum roll amplitude requires an accurate estimation of the PDF of the roll amplitude down to the tail part. Therefore, when consider the PDF of the maximum roll amplitude in Section 6 using theoretical methods, the results obtained by the moment method, which are successful in quantitatively estimating the tail, are considered more suitable for use in Eq. 5.

5 Correction of GM variation

As shown in Figs. 6 and 7, there was a discrepancy between the experimental PDF and the PDF obtained by MCS and each theoretical calculation. We attempt to reduce this discrepancy in order to improve the estimation accuracy of the PDF of the parametric roll amplitude. In this study, GM variation is calculated by the equation relating wave displacement and GM variation. We consider only the wave component based on the Froude-Krylov assumption under quasi-statically balancing heave and pitch in waves, and introduce Grim’s effective wave theory Grim1961 for irregular waves as explained in Section 4.2 Maruyama2022 . The specific calculation procedure is described in Appendix Appendix A. However, the results of captive model tests in regular following waves by Hashimoto et al. Hashimoto2004 suggests that the estimation method of the roll restoring variation considering only the wave component based on the Froude-Krylov assumption under quasi-statically balancing heave and pitch in waves might overestimate the risk of parametric rolling due to the large roll restoring variation during wave passage. Furthermore, Hashimoto et al. Hashimoto2006 measured the roll restoring variation in a captive model test in irregular waves and compared it with the results of the calculation of the roll restoring variation using Grim’s effective wave. The results suggest that the accuracy of the roll restoring variation using Grim’s effective wave remain a problem. Therefore, we consider the accuracy of the estimation of the roll restoring variation to be one of the reasons for the discrepancy between the experimental and theoretical results for the PDF of the parametric roll amplitude in irregular waves. Hence, captive model tests were conducted to evaluate the estimation accuracy of the theoretical equation for GM variation.

5.1 Captive model test

The experimental method for the captive model testing is explained in the following. The model ship was free to heave and pitch, and was attached to the towing vehicle via a quarter force gauge. The fore-and-aft force (FXsubscript𝐹𝑋F_{X}italic_F start_POSTSUBSCRIPT italic_X end_POSTSUBSCRIPT), lateral force (FYsubscript𝐹𝑌F_{Y}italic_F start_POSTSUBSCRIPT italic_Y end_POSTSUBSCRIPT), turning moment (N𝑁Nitalic_N), lateral tilting moment (K𝐾Kitalic_K), heaving, and pitching on the hull were measured. In addition, irregular waves were generated by a flap wave generator at the end of the tank, and the generated irregular wave forms were measured using a capacitance type wave height meter installed between the wave generator and the model ship. Three wave conditions were used in this experiment, see Table 3. For each wave condition, the number of realizations is 10101010.

Table 3: Wave condition at full scale on the captive model test
Items Value
H1/3(m)subscript𝐻13mH_{1/3}~{}\left({\mathrm{m}}\right)italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT ( roman_m ) 5.0 5.0 7.0
T01(s)subscript𝑇01sT_{01}~{}\left({\mathrm{s}}\right)italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT ( roman_s ) 10.0 12.37 10.0
Num. of realizations 10 10 10

In the experiment of this study, the fixed roll angles ϕitalic-ϕ\phiitalic_ϕ were 0, 5, 10, and 15 degrees. The GM variation ΔGMΔGM\Delta\mathrm{GM}roman_Δ roman_GM was obtained from the difference between the roll moment at ϕ=0degitalic-ϕ0deg\phi=0~{}\mathrm{deg}italic_ϕ = 0 roman_deg.

Refer to caption
Figure 8: Condition during captive model test

5.2 Experimental result and Method of correcting GM variation

The theoretical PDFPDF\mathrm{PDF}roman_PDF is calculated based on ΔGM(t)ΔGM𝑡\Delta\mathrm{GM}(t)roman_Δ roman_GM ( italic_t ) created in Appendix Appendix A. The experimental PDFPDF\mathrm{PDF}roman_PDF and theoretical PDFs are shown in Figs. 9 and 10. The black line is the experimental value and the blue line is the theoretical value based on the relation Eq. A.1. By comparison, the theoretical PDF has a wider range of ΔGMΔGM\Delta\mathrm{GM}roman_Δ roman_GM than the experimental PDF. Therefore, it is found that the conventional theoretical PDFs are overestimated compared to the experimental PDFs. Now, assuming that ΔGMΔGM\Delta\mathrm{GM}roman_Δ roman_GM is normally distributed, we correct Eq. A.1 as Eq. 1 using the ratio of the standard deviation σGMobssubscriptsuperscript𝜎obsGM\sigma^{\mathrm{obs}}_{\mathrm{GM}}italic_σ start_POSTSUPERSCRIPT roman_obs end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT of the experimental values to the standard deviation σGMthrsubscriptsuperscript𝜎thrGM\sigma^{\mathrm{thr}}_{\mathrm{GM}}italic_σ start_POSTSUPERSCRIPT roman_thr end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT of the theoretical values. In Figs. 9 and 10, the red line is the theoretical value based on Eq. 1.

ΔGM(ζmid)=σGMobsσGMthrk=06CkζmidkΔGMsubscript𝜁midsubscriptsuperscript𝜎obsGMsubscriptsuperscript𝜎thrGMsuperscriptsubscript𝑘06subscriptCksuperscriptsubscript𝜁midk\Delta\mathrm{GM}(\zeta_{\mathrm{mid}})=\frac{\sigma^{\mathrm{obs}}_{\mathrm{% GM}}}{\sigma^{\mathrm{thr}}_{\mathrm{GM}}}\sum_{k=0}^{6}\mathrm{C}_{\mathrm{k}% }\zeta_{\mathrm{mid}}^{\mathrm{k}}roman_Δ roman_GM ( italic_ζ start_POSTSUBSCRIPT roman_mid end_POSTSUBSCRIPT ) = divide start_ARG italic_σ start_POSTSUPERSCRIPT roman_obs end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT end_ARG start_ARG italic_σ start_POSTSUPERSCRIPT roman_thr end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT end_ARG ∑ start_POSTSUBSCRIPT italic_k = 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 6 end_POSTSUPERSCRIPT roman_C start_POSTSUBSCRIPT roman_k end_POSTSUBSCRIPT italic_ζ start_POSTSUBSCRIPT roman_mid end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_k end_POSTSUPERSCRIPT (1)
Refer to caption
Figure 9: PDFPDF\mathrm{PDF}roman_PDF of ΔGMΔGM\Delta\mathrm{GM}~{}roman_Δ roman_GMwithH1/3=0.07msubscriptH130.07m~{}\mathrm{H}_{1/3}=0.07~{}\mathrm{m}roman_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 0.07 roman_m,T01=1.0ssubscriptT011.0s~{}\mathrm{T}_{01}=1.0~{}\mathrm{s}roman_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 1.0 roman_s, linear scale
Refer to caption
Figure 10: PDFPDF\mathrm{PDF}roman_PDF of ΔGMΔGM\Delta\mathrm{GM}~{}roman_Δ roman_GMwithH1/3=0.07msubscriptH130.07m~{}\mathrm{H}_{1/3}=0.07~{}\mathrm{m}roman_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 0.07 roman_m,T01=1.0ssubscriptT011.0s~{}\mathrm{T}_{01}=1.0~{}\mathrm{s}roman_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 1.0 roman_s, logarithmic scale)

5.3 Calculation results of parametric rolling amplitude

The conventional equation of motion for one-degree-of-freedom roll motion is shown in Eq. 2, where ϕitalic-ϕ\phiitalic_ϕ is the roll angle. b1subscript𝑏1b_{1}italic_b start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT, b2subscript𝑏2b_{2}italic_b start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT, and b3subscript𝑏3b_{3}italic_b start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT are the linear, quadratic, and cubic damping coefficients divided by Ixxsubscript𝐼𝑥𝑥I_{xx}italic_I start_POSTSUBSCRIPT italic_x italic_x end_POSTSUBSCRIPT, where Ixxsubscript𝐼𝑥𝑥I_{xx}italic_I start_POSTSUBSCRIPT italic_x italic_x end_POSTSUBSCRIPT is the moment of inertia in roll including the corresponding added moment of inertia. The coefficients of the polynomial approximation of the GZ curve are ai(i=1,3,5,7,and9)subscript𝑎𝑖𝑖1357and9a_{i}~{}(i=1,~{}3,~{}5,~{}7,~{}\mathrm{and}~{}9)italic_a start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_i = 1 , 3 , 5 , 7 , roman_and 9 ), Mw(t)subscript𝑀w𝑡{M}_{\mathrm{w}}({t})italic_M start_POSTSUBSCRIPT roman_w end_POSTSUBSCRIPT ( italic_t ) is the moment related to waves. The GM variation term is P(t)𝑃𝑡{P}({t})italic_P ( italic_t ), which can be expressed by Eq. 3. Moreover, ω0subscript𝜔0\omega_{0}italic_ω start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT is natural roll frequency, and GM0subscriptGM0\mathrm{GM}_{0}roman_GM start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT is the metacentric height in still water.

d2ϕdt2+b1dϕdt+b2dϕdt|dϕdt|+b3(dϕdt)3+n=15a2n1ϕ2n1+P(t)ϕ=Mw(t)superscriptd2italic-ϕdsuperscript𝑡2subscript𝑏1ditalic-ϕd𝑡subscript𝑏2ditalic-ϕd𝑡ditalic-ϕd𝑡subscript𝑏3superscriptditalic-ϕd𝑡3superscriptsubscript𝑛15subscript𝑎2𝑛1superscriptitalic-ϕ2𝑛1𝑃𝑡italic-ϕsubscript𝑀w𝑡\begin{split}&\frac{\mathrm{d}^{2}\phi}{\mathrm{d}{t}^{2}}+{b}_{1}\frac{% \mathrm{d}\phi}{\mathrm{d}{t}}+{b}_{2}\frac{\mathrm{d}\phi}{\mathrm{d}{t}}% \left|\frac{\mathrm{d}\phi}{\mathrm{d}{t}}\right|+{b}_{3}\left(\frac{\mathrm{d% }\phi}{\mathrm{d}{t}}\right)^{3}\\ &+\sum_{n=1}^{5}{a}_{2n-1}\phi^{2n-1}+{P}({t})\phi={M}_{\mathrm{w}}({t})\end{split}start_ROW start_CELL end_CELL start_CELL divide start_ARG roman_d start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT italic_ϕ end_ARG start_ARG roman_d italic_t start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_ARG + italic_b start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d italic_t end_ARG + italic_b start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d italic_t end_ARG | divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d italic_t end_ARG | + italic_b start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT ( divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d italic_t end_ARG ) start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL + ∑ start_POSTSUBSCRIPT italic_n = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 5 end_POSTSUPERSCRIPT italic_a start_POSTSUBSCRIPT 2 italic_n - 1 end_POSTSUBSCRIPT italic_ϕ start_POSTSUPERSCRIPT 2 italic_n - 1 end_POSTSUPERSCRIPT + italic_P ( italic_t ) italic_ϕ = italic_M start_POSTSUBSCRIPT roman_w end_POSTSUBSCRIPT ( italic_t ) end_CELL end_ROW (2)
P(t)=ω02GM0GM(t)𝑃𝑡superscriptsubscript𝜔02subscriptGM0GM𝑡{P}({t})=\frac{\omega_{0}^{2}}{\textup{GM}_{0}}{\textup{GM}}({t})italic_P ( italic_t ) = divide start_ARG italic_ω start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_ARG start_ARG GM start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT end_ARG GM ( italic_t ) (3)

According to Section 5.2, we multiplied the modification coefficient σGMobs/σGMthrsubscriptsuperscript𝜎obsGMsubscriptsuperscript𝜎thrGM{\sigma^{\mathrm{obs}}_{\mathrm{GM}}}/{\sigma^{\mathrm{thr}}_{\mathrm{GM}}}italic_σ start_POSTSUPERSCRIPT roman_obs end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT / italic_σ start_POSTSUPERSCRIPT roman_thr end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT by P(t)𝑃𝑡{P}({t})italic_P ( italic_t ). The final equation of motion is given by

d2ϕdt2+b1dϕdt+b2dϕdt¯|dϕdt|+b3(dϕdt)3+n=15a2n1ϕ2n1+σGMobsσGMthrP(t)ϕ=Mw(t).superscriptd2italic-ϕdsuperscript𝑡2subscript𝑏1ditalic-ϕd𝑡subscript𝑏2ditalic-ϕd¯𝑡ditalic-ϕd𝑡subscript𝑏3superscriptditalic-ϕd𝑡3superscriptsubscript𝑛15subscript𝑎2𝑛1superscriptitalic-ϕ2𝑛1subscriptsuperscript𝜎obsGMsubscriptsuperscript𝜎thrGM𝑃𝑡italic-ϕsubscript𝑀w𝑡\begin{split}&\frac{\mathrm{d}^{2}\phi}{\mathrm{d}{t}^{2}}+{b}_{1}\frac{% \mathrm{d}\phi}{\mathrm{d}{t}}+{b}_{2}\frac{\mathrm{d}\phi}{\mathrm{d}\bar{t}}% \left|\frac{\mathrm{d}\phi}{\mathrm{d}{t}}\right|+{b}_{3}\left(\frac{\mathrm{d% }\phi}{\mathrm{d}{t}}\right)^{3}\\ &+\sum_{n=1}^{5}{a}_{2n-1}\phi^{2n-1}+\frac{\sigma^{\mathrm{obs}}_{\mathrm{GM}% }}{\sigma^{\mathrm{thr}}_{\mathrm{GM}}}{P}({t})\phi={M}_{\mathrm{w}}({t}).\end% {split}start_ROW start_CELL end_CELL start_CELL divide start_ARG roman_d start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT italic_ϕ end_ARG start_ARG roman_d italic_t start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_ARG + italic_b start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d italic_t end_ARG + italic_b start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d over¯ start_ARG italic_t end_ARG end_ARG | divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d italic_t end_ARG | + italic_b start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT ( divide start_ARG roman_d italic_ϕ end_ARG start_ARG roman_d italic_t end_ARG ) start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL + ∑ start_POSTSUBSCRIPT italic_n = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 5 end_POSTSUPERSCRIPT italic_a start_POSTSUBSCRIPT 2 italic_n - 1 end_POSTSUBSCRIPT italic_ϕ start_POSTSUPERSCRIPT 2 italic_n - 1 end_POSTSUPERSCRIPT + divide start_ARG italic_σ start_POSTSUPERSCRIPT roman_obs end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT end_ARG start_ARG italic_σ start_POSTSUPERSCRIPT roman_thr end_POSTSUPERSCRIPT start_POSTSUBSCRIPT roman_GM end_POSTSUBSCRIPT end_ARG italic_P ( italic_t ) italic_ϕ = italic_M start_POSTSUBSCRIPT roman_w end_POSTSUBSCRIPT ( italic_t ) . end_CELL end_ROW (4)

Applying Roberts stochastic averaging method Roberts1982 and the moment method Maruyama2023_momentamp to the equations of motion  Eq. 2 and Eq. 4 before and after correction, respectively, the PDFs of the roll amplitude are obtained. These theoretical PDFs are compared with the PDF obtained from the experiment in Figs. 11 and 12. The black dots denote the experimental result, the blue line and the red line are the theoretical results obtained by Roberts’ stochastic averaging method before and after correction, and the purple line and the green line are the theoretical results obtained by the moment method before and after correction.

Refer to caption
Figure 11: Comparison of the PDFs of the roll amplitude before and after correction, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, linear scale
Refer to caption
Figure 12: Comparison of PDF of roll amplitude before and after correction, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, logarithmic scale

The results of Figs. 11 and 12 show that the calculation results after the correction of the GM variation are closer to the experimental results than before the correction, and the accuracy of the quantitative estimation of the roll amplitude is successfully improved. However, there is still a large difference between the experimental and theoretical results, mainly due to the discrepancy in the region of small roll angles. Our theoretical method using the corrected equation of motion could not completely represent an intermittent rolling behavior as in the experiments. Therefore, further improvement of used model equations and the theoretical calculation methods is desirable.

6 PDF of the maximum roll amplitude

In order to estimate the distribution of the maximum values of the parametric roll amplitude, a theoretical expression for the PDF of the maximum roll amplitude is derived. The PDF of the maximum roll amplitude is then estimated based on the PDF obtained by the moment method, which is in quantitative agreement with the MCS results in Section 4. Moreover, the PDF estimated by the theoretical method is compared with the results of MCS in order to confirm the accuracy of the estimation of the PDF of the maximum roll amplitude.

6.1 Derivation of a theoretical formula to estimate the PDF of the maximum roll amplitude

Suppose now that N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT amplitudes are extracted from the population of roll amplitudes and the highest value among them is 𝒜Msubscript𝒜M\mathcal{A}_{\mathrm{M}}caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT, where 𝒜Msubscript𝒜M\mathcal{A}_{\mathrm{M}}caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT is a dimensionless value. Denoting the PDF of 𝒜Msubscript𝒜M\mathcal{A}_{\mathrm{M}}caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT as 𝒫(𝒜M)superscript𝒫subscript𝒜M\mathcal{P}^{*}\left(\mathcal{A}_{\mathrm{M}}\right)caligraphic_P start_POSTSUPERSCRIPT ∗ end_POSTSUPERSCRIPT ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ), the probability that the maximum value of the amplitude 𝒜𝒜\mathcal{A}caligraphic_A is in [𝒜M,𝒜M+d𝒜M]subscript𝒜Msubscript𝒜Mdsubscript𝒜M[\mathcal{A}_{\mathrm{M}},~{}\mathcal{A}_{\mathrm{M}}+\mathrm{d}\mathcal{A}_{% \mathrm{M}}][ caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT , caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT + roman_d caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ] is 𝒫(𝒜M)d𝒜Msuperscript𝒫subscript𝒜Mdsubscript𝒜M\mathcal{P}^{*}\left(\mathcal{A}_{\mathrm{M}}\right)\mathrm{d}\mathcal{A}_{% \mathrm{M}}caligraphic_P start_POSTSUPERSCRIPT ∗ end_POSTSUPERSCRIPT ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ) roman_d caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT.

This is the probability that only one of the N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT amplitudes is between 𝒜Msubscript𝒜M\mathcal{A}_{\mathrm{M}}caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT and 𝒜M+d𝒜Msubscript𝒜Mdsubscript𝒜M\mathcal{A}_{\mathrm{M}}+\mathrm{d}\mathcal{A}_{\mathrm{M}}caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT + roman_d caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT and the remaining (N01)subscript𝑁01\left(N_{0}-1\right)( italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - 1 ) amplitudes are less than 𝒜Msubscript𝒜M\mathcal{A}_{\mathrm{M}}caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT. Thus 𝒫(𝒜M)superscript𝒫subscript𝒜M\mathcal{P}^{*}\left(\mathcal{A}_{\mathrm{M}}\right)caligraphic_P start_POSTSUPERSCRIPT ∗ end_POSTSUPERSCRIPT ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ) can be expressed as follows Eq. 5 Goda2000

𝒫(𝒜M)d𝒜M=N0[1𝒜M𝒫(𝒜)d𝒜]N01𝒫(𝒜M)d𝒜Msuperscript𝒫subscript𝒜Mdsubscript𝒜Msubscript𝑁0superscriptdelimited-[]1superscriptsubscriptsubscript𝒜M𝒫𝒜differential-d𝒜subscript𝑁01𝒫subscript𝒜Mdsubscript𝒜M\begin{split}\mathcal{P}^{*}\left(\mathcal{A}_{\mathrm{M}}\right)\mathrm{d}% \mathcal{A}_{\mathrm{M}}=&N_{0}{\left[1-\int_{\mathcal{A}_{\mathrm{M}}}^{% \infty}\mathcal{P}\left(\mathcal{A}\right)\mathrm{d}\mathcal{A}\right]}^{N_{0}% -1}\\ &\mathcal{P}\left(\mathcal{A}_{\mathrm{M}}\right)\mathrm{d}\mathcal{A}_{% \mathrm{M}}\end{split}start_ROW start_CELL caligraphic_P start_POSTSUPERSCRIPT ∗ end_POSTSUPERSCRIPT ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ) roman_d caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT = end_CELL start_CELL italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT [ 1 - ∫ start_POSTSUBSCRIPT caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ∞ end_POSTSUPERSCRIPT caligraphic_P ( caligraphic_A ) roman_d caligraphic_A ] start_POSTSUPERSCRIPT italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - 1 end_POSTSUPERSCRIPT end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL caligraphic_P ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ) roman_d caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_CELL end_ROW (5)

where 𝒫(𝒜M)𝒫subscript𝒜M\mathcal{P}\left(\mathcal{A}_{\mathrm{M}}\right)caligraphic_P ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ) is the value of PDF when the amplitude 𝒜𝒜\mathcal{A}caligraphic_A is 𝒜Msubscript𝒜M\mathcal{A}_{\mathrm{M}}caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT.

From the above, 𝒫(𝒜M)superscript𝒫subscript𝒜M\mathcal{P}^{*}\left(\mathcal{A}_{\mathrm{M}}\right)caligraphic_P start_POSTSUPERSCRIPT ∗ end_POSTSUPERSCRIPT ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ), the PDF of the maximum roll amplitude, is obtained.

Also, when N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT is large enough, i.e. N0subscript𝑁0N_{0}\rightarrow\inftyitalic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT → ∞, then

ξ=N0𝒜M𝒫(𝒜)d𝒜𝜉subscript𝑁0superscriptsubscriptsubscript𝒜M𝒫𝒜differential-d𝒜\xi=N_{0}\int_{\mathcal{A}_{\mathrm{M}}}^{\infty}\mathcal{P}\left(\mathcal{A}% \right)\mathrm{d}\mathcal{A}italic_ξ = italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ∫ start_POSTSUBSCRIPT caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ∞ end_POSTSUPERSCRIPT caligraphic_P ( caligraphic_A ) roman_d caligraphic_A (6)
limN0[1𝒜M𝒫(𝒜)d𝒜]N0=limN0[1ξN0]N0=eξsubscriptsubscript𝑁0superscriptdelimited-[]1superscriptsubscriptsubscript𝒜M𝒫𝒜differential-d𝒜subscript𝑁0subscriptsubscript𝑁0superscriptdelimited-[]1𝜉subscript𝑁0subscript𝑁0superscript𝑒𝜉\begin{split}\lim_{N_{0}\to\infty}{\left[1-\int_{\mathcal{A}_{\mathrm{M}}}^{% \infty}\mathcal{P}\left(\mathcal{A}\right)\mathrm{d}\mathcal{A}\right]}^{N_{0}% }&=\lim_{N_{0}\to\infty}{\left[1-\frac{\xi}{N_{0}}\right]}^{N_{0}}\\ &=e^{-\xi}\end{split}start_ROW start_CELL roman_lim start_POSTSUBSCRIPT italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT → ∞ end_POSTSUBSCRIPT [ 1 - ∫ start_POSTSUBSCRIPT caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ∞ end_POSTSUPERSCRIPT caligraphic_P ( caligraphic_A ) roman_d caligraphic_A ] start_POSTSUPERSCRIPT italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT end_POSTSUPERSCRIPT end_CELL start_CELL = roman_lim start_POSTSUBSCRIPT italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT → ∞ end_POSTSUBSCRIPT [ 1 - divide start_ARG italic_ξ end_ARG start_ARG italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT end_ARG ] start_POSTSUPERSCRIPT italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT end_POSTSUPERSCRIPT end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL = italic_e start_POSTSUPERSCRIPT - italic_ξ end_POSTSUPERSCRIPT end_CELL end_ROW (7)

Therefore, Eq. 5 can be expressed in terms of Eq. 8.

𝒫(𝒜M)d𝒜M=ξeξ[𝒜M𝒫(𝒜)d𝒜{1𝒜M𝒫(𝒜)d𝒜}]1𝒫(𝒜M)d𝒜Msuperscript𝒫subscript𝒜Mdsubscript𝒜M𝜉superscript𝑒𝜉superscriptdelimited-[]superscriptsubscriptsubscript𝒜M𝒫𝒜differential-d𝒜1superscriptsubscriptsubscript𝒜M𝒫𝒜differential-d𝒜1𝒫subscript𝒜Mdsubscript𝒜M\begin{split}&\mathcal{P}^{*}\left(\mathcal{A}_{\mathrm{M}}\right)\mathrm{d}% \mathcal{A}_{\mathrm{M}}\\ =&\xi e^{-\xi}{\left[\int_{\mathcal{A}_{\mathrm{M}}}^{\infty}\mathcal{P}\left(% \mathcal{A}\right)\mathrm{d}\mathcal{A}\left\{1-\int_{\mathcal{A}_{\mathrm{M}}% }^{\infty}\mathcal{P}\left(\mathcal{A}\right)\mathrm{d}\mathcal{A}\right\}% \right]}^{-1}\\ &\mathcal{P}\left(\mathcal{A}_{\mathrm{M}}\right)\mathrm{d}\mathcal{A}_{% \mathrm{M}}\end{split}start_ROW start_CELL end_CELL start_CELL caligraphic_P start_POSTSUPERSCRIPT ∗ end_POSTSUPERSCRIPT ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ) roman_d caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_CELL end_ROW start_ROW start_CELL = end_CELL start_CELL italic_ξ italic_e start_POSTSUPERSCRIPT - italic_ξ end_POSTSUPERSCRIPT [ ∫ start_POSTSUBSCRIPT caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ∞ end_POSTSUPERSCRIPT caligraphic_P ( caligraphic_A ) roman_d caligraphic_A { 1 - ∫ start_POSTSUBSCRIPT caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ∞ end_POSTSUPERSCRIPT caligraphic_P ( caligraphic_A ) roman_d caligraphic_A } ] start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL caligraphic_P ( caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT ) roman_d caligraphic_A start_POSTSUBSCRIPT roman_M end_POSTSUBSCRIPT end_CELL end_ROW (8)

6.2 Calculation Method

The specific procedure for the calculation of the PDF of the maximum roll amplitude using MCS and theoretical method based on Eq. 5 is explained as follows. The initial value of the roll angle on MCS is set to 5deg5deg5~{}\mathrm{deg}5 roman_deg, and the first N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT amplitudes are extracted from the time series data of the roll angle after 500s500s500~{}\mathrm{s}500 roman_s from the start of the simulation. The maximum value among the N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT amplitudes is recorded as the maximum value in one trial. By repeating this N𝑁Nitalic_N times, N𝑁Nitalic_N maximum values are obtained in total. Based on these N𝑁Nitalic_N maxima, the PDF of the maximum roll amplitude is calculated. Then, the shapes of each PDF with N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT as the variable are compared.

On the other hand, 𝒫(𝒜)𝒫𝒜\mathcal{P}\left(\mathcal{A}\right)caligraphic_P ( caligraphic_A ), the PDF of the roll amplitude, is calculated using the moment method and Roberts’ stochastic averaging method. Using this, the PDF of the maximum roll amplitude is theoretically derived by using Eq. 5. Then, the PDF of the maximum roll amplitude calculated by the theoretical method is compared with that obtained from the MCS results.

6.3 Numerical Results

With N0=20,50,102,103,104subscript𝑁02050superscript102superscript103superscript104N_{0}=20,~{}50,~{}10^{2},~{}10^{3},~{}10^{4}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 20 , 50 , 10 start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT , 10 start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT , 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT, the PDFs of the maximum roll amplitude calculated by MCS are shown in Figs. 13 and 14. When N0=20,50,102,103subscript𝑁02050superscript102superscript103N_{0}=20,~{}50,~{}10^{2},~{}10^{3}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 20 , 50 , 10 start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT , 10 start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT, the number of MCS trials N𝑁Nitalic_N is 104superscript10410^{4}10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT. When N0=104subscript𝑁0superscript104N_{0}=10^{4}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT, N𝑁Nitalic_N is 5×1035superscript1035\times 10^{3}5 × 10 start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT.

Refer to caption
Figure 13: PDF of the maximum of roll amplitudes obtained by MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, linear scale
Refer to caption
Figure 14: PDF of the maximum of roll amplitudes obtained by MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, logarithmic scale

The results calculated by Eq. 5 using the moment method and Roberts’ stochastic averaging method are compared with that obtained from MCS and are shown in Figs. 15 and 16, and Figs. 17 and 18, and Figs. 19 and 20.

Refer to caption
Figure 15: Comparison of PDF of the maximum of roll amplitudes obtained theoretical method and MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, N0=102subscript𝑁0superscript102N_{0}=10^{2}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT, linear scale
Refer to caption
Figure 16: Comparison of PDF of the maximum roll amplitudes obtained by the theoretical methods and MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, N0=102subscript𝑁0superscript102N_{0}=10^{2}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT, logarithmic scale
Refer to caption
Figure 17: Comparison of PDF of the maximum roll amplitudes obtained by the theoretical methods and MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, N0=103subscript𝑁0superscript103N_{0}=10^{3}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT, linear scale
Refer to caption
Figure 18: Comparison of PDF of the maximum roll amplitudes obtained by the theoretical methods and MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, N0=103subscript𝑁0superscript103N_{0}=10^{3}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT, logarithmic scale
Refer to caption
Figure 19: Comparison of PDF of the maximum roll amplitudes obtained by the theoretical methods and MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, N0=104subscript𝑁0superscript104N_{0}=10^{4}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT, linear scale
Refer to caption
Figure 20: Comparison of PDF of the maximum roll amplitudes obtained by the theoretical methods and MCS, H1/3=7.0msubscript𝐻137.0mH_{1/3}=7.0~{}\mathrm{m}italic_H start_POSTSUBSCRIPT 1 / 3 end_POSTSUBSCRIPT = 7.0 roman_m, T01=10.0ssubscript𝑇0110.0sT_{01}=10.0~{}\mathrm{s}italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT = 10.0 roman_s, N0=104subscript𝑁0superscript104N_{0}=10^{4}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT, logarithmic scale

6.4 Discussion and Limitations

From Figs. 13 and 14, it can be seen that as N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT increases, the peak of the PDF moves to the right and the range becomes narrower. As N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT increases, the number of samples for searching the maximum value per trial increases. Therefore, a larger value is counted as the maximum value per trial, and the peak moves to the right. The fact that the shape of the PDF changes with N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT suggests that we need to pay attention to the number of samples N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT per trial when estimating the maximum value.

Also, looking at the Figs. 15 and 16, in terms of the position of the peak, the result obtained by the moment method is in better agreement with the result of MCS than that by Roberts’ stochastic averaging method. This may be due to the accuracy of the estimation of the tail of the PDF of the roll amplitude by the theoretical calculation, as mentioned in Section 4.4. Therefore, it can be concluded that the moment method provides a better qualitative result for the PDF of the maximum roll amplitude. On the other hand, under conditions where N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT is large, the deviation between the theoretical value and the MCS result becomes large, and the estimation accuracy remains somewhat problematic. In the future, it may be necessary to improve the theoretical estimation method of the PDF of the roll amplitude.

7 Conclusion

The accuracy of the PDF estimation of the parametric roll amplitude in long-crested irregular waves was validated by comparing the results of each theoretical method with the experimental results. In this study, the experimental results were generally smaller, especially after 20deg20deg20~{}\mathrm{deg}20 roman_deg, and differed from the results of each theoretical method. In the tail section of the PDF, the PDF result of the moment method is smaller than those obtained by Roberts’ stochastic averaging method and ESAM. Therefore, in the considered cases, the moment method is the most suitable for extreme value estimation. The experimental PDF of the roll amplitude had a shape with humps and hollows. We believe that this is due to the change in the yaw angle in the experiment. Further study is needed to develop a model which introduces an additional yaw disturbance and also to improve the experimental method.

Then, the PDFs of the maximum roll amplitude were calculated using the moment method and Roberts’ stochastic averaging method, and the results were compared with that using the MCS. The results showed that the moment method gave a better agreement with the MCS results than Roberts’ stochastic averaging method, and that the maximum value could be estimated quantitatively. However, it was also confirmed that the deviation from the MCS results becomes larger when N0subscript𝑁0N_{0}italic_N start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT is large. We conclude that more accurate estimation of the PDF of the parametric roll amplitude is needed to resolve this issue.

Acknowledgements.
This study was supported by a Grant-in-Aid for Scientific Research from the Japan Society for the Promotion of Science (JSPS KAKENHI Grant #22H01701).

References

  • (1) W.N. France, M. Levadou, T.W. Treakle, J.R. Paulling, R.K. Michel, C. Moore, An investigation of head-sea parametric rolling and its influence on container lashing systems, Marine technology and SNAME news 40(1), 1 (2003)
  • (2) A. Rosén, M. Huss, M. Palmquist, Experience from parametric rolling of ships, Parametric Resonance in Dynamical Systems pp. 147–165 (2012)
  • (3) J.E. Kerwin, Notes on rolling in longitudinal waves, International Shipbuilding Progress 2(16), 597 (1955)
  • (4) L.D. Zavodney, A. Nayfeh, N. Sanchez, The response of a single-degree-of-freedom system with quadratic and cubic non-linearities to a principal parametric resonance, Journal of Sound and Vibration 129(3), 417 (1989)
  • (5) A. Francescutto, An experimental investigation of parametric rolling in head waves, J. Offshore Mech. Arct. Eng. 123(2), 65 (2001)
  • (6) G. Bulian, Approximate analytical response curve for a parametrically excited highly nonlinear 1-dof system with an application to ship roll motion prediction, Nonlinear Analysis: Real World Applications 5(4), 725 (2004)
  • (7) K. Spyrou, Design criteria for parametric rolling, Oceanic Engineering International 9(1), 11 (2005)
  • (8) N. Umeda, H. Hashimoto, D. Vassalos, S. Urano, K. Okou, Nonlinear dynamics on parametric roll resonance with realistic numerical modelling, International Shipbuilding Progress 51(2-3), 205 (2004)
  • (9) A. Maki, N. Umeda, S. Shiotani, E.��Kobayashi, Parametric rolling prediction in irregular seas using combination of deterministic ship dynamics and probabilistic wave theory, Journal of Marine Science and Technology 16(3), 294 (2011)
  • (10) M. Sakai, N. Umeda, T. Yano, A. Maki, N. Yamashita, A. Matsuda, D. Terada, Averaging methods for estimating parametric roll in longitudinal and oblique waves, Journal of Marine Science and Technology 23(3), 413 (2018)
  • (11) L. Dostal, E. Kreuzer, N. Sri Namachchivaya, Non-standard stochastic averaging of large-amplitude ship rolling in random seas, Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 468(2148), 4146 (2012)
  • (12) A. Maki, L. Dostal, Y. Liu, L. Dostal, Comparison of stochastic stability boundaries for parametrically forced systems with application to ship rolling motion, Journal of Marine Science and Technology (Under Review) (2023)
  • (13) J.B. Roberts, Effect of parametric excitation on ship rolling motion in random waves, Journal of Ship Research 26(4), 246 (1982)
  • (14) J.B. Roberts, P.D. Spanos, Stochastic averaging: An approximate method of solving random vibration problems, International Journal of Non-Linear Mechanics 21(2), 111 (1986)
  • (15) J.B.Roberts, M.Vasta, Markov modelling and stochastic identification for nonlinear ship rolling in random waves, Philosophical Transaction of Royal Society of London A 358, 1917 (2000)
  • (16) L. Dostal, E. Kreuzer, Probabilistic approach to large amplitude ship rolling in random seas, Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science 225(10), 2464 (2011)
  • (17) Y. Maruyama, A. Maki, L. Dostal, N. Umeda, Improved stochastic averaging method using hamiltonian for parametric rolling in irregular longitudinal waves., Journal of Marine Science and Technology 27(1), 186 (2022)
  • (18) Y. Maruyama, A. Maki, L. Dostal, N. Umeda, Application of linear filter and moment equation for parametric rolling in irregular longitudinal waves, Journal of Marine Science and Technology 27, 1252 (2022)
  • (19) D.C.C. Bover, Moment equation methods for nonlinear stochastic systems, Journal of Mathematical Analysis and Applications 65(2), 306 (1978)
  • (20) Y. Maruyama, A. Maki, L. Dostal, Probability density function of roll amplitude for parametric rolling using moment equation (under review), Journal of Marine Science and Technology (2023)
  • (21) O. Grim, Beitrag zu dem problem der sicherheit des schiffes im seegang, Schiff und Hafen (in German) 6, 490 (1961)
  • (22) H. Hashimoto, N. Umeda, Nonlinear analysis of parametric rolling in longitudinal and quartering seas with realistic modeling of roll-restoring moment, Journal of Marine Science and Technology 9, 117 (2004)
  • (23) L. Yu, N. Ma, S. Wang, Parametric roll prediction of the kcs containership in head waves with emphasis on the roll damping and nonlinear restoring moment, Ocean Engineering 188, 106 (2019)
  • (24) H. Hashimoto, N. Umeda, G. Sakamoto, G. Bulian, Estimation of roll restoring moment in long-crested irregular waves, Conference Proceedings The Japan Society of Naval Architects and Ocean Engineers 3, 201 (2006)
  • (25) A. Maki, Y. Maruyama, U. Umeda, Y. Miino, T. Katayama, M. Sakai, T. Ueta, A perspective on theoretical estimation of stochastic nonlinear rolling, Proceedings of the 17th International Ship Stability Workshop, Helsinki, Finland pp. 39–46 (2019)
  • (26) Y. Goda, Random seas and design of maritime structures (Advanced Series on Ocean Engineering, vol 15. World Scientific Pub Co Inc, Singapore, 2000)

Appendix Appendix A. Method of the estimation of GM variation

A.1 Procedure for calculating GM variation

The following procedure is used to calculate the GM variation (ΔGMΔGM\Delta\mathrm{GM}roman_Δ roman_GM) in irregular waves.

  1. 1.

    Generation of time series for Grim’s effective wave

  2. 2.

    Equation relating wave elevation and GM Variation

  3. 3.

    Generation of time series of GM variation

Each calculation method is explained following.

A.1.1 Grim’s effective wave concept

This is the calculation method performed in Item 3. Grim’s effective wave is the replacement of a spatially irregular waveform around a ship by a single regular wave using the least-squares method, and this regular wave is called the effective waveGrim1961 . As in Fig. A.1, this effective wave is assumed to have a wavelength to ship length ratio of 1, and the crest or trough of the wave is in the center of the hull. The amplitude of this effective wave has a linear relationship with ocean wave displacement, while the effect of the wave on the righting lever has a nonlinear and non-memory relationship. As a result, there are no major obstacles to the statistical treatment of the GM variation in irregular waves. The time series of the effective wave displacement ζGrim(t)subscript𝜁Grim𝑡\zeta_{\mathrm{Grim}}(t)italic_ζ start_POSTSUBSCRIPT roman_Grim end_POSTSUBSCRIPT ( italic_t ) is calculated by using the ocean wave spectrum S(ω)𝑆𝜔S(\omega)italic_S ( italic_ω ). Two methods for generating the time series of effective wave is presented by Maruyama et al Maruyama2022_moment . One method is based on the superposition of component waves [Method 1] Grim1961 , and the other is based on the application of a linear filter [Method 2] Maruyama2022_moment . In this study, in Section 4 and Section 6, Method 1 was used on Roberts’ stochastic averaging method and ESAM, and Method 2 was used on MCS and the moment method. In Section 5, Method 1 was used.

Refer to caption
Figure A.1: The schematic view of Grim’s effective wave

A.1.2 Method of non-memory transformation

This is the calculation method performed by item 2. Considering only the wave component based on the Froude-Krylov assumption under quasi-statically balancing heave and pitch in waves, we obtain the equation relating wave displacement and GM variation Maruyama2022 in regular waves with the ratio of the wavelength to the ship length of 1. The wave peaks and troughs are assumed to always exist in the center of the hull. Using the numerical simulation, we calculate the GM variation ΔGMϕ,calthrΔsubscriptsuperscriptGMthritalic-ϕcal\Delta\mathrm{GM^{thr}_{\phi,cal}}roman_Δ roman_GM start_POSTSUPERSCRIPT roman_thr end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_ϕ , roman_cal end_POSTSUBSCRIPT when the wave peaks and troughs are located in the center of the hull and the heel angle of the hull is ϕ[deg]italic-ϕdelimited-[]deg\phi~{}[\mathrm{deg}]italic_ϕ [ roman_deg ]. In doing so, calculations are performed for several wave amplitudes ζmidsubscript𝜁mid\zeta_{\mathrm{mid}}italic_ζ start_POSTSUBSCRIPT roman_mid end_POSTSUBSCRIPT. As a result, the relationship diagram in Fig. A.2 is obtained. Where the wave displacement is positive when the wave trough is located in the center of the hull and it is negative when the wave crest is located in the center of the hull. The original data is plotted by the black dotted line in Fig. A.2. The red dotted line in Fig. A.3 is a polynomial approximation of the original data as in Eq. A.1, where Np=6subscript𝑁p6N_{\mathrm{p}}=6italic_N start_POSTSUBSCRIPT roman_p end_POSTSUBSCRIPT = 6. Using this relation given by Eq. A.1 between the wave displacement at the center of the hull and the GM variation, the time series of the effective wave displacement calculated in item 1 is transformed into the time series of the GM variation.

ΔGMϕ,fitthr(ζmid)=k=0NpCkζmidkΔsubscriptsuperscriptGMthritalic-ϕfitsubscript𝜁midsuperscriptsubscriptk0subscript𝑁psubscriptCksuperscriptsubscript𝜁midk\Delta\mathrm{GM^{thr}_{\phi,fit}}(\zeta_{\mathrm{mid}})=\sum_{\mathrm{k}=0}^{% N_{\mathrm{p}}}\mathrm{C}_{\mathrm{k}}\zeta_{\mathrm{mid}}^{\mathrm{k}}roman_Δ roman_GM start_POSTSUPERSCRIPT roman_thr end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_ϕ , roman_fit end_POSTSUBSCRIPT ( italic_ζ start_POSTSUBSCRIPT roman_mid end_POSTSUBSCRIPT ) = ∑ start_POSTSUBSCRIPT roman_k = 0 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N start_POSTSUBSCRIPT roman_p end_POSTSUBSCRIPT end_POSTSUPERSCRIPT roman_C start_POSTSUBSCRIPT roman_k end_POSTSUBSCRIPT italic_ζ start_POSTSUBSCRIPT roman_mid end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_k end_POSTSUPERSCRIPT (A.1)
Refer to caption
Figure A.2: The relationship between ΔGMΔGM\Delta\mathrm{GM}roman_Δ roman_GM and ζmidsubscript𝜁mid\zeta_{\mathrm{mid}}italic_ζ start_POSTSUBSCRIPT roman_mid end_POSTSUBSCRIPT
Refer to caption
Figure A.3: The polynomial approximation of the relationship between ΔGMΔGM\Delta\mathrm{GM}roman_Δ roman_GM and ζmidsubscript𝜁mid\zeta_{\mathrm{mid}}italic_ζ start_POSTSUBSCRIPT roman_mid end_POSTSUBSCRIPT