舰船科学技术  2026, Vol. 48 Issue (7): 18-24    DOI: 10.3404/j.issn.1672-7649.2026.07.004   PDF    
一种基于马尔可夫过程的船舶优化设计方法
董新硕, 刘可峰     
江苏科技大学 船舶与海洋工程学院,江苏 镇江 212100
摘要: 针对船舶方案设计阶段难以评估随机环境条件下船舶倾覆风险的问题,本文提出一种基于马尔可夫过程的船舶主尺度优化设计方法。该方法充分考虑阻尼力矩与复原力矩的非线性特性以及风浪的随机性,建立随机风浪中船舶非线性横摇运动的伊藤随机微分方程模型;采用二阶随机龙格-库塔法求解方程,结合蒙特卡洛模拟计算船舶在随机环境条件下的倾覆概率;进而,将倾覆概率作为安全性指标融入优化过程,构建以主尺度为变量的优化模型。算例验证了方法的有效性,并进行了参数敏感性分析。研究表明,该方法计算高效、易于实施,能在方案设计阶段有效评估船舶倾覆风险,为船型方案优选提供关键决策依据。
关键词: 非线性横摇     马尔可夫过程     倾覆风险评估     船舶优化设计    
A Markov process-based method for optimal ship design
DONG Xinshuo, LIU Kefeng     
School of Naval Architecture and Ocean Engineering, Jiangsu University of Science and Technology, Zhenjiang 212100, China
Abstract: To address the challenge of evaluating ship capsize risk under stochastic environmental conditions during the preliminary design phase, this paper proposes a Markov process-based method for optimal ship design focusing on principal dimensions. The approach comprehensively accounts for the nonlinear characteristics of damping and restoring moments and the stochastic nature of wind waves, establishing a model of Itô stochastic differential equations for nonlinear ship rolling in random seas. The SRI-2 stochastic Runge-Kutta method is employed to solve the equations, combined with Monte Carlo simulations to compute ship capsize probability under stochastic environmental conditions. Subsequently, the capsize probability is incorporated as a safety metric into the optimization process, constructing an optimization model with principal dimensions as variables. Case studies validate the method's effectiveness and include parametric sensitivity analyses. Research demonstrates that this approach achieves high computational efficiency and ease of implementation, enabling effective capsize risk assessment during preliminary design while providing critical decision-making basis for hull form selection.
Key words: nonlinear roll     Markov process     capsize risk assessment     optimal ship design    
0 引 言

船舶倾覆与船上人员和货物的安全密切相关,对运营中的船舶来说是十分严重的恶性事故。船舶倾覆一旦发生,在数秒钟内就会发展到完全倾覆的状态,在如此短暂的时间内对船上人员和物资进行抢救几乎是不可能的。因此,建立船舶倾覆的理论模型并据此优化船型方案,对提升航行安全性具有重要工程价值。

船舶在风浪载荷作用下的倾覆是一个典型的具有“一果多因”特点的反问题,即船舶倾覆的发生机理并不唯一。国际海事组织(IMO)于 2020 年正式批准第二代完整稳性衡准[1]。新衡准将瘫船、参数横摇、纯稳性丧失、骑浪/横甩和过度加速度这 5 种波浪中的稳性失效模式列为船舶设计中需要核准的危险状态。处于瘫船状态时,船舶在横风横浪联合作用下易产生大幅横摇甚至倾覆,因此瘫船稳性是公认相对危险的稳性失效模式之一。本文针对随机风浪中的瘫船稳性问题,提出一种以倾覆概率为安全性指标的船舶主尺度优化设计方法。

目前,国内外不少学者对随机载荷激励下的船舶非线性横摇运动进行了研究,主要采用的研究方法有时域仿真法、Melnikov法及FPK方程法等。Kat等[2]提出了一种对船舶进行倾覆概率评估的方法,用时域仿真得到了极端横摇角的短期和长期统计分布。Belenky等[3]应用分段线性化方法对横摇运动微分方程进行简化与求解。朱杰等[4]结合线性滤波器技术,计算得到随机海浪下船舶横摇运动响应时历,运用平均条件超越率(ACER)法对时历数据进行了分析。曾柯等[5]研究了 GZ 曲线分段线性方法和不同有效波倾系数对概率计算方法的影响,计算了不同风浪环境条件下船舶的倾覆概率。唐友刚等[6]应用Melnikov函数和相空间转移率研究船舶在随机波浪中的强非线性横摇运动及其倾覆问题。

当一个线性或非线性动力学系统受到高斯白噪声激励时,响应是马尔可夫扩散过程,其转移概率密度满足FPK方程,可通过求解FPK方程得到动力学系统响应的概率密度。Nekrasov[7]运用矩函数方程研究了随机风浪作用下的船舶倾覆问题。Roberts等[8]将横摇响应的能量包络线处理为马尔可夫过程。沈栋等[9]运用马尔可夫过程理论和首次通过概念研究了船舶倾覆前的持续时间。施兴华等[10]以二维路径积分法为基础,采用 Gauss-Legendre公式计算了白噪声随机扰动和定常风倾力矩作用时横摇角概率密度函数随时间的演变。Chai等[11]提出了一种四维路径积分法来研究随机横浪中船舶的随机横摇响应和可靠性。现有基于马尔可夫过程理论的研究多是直接数值求解FPK方程,计算资源需求较高,因而鲜少应用于强调快速迭代的船舶方案设计阶段。

传统的船舶设计方法是一种基于规则的设计方法。以完整稳性为例,设计方案需满足规范要求,但无法评估最终方案的安全程度。在这样的背景下,基于风险的船舶设计方法应运而生:安全作为设计目标,而不是通过规则进行制约[12-13]。于是,安全可以被度量,更多可行方案能在设计阶段得到充分比选。当前,船舶优化设计领域,特别是针对阻力、耐波性等确定性性能指标的优化,已广泛采用遗传算法等智能优化算法,并取得了显著成效。刘鑫旺等[14]通过自由变形方法,开展了豪华游轮的型线优化设计,结合优化算法得到了不同航速下兴波阻力系数最优的船型。章瑾等[15]运用了改进的粒子群算法,完成了船型多学科优化设计,重点研究了艏部和螺旋桨的性能优化。然而,当优化目标转向随机动力可靠性这类本身即需大量计算(如蒙特卡洛模拟)的概率指标时,若直接将高耗时的概率评估方法与需要大量迭代的智能算法结合,总体计算成本将快速增长。

本文基于马尔可夫过程理论,通过直接数值求解伊藤随机微分方程,旨在降低随机风浪中船舶非线性横摇倾覆概率的计算复杂度。进而,结合试验设计法构建了以倾覆概率为安全性指标的优化流程,实现了在船舶方案设计阶段对倾覆风险的量化评估与安全性的主动提升。

1 理论基础与优化设计框架 1.1 随机风浪中船舶非线性横摇运动微分方程

风扰动力矩可以看作是平均值$ {M}_{a} $与脉动部分$ {M}_{d}\left(t\right) $两部分组成[10],其中$ {M}_{a} $是力矩取值的主要部分,即

$ {M}_{a}=\frac{1}{2}\rho A\Delta z{C}_{m}{\upsilon }_{a}{\left(z\right)}^{2}。$ (1)

式中:$ \rho $为空气密度;$ A $为船舶水线以上部分侧投影面积;$ \Delta z $为船舶水线以上部分侧投影面积形心至水压力作用点的距离;$ {\upsilon }_{a}\left(z\right) $为平均风速;$ {C}_{m} $为风压倾斜力矩系数。

$ {C}_{m}=1.4365\frac{z{L}_{\text{OA}}}{A}+0.1476\frac{B{L}_{\text{OA}}}{A} 。$ (2)

式中:$ {L}_{\text{OA}} $为船舶总长;$ B $为船宽;z为船舶水线以上部分侧投影面积形心至水线的距离。

式(2)中$ {C}_{m} $的变化范围通常在0.955~1.418之间,其取值充分反映了船舶主尺度对风扰动力矩的影响,通过选取合适的主尺度可有效降低风扰动力矩。

船舶在随机风浪中的运动理论研究的主要困难之一在于随机外载荷的处理。为了简化分析,波浪扰动力矩$ {M}_{w}\left(t\right) $计算为:

$ {M}_{w}\left(t\right)=D\cdot h\cdot {K}_{\phi }\cdot k\cdot \zeta \left(t\right) 。$ (3)

式中:$ D $为船舶排水重量;$ h $为初稳性高;$ {K}_{\phi } $为有效波倾系数,按$ {K}_{\phi }=0.13+0.6{z}_{g}/d $计算[16],其中,$ {z}_{g} $为船舶重心垂向坐标,$ d $为吃水;$ k $为波数;$ \zeta \left(t\right) $为波面升高。

大量海上观测证实,船舶工程中波面升高$ \zeta \left(t\right) $可近似看作平稳正态随机过程。当风浪谱密度具有有理谱特性时,依据马尔可夫过程理论中的Doob第二定理,$ \zeta \left(t\right) $等价于多维马尔可夫过程的一个分量。此时,$ \zeta \left(t\right) $可表征为谱密度为$ 1/\left(2\text{π} \right) $的单位高斯白噪声过程通过某线性时不变系统的响应输出。显然,此时$ {M}_{w}\left(t\right) $也为平稳正态随机过程。

建立运动方程时,选取随船移动的平衡坐标系Oxyz:当船在静水中以航速v航行时,该坐标系随船同速前进,Oxy位于静水面上,Ox正向与航速同向。为进一步简化且不妨碍所研究问题的普遍性,忽略船舶重心G偏离参考点O对运动方程的影响。船舶在随机波浪中的横摇运动受复原力矩、阻尼力矩、惯性力矩及环境载荷的作用,其可用以下随机微分方程描述。

$ \begin{split}&\left({I}_{xx}+{A}_{44}\right)\ddot{\phi }\left(t\right)+{B}_{44}\dot{\phi }\left(t\right)+{B}_{44q}\dot{\phi }\left(t\right)\left| \dot{\phi }\left(t\right)\right| +\\&\qquad D\cdot l\left(\phi \left(t\right)\right)={M}_{a}+{M}_{w}\left(t\right)。\end{split} $ (4)

式中:$ \phi \left(t\right) $ 为横摇角时历;$ {I}_{xx} $为船舶自身对Ox轴的惯性矩;$ {A}_{44} $为附加惯性矩;$ {B}_{44} $为横摇线性阻尼系数;$ {B}_{44q} $为平方项阻尼系数;$ l\left(\phi \left(t\right)\right) $为复原力臂。

静稳性曲线的非线性效应不能忽略,采用下式拟合$ l\left(\phi \left(t\right)\right) $

$ l\left(\phi \left(t\right)\right)=h\cdot \phi \left(t\right)+{C}_{3}{\phi }^{3}\left(t\right) 。$ (5)

式中:$ {C}_{3} $为三次复原力系数。

将式(4)两边同时除以$ ({I}_{xx}+{A}_{44}) $可得到:

$ \ddot{\phi }\left(t\right)+2\nu \dot{\phi }\left(t\right)+{\nu }_{3}\dot{\phi }\left(t\right)\left| \dot{\phi }\left(t\right)\right| +\omega _{\phi }^{2}\phi \left(t\right)+{c}_{3}{\phi }^{3}\left(t\right)={X}_{a}+Y\left(t\right)。$ (6)

式中:$ 2\nu = {B}_{44}/\left({I}_{xx} + {A}_{44}\right) $为横摇衰减系数;$ {\nu }_{3} = {B}_{44q}/ \left({I}_{xx} + {A}_{44} \right) $为阻尼力矩非线性项系数;$ {\omega }_{\phi } = \sqrt{D \cdot h / \left( {I}_{xx} + {A}_{44} \right)} $为横摇固有频率;$ {c}_{3}=D\cdot {C}_{3}/ \left({I}_{xx}+{A}_{44}\right) $为回复力矩非线性项系数;$ {X}_{a} = {M}_{a}/ \left({I}_{xx} + {A}_{44}\right) $$ Y\left(t\right) = {M}_{w}\left(t\right)/ \left({I}_{xx} + {A}_{44}\right) $

1.2 船舶倾覆概率计算方法

使用滤波系统将随机波浪扰动力矩处理为有色噪声,这会显著增加随机微分方程的维数,从而加大数值求解的难度。当研究大浪(5级风浪)中的共振横摇时,风浪谱是相对较宽的,在方案设计阶段可将方程(6)中的$ Y\left(t\right) $近似为具有以下谱密度的高斯白噪声[7]

$ {S}_{w}={\left(D\cdot h\right)}^{2}\frac{\omega _{\phi }^{4}}{{g}^{2}}{K}_{\phi }{S}_{\zeta }\left({\omega }_{\phi }\right)\frac{1}{{\left({I}_{xx}+{A}_{44}\right)}^{2}}。$ (7)

式中:$ {S}_{\zeta }\left({\omega }_{\phi }\right) $为风浪谱密度取横摇固有频率$ {\omega }_{\phi } $时的值,这将降低方程(6)的求解难度。

将方程(6)转化为一阶状态方程组。

$ \left\{\begin{aligned} &{\dot{x}}_{1}\left(t\right)={x}_{2}\left(t\right),\\ &{\dot{x}}_{2}\left(t\right)=-2\nu {x}_{2}\left(t\right)-{\nu }_{3}{x}_{2}\left(t\right)\left| {x}_{2}\left(t\right)\right| -\\&\qquad \omega _{\phi }^{2}{x}_{1}\left(t\right)-{c}_{3}x_{1}^{3}\left(t\right)+{X}_{a}+Y\left(t\right)。\end{aligned}\right. $ (8)

其中,$ {x}_{1}\left(t\right)=\phi \left(t\right) $$ {x}_{2}\left(t\right)=\dot{\phi }\left(t\right) $

式(8)给出了均方意义下的随机微分方程组。由于输入端的高斯白噪声在均方收敛框架下无定义,需将其表征为维纳过程的形式导数。由此,式(8)可改写为以下伊藤随机微分方程组。

$ \left\{\begin{aligned} &{\rm{d}}{x}_{1}\left(t\right)={x}_{2}\left(t\right){\mathrm{d}}t,\\ &{\rm{d}}{x}_{2}\left(t\right)=-2\nu {x}_{2}\left(t\right){\mathrm{d}}t-{\nu }_{3}{x}_{2}\left(t\right)\left| {x}_{2}\left(t\right)\right| {\mathrm{d}}t-\\ &\qquad \omega _{\phi }^{2}{x}_{1}\left(t\right){\mathrm{d}}t-{c}_{3}x_{1}^{3}\left(t\right){\mathrm{d}}t+{X}_{a}{\mathrm{d}}t+{\sigma} W\left(t\right)。\end{aligned} \right.$ (9)

式中:$ W\left(t\right) $为单位维纳过程;$ {\sigma }^{2} $为强度系数。

二阶随机龙格-库塔法是求解伊藤随机微分方程组的有效数值方法[17]。将方程组(9)用下标形式改写后可得:

$ {\rm{d}}\left[\begin{array}{l} {\phi }_{t}\\ {\dot{\phi }}_{t} \end{array}\right]=\underset{漂移项a\left(\cdot \right)}{\underbrace{\left[\begin{array}{l} {\dot{\phi }}_{t}\\ f\left({\phi }_{t},{\dot{\phi }}_{t}\right) \end{array}\right]} }{\rm{d}}t+\underset{扩散项b\left(\cdot \right)}{\underbrace{\left[\begin{array}{l} 0\\ g\left({\phi }_{t},{\dot{\phi }}_{t}\right) \end{array}\right]} }{\rm{d}}{W}_{t} 。$ (10)

式中:$ f\left(\cdot \right) $包含阻尼与恢复力矩;$ g\left(\cdot \right) $表征环境载荷强度。

迭代格式为:

1)第一阶段计算

$ \left\{\begin{aligned} & k_{1}^{p}={\dot{\phi }}_{n}\Delta t,\\ &k_{1}^{v}=f\left({\phi }_{n},{\dot{\phi }}_{n}\right)\Delta t+g\left({\phi }_{n},{\dot{\phi }}_{n}\right)\Delta {W}_{n}。\end{aligned} \right.$ (11)

式中:$ \Delta {W}_{n}={W}_{{{t}_{n+1}}}-{W}_{{{t}_{n}}}\sim N(0,\Delta t) $,通过生成随机数实现。

2)中间状态计算

$ \left\{\begin{aligned} &{\phi }_{m}={\phi }_{n}+0.5k_{1}^{p},\\ &{\dot{\phi }}_{m}={\dot{\phi }}_{n}+0.5k_{1}^{v}。\end{aligned} \right.$ (12)

3)第二阶段计算

$ \left\{\begin{aligned} & k_{2}^{p}={\dot{\phi }}_{m}\Delta t,\\ &k_{2}^{v}=f\left({\phi }_{m},{\dot{\phi }}_{m}\right)\Delta t+g\left({\phi }_{m},{\dot{\phi }}_{m}\right)\Delta {W}_{n}。\end{aligned}\right. $ (13)

4)状态更新

$ \left\{\begin{aligned} &{\phi }_{n+1}={\phi }_{n}+k_{2}^{p}+0.5g\left({\phi }_{n},{\dot{\phi }}_{n}\right)\Delta {W}_{n}\left(\Delta t\right),\\ &{\dot{\phi }}_{n+1}={\dot{\phi }}_{n}+k_{2}^{v}+g\left({\phi }_{n},{\dot{\phi }}_{n}\right)\Delta {W}_{n}\Delta t^{2}。\end{aligned}\right. $ (14)

本文通过如下步骤计算随机风浪中船舶倾覆概率:

1) 采用二阶随机龙格-库塔法求解描述船舶横摇运动的伊藤随机微分方程组(9),模拟风浪耦合激励下船舶横摇角时历;

2)设定倾覆阈值角$ {\theta }_{a} $,当模拟中横摇角超越$ {\theta }_{a} $时判定船舶倾覆;

3)基于蒙特卡洛法进行N次独立数值试验,统计倾覆事件发生次数$ N_{\rm{capsize}} $

4)通过频率估计法计算倾覆概率,即$ P_{\rm{capsize}}= N_{\rm{capsize}}/N $

1.3 船舶主尺度优化方法

初稳性高$ h $的估算公式为:

$ h={a}_{1}d+{a}_{2}\frac{{B}^{2}}{d}-\xi H。$ (15)

式中:$ {a}_{1} $$ {a}_{2} $均为系数,与船型系数及型线有关,可按近似公式选取[18]H为型深;$ \xi $为系数,按所研究船型选取。

船舶自身对Ox轴的惯性矩$ {I}_{xx} $与附加惯性矩$ {A}_{44} $可依据经验公式计算[16]

$ {I}_{xx}+{A}_{44}=\frac{D}{g}{\left(K_{xx}^{\prime}\right)}^{2}。$ (16)

其中,$ K_{xx}^{\prime} $为船体质量连同附加质量的惯性半径。计算时可取$ K_{xx}^{\prime}=CB $,系数C取决于船舶类型及装载情况。

横摇线性阻尼系数$ {B}_{44} $的计算公式为[19]

$ {B}_{44}=2\mu \sqrt{\left({I}_{xx}+{A}_{44}\right)Dh} 。$ (17)

其中,$ \mu =0.06 $为无量纲衰减系数。

利用能量关系可以证明,平方项阻尼系数$ {B}_{44q}=\displaystyle\frac{3\left({I}_{xx}+{A}_{44}\right){B}_{\phi }}{4} $ ,其中衰减系数$ {B}_{\phi } $可按渡边公式计算[16]

$ \begin{split}&{B}_{\phi }=\frac{Ld}{DhT_{\phi }^{2}}\left[\left(0.02+1.1{C}_{B}\frac{d}{L}+{\sigma }_{0}\frac{{A}_{b}}{{L}^{2}}\right)\right]\\&\qquad\left\{l_{h}^{3}\left[1+\frac{1}{4}{\left(\frac{d}{{l}_{h}}\right)}^{2}\right]+\frac{f\left({C}_{w}\right){B}^{4}}{64d}\right\}。\end{split} $ (18)

式中:L为船长(垂线间长);$ {T}_{\phi } $为横摇固有周期;$ {C}_{B} $为方形系数;$ {l}_{h} $为重心到吃水之半的距离;$ {A}_{b} $为单边舭龙骨的面积,可按$ {A}_{b}=0.03LB $估算;系数$ {\sigma }_{0} $$ f\left({C}_{w}\right) $可查表获得。

对船舶静稳性曲线进行拟合时,常采用线性项加立方项的组合,其中三次复原力系数$ {C}_{3} $取决于所研究船型及装载情况。可将$ {C}_{3} $表示为初稳性高$ h $的线性函数,同时需考虑船舶干舷f的影响。本文提出按下式对$ {C}_{3} $进行估算。

$ {C}_{3}=-\frac{{K}_{c3}}{f}h。$ (19)

式中:$ {K}_{c3} $为取决于船型的正系数。

船舶倾覆概率不仅受环境因素影响,还显著取决于其几何形状及重心位置。主尺度是描述船舶几何形状的一些最基本的特征数据,对航海性能具有决定性影响,这从式(2)、式(15)~式(19)中可直接体现。因此,在方案设计阶段,合理选取主尺度是核心任务。本文以船长L、船宽B、吃水d及干舷f为研究对象,旨在从最小倾覆概率的角度出发,确定最优的主尺度组合。针对每组船舶主尺度组合,计算其倾覆概率需要数值求解随机微分方程(9)并结合蒙特卡洛模拟,计算成本较高。为此,本文采用试验设计中的拉丁超立方抽样方法,在保证精度的同时有效降低了计算量。

在实际船舶设计中,船型方案选取需考虑多种约束条件。然而,在方案设计初始阶段,由于船舶型线尚未确定,部分约束条件难以充分评估。针对这一情况,本文提出:通过拉丁超立方抽样生成的船舶主尺度方案,应满足以下约束条件[18]。不满足约束的方案将被移除。

1)重力与浮力相平衡

满足重力与浮力的静水平衡是任何可行设计方案必须遵循的基础物理定律和首要约束条件。要求船舶各部分所受重力(重量)之和与浮力之间的相对误差不超过给定值ε,即:

$ \frac{{F}_{浮}-\sum{W}_{i}}{{F}_{浮}}\leqslant \varepsilon 。$ (20)

式中:船舶所受浮力$ {F}_{浮}=\gamma {k}_{0}LBd{C}_{B} $$ \gamma $为水的重量密度,$ {k}_{0} $为附体体积系数;$ \sum{W}_{i} $为船舶各部分重量之和,其中,船体钢料重量$ {W}_{h}={c}_{h}LBH $,舾装设备重量$ {W}_{f}={c}_{f}{\left(LBH\right)}^{2/3} $,机电设备重量$ {W}_{m}={c}_{m}{P}_{s} $,其中$ {c}_{h}、{c}_{f}、{c}_{m} $均为重量系数,按母型船选取;主机功率$ {P}_{s} $及载重量各项根据任务书及经验公式进行估算[18]

2)初稳性高约束

确保满足最低初稳性高是至关重要的强制性安全约束,该约束旨在保证船舶在小幅扰动下具备足够的复原力矩,防止倾覆,是船舶稳性设计的基础。同时,从横摇缓和性方面考虑,初稳性高应不大于其上限值,即

$ {h}_{\min }\leqslant h\leqslant {h}_{\max }。$ (21)

式中:$ {h}_{\min } $为基于稳性要求而设定的初稳性高下限,按式(10)计算;$ {h}_{\max } $为初稳性高上限,按下式计算[20]

$ {T}_{\phi }=0.58\sqrt{\frac{{B}^{2}+4z_{g}^{2}}{h}}。$ (22)

3)航速约束

航速校核的目的是初步估算船舶在给定主机功率下的设计航速,该航速需大于任务书中对航速的要求。

4)干舷约束

船舶干舷的大小对船的耐波性、抗沉性及稳性都有重大影响。对民用船舶来说,最小干舷受技术法规的约束,设计船的干舷应当符合其要求。

5)其他校核

根据设计船具体情况的不同,有些船舶还需要进行一些其他校核,比如容量校核、回转半径校核及振动频率校核等。

2 优化实例与分析 2.1 小型拖网渔船的最优主尺度选取

作为专用于拖网捕捞的船型,小型拖网渔船在渔业中扮演着重要的角色,在全球渔船中占有相当大的比例。该类型船舶在海上作业时间较长,相比其他船舶具有更大的概率遭遇恶劣的海况,造成失稳从而发生安全事故。本例以小型拖网渔船为研究对象,选定满载返航载况为典型作业工况,其主尺度取值范围如表1所示。为确保船员作业安全并留足一定安全裕度,将倾覆阈值角$ {\theta }_{\alpha } $设定为50°,其他主要计算参数如表2所示。

表 1 设计变量取值范围 Tab.1 Design variable ranges

表 2 主要计算参数 Tab.2 Main calculation parameters

约束方面:1)通过分项估算法对船舶各部分所受重力进行估算,要求总重力与浮力之间的相对误差不超过5%;2)依据国内海洋渔船法定检验技术规则(2019)[20],小型拖网渔船初稳性高应大于0.35 m,横摇固有周期按大于4.5 s计算;3)本例中不对设计航速进行校核;4)干舷约束为$ f\geqslant 0.035{L}^{2}+ 3.5L+ 190\;\text{mm} $;5)由于尺度比与船舶航海性能密切相关且存在典型范围,本文依据表3限制尺度比,以约束船型、容积及回转半径。

表 3 小型拖网渔船尺度比范围 Tab.3 Small trawler scale ratio range

进行蒙特卡洛模拟时,设置样本数量$ N=100\;000 $;模拟时间$ t=300\;{\mathrm{s}} $;时间步长$ \Delta t=0.02\;{\mathrm{s}} $;置信水平95%;初始横摇角及角速度均为0;拉丁超立方的抽样个数为100。

基于自主开发的 Matlab 程序,分别采用马尔可夫过程方法与谱方法(Pierson-Moskowitz 谱)[16],对一小型拖网渔船(L = 17.56 m, B = 5.5 m, d = 1.88 m, f = 0.52 m, CB = 0.652)在 5 级风浪下的横摇运动进行计算。其中,马尔可夫过程方法取风浪谱的极值点作为白噪声的谱密度值[7],这一处理在物理上等效于考虑了最恶劣的载荷工况。计算结果表明,马尔可夫过程方法所得横摇角的最大值与标准差(32.80°,10.15°)均大于谱方法结果(21.17°,7.04°)。此差异源于马尔可夫过程方法所采用的保守载荷假设,虽在数值上会高估响应,却为设计提供了必要的安全裕度。2种方法所揭示的横摇运动统计规律与响应时历的变化趋势(见图1)吻合,这验证了本文所采用的马尔可夫过程方法在机理与工程应用上的有效性。

图 1 2种方法的横摇角时历对比 Fig. 1 Comparison of roll angle time histories between two methods

以最小倾覆概率作为安全性指标,确定了小型拖网渔船的最优船型方案(计算结果见表4):当主尺度L = 18.88 m, B = 5.71 m, d = 1.81 m, f = 0.76 m, CB = 0.652时,船舶在选定海况下的倾覆风险最低。表4中的“平均最大横摇角$ {\overline{\theta }}_{\max } $”定义为基于蒙特卡洛法进行N次独立数值试验所得横摇角最大值的平均值。该值越小,表明对应方案的倾覆概率越低。

表 4 最小倾覆概率排序前10位的船型方案 Tab.4 Top 10 ship design solutions ranked by minimum capsizing probability
2.2 计算结果分析

优化方案同母型船的计算结果对比见表5。可以看出,优化方案在满足全部约束的前提下,适当地增加了船宽、干舷及船长值,减少了吃水值,从而有效降低了船舶在随机环境条件下的倾覆风险。

表 5 优化方案同母型船的对比 Tab.5 Comparison between optimized solution and parent ship

为更清晰地反映各主尺度变量对船舶倾覆风险的影响,本文采用单变量分析法(其他主尺度取母型船值)进行进一步研究。选取平均最大横摇角作为表征船舶倾覆风险的指标(图2),以使图例表达更为清晰。

图 2 不同参数对平均最大横摇角的影响 Fig. 2 Effect of different parameters on average maximum roll angle

图2(a)可知,船舶在横风横浪中的倾覆风险随船长L的增加呈现有限降低趋势。然而,表4的数据表明这种趋势并非绝对线性或单调;例如,最优方案中L = 18.88 m在表4中仅位列第三。在船舶设计中,船长L的选择需综合考量多种因素,如静水阻力性能、总布置可行性、操纵性及纵摇垂荡等。因此,在评估横摇倾覆风险时,船长通常并非首要或核心考量因素。

图2(b)可知,船宽B对船舶平均最大横摇角$ {\overline{\theta }}_{\max } $具有显著影响:当船宽为5.8 m时,$ {\overline{\theta }}_{\max } $为24.7°,倾覆概率仅为0.04%;然而,当船宽减小至4.7 m时,$ {\overline{\theta }}_{\max } $急剧上升至51°,倾覆概率激增至63%。这一强烈的负相关趋势在表4的数据中得到了进一步证实。其原因在于:增大船宽显著提高了初稳性高h,并使复原力臂也增大,从而增强了船舶抵抗横倾的能力。然而,船宽增加也带来一定的不利影响:一方面,它会降低船舶的进水角;另一方面,会导致横摇固有周期缩短,使船舶在风浪中产生剧烈的摇摆,影响人员舒适性和操作安全。

图2(c)可知,船舶吃水d在较宽范围内变化时,平均最大横摇角$ {\overline{\theta }}_{\max } $变化并不显著(仅从29.5°增至34°);表4数据亦表明,吃水变化并非倾覆概率改变的主因。值得注意的是,船舶设计中吃水数值的选取通常基于港口与航道水深及螺旋桨效率等因素,而非倾覆风险考量。

图2(d)可知,干舷增加仅使平均最大横摇角$ {\overline{\theta }}_{\max } $小幅下降。结合式(19)分析可知:干舷增大主要提升船舶甲板边缘入水后(静稳性曲线后段)的复原性能,而对曲线前段影响有限,故对抑制$ {\overline{\theta }}_{\max } $效果较弱。然而,船舶倾覆概率与静稳性曲线后段形状密切相关,因此受干舷数值影响显著,这与表4数据趋势一致。需指出,在实际设计中增加干舷虽可提升抗倾覆能力,但需权衡其负面效应:上层建筑增高可能引起重心上升,同时船舶受风面积亦会增大。

3 结 语

本文提出了一种基于马尔可夫过程的船舶主尺度优化设计方法。该方法将随机波浪扰动力矩近似为高斯白噪声,并以最小倾覆概率作为安全性指标。通过理论分析和实例验证,得到以下主要结论:

1)基于马尔可夫过程的倾覆概率计算具有保守性特征。据此获得的最优船型方案存在安全裕度,在方案设计阶段具备可接受的容错空间。

2)对于小型拖网渔船一类尺寸较小的船舶,其他主要要素及重心不变的条件下,增加船宽是降低横风横浪中船舶倾覆概率的最有效途径。然而,需权衡其带来的不利影响,如进水角降低和横摇固有周期缩短。增加干舷对降低倾覆概率有积极影响,但效果弱于船宽,且需充分考虑由此可能导致的重心上升和受风面积增大。在倾覆概率控制策略中,船长和吃水通常不属于敏感性设计变量。

3)当前研究聚焦于倾覆概率单目标优化,更完善的方案需融合纵摇垂荡、快速性、经济性等多目标优化,并探索智能优化算法的应用。采用高斯白噪声近似实际风浪(有色噪声)过程存在固有误差。引入滤波系统虽能提升模型精度,但将导致求解高维伊藤随机微分方程的挑战。因此,发展高效求解高维伊藤方程并将其融入船型优化框架,是极具价值的未来研究方向。

参考文献
[1]
IMO. Finalization of second generation intact stability criteria: SDC 6/WP. 6[R]. London: International Maritime Organization, 2019.
[2]
KAT J, PINKSTER D, MCTAGGART K. Random waves and capsize probability based on large amplitude motion analysis[C]//Proceedings of the International Conference on Offshore Mechanics and Arctic Engineering (OMAE). Oslo, Norway, 2002, 4: 685−694.
[3]
BELENKY V, REED A, WEEMS K. Probability of capsizing in beam seas with piecewise linear stochastic GZ curve[C]//Contemporary Ideas on Ship Stability and Capsizing in Waves, New York, 2011.
[4]
朱杰, 刘在良, 林艳, 等. 随机海浪下船舶横摇运动响应极值预报研究[J]. 中国舰船研究, 2025, 20(2): 196-202.
ZHU J, LIU Z L, LIN Y, et al. Extreme value prediction of the roll motion under random seas[J]. Chinese Journal of Ship Research, 2025, 20(2): 196-202. DOI:10.19693/j.issn.1673-3185.03589
[5]
曾柯, 顾民, 鲁江, 等. IMO船舶瘫船稳性倾覆概率研究[J]. 中国造船. 2015, 56 (4): 17−24
ZENG K, GU M, LU J, et al. Study on IMO capsizing probability under dead ship condition[J]. Shipbuilding of China. 2015, 56 (4): 17−24.
[6]
唐友刚, 谷家扬, 郑宏宇, 等. 用Melnikov方法研究船舶在随机横浪中的倾覆[J]. 船舶力学, 2004, 8(5): 27-34.
TANG Y G, GU J Y, ZHENG H Y, et al. Study on the ship capsize in random beam seas using Melnikov method[J]. Journal of Ship Mechanics, 2004, 8(5): 27-34. DOI:10.3969/j.issn.1007-7294.2004.05.004
[7]
NEKRASOV V. A. mean-square non-local stability of ship in storm conditions of operation[J]. Polish Maritime Research, 2019, 26(4): 6-15. DOI:10.2478/pomr-2019-0061
[8]
ROBERTS J. , VASTA M. Markov modelling and stochastic identification for nonlinear ship rolling in random waves[J]. Philosophical Transactions of The Royal Society A Mathematical Physical and Engineering Sciences, 2000, 358(1771): 1917-1941. DOI:10.1098/rsta.2000.0621
[9]
沈栋, 黄祥鹿. 随机波浪作用下船舶倾覆前持续时间的研究[J]. 中国造船, 2000, 41(3): 14-22.
SHEN D, HUANG X L. Lasting time before capsize of ship in random beam waves[J]. Shipbuilding of China, 2000, 41(3): 14-22.
[10]
施兴华, 张婧, 王善. 随机风浪中舰船横摇倾覆概率分析[J]. 船舶力学, 2011, 15(5): 473-479.
SHI X H, ZHANG J, WANG S. Analysis of rolling capsizing probability of warship under random wind and beam seas[J]. Journal of Ship Mechanics, 2011, 15(5): 473-479. DOI:10.3969/j.issn.1007-7294.2011.05.005
[11]
CHAI W, NAESS A, LEIRA B. Stochastic dynamic analysis and reliability of a vessel rolling in random beam seas[J]. Journal of Ship Research, 2015, 59(2): 113−131.
[12]
希)帕帕尼古拉夫, 丁毅, 等译. 基于风险的船舶设计[M]. 上海: 上海交通大学出版社, 2011.
[13]
汪敏, 娄月华, 王丽铮, 等. 基于完整稳性倾覆概率风险的内河客船主尺度优化[J]. 船舶工程, 2021, 43(10): 49-54+61.
WANG M, LOU Y H, WANG L Z, et al. Principal dimension optimization of inland river passenger ship based on intact stability capsizal probability risk[J]. Ship Engineering, 2021, 43(10): 49-54+61. DOI:10.13788/j.cnki.cbgc.2021.10.08
[14]
刘鑫旺, 万德成. 豪华邮轮多航速兴波阻力的船型优化[J]. 中国舰船研究, 2020, 15(5): 1-10+40.
LIU X W, WAN D C. Hull form optimization of wave-making resistance in different speeds for a luxury cruise ship[J]. Chinese Journal of Ship Research, 2020, 15(5): 1-10+40. DOI:10.19693/j.issn.1673-3185.01845
[15]
章瑾, 叶杨, 朱婷. 多学科设计优化在复杂船型开发中的应用[J]. 舰船科学技术, 2025, 47(7): 59-63.
ZHANG J, YE Y, ZHU T. Application of multidisciplinary design optimization in the development of complex ship forms[J]. Ship Science and Technology, 2025, 47(7): 59-63. DOI:10.3404/j.issn.1672-7649.2025.07.012
[16]
盛振邦. 船舶原理(下册)[M]. 上海: 上海交通大学出版社, 2019.
[17]
KLOEDEN P E, PLATEN E. Numerical solution of stochastic differential equations[M]. Berlin: Springer, 1992.
[18]
林焰. 船舶设计原理[M]. 大连: 大连理工大学出版社, 2016.
[19]
孙树政, 李积德, 胡开业. 船舶耐波性[M]. 哈尔滨: 哈尔滨工程大学出版社, 2022.
[20]
中华人民共和国海事局. 国内海洋渔船法定检验技术规则(2019)[S]. 北京: 人民交通出版社, 2019.