福建近海漂移物漂移轨迹分析及模拟

曾银东,郭民权,方许闻,江毓武,郑祥靖

(1.中国海洋大学海洋与大气学院,山东 青岛 266100;2.福建省海洋预报台,福建 福州 350003;3.厦门大学海洋与地球学院,福建 厦门 361102)

台湾海峡及其邻近海域由于自然条件及生产活动的原因,船只事故时有发生,造成人员落水、海上溢油,加上近年频发的有毒赤潮,这些事件对近海海洋经济发展和环境安全造成严重危害,并导致巨大的经济损失.因此,在台湾海峡及其邻近海域建立漂移轨迹预报模式,对落水人员、溢油、赤潮等漂浮物的运动轨迹进行及时可靠的模拟预测,能为开展有效地搜救、防控、拦截以及减轻生态灾害损失、保障海洋生态环境安全等提供决策支持.

影响漂移运动的海表动力环境主要包括海面的风、浪、流等因素,漂移物主要随海流运动.浮子漂移实验能够捕捉研究海域表面海流的特征,基于GPS的浮子具有高的精度(2～50 m)与采样频率等特点,常被运用于近岸海域的浮子漂移实验[1].Lin等[2]于2018年1月在大亚湾开展了SUCE(surface current experiment)浮子漂移实验,通过比较浮子漂移实验计算的流速与广东沿海已有的再分析数据,认为浮子漂移实验结果能够反映表层海流的特征.2020年Gao等[3]分析部署在大亚湾海域的120个浮子的实测资料,认为漂流浮子总体上能够捕捉到表面流(潮流与余流)特征.随着互联网以及通信定位技术的发展,漂流浮子由于成本更低,体积更小,已成为观测近岸海域表层海流的重要工具.

目前国外海事发达国家的目标漂移轨迹预测模型基本已形成较成熟的业务化预报体系,如英国的SARIS[4]、美国的HACSALV[5]、挪威的LEEWAY[6]等.国内也开展了不少相关的研究,如胡志武等[7]基于实测海流和水文图表,以遇难船舶为对象,研究了在风、浪、流共同作用下船舶的漂移轨迹,认为考虑波浪作用后并未提高预测精度.近年来,部分研究均采用数值模式,如肖文军等[8]利用中尺度气象研究与预报(weather research and forecasting,WRF)模式、无结构网格近岸海洋模式(finite-volume coastal and ocean model,FVCOM)和基于蒙特卡洛随机统计理论的海上搜救目标漂移轨迹模式LEEWAY,建立了长江口及邻近海域海上搜救物漂移轨迹预测模型系统.尽管大部分研究工作都不考虑波浪对物体的作用,但是部分学者认为波浪的作用不可忽略.徐江玲等[9]基于海上大型浮标(直径10 m)脱离锚定(跑位)后的漂移轨迹实测数据,对比了只考虑风和流的漂移轨迹和增加波浪作用后的漂移轨迹,认为波浪作用引起的漂移速度随波陡或浮在海面上的体积增大而增加.

目前关于台湾海峡浮子漂移实验研究的报道较少,且局限于平潭等部分海域[10],本研究通过覆盖福建近海的浮子漂移实验,分析了研究海域的表面流特征,讨论了FES潮流模型在福建近海的适用性;并采用已在台湾海峡业务化运行的具有港湾分辨率的三维海洋动力模型,在对港湾海流模拟的基础上,结合拉格朗日粒子追踪方法,建立了福建近海海上目标物漂移轨迹预测模式,利用浮子漂移实验对漂移模型在福建近海的适用性进行了较全面的验证,并分析了不同海域漂移轨迹预测误差来源.

1 漂移物漂移轨迹预测模式

1.1 拉格朗日粒子追踪法

漂移是一个连续变化的过程,但在实际处理中认为漂移物在某时间段内的速度不变.假定漂移物的初始位置为S0,经过一定的时间步长Δt后漂移到新的位置S,根据拉格朗日粒子追踪法可得[11]:

(1)

式中:t0为初始时刻;Vt为t时刻漂移物的漂移速度,是各环境因子综合作用所产生的速度;Δα为湍流涡动造成的随机运动距离(α代表x、y或z方向).

本研究认为海上漂移物受海面风和流的作用,即Vt主要由风导漂移和海表流场两个分速度构成,风与流产生的漂移物漂流速度可采用下式计算[10]:

Vt=k×Vwind+Vcurrent.

(2)

式中:Vwind为海面风速;Vcurrent为海表流速;k为风导系数,其值受漂移物种类、漂移物沉浸状态等因素影响,当漂移物为船舶和浮标时,k=0.01～0.05,对于浮子本研究取k=0.02.

漂移物在湍流涡动作用下的随机运动距离可以由下式计算[12]:

(3)

式中,R为[-1,1]区间的均匀分布随机数,Kα为α方向上的扰动系数.

1.2 流场模型

ROMS(regional oceanic modeling system)模型是一个开源的三维非线性斜压海洋模型,具有自由表面、地形跟踪、高分辨率等特点[13],此模型可以预测不同尺度的运动.本研究在已研发的大小网格嵌套精细化三维温、盐、流数值预报模式[14]基础上,采用最新的水深、地形等基础数据,利用全球模型提供边界条件,采用大小网格嵌套的曲线正交网格方案,其中西北太平洋区为大网格区,厦门湾、兴化湾等重点港湾为小网格区;通过不断优化完善,建立了具有港湾分辨率的重点海域三维海洋动力模型.

1.2.1 嵌套模型与网格分布

本模型用大小区域单向嵌套方法进行,大区域网格计算的结果进行插值后作为小网格的边界条件;小区域网格的环流参数,如水位、温、盐、环流由大网格提供,潮汐参数(潮位、潮流)由调和常数计算.

水平网格(图1):大网格区为95°～147.9° E,9° S～44.01° N.大区域网格在开边界处网格距为15 km,在台湾海峡西岸网格距为2 km.小区域网格在厦门湾口、兴化湾等海域达到100 m的分辨率.在模型计算中,与大洋资料相一致,取最大水深为5 500 m.

(a)大网格区,一个网格表示15个实际网格;(b)厦门湾网格,一个网格表示60个实际网格.图1 模型计算区域和水深分布Fig.1 Model calculation area and water depth distribution

垂向网格:大区域网格垂向均取30层,最小水深为10 m.小区域网格垂向取10层,并考虑潮汐漫滩过程,临界水深为0.5 m.本模型采用对数形式进行垂向分层,不同水深区域和不同层次的间隔不同.在海表分层较密,以反映海表混合层及温、盐跃层的复杂直结构;在海底分层较疏.

1.2.2 开边界控制条件

大网格水动力的开边界采用欧洲MyOcean的业务化产品,包括海表高度、温度、盐度和流速(https:∥www.copernicus.eu/en/myocean).海表的净热通量、水气通量和动量通量则采用福建省海洋预报台每日业务化运行的WRF模式的预报产品进行计算.模型计算中的江河输入均以《中国海湾志》[15]的主要河流各月多年平均流量代入模型计算,入海盐度取0.开边界潮汐所引起的水位与潮流过程采用FES2014潮汐模型(https:∥www.aviso.altimetry.fr/en/data/products/auxiliary-products/global-tide-fes.html)的天文潮调和常数,包含M2、S2、K1、O1、N2、P1、K2、Q1、M1、J1、2N2、Mu2、Nu2、L2和T2,共计15个分潮.模型设定详细方案可参考文献[14].小网格开边界来源于大网格的滤潮计算结果,以及潮汐调和常数所计算的潮位及潮流.

已有研究[14]利用福建近岸海洋监测系统数据评估模型性能,并以平均偏差、相关系数、均方根误差等指标比较了模型结果与观测数据,认为模型能较准确地模拟台湾海峡海流状态,本研究不再详述.

2 结果与讨论

2.1 浮子漂移实验与时段选取

本研究中浮子漂移实验的浮子内置GPS定位模块(含GSM通讯模块),能定时将浮子的运动轨迹坐标信息发回服务器,内置锂电池的电量理论上在连续运动状态可供5 d以上续航时间,外部采用密闭性塑料瓶安装保证不进水.同时会在浮子内部装有适量砂石,使浮子在海面上的部分不超过2 cm,以减少风对浮子运动的直接影响.浮子结构见图2.

图2 浮子结构图(单位:cm)Fig.2 Structure diagram of float (unit:cm)

2019年6月至2021年6月,福建省海洋预报台在福建近海共开展了17个航次调查,期间共释放245个漂移浮子,获取漂移点信息42万个,累计1 600多天,累计漂移距离4.2万km.

浮子在漂移过程中存在搁浅、信号中断、定位数据异常、搁浅后重新漂移等情况,因此通过判读漂移数据的连续性,选取连续漂移时段>24 h的漂移过程,最终获得273个漂移过程,图3为遴选出的部分浮子漂移轨迹图.可以看出:湾内的浮子(如厦门湾与兴化湾海域)主要受潮流作用,在湾内呈往复运动,并有流出湾外的趋势.湾外的浮子则主要沿着福建海岸运动;在春、秋季东北风盛行时节,在各河口淡水形成的浮力推动下[16],沿福建近海南下的沿岸流势力较强,因此浮子主要是沿岸向南运动,并具沿等深线运动的趋势;而在夏季西南风盛行时节,同时在南高北低的水位条件下[17],沿福建近海以北向流为主,并有跨等深线向外海运动的趋势,这主要由与上升流有关的离岸运动造成.

图3 福建近海部分浮子漂移轨迹Fig.3 Part of float drift trajectory in Fujian offshore

2.2 浮子漂移实验结果分析

综合统计分析福建近海预测时段内浮子的漂移,主要有3种形态:第1种,受涨落潮作用,浮子在一段时间内在某一海域往复运动,这主要发生在福建沿岸的强潮海湾内;第2种,受长时间的季风和海流的作用,浮子在一段时间内沿一个方向持续运动,这主要在沿岸近海区域,也会受离岸流的作用跨岸迁移;第3种是前两种的结合,主要在河口区,如九龙江河口区等.

针对不同漂移形态,本研究从北至南选取3种特征浮子进行分析(图4):1号浮子在兴化湾内往复运动,运动时间段为2021年3月30日18:00至2021年4月1日6:00;2号浮子沿福建近岸海域往南运动,运动时间段为2019年10月24日8:00至2019年10月30日0:00;3号浮子在厦门湾内往复运动,运动时间段为2021年3月28日8:00至2021年3月29日20:00.

图4 特征浮子漂移轨迹Fig.4 Typical float drift trajectory

2.2.1 海湾内浮子运动轨迹分析

1) 兴化湾1号浮子

兴化湾海域整体呈长方形,无大河流直接流入,仅在湾顶附近有木兰溪和荻芦溪两条小溪,与湾外台湾海峡相连的有南日水道和兴化水道.兴化湾的潮流为正规半日潮,湾内涨落潮流由南日水道和兴化水道两股水流形成,流路稳定,呈现往复流,水流流向与深槽轴线基本一致[18].

从图5可以看出,1号浮子运动海域为兴化湾海域南日水道与兴化水道的交汇处,浮子主要在湾内兴化水道呈往复运动,流向与深槽轴线基本一致.浮子在海域最大漂移速度为0.87 m/s,浮子在转向时流速最小.流速大小与流速方向均呈周期性变化.

u、v分别表示流速的东分量和北分量.图5 1号浮子漂移轨迹Fig.5 Drift trajectory of float No.1

FES是基于有限元流体动力模型的同化潮流模型,FES2014是FES系列潮流模型的最新版本,空间分别率为1/16°×1/16°[19].而福建沿海海域为正规半日潮,S2分潮振幅约为M2分潮的一半,N2和K2分潮则要小得多[20],M2分潮占主导.本研究从FES模型获取15个分潮的调和常数,图5插图为根据M2分潮调和常数计算得到的M2分潮潮流椭圆图,其最大流速为0.45 m/s,与郭玉臣等[21]利用实测数据调和分析得出的结果(0.47 m/s)基本一致,可以看出浮子的运动方向与M2分潮最大流速方向基本一致.

余流是海流分离出潮流后的水体运动,具有单向流动特性,因此与近海的物质输运密切相关.在余流非线性作用比较弱的海域,拉格朗日余流可以代表长期平均水量传输速度.拉格朗日余流可以使用Awaji等[22]提出的近似算法进行计算.假设Xt是浮子在t时刻的位置,Xt+T是在一个潮周期T(本研究取12 h)之后的位置,则拉格朗日余流流速VL可以用下式进行计算:

(4)

利用上述公式,将实测海流分离成余流与潮流,同时利用从FES模型获取的15个分潮的调和常数计算潮流,将模型结果与实测潮流结果对比.计算结果如图6所示,兴化湾内的海流主要以潮流为主,拉格朗日余流较小,最大值为0.13 m/s,略大于郭玉臣等[21]利用实测数据调和分析得出的兴化湾最大余流(0.10 m/s).基于FES模型调和常数计算得到的潮流与浮子轨迹计算得到的潮流存在一定差异,前者流速比后者小0.20 m/s左右,且转向时间约延迟1 h.这是由于兴化湾海域虽然是M2分潮占主导,但是浅海分潮影响明显[20],同时FES2014作为全球潮汐数据可能在沿岸港湾的准确性不足,所以振幅与相位存在一定偏差.

图6 调和常数计算潮流流速与1号浮子轨迹计算流速Fig.6 Tidal velocity calculated from harmonic constant and velocity calculated from the trajectory of float No.1

综上,兴化湾内径流影响较小,浮子在兴化湾内的运动主要受潮流作用,浮子的运动轨迹也与该区域的潮流运动一致.

2) 厦门湾3号浮子

厦门湾位于台湾海峡西岸南口,是一个半封闭型海湾,潮流属半日潮流,往复流.厦门湾海域主要受河流径流与海洋潮流的相互作用.九龙江是流入厦门湾海域的最大河流,年平均径流量为150亿m3[23].

如图7所示,3号浮子运动海域为九龙江口海域以及青屿水道,浮子有时也会发生偏移,流向厦门东侧水道.浮子在海域最大漂移速度为1.39 m/s,出现在九龙江河口处,从整体运动轨迹上看,浮子呈往复运动,涨落行程约15 km.

图7 3号浮子漂移轨迹及流速Fig.7 Drift trajectory and velocity of float No.3

将实测海流分离成余流与潮流,同时利用从FES模型获取的15个分潮的调和常数计算潮流,结果如图8所示.FES模型M2分潮的最大流速为0.11 m/s,远小于靖春生等[24]和綦梦楠[25]在厦门湾的研究结果(最大流速约0.48 m/s).同时从图8可以看出,基于FES模型调和常数计算得到的潮流与浮子轨迹计算得到的潮流存在很大差异,基于FES模型调和常数计算得到的潮流流速均小于0.2 m/s,远小于浮子轨迹计算得到的潮流流速,因此对于九龙江口-厦门湾海域,FES模型数据并不适用.图8显示研究海域的拉格朗日余流最大值为0.25 m/s,出现在九龙江口海域,平均为0.15 m/s,方向为流出河口区,这主要体现的是九龙江径流的作用.此外,采用3号浮子轨迹计算的潮流并不像1号浮子一样呈现周期性变化,这是由于3号浮子运动所属海域在河口区及厦门港交汇区,流场较为复杂,九龙江河口区的转流时刻落后于高、低潮,而厦门港海域转流时刻却比高、低潮超前20～40 min[26].

图8 调和常数计算潮流流速与3号浮子轨迹计算流速Fig.8 Tidal velocity calculated from harmonic constant and velocity calculated from the trajectory of float No.3

2.2.2 海湾外浮子运动轨迹分析

2号浮子运动轨迹分布在福建近海,从图9可以看出,浮子漂移轨迹受沿岸流与潮流的共同作用:沿岸流导致的沿岸向南运动以及潮流导致的周期性往复运动.浮子起点为兴化湾湾口附近,沿着海岸向南漂移,每移动一段时间后,浮子漂移方向会发生改变,向岸移动一段距离;在到达厦金海域附近时,浮子会改变路径,沿岸向北运动一段距离,回旋形成一个椭圆,然后接着继续往南运动.根据图9中的近岸海域基于FES模型数据计算的M2分潮潮流椭圆分布可以看出,在兴化湾南部海域,潮流椭圆的椭率接近1,为旋转流特征,因此浮子漂移运动表现为垂直向岸、离岸运动与沿岸向南运动的组合.到了泉州和厦门附近海域,潮流椭圆长轴方向为沿岸方向,涨落潮流较大,因此浮子轨迹受潮流影响,可以回旋形成一个椭圆,浮子的漂移轨迹与潮流椭圆有一致性.

图9 2号浮子漂移轨迹与M2潮流椭圆Fig.9 Drift trajectory of float No.2 and M2 tidal current ellipse

将实测海流分离成余流与潮流,同时利用从FES模型获取的15个分潮的调和常数计算潮流,将模型结果与实测潮流结果对比.从图10可以看出,福建近海的拉格朗日余流较大,最大值为0.56 m/s,平均值为0.31 m/s,方向主要沿岸向南;10月28日余流流速减小为0.2 m/s左右,此时浮子位于厦金海域附近,余流的减小使得潮流作用效果明显,因此该海域的浮子回旋距离较远,这与前文浮子轨迹基本一致.根据FES模型调和常数计算的潮流与浮子轨迹计算的潮流基本一致,说明在外海开阔海域,浮子的实测数据分离出的潮流数据能较好地体现潮流运动,同时FES模型能较好地模拟近海潮波运动.从图11的余流场图与风场图可以看出:浮子运动这段时间,10月东北风盛行,南下的沿岸流增强.将余流流速与风场速度投影到30 m等深线方向,并计算得两组数据的相关系数为0.37,判断浮子的余流主要受沿岸流与风的作用.

图10 调和常数计算潮流流速与2号浮子轨迹计算流速Fig.10 Tidal velocity calculated from harmonic constant and velocity calculated from the trajectory of float No.2

图11 2号浮子余流场与风场图Fig.11 Residual flow of float No.2 and corresponding wind

综上,秋季浮子在福建近海的运动主要受潮流、风与沿岸流的共同作用,浮子的运动轨迹也显现出潮流作用的周期性以及风与沿岸流作用的均一性.当沿岸流较弱且潮流呈现较强往复流时,浮子可呈现回旋轨迹;反之浮子具有沿等深线南下趋势.

2.3 预测模式计算结果与分析

本研究用漂移轨迹距离平均误差随时间和漂移距离的变化两种形式对漂移轨迹预测的准确性进行评判(图12).从漂移浮子中选取273个有效漂移过程进行评估,累计预测天数为570 d,进行24和48 h的漂移预测次数分别为273和181.根据误差统计,得到24和48 h距离的平均相对误差分别为17.3%和12.5%,而在120 km的模拟范围内距离平均相对误差大多小于20%,说明本预测模式在台湾海峡海域及港湾区总体具有较高的准确性.

图12 距离误差及其平均误差随漂移时间(a)和漂移距离(b)的变化Fig.12 Change of distance errors and their average errors with drift time (a) and drift distance (b)

图13为湾内(受潮流作用为主)、湾外(受潮流、环流共同作用)的距离平均误差,可以看出湾外浮子的距离平均误差小于湾内浮子.因此,相对于港湾区域,本轨迹预测模式对于福建近岸海域环流的模拟有更高精确度.选取前文的3个特征浮子进行分析,从图14可以看出,模式预测轨迹总体与实测浮子轨迹趋势基本一致.对于1号浮子与3号浮子,预测浮子在湾内呈往复运动,且最大往复距离基本一致;对于2号浮子,受潮流与沿岸流的共同作用,预测浮子沿福建沿岸南下运动并在离岸线方向往复漂移,虽然在垂直岸线方向的距离存在偏差,但是总的迁移方向和距离基本一致.因此模拟结果能够再现浮子实际漂移特点.

图13 距离平均误差随漂移距离的变化Fig.13 Change of distance average error with drift distance

图14 特征浮子实测与预测轨迹对比Fig.14 Comparison of measured and predicted trajectories of characteristic floats

图15为特征浮子预测轨迹与实测轨迹的距离平均误差,从模拟的距离平均误差来看,湾内主要受潮流作用,厦门湾与兴化湾湾内潮流的往复流特征明显,因此模拟轨迹漂移距离的平均误差并非单调增加.1号浮子最大漂移距离平均误差约11.6 km,这是由于浮子在转向时间上存在差异,导致出现较大误差,所以对于湾内主要受潮流作用的海域,潮流模拟的准确性会极大影响预测模式的结果.2号浮子主要受潮流与沿岸流作用,与湾内浮子相比,距离平均误差反复性变化较弱,其最大漂移距离平均误差约7.2 km.

图15 特征浮子预测与实测轨迹的距离平均误差Fig.15 Distance mean errors between predicted and measured trajectories of characteristic floats

3 结 论

本研究采用风场数值模型和业务化海流数值模型,结合拉格朗日粒子追踪法,建立了福建近海海上漂移物漂移轨迹预测模式,模拟并分析实测浮子轨迹数据,得到以下结论:

兴化湾内浮子的运动主要受潮流作用,其运动轨迹也与该区域的潮流运动一致,呈现往复运动的趋势;厦门湾内的浮子运动主要受潮流与径流的共同作用;而湾外近海的浮子运动主要受潮流与沿岸流的作用,在外流环流的驱动下有沿等深线运动的趋势,在潮流作用下可能局部形成回旋轨迹.分析结果可以看出浮子漂移实验能够捕捉到研究海域表面海流的潮流、沿岸流及径流共同作用的特征.结果也表明FES2014全球潮汐数据在湾外近海区有较高的准确性,而在港湾区准确性差.福建近海海上浮子漂移实验表明,本模式模拟结果与浮子漂移结果基本一致;相对于湾内的轨迹预测,湾外的模拟有更高的可信度.这与浮子轨迹得到的潮流和余流(环流)分析结果是一致的.未来可通过提高潮流模拟的准确度来进一步提高湾内漂移物漂移轨迹模拟的准确度.