1、收 稿日期： 2010-10-18； 修回日期： 2010-11-26基金项目： 国家创新研究群体科学基金（ 5092100）作者简介： 姜华林（ 1986-），男，河南省人，硕士研究生，主要从事海岸和近海工程研究 。Biography： JIANG Hua-lin（ 1986-）， male， master student.基于蒙特卡洛方法的海上搜寻区域确定模型姜华林1，孙昭晨1，李 黎2，刘 波2，吴俊萍2（ 1.大连理工大学 海岸和近海工程国家重点实验室，大连 116023；2.天津市西戈数字科技有限公司，天津 300384）摘 要： 为了提高海上搜寻与救助的成功率，采用 FVCOM 模
2、型进行遇难海域表面流场的数值预报，并与实测结果进行比较 。采用经验公式法确定不同漂浮物种类风致漂移速度与风速的关系 。引入蒙特卡洛方法，通过量化影响漂移的不确定性因素，确定漂浮物位置随时间变化的概率密度函数，并以统计区域作为搜寻区域，编制程序建立搜寻区域确定模型 。模型实例验证效果较好，可为搜救工作的开展提供依据 。关键词： 搜寻区域；风致漂移；蒙特卡洛方法；数值模拟中图分类号： P 73123； O 242.1 文献标识码： A 文章编号： 1005-8443（ 2011） 04-0285-06海上搜救是通过量化一系列不确定因素，包括搜寻目标最后的已知位置 、搜救对象类型以及影响漂移的风况
3、、海况等，计算搜寻区域以及搜寻区域随时间的演变，并且快速调动搜救力量在确定的搜寻区域进行救援 。海上搜寻区域的确定很大程度上决定着搜救行动的成败，因此建立快速的海上搜寻区域确定模型具有重要意义 。本文旨在建立一个操作性较好 、计算快速的海上搜寻区域确定模型，为搜救工作的开展提供依据 。海运业发达的国家对海上搜寻理论进行了深入的研究，开发了一些搜救的计算软件系统，如加拿大的CANSARP 系统 、美国的 HACSALV 系统，以及英国的 SARIS 系统 。我国学者近年来在这方面也进行了一些研究，文献 1-2介绍了海上目标漂移区域计算方法，但这些计算基于实测海流和水文图表，在一定程度上缺乏可操作
4、性；文献 3研究了遇难船舶漂移，未涉及落水人员 、救生艇筏等目标漂移；文献 4在基于流场数值模拟的基础上研究了遇险目标漂移，但对风致漂移作用不够明确 。漂浮物在海面的运动是风 、浪 、流联合作用的结果，为了计算方便，模型中引入简单化处理，同时也将产生误差，包括风 、流的观测或计算误差 、风致漂移误差 、物体最后已知位置误差等 。本文在表面流场数值预报的基础上，采用概率公式处理各种漂移误差，并且引入蒙特卡洛计算方法，在仅考虑流场和风的作用下，建立针对遇难船舶 、落水人员及救生艇筏等的海上搜寻区域确定模型 。1 影响漂浮物漂移的主要环境因素海面物体所受作用力可用下式计算（ m+m）dVdt=F （
5、 1）式中： V 为物体速度； F 为所受合力；附加质量 m是由船体周围水粒子的加速度引起的 。研究表明救生筏在强风（ 20 m/s）条件下大约 20 s 能达到极限速度，这就意味着小体积物体（通常长度小于 30 m）在海面加速很快，因此假设漂浮物在单个时间步长内具有无限加速度和均匀速度 。水 道 港 口Journal of Waterway and Harbor第 32 卷第 4 期2011 年 8 月Vol.32 No.4Aug. 2011水 道 港 口 第 32 卷第 4 期1.1 风和风致漂移风致漂移是由漂浮物暴露于空气中的结构受到风力作用引起的 。由于漂移物体形状不对称，存在一个侧向
6、分力使物体的风致漂移方向与下风方向存在一定夹角，而同一角度出现在下风方向左右的概率相同 。文献 5-6对多种不同的漂移物体（落水人员 、救生艇筏 、遇难船舶等）的风致漂移进行了实验研究，给出了经验系数 。风压作用可分解为顺风方向分量 Ld（ DWL）和侧风方向的分量Lc（ CWL）（图 1），前人研究表明，顺风分量和侧风分量与海面以上10 m 风速之间存在较好的线性关系，利用最小二乘法即可确定该线性关系 。Ld=adW10+bd（ 2）式中： ad和 bd为通过拟合实验数据所得线性系数，与漂浮物的类型有关； Ld为估算的最可能的顺风分量 。侧风分量 Lc也有类似的计算方法 。1.2 流场和随流
7、漂移随流漂移是漂移物体在周围水体的作用下的运动，不考虑风的情况下，假设漂移物体的运动速度与表面水流运动一致，当水流速度改变时，物体的运动也随之改变（无限加速度） 。本文采用 FVCOM 模型进行表面流场的数值预报 7。FVCOM 模型的全称是非结构网格有限体积法海洋模式 。采用非结构三角形网格，能够更精确的拟合复杂曲率的岸界；采用有限体积法，在整个计算区域内更好的保证了各物理量的守恒 。该水动力模型采用垂向静压假定，在水平方向上采用非结构化三角形网格，垂直方向采用 坐标变换，沿重力方向分层 。考虑到快速计算，模型中不考虑温度 、盐度 、密度等影响 。动量方程和连续方程如式（ 3） 式（ 5）坠
8、uD坠t+坠u2D坠x+坠uvD坠y+坠u坠-fvD=-gD坠坠x-gDe坠坠xD0乙d乙 乙+坠D坠x +1D坠坠Km坠u坠 乙+DFx （ 3）坠vD坠t+坠v2D坠y+坠uvD坠x+坠v坠-fuD=-gD坠坠y-gDe坠坠yD0乙d 乙+坠D坠y +1D坠坠Km坠v坠 乙+DFy （ 4）坠坠t+坠Du坠x+坠Dv坠y+坠坠=0 （ 5）2 漂移模型的建立与软件实现经过上述的假设和简化处理，漂移模型简化为随流漂移和风致漂移的叠加 。x（ t） -x0=t0乙V（ t） dt=t0乙L（ t） +uW（ t ） dt （ 6）采用二阶龙哥库塔法进行计算 。然而计算过程中存在着很多的不确定因
9、素和误差，将导致计算的漂移轨迹与真实漂移轨迹随着时间推移不断偏离，因此该研究引入蒙特卡洛算法 8。在漂移模型中，其基本思想是将漂浮物定义为单个粒子，粒子的移动轨迹是一个满足一定概率分布的随机过程，在同一规则下，大量复制这些粒子，并设这些粒子的运动互不相关，那么这些粒子在一段时间后的统计区域即是漂浮物漂移的概率区域 。本文假设漂浮物在海面上的位置满足马尔科夫过程，即漂浮物未来状态的概率密度函数仅依赖于漂浮物当前所处的外力状态，而与物体到达当前状态的特殊过程无关 。2862011 年 8 月p（ xt+1xt， xt-1， xt-2， .， x1） =p（ xt+1xt） （ 7）随机漂移过程用下
10、式表达dx=V（ x， t） dt+d （ 8）式中： x 为漂浮物的位置； V 为漂浮物在外力作用下的速度函数； d 为具有零均值和已知方差的随机扰动项，代表漂浮物自身特性以及外力（风和流）的不确定性 。水流扰动对搜寻区域的确定影响较小，随机扰动项 d仅考虑风致漂移的作用 。2.1 风致漂移系数扰动漂浮物的漂移方向（风向的左 、右）是平均分布的，并且假定在漂移过程中是固定不变的，即假设最初向顺风的左侧漂移，漂浮物在模拟的过程中将一直沿左侧漂移 9。如前所述，风致漂移与风速之间存在线性关系，由于存在偶然因素，线性关系的斜率和截距分别加上一个满足正态分布 N（ 0， ）的随机扰动 。an=a+n
11、/20bn=b+n/2当 W10=10 m/s 时，令斜率和截距具有相同的扰动，则L=（ a+n/20） W10+（ b+n/2） （ 9）2.2 风的扰动假设扰动风场满足圆正态分布un（ un， vn） N（ 0， W）Wn= W10+un整体的风致漂移如下式Ln=anWn+bn， n=1， .， N （ 10）2.3 搜寻区域确定模型的软件实现根据前文的分析，在只考虑风场和流场对漂浮物影响的情况下，对遇难海域进行流场的数值预报，用Matlab 语言编写一个演示性的海上搜寻区域确定软件 。只要给出船舶失事地点的初始位置 、失事时间 、预计漂流时间 、失事海域的风况及船舶类型等资料，就可预报在
12、此时间段内失事船舶的漂移轨迹和最终的搜寻区域 。软件的主界面如图 6 所示 。3 算例验证本案例船舶为 “辽盖渔 5088”轮 10，当时正在进行海上测量任务，后来由于风浪较大，锚链被拉断，且船舶失去动力 。“辽盖渔 5088”于 2004 年 10 月 25 日 3 时发出求救信号，通信在初始要求救助后中断， 26 日 7时 29 分，被 “辽盖渔 4028”发现 。事故位置为 4029.7N/12156.1E；发现位置为 4023.1N/12148.1E。3.1 流场计算3.1.1 模型计算区域及网格设置根据船舶遇难位置，确定建模计算区域为辽东湾：南起 39.5N、北至 41N、西起 11
13、9.5E、东至 122.3E，利用海图获取计算海域的岸线和水深数据 。模型采用无结构的三角形网格系统，计算区域内共有三角形网格节点 7 180 个，三角形单元 13 807 个，时间步长为 6.0 s。3.1.2 边界条件潮流模式开边界的强迫潮位值，通过给出秦皇岛（ 3940N， 12121E）和长兴岛（ 3939N， 12128E）的调和常数，以及模拟的时间，经过潮位调和计算给出，开边界其他网格点的潮位由这两点经线性插值得到 。3.1.3 计算及结果依据上述的计算模式进行运算，预报事故海域内的潮位 、流速等信息，本算例计算时间为 2004 年 10 月20 日 30 日 。图 2、图 3 给
14、出了葫芦岛与营口 2 个验潮站的潮位实测值与计算值对比，图 4 给出了测点（ 4020.5N/12041.25E）表层流速实测值与计算值的对比，总体来看吻合情况较好 。图 5 给出了 2004 年 10 月 25姜华林，等 基于蒙特卡洛方法的海上搜寻区域确定模型 287水 道 港 口 第 32 卷第 4 期日 18 时的流场图 。3.2 风致漂移本案例遇难船舶类型为渔船，根据文献 6，风致漂移经验系数为 ad=1.87， bd=0， ac=2.09， bc=0，故该事故图 4 测站表层流速实测与计算值对比图Fig.4 Comparison between simulated and obser
15、ved data on the surface layer of Liaodong Bay4-a 流速 4-b 流向时间（ h）2882011 年 8 月船舶的风致漂移经验公式为LDWL=（ 1.87+n/20） Wn+n/2 nN（ 0， 8.3） （ 11）LCWL=（ 2.09+n/20） Wn+n/2 nN（ 0， 6.7） （ 12）图 6 给出了由于随机扰动产生的风致漂移顺风和侧风分量集合 O（ 100） 。3.3 计算结果与分析根据前文计算的风和流场信息，将漂移模型的计算结果与实测值相比较 。模型计算结果与真实情况基本一致，漂移经过 28.5 h，最终船舶发现位置在计算所得的搜寻
16、区域，距离搜寻中心约 3 km，说明模型计算基本能够满足海上搜救的需要 。4 结语本文对海上搜救行动中搜寻区域的确定进行了数值模型研究，在仅考虑风和流对漂移影响的前提下，首先利用 FVCOM 模型对流场进行数值预报，利用经验公式的方法研究风与风致漂移的关系，采用蒙特卡洛算法经过多次模拟计算得出搜寻区域以及搜寻区域随时间的演变 。实际案例验证结果表明，该模型计算快速 、准确，操作性较好，能够为搜救工作提供依据 。但模型中并未考虑波浪对漂移的影响，还需进行进一步分析和研究 。参考文献： 1于卫红，贾传荧 .海上搜救中搜寻区域确定方法研究 J .中国航海， 2006（ 2）： 34-37.YU W
17、H， JIA C Y.Methods of determining search area for SAR at sea J .Navigation of China， 2006（ 2）： 34-37.姜华林，等 基于蒙特卡洛方法的海上搜寻区域确定模型图 6 事故船舶风致漂移顺风分量和侧风分量Fig.6 Downwind and crosswind leeway component of accident ship6-a 风致漂移顺风分量020406080100-20风致漂移顺风分量（cm/s）2520151050海面上 10 m 处风速（ m/s）289水 道 港 口 第 32 卷第 4 期
18、Determining maritime search area model based on Monte Carlo methodJIANG Hua-lin1， SUN Zhao-chen1， LI Li2， LIU Bo2， WU Jun-ping2（ 1.State Key Laboratory of Coastal and Offshore Engineering， Dalian University of Technology， Dalian 116023，China； 2.Tianjin CG Digital Science and Technology Co. Ltd.， Tia
19、njin 300384， China）Abstract： In order to improve the success rate of maritime search and rescue， FVCOM model was adoptedfor numerical forecasting of the surface flow field in the accident sea area， and the calculated results werecompared with experimental results. Empirical formula method for the re
20、lation between wind and the motion ofthe drifting object was employed. A Monte Carlo technique was introduced to generate an ensemble thataccounting for the uncertainties. The ensemble yielded an estimate of the time-evolving probability densityfunction of the search object， and its envelope defined
21、 the search area. This model was validated using actualcase. The results show that it agrees well with the actual data， and it can provide the basis for carrying out searchand rescue operations.Key words： search area； drift associated with wind forces； Monte Carlo technique； numerical simulation 2胡志
22、武，毕曙光，耿鹤军 .海上搜寻区域确定的计算机辅助方法 J .上海海事大学学报， 2008， 29（ 2）： 11-16.HU Z W， BI S G， GENG H J.Computer-aided method for determining maritime search area J .Journal of Shanghai Maritime Univer-sity， 2008， 29（ 2）： 11-16. 3胡志武，张秋荣，顾维国 .遇难船舶漂流轨迹预测技术 J 航海技术， 2007（ 3）： 18-21.HU Z W， ZHANG Q R， GU W G. Prediction
23、on technique of drift course of distressed ships J .Navigation Technology， 2007（ 3）：18-21 4马文耀，陈达森，毕修颖 .基于流场数值模拟的遇险目标漂移计算研究 J .中国航海， 2009（ 2）： 45-48.MA W Y， CHEN D S， BI X Y.Study of distress target drift based on numerical simulation of flow field J .Navigation of China，2009（ 2）： 45-48. 5 Allen A，
24、Plourde J V.Review of leeway： field experiments and implementation R .USA： US Coast Guard Research and Develop-ment Center， 1999. 6 Allen A.Leeway divergence R .USA： US Coast Guard Research and Development Center， 2005. 7刘广生 .非结构网格下河流对近岸水环境影响的数值模拟 D .大连：大连理工大学， 2010. 8邱平 .海上最佳搜寻区域与搜寻方法的研究 D .大连：大连海事
25、大学， 2006. 9 Breivik O， Allen A.An operational search and rescue model for Norwegian Sea and the North Sea J . Journal of Marine Systems，2008， 69： 99-113. 10欧阳 .海上失控船舶漂移模型研究 D .大连：大连海事大学， 2008.福建 “十二五 ”将投资四千亿打造综合交通运输体系据 福建省交通运输 “十二五 ”发展规划 ，“十二五 ”期间，福建省将以 3 倍于 “十一五 ”规划投资的力度，投资 4 000 亿元，着力推进综合交通运输体系建设，实现到 2015 年，福州港 、湄洲湾港 、厦门港三大港口群基本形成，对接两洲 、拓展中西部 、服务周边地区发展新的对外开放综合通道全面构建，具有区域性集散效应的东南沿海物流中心和国家综合运输枢纽初步建成 。到 2015 年末，全省沿海港口通过能力达到 5.5 亿 t，比“十一五 ”末增加 2.3 亿 t，适应度达到 1： 1；沿海港口深水泊位约 200 个，深水泊位比率约 40%，沿海港口专业化泊位通过能力比例达到 85%。（殷缶，梅深）290
