2008年黄海绿潮路径的数值模拟

郑向阳,邢前国,李 丽,施 平

(1.中国科学院 烟台海岸带研究所,山东 烟台 264003; 2.中国科学院 研究生院,北京 100049; 3.中国科学院 海岸带环境过程重点实验室,山东 烟台 264003)

从2007年开始,山东半岛沿岸的威海、青岛、日照连续 4年爆发大规模的绿潮灾害。造成绿潮的是一种叫做浒苔的大型石莼科(Ulvaceae)藻类[1],浒苔本身没有毒性,但大规模的绿潮发生时会覆盖近岸水体和沙滩,对旅游业和航运造成影响。浒苔大量死亡后,会消耗海水中的溶解氧,导致其他海洋生物的死亡,对近海养殖业造成致命打击[2]。2008年登陆青岛的大规模绿潮对当时即将开幕的奥运会帆船比赛的顺利进行造成了重大影响。2009年山东省因浒苔灾害,致直接经济损失 6.41亿元,占该省海洋灾害总损失的68%[3]。因此追溯浒苔来源、弄清浒苔运移路径对于从根本上防止、减少浒苔灾害具有重要意义。

最初有学者认为浒苔来源于青岛近海的富营养化水体,但后来的研究者通过DNA测序、生物形态学、遥感、物理海洋、实地观测等手段发现浒苔并不是在青岛近海产生的,而是最早出现在南黄海的江苏近岸,在表层风海流的作用下,漂至山东半岛沿岸,目前这一观点基本得到认可[4-7]。遥感在监测浒苔分布和分析浒苔来源方面具有重要的作用,在宽阔的海面上,只有利用卫星遥感图片才能获得可靠而及时的浒苔信息,但遥感观测也有局限性,当被监测海面上空云量较多时,可能会出现连续几天不能获得浒苔信息的问题,而且遥感不能对浒苔的路径做出定量的预测; 数值模拟可以对海洋现象进行连续模拟,在有良好资料的基础上,利用数值模拟来预测浒苔路径可以补充遥感手段的不足。本文利用 2008年黄海 MODIS遥感数据作为验证,模拟了2008年黄海绿潮的输运路径。

1 模式介绍

本次模拟使用的海洋模式是 FVCOM (Finite Volume Coast and Ocean Model)[8]。FVCOM 是由美国麻州大学(The University of Massachusetts)陈长胜博士研究组建立的一个三角形网格、有限体积三维原始方程组海洋模式。它包括质量守恒方程、动量守恒方程、温度方程、盐度方程和一个湍流闭合方程。目前,该模型已经成功地应用于许多河口、近岸海区和某些内陆湖泊,如美国的乔治浅滩、五大湖,中国的渤海、长江口、南海等海区。它代表数值模式新的发展方向,具有广阔的应用前景[9]。

相对于深海大洋,河口、近岸海区与人类的生产、生活的联系更为密切。近岸海区有曲折的海岸线,有的海区有非常宽阔的潮间带。这些几何形状不规则的河口海岸系统给近海水动力模拟精度的提高带来了很大困难。FVCOM使用三角形网格,有限体积方法对海洋原始方程进行离散,综合了传统的有限差分法和有限元法的优点,不仅原理简单、计算速度快,而且能够很好地拟合岸线,特别适合近海水动力模拟。FVCOM采用模式分离技术,将与表面重力波相关的快速运动和与密度相关的较慢的运动分开计算,根据 CFL条件分别选择时间步长,可以大大节省计算时间; 使用干湿点处理技术,能够较好地模拟具有宽广潮间带的海域[10]。

在水动力模拟的基础上,FVCOM可以进行生态模拟、泥沙输运模拟、水质模拟和拉格朗日粒子模拟(3-D Lagrangian Particle Tracking)。本文进行的浒苔轨迹模拟就是利用拉格朗日粒子跟踪模块,此模块的方程如下：

这里x表示粒子在时刻t的位移,dx/dt表示在此时刻粒子位移的时间变化率,v( x,t)是由水动力模块生成的三维速度场。FVCOM中利用显式龙格-库塔多步法来计算此常微分方程。

拉格朗日粒子跟踪模拟是在某一时刻把某一粒子投入到模拟区域中,跟踪过程中不考虑粒子自身的变化,只考虑粒子在水流的作用下运动。浒苔虽然在漂移过程中会不断地生长,但它们基本是作为一个整体来运动的,不会像营养盐或浮游植物扩散到海水中,所以用拉格朗日粒子跟踪模型来模拟浒苔路径是可行的[11]。

2 模型的配置

研究区域主要是黄海中部。为了减小边界对研究区域的影响,将开边界取在长江口北岸启东嘴至朝鲜半岛西南角一线,计算海区包括整个黄海和渤海,117.5°～126.8°E,32°～41°N,水深(单位：m)数据采用分辨率为 1′的 ETOPO1地形数据,如图 1所示。模式水平方向采用三角形网格,在岸线曲折的海区适当加密,包含8 603个网格节点和16 543个三角形单元,最小网格步长为6 300 m。垂向有6个sigma层。内模式时间步长为 60 s。模拟时间从 2008年5月15日到7月1日。风场采用quickscat日平均的资料,底摩擦系数取0.002 5。

浒苔主要是因表层风海流的作用而漂移,所以本文只考虑正压模式,不考虑热盐环流,温盐初始条件在整个海区是均匀的,也没有考虑海洋和大气的热交换。开边界上选取此海区四个主要分潮M2,S2,K1,O1作为潮汐强迫。模式采用零初始条件,为了保证计算稳定性,开边界潮位强迫从0逐渐增加,经过两个潮周期达到正常变化[12],前三天只加潮位不加风,三天后潮流稳定以后再加上风,产生风海流。

首先来验证此模型在计算海区的可行性。只在开边界处加M2,S2,K1,O1四个主要分潮,将模式运行33 d,输出后30 d的每隔一小时的潮位数据,用最小二乘法做调和分析,求出各个网格点上主要的半日分潮 M2分潮和主要的全日分潮 K1分潮的振幅和迟角[13],画出这两个分潮的等振幅图和等迟角图(见图2)。M2分潮在研究海区共4个无潮点,模式仅在秦皇岛外海的M2无潮点没有模拟出来,可能的原因是渤海的底摩擦系数比黄海小很多,但为了计算方便,模式取了一个统一值 0.002 5,而在黄海区域M2分潮无潮点的位置与振幅的分布与图集基本一致[14]; K1分潮有2个无潮点,无潮点的位置和等振幅线的分布与图集比较符合,所以模拟结果基本能够刻画研究区域的潮波系统。

图2 计算得出的M2,K1分潮的同潮图Fig.2 Computed co-tide lines of M2 and K1 constituents

浒苔的漂移主要受表层流场的影响,表层潮流场是表层流场的重要组成部分,本文对表层潮流场进行了调和分析,计算了 M2,K1分潮流表层潮流椭圆要素,图 3给出了这两个分潮流的潮流椭圆图。M2分潮流在辽东湾、渤海湾、渤海海峡以及山东半岛沿岸的椭圆长轴几乎与岸线平行,并且远远大于短轴,为往复流,在离岸较远处则表现出旋转流的特征; 烟台、秦皇岛近岸以及青岛以东的黄海中部,M2分潮最大潮流很小,是M2分潮的圆流点; 强流区主要分布在西朝鲜湾、江华湾、朝鲜半岛西南、苏北近海,这些区域的M2分潮最大流速超过1 m/s,渤海的渤海湾、辽东湾、渤海海峡最大流速也较大,可达0.6～0.8 m/s。K1分潮流跟M2分潮流一样,在近岸表现为往复流,离岸较远处为旋转流; 在渤海海峡可达 0.5～0.6 m/s,在渤海海峡两侧的黄海和渤海也可以达到0.3 m/s,整个黄海K1分潮流都比较小,除了在靠近黄海南部边界处超过0.2 m/s外,其他海区都在0.1 m/s以下。

图3 计算得出的M2,K1分潮流表层潮流椭圆图Fig.3 Computed tidal current ellipses of M2 and K1 on sea surface

综上分析,主要分潮和分潮流的特征都符合实际情况,可以用此模型来研究黄海浒苔漂移路径。

3 粒子跟踪模块配置

通过遥感图像能够清晰监测到浒苔的日期是2008年5月15日,之前或是因为遥感的空间分辨率不够或是因为云层的覆盖,不能监测到浒苔,所以本次模拟将粒子释放的时间定在5月15日,将粒子的初始位置选择在当天所监测到浒苔的位置(121.45°E,34.70°N),粒子位于海表层,如图 4。粒子初始位置对于粒子的漂移轨迹影响很大,考虑到浒苔呈片状分布,选取了1组共15个粒子,利用粒子的整体分布特征来模拟浒苔的分布。

图4 粒子初始位置和初始时间示意图Fig.4 Initial positions and times of particles

4 结果分析

本次模拟没有考虑浒苔的生长率,只是模拟浒苔随表面水流的运动轨迹,而粒子的个数并不随时间改变,粒子只是表示一个位置,粒子的大小也不会随时间而改变。而实际上浒苔在运动过程中不断生长和死亡,浒苔数量和覆盖面积一直在变化。粒子跟踪模拟虽然不能表示完整的绿潮爆发过程,但是从粒子所覆盖的范围与运动趋势上来看,粒子的运动轨迹和浒苔的输运路径还是比较一致的。由于黄海上空云的影响,能够监测到浒苔的天数并不多,我们从中选取了 4 天的浒苔遥感图像与模拟结果进行对比,如图5。5月20日浒苔的中心位置与初始位置几乎没有变化,只是覆盖面积扩大了数倍,这跟数值模拟的结果几乎是一样的; 5月30日,浒苔中心位置向北移动到青岛东南120 km处,在同一天15个粒子中的大部分散布在中心位置周围。之后一直到 6月 18日,由于云层的阻挡,这期间遥感图像没有能够监测到浒苔; 6月 18日,大规模的浒苔向西漂移,边缘的浒苔已经接近胶南和青岛,并且向日照方向漂移,中心位置在青岛以南65 km,数值模拟也基本反映了这一现象。6月27日,数值模拟的结果显示,到达青岛的粒子经青岛近海往东北方向漂移,同一时期的遥感影像也显示,青岛东北方向的乳山南部海域已经有分布较为分散的浒苔。然而在青岛近岸海域,遥感影像显示出大量的浒苔在此富集,数值模拟的结果却没有再现浒苔在青岛近岸的滞留过程。可能的原因如下：浒苔靠近近岸水域之后可以广泛附着于沿岸中低潮区的砂砾、岩石、滩涂和石沼中。另一方面,稳定下来的浒苔经过大量的生长,聚集成团后,其质量增大,受表层海流的影响变小,这在浅滩区域是必须考虑的因素。数值模拟仅仅考虑水动力对浒苔漂移的影响,并没有考虑浒苔在近岸的附着作用。由于遥感图像分辨率的限制,6月27日的遥感图像没有显示浒苔到达日照,但根据当地新闻媒体的报道[15],日照辖区北部海域出现覆盖面积接近100 km2浒苔,这与模拟结果是一致的。

数值模拟结果和遥感图像良好的一致性反映了水动力场对浒苔输运的重要作用。利用模拟结果,绘出粒子的运动轨迹,如图 6,粒子的运动趋势是首先向青岛方向漂移,然后沿海岸线向日照和乳山方向运动,这与遥感分析结果也是一致的[16],所以在不考虑浒苔自身生长率的情况下可以利用 FVCOM 的粒子追踪模型进行初步的绿潮路径模拟。这对于深海处面积较小的浒苔块是适用的,浒苔到达近岸后,容易聚集成团,又由于近岸水深较浅,浒苔的水下部分使其漂移方向与表层海流的相关性降低。今后的工作将重点放在对浒苔到达近岸后的生消模拟,这需要对浒苔的生理特性有充分的了解。

5 总结

本文的粒子追踪模拟是对绿潮生态动力学的初步研究,但2008年浒苔大规模的爆发是受多方面因素影响的,总体可以分为直接因素和间接因素[17]。直接因素包括光照强度、温度、海水透明度、营养盐的理化性质,它们直接影响浒苔的生长率; 间接因素包括浒苔源头的位置和表层流场,它们主要影响浒苔的漂移路径。目前我国对绿潮生态学的研究还远远不够,需要综合现场观测、实验数据来增强光照、温度、营养盐、海水pH值、透明度、摄食压力对浒苔生长率影响的认识,建立适合某一海区的经验公式[18],为绿潮生态模型的建立提供理论依据。本文的模拟是研究绿潮生态动力学的初步阶段,仅考虑流场对浒苔的影响,对于浒苔爆发的模拟还不完整。风是决定浒苔路径的关键因素,目前还无法做到在春季对夏季风场做出预报,风场资料的精度和分辨率也远远不够,这也影响了浒苔路径的预报。浒苔初始位置的选取也会影响浒苔的未来路径,精确地确定浒苔的初始位置和初始时间非常重要; 苏北浅滩的浑浊水体以及遥感图像分辨率的不足增大了对浒苔溯源的难度,这既需要提高遥感监测的技术,又需要加强和数值模拟的有机结合,来弥补遥感的不足。

图5 浒苔位置模拟值和遥感反演对比Fig.5 Comparison of the results of numerical simulation and remote sensing

图6 粒子轨迹图Fig.6 Tracks of particles

[1]梁宗英,林祥志,马牧,等.浒苔漂流聚集现象的初步分析[J].中国海洋大学学报,2008,38(4)：601-604.

[2]王文娟,赵宏,米锴,等.大型绿藻浒苔属植物研究进展[J].齐鲁渔业,2009,26(10)：3-6.

[3]国家海洋局.2009年中国海洋灾害公报[EB/OL].[2010-03-05]http：//www.soa.gov.cn/soa/hygb/zhgb/webinfo/2010/03/1271382648968147.htm.

[4]Hu Chuanmin,Li Daqiu,Chen Changsheng,et al.On the recurrent Ulva prolifera blooms in the Yellow Sea and East China Sea[J].Journal of Geophysical Research,2010,115(C5)：C05017.

[5]邢前国,郑向阳,施平,等.基于多源、多时相遥感影像的黄、东海绿潮影响区监测[J].光谱学与光谱分析,2011,31(6)：1644-1647.

[6]Liu Dongyan,John K K,Xing Qianguo,et al.World’s largest macroalgal bloom caused by expansion of seaweed aquaculture in China[J].Marine Pollution Bulletin,2009,58(6)：888-895.

[7]夏斌,马绍赛,崔毅,等.黄海绿潮(浒苔)暴发区温盐、溶解氧和营养盐的分布特征及其与绿潮发生的关系[J].渔业科学进展,2009,30(5)：94-101.

[8]Chen Changsheng,Liu Hedong,Robert C B.An unstructured grid,finite-volume,three-dimensional,primitive equations ocean model：Application to coastal ocean and estuaries[J].Journal of Atmospheric and Oceanic Technology,2003,20：159-186.

[9]朱建荣.海洋数值计算方法和数值模拟[M].北京：海洋出版社,2003：128.

[10]陈长胜.海洋生态系统动力学与模型[M].北京：高等教育出版社,2003：312.

[11]李宁.近岸水质的遥感监测和数值模拟研究[D].天津：天津大学建筑工程学院,2009.

[12]宋文鹏.渤海冬、夏季温盐场结构及其海流特征分析[D].青岛：中国海洋大学海洋环境学院,2009.

[13]黄祖珂,黄磊.潮汐原理与计算[M].青岛：中国海洋出版社,2005：118-119.

[14]海洋图集编委会.渤海、黄海、东海海洋图集——水文分册[M].北京：海洋出版社,1992：429.

[15]王霞.日照海事局查看浒苔情况[N].黄海晨报,2008-06-27(A3).

[16]李三妹,李亚君,董海鹰.浅析卫星遥感在黄海浒苔监测中的应用[J].应用气象学报,2010,21(1)：76-82.

[17]徐兆礼,叶属峰,徐韧.2008年中国浒苔灾害成因条件和过程推测[J].水产学报,2009,33(3)：430-437.

[18]刘桂梅,李海,王辉,等.我国海洋绿潮生态动力学研究进展[J].地球科学进展,2010,25(2)：147-153.