基于行走时间最短的机场停机位分配问题的一种混合算法:遗传算法加涟漪扩散模型 A Study on Marine Vessels’ Path Optimization under Typhoon Scenarios Li Xu1,2, Xiaobing Hu1,2 1 State Key Laboratory of Earth Surface Processes and Resource Ecology, Beijing Normal University, Beijing 100875, China 2Academy of Disaster Reduction and Emergence Management, Beijing Normal University, Beijing 100875, China Abstract Typhoon may cause a huge impact on the safety and costs of marine vessels on the voyage. In order to effectively solve the path optimization problem for marine vessels under typhoon scenarios, this paper proposes a hybrid algorithm integrating Genetic Algorithm (GA) with Receding Horizon Control (RHC), which considers real-time typhoon data, marine environment data including ocean currents and reefs etc., as well as safety operation requirements for marine voyage. Then, simulation experiments are conducted to test the proposed method. The experimental results indicate that the new algorithm is effective and efficient to optimize voyage paths for marine vessels under typhoon scenarios, i.e., when compared with global path optimization, the new algorithm can improve the cost-efficiency of voyage path with safety guarantee. Keywords: Marine Vessel; Real-Time Path Optimization; Typhoon; Genetic Algorithm; Receding Horizon Control 船舶海上躲避台风航线优化研究 徐里 1,2 ,胡小兵 1,2 1.北京师范大学 地表过程与资源生态国家重点实验室,北京 100875 2.北京师范大学 减灾与应急管理研究院,北京 100875 摘要:台风对海上航行船舶的安全和效益具有巨大的影响。为了更好地解决船舶海上航行避台路径优化问题, 以达到兼顾船舶安全性及经济性运营,本文提出了一种遗传算法(GA)加滑动地平线策略(RHC)的混合 算法;并结合台风实时数据、海洋洋流及岛礁等环境数据、船舶避台安全运营数据,开展了仿真实验。结果 表明,该算法在船舶海上航行躲避台风问题上有效可行,并相对于传统的全局路径规划算法,在保证船舶安 全航行的前提下,能更好的改善船舶运行的成本效益。 关键词:船舶路径优化;台风;应急管理;混合算法;滑动地平线策略 1 引言 自然灾害系统是由致灾因子、孕灾环境和承灾 体共同组成的复杂系统 [1] 。灾害风险往往是天、地、 人综合作用的结果 [2] 。对于台风这种自然灾害,其 致灾因子有大风、强降雨、巨浪和风暴潮:台风的 大风具有极强的摧毁力,可拔树、倒屋、翻船、甚 至摧毁城镇;台风过境时所形成的强降雨,在短时 间内形成大量的雨水汇集,可导致山洪,城市内涝, 农田淹没等灾害,还会导致山体滑坡,泥石流等次 生灾害的发生;巨浪和风暴潮与天文大潮共同作用 会在沿海地区形成多达几米的增水,造成沿海城市 海水倒灌、港口基础设施损毁、沿海养殖业毁灭性 损失等严重后果 [3] 。就航运领域而言,台风对海上 航行船舶的安全和效益具有巨大的影响。 在海上,强大台风中心附近的最大风速可达到 100~120 米/秒,形成的海浪可达到十几米。在风 速和巨浪的影响下,船舶承受着巨浪的撞击。当船 长小于涌浪波长时,船舶在巨浪的涌动过程中,船 舶螺旋桨会露出水面,使螺旋桨空转,极易造成螺 Received 7 April 2017 Accepted 7 May 2017 Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 72 Copyright © 2017, the Authors. Published by Atlantis Press. This is an open access article under the CC BY-NC license (http://creativecommons.org/licenses/by-nc/4.0/). 旋桨掉落,尾轴断裂等事故;当船长接近于涌浪波 长时,船舶可能会同时跨一个或两个波峰,导致船 舶发生中拱或中垂,船体结构变形,甚至船体断裂 等事故;另外船舶顶风航行时,船舶的机电设备长 时间高负荷运转,极易导致机械故障,而一旦发生 机械故障,船舶的动力受到影响,船舶将随风飘荡 而发生触礁、搁浅甚至倾覆的事故,尤其是老旧船 舶,在遭遇台风时,由于船体结构和机电设备老化, 抗台能力大大下降,更容易发生损失 [4] 。由于船舶 在抗台过程中面临各种危险,因此,各船舶管理公 司在船舶应对台风问题上均极为重视,但主要采取 被动措施,例如锚泊避风、长距离绕行,或较长时 间的滞航,其结果是虽然保证了船舶的安全,却造 成船舶燃油消耗的增加,并导致航期延误,无论是 对货主还是船务公司或租船人而言,都产生了不必 要的经济损失。由于台风对于海上航行船舶造成的 巨大影响,因此,船舶在海上航行应对台风时,寻 找一条兼顾安全性和经济性的航行路径是一个值 得研究的课题。本文提出一种将滑动地平线控制 (RHC)策略和遗传算法(GA)相结合的混合算 法,应用于船舶海上航行躲避台风航线优化问题, 以期达到在保证船舶安全航行的前提下,减少因绕 行和延期所带来的经济损失,为船舶企业运营带来 更好的成本效益。本文的余下部分安排如下。第 2 节给出了船舶在海上自由海域的航线实时优化问 题的数学描述。第 3 节给出了遗传算法(GA)加 滑动地平线控制(RHC)策略的混合算法的设计细 节和主要技术参数设定。第 4 节通过仿真实验证实 此算法的有效性。第 5 节给出了本研究的主要结 论。 2 船舶在海上自由海域航行时的路径实时优化问 题 路径优化和障碍规避是交通研究邻域里的两 个重要问题 [5-10] 。传统的路径优化问题通常都是基 于特定的路径网络进行求解 [5,6] ,而路径网络一般 是基于有限数目的结点和联结结点的链接来定义 的。有了结点和结点间的链接,路径优化的目的就 是找出一组结点和链接的子集,从而连通给定的起 点和终点,并使得给定的路径指标最优化。自由区 域路径优化问题则不同,因为在自由区域路径优化 问题中不存在路径网络,而是可以在除障碍区之外 的整个区域内自由规划路径(因而可以视为有无穷 多的结点和链接)。障碍区可以是固定于某位置的 障碍物,也可以是动态变化的特殊区域,例如台风 区域、其它运动物体等等。因此,自由区域路径优 化问题通常要求针对动态障碍区进行实时路径优 化。带动态障碍区的自由区域路径优化问题具有十 分广泛的应用背景。船舶在公海自由航行时,需要 根据暗礁区域、台风区域、及其它船舶的位置等障 碍区信息进行实时的航路规划。 目前,国内外学者已开展了许多关于船舶在海 上躲避台风航线优化这类自由区域路径优化问题 的研究,提出了扇形避台法 [11] ,目标圆法 [12,13] ,这 些方法对于船舶躲避台风优化路径选择的研究,多 凭借航海经验进行避台决策制定,可以肯定的是上 述研究对于船舶避台的安全性均有显著的效果,但 难以实现船舶运营的最佳的成本效益。文献 [14] 提出 了一种改进的遗传算法用以求解带有动态障碍区 的自由飞空域中的飞行路径优化问题。然而,由于 该改进算法采用了全局优化策略,所以很难保证在 动态环境中应用的时效性。 为了更好地解决船舶在海上躲避台风航线优 化问题,本文提出一种结合遗传算法(GA)加滑动 地平线控制(RHC)策略的混合算法。RHC 策略又 称为模型预测控制(MPC),专门用于解决在动态 环境下的实时优化问题。简单地讲,RHC 策略是提 前 N 步进行的在线优化策略。在每个时间间隔起 点,基于当前的可用信息,RHC 策略优化接下来最 近的 N 段时间间隔的路径,同时只执行当前的一 个时间间隔的优化路径。到下一个时间间隔,RHC 将基于更新的信息,重复计算接下来的 N 段时间 间隔的优化路径,以此类推。RHC 策略现在已经被 控制工程领域广泛接受,与其他控制方法相比, RHC 策略具有许多优点 [15] 。最近,RHC 策略在运筹 学领域的应用研究也受到了很多重视。例如,在 文献 [16] 讨论了关于将 RHC 策略应用于各类离散事 件系统的理论研究工作。文献 [17] 总结了许多 RHC 策略在管理学领域的应用研究。文献 [18] 尝试了应 用 RHC 策略进行有效的动态实时飞机进近排序研 究。 2.1 可选自由航线 在船舶可自由航行的海域进行路径优化,有无 数可选的子航线。因此,有必要合理并适当地简化 初始问题。使用传统的结构化的路径网络就是简化 此问题的一个简单方式。为了能充分利用好船舶可 自由航行的海域,我们可定义一个有限的可选方向 Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 73 的集合,并引入“时间段”概念,以此将船舶可自 由航行的海域转换成一个动态路径网络。这样就可 以实现在简化问题和使用的恰当算法找到的最优 方案之间来进行合理的折中。 图 1 船舶离散化自由海域航线 Fig. 1. Discretization of free navigation area for marine vessel’s path planning 使用一个有限的数值的集合来代表船舶可选 航线方向,以代替原来的无限的船舶航线方向的集 合,是离散化船舶可自由航海域的关键技术之一。 本文中,有限的船舶可选航线方向集合假定如下: Ω = [   50,,20,10,0 , direθ ] (1) 其中 dire θ 表示从当前船舶位置点到终点的直接方 向。所以,每次需要确定船舶可选航线方向时,仅 有 37 个值可供选用。 关于“时间段”的概念,我们假设有一个信息 发布系统将周期性地发布船舶航行海域的环境数 据,即台风七级风圈影响海域数据,这个周期叫“时 间段”。船舶使用最新更新的信息以便在下一个时 间段开始来优化剩余的航线。一个可选的航线包括 一系列和“时间段”相关的子航线。用于当前时间 段的子航线是由之前的优化来决定的。图 2 示例了 船舶离散化自由海域航线选择的过程。如果时间段 太长,而在 Ω 中的可选航线方向太少,则离散化 后的航线网络可能和传统的结构化的航线网络类 似,从而几乎无法包括全局最优的航线。另一方面, 如果时间段太短,而可选航线方向太多,将会有许 多可选的航线,从而导致找到最优航线的计算时间 大大增加,实时执行可能性就大大降低了。 2.2 航线优化的性能指标 如 2.1 节所述,一个可选的船舶航线包括一 系列和时间段相关的子航线。每个子航线(可选航 线的最后一个子航线除外)的航行时间都可以被看 做一个时间段的长度。然后,可选航线的总航行 时间由其所包括的子航线的数值来决定。因此, 虽然指标是总航行时间,但现在优化的基本变量 是子航线起始点和终点的坐标。这些基本变量和 一些重要参数如示意图 2 所示,(x,y)是一个点的 坐标,SAB 是 A 和 B 两点的距离,( , vϕ )是一个 点的洋流运动方向和速度,θ 表示船舶航向(本文 中所有方向都是相对于正北方向而定义的)。严格 讲,计算一个子线段的终点坐标,即(xB,yB),是 不可能的,因为(xB,yB)和( BB v,ϕ )互为前提条 件。但是,如果时间段定义得足够短,通常可以 假设子线段终点处的洋流运动参数与起始点处一 致,即( BB v,ϕ )=( AA v,ϕ )。根据此假设,计 算(xB,yB)就不再需要( BB v,ϕ )。 (a) 目标点 (b) 变量 图 2 计算船舶航行的一些重要参数 Fig. 2. Key parameters for calculating the movement of vessel θR 点 A 点 B ABS BAθ 坐标:( AA yx , ),( BB yx , ) 点 A 处的洋流参数: ( AA v,ϕ ) 点 B 处的洋流参数: ( BB v,ϕ ) 正北方向 船舶绝对速度: vAbs 洋流运动速度:vE 船 舶 相 对 速 度:vR θAbs θE 正北方向 船 舶 当 前 位 置 目 的 港 口 船舶自由 航行海域 台风七级风圈 或暗礁、其它船舶 可能的子航线 优化后船舶航线 Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 74 一条子航线的终点坐标(xB,yB),可通过如 下计算。 cosB A AB BAx x S θ= + (2) sinB A AB BAy y S θ= + () 其中 tsAbsAB TvS = , AbsBA θθ = (4) )cos(222 REREREAbs vvvvv θθ −++= (5) )/)sin((sin 1 AbsREERAbs vv θθθθ −+= − (6) AE ϕθ = , AE vv = W Av v= (7) 其中 Tts 是一个时间段的长度。坐标(xB,yB)随后 将被用于下个子航线的起始点。然后,基于气象 服务部门发布的洋流运动参数,新子航线的终点 坐标也可通过相同的方式进行计算。子航线的计 算可以一直进行到到达目的地港口为止。 对于一个可选航线最后的子航线,其终点就 是目的港,因此,(xB,yB) 是已知的,不再需要按 公式(1)和(2)计算。但是,最后子航线的航行 时间一般不再等于一个时间段的长度,而是需要 进行计算。假定图  中的 B 点就是目的港,则最后 子航线的航行时间可以通过如下计算: AbsABlast vSt /= (8) 其中 22 )()(),( BABABAAB yyxxPPdisS −+−== (9) ))(),(((2tan90 ABABBAAbs xxyya −−−== θθ (10) )/)sin((sin 1 REAbsEAbsR vv θθθθ −−= − (11) )cos()cos( EAbsERAbsRAbs vvv θθθθ −+−= (12) 而“ )(2tan ⋅a ”是计算四象限反正切的函数。 假设,排除最后的子航线,在一条可选航线 里有 N 条子航线。那相应的航行时间成本为: lastts tNTJ +=1 (1) 3 遗传算法(GA)加滑动地平线控制(RHC)策略 的混合算法 应用于船舶在海上自由海域航线实时优化问 题的现存方法有一个共同点:即在每一个时间段 ,优化的对象都是从当前子航线的结束点到目的 港的剩余航程。因此,这些方法总会面临两个常 见问题。一个问题是,在受台风风圈影响的自由 海域中,要在一个很短的时间段内完成对于长途 航线优化,通常很难实现。另一个问题是,在受 台风风圈影响的自由海域环境中,由于未来不确 定信息的原因,从目的港反推到当前子航线结束 点的优化方法所计算出的优化航线,其性能和可 靠性都会大大降低。 3.1 滑动地平线策略(RHC) 本文提出使用滑动地平线策略(RHC)来克服 现存方法中的上述问题。在每个阶段,即在每个 时间段,RHC 策略将只优化近期 N 个时间段内的航 线。这样,无论到目的港距离有多远,每次优化 所需的计算时间都有一个上限,这个计算时间的 上限主要取决于 N,即滑动地平线的长度。而且, 适当地选择滑动地平线的长度能像过滤器一样滤 除未来的不可靠信息。图 3 直观地显示了 RHC 策略 以及其相对于传统全局路径动态优化策略的潜在 优势。 图  RHC 策略与传统全局路径优化策略在船舶避台航线 优化问题的比较 Fig. . Comparison between RHC and traditional global dynamic optimization strategy in the route optimization problem for marine vessels under typhoon scenarios 公式(1)中所定义的 J1 是根据传统全局路 径动态优化策略而设计的优化性能指标。对于 RHC 策略,在每个时间段,仅需优化滑动地平线 内的航线,这只与 N 的大小有关,而与到目的港 的距离无关。因此,最小航行时间看起来对 RHC 策略没有什么意义。但事实是,在船舶可自由航 船 舶 1 和 船 舶 2 目 的 港 目 的 港 船舶 1(RHC 策略) 船舶 2(全路 径优化策略) (RHC 策略) (全路径优化策略) Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 75 行海域内,与 N 个时间段对应的多数可选航线都 是 Z 字型(因为如果每个时间段都可以随机选取方 向的话,正如后文中遗传算法的做法,那么所得 到的 N 个时间段的总航线很少可能会是直线),这 意味着捷径通常存在于 Z 字型路径的子航线之间 ,如图 4 所示。显然,即使原来计划的可选 Z 字型 航线是基于 N 个时间段的固定长度的滑动地平线 ,但走捷径后,除了实际航线长度缩短了,沿捷 径从 Z 字型航线起点到其终点的航行时间也会缩 短,不过具体缩短了多少并不确定。因此,通过 找捷径降低基于固定长度的滑动地平线的航行时 间仍是有意义的。事实上为了找到最短航线,传 统的全局路径动态优化策略也需要找捷径。 假定,在第 k 个时间段,在走捷径后,原来的 可选 Z 字型航线变为了 M(k)个时间段,其中 0 ≤ M(k) ≤ N 是一个实数,而 M(k)的小数部 分等于通过对应捷径上最后一个子航线的航行时 间除以 Tts,即有: ( ) ( ( )) /last tsM k floor M k t T− = (14) 其中“floor”使 M(k)四舍五入向负无穷大就近 取整。 图 4 N 个时间长度的 Z 字型航线和捷径 Fig. 4. Zigzag route of N-time-slice-long and shortcut 这样一来,RHC 策略采用的性能指标如下: J2(k)=M(k)Tts+Wterm(k) (15) 其中 Wterm(k)是一个与最后一条子航线和目的港 相关的终端加权函数。稍后我们会提供更多关于 Wterm(k)细节的论述。然后,在一个船舶可自由 航行海域环境下优化航线的 RHC 策略可以按如下 描述: 步骤 1:当船舶从起点出发,让 k =0。 步骤 2:接收气象部门发布的更新台风数据 ,以当前运动物体所在点 P(k)作为开始航线优 化的起始点,然后解决如下的最小化问题 )(min 2)|(),...,|2(),|1( kJkNkPkkPkkP +++ (16) 其中 P(k+i|k), 1, ,i N=  是与第 k 时刻所对应 的滑动地平线的一个可选 Z 字型航线中的第 i 子航 线的终点,这些点不能位于台风七级风圈等船舶 不 可 航 行 海 域 之 内 。 假 设 最 优 解 为 )]|(),...,|1([ ** kNkPkkP ++ ,而相应的 捷 径 为 )]|))(((),...,|1([ ** kkMceilkPkkP ff ++ ,其中“ceil”将 M(k)四舍五入向正无穷大就近 取整。 步骤 :当船舶到达 P(k),设定 )|1()1( * kkPkP f +=+ (17) 则船舶在当前时间段按照[P(k), P(k+1)]所确 定的子航线航行。 步骤 4:如果 P(k+1)不是目的港,让 k=k+1 ,然后去步骤 2。否则,计算结束。 3.2 基于遗传算法(GA)的优化程序 第 2 节里的数学模型提供了离散化的船舶可自 由航行海域的航线,求解最小化问题(16)首先就 需要在离散化的船舶可自由航行海域的航线中快 速高效地找到可行的船舶航线。众所周知,遗传算 法(GA)是一个大规模平行随机搜索的优化算法, 它能很好地解决最小化问题(16)。本文则提出基 于 RHC 策略的 GA 算法,用以在船舶可自由航行 海域的实时航线优化问题。在 GA 中的染色体是根 据原始的 Z 字型航线或捷径航线中的子航线端点 来构造的。由于 M(k)是一个不确定的有限实数, 不同的染色体可以有不同的长度。换句话说,本文 中的遗传算法采用的是可变长度的染色体设计。因 此,一个染色体的构造可以是:第一个基因记录 M (k)的值,“ceil(M(k))”是相应船舶航线中子 航线端点的数目,而接下来的基因按顺序记录下这 些端点的坐标,如图 5 所示。 基于一个染色体上记录的信息,相应的船舶航 线的 J2(k)可以根据公式(2)到公式(15)计算 出来。在第 k 时刻,假定 GA 中每代群体有 n 个染 色体,第 i 个染色体的 J2(k)值是 qi(k),而 qmax (k)和 qmin(k)表示 J2(k)在此代中的最大和最 小值。那么,第 i 个染色体的适应性由式(18)定 船舶 与有 N 个时间段的滑动地 平线对应的 Z 字型航线 与 Z 字型航线对应的捷 径 Z 字型航线中的子航线 的终点 捷径中的子航线的终点 Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 76 义。 图 5 GA 中的染色体的结构和含义 Fig. 5. The structure and meaning of chromosomes in GA ( ) ( ) ( ( ) ( )) / , ( ) ( ) max max min max min( ) ( ) ( ) ( ( ) ( )) / ( ), ( ) ( ) max max min max max min q k q k q k q k n q k q k iF k i q k q k q k q k n q k q k q k i − + − ≠ = − + − + =    (18) 在文献 [15] 中,提出了一些改进 GA 性能的有 效技术,如自适应交叉和变异概率,以及启发法 规则等。本文中的 GA 将采用这些技术,但限于篇 幅原因,不多赘述。 3.3 滑动地平线的长度和终端加权 滑动地平线的长度,即 N 的选择,是至关重 要的。RHC 策略每轮在线计算时间有一个上限, 主要取决于 N 的大小,通常可以通过模拟仿真进 行估计。因此,只要 RHC 策略中的时间段大于这 个计算时间的上限,则无论全程距离有多远,我 们总能保证算法的实时有效性。另外,适当地选 择滑动地平线的长度,可以像过滤器一样滤除长 远未来的不可靠信息。如果 N 太大,RHC 策略将 面临和现存的全局路径优化方法一样的关于实时 计算和台风等环境信息动态变化的问题。反之, 如果 N 太小,RHC 策略将变得短视,而且性能将 极大的降低。因此,需要合理选择 N 的大小,以 使 RHC 策略在在线计算时间和算法的稳定性能之 间得到最好的平衡。 但是,滑动地平线的性质让 RHC 策略在某种 意义上不可避免地变得短视,特别是和静态环境 中的传统全局路径优化方法相比较时。因此,我 们需要精心设计 J2(k)中的终端加权 Wterm(k) ,以求能有效降低 RHC 策略的短视性。 在控制工程领域,无终端加权的 RHC 策略曾 被广泛应用,然而,这种 RHC 策略经常导致系统 变得不稳定。为了解决这个问题,人们专门引进 了终端加权的技术 [15] 。现在,终端加权已成为保 证 RHC 策略的稳定性的一个关键技术。本文针对 于船舶在自由海域航行遭遇台风的情况 假设 Wterm(k)从 J2(k)中移除,即 0term k =W ( ) (19) 那么,如果目的地不在第 k 个时间段的范围内的话, J2(k)将没有目的港的信息。这种情况下,航线优 化的结果将导致一个随机的船舶航线,这很可能永 远不会指向目的港,如图 6 中的虚线所示。这是一 个由于在 J2(k)中没有终端加权所造成的不稳定 情况。 将目的地港信息加到 J2(k)中的一个简单方 法是使用如下的终端加权计算: AbsADlastTerm vPkPdiskW /)),(()( ..= (20) 其中 Plast(k)是可选航线中的最后一条子航线的终 点,PD.A.是目的港,而地上速率 vAbs 和函数“dis” 在 公式(9)和公式(12)中分别给出了。在公式 (20)中的 Wterm(k)能有效避免由于公式(19) 而产生的随机船舶航线,并可以在许多条件下让船 舶抵达目的港。但是,公式(20)中没有台风等不 可航海域的信息,有时会产生一个新的问题,即船 舶被困在一个海域内,RHC 策略几乎无法使船舶 驶离该海域,如图 6 上的点虚线所示。 图 6 各种终端加权的效果 Fig. 6. Effects under different terminal weighting terms 为了避免被困的现象,应在终端加权中加入一 些必要的障碍区域的信息。基本上,在 PD.A.(目的 地)和 Plast(k)(可选航线中的最后一条子航线的 终点)之间存在的七级台风风圈、附近航行船舶、 无目的港的 Wterm 可能被困的 Wterm 最 优 化 的 可能绕远的 Wterm 船舶 起始港 目的港 GA 中的一条染色体: (x1, y1) … M(k)的 数值 表达航线的第一 条子航线的终点 表达航线的最后一 条子航线的终点 M(k) (xceil(M(k)),yceil(M(k))) (x2, y2) (x1, y1) 4. (x2, y2) (x, y) (x4, y4) Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 77 岛礁等不可航行海域是导致船舶被困问题的主要 原因。为了方便考虑问题,下文中将这些在 PD.A. 和 Plast(k)之间的海域叫做 IW 不可航行海域,其 它海域叫 OW 不可航行海域。如果没有 IW 不可航 行海域,那么 Wterm(k)就如公式(20)。否则, 到 Plast(k)最近的 IW 不可航行海域(可能包括相 互重叠的一些台风风圈、岛礁等不可航海域)可以 用图 7 进行分析: 图 7 静态环境下第一种终端加权说明 Fig. 7. Illustration of the first terminal weighting term in static environment 其中 1θ 和 2θ 是显示在图 7 上的角度。Wterm(k)能 防止运动物体被困在某一区域,因为在潜在的受困 区域, ),max(/),min( 2121 θθθθ 接近于 1,将导致 很大的终端加权函数值。 然而,按照图 7 中的 Wterm(k)可能会导致较长 的航行时间。如图 6 中的双点虚线所示,为了避免 被困,运动物体可能会偏离直接航向 direθ 很远。为 了使 RHC 算法能更有效地找到最优航线,而不只 是可行的航线,需要对终端加权做更多修改。假设 Pprev(k)是遵照可选航线的点 Plast(k)而延伸出 来的点。如果 IW 不可航行海域的数目不为零,那 么一个更有效的终端加权可通过图 8 进行分析。图 中 θ 和 4θ 是图 8 所示的角度。  0θ > 表示可选路 径中的最后一条子线段在向外偏离目的地。相反,  0θ < 表示在向内偏离。无论哪种情况,都会通过 终端加权进行调整,保证的 Wterm(k)是一个很有 效的定义,其效果如图 6 中的实线所示。 然而,船舶在海上航行遭遇台风的情况,IW 不可航行海域是动态变化的,面积的大小、形状、 方向都随时发生改变甚至消失。为了方便处理问 题, 图 8 静态环境下第二种终端加权说明 Fig. 8. Illustration of the second terminal weighting term in static environment 不可航行海域的动态特征也可以简单地包括在 Wterm(k)。某种程度上,利用不可航行海域的动态 特征的一个简单方式就是考虑到 Plast(k)最近的 IW 不可航行海域的移动方向,即: AbsADlastTerm vPkPdiskW /)),(()1/||)1(()( ..4 ++= θθbr (21) 6 5 6IWsign signr g θ θ θ θ= − −( ) ( ) (22) 其中 θ 和 4θ 如图 8 所示, 0α > 是一个系数,而 0b > 是一个系数, 5θ , 6θ 和 IWθ 是相对于正北 方向的顺时针方向旋转角,如图 9 所示, IWθ 是到 Plast(k)最近的 IW 不可航行海域的移动方向, 0g > 是一个转向参数,而“sign”是一个符号函 数。 目前,本文的终端加权计算仅使用了到 Plast(k) 最近的 IW 不航行海域的信息。而关于如何使用其 它的 OW 不可航行海域信息可以进行进一步的研 究。例如,结合具体某一次航行任务中的不可航行 海域的空间分布特征和动态变化特性开展研究,但 已就不在本文的讨论范围内了。 图 9 动态环境下终端加权说明 Fig. 9. Illustration of terminal weighting term in Dynamic environment θ ( 此 图 中  0θ < ) 4θ 起始港 Pprev(k) PD.A. 船舶 已走航行过的航线 可能的计划航线 OW 不可航行海域 IW 不可航行海域 Plast(k) 5θ θ Pprev(k) PD.A. 船舶 可能的计 划航线 OW 不可 航行海域 IW 不可 航行海域 Plast(k) 6θ 1θ 2θ 起始港 Plast(k) PD.A. 船舶 已航行过的航线 可能的计划航线 OW 不可航行海域 IW 不可航行海域 Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 78 4 仿真结果 为了评估本文所提出的遗传算法(GA)加滑 动地平线控制(RHC)策略的混合算法对解决船舶 海上躲避台风的航线优化问题的有效性,本文建立 了的一套仿真系统来模拟船舶在海上自由海域航 行的环境,并以 2016 年 14 号“莫兰蒂”台风为例, 开展了一个仿真案例。图 10 是“莫兰蒂”台风的 实际路径,为了更准确的验证算法的效果,仿真案 例中的台风相关数据设置均以实际莫兰蒂台风的 实际数据为基础进行设定。 在航线设计方面,假定一艘集装箱船舶正由香 港港口以 90 度角的方向台湾航行,航速的设定以 参照大型集装箱船舶的运营航速,从船舶航行的安 全性及实验的数据误差两方面考虑,本文设定船舶 受台风影响的不可航行海域是以台风的七级风圈 影响范围作为数据基础。 台风所在位置的设计,假定仿真实验开始时, “莫兰蒂”台风正由东南向西北方向迅速移动,台 风中心在台湾高雄以南附近海域。仿真案例的具体 参数设定如下: 参数 船舶速 度 (km) 船舶航 向 (度) 台风移 动速度 (km) 台风移 动方向 (度) 台风七 级风圈 半径 (km) 数值 7 90 18 25 200 (备注:本文中所有方向都是相对于正北方向而定义的) 图 10 “莫兰蒂“台风实时路径 Fig. 10. The“MERANTI”Typhoon real-time path 为了对比算法的效果,本文特地引用了基于 GA 的传统全路径动态优化算法。为了便于区别, 此后使用了 RHC 策略的混合算法简称为 RHC,而 传统全路径优化的遗传算法为 CDO。由于它们都 使用了相同的 GA 作为航线优化算子,因此将 RHC 和 CDO 进行对比就很公平。关于 GA 优化算子的 更多细节本文不加以赘述。在仿真试验中,除非特 别指定,仿真试验中一个时间段是 10 分钟长,滑 动地平线的长度为 N=12(表 5 对应的试验除外, 因为表 5 专门研究不同 N 的影响),或 2 小时长, RHC 使用公式(21)中定义的终端加权计算 Wterm (k)(表 6 对应的试验除外,因为表 6 专门研究不 同终端加权的影响)。 在表 1 中定义了 4 个仿真案例,其自由航行海 域环境的复杂程度都各不相同,DD 表示从始发港 到目的地港的直接距离,UR 表示不可航行区域。 在案例 1-2 中,UR 是静态的,而在案例 -4 中, 其随着时间而产生变化,换言之,他们可以移动, 面积可以发生改变,消失或可能随机出现新的 UR。 此仿真主要用于比较 RHC 和 CDO 的在线计算时间 (OCT)及性能,即从始发港到目地港的实际航行 时间(AST)由于篇幅限制,我们仅采用并显示了 整个航行过程中 8 个时间段相关的结果。数值仿真 结果见表 2-表 6,其中对于每个静态案例,RHC 或 CDO 进行了 10 次仿真模拟,而对于每个动态方案 则进行了 200 次仿真模拟。 表 1 仿真试验中的 4 个仿真案例 Table 1. Four simluation cases 静态环境 动态环境 Case 1 Case 2 Case  Case 4 DD (km) 500 1000 500 1000 UR 的数目 1 6 1 6 虽然 RHC 主要应用于带台风等不可航海域的 动态方案,但在静态方案中它也需要能正常运转。 表 2 给出了根据不同算法案例 1-2 的仿真模拟结果, 从中可以看出,在两个案例中 CDO 实现了最好的优 化性能,即取得了最小的 AST。这是因为从理论上 说,在静态案例中像 CDO 这样的传统动态全路径优 化策略本就应该达到最好的性能。表 2 也显示了 RHC 的性能非常接近 CDO,这表示 RHC 在静态案 例中的运行效果不错。关于 OCT,RHC 很显然比 CDO 更有效率。由于一个时间段是 10 分钟长,我 们可以看出 RHC 在实时运行方面应该没有问题,而 在一些案例中 CDO 却很难实现在线运行。 我们主要的关注点是带台风等不可航海域的动 态案例,相关的仿真模拟结果见表 。在性能方面, Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 79 对于像案例  和 4 这样相对简单的案例,CDO 和 RHC 的 OCT 虽然相差不大,但可以看出 RHC 的性 能要优于 CDO。原因已在第 2 节和第  节中详细阐 述过了。因此,RHC 能提供比 CDO 更加可靠的实 时运算性能。 表 2 静态案例的仿真试验结果 Table 2. Simuluation results in static cases (单位: 分) CDO RHC 0.1590 1.0489 0.218 0.6079 平均 OCT 871.84 1628.28 871.97 161.40 平均 AST 0.5900 4.0975 0.68 0.8099 最大 OCT 872.1 164.5 872.52 1644.5 最大 AST 0.1590 1.0489 0.218 0.6079 如表 1 所示,从案例  到案例 4,DD 和 UR 的数值都增加了。那么是 DD 还是 UR 的数值更能 影响 RHC 的 OCT 呢?表 4 回答了这个问题,其中 DD 在[500,1000]之间变更,UR 数值在[1,6]之间变 更,而所有案例都是动态的。从表 4 中可以看出, RHC 的 OCT 主要取决于 UR(因为 UR 的数值极 大地影响了寻找可选飞行路径和计算终端加权的 GA 优化算子的计算负荷),和 DD 关系不大(因为, 对于 RHC,不是由 DD 而是由 N 来决定一条可选 飞行路径的最长飞行时间)。 表 3 动态态案例的仿真试验结果 Table 3. Simuluation results in dynamic cases (单位: 分) CDO RHC Case  Case 4 Case  Case 4 平均 OCT 0.1206 1.1844 0.125 0.4816 平均 AST 928.21 164.0 928.12 168.61 最大 OCT 0.581 4.260 0.6449 0.6908 最大 AST 928.6 1866.8 928.45 1757.65 表 4 DD 和 UR’s 对 RHC 性能的影响 Table 4. Influence of DD and UR on the OCT of the RHC OCT (分) DD=500 (km) DD=1000 (km) Ave. Max. Ave. Max. 1 个 UR 0.125 0.6449 0.77 0.7072 6 个 UR 0.495 0.7455 0.4816 0.6908 表 4 很清楚地表明,应合理地选择 N,即滑动 地平线的长度。如果 N 太小,RHC 的优化性能就 差,如表 5 中 N=1 和 N= 的案例。而如果 N 太大, OCT 增加了,但性能并没有进一步提高。相反,在 动态案例中,N 即太大时 RHC 的优化性能还可能 降低,如表 5 中 N=9 时所示。 表 5 滑动地平线的长度 N 对 RHC 性能的影响 Table 5. Influence of N on the RHC (单位:分) 静态环境 动态环境 Case 1 Case 2 Case  Case 4 N=1 OCT 0.1045 0.1174 0.0920 0.1061 AST 880.8 541.76 928.89 175.47 N= OCT 0.160 0.2445 0.1618 0.182 AST 871.71 1717.02 929.20 1644.8 N=6 OCT 0.218 0.6079 0.125 0.4816 AST 871.97 161.40 928.12 168.61 N=9 OCT 0.5799 1.290 0.515 1.0750 AST 871.91 1628.5 928.19 1640.9 表 6 显示了 RHC 中终端加权对优化性能的影 响。因为公式(19)中定义的 Wterm(k)使算法变 得不稳定,所以表 6 中没有给出相关的结果。基本 上,在 OCT 保持在相同水平时,可以看出在分别 使用公式(20),公式(21)中定义的 Wterm(k) 后,RHC 的性能一步步地得到了改善。原因已在 . 节中详细阐述了。 表 6 终端加权对 RHC 性能的影响 Table 6. Influence of terminal weighting on the RHC Wterm(k) 静态环境 动态环境 Case 1 Case 2 Case  Case 4 公式 (20) OCT 2.8102 4.9121 2.8994 .845 AST 904.67 164.44 978.54 1647.80 公式 (21) OCT 2.5675 4.8498 2.490 .8419 AST 871.97 161.40 928.12 168.61 下面的图 11 截取了带有单个台风动态不可航行 区域的仿真案例中,船舶避台航线规划过程及结果, 共有 8 张图片展示船舶避台航行的航线优化的过程。 5 结论 船舶海上躲避台风航线优化问题与传统的路径 优化问题不同,船舶在海上躲避台风时的不可航海 域是实时变化的,不可航海域可以是固定的岛礁、 浅滩,也可以是动态变化的台风风圈和其它船舶在 航海域等等,因此,船舶海上躲避台风航线优化问 题要求针对动态变化的不可航海域进行实时航线优 化。本文提出了一种基于遗传算法(GA)和滑动地 平线控制(RHC)策略的混和算法来解决此问题。 Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 80 为此目的,本文首先建立了带动态变化不可航海域 备注:图中的横纵坐标代表仿真实验航行海域范围,单 位:km) 图 11 RHC 仿真案例 Fig. 11. An example of RHC simulation 的自由海域路径优化问题的数学模型后,然后详细 介绍了混和算法的设计思路,全面探讨了此算法在 如何选择 RHC 的滑动地平线的长度,以及如何使 用终端加权等主要技术问题。仿真试验结果显示, 在缺少台风等动态变化环境因素时,混和算法可以 和现存算法达到一样好的性能,而在带有台风等动 态变化的环境中,混和算法的性能(不论求解质量, 还是求解速度)则优于现存算法,混合算法在解决 船舶海上躲避台风航线优化问题能够做到在保证 船舶安全航行的前提下,更好的实现船舶运营的成 本效益。 参考文献: [1] 史培军. 论灾害研究的理论与实践. 南京大学学报 (自然科学版). 自然灾害研究专辑, 1991(11): 7-42. [2] 张兰生,史培军,方修琦. 我国农业自然灾害灾情分 析 . 北京师范大学学报 ( 自然科学版 ). 1990(0): 94-100. [] 梁必骐,梁经萍,温之平. 中国台风灾害及其影响的 研究. 自然灾害学报. 1995(01): 84-91. [4] 夏施敏. 台风对沿海船舶的影响及防抗措施研究. 中 国水运. 2008(09): 50-51. [5] CORMEN T H, LEISERSON C E, RIVEST R L, et al. Introduction to Algorithms, Second Edition. Cambridge: MIT Press and McGraw-Hill, 2001. [6] 韩雪, 刘英舜, 郭唐仪. 城市轨道交通网络中断下的 有效路径搜索模型. 公路交通科技, 2015, 2(10): 97-101. HAN Xue, LIU Ying-shun, GUO Tang-yi. A Valid Path Searching Model for Interrupted Urban Rail Transit Network. Journal of Highway and Transportation Research and Development, 2015, 2(10): 97-101. [7] WICKENS C D, MAVOR AS, PARASURAMAN R,et al. Airspace System Integration-The Concept of Free Flight. Washington: National Academy Press, 1998. [8] KAHNE S. Research Issues in the Transition to Free Flight. Annual Reviews in Control, 2000, 24(1): 21-29. [9] MCLEAN D. Operational Requirements for Fully Automatic Flight. Aircraft Engineering and Aerospace Technology, 200, 75(6): 570-574. [10] 胡晓明, 李万莉, 朱为国, 等. 半挂液罐车的避障运 动仿真分析. 公路交通科技 201, 0(7): 151-158. HU Xiao-ming, LI Wan-li, ZHU Wei-guo, et al. Simulation of Obstacle Avoidance Motion of Liquid Tank Semi-trailer. Journal of Highway and Transportation Research and Development, 201, 0(7): 151-158. [11] 刘大刚, 吴兆麟, 李国平. 船舶绕避热带气旋技术方 法述评. 航海技术, 2004(5): 9-11. [12] 刘涛,刘大刚. 目标圆方法在船舶导航中的应用. 气 象科技. 2007, 5(6): 867-87. [1] 林明智, 庄丽, 张楠, 等. 船舶避台技术方法新探. 气 象科技. 1999(4): 49-5. [14] HU X B, WU S F, JIANG J. On-Line Free-Flight Path Optimization Based on Improved Genetic Algorithms. Engineering Applications of Artificial Intelligence, 2004, 17(8): 897-907. [15] CLARKE D W. Advances in Model-based Predictive Control. NEW YORK: Oxford University Press, 1994. [16] BART D S, TON V D B. Model Predictive Control for Max-Plus-Linear Discrete Event Systems. Automatica, 2001, 7(7): 1049-1056. [17] CHAND S, HSU V N, SETHI S. Forecast, Solution, and Rolling Horizons in Operations Management Problems: A Classified Bibliography. Manufacturing & Service Operations Management, 2002, 4(4): 25-4. [18] HU X B, CHEN W H. Receding Horizon Control for Aircraft Arrival Sequencing and Scheduling. IEEE Transaction on Intelligent Transportation System, 2005, 6(2): 189-197. Journal of Risk Analysis and Crisis Response, Vol. 7, No. 2 (July 2017) 72–81 ___________________________________________________________________________________________________________ 81 http://en.wikipedia.org/wiki/Introduction_to_Algorithms A Study on Marine Vessels’ Path Optimization under Typhoon Scenarios 1 引言 2 船舶在海上自由海域航行时的路径实时优化问题 2.1 可选自由航线 2.2 航线优化的性能指标 3.2 基于遗传算法(GA)的优化程序 3.3 滑动地平线的长度和终端加权 (单位:分) Wterm(k) [12] 刘涛,刘大刚. 目标圆方法在船舶导航中的应用. 气象科技. 2007, 35(6): 867-873.