陳勤思，胡松，魏永亮

（1.上海海洋大學(xué) 海洋生態(tài)系統(tǒng)與環(huán)境實(shí)驗(yàn)室，上海 201306；2.上海海洋大學(xué) 海洋科學(xué)學(xué)院，上海 201306）

渤海是我國(guó)重要的石油產(chǎn)區(qū)，也是重要的海上運(yùn)輸通道。隨著海洋油氣業(yè)和海上運(yùn)輸業(yè)的發(fā)展，渤海逐漸成為溢油污染事故的高發(fā)區(qū)。渤海海域分布著眾多海上石油平臺(tái)，石油開采過程中的泄漏和船舶事故引起的溢油污染事件時(shí)有發(fā)生，給區(qū)域海洋生態(tài)環(huán)境造成損害。柯麗娜等總結(jié)了渤海海域中來自船舶和石油平臺(tái)的溢油風(fēng)險(xiǎn)源分布情況（圖1）[1]。此外，一些可能來源于石油管道滲油、沉船和非法傾倒等多種原因引發(fā)的溢油，常常因其來源不明而形成無主漂油，這類溢油的污染危害同樣不可小覷，需要對(duì)此進(jìn)行風(fēng)險(xiǎn)分析從而及時(shí)采取防污染的措施。

圖1 渤海海域溢油風(fēng)險(xiǎn)源分布

山東長(zhǎng)島國(guó)家級(jí)自然保護(hù)區(qū)（以下簡(jiǎn)稱長(zhǎng)島保護(hù)區(qū)）是位于渤海海峽南部的國(guó)家級(jí)自然保護(hù)區(qū)，具有重要的生態(tài)和旅游價(jià)值，是鳥類和太平洋斑海豹生存和繁殖的重要區(qū)域。但是，保護(hù)區(qū)海域時(shí)常受到來源不明的溢油污染物的影響[2-5]，不僅給當(dāng)?shù)氐臐O業(yè)、旅游業(yè)等產(chǎn)業(yè)造成巨大損失，而且危及鳥類和斑海豹等保護(hù)物種的生存。溢油可能會(huì)隨著潮流、風(fēng)等到達(dá)保護(hù)區(qū)海域，也可能是海域內(nèi)的船舶泄漏造成，因此很難從溢油來源進(jìn)行污染預(yù)測(cè)，給風(fēng)險(xiǎn)分析帶來挑戰(zhàn)。

以往研究中，海洋溢油溯源數(shù)值模擬方法常常被用于解決一類污染受體已知、溢油源未知（或過多）情景的災(zāi)害預(yù)防問題，溯源方法在這類情景中所需要花費(fèi)的計(jì)算時(shí)間少[6]，并且可以提供更為全面的潛在溢油泄漏情景和污染物擴(kuò)散軌跡信息[7-9]。不同于一些隨機(jī)和典型情景的溢油事故模擬方法[10-11]，溯源方法從污染受體位置出發(fā)開始，采用溯源軌跡對(duì)溢油的可能來源進(jìn)行追蹤，判斷可能對(duì)污染受體造成危害的溢油來源區(qū)域和概率，從而對(duì)高概率區(qū)域和低概率區(qū)域執(zhí)行相應(yīng)等級(jí)的監(jiān)測(cè)和預(yù)防措施[12]。Ciappa 等將該方法應(yīng)用在地中海的海洋保護(hù)區(qū)溢油風(fēng)險(xiǎn)分析中，評(píng)估了該區(qū)域來自船舶的潛在污染風(fēng)險(xiǎn)和溢油監(jiān)測(cè)的重點(diǎn)海域[8]。高松等基于該方法對(duì)秦皇島海域無主漂油的可能來源區(qū)域進(jìn)行了估計(jì)，為溢油來源的排查和應(yīng)急處置提供了幫助[3,13]。位于渤海的長(zhǎng)島保護(hù)區(qū)海域位置特殊，受到周期性變化的強(qiáng)潮流和季風(fēng)的共同影響，在進(jìn)行溢油污染風(fēng)險(xiǎn)分析時(shí)應(yīng)當(dāng)盡可能全面地考慮這些狀況。

本文采用海洋溢油溯源數(shù)值模擬方法對(duì)長(zhǎng)島保護(hù)區(qū)海域進(jìn)行溢油災(zāi)害風(fēng)險(xiǎn)評(píng)估和分析。研究使用FVCOM模擬的潮流數(shù)據(jù)和ECMWF提供的ERA5海面10 m 風(fēng)數(shù)據(jù)，在溢油模型GNOME （General NOAA Oil Modeling Environment）中進(jìn)行溯源模擬。首先，通過一次真實(shí)溢油事故的后報(bào)驗(yàn)證了溯源模擬方法的有效性，然后將這套資料和方法應(yīng)用于渤海長(zhǎng)島保護(hù)區(qū)溢油災(zāi)害風(fēng)險(xiǎn)的分析，評(píng)估了月平均統(tǒng)計(jì)條件下保護(hù)區(qū)溢油污染的來源分布和危害風(fēng)險(xiǎn)，探討了該方法的適用性，以期為我國(guó)海洋溢油風(fēng)險(xiǎn)評(píng)估相關(guān)的研究工作提供參考。

1 研究區(qū)域和數(shù)據(jù)

1.1 研究區(qū)域

渤海海峽是我國(guó)渤海的唯一出口，地理位置特殊。長(zhǎng)島保護(hù)區(qū)海域位于渤海海峽南側(cè)，自北向南與老鐵山水道、北砣磯水道和長(zhǎng)山水道等主要航道相鄰。長(zhǎng)島保護(hù)區(qū)海域受強(qiáng)潮流影響，海流多為東-西向往復(fù)流。區(qū)域盛行風(fēng)向隨季節(jié)轉(zhuǎn)換而有明顯變化，全年中春、秋、冬季的季風(fēng)強(qiáng)度較大且持續(xù)時(shí)間較長(zhǎng)，夏季大風(fēng)天氣較少。

在長(zhǎng)島保護(hù)區(qū)海域的三處位置周圍選取了66個(gè)位置并設(shè)置為污染受體地點(diǎn)，即圖2（b）中A、B和C 三處位置所示的受體點(diǎn)。對(duì)于A、B 和C 處的每一部分，相鄰受體點(diǎn)之間的距離是相等的。研究不設(shè)置潛在溢油事故的發(fā)生地點(diǎn)，油粒子從受體點(diǎn)位置釋放，隨時(shí)間逆向回溯其潛在的來源地點(diǎn)和粒子溯源軌跡。

圖2 研究區(qū)域

1.2 數(shù)據(jù)

海面10 m 高的風(fēng)場(chǎng)數(shù)據(jù)來自ECMWF 提供的第五代大氣再分析產(chǎn)品ERA5再分析數(shù)據(jù)（https://cds.climate.copernicus.eu/#!/home）。該產(chǎn)品同化了衛(wèi)星數(shù)據(jù)、現(xiàn)場(chǎng)觀測(cè)數(shù)據(jù)等多種數(shù)據(jù)，具有較高的精度。本文采用2000－2019年20年間的逐時(shí)風(fēng)場(chǎng)資料，空間分辨率為1/4°×1/4°。在風(fēng)險(xiǎn)分析的研究部分，為體現(xiàn)保護(hù)區(qū)海域風(fēng)場(chǎng)季節(jié)性差異和多年平均狀態(tài)，本文對(duì)數(shù)據(jù)進(jìn)行了月平均的預(yù)處理。

蓬萊19-3 溢油事故的衛(wèi)星遙感圖像數(shù)據(jù)來自歐空局ENVISAT 衛(wèi)星ASAR（Advanced Synthetic Aperture Radar）圖像，并采用TCNNA（Texture-Classifying Neural Network Algorithm）算法[14]提取了圖像中油膜的范圍信息，與前人的研究結(jié)果較為一致[15-16]。為了驗(yàn)證模型的可靠性，本文采用同一次溢油不同時(shí)刻的兩張衛(wèi)星圖像，以較晚時(shí)刻的圖像中獲得的溢油時(shí)間和位置信息作為溯源模擬的初始條件，以較早時(shí)刻的圖像中獲得的溢油時(shí)間和位置信息作為溯源模擬結(jié)果的驗(yàn)證。

2 分析方法和模型設(shè)置

2.1 分析方法

受體模式（receptor mode）由Galt 等[17]提出并較早應(yīng)用于溢油溯源模擬中。受體模式基于拉格朗日觀點(diǎn)，將溢油視為以油粒子的形式運(yùn)動(dòng)，在溯源模擬中會(huì)根據(jù)平流方程，逆著洋流和海面風(fēng)作用方向，在時(shí)間上逆向反推油粒子的來源位置和溯源軌跡。在受體模式中，溢油受海流影響的輸運(yùn)部分被簡(jiǎn)化為位移，即油粒子代表的浮油隨著水體以拉格朗日方式運(yùn)動(dòng)：

式中：x代表油粒的位置；為表層海流速度。根據(jù)Galt 等的描述，風(fēng)的效應(yīng)和溢油粒子的擴(kuò)散通過額外的方法模擬，溢油輸運(yùn)的最終結(jié)果是矢量疊加的效果。Yu等[18]將這個(gè)過程表述為：

式（2）中：p(t)為監(jiān)測(cè)到的溢油粒子位置；p(t-1)表示僅受海流和風(fēng)影響下的前一時(shí)刻的位置；和uc(t-1)表示每個(gè)時(shí)間步長(zhǎng)中由風(fēng)和海流引起的速度。式（3）中：R(t)表示隨機(jī)行走產(chǎn)生的距離半徑；uR()t-1為隨機(jī)行走的速度。式（4）中：St-1則表示由隨機(jī)行走造成的溢油粒子的可能影響范圍。通過式（2）和（4）可以得到溢油來源的高概率區(qū)域。

由于油粒子的隨機(jī)行走過程在溯源過程中是不可逆的，因此應(yīng)根據(jù)實(shí)際情況判斷回溯模擬的參數(shù)以滿足對(duì)當(dāng)?shù)丨h(huán)境特征的適用性。關(guān)于溯源模擬中的隨機(jī)行走帶來的問題可以參考以往研究[6,19]。

本文根據(jù)溢油溯源軌跡分析了可能對(duì)長(zhǎng)島保護(hù)區(qū)造成污染危害的來源區(qū)域，以及來源區(qū)域的溢油到達(dá)保護(hù)區(qū)的概率分布。首先，在長(zhǎng)島保護(hù)區(qū)海域的66 個(gè)受體點(diǎn)（圖2（b））釋放了用于溯源的油粒子，采用月平均的風(fēng)數(shù)據(jù)和隨時(shí)間變化的潮流數(shù)據(jù)模擬得到研究區(qū)域的溢油溯源軌跡。然后，從這些軌跡中分析可能對(duì)保護(hù)區(qū)造成污染影響的溢油來源位置，將這些信息按月合并在一起，最終得到了研究區(qū)域整體的概率圖。因此，概率圖展示了保護(hù)區(qū)周圍海域某一位置的溢油到達(dá)受體點(diǎn)的概率，其由穿過統(tǒng)計(jì)單元的油粒子軌跡數(shù)和由受體點(diǎn)釋放的油粒子間的比例組成[7]。本文標(biāo)記了概率圖中置信度為95%的區(qū)域。

2.2 動(dòng)力模型及驗(yàn)證

通過FVCOM 海洋模型模擬得到了渤海海域的潮流數(shù)據(jù)。FVCOM 是基于三維海洋原始方程的非結(jié)構(gòu)網(wǎng)格、有限體積的海洋模式[20]，結(jié)合了有限元法擬合邊界、局部加密和有限差分法計(jì)算高效等優(yōu)點(diǎn)。本文建立了渤海三維正壓潮汐模型，模型共生成32 237 個(gè)三角元與17 173 個(gè)節(jié)點(diǎn)。為保證長(zhǎng)島保護(hù)區(qū)海域網(wǎng)格質(zhì)量，對(duì)該海域進(jìn)行了網(wǎng)格加密，使海域網(wǎng)格最高空間分辨率在100 m 以內(nèi)。模型開邊界由1/30°分辨率的TPXO9主要分潮的潮汐水位驅(qū)動(dòng)[21-22]，運(yùn)行時(shí)間為2011 年4 月1 日至7 月31 日。此外，在風(fēng)險(xiǎn)分析的研究部分采用模型模擬的5 月1 日至7 月1 日的潮流數(shù)據(jù)，從而體現(xiàn)出研究區(qū)域漲落潮、大小潮等潮汐狀態(tài)下溢油溯源模擬的結(jié)果。

通過事故同期的潮汐資料對(duì)模型水動(dòng)力模擬結(jié)果進(jìn)行驗(yàn)證，結(jié)果如圖3 所示。本文驗(yàn)證了2011年6月7－12日P1位置潮位資料與2011年5月18日P2位置潮流資料，潮位與潮流的模式計(jì)算結(jié)果與事故同期實(shí)測(cè)資料的比對(duì)結(jié)果良好。

圖3 FVCOM 模式計(jì)算網(wǎng)格及驗(yàn)證

2.3 溢油模型

GNOME 是美國(guó)國(guó)家海洋和大氣管理局（National Oceanic and Atmospheric Administration，NOAA）開發(fā)的溢油建模工具，可對(duì)溢油漂移和擴(kuò)散軌跡進(jìn)行預(yù)測(cè)[23]。GNOME 在診斷模式（diagnostic mode）下可輸出“最佳預(yù)測(cè)軌跡”和“最小遺憾軌跡”。前者默認(rèn)輸入信息完全正確，代表了對(duì)溢油去向的最佳判斷，后者在模型中包含不確定性，顯示可能受到溢油影響的更大區(qū)域。本文采用了“最小遺憾軌跡”模擬結(jié)果進(jìn)行分析。

GNOME 自 早 期 版 本OSSM （On-Scene Spill Model）開始采用受體模式，并支持用于溢油溯源的“Run Backwards”模塊。在本文中，用于溯源模擬的油粒子是在圖2（b）所示的受體點(diǎn)位置釋放的。模型運(yùn)行中設(shè)置的時(shí)間步長(zhǎng)為30 min，每個(gè)受體點(diǎn)在一個(gè)時(shí)間步長(zhǎng)釋放3 個(gè)油粒子，并在一個(gè)月中連續(xù)釋放15 日。不考慮油粒子風(fēng)化作用和歸宿，單個(gè)粒子在模擬中設(shè)置的有效時(shí)長(zhǎng)為30日，即模擬時(shí)長(zhǎng)超過30 日的粒子會(huì)被剔除。因此，模擬過程的整體有效模擬周期為45日。

3 結(jié)果與討論

3.1 模擬結(jié)果驗(yàn)證

本研究首先對(duì)蓬萊19-3 油田溢油事故中溢油真實(shí)的漂移運(yùn)動(dòng)情況進(jìn)行了后報(bào)，其目的是驗(yàn)證數(shù)據(jù)和溢油溯源模擬方法的有效性，并確定符合本文研究區(qū)域環(huán)境特征的溢油溯源模擬的適用參數(shù)。蓬萊19-3 溢油事故最早發(fā)生于2011 年6 月4日，持續(xù)泄漏時(shí)間在2個(gè)月以上，約700桶原油和2 500 桶含油泥漿泄漏到海水中。6 月19 日，兩張不同時(shí)刻的衛(wèi)星圖像記錄了事故中溢油在海面的位置和分布情況，圖4（b）和圖4（c）展示了提取的溢油位置和范圍信息。

圖4 蓬萊19-3事故中不同時(shí)刻的溢油位置

以較晚時(shí)刻的圖像中獲得的溢油時(shí)間和位置信息作為模擬的初始條件（圖5（a）），溯源模擬的開始時(shí)刻為2011 年6 月19 日9 時(shí)0 分。以較早時(shí)刻的圖像中獲得的溢油時(shí)間和位置信息作為溯源模擬結(jié)果的驗(yàn)證，模擬結(jié)束的時(shí)刻為2011 年6 月19 日2 時(shí)21 分。圖5（b）展示了本文采用的溢油溯源模擬參數(shù)方案在事故后報(bào)中得到的模擬結(jié)果，且與真實(shí)位置接近。結(jié)果表明本文的資料和方法可以在短期溯源模擬中得到較好的結(jié)果。

圖5 溢油事故溢油溯源模擬結(jié)果

由于污染風(fēng)險(xiǎn)分析采用溯源方法評(píng)價(jià)的研究尚不多見，對(duì)于溯源模擬結(jié)果的精度評(píng)價(jià)案例相對(duì)較少。因此，本文在表1 中給出了通過此次模擬所優(yōu)化后的相關(guān)參數(shù)取值。在GNOME 中，各種輸入項(xiàng)的不確定性部分可以設(shè)置為一個(gè)區(qū)間內(nèi)的隨機(jī)值。比如，風(fēng)漂移系數(shù)默認(rèn)在1%～4%中取值，并且該值會(huì)每15 min重置一次。

表1 溯源模擬中相關(guān)參數(shù)取值

3.2 溢油概率分布模擬

GNOME 以月平均的風(fēng)數(shù)據(jù)和隨時(shí)間變化的潮流數(shù)據(jù)作為驅(qū)動(dòng)，在圖2（b）所示的受體點(diǎn)位置釋放了用于溯源的油粒子，每個(gè)受體點(diǎn)每日的釋放頻率為144 個(gè)，一個(gè)月中連續(xù)釋放時(shí)長(zhǎng)為15 日，每個(gè)粒子的有效模擬時(shí)長(zhǎng)為30 日。然后，從66×144×45 條溯源軌跡中得到可能對(duì)保護(hù)區(qū)造成污染影響的溢油來源位置，這些信息按月合并在一起，最終得到了研究區(qū)域整體的概率圖。

溯源軌跡的概率圖顯示了保護(hù)區(qū)周圍海域某一位置的溢油到達(dá)保護(hù)區(qū)造成污染的概率。圖6不同顏色的區(qū)域表示了對(duì)應(yīng)位置溢油的污染概率，黑色實(shí)線標(biāo)記了圖中置信度為95%的區(qū)域。結(jié)果表明，①受強(qiáng)潮流作用影響，區(qū)域內(nèi)東-西走向的往復(fù)流運(yùn)動(dòng)使東岸和西岸的區(qū)域周圍常年表現(xiàn)為高概率區(qū)域，也使保護(hù)區(qū)內(nèi)島嶼的東西兩岸更容易受到溢油污染物的影響。②受季節(jié)性變化的典型風(fēng)場(chǎng)影響，高概率區(qū)域的主要分布態(tài)勢(shì)呈現(xiàn)出較為明顯的季節(jié)性差異。高概率區(qū)域的分布形態(tài)在夏季向東南方向延伸，對(duì)應(yīng)季節(jié)中盛行東南風(fēng)，在冬季向西北方向延伸，對(duì)應(yīng)季節(jié)中盛行西北風(fēng)。春夏季節(jié)，保護(hù)區(qū)內(nèi)島嶼更容易受到東部和南部海域溢油污染的影響，在對(duì)應(yīng)季節(jié)中溢油可能在島嶼的東部和南部更早被發(fā)現(xiàn)，而在溢油污染事件多發(fā)的秋冬季節(jié)，保護(hù)區(qū)內(nèi)島嶼更容易受到來自西北方向溢油的污染影響，季風(fēng)會(huì)導(dǎo)致來自渤海中部及北部海域的油污更容易到達(dá)保護(hù)區(qū)，這可能是在保護(hù)區(qū)內(nèi)島嶼的西海岸往往更早發(fā)現(xiàn)溢油的原因。③在表2 中給出了置信度為95%區(qū)域的風(fēng)險(xiǎn)分析表，從區(qū)域的分布來看，受主導(dǎo)風(fēng)向、往復(fù)流的共同影響，在1－2 月以及9－12 月，保護(hù)區(qū)周圍的高概率區(qū)域向保護(hù)區(qū)北部延伸，分布面積較大，并在12 月達(dá)到最大。而在4－8 月保護(hù)區(qū)南部的近岸部分成為高概率區(qū)域，分布面積較小。高概率區(qū)域的覆蓋面積表現(xiàn)為春夏季節(jié)較小、秋冬季節(jié)較大的特征。

表2 溢油污染風(fēng)險(xiǎn)分析

圖6 月平均條件下長(zhǎng)島保護(hù)區(qū)海域溢油溯源模擬結(jié)果的概率圖

4 結(jié)論

本文基于海洋溢油溯源數(shù)值模擬方法，評(píng)估月平均統(tǒng)計(jì)條件下長(zhǎng)島保護(hù)區(qū)海域整體污染危害風(fēng)險(xiǎn)，分析可能對(duì)保護(hù)區(qū)造成危害的來源區(qū)域及其溢油對(duì)保護(hù)區(qū)造成污染影響的概率等信息。模擬結(jié)果表明，①受到長(zhǎng)島保護(hù)區(qū)海域強(qiáng)潮流作用和主導(dǎo)風(fēng)向的顯著影響，高概率區(qū)域的分布常年表現(xiàn)為東-西走向的分布態(tài)勢(shì)，并在夏季向東南方向延伸，最遠(yuǎn)可延伸至山東煙臺(tái)市北岸；在冬季向西北方向延伸，最遠(yuǎn)可延伸至渤海中部和北部。②高概率區(qū)域的區(qū)域面積表現(xiàn)為春夏季節(jié)較小、秋冬季節(jié)較大的特征。溢油污染的最大風(fēng)險(xiǎn)出現(xiàn)在冬季。