一种计算随机海况下船舶参‑强激励横摇响应概率的方法与流程

本发明属于运输、信息等技术领域，涉及船舶在真实随机海况下的运动响应预报。

背景技术：

目前，随着海上运输船舶越来越向着大型化、高技术化、高性能化的方向发展，完整稳性要求成为确保船舶安全性的重要基础，也成为船舶设计建造的强制性法定要求。众所周知，横风及横浪是船舶航行的极端危险海况，船舶在该海况下的摇荡运动是造成船舶倾覆的主要原因。一直以来人们都重点研究关注了这一问题，到目前为止，已经针对其取得了大量而深入的研究成果。但是，通过对完整船舶海难事故的观察与统计，结果表明，完整船舶的海难事故竟有一半是发生在船舶处于纵浪或斜浪的航行状态，而且长久以来，船舶在遭遇恶劣的大风浪海况时，都会尽量避免直接承受横向的风浪载荷，而选择把船舶调整为纵浪航行或斜浪航行，这也就使得船舶在纵浪和斜浪中的稳性损失问题愈加凸显。根据目前研究进展可知，即使船舶在静水中的稳性满足稳性衡准要求，船舶在纵浪或斜浪航行时，仍发生大角度横摇甚至倾覆的原因之一，就是船舶发生了参数横摇。

作为参数横摇现象的一个典型案例，1998年10月，巴拿马型c11级集装箱船“aplchina”号在从高雄驶往西雅图的途中遭遇到了恶劣风暴天气，在船长按照常规选择减速迎浪航行后还是产生了极其剧烈的横摇运动，横摇角甚至达到了三十五度至四十度，事故直接导致了三分之一的集装箱丢失，三分之一的集装箱损坏。

特别是斜浪中航行的船舶，该船舶除了可能受到参数激励外，由于航向角的存在，其各横剖面左右两舷会存在一定的波面高度差，因此该船舶还会受到横向的波浪强迫激励，即同时受到参-强联合激励的作用，这时船舶的横摇现象也就变得更为复杂。

对于规则波中斜浪航行的船舶，可以通过数值模拟的方法，较精确的求得船舶在某一特定波浪条件(波长l、波高h)下的稳定横摇角。但是对于真实随机海况，我们不能确定船舶出海后遭遇到的具体波浪情况，因此单纯数值模拟船舶在某一组随机种子数下的参-强激励横摇响应的时间历程曲线，并不具有太大的现实意义，发明人更重视船舶在随机参-强激励下的横摇响应的统计特征，也就是对船舶的响应概率分布进行预报。

技术实现要素：

本发明提供一种计算随机海况下船舶参-强激励横摇响应概率的方法。本发明的方法可以计算船舶在真实随机海况下的参-强激励横摇运动响应的概率分布。本发明采用如下技术方案：

一种计算随机海况下船舶参-强激励横摇响应概率的方法，包括下列步骤：

一、建立船舶参-强激励单自由度横摇运动方程，即数学振动方程：

其中，φ为横摇角；ixx为横摇惯性矩，axx(ωn)为横摇附加惯性矩；b1与b3分别为线性阻尼系数与立方非线性阻尼系数；δ为排水量；g为重力加速度；gz为船舶复原力臂，由φ、η和ψ决定；m为波浪强迫激励力矩；

固定坐标系下，沿x轴传播的波长等于船长的长峰波在某一时刻的波形表达式为ηcos(2πx/l+ψ)，其中η为该长峰波的波动幅度，η＞0为波峰，η＜0为波谷，ψ为船和波的相对位置，取值范围为[0-2π]。

二、建立基于切片理论求解船舶复原力臂函数，用于得到数学振动方程中的船舶复原力臂gz(φ，η，ψ)：

1)应用gauss积分，得到船舶横摇复原力矩为：

式中，与为各浸水横剖面形心b0在参考坐标系下的坐标，s(x′)为各浸水横剖面面积，ρ为海水密度，l为船舶垂线间长，β为航向角，f(x′)为各横剖面的压力梯度系数，其表达式为：

式中，为船舶静水中直立状态时各横剖面的水线宽，d(x′)为船舶静水中直立状态时各横剖面的吃水，k为波数。

2)将gz(φ，η，ψ)展开为关于横摇角φ的多项式和傅里叶级数的组合形式，采用最小二乘法对gz(φ，η，ψ)进行四维拟合，得到拟合参数qi(i＝1，2，3)。

gzapp(φ，η，ψ)＝q1φ+q2φ³+q3ηcφ

三、利用波面角求解随机波浪强迫激励力矩，用于得到数学振动方程中的波浪强迫激励力矩m。

其中，q0为船舶的初稳性高gm。

四、针对船舶参-强激励横摇运动方程进行无因次处理，得到无因次的简化的数学振动方程。

将gzapp(φ，η，ψ)与m代入船舶参-强激励横摇运动方程，令ηc＝ξt，并对方程进行无因次处理，得到：

其中，

五、应用谱分析方法处理随机海况，用于得到数学振动方程中的随机激励项ηc＝ξt：

1)将随机海浪波面升高处理为窄带随机过程，根据窄带随机过程的莱斯表示法，真实海况的有效波面函数zeff(x，t)存在如下解析表达式：

zeff(x，t)＝ηc(t)cos(2πx/l)-ηs(t)sin(2πx/l)

其中ηs(t)和ηc(t)均为随机过程；

2)采用二阶受控自回归滑动平均模型准确近似随机海况，即应用carma(2，1)过程的谱密度函数近似已知的海浪谱，从而得到描述随机海况的另一重要特征参数——自相关函数。

六、应用能量包线随机平均法求解数学振动方程，得到船舶横摇响应的概率分布。

本发明可以应用在以下两个方面：

1)船舶设计阶段，由设计船舶的计划航行海域的海浪谱，可以预报该设计船舶在计划海域航行时，横摇角度超过某一极限值的概率，如果该概率过大，即可认为该船舶容易发生参数横摇，需要重新调整设计船型。

2)船舶航行阶段，可协助船长做出操船决策。如果船舶突遭恶劣海况，船长在决定减速迎浪航行前，可通过本发明迅速预报该恶劣海况下船舶横摇角度超过某一极限值的概率，再通过适当调整航速、航向，避开容易发生参数横摇现象的遭遇频率范围。

本发明的技术效果如下：

1)快速性

预报船舶随机参-强激励横摇响应概率分布，常规做法是应用蒙特卡洛原理，大量统计不同随机种子数下，数值模拟得到的船舶长时间的横摇响应时间历程曲线。通过对响应样本值进行统计，得到横摇响应各横摇角度区间的概率。该方法的精确性，要以统计数据量为基础，即会不可避免地牺牲时间效率。

而本发明对随机海况的处理是基于海浪谱的表达式，而非随机种子数，任意特定海况对应唯一的特征海浪谱，因此无需重复计算，大大提高计算速度。同时还可避免因随机种子数不足(无法完全描述真实海况)而产生的误差。

2)准确性

本发明基于非线性随机动力学理论。不同于数值模拟方法，本发明是基于严谨的数学推导过程给出的计算方法，计算结果具备更高的准确性。

3)普遍性

本发明求解船舶响应概率的本质，即解析求解其数学振动方程。船型与海况的变化并不会影响方程求解的理论推导过程，因此，本发明适用于任何船型在任何随机海况下的参-强激励横摇响应概率求解，且无需进行任何特殊处理，通用性强。

附图说明

图1数值求解复原力臂函数gz(φ，η，ψ)的流程图

图2求解平稳概率密度函数pst(h)的流程图

图3c11集装箱船型线图

图4复原力臂曲面(ψ＝0)

图5复原力臂函数的三维拟合结果和数值计算结果对比图，(a)关于(φ，ψ)的对比结果(η＝10m)；(b)关于(φ，η)的对比结果(ψ＝0rad)

图6海浪谱密度分解

图7随机过程ηc(t)的谱密度函数与carma(2，1)的拟合效果

图8关于能量h的概率密度函数

图9关于特定能量值h下船舶最大横摇角b的概率密度函数

具体实施方式

下面对本发明进行详细说明

一、建立船舶参-强激励单自由度横摇运动方程，即数学振动方程。

其中，φ为横摇角；ixx为横摇惯性矩，axx(ωn)为横摇附加惯性矩；b1与b3分别为线性阻尼系数与立方非线性阻尼系数；δ为排水量；g为重力加速度；gz为船舶复原力臂，由φ、η和ψ决定；m为波浪强迫激励力矩。

固定坐标系下，沿x轴传播的波长等于船长的长峰波在某一时刻的波形表达式为ηcos(2πx/l+ψ)，其中η为该长峰波的波动幅度(η＞0为波峰，η＜0为波谷)，ψ为船和波的相对位置，取值范围为[0-2π]。

二、基于切片理论求解船舶复原力臂函数，用于得到数学振动方程中的复原力项gz(φ，η，ψ)。

1)应用gauss积分，得到船舶横摇复原力矩为：

式中，与为各浸水横剖面形心b0在参考坐标系下的坐标，s(x′)为各浸水横剖面面积，ρ为海水密度，l为船舶垂线间长，β为航向角，f(x′)为各横剖面的压力梯度系数，其表达式为：

式中，为船舶静水中直立状态时各横剖面的水线宽，d(x′)为船舶静水中直立状态时各横剖面的吃水，k为波数。

2)将gz(φ，η，ψ)展开为关于横摇角φ的多项式和傅里叶级数的组合形式，采用最小二乘法对复原力臂函数gz(φ，η，ψ)进行四维拟合，得到拟合参数qi(i＝1，2，3)。

gzapp(φ，η，ψ)＝q1φ+q2φ³+q3ηcφ

三、利用波面角求解随机波浪强迫激励力矩，用于得到数学振动方程中的强迫力项m。

其中，q0为船舶的初稳性高gm。

四、针对船舶参-强激励横摇运动方程进行无因次处理，得到无因次的简化的数学振动方程。

将gzapp(φ，η，ψ)与m代入船舶参-强激励横摇运动方程，令ηc＝ξt，并对方程进行无因次处理，得到：

其中，

六、应用谱分析方法处理随机海况，用于得到数学振动方程中的随机激励项ηc＝ζt(与gz(φ，η，ψ和m有关)。

1)将随机海浪波面升高处理为窄带随机过程，由窄带随机过程的莱斯表示法可知，真实海况的有效波面函数zeff(x，t)存在如下解析表达式：

zeff(x，t)＝ηc(t)cos(2πx/l)-ηs(t)sin(2πx/l)

其中ηs(t)和ηc(t)均为随机过程。不使用谐波叠加法的表达形式，从而大大简化了后续的分析过程；

2)采用数学方法描述物理随机过程。本发明采用概率模型(二阶受控自回归滑动平均模型)准确近似随机海况，即应用carma(2，1)过程的谱密度函数近似已知的海浪谱，从而得到描述随机海况的另一重要特征参数——自相关函数。

六、应用能量包线随机平均法求解数学振动方程，得到船舶横摇响应的概率分布。

基于流程图2，可得到船舶随机参-强激励下的横摇响应的平稳概率密度函数pst(h)、pst(b)，将平稳概率密度函数在某一角度范围内积分，可进一步得到b(h)在该角度范围内的概率。

下面结合实施例进行说明。

c11级集装箱船属于容易发生参数横摇现象的“非常规”船型，以c11集装箱船为例，具体研究其在特征波高3m、特征波长262m的ittc海浪谱和航向角30°、航速1.43m/s的船舶航行参数下的横摇响应概率密度函数。已知该船的主尺度与型线图如图3。

表1c11集装箱船主尺度

在φ∈[-0.60.6]rad、η∈[-1010]m、ψ∈[02π]rad范围内数值模拟c11集装箱船的复原力臂函数gz(φ，η，ψ)，得到不同φ、η和ψ时的复原力臂函数gz(φ，η，ψ)，图4列出了ψ＝0rad时复原力臂函数关于变量φ和η的变化曲面。

对上述复原力臂函数gz(φ，η，ψ)采用最小二乘法求取拟合参数qi(i＝1，2，3)，得到q1＝2.164、q2＝-1.7753、q3＝-0.1061。c11集装箱船的部分三维拟合图像如图5所示。

由拟合参数qi(i＝1，2，3)与c11集装箱船主尺度(表1)，得到无因次数学振动方程的各项系数分别为其中，和特征波长及航向角有关。

选用1978年第八届国际拖曳水池会议(ittc)海浪谱作为外部海洋环境，当特征波高3m、特征波长262m、船舶航速为1.43m/s、航向角为30°时，可得到：

2)遭遇频率下的海浪谱s(ωe)与其分解出的随机过程ηs(t)和ηc(t)的谱密度函数sηs(ωe)和sηc(ωe)

3)应用carma(2，1)过程拟合随机过程ηc(t)，得到拟合效果图与拟合系数：

表2carma(2，1)拟合系数

由系统参数a1、a2、b1，可得到随机过程的自相关函数，其只与时间间隔有关，具体表达式如下：

其中，τ为时间间隔，μ1，2的表达式为：

由能量包线随机平均法，并根据流程图2，得到c11集装箱船在随机参-强激励下的横摇响应概率密度函数pst(h)、pst(b)如图8所示。

将概率密度函数在某一角度范围内积分，可最终得到b(h)在该角度范围内的概率。