高 松 ,黄 娟 ,白 涛 ,曹雅静 ,徐江玲
(1.国家海洋局北海预报中心,山东 青岛 266033;2.山东省海洋生态环境与防灾减灾重点实验室,山东 青岛 266033)
随着渤海海上油田迅速发展以及一些输油管线进入老化期,海上溢油渗油事件和海上无主漂油时有发生,如2002年河北黄骅附近海域发生大面积原油污染事件,由于无法确认事故责任者,使受损失渔民无法得到有效赔偿。2006年,山东省长岛、龙口、蓬莱岸段及其附近海域发生大面积原油污染,由于无法直接认定溢油源,在确认事故责任者时颇费周折。2010年秦皇岛市某浴场沙滩发现黑色油块,严重影响浴场的正常开放。因此,发展有效、准确的溢油源鉴别预测技术,对责任事故的认定和溢油污染的防治具有非常重要的意义。
海洋溢油漂移预测数值模型的开发研究已有几十年的历史。最早的溢油漂移数值模型是连续介质模型[1],目前油粒子模型占据了主导地位[2]。起初的粒子模型主要是模拟溢油水平漂移路径,后来进一步发展成能够模拟溢油三维扩散运动的模型[3]。这些预测模型尽管还有许多值得进一步完善的方面,但已经取得了较为广泛的应用[4-7]。溢油溯源模型的研究发展相对较晚,但其应用前景已经取得了关注[8]。海洋动力环境的复杂性和数值模型自身的缺陷常常造成较大的误差,溢油溯源模型开发的技术关键在于逆向追溯的精确性。
集合预报方法在气象预报方面应用较早,且已比较成熟。其基本方法是通过一定的数学方法结合天气学原理,构造一组有限数量的预报样本,把不同的预报样本进行“集合”,形成集合预报。简单地说,就是从一组初值出发,得到一组数值预报。因为得到这组数值预报所用的初值和模型过程彼此之间有一定的差别,所以就反映大气的不确定性。从而把过去传统意义上单一的确定性天气预报,变成不确定性预报。集合预报的目的就是尽可能使得到的这组不确定性预报包含未来大气可能出现的所有状态,从而达到提高预报水平的目的。集合预报是数值预报的一次革命,即人们从“确定论”向“随机论”转变的思维模式的革命。
油粒子在风力、波浪、潮流和环流的联合作用下的运动,可采用拉格朗日粒子追踪方程模拟。通常,水动力数学模型基于欧拉场建立,而要描述质点的运动,需采用拉格朗日的观点,这就涉及如何将欧拉场中的结果转换为拉格朗日质点位移。
在欧拉场中,对平面二维问题,任意空间点的速度V可由公式(1)表示,其中,x和y表示空间点的位置,t表示时间。
采用拉格朗日观点,任意质点的速度 VL可由公式(2)表示,其中,xL和 yL表示任意质点的空间位置,X表示质点位移。
公式(2)实际上建立了求解质点位移的一阶常微分方程。改写公式(2),质点的运动轨迹可通过公式(3)中的积分求得。
若每一时刻的 VL已知,可通过数值积分的方法由上式求出质点的运动轨迹。
粒子的漂移速度 VL计算公式为: 由公式(4)计算,其中,Vw表示风力和波浪作用产生的速度分量,Vt表示潮流作用产生的速度分量,Vr表示潮致余流作用产生的速度分量,Vh表示环流(包括: 风海流和密度流)作用产生的速度分量,VR表示由波浪等海洋环境引起的粒子不确定运动。
溢油溯源是在已知溢油被发现的时间和地点以及海洋环境背景场的前提下,逆向(时间反方向)积分溢油漂移的轨迹,得到在过去任意时刻溢油可能所处的位置。溢油溯源数值模型是在溢油漂移预测模型的基础上,采用逆向积分方法发展而来。溢油粒子逆向溯源过程中的随机运行,通过“蒙特卡罗”方法[9]进行描述。
目前,受观测手段和预报技术的限制,所获得的海洋背景环境要素(主要是海面风场和流场)都存在一定的误差。在溢油溯源模拟过程中,通过对海洋背景环境要素数据在一定误差范围内的修正设计,进行多种情况溯源模拟,将获得的一组溯源模拟结果,进行“集合”,统计可能溢油源范围和概率。如:获得风速为10 m/s,为误差20%,可设计修正数据5组,分别为12、11、10、9、8 m/s。由于无主溢油的溢油发生时间无法准确的确定,则将一定的时间段(如24 h)设为油粒子的统计时间域,累计这一时间域内每个统计单元内的粒子数,经过归一化处理后作为该单元可能是溢油源的概率。
最终每个统计网格溢油源的可能概率G由公式(5)表示。其中,S为某个统计单元的累计粒子数,M为某种溯源试验抛放的粒子总数,N为参与集合统计的试验个数。
2010年 6月底,在秦皇岛某浴场沙滩发现黑色油块,根据之后的鉴定与调查情况,确认此次发现油污为渤海中西部某平台 5月中旬泄露原油所致。利用本文研发的溢油溯源模型,在秦皇岛该浴场附近海域进行粒子抛放,模拟粒子在海面风和海流的共同作用下,在海面的逆向漂移轨迹。溢油溯源模式计算范围为渤海区域(37°~41°N,117.5°~122°E),水平分辨率为1/90°(见图1,其中1个网格代表模式中10×10个计算网格)。溢油溯源模式采用的渤海气象和流场数据,为国家海洋局北海预报中心业务化数值预报结果。气象模式采用 Weather Research and Forecasting Model(WRF)研发,采用两重嵌套方式,大区计算范围为渤黄东海,空间分辨率为27 km;小区计算范围为渤黄海,水平分辨率为 9 km。海流模式采用Regional Ocean Model System (ROMS)研发,采用两重嵌套方式,大区为渤黄东海,空间分辨率为 1/30°;小区为渤黄海,水平分辨率为 1/90°,垂向分为6层,进行潮流和环流的耦合模拟。
图1 溢油模式计算网格及地形(m)图Fig.1 The computing grid and topographic map of oil spill model (m)
图2a是采用 100%风速情况,粒子溯源轨迹及分布图。由图2a可以看出: 油粒子6月30日~5月13日由发现油污的浴场(YC)开始,向偏南方向移动,在5月13日前后,大部分粒子位于某石油平台(PT)的西侧海域。各特征时刻粒子形成不同的分布形态,总体分布范围不断扩大,这也符合粒子随机走动的一般规律。
图2b是采用 100%风速情况,粒子溯源模拟概率分布图。概率分布的统计粒子样本选取5月13日00~23时,24 h内的累加粒子。粒子统计范围为(119.205°~119.745°E,37.925°~ 38.465°N),网格设计间距为 0.03°×0.03°,统计网格溢油源的可能概率为每一个网格中的粒子数占统计样本的百分比。从图2b可以看出,粒子的分布概率呈现由中心向四周逐渐减少的趋势,PT位于高概率区域的东侧,概率约为1.5%。
图2 油粒子溯源轨迹及概率分布图Fig.2 Thetracing tracks and possible distribution of oil particle
图3 不同集合预测的油粒子概率分布比较图Fig.3 Comparison of possible oil particle distribution with different ensemble sampling
图 3为不同集合预测情况下的油粒子概率分布图,图 3a参与集合预测的统计样本是 105%、100%和 95%三种风速情况,图 3b参与集合预测的统计样本是110%、105%和100%三种风速情况。可以看出,这两种集合预测概率统计结果,均优于100%风速概率统计结果(图 2b)。与 100%风速相比,两种集合预测的PT所在网格的概率值均有增大,分别为 2.5%和3.2%,增加比值为 66%和 110%,且预测平台基本位于概率统计值较高的中心区域,可能溢油源的分布区域略有增加。此外,图 3b的集合预测概率结果明显优于的结果图3a,图3b中PT所在网格的概率值比图 3a增大28%,且PT位于可能溢油源的分布区域的中心位置。因此,无主漂油源集合预测可以通过对环境背景场和溢油溯源模型的系统误差科学分析,合理设计参与集合概率统计的模拟结果,进一步提高溯源预测的准确度,缩小高概率区的分布范围。
本文采用基于“拉格朗日”粒子追踪方法的油粒子数值模拟溯源模型,引入天气预报中的“集合预报”方法,考虑气象、海洋环境背景场误差因素,建立了无主溢油源的集合预测方法。通过对2010年溢油溯源事件的试验模拟,说明集合预测方法的科学性和正确性,为进一步解决溢油源的判定问题提供了新的思路和方法。
[1]Fay J A.Physical processes in the spread of oil on a water surface[J].International Oil Spill Conference Proceedings,1971,1971(1): 463-468.
[2]张存智,窦振兴,韩康,等.三维溢油动态预报模式[J].海洋环境科学,1997,16(1): 22-29.
[3]娄安刚,吴德星,王学昌,等.三维海洋溢油预测模型的建立[J].青岛海洋大学学报(自然科学版),2001,31(4): 473-479.
[4]DHI.MIKE 21/3 oil spill[R].Shanghai: DHI,2012.
[5]Delft.Simulationof mid-field water quality and oil spills using particle tracking[R].Netherlands: Delft,2005.
[6]汪守东,沈永明,郑永红.海上溢油迁移转化的双层数学模型[J].力学学报,2006,38(4): 452-460.
[7]李彤,谢志宜.水上事故溢油漂移轨迹预测模型研究与应用[J].环境科学与管理,2013,38(7): 56-61.
[8]赵如箱.浅谈溢油模型的发展及其应用设想[J].交通环保,2000,21(4): 15-17.
[9]娄安刚,王学昌,于宜法,等.蒙特卡罗方法在海洋溢油扩展预测中的应用研究[J].海洋科学,2000,24(5): 7-9.