大风浪环境下的船舶气象航线优化设计

崔康靖 郑元洲 陈国成 胡卫东

(武汉理工大学航运学院1) 武汉 430063) (内河航运技术湖北省重点实验室2) 武汉 430063)(湄洲湾港引航站3) 莆田 351100)

0 引 言

随着经济全球化，海运量已占据全球商品贸易的80%以上[1]，而由气象因素导致的海难事故时有发生．在气象航线及其优化方面国内外学者已有众多研究．Park等[2]通过分阶段来优化油耗问题，在航速恒定下用A*算法求解转向点位置，结合几何规划法重新分配各航路航速．Zaccone等[3]提出了一种基于三维动态规划的船舶航线优化方法，通过分配航路点寻求最优航线．张进峰等[4]针对我国近海船舶进行了避台航线优化设计，基于动态规划算法确定各航路点．在航线中船舶各航段的航速优化方面，李明峰等[5]以浪高为约束条件构建气象航线模型，引入惩罚函数应用于遗传进化算法进行优化．刘洋等[6-7]拟定危险风区为障航物，提出一种基于向量的多步滚动窗口优化算法进行航线设计．

上述研究在模拟海上气象时考虑的气象条件都比较单一，应用算法在规避恶劣气象区域时也不易设置合适的安全距离．针对大风浪问题，文中提出了一种结合人工势场法和模拟退火算法的优化方法．该方法在考虑风速和浪高的双重气象因素下对各航段航速影响，并以航程和航行时间为优化目标，利用人工势场法能较好的避让大风浪区，并通过模拟退火算法替代原有的负梯度下降法对相应的航路点位置进行调整，增加可选航路点位置来优化人工势场法航线，以寻求最优气象航线．

1 船舶气象航线设计的数学模型

1.1 航线模型定义

在考虑船舶航线设计时设定为整条航线由多个航路点组成，即将整条航线分为多个航路段，每段航路以恒向线方式表示，包含该段航路端点的航行速度、经纬度位置和航向角等信息．

设定航线中有n个航路点，用向量R表示．各个航路点的位置坐标为

R={(Lat1,Long1),…,(Lati,Longi),…,

(Latn,Longn)}

式中：Lati为第i个航路点坐标的纬度；Longi为第i个航路点坐标经度．

船舶在每段航路间将以恒向方式航行，航线中各航段的航向变量φ为

φ=[φ1,…,φi,…,φn]

(1)

式中：φi为第i段航路的航向z，由第i个和第i+1个航路点经纬度坐标决定，φi∈[0°,360°]，以正北方向为0°，顺时针方向递增．

V0为船舶在静水中航行速度，而在实际航行中船舶容易受到风浪流等外界因素的影响而导致速度发生变化，因此用向量V表示各航路段航速，即

V=[V1,…,Vi,…,Vn]

(2)

式中：Vi为第i段航路的航速．

1.2 风浪对船舶速度的影响

船舶在海上航行时，受到风、浪、涌等外界环境因素影响，船舶行驶受到的阻力增加，在主机输出功率不变时，船舶实际航速低于静水时航速，这种现象称为船舶失速．采用文献[8]以最小二乘法为基础迭代计算船舶失速公式.

Va=V0-(1.08h-0.126qh+2.77×10-3×

Fcosα)(1-2.33×10-7DV0)

(3)

式中：Va为船舶在风浪中实际速度，kn；V0为所给定的船舶静水速度，kn;F为风速大小，m/s；D为船舶实际排水量，t;h为有义波高，m;q为船舶航向与浪向间的夹角，rad；α为船舶航向与风向之间的相对夹角，rad．

1.3 大风浪区的界定

在海上航行时，不同船舶所能应对的风浪情况并不相同．根据船舶的船型、吨位、结构强度、操作性能，以及装载货物等，风浪对船舶造成的影响会有很大差异．因此，在大风浪区的具体界定时，应当考虑该船的现实情况，只要某片海域的风浪环境超出该船的承受能力，便能将其定义为大风浪区域．船舶抗风浪等级的公式为[9]

(4)

h=2(M0/g+K1DLw1-K2W0Lw1-

(5)

式中：V10为船舶在改装载状态下可安全航行的最大风速；D为船舶排水量；lq为最小倾覆力臂；Af为船舶受风面积；K为稳性衡准数；Zf为风作用力臂；h为船舶船体强度所能承受的最大浪高；M0为船舶的船中弯矩；g为重力加速度；K1为平水排水量力矩系数；Lw1为船舶水线长；K2为空船重量力矩系数；W0为空船重量；∑PiXi为所有载荷对船中的力矩；K3为波形排水量力矩修正系数；B为船宽．

2 动态航线规划模型

2.1 目标函数

1) 航线航程 根据前文中航线模型的设定，航线总航程应为各个航路段相加之和，为

(6)

式中：L为总航线航程；Li为各分段航程．

在求航线航程时，需要根据坐标点的经纬度信息求实际航行距离，由前文介绍在各航路段中航线表示为恒向线，选用墨卡托算法计算航程，为

(7)

L=(Lat2-Lat1)×secφ

(8)

式中：Lat1,Long1分别为第一个点的纬、经度坐标；Lat2,Long2分别为第二个点的纬、经度坐标；φ为恒向线方向；L为两点间航程；e为地球偏心率．

2) 航线航时 从起点到终点的总航行时间可由各航路段航行时间求和得到，为

(9)

式中：Ttotal为总的航行时间；Vi为第i段的实际航速．

2.2 人工势场算法

文中在模型仿真中采用的是人工势场算法，其原理为人通过引力场和斥力场共同作用使物体能够沿梯度下降最快的方向运动，其中目标点对物体产生引力，引导物体朝向其运动．障碍物对物体产生斥力，避免物体与之发生碰撞．物体在路径上每一点所受的合力等于该点所有斥力和引力之和．具体的算法模型如下.

常见的引力场：

(10)

式中：ξ为引力尺度因子；ρ(q,qgoal)为物体当前状态与目标的距离．

引力即为引力场对距离的导数：

Fatt(q)=-Uatt(q)=ρ(qgoal-q)

(11)

斥力场：

(12)

式中：η为斥力尺度因子；ρ(q,qobs)为物体与障碍物之间的距离；ρ0为每个障碍物的影响半径．

斥力为斥力场的导数：

(13)

引力和斥力尺度因子的取值对人工势场算法航路规划结果影响较大，如果取值不当，会导致航路不可行或航路点冗余[10]．取ξ=10且保持不变，η分别取5、8、10、12、15、18时人工势场算法规划结果见图1．

图1 η取不同值时对航线规划的影响

由图1可知，当η取5、8、10时，会出现航路途径大风浪障碍区的情形，规划航线不可取；而相比η=12，当η取15和18时，发现航路点有明显冗余现象，因此本文参数选择取ξ=10，η=12，此时人工势场算法得到的规划航路最佳．

2.3 模拟退火算法

经由人工势场法得到航线并非是最优航线，船舶驶于航线中途的航路点时，一方面可能会面临引力和斥力的合力矢量为零陷入局部最小而终止路径规划．另一方面人工势场法可能受引、斥力增益系数的影响造成迭代步伐过大，引起抖动现象，所以需采用模拟退火算法进行动态航线优化．

在人工势场法航线中加入模拟退火算法，可以通过增加随机扰动值调整新航路点，避免陷入局部最小值，还可以搜索大风浪区附近的航路点位置，替代原有的负梯度下降受力控制的方法，避免人工势场法中参数系数对航路点选取的影响，可有效缓解抖动现象．

每次经过调整之后的航路点组合便构成了一条新航线，可以利用目标函数值比较以及劣航线的接受准则进行处理．将新航线的目标函数值设为f(A)，当前最优航线目标函数值设为f(B)．若f(A)

Metropolis准则判断具体方法为

设当前航线为xi，新航线为xj，其目标函数值分别记为f(xi)和f(xj)，新航线被接受可以作为最优航线的概率Pt为

(14)

在判断过程中，系统将产生一个随机数ε∈(0,1)，并与Pt进行比较，若Pt>ε，则将新航线设为当前最优航线；否则，舍弃新航线，继续使用当前航线作为最优航线．从式(14)可见，温度越低，接受概率Pt的值也就越小，即劣航线被接受的概率会越来越小，当退火到最小值时，劣航线基本被全部淘汰．

经过多次实验，选用了以下参数：初始温度t0=100；衰减系数α=0.7；终止温度tf=10；迭代次数L=500．

模拟退火算法的具体流程为

步骤1初始参数设置 初始温度t0，衰减系数α，终止温度tf，迭代总次数L．

步骤2初始值设置 初始航线x0，当前航线xcurrent，新航线为xnew，最优航线为xbest，初始目标值设为E0，当前航线航程为Ecurrent，新航线航程为Enew，最优航线航程为Ebest，当前迭代次数l=0．

步骤3调整航路点 生成随机数n和n个2维扰动变量ΔLat,ΔLon，随机选择当前航线上n个航路点(除始点和终点)，令Lat′=Lat+ΔLat,Lon′=Lon+ΔLon,得到新航线xnew，计算其目标值Enew，迭代次数l=l+1．

步骤4计算新航线xnew，与最优航线xbest之间的航程差值ΔE，ΔE=Enew-Ecurrent，如果ΔE<0，那么令xbest=xcurrent=xnew，Ebest=Ecurrent=Enew，转至步骤6．否则转至步骤5．

步骤5令ΔE=Enew-Ecurrent，若exp(-ΔE/ti)>ε，则xcurrent=xnew，Ecurrent=Enew，其中ε为(0，1)内的随机数；否则xcurrent和Ecurrent的值不变．

步骤6若l≥L，令ti+1=αti、l=0，转至步骤7；否则返回步骤3．

步骤7若ti+1≤tf，算法终止；否则返回步骤3．

3 仿真结果及分析

3.1 气象资料获取和处理

研究采用欧洲中期天气预报中心(european center for medium-range weather forecasts，ECMWF)提供的气象数据．ECMWF主要提供中长期天气预报如实时的天气、海况信息来为远洋航行船舶提供服务．由于本文在模型中主要考虑风浪对船舶的影响，所以采用ECMWF提供海平面10 m高处的风速与风向数据．气象数据文件格式采用网络通用数据格式NetCDF(Network Common Data Format)，数据更新间隔为6 h，数据精度的网格大小为1°×1°．利用Matlab可以将气象数据可视化，更直观显示风浪分布情况[11]．图2为2019年6月太平洋区域第8日12:00时刻的风速数据信息．图3为2019年6月太平洋区域第8日中12:00时刻的浪高数据信息．

图2 风速与风向数据信息

图3 浪高数据信息

3.2 结果及分析

选用大海域太平洋作为航行水域，起点为日本横滨港(35°N,141°E)，终点为美国旧金山港(37°N,123°W)．仿真选用的船舶为标准排水量23 740 t的集装箱船S175，船舶静水速度设为15 kn，浪高与风级合并考虑，仅将船舶与风浪的遭遇角作为海况考虑[12]．选用2019年5月30日和6月8日的大风浪数据在仿真中展示效果图(见图4)，航线为大圆航线．

图4 航行途径大风浪区示意图

大圆航线为球面上两点之间距离最短的航线，即在理想环境下船舶走大圆航线将是最优选择，但在现实中船舶并不容易以大圆航线航行，如考虑危险气象环境时，航线必将有所调整．因此本文在进行航线规划时将以大圆航线为基准并利用人工势场算法进行动态调整，在经过人工势场法避障之后的航线其航程里数并不是最优解，因此再利用模拟退火算法对该航线进行优化，优化结果见图5，几种航线的航程与航时对比见表1．

图5 优化航线图

表1 三种航线结果对比

由图5和表1可知：单纯的人工势场算法所规划航线可以很好的避开大风浪区域，且相较大圆航线船舶平均航速有所提高，但航行距离和航行时间有明显增加．在人工-模拟法航线中则有较好的改进，虽总航程相较大圆航线有所增加，但总航时明显降低而且平均航速也有明显提高．

4 结 束 语

文中提出了一种结合人工势场法和模拟退火算法的气象航线优化方法，并用Matlab进行模拟仿真，在考虑风速和浪高气象因素下，调整风浪中各航段航速，优化缩短航程和航时．该方法能根据船舶自身条件设置安全距离，有效避过大风浪区，人工-模拟退火算法对航线优化效果良好，可以显著提高船舶平均航速，缩短航行时间．

文中的气象因素仅考虑了风速和浪高，为更准确模拟海上环境，后续研究中还应考虑云雾、暴雨等更多海上气象要素的影响，并从经济层面考虑，构建如油耗等模型以进一步完善优化算法．