牟林，武雙全，宋軍，李歡，劉首華，李琰，高佳

（國家海洋信息中心，天津 300171）

渤海海域溢油應(yīng)急預(yù)測預(yù)警系統(tǒng)研究I. 海洋動(dòng)力要素預(yù)測技術(shù)研究

牟林，武雙全，宋軍，李歡，劉首華，李琰，高佳

（國家海洋信息中心，天津 300171）

采用“油粒子”模型，模擬溢油在海洋中的漂移擴(kuò)散運(yùn)動(dòng)，將海流、波浪和風(fēng)場的模擬結(jié)果作為溢油漂移計(jì)算的外強(qiáng)迫，來驅(qū)動(dòng)溢油漂移模型，計(jì)算其漂移路徑。把計(jì)算的溢油漂移結(jié)果與地理信息系統(tǒng)等平臺(tái)相結(jié)合，建立了渤海海域應(yīng)急預(yù)測預(yù)警系統(tǒng)。本文作為該系統(tǒng)研制中的一部分，主要介紹了海流數(shù)值計(jì)算模型、波浪數(shù)值計(jì)算模型、風(fēng)場數(shù)值計(jì)算模型和溢油預(yù)測模型的研究工作。

油粒子；海流；波浪；溢油預(yù)測

1 引 言

隨著全球經(jīng)濟(jì)和社會(huì)的發(fā)展，需油量不斷增加，石油運(yùn)輸業(yè)和開采業(yè)的規(guī)模不斷擴(kuò)大，海上溢油已經(jīng)帶來了多起災(zāi)難性事故。溢油事故發(fā)生后，油膜在海面的漂移主要受制于表層流的驅(qū)動(dòng)，而其上邊界直接受到風(fēng)應(yīng)力的作用。海面某水體微團(tuán)，在水中的運(yùn)動(dòng)由三部分組成，在不計(jì)彼此間相互作用的前提下，包括風(fēng)生流、非風(fēng)生流和波余流[1]。風(fēng)生流與非風(fēng)生流是在遠(yuǎn)大于油膜尺度的驅(qū)動(dòng)力如引潮力、壓強(qiáng)梯度力和風(fēng)應(yīng)力等作用下形成的海水運(yùn)動(dòng)，它們并不會(huì)因油膜的存在而產(chǎn)生較大的變化。波余流則不同，根據(jù)Stokes理論，波余流的量值可以達(dá)到風(fēng)速的 2%。但是油膜的存在使得表面張力增加，海面趨于平坦，波浪的非線性作用大為削弱，因此，事實(shí)上波余流是可以忽略的[2]。因此，在近海海域，潮流和風(fēng)生流是決定溢油漂移的最重要的因素。

實(shí)際觀測表明，在開闊海域，溢油的漂移速度主要由風(fēng)的作用決定；而在近海或沿岸時(shí)，潮流的作用就不可忽視；尤其在港灣或碼頭，潮流的作用更加重要。風(fēng)對油膜的影響一般按風(fēng)速的 3%～3.5%計(jì)算，而在封閉或者半封閉的海域，則取更小值。如果溢油事故發(fā)生海域受沿岸流和深海流的影響不大，用上述計(jì)算方法可獲得較好的結(jié)果。

溢油預(yù)報(bào)需要高分辨率高精度的海面風(fēng)、表層流的數(shù)值預(yù)報(bào)和完善的溢油漂移擴(kuò)散數(shù)值預(yù)報(bào)模式[3]。為此，國家海洋信息中心開展了渤海海域溢油污染預(yù)測預(yù)警技術(shù)研究工作，并基于FVCOM、WRF等高分辨率數(shù)值預(yù)報(bào)結(jié)果實(shí)現(xiàn)了對溢油污染物快速準(zhǔn)確的模擬預(yù)報(bào)及溯源回推。本文就預(yù)報(bào)預(yù)警系統(tǒng)中表層流、風(fēng)場、波浪等海洋要素的研究工作進(jìn)行了介紹、分析。

2 預(yù)測模式簡介

溢油漂移路徑數(shù)值模擬需要高精度高分辨率的海面風(fēng)場和表層流場動(dòng)力因素，溢油在海洋環(huán)境條件下的生物降解和化學(xué)過程等作用，是溢油環(huán)境影響評估的一個(gè)重要問題，但是在短時(shí)間的溢油漂移模擬問題中是可忽略不計(jì)的。因此，溢油漂移模擬的動(dòng)力要素主要包括高分辨率海面風(fēng)場、表層流和波浪等[4]。

2.1 海流模塊

2.1.1 FVCOM 模型介紹 海流模塊采用了三維原始方程組海洋數(shù)值模型FVCOM。該模型由美國麻省大學(xué)海洋科學(xué)技術(shù)學(xué)院陳長勝博士領(lǐng)導(dǎo)的研究小組與美國伍茲霍爾海洋學(xué)協(xié)會(huì)的羅伯特·C·比爾茲利博士針對近海和河口潮汐環(huán)流合作開發(fā)，其最大的特色和優(yōu)點(diǎn)是結(jié)合了有限元法易擬合邊界、可進(jìn)行局部加密以及利用有限差分法便于海洋原始方程組的離散計(jì)算。有限元法采用非結(jié)構(gòu)化的三角形網(wǎng)格，給出線性無關(guān)的基函數(shù)，求其待定系數(shù)，特點(diǎn)是三角形網(wǎng)格易擬合邊界、易進(jìn)行局部加密，可以比較完美的擬合各種復(fù)雜的岸線；而有限差分法直接離散差分海洋原始方程組，特點(diǎn)是動(dòng)力學(xué)基礎(chǔ)明確、差分直觀、計(jì)算高效[5]。FVCOM 兼有兩者的優(yōu)點(diǎn)，數(shù)值計(jì)算采用的方程的積分形式和更好的計(jì)算格式，使動(dòng)量、能量和質(zhì)量具有了更好的守恒性，利用干濕判斷法處理潮灘移動(dòng)邊界，應(yīng)用Mellor和Yamada的2.5階湍封閉子模型使整個(gè)模型在物理和數(shù)學(xué)上閉合，垂向采用 變換來體現(xiàn)不規(guī)則的底邊界，以及外模和內(nèi)模分裂以節(jié)省計(jì)算時(shí)間[6]。

2.1.2 模型配置及結(jié)果檢驗(yàn)

（1）模型配置

模型采用的溫鹽初始場為歷史觀測資料得到的多年平均溫度和鹽度場；開邊界采用潮汐調(diào)和常數(shù)預(yù)報(bào)的水位進(jìn)行強(qiáng)迫，采用6個(gè)主要半日分潮和全日分潮的調(diào)和常數(shù) M2、S2、N2、K1、O1、P1進(jìn)行預(yù)報(bào)，所用的潮汐調(diào)和常數(shù)來自于伴隨同化的渤、黃、東海潮汐模型模擬結(jié)果。模型在垂向分為11層，采用 坐標(biāo)，并且表層具有較高的分辨率。模型的外模的時(shí)間步長為 2 s，內(nèi)模的時(shí)間步長為12 s。模型每天凌晨從零點(diǎn)開始自動(dòng)預(yù)報(bào)48 h，每小時(shí)輸出1次計(jì)算結(jié)果。其中模型所采用的地形來自于山東海事局所提供的不同分辨率的海圖（如圖1）。

圖 1 處理后的高精度渤黃海地形Fig. 1 Processed topography of Bohai Sea and Yellow Sea with high precision

模型計(jì)算的海區(qū)為渤海和黃海區(qū)域，其范圍為：117.42°～126.25°E，37°～41°N，在膠州灣及青島近海、成山頭定線制水域和渤海海峽等重點(diǎn)區(qū)域進(jìn)行了網(wǎng)格加密，加密區(qū)分為3個(gè)等級：一級加密區(qū)網(wǎng)格精度為100 m，包括以下3個(gè)海域：膠州灣及青島近海、成山頭定線制水域、渤海海峽；二級加密區(qū)網(wǎng)格精度為300 m，包括以下海區(qū)：整個(gè)渤海及124°E以西海域；三級加密區(qū)為剩余的海區(qū)，網(wǎng)格分辨率為1分（如圖2所示）。分辨率在渤海海峽處的島嶼附近可達(dá)到百米的量級，基本上所有的島嶼都能保留和分辨出，避免了在以往的預(yù)報(bào)過程中，出現(xiàn)海流穿島而過的情況，而在非重點(diǎn)區(qū)域又降低了網(wǎng)格的分辨率，既保證了重點(diǎn)區(qū)域的高分辨率模擬，又控制了計(jì)算成本，使海流模塊的預(yù)報(bào)效率達(dá)到了最優(yōu)化。

圖 2 渤黃海網(wǎng)格分布圖，顏色代表水深Fig. 2 Grid distribution in Bohai Sea and Yellow Sea, where the water depth is represented by the color

（2）模型檢驗(yàn)

表 1 模型模擬的M2分潮調(diào)和常數(shù)同驗(yàn)潮站調(diào)和常數(shù)對比Tab. 1 The comparison of the M2 tide harmonic constant counted by tide station and simulated by model

調(diào)和常數(shù)檢驗(yàn)

表1為模型模擬的M2分潮調(diào)和常數(shù)與環(huán)渤、黃海各驗(yàn)潮站觀測的調(diào)和常數(shù)對比，19個(gè)驗(yàn)潮站模擬的平均振幅誤差只有2.8 cm，平均遲角誤差只有3.7°，計(jì)算穩(wěn)定后所得到的調(diào)和常數(shù)與觀測結(jié)果符合較好。它們之間的差異，一方面可能是由網(wǎng)格點(diǎn)與驗(yàn)潮站之間的位置差異所造成的，因?yàn)槟Ｐ偷乃骄W(wǎng)格距都在數(shù)千米間，這樣的網(wǎng)格對于海岸線的分辨率仍然較差；另一方面可能是由于岸界變遷以及地形和水深的誤差所致，這在黃河口和萊州灣表現(xiàn)得尤為明顯；另外還有可能與邊界處各分潮的調(diào)和常數(shù)精度等條件有關(guān)。

潮流檢驗(yàn)

選取渤海、黃海中的 A 點(diǎn)（120°43′E，39°52′N）、B 點(diǎn)（121°21′E，40°30′N）、C 點(diǎn)(117°57′E，38°39′N)、D 點(diǎn)（121°24 E，40°38′N）、E 點(diǎn)（121°24′E，40°38′N）、F 點(diǎn)（119°25′E，37°59′N）共 13個(gè)過程進(jìn)行了潮流比較，比較結(jié)果如表2所示?？梢钥闯瞿M的潮流和實(shí)際觀測的潮流較接近，二者流速、流向、高低潮流、轉(zhuǎn)流等都基本一致。

表2 渤、黃海潮流實(shí)測結(jié)果與模擬結(jié)果的對比Tab. 2 Comparison of measured tidal current data and the simulated result of Bohai Sea and Yellow Sea

通過上述調(diào)和常數(shù)、潮流、水位等結(jié)果的比較可以發(fā)現(xiàn)，F(xiàn)VCOM模型能較好地模擬黃、渤海的潮汐、潮流特征，這說明該模型在黃、渤海的應(yīng)用較為可靠，流場模擬結(jié)果將會(huì)較為真實(shí)的反映預(yù)報(bào)海域的流場特征。

2.2 海浪模塊

波浪對于溢油的影響，主要表現(xiàn)為波浪的擾動(dòng)（尤其是破碎波）。波浪的擾動(dòng)影響著溢油的破碎、分散和乳化過程，會(huì)使溢油以油粒子的形態(tài)進(jìn)入水體，因此波浪是油粒子入水率的控制因素，也是油膜發(fā)生破裂的環(huán)境動(dòng)力之一。此外，由風(fēng)浪、涌浪等非線性波產(chǎn)生的波余流對油膜的漂移過程也有一定的影響。但是，由于油膜的存在會(huì)使表面張力增加，海面趨于平滑，大大減弱了風(fēng)浪、涌浪等的非線性效應(yīng)，因此在模型計(jì)算中一般可以忽略波余流對于油膜漂移的影響。

在本系統(tǒng)中，海浪模塊為溢油行為與歸宿預(yù)測模塊提供必要的波浪要素場，以對溢油入水的過程進(jìn)行精確的模擬。

2.2.1 SWAN模型介紹 海浪模塊所采用的模式為國際上非常先進(jìn)的第三代淺海海浪模式SWAN。SWAN模式考慮了較多的物理過程，包含了當(dāng)前海浪預(yù)報(bào)研究的最新成果，其考慮的因素有：

（1）波浪的傳播過程，包括：非平穩(wěn)的水深變化引起的折射；由水底和流的變化引起的變淺作用；逆流傳播時(shí)的阻礙和反射；波浪在幾何空間的傳播；次網(wǎng)格障礙物對波浪的阻礙和波浪通過次網(wǎng)格障礙物的傳播；波生增水。

（2）波浪的產(chǎn)生和耗散，包括：風(fēng)輸入；白冠破碎；水深變淺引起的破碎；水底摩擦；三波和四波的非線性相互作用。

2.2.2 模型配置和結(jié)果檢驗(yàn) 波浪計(jì)算采用了最新的具備三角網(wǎng)格的海浪模式——SWAN，同樣對重點(diǎn)關(guān)注區(qū)域進(jìn)行網(wǎng)格加密，便于與FVCOM耦合。兩者采用海洋環(huán)境要素交換的方式進(jìn)行耦合，F(xiàn)VCOM給SWAN提供流場，SWAN給FVCOM提供波浪要素場，兩者在每個(gè)時(shí)間步長上交換數(shù)據(jù)進(jìn)行耦合，最后為溢油行為與歸宿預(yù)測模塊提供更精確的流場和波浪要素場。

利用SWAN模型對2009年10月11－12日秦皇島附近的大浪進(jìn)行了數(shù)值模擬，并將計(jì)算結(jié)果與秦皇島站2009年10月11日的波浪觀測值進(jìn)行對比，如圖3所示。可以看出，模擬結(jié)果與實(shí)測數(shù)據(jù)吻合較好。

圖 3 2009年10月11日秦皇島站實(shí)測值與模式計(jì)算值的對比Fig. 3 Comparison of measured data and the data computed by model at Qinhuangdao Station, Oct.11th ,2009

2.3 大氣模塊

2.3.1 模型介紹 海面風(fēng)場預(yù)報(bào)系統(tǒng)采用了動(dòng)力模塊十分成熟的大氣中尺度模型——WRF模型。

WRF模型（Weather Research Forecast）是由多個(gè)美國研究部門及大學(xué)的科學(xué)家共同研發(fā)的新一代中尺度預(yù)報(bào)模型。該模型為一個(gè)完全可壓、非靜力的模型，控制方程組都為通量形式。垂直坐標(biāo)采取地形追隨的靜力氣壓坐標(biāo)。采用Arakawa C型水平和垂直交錯(cuò)網(wǎng)格，有利于在高分辨率模擬中提高準(zhǔn)確性。時(shí)間積分采用完全時(shí)間分裂格式，外循環(huán)Runge-Kutta技術(shù)的較大時(shí)間步長，內(nèi)循環(huán)為聲波時(shí)間積分，可以允許較大的時(shí)間步長，在保證積分穩(wěn)定性的情況下縮短了計(jì)算時(shí)間。為了滿足模擬實(shí)際天氣的需要，該模型有一套物理過程及參數(shù)化過程，包括：云微物理過程，積云對流參數(shù)化，長波輻射，短波輻射，邊界層湍流，近地面層，陸面參數(shù)化，以及次網(wǎng)格湍流擴(kuò)散。該模型應(yīng)用了繼承式軟件設(shè)計(jì)、多級并行分解算法、選擇式軟件管理工具、中間軟件包結(jié)構(gòu)。目前WRF模型水平網(wǎng)格精度現(xiàn)已能精確到1 km甚至更高；使得WRF模型成為改進(jìn)從云尺度到天氣尺度等不同尺度重要天氣特征預(yù)報(bào)精度的工具。同時(shí)WRF模型擁有先進(jìn)的資料3維變分同化技術(shù)（3DVAR），能充分有效地將各種資料信息同化到模型初始場中，為模型提供了更高質(zhì)量的初值，達(dá)到了明顯改進(jìn)數(shù)值預(yù)報(bào)質(zhì)量的目的。

2.3.2 模型配置

圖 4 WRF模擬區(qū)域設(shè)置（D1：30 km；D2：10 km）Fig. 4 The simulation field set by WRF (D1:30 km; D2:10 km)

表 3 WRF模型設(shè)置Tab. 3 Setting of WRF model

考慮到天氣系統(tǒng)具有遷移性，所以實(shí)際的模擬區(qū)域遠(yuǎn)大于渤海與北黃海區(qū)域。綜合考慮區(qū)域與分辨率，模擬采用了雙重嵌套技術(shù)，分辨率較粗的區(qū)域（30 km）能夠抓住整個(gè)影響東中國近海的天氣系統(tǒng)，而分辨率較高的區(qū)域（10 km）針對渤海與北黃海區(qū)域的天氣系統(tǒng)進(jìn)行高分辨率的模擬再現(xiàn)。區(qū)域設(shè)置見圖4。

為了與水平分辨率相匹配，模型的垂直分辨率做了很大程度的調(diào)整，一共設(shè)計(jì)了 44個(gè) 層（一般的模擬應(yīng)用為27層）。由于WRF模型本身提供了眾多的物理選項(xiàng)，針對研究海域確定了合適的方案，具體見表 3。對于海面風(fēng)來說，最關(guān)鍵的是大氣邊界層湍流方案，本模型采用了YSU（韓國延世大學(xué)）開發(fā)的邊界層方案。經(jīng)過大量的數(shù)值試驗(yàn)證實(shí)，其非常適合海面的數(shù)值模擬。積云參數(shù)化方案主要與模型分辨率有關(guān)，針對 10 km分辨率，Kain-Fritsch方案是比較合適的。運(yùn)行WRF模式，為保證模式積分穩(wěn)定性，取時(shí)間步長為120 s，預(yù)報(bào)未來48 h，運(yùn)行完畢需要2～3 h左右；將計(jì)算出來的預(yù)報(bào)風(fēng)場用作流場的驅(qū)動(dòng)場預(yù)報(bào)未來 48 h的流場和波浪場，大約需2～3 h，一旦有溢油發(fā)生，可以馬上投入預(yù)報(bào)預(yù)警。此外，利用該系統(tǒng)建立6個(gè)月的模擬結(jié)果數(shù)據(jù)庫亦可實(shí)現(xiàn)對溢油的溯源回推。

2.3.3 模型檢驗(yàn) 利用環(huán)渤海與北黃海的岸邊氣象地面觀測站、海上船舶及島嶼觀測的風(fēng)速數(shù)據(jù)，我們對預(yù)報(bào)的風(fēng)場進(jìn)行了統(tǒng)計(jì)分析，時(shí)間為 3月10日－5月10日，共2個(gè)月。統(tǒng)計(jì)結(jié)果見圖5，結(jié)果表明，風(fēng)速均方根誤差在2.2～3.2 m/s之間，風(fēng)向平均誤差在24～36°之間，24 h之內(nèi)的風(fēng)速均方根誤差小于2.5 m/s，風(fēng)向平均誤差小于 28°。預(yù)報(bào)結(jié)果較為理想。風(fēng)場的預(yù)報(bào)結(jié)果見圖6。

與過去求解對流擴(kuò)散方程的算法相比，應(yīng)用油粒子方法不但較好地解決了油膜在環(huán)境動(dòng)力作用下的變形和破碎過程，并能準(zhǔn)確地預(yù)報(bào)油膜邊緣的擴(kuò)展過程和油膜形狀在風(fēng)向上的明顯拉伸現(xiàn)象，可以切和實(shí)際地模擬油膜的不規(guī)則形狀和漂移軌跡，而且可以有效地消除數(shù)值發(fā)散問題。溢油在水表面的漂移和風(fēng)化過程可通過表面擴(kuò)散、平移、輸送、乳化和蒸發(fā)的算法來模擬計(jì)算。

油粒子方法實(shí)際上就是追蹤油粒子微團(tuán)隨周圍水體平流和湍流擴(kuò)散運(yùn)動(dòng)的過程，因?yàn)椴捎昧死窭嗜辗椒ㄟ@種確定性方法模擬平流過程，并采用了隨機(jī)走動(dòng)法這種隨機(jī)性方法模擬擴(kuò)散過程，所以也可以說該方法是確定性方法和隨機(jī)性方法的結(jié)合。

圖 5 風(fēng)速與風(fēng)向的統(tǒng)計(jì)結(jié)果Fig. 5 Statistical result of wind speed and direction

圖 6 2010年4月14日海面風(fēng)場預(yù)報(bào)結(jié)果示意（色標(biāo)表示不同風(fēng)速等級，單位：m/s）Fig. 6 Sketch of the sea surface wind field prediction in Apr. 14th 2010 (The colorbar represents for the different wind speed scales, Unit: m/s)

3 結(jié) 語

考慮到溢油預(yù)報(bào)預(yù)警需要高分辨率的海流、風(fēng)場以及波浪的數(shù)值預(yù)報(bào)和多種過程較為完善的溢油漂移擴(kuò)散數(shù)值預(yù)報(bào)模式，因此，本文主要介紹了渤海海域溢油應(yīng)急預(yù)測預(yù)警系統(tǒng)中溢油預(yù)測的動(dòng)力要素，主要包括海流數(shù)值預(yù)報(bào)模型、波浪數(shù)值預(yù)報(bào)模型、海面風(fēng)場數(shù)值預(yù)報(bào)模型和溢油漂移擴(kuò)散數(shù)值預(yù)報(bào)模型，通過對四個(gè)模型的研究和開發(fā)，為溢油應(yīng)急預(yù)測預(yù)警系統(tǒng)的研制提供了基本的要素。

[1] 王耀華. 海上溢油的歸宿和形態(tài)研究 [D]. 大連: 大連海事大學(xué),2010: 14-28.

[2] 趙冬至, 張存智, 徐恒振. 海洋溢油災(zāi)害應(yīng)急響應(yīng)技術(shù)研究 [M].北京: 海洋出版社, 2006: 184-185.

[3] 劉欽政, 張存智. 渤海溢油數(shù)值預(yù)報(bào)研究 [J]. 海洋預(yù)報(bào), 2005,22(增刊): 70-76.

[4] 牟林, 鄒和平, 武雙全, 等. 海上溢油數(shù)值模型研究進(jìn)展 [J]. 海洋通報(bào), 2011, 30(4): 473-480.

[5] 朱建榮. 海洋數(shù)值計(jì)算方法和數(shù)值模式 [M]. 北京: 海洋出版社,2003: 128-129.

[6] Mellor G L, Yamada T. Development of a turbulence closure model for geophysical fluid problem [J]. Reviews of Geophysics, 1982,24(4): 851-875.

Numerical model research on Emergency Warning & Predicting system of ocean oil spill: I. Research on predicting of ocean dynamical factors

MU Lin, WU Shuang-quan, SONG Jun, LI Huan, LIU Shou-hua, LI Yan, GAO Jia

（National Marine Data and Information Service，Tianjin 300171，China）

The “oil particle” model is adopted to simulate the drift and diffusion of the oil spill in ocean. The simulative results of the current, wave and wind field are taken as the outside force for the computation to drive the oil spill drift model, and the drift trace is consequently calculated. The Emergency Warning& Predicting system of Bohai Sea is established by combining those computed oil spill results and several platforms like GIS. As a part of this system,this paper introduces the research status of the numerical computing models of current, wave and wind field, and also the oil spill predicting model.

Oil Particle；current；wave；Oil Spill Weathering Model

P731.2; X834

A

1001-6932(2011)05-0502-07

2011-06-15；

2011-08-02

國家發(fā)展與改革委員會(huì)高技術(shù)產(chǎn)業(yè)發(fā)展項(xiàng)目( 2008-311-000-035 )；海洋公益性行業(yè)科研專項(xiàng)（200905001）。

牟林 （1977－），男，博士，副研究員。電子郵箱：moulin1977@hotmail.com。