NONLINEAR ROLLING STABILITY AND CHAOS RESEARCH OF TRIMARAN VESSEL WITH VARIABLE LAY-OUTS IN REGULAR AND IRREGULAR WAVES UNDER WIND LOAD

Summary The trimaran vessel rolls strongly at low forward speed and may capsize in high sea conditions due to chaos and loss of stability, which is not usually considered in conventional limit-based criteria. In order to perfect the method of measuring roll performance of trimaran, a set of nonlinear roll motion stability analysis method based on Lyapunov and Melnikov theory was established. The nonlinear roll motion equation was constructed by CFD and high-order polynomial fitting method. The wave force threshold of rolling chaos in regular waves is calculated by Gauss-Legendre numerical integration method. The limited significant wave height of rolling chaos in random sea conditions is deduced by the phase space transfer rate, and the complex effect of wind load is superposed in the calculation. The influence of trimaran configuration on the roll system is analyzed through the state differentiation of homoclinic and heteroclinic orbit in phase portrait. The calculation of the maximum Lyapunov exponent further verified the applicability of Melnikov method, and the topological structure change of gradual failure of the rolling system is analyzed by the erosion of safe basin. The complex changes of the nonlinear damping coefficient and the nonlinear restoring moment coefficient caused by the change of the transverse lay-outs between the main hull and side hull have a significant influence on chaos and stability, and the existence of wind load has a certain weakening effect on the stability and symmetry of the system. The conclusion also further indicates the importance of the lay-outs to the dynamic stability of the trimaran vessel, which is significant for its seakeeping design.


Introduction
In Naval Architecture, "Stability" has a wide meaning, which usually involves static and dynamic stability, and its fundamentals have wider implications for the design and operation of ships [1].In recent years, with the continuous development of second-generation stability, the research on dynamic stability has been increasing for mono-hull ships.However, the method Melnikov process, Li studied the rolling motion stability of a mono-hull ship in random sea conditions and analyzed the influence of parameters on the safe basin [28].Frey and Jiang first proposed to judge the stability of random rolling based on the random Melnikov theory through the phase space transfer rate, and he defined the phase space transfer rate as the area from the safe domain to the unsafe domain within a unit period [29][30].Subsequently, it has been widely used in the study of random roll stability of ships, but mainly concentrated on mono-hull ships.Liu studied the phase -space transfer rate of random wave ship nonlinear rolling under green water [31].Zhang studied the change of phase space transfer rate with wave direction and forward speed [32].
At present, there are few studies on the roll stability of multi-hulled ships, such as trimaran vessel, in either regular wave or random sea conditions.The trimaran vessel, as a popular multihull ship type in the last two decades, its configuration is more flexible than single-hull ship, and the designer can get a better scheme by changing the transverse and longitudinal lay-outs between the main hull and the side hull [33][34][35].As we all know, with the increasing of transverse spacing, the trimaran vessel has a larger radius of transverse meta-centering and better stability than mono-hull ship in static stability, but the roll period is reduced quickly.The rolling system is a nonlinear system by the function of the nonlinear of rolling damping and restoring moment.In addition to the conventional component such as friction damping, wavemaking damping, vortex damping, and appendage damping, nonlinear damping caused by flow interference between the main hull and the side hull should be considered for the trimaran vessel [36].The nonlinearity of the restoring moment coefficient is mainly related to the irregularity of the hull, which is also more prominent for the trimaran vessel.Under severe sea conditions and low forward speed, this kind of ship will have a large nonlinear roll motion [37][38].The change of the coefficient caused by the change of the lay-outs will cause some uncertainties to the stability.So, it is necessary to establish the method to study the rolling stability for trimaran vessels.This paper focuses on the establishment of roll motion stability prediction methods for trimaran vessel in complex environments with regular and irregular waves and wind loads under the variation of transverse layouts.

Principal Dimensions
A 600-ton high-speed transportation trimaran vessel is adopted as the research object.Table 1 shows the principal dimensions.The characteristics of different transverse spacing from 4.5m to 10m will be studied.When the spacing changes, the draft, displacement, and center of gravity remain unchanged, but the hydrodynamic coefficients are changing, and the detailed data is shown in table 1 and subsequent studies.Fig. 1 shows the whole model and section lines when transverse spacing is 7.0m.
The ship roll motion equation in waves can be expressed as equation (1).Where,   , ∆  are the transverse moment of inertia and the additional moment of inertia, respectively.(̇) is damping moment, () is restoring moment, () is wave force, () is wind force.

Nonlinear roll damping model
Roll damping is always obtained by roll decay motion in still water which can be simulated by CFD or model experiments [39][40].Due to the limited test conditions, one roll decay model test of the transverse spacing CL=7m is carried out as shown in Fig. 2, and other conditions are simulated by CFD method.The implicit inconstant method is applied to solve the Reynolds Average Equations (RANS)，while  − ω model is applied as the turbulence model.The overlapping mesh technique is used as shown in Fig. 3.The CFD results are compared with experiment.It can be seen from Fig. 4 that the accuracy is well verified.Fig. 6 shows the roll decay curves of different CLs.Preliminary observation shows that trimarans with large transverse spacing decay more quickly.
The selection of roll damping model may be different for different ship types.Wassermann [41] studied a post panama ship with different damping models and found that the linear plus square plus cubic model was more suitable.Bikdash et al. [42] studied the effect of damping models of a mono-hull ship, they derived an equivalent damping condition under which the two models yields the same Melnikov predictions.For trimaran vessels, Zhang [43] studied the roll damping in detail, and LPQD damping model was discussed at low speed while the LPCD damping model was not studied in his research.As we all know, the roll damping coefficient is usually obtained by using the   − ∆ curve according to the peak points of the decay curve, where,   and ∆ are the mean values of the adjacent amplitudes and the difference of the adjacent amplitudes of the decay curve.In this paper, the two damping models were used to fit the   − ∆ curve in Fig. 5, respectively, and it was found that the error between LPQD and LPCD models is very small within 0.5%.Therefore, any of the two models can be used for the trimaran vessel.However, the LPQD model has the absolute item || ̇, in order to deal with the equations conveniently for stability and chaotic analysis, LPCD model is adopted in following studies.

Restoring moment model
The direction of restoring moment is always opposite to the roll motion.In the case of small motion, the restoring moment has a linear relationship with the roll angle, but in the case of large roll motion, this relationship turns to nonlinear.In this paper, linear plus cubic plus quintic model is adopted to describe the nonlinear restoring moment, as shown in equation (4).
() =  1  +  3  3 +  5  5 (4) It could also be written as () = ( 1  +  3  3 +  5  5 ) (5) Where,  1 ,  3 ,  5 are the lengths of restoring moment arm respectively, and can be obtained by fitting the GZ curve of a real ship.The GZ curves are obtained based on ship statics and corrected by free liquid level. is the gravitational force or buoyancy.The fitting curve of restoring moment arm when the transverse spacing is 5m and the transverse spacing is 8m are shown in Fig. 7 and Fig. 8, indicating that the above nonlinear restoring moment model has a high precision.
0 is mainly affected by the wind speed and shape above the waterline of the trimaran as shown in equation (8).The CFD method is used for calculation of the wind force coefficient   .Fig. 9 is a scalar diagram of the wind field.The calculated wind load coefficient   is 0.805 for this trimaran vessel.Then, according to different wind speeds,  0 can be obtained.
Where,  is the velocity of wind,  is windward area,  is the air density,  is the arm of the wind force.Table 2 shows the detailed parameters in equation (10) at different transverse spacing.To study the transverse spacing as a variable, the center of mass and displacement are kept constant.As we can see, different CL has different nonlinear restoring moment coefficients  1 ,  3 ,  5 , which are fitted by corresponding GZ curves.With the increase of CL, the damping generally presents an increasing trend, with the linear damping increasing all the time, while the nonlinear damping increases first and then decreases.When the CL increases to 8m, the damping coefficient can only be fitted as linear damping, which may be due to the reduction of nonlinear interference caused by too large transverse spacing between main hull and side hull.

Phase portraits analysis
The equilibrium point calculation and phase portraits analysis in the nonlinear dynamic system are usually taken as the key stages in whole stability analysis process.The stability of the autonomous system near the equilibrium point can be qualitatively analyzed by phase portrait of homoclinic or heteroclinic orbit, which is the basement and precondition for stability judgement of the non-autonomous system.Nayfeh and Balachandran [14] has provided the homoclinic and heteroclinic orbits and phase portraits research of some typical quantic nonlinear systems, which is useful for nonlinear roll system of the trimaran vessel.
We study the equation of roll motion under the action of regular beam waves and transverse wind： Wave force amplitude  is usually expressed as  = ℎ. 0 [44], where ℎ is transverse metacentric height,  0 is the effective wave slope angle, which is expressed as  0 =    0 ,   is the effective wave slope angle correction which consider the influence of ship breadth and draught limitation [45].The wave slope angle  0 =   ，  is the wave number,   is the wave amplitude.

Symmetric system with no wind
Firstly, we analyze the situation when there is no wind, for the roll motion equation (10), Set  1 = 2  ,  3 =  3 , then equation (10) becomes： Yihan Zhang, Ping Wang, Yachong Liu, Nonlinear rolling stability and chaos research of trimaran vessel Jingfeng Hu.
with variable lay-outs in regular and irregular waves under wind load 104 Set  = ,  = ̇, then, According to the nonlinear dynamic method, the small parameter ε is used for the external force，so，  1 =  1 ,  3 =  3 ,  =  (14) and equation ( 13) turns into: When  = 0, the system degrades to a non-interference autonomous system, From the perspective of energy, the restoring moment corresponds to the potential energy, which is expressed as follows, For a system without damping and external excitation, it can be regarded as a Hamilton system.In order to solve the homoclinic or heteroclinic orbits, the Hamilton quantity () is used.
According to the initial condition, ̇= 0, ̇= 0, five equilibrium points can be calculated for the system ( 16), The characteristics of five equilibrium points are closely related to the stability of the system.The characteristic equation of system (16) is shown in equation (20).
We calculated the equilibrium points and eigenvalues of all layouts, as shown in table 3，  1,2 are the characteristic values corresponding to  2 ,  3 , and  3,4 are the characteristic values corresponding to  4,  5 .
Fig. 10 shows the phase portraits of the non-interference autonomous system at different CLs.In general, for the rolling system of the trimaran vessel, according to the classification method of homoclinic and heteroclinic orbits [14], stable attractor domain consists of two heteroclinic orbits Γ 1 and Γ 2 , shown in Fig. 10-(a).
Nonlinear rolling stability and chaos research of trimaran vessel Yihan Zhang, Ping Wang, Yachong Liu, with variable lay-outs in regular and irregular waves under wind load Jingfeng Hu.
105  The characteristic value of  1 = (0,0) is a pure virtual root; therefore, it is impossible to determine its properties according to the linearization method.However, the Hamiltonian, as Lyapunov function, is greater than zero in the vicinity of  1 = (0,0)，so,  1 = (0,0) can be Yihan Zhang, Ping Wang, Yachong Liu, Nonlinear rolling stability and chaos research of trimaran vessel Jingfeng Hu.
with variable lay-outs in regular and irregular waves under wind load 106 judged as the center point which is at stable equilibrium state.When the transverse space CL is less than 8m, the characteristic values corresponding to 2 ,  3 are real numbers, and it can be determined that the equilibrium points corresponding to  2 ,  3 are the saddle points where two heteroclinic orbits intersect, which is exactly the vanishing angle of the ship's stability.The characteristic value of  4,  5 are pure imaginary roots, similar to  1 , the corresponding type of equilibrium points are also the central point, which are the unstable attractors for ship roll motion.The left orbit passing through the left saddle point and the right orbit passing through the right saddle point are homoclinic orbits, as shown in Fig. 10 However, the homoclinic orbits has exceeded the motion bounds( =  , ), in reality, ship roll motion loses its practical significance, therefore, it has no significance for discussion.It can also be seen from Fig. 10 and Table 3 that, with the increase of the transverse spacing, the abscissa of the saddle point turns away from the center point C, and the attractor area of the intersection of the two heteroclinic rails becomes larger and larger.In a physical sense, that is, the larger the stability range of the corresponding Hamilton system, and the vanishing angle of stability is infinitely close to 90 degrees, leading to more stable state for the trimaran vessel.On the other hand, with the increase of CL, the position of the saddle point is further and further away from  1 = (0,0), while the positions of the other two unstable center points are closer and closer to  1 = (0,0). 2 ,  3 and  4,  5 get closer, and finally disappears.The attractor areas of  4 and  5 bounded by homoclinic orbits is shrinking and finally disappearing as the saddle point disappears.When CL≥8m, the system is left with only one central point, which is in a globally stable stat.Physically, when other conditions of the trimaran vessel remain unchanged, and the transverse spacing CL reaches a certain distance, the vanishing Angle of stability will not exist in the quadrants of 0~ 2 , no matter how much angle it turns, it can be restored to the original state.Looking at the potential function, a potential well can be found both at the center point and the saddle point.As the CL increases, the potential energy generally tends to increase in the phase plane of roll motion.

Asymmetric system with wind
When the wind speed is constant, the wind load action to the trimaran is related to the heeling angle, which is similar to the restoring force.So, the wind load can be taken as a part of the undisturbed Hamilton system, and the nonlinear rolling motion system of the trimaran can be expressed as: When,  = 0，rolling system become the undisturbed Hamilton system.
The potential energy equation is shown as below， It is similar to the state of no wind, where Hamilton quantity is: The characteristic equation of the undisturbed autonomous system is used to solve the equilibrium point, and the Hamiltonian of the equilibrium point is used to solve the homoclinic or heteroclinic orbits, which forms the phase plane of the rolling motion.For the trimaran ships 107 with different transverse spacing, the wind speed are 0m/s, 20m/s,30m/s and 40m/s.Fig. 11 shows the phase plane and potential energy function when CL=4.5m under different wind speed.When the wind load is not considered, the undisturbed system is a symmetric system, with a center point C and two saddle points ml and mr.The saddle point is the intersection point of two heteroclinic orbits, and the potential function is symmetric about the Y-axis.With the increase of wind speed, the wind force moment gradually increases, and the symmetry of the system is destroyed.The central equilibrium point C of the system is shifted to the right, and the right saddle point develops into the intersection point of a homoclinic orbit.With the increase of wind speed, the safe basin surrounded by the homoclinic orbit becomes smaller, making the ship more vulnerable to loss of stability.It is also worth noting that the potential energy function is also asymmetrical, with the potential energy on the right plane decreasing and that on the left plane increasing.Fig. 12 shows the phase plane and potential energy function of trimaran at different transverse spacing.It can be intuitively found that, with the increase of spacing, the center point moves to the origin, and the coordinates of the saddle point and the area of the safe basin increases gradually, indicating that nonlinear system become more stable.The asymmetry of the potential function is decreasing.For the unstable left saddle point ml, the intersection point of the homoclinic orbit gradually develops to the homoclinic orbit.It is worth noting that when CL=8, the saddle point disappears and the system has only one central equilibrium point.For the ship moving plane [−1.57,1.57], the state of global stability will be reached.

Melnikov function
Melnikov function is a typical analytic method to solve the external excitation threshold of a chaotic system, although it returns the conservative side results in some engineering applications.According to the theory of Melnikov's function, when the non-autonomous system connects the unstable shape flow and stable shape flow of homoclinic or heteroclinic saddle points, it can be determined that the system has chaos in the sense of Smale horseshoe transformation [14].The meaning of Melnikov function is the distance between stable form flow and unstable form flow, which is always applied to calculate the threshold of external excitation when the system is entering the chaotic state.
For nonlinear roll motion model of the trimaran vessel, there is only heteroclinic orbits in the real motion plane as discussed above.The distance between stable and unstable form flow through saddle point  = 0 is expressed as: Where ̅ 0 () is the heteroclinic orbital equation through 0, and ((̅ 0 ())) is the trace of the matrix.
By further simplification, Where, ,， are the x coordinates of the left and right boundaries of safe basin in the phase plane.The above integrals are solved by numerical integration method of Gauss-Legendre.When ( 0 ) = 0, Melnikov function appears simple zero root, and the system loses stability due to chaos.

𝑓 ̅
Multiply  for both sides of this equation, the critical value of external excitation for the real ship is (31) For the rolling system,   is the threshold of wave force when chaos happens.According to the quantity of Hamilton,  = ±√2(, ) − ( 1  2 +  3  4 /2 +  5  6 /3) +  0 ( −  3   3 Because the above equation involves the  6 term, it is impossible to get the analytic expression of heteroclinic orbit, which usually needs to be solved by numerical integration.
From  = ̇= /，after integrating both sides, The three points integral formula of Gauss-Legendre integral is used to calculate the integral 1 ,  2 ,  3 ,  4 , and the gaussian points and weight coefficients can refer to the numerical integration textbook.
The above method is applied to calculate the wave force threshold under different wave frequencies for the trimaran with different transverse spacing.As for the trimaran itself, with the increase of wave excitation frequency, the threshold also increases.Fig. 13 shows the wave threshold with no wind force.From the perspective of layout variation, the threshold value generally increases with the increase of CL, which is identical to the results of phase portrait analysis of Hamilton system.However, considering the nonlinear roll damping, when the roll motion is violent, nonlinear roll damping will produce great effect to the roll system, which can be found from table 2, when spacing increasing after CL=5.5m, nonlinear damping shows the tendency of decrease, so for CL=6.0 m, CL=6.5 m and CL=7.0 m, their threshold values are less than some small spacing, especially in the high frequency obviously.Theoretically, when the wave force excitation exceeds the threshold, chaos will occur in the system, causing unstable roll motion and capsizing.Fig. 13 Wave force threshold at different CLs of the trimaran vessel (Vwind=0m/s) Fig. 14 shows the threshold values of the distance CL=4.5m for different wind speeds.It can be seen that the wave threshold gradually decreases with the increase of wind speed, which is consistent with the physical phenomenon.Fig. 15 shows the wave threshold values at different CL.With the increase of CL, the system stability is enhanced, which is consistent with the conclusion of phase portrait analysis.

Numerical Verification
The threshold of chaos can be calculated based on Melnikov function, but it is not a sufficient condition for chaos to occur in the system.In order to verify the accuracy of the wave force threshold for trimaran roll system, other methods are needed for numerical verification.In the field of nonlinear dynamics, there are some methods for numerical verification of chaos, including Lyapunov exponent, power spectrum, fractal dimension, etc.In this paper, Lyapunov exponent is used to judge the chaotic characteristics of a certain threshold, and the global stability of rolling system is studied by the erosion of safe basin.

Lyapunov exponent
To the equilibrium point, the stability can be judged by the real part of characteristic value of the system's Jacobi matrix or the characteristic exponent of the periodic motion perturbation equation.Lyapunov exponent is a generalization of characteristic value and characteristic exponent, and gives a measure of the average divergence or average convergence of any adjacent orbits of the system, which is the most reliable measure to judge chaos.
If the maximum Lyapunov exponent is less than zero, the system is stable.If the maximum Lyapunov exponent is greater than zero, chaos will appear in the system.For the nonautonomous two-dimensional roll motion system, there are only two Lyapunov exponents.The most common solution of Lyapunov exponents is RHR algorithm [46].
For two dimensional non-autonomous systems as follows.
Obviously, the two Lyapunov exponents of the two-dimensional non-autonomous system are the same as the first two Lyapunov exponents of the corresponding three-dimensional autonomous system.
Carry on the division to autonomous system, and  ̇() = ()() , (0) =  3 is the identity matrix，() is the Jacobian matrix.According to (), the third-row element of Y(t) is known, and Y(t) is a non-singular matrix, so QR decomposition can be performed as follows. [ Thus, the solution system of the Lyapunov exponent can be expressed as,， Where, The two Lyapunov exponents of two-dimensional non-autonomous system were further solved as： The parameters of the roll motion system are substituted into the equations ( 39), we can get the final calculation system.
)]cos 2  (42) In this paper, the 4-order Runge-Kutta method is applied to solve the equations ( 22) and (42), and the Lyapunov exponent can be obtained by combining equations ( 40) and (41).For the rolling system of trimaran vessel with different transverse layouts, it is verified according to the wave force excitation threshold obtained by Melnikov function.Fig. 16 and Fig. 17 show the Lyapunov exponent of two CLs when the wave frequency is 0.4 rad/s, and the ft refers to wave force threshold.For the selected three lay-outs, when the excitation value is greater than the threshold, the maximum Lyapunov exponent is greater than zero, and the system is in a chaotic state, however, when the excitation is less than the threshold, the maximum Lyapunov exponent is less than zero, and the system is in a stable state.This result effectively verifies the wave force threshold calculated by Melnikov function.

Safe basin analysis
For the roll motion system, the safety region can be defined as follows: on the phase plane, the numerical simulation of roll motion is carried out with different rolling angle and angular velocity as initial values.When the motion is not divergent, the rolling system is stable.By representing the distribution of the stable points and unstable points on the phase plane, the safe basin diagram of roll motion can be obtained [28,[47][48].It can also be applied for the further verification of the wave force threshold obtained by Melnikov function.Take three kinds of transverse lay-outs for example, Fig. 18 show the safe basin of CL=4.5m,CL=5.0m and CL=6.0m, and keep the wave frequency as 0.4 rad/s as the same as Lyapunov calculation above.As we can see, when the external excitation increases from zero to the same threshold obtained from Melnikov method, the safe basin is relatively complete, and when the excitation exceeds the threshold, the stability domain begins to break, and chaos will be generated in the system.The wave force threshold obtained from Melnikov method is well validated by safe basin analysis.Fig. 19 shows the area of the safe basin.The beginning areas of the three CLs are different, and the larger the CL is, the larger the area will be, which is related to the restoring force coefficient [49].However, with the increase of wave force excitation, when the threshold is Nonlinear rolling stability and chaos research of trimaran vessel Yihan Zhang, Ping Wang, Yachong Liu, with variable lay-outs in regular and irregular waves under wind load Jingfeng Hu.
115 reached, the area of safe basin starts to break.The trimaran vessel with small CL first begins to break down, and the rolling system becomes more unstable than the ship with big CL.From the topological analysis, when the wave force excitation reaches the threshold, the system is in an unstable and chaotic state, and the phenomenon of multi-value and bifurcation of rolling will appear.We can track the trail of a point in motion plane to analyze this phenomenon.Point a is selected for example, when the transverse spacing CL of the trimaran vessel is 4.5m, and its initial state is (1.45, -0.2), as shown in the Fig. 18-(c).The trail of point a are obtained by solving the nonlinear motion equation (1) with numerical method.
As shown in Fig. 20-(a), when the wave force excitation is less than the threshold, the motion is attracted to the equilibrium point  1 from the initial point a and makes stable resonant motion.when the wave force excitation exceeds the threshold value, as shown in the Fig. 20-(b), it can be found that the trail starting from point a escapes from the stable safety region, and a large roll motion and capsizing occurs.Obviously, the trial is attracted by the attractor corresponding to the equilibrium point  4 in the phase plane, and the detailed value of  4 is listed in Table3.In addition, to continue to increase the wave force excitation, the trail starting from point a has already lost the significance of discussion in the actual ship concept, and the trimaran vessel has already capsized.However, from the perspective of theoretical research， based on the bifurcation topological theory, we can find in Fig. 20-(c) that when the excitation exceeds a certain level, the motion begins to be bifurcated and multi-valued.The ship will move back and forth between the two attractive fields of  4 and  5 , and the rolling system will be in a strong chaotic situation.

Chaotic analysis in irregular waves
In the actual sea environment, the waves are irregular and random, so the deterministic method under regular wave state can no longer be used.Based on the Melnikov function theory of deterministic systems, stochastic dynamics can be studied from the perspective of statistical with variable lay-outs in regular and irregular waves under wind load 116 analysis.Melnikov random theory is a statistical extension of the simple zero-point solution of Melnikov function for deterministic systems.Combined with phase flow function and phase space transfer rate, the global stability of stochastic systems can be analyzed.

Stochastic Melnikov process and phase space transfer rate
For Hamilton single degree of freedom system with bounded weak disturbance, using small parameter  , the motion equation can be expressed as: The Hamilton Value is: Where, () is the bounded stochastic excitation.The Hamilton system has hyperbolic saddle points and homoclinic or heteroclinic orbits. 0 (),  0 ().
Based on the deterministic Melnikov function of equation ( 26), the stochastic Melnikov process can be expressed as： According to the previous research, for a deterministic system, the chaotic threshold of external excitation can be calculated by solving the first-order zero of Melnikov function, but the first-order zero of stochastic Melnikov process cannot be solved.In order to quantitatively reflect the damage to the stability of the stochastic system, Frey and Jiang further extended the Melnikov method and introduced the phase flow function and phase space transfer rate to describe the instability of the random system [29][30].Phase flow function  is the flow from stable domain into unstable domain in the phase space.The phase space transfer rate refers to the area where the phase flow function of the system flows from the safe domain to the unsafe domain within a unit period of motion or excitation period.When the phase flow function  > 0, system lost stability into chaotic state.Therefore, the external excitation threshold of the stochastic system can be solved by solving the zero of the phase flow function.
The phase flow function is defined as： Where,  is average value of Melnikov function.According to equation (45), the Melnikov process of the random system can be expressed as: Then, the phase space transfer rate of the system is: Define the function: According to the stationary random characteristics of ergodic states of random waves, the mean time results can be described by statistical means.) −   = 0, it can be considered as the zero of the phase space transfer rate of the stochastic system, and then the threshold of external excitation can be solved.

Extreme significant wave height
The nonlinear roll motion equation of the trimaran vessel in crosswind and random beam waves can be expressed as： Where, () is the wave height at time t， 0 =  0 /( + Δ),  0 is wave force of unit wave amplitude.
According to equation ( 47), the random Melnikov process of the roll system is:   (Ω) is the power spectral function of the roll angular velocity () along a homoclinic or heteroclinic orbit.For the high-order nonlinear rolling system of the trimaran vessel, () cannot be solved analytically, so the discrete solution can only be obtained by numerical calculation.The Fig. 22 shows the discrete orbit equation of () corresponding to different wind speeds obtained by the 4-order Runge Kutta method.The power spectrum of () can be obtained by Fourier transform.As shown in Fig. 23, the power spectrum of different wind speeds in the low-frequency region is different obviously, and the higher the wind speed is, the smaller the power spectrum is.Then, according to equations ( 63) and (64), the   corresponding to different wave frequencies can be obtained.As shown in Fig. 24, with the increase of wind speed, the   when entering the rolling chaos state in the irregular wave gradually decreases.Moreover, with the change of wave period, the   tends to decrease first and then increase.For the trimaran vessel with different transverse spacing CL, when the wind speed is 30m/s, the () is obtained according to the above method.As shown in the Fig. 25, the trimaran vessel with large transverse spacing has lager area bounded by the homoclinic orbit.Fig. 26 shows power spectrum of () corresponding to different transverse spacing, and the Fig. 27 shows   of different transverse spacing.It can be clearly found that the trimaran vessel with large transverse spacing has a lager   , the stronger the stability of the rolling motion in the irregular wave is, which is consistent with the conclusion in the regular wave.It also means that the roll motion in irregular waves is more stable, which is consistent with the conclusion in regular waves.Combined with the actual sea conditions, the trimaran ship studied in this paper is prone to unstable rolling in the region with a small wave period.In the condition of 30m/s of the cross wind, the minimum   of a trimaran vessel with CL= 7.0m is about 12m, so it will not encounter such sea conditions in its sailing sea area, but unstable rolling may occur in severe weather such as typhoon, and it's even more dangerous for trimarans of small CL.From the point of view of limitation, a single-hull ship can be considered as a trimaran with a very small transverse spacing.Therefore, the rolling stability of trimaran vessels with the same tonnage and the same center of gravity is obviously better than that of mono-hull ships in theory.

Conclusions
It is significant to study the stability and chaos of roll motion for the trimaran vessel with the variation of layout, which is the main parameter in the design of the trimaran vessel.In  Nonlinear rolling stability and chaos research of trimaran vessel Yihan Zhang, Ping Wang, Yachong Liu, with variable lay-outs in regular and irregular waves under wind load Jingfeng Hu.
121 rolling system, both damping coefficients and restoring moment coefficients are important for the stability analysis.With the CFD method and high-order polynomial fitting of GZ curves, the characteristics of the damping coefficients and restoring moment coefficients are found, which has the deep relationship with the chaotic analysis.
In the phase portrait analysis, with the increase of the transverse spacing, the stable attraction domain gradually expands, and the saddle point disappears after reaching a certain distance, and the system reaches the state of global stability.After the wind load is added, the system has the asymmetric migration, and the heteroclinic orbits surrounding the stability region change topologically, from heteroclinic orbits to homoclinic orbits, and the equilibrium point also shifts accordingly with the decrease of the stability.Based on the Lyapunov stability theory and Melnikov method, the wave force threshold of rolling chaos for trimaran vessel is calculated in regular waves.The wave force threshold increases gradually with the increase of transverse spacing, indicating the enhancement of system stability.We also find that the nonlinear roll damping is sensitive to wave force threshold.The increase of the restoring force will increase the area of the safety region and increase the external excitation of chaos.The Lyapunov exponent and safe basin analysis verifies the threshold by Melnikov method.The extreme significant wave height of the trimaran vessel is calculated by using the phase space transfer rate in irregular waves, which can be the reference to the design of the trimaran vessel.

Fig. 19
Fig.19 The change curves of the safe basin area Yihan Zhang, Ping Wang, Yachong Liu, Nonlinear rolling stability and chaos research of trimaran vessel Jingfeng Hu.
Based on RHR algorithm ，the orthogonal matrix  ̅ () is transformed as follows.