綠華島周邊海域溢油擴散風險數值模擬研究

顧一凡，王彬諭，倪云林，陳 維,2*

(1.浙江海洋大學 海洋工程裝備學院，浙江 舟山 316022；2.自然資源部 海洋空間資源管理技術重點實驗室，浙江 杭州 310012；3.廣西交通職業(yè)技術學院，廣西 南寧 530023)

引 言

工程施工過程中可能發(fā)生的海洋溢油事故會給海洋生態(tài)環(huán)境帶來嚴重的影響[1]。溢油在海面上擴散形成大面積的油膜覆蓋于海面，隔絕了水體與大氣的氣體交換進程，導致水體缺氧，降低海洋的生產力；同時油膜在風和潮流等作用下擴散到海岸，將破壞海產養(yǎng)殖和鹽田生產，污染濱海旅游區(qū)等，造成嚴重的經濟損失。另一方面，海上溢油中的芳香烴化合物能在生物體內長期積累，最終將危害人類健康[2]。本文以綠華—花鳥海底電纜工程項目(圖1)為例，分析工程中可能發(fā)生的溢油擴散情況。

圖1 項目位置圖

綠華島及其周邊海域有較多數量的環(huán)境敏感點如：嵊泗縣馬鞍列島海洋牧場示范區(qū)(圖2)，嵊泗綠華養(yǎng)殖區(qū)(圖2中編號為①～⑦的區(qū)域)以及花鳥島海水養(yǎng)殖區(qū)(圖2編號為A—F的區(qū)域)。為了在事故發(fā)生后及時采取相應的措施減少溢油對環(huán)境的影響，研究溢油在水體中的運動軌跡，建立溢油模型對溢油事故進行分析，制定合理的預防、治理方案是必要的。

圖2 工程海域環(huán)境敏感點及測點分布

溢油擴散是一個非常復雜的過程，受油品的物理、化學性質及事故發(fā)生地的氣象條件和水動力特性等共同影響，進入水體的石油將發(fā)生擴展、漂移、擴散、蒸發(fā)、溶解、乳化、沉降、擱淺以及在生物作用下的降解等變化過程[3]。同時油膜與水體和大氣之間發(fā)生復雜的物理化學變化，表面張力等屬性也隨著組分和溫度的變化發(fā)生不斷變化[4]。在90年代，Cuesta[5]等運用FORTRAN語言編寫建立二維溢油模型對浮油的漂移情況進行預測。Sugioka[6]等建立了3-D溢油模型，分析了油粒子在水平和垂直空間內的運動軌跡。隨著時間的推移，更多的溢油模型被開發(fā)，如GNOME模型[7]、OILMAP模型[8]、MS4模型[9]、Navy模型[10]，Delft3D-PART模塊[11]以及MIKE 21模型[12]等。后來Marques[13]等開始將Telemac3D和Ecos模型耦合起來進行溢油事故模擬。郭健[14]等通過運用MIKE-SA模塊分析了各種海況下油粒子的漂移軌跡、油膜掃海面積?？锎淦糩15]等運用Delft3D水動力模型的粒子追蹤模塊，通過釋放大量粒子云團來模擬溢油的運移與歸宿。

本文基于綠華島海底電纜敷設工程，建立綠華島及周邊海域的MIKE二維水動力潮流模型，并根據實測數據進行驗證，在此基礎上建立溢油事故預測分析模型，通過釋放大量粒子模擬溢油的運移與歸宿，分析對周邊環(huán)境敏感區(qū)域的影響。

1 潮流數值模型

1.1 模型簡介

MIKE軟件是丹麥水力研究(DHI)開發(fā)的數值模型，包含多種模塊，用于模擬各種情況下的水動力環(huán)境，有著十分強大的前后處理能力[16]。邊界處可以根據具體需要自由定義，如岸線、海洋、島嶼、水利建筑物等。

粒子追蹤模塊(PT)是MIKE21所包含的用于水質和環(huán)境評價的模塊之一，主要用于模擬懸浮物質的擴散軌跡和狀態(tài)，可用于預報溢油問題、評估意外溢油預案分析。

1.2 模型區(qū)域和網格

模型計算范圍西起倉前，北至蘆潮港以北，南臨象山，東到東經124°，包含了杭州灣、舟山群島海域，計算域的橫向寬約為378 km，縱向長度為216 km，計算面積約為81 648 km2。由于計算海域所處舟山群島，島嶼眾多，擁有漫長曲折的岸線，近岸海域地形復雜，模型采用三角形網格，使模型更能真實地模擬岸線情況。為了提升模型的計算效率，遠離工程區(qū)的網格分辨率在200～4 000 m之間。工程區(qū)附近海區(qū)島嶼眾多、地形和潮流結構復雜，潮動力對環(huán)境的影響較強。為了獲得足夠的計算精度，對工程區(qū)網格進行加密處理，網格尺度最小為10 m，能夠對項目海域水下地形及岸線進行較好地刻畫。模型共有22 147個節(jié)點，42 088個單元，計算網格見圖3。

圖3 計算區(qū)域和網格

水下地形資料的精確性對模型計算有著極其重要的影響。計算域內大范圍水下地形由海軍航保部海圖通過GIS數字化得到，工程附近海域采用最新測量的水下地形CAD數據，所有數據基面均統(tǒng)一國家八五高程基準面。

1.3 相關參數設定

模型的時間步長設置由模型自動調節(jié)，設置的范圍為0.000 1～30 s，柯朗數限值為0.8，床阻糙率采用曼寧系數，數值范圍為0.012～0.014，河流邊界采用流量控制，外海開邊界采用東中國海大模型[17]計算得到的潮位控制，大模型的開邊界僅考慮外海的開邊界，外海開邊界采用全球潮汐模型中提取的潮位數據，模型考慮的分潮包含M2，S2，K2，N2，S1，K1，O1，P1，Q1和M4，共計10個分潮。陸地及島嶼設為不透水，水質點沿岸線的切線方向自由滑移的岸線邊界。在近岸區(qū)域內存在淺灘，隨著潮水漲落會出現淹沒和露灘交替的現象，于是采用干濕網格法進行處理，取干點臨界水深0.005 m，濕點臨界水深0.05 m。

1.4 模型驗證

潮位觀測布設1個臨時驗潮站，觀測為期12 d，覆蓋了大、中、小潮，臨時水尺高程精度不低于圖根水準要求，采用1985國家高程基準。潮流觀測布設6個站點，采用船舶觀測方式分別在大、小潮期間連續(xù)觀測27 h。通過分析海域的潮流特征和實測數據對模型的可靠性進行驗證、評估。驗證資料所采用的潮流站和潮位站的具體位置如圖2所示。

實測潮位變化具有明顯規(guī)律，測區(qū)中實測潮位的變化過程多呈現為正規(guī)半日潮的特征。測區(qū)中潮流性質歸屬于非正規(guī)半日淺海潮流的類型。

圖4為潮位驗證過程，大、中潮期間計算與實測值有些微小的相位差，計算值的振幅在整個潮周期上稍大于實測值；圖5和圖6所示的潮流計算值與實測值基本吻合。

圖4 綠華島臨時驗潮站潮位驗證

圖5 綠華島測潮站大潮流速流向驗證

圖6 綠華島測潮站小潮流速流向驗證

表1 模型效率系數

1.5 潮流現狀模擬及分析

通過模擬結果可知，維系大區(qū)域模型流場的主要動力是來自東海并由東南和東方向傳入的潮波。區(qū)域模型的流速具有自東向西逐漸增大的特征。圖7為工程區(qū)域的漲急、落急流矢圖。由圖可見，位于綠華-花鳥島之間的工程區(qū)域，其東西兩側基本為開闊海域，受狹道地形制約，漲落潮基本呈往復運動，近岸區(qū)漲落潮流流向大致與岸線走向平行，且主槽的流速明顯大于邊灘；漲落潮路徑基本一致。

圖7 工程區(qū)域附近漲落潮流場圖

2 溢油模型

溢油事故預測采用Johansen[19]等提出的“油粒子”模式，認為海面上的油膜是由大量油粒子組成，每個油粒子代表一定的油量，油粒子之間彼此互相獨立、互不干擾，其運動過程分為平流和擴散過程。它們匯聚所構成的“云團”就是油膜，油膜的厚度則通過一定水面面積內的油粒子量來體現，它在宏觀上體現著油粒子的擴散運動。對于小規(guī)模溢油，油膜自身擴展過程時間較短，在水流的紊動作用下開始分散，其運動方式以擴散為主，采用“油粒子”模型顯得更為合理。

“油粒子”模型中油膜的運動和形變都是通過油粒子來體現。在運動過程中，用確定性方法——拉格朗日法模擬平移過程，擴散過程則運用不確定方法——隨機走動法來模擬[20]。油粒子可通過不同的擴散系數來模擬油膜自身擴展過程與油膜在流場和風場影響下的紊流擴展過程[21]，并能確切地模擬油膜邊緣擴展過程和油膜形狀在風向上的明顯拉伸過程[22]。

在風和流的共同作用下，油粒子群的每一個油粒子的運動可用下式表示：

(1)

2.1 工況選取

根據施工船4個油艙總容量為100 m3，選取最大不利情況下即4個油倉全破情況下，根據比重可計算出發(fā)生溢油80 t，釋放的總粒子數為90 000個，顯示頻率為100。油粒子將在溢油開始1 h內釋放完畢。

在航道處發(fā)生碰撞后發(fā)生泄漏，泄漏時間在平潮后2小時(流速最小時刻)，預測年常風向為SE與NW，常風向年平均風速7.1 m/s。本文對靜風、SE與NW三種風況進行了模擬。

2.2 風拖曳系數選取

風海流采用如下計算公式：U=α·W10Δf(θ)，式中α為風拖曳系數，f(θ)為科氏力引起的偏轉角的函數，θ為偏轉角。

風拖曳系數是海洋大氣物理學中的重要參數，本文中采用WuJin公式，即：

α=CaW10

(2)

α=CbW10>Wb

式中，Ca=1.255×10-3，Cb=2.425×10-3，Wa=7 m/s，Wb=25 m/s 。

在北半球，風海流向右偏轉于風向，本報告中偏轉角取15°，

3 計算結果分析

根據溢油計算結果，分別統(tǒng)計了48 h內各風況下油粒子在漲、落潮時段泄露后的散布情況；及6 h、12 h、24 h、48 h油膜的面積及油粒子的掃海面積，見表2。具體分析如下：

表2 各工況下油膜面積及油粒子掃海面積

3.1 漲潮階段溢油

溢油在漲潮發(fā)生后，油粒子受潮流的影響，油膜中心向西漂移。3 h后，由于潮流的作用，油粒子開始向東北轉向，在模擬的各種風況下，敏感區(qū)域均未受溢油影響。6 h后油粒子在落潮潮流的影響下回返運動，油膜面積約為12 km2，影響范圍約為21 km2。

靜風情況下，12 h后，油粒子隨潮流向東南方向漂移，部分油粒子在經過東綠華島后向綠華島南部漂移，部分停留在綠華島西部；24 h后，北部的油粒子向東擴散，南部的油粒子持續(xù)向南擴散，油膜面積擴展至35.2 km2，48 h后完全覆蓋綠華島及花鳥島周邊海域，油膜面積為58.57 km2，影響242.23 km2海域，見圖8(d)、圖8(e)和圖8(f)。

圖8 靜風漲潮階段溢油后油粒子的分布

SE風況下油粒子隨潮流的漂移受阻，持續(xù)影響區(qū)域D和E，見圖9(d)；24 h后，油粒子向綠華島的漂移受阻，油膜主要覆蓋了花鳥島東側以外的海域，面積為33.61 km2，影響了綠華島東部養(yǎng)殖區(qū)及海域牧場示范區(qū)部分海域，見圖9(e)，119.18 km2海域受到溢油影響；48 h后油膜面積為57.06 km2，完全覆蓋花鳥島及其臨近海域，對綠華島東北部約252.54 km2海域造成影響。

圖9 SE風況漲潮階段溢油后油粒子的分布

NW風況下，12 h后油粒子加速向東南方向漂移，并影響到海洋牧場示范區(qū)部分區(qū)域，見圖10(d)。24 h后，油膜面積擴展至36.21 km2，擴散至東西綠華島和花鳥島附近的大部分海域，并在48 h后油膜面積為56.97km2，完全覆蓋綠華島及花鳥島附近海域，漂移至綠華島東南部并在潮流和風的影響了綠華島南部大片海域，見圖10(e)與圖10(f)，影響247.49 km2海域。

圖10 NW風況漲潮階段溢油后油粒子的分布

3.2 落潮階段溢油

落潮時刻溢油后，油粒子受潮流的影響，向東南方向擴散。3 h后，油粒子在潮流的影響下開始向南漂移并對部分敏感區(qū)域造成影響，見圖11(b)。6 h后，在潮流的作用下向西漂移動，見圖11(c)，覆蓋面積為10.9 km2，23 km2的海域受到溢油的影響；12 h后油粒子向西漂移，油膜擴展至19.18 km2，影響綠華島南側全部養(yǎng)殖區(qū)，見圖11(d)，SE風況下油膜順流加速擴散，油膜面積擴大，NW風況下油膜擴展受限，影響范圍減小。

在靜風情況下，24 h后，油膜擴展至34.31km2，油粒子在主要集中在西綠華島南岸并向綠華島北部漂移，見圖11(e)，影響海域為96.84 km2；48 h后油粒子主要在綠華島南部東西向反復漂移，油膜擴展至53.32 km2，約171.87 km2海域受到溢油影響，少數油粒子漂移并影響綠華島養(yǎng)殖區(qū)D。

圖11 靜風落潮階段溢油后油粒子的分布

在SE風況下，24 h后少數粒子漂移至花鳥島附近，見圖12(e)，油膜進一步擴展至38.87 km2，95.56 km2海域受到溢油影響；在48 h后覆蓋范圍達到54.39km2，影響綠華島南部及其西北部大片海域，面積達185.84 km2并對花鳥島養(yǎng)殖區(qū)D產生一定影響，見圖12(f)。

圖12 SE風況落潮階段溢油后油粒子的分布

在NW風況下，24 h后大部分油粒子向東南方向漂移并遠離綠華島，油膜面積為32.36 km2，掃海面積為105.4 km2，見圖13(e)，在48 h后油膜擴展至53.42 km2，油粒子整體漂移至嵊泗島東部海域，對綠華至嵊泗之間的大片海域造成影響，影響范圍明顯擴大，受影響海域面積達204.38 km2，整體呈放射狀擴散，見圖13(f)。

4 結論

本文運用MIKE模型，用大量的油粒子模擬溢油在海面上形成的油膜，計算模擬其運動軌跡及歸宿。分析油粒子對溢油點周邊海域內環(huán)境敏感點的影響。模擬結果顯示：

(1) 在溢油初期油粒子聚集在一起，潮流對油粒子的運動軌跡起決定性作用，但隨著油粒子自身的擴展，油膜厚度減小，在大約6 h后油粒子的漂移軌跡在受風影響與靜風情況下的漂移軌跡產生差別。若溢油發(fā)生在漲潮時刻，油膜擴散范圍和油粒子的掃海面積均明顯大于在落潮時刻發(fā)生的溢油情況，見表2。

(2) 若溢油發(fā)生在漲潮時，在6 h內油粒子不會對敏感區(qū)域造成影響，12 h后溢散的油粒子主要向花鳥島移動并影響其環(huán)境敏感區(qū)域，NW風況下同時對綠華島及其周邊區(qū)域造成威脅；若溢油發(fā)生在落潮時，油粒子將在3 h后影響海洋牧場示范區(qū)，12 h內對綠華島養(yǎng)殖區(qū)造成較大影響，但在48 h內不會威脅到花鳥島養(yǎng)殖區(qū)。