苗琪,徐福敏,俞茂玲
( 1. 中國電建集團(tuán) 華東勘測設(shè)計研究院有限公司,浙江 杭州 311122;2. 河海大學(xué) 港口海岸與近海工程學(xué)院,江蘇 南京210098)
20 世紀(jì)50年代以來,作為北冰洋一部分的波弗特海(Beaufort Sea)近海大陸架被發(fā)現(xiàn)具有豐富的石油、天然氣儲量。隨著全球變暖,北極海冰儲量大量減少[1–3],北極航道潛力巨大。北極地區(qū)開采所得的礦產(chǎn)利用北極航道經(jīng)由波弗特海運(yùn)往世界各地,較傳統(tǒng)方式可節(jié)約近一半時間,因此近年來隨著氣候變暖和亞洲天然氣需求激增,各國均對該海域的極地研究展現(xiàn)出極大關(guān)注。對北極風(fēng)暴下波弗特海?馬更些河(Mackenzie River)河口水域的海浪特性進(jìn)行分析研究,是油氣開采以及北極航道開發(fā)中的關(guān)鍵一環(huán)。
風(fēng)暴下短期內(nèi)海冰密集度、海冰厚度、冰緣線位置均劇烈變化的動態(tài)海冰,往往會給北極地區(qū)海浪的精確模擬帶來巨大困難。國外早期針對波弗特海海域的波浪模擬及海浪特性研究中[4–7],在海冰的處理上多采用的是未經(jīng)可靠驗(yàn)證的模型或經(jīng)驗(yàn)公式,或著眼于夏季無海冰覆蓋的近岸區(qū)域,對海冰?海浪間的相互作用考慮較為粗糙,且關(guān)于北極風(fēng)暴過程中海冰的短期運(yùn)動變化、不同風(fēng)暴下的海浪分布特征的研究依然較少,因此對海冰影響下波弗特海至馬更些河河口海域的海浪進(jìn)行精確數(shù)值模擬,以及對如何在海浪數(shù)值模擬中正確計量海冰的影響開展研究具有一定的研究意義。
海浪模式是現(xiàn)今較為成熟的海浪計算及預(yù)報手段。第三代海浪模式WAVEWATCH Ⅲ在控制方程、程序結(jié)構(gòu)、數(shù)值方法等很多方面都作出改進(jìn),因此近年來在海浪業(yè)務(wù)預(yù)報上得到了廣泛的應(yīng)用[8–11]。在版本V5.16 中,WAVEWATCH Ⅲ模式已具備多種海冰對海浪的能量耗散源項(xiàng)(IC1?IC4,IS1?IS2),能在極地地區(qū)海浪數(shù)值模擬中相對精確地考慮海冰的影響,因此近年來,一些學(xué)者開始利用WAVEWATCH Ⅲ模式開展波弗特海海域的海冰?海浪相互作用研究,如Erick 等[12]利用WAVEWATCH Ⅲ模式中的IC3源項(xiàng)對秋季波弗特海的波浪進(jìn)行后報,并利用后報結(jié)果和浮標(biāo)實(shí)測數(shù)據(jù)分析了不同冰種對海冰的能量耗散能力;Cheng 等[13]利用實(shí)測冰情、浪情率定了WAVEWATCH Ⅲ模式中IC3 源項(xiàng)的參數(shù),并對波弗特海海域的海浪進(jìn)行模擬。這些研究多針對特定海冰源項(xiàng),而對于WAVEWATCH Ⅲ模式中的多種海冰源項(xiàng)缺乏詳細(xì)而完整的評述。據(jù)此,本文利用WAVEWATCH Ⅲ模式中的不同海冰源項(xiàng)分別對海冰影響下波弗特海夏季開冰期的暴風(fēng)浪進(jìn)行數(shù)值模擬,以期為極地地區(qū)海洋工程建設(shè)、海洋防災(zāi)減災(zāi)等提供一定的參考。
WAVEWATCH Ⅲ模式中,球坐標(biāo)系下關(guān)于動譜密度的控制方程通常表示為
式中,N為動譜密度;t為時間; ?x為水平拉普拉斯算子;k為波數(shù)矢量;θ為波向;S為源函數(shù);σ為圓頻率;Cg為群速度;U為環(huán)境流速;d為水深;s為與波向相同的坐標(biāo);m為與s垂直的另一坐標(biāo)。
WAVEWATCH Ⅲ模式中海冰損耗作用最傳統(tǒng)的處理方式是采用經(jīng)驗(yàn)方法[14],近似地將海冰密集度大于某臨界值(由用戶指定)的網(wǎng)格點(diǎn)視作陸地,小于某臨界值的視為不受海冰影響的開闊水域,兩者之間進(jìn)行線性插值,并通過修改控制方程使波能在每個海冰網(wǎng)格點(diǎn)處進(jìn)行百分比傳遞,該方法記為IC0。WAVEWATCHⅢ V5.16 版本中引入了4 種海冰損耗源項(xiàng),分別以IC1?IC4 來表示:第一種海冰損耗源項(xiàng)(IC1)對所有的波浪成分都采用相同的耗散率,是一種不依賴于頻率的簡單損耗過程,在計算時需要輸入海冰的密集度以及損耗率;第二種海冰損耗源項(xiàng)(IC2)將海冰視作一個連續(xù)的彈性薄盤,能量損耗由冰下邊界層的摩擦產(chǎn)生[15–17],在計算時需要輸入海冰密集度、海冰有效剪切模量、黏度;第三種海冰損耗源項(xiàng)(IC3)認(rèn)為海冰是一個黏滯彈性層[18],對損耗的計算需要輸入海冰的密集度、厚度、有效黏度系數(shù)、密度和有效剪切模量;第四種海冰損耗源項(xiàng)(IC4)對能量損耗給出了幾組簡單的依賴于頻率的經(jīng)驗(yàn)或參數(shù)化的公式[19–24],從而使得源項(xiàng)能夠簡單靈活卻有效地再現(xiàn)冰?浪相互作用過程中海浪能量的衰減。
散射源項(xiàng)有兩種,分別以IS1 和IS2 表示:第一種散射源項(xiàng)(IS1)是比較簡單的海冰對海浪的守恒作用,僅考慮海冰密集度,海冰只散射海浪能量,但是不耗散海浪的能量;第二種散射源項(xiàng)(IS2)依賴于浮冰尺寸,并且考慮了海浪作用下海冰的破碎從而能夠不斷更新最大浮冰尺寸[24]。
2014年7月 末,陸 緣 海 冰 帶(Marginal Ice Zone,MIZ)項(xiàng)目由美國海軍在波弗特海/楚科奇海開展(http://www.apl.washington.edu/project/project.php?id=miz),范圍為72°~79°N,137°~170°W。項(xiàng)目過程中共投放3 個配備定位系統(tǒng)的SWIFT(Surface Wave Instrument Floats with Tracking)漂移浮標(biāo),并于9月底回收。SWIFT10 和SWIFT11 浮標(biāo)于7月27 日在阿拉斯加北部約100km 的近岸處從RV“Ukpik”號極地考察船上投放,向西北方向漂移。在9月1 日之前兩者漂移路徑大致相同,然后SWIFT10 浮標(biāo)開始進(jìn)入高密集度海冰區(qū)域,直至9月15 日。SWIFT15 浮標(biāo)在8月5?17 日期間由R/V “Araon”號極地考察船投放并回收,且全程被海冰環(huán)繞。各浮標(biāo)漂移路徑見圖1。
由于SWIFT15 浮標(biāo)全程被高密集度海冰包圍,幾乎處于擱淺狀態(tài),大部分時間測得有效波高均為0m,最大有效波高不足0.7m,故本次研究僅選取SWIFT10和SWIFT11 兩個浮標(biāo)用作后續(xù)驗(yàn)證。
圖1 SWIFT 浮標(biāo)漂移路徑Fig. 1 The track of SWIFT drifting buoys
2014年8?9月有多場風(fēng)暴襲擊波弗特海至馬更些河河口水域,且影響范圍較廣。研究區(qū)域過大,通常會導(dǎo)致計算效率過于低下,故研究采用自嵌套的方案,即先用分辨率較低的網(wǎng)格對較大的海域進(jìn)行計算,然后將計算所得波浪譜作為邊界輸入到內(nèi)層分辨率更高、范圍更小的模型中進(jìn)行計算,以提高計算效率。嵌套外層包括整個東西伯利亞海、楚科奇海、波弗特海和部分位于西北太平洋的白令海,內(nèi)層為波弗特海至馬更些河口外水域,據(jù)此建立自波弗特海至馬更些河口的WAVEWATCH Ⅲ兩級嵌套模型,由外層的計算獲取傳播至波弗特海的所有低頻涌浪成份,輸入到第二層中,然后在存在大量海冰的波弗特海進(jìn)行風(fēng)?冰?浪相互作用模擬研究,并探究WAVEWATCH Ⅲ對海冰影響下的海浪場模擬的準(zhǔn)確性。
嵌套計算時,中低頻截斷頻率取為0.04118 Hz,高頻截斷頻率取為0.7904 Hz,以對數(shù)分布fn+1=1.1fn劃分為32 個頻率網(wǎng)格,外層譜方向劃分為30 個波向,分辨率為12°,內(nèi)層譜方向劃分為36 個波向,分辨率為10°。
3.2.1 計算區(qū)域和水深
本文建立的WAVEWATCH Ⅲ兩級自嵌套海浪模型網(wǎng)格采用球面坐標(biāo),各層范圍為:(1)嵌套Ⅰ區(qū):整個東西伯利亞海、波弗特海和部分位于西北太平洋的白令海,范圍為58°~78°N,140°E~110°W,空間分辨率為0.25°×0.25°,網(wǎng)格節(jié)點(diǎn)數(shù)為441×81;(2)嵌套Ⅱ區(qū):波弗特海至馬更些河河口外水域,范圍為65°~75°N,125°~175°W,空間分辨率為5′×5′,網(wǎng)格節(jié)點(diǎn)數(shù)為601×121(圖2)。
圖2 兩級嵌套模型研究區(qū)域示意圖Fig. 2 Two-level nested computational domain
模型水深數(shù)據(jù)采用美國國家海洋與大氣管理局(NOAA)的ETOPO1 全球地形數(shù)據(jù)集,其空間分辨率為1′×1′,在各級嵌套區(qū)域分別插值至各網(wǎng)格點(diǎn)處。
3.2.2 驅(qū)動風(fēng)場
模型風(fēng)場數(shù)據(jù)采用CCMP(Cross-Calibrated Multi-Platform)V2.0 風(fēng)場,它是CCMP 風(fēng)場的升級版本,在繼承CCMP 風(fēng)場高時間、空間分辨率的基礎(chǔ)上(時間分辨率為6 h,空間分辨率為0.25°×0.25°),又將數(shù)據(jù)覆蓋時間從2011年12月延長至今。對2014年8月1 日00 時 至9月30 日18 時 的CCMP V2.0 風(fēng) 場 進(jìn) 行空間線性插值,作為模型驅(qū)動風(fēng)場。
鑒于CCMP V2.0 風(fēng)場在相關(guān)研究中應(yīng)用暫時較少,將CCMP V2.0 風(fēng)場的風(fēng)速與SWIFT 漂移浮標(biāo)的實(shí)測風(fēng)速進(jìn)行對比,以驗(yàn)證其可靠性。將CCMP 風(fēng)速分別插值至每個時刻兩個浮標(biāo)所在位置處,對比結(jié)果如圖3 所示??梢钥闯?,兩條風(fēng)速曲線趨勢和變化幾乎完全一致,數(shù)據(jù)吻合較好,證明研究所用的CCMP V2.0 風(fēng)場風(fēng)速與波弗特海的風(fēng)速吻合較好,能夠較好地反映北極風(fēng)暴過程。
3.2.3 海冰資料
海冰密集度數(shù)據(jù)資料來源有兩種。第一種為NOAA 最優(yōu)插值海面溫度分析數(shù)據(jù)集(NOAA Optimum Interpolation 1/4 Degree Daily Sea Surface Temperature Analysis, Version 2),被用于嵌套模型的外層,即從白令海峽至波弗特海的整個海域。該數(shù)據(jù)集由NOAA/NCDC 提供,空間上覆蓋全球,分辨率為0.25°×0.25°,時間上從1981年9月1 日至今,分辨率為1 d。第二種為MASAM2 逐日4km 北極海冰密集度數(shù)據(jù)集(MASAM2: Daily 4km Arctic Sea Ice Concentration,Version 1),被用于模型的第二層,即波弗特海至馬更些河河口。該數(shù)據(jù)集時間覆蓋范圍為2012年7月至今,分辨率為1 d,空間覆蓋范圍為29.08°~90°N,環(huán)全球經(jīng)度,空間分辨率為4km×4km。
海冰厚度數(shù)據(jù)采用美國國家環(huán)境預(yù)報中心(NCEP)氣候預(yù)報系統(tǒng)2.0 版本(CFSV2)的逐時時間序列產(chǎn)品(NCEP Climate Forecast System Version 2(CFSV2)Selected Hourly Time-Series Products),范圍覆蓋全球,有0.2°×0.2°、0.5°×0.5°、1.0°×1.0°和2.5°×2.5°的空間分辨率可供選擇,時間覆蓋為2011年1月1 日至今,分辨率為6 h。CFSV2 海冰厚度數(shù)據(jù)相比實(shí)測數(shù)據(jù)往往偏大[25–26],故在本次模擬中對其乘以0.6 的折減。
3.2.4 參數(shù)設(shè)置
計算時,對于海冰的散射作用,由于計算海域中浮冰直徑相對波長來說較?。蛇_(dá)2m,但通常小于1m),能帶來的散射作用有限,且WAVEWATCH Ⅲ模式中的海冰散射模塊尚不能很好地體現(xiàn)海冰對海浪能量的衰減作用,故本次計算不考慮海冰的散射作用。對于海冰對海浪能量的衰減作用,分別采用IC0、IC1、IC2、IC34 種機(jī)理進(jìn)行考慮,其中IC0 額外設(shè)置海冰損耗上下限為50%~50%及25%~75%兩種方案。其中50%~50%方案為模型的缺省設(shè)置方案,而25%~75%方案為一些學(xué)者采用過的改進(jìn)方案。
圖3 CCMP V2.0 風(fēng)場風(fēng)速驗(yàn)證Fig. 3 Validation of CCMP V2.0 wind speed
由于應(yīng)用于該大范圍海域時,缺乏海冰的剪切模量、黏滯系數(shù)等實(shí)測屬性參數(shù),故除海冰密集度、厚度由衛(wèi)星資料提供外,其余海冰相關(guān)的模式參數(shù)均取為模型缺省值。其他源項(xiàng)設(shè)置還包括風(fēng)能輸入和白浪耗散項(xiàng)選取Ardhuin 等[27]的WAM4 方案,四波非線性相互作用項(xiàng)選取DIA 方法[28],三波非線性相互作用項(xiàng)選取LTA 方法[29],底摩擦項(xiàng)選取JONSWAP 底摩擦公式[30],水深變淺引起的波浪破碎項(xiàng)選取Battjes-Janssen 的方法[31],參數(shù)均為默認(rèn)值。
由于SWIFT 浮標(biāo)屬于漂移浮標(biāo),所處位置隨時間在不斷變化,需利用WAVEWATCHⅢ模式中的沿軌跡輸出功能輸出兩個SWIFT 浮標(biāo)的軌跡在對應(yīng)時刻的波浪要素以進(jìn)行對比驗(yàn)證。WAVEWATCH Ⅲ模式的沿軌跡輸出功能僅能輸出一維海浪譜,利用公式編程由離散的海浪譜換算有效波高(其中m0為海浪譜的零階譜矩),從而與實(shí)測數(shù)據(jù)進(jìn)行驗(yàn)證比較。圖4 為不同海冰方案下兩個SWIFT 浮標(biāo)處的有效波高模擬結(jié)果及浮標(biāo)實(shí)測結(jié)果。
由圖4 可見,在9月1 日之前,由于SWIFT10 浮標(biāo)和SWIFT11 浮標(biāo)均處于開闊水域中,海冰對兩浮標(biāo)唯一的影響為當(dāng)風(fēng)向背離海冰區(qū)時風(fēng)區(qū)長度的不同,但該影響較為微弱,故此時不同的海冰耗散機(jī)理對浮標(biāo)所在位置的海浪模擬并無明顯影響,5 種機(jī)理結(jié)果并無明顯區(qū)別,部分方案結(jié)果幾乎完全一致,因此部分方案的結(jié)果被遮擋,但各方案結(jié)果均較為理想。
9月1?14 日,SWIFT10 浮標(biāo)實(shí)測有效波高開始迅速衰減,從浮標(biāo)掛載的攝像頭錄像可以發(fā)現(xiàn),此時SWIFT10 浮標(biāo)已被海冰所包圍(圖5),海面上大面積覆蓋的海冰對波浪能量產(chǎn)生劇烈的衰減作用。這段時間內(nèi),IC0(25%)、IC2、IC3 有效波高模擬結(jié)果的趨勢與實(shí)測數(shù)據(jù)相同,但數(shù)值明顯偏高,IC1 趨勢與實(shí)測相同,且數(shù)值相對較為接近。IC2 模擬結(jié)果與IC1 非常接近,幾乎被完全遮擋(從下文誤差分析中也可以看出兩方案差距較小,SWIFT10 浮標(biāo)處IC2 誤差略微高于IC1,而在SWIFT11 浮標(biāo)處略低于IC1)。IC0(50%)在該時段內(nèi)波要素始終為0。
9月14 日以后,SWIFT10 浮標(biāo)重新回到開闊水域,IC0(50%)方案仍圍繞觀測值出現(xiàn)劇烈上下波動,而其余方案均較為理想。與SWIFT10 浮標(biāo)不同,SWIFT11 浮標(biāo)全程位于冰密集度較低的開闊海域內(nèi),從圖4b 可以看出該浮標(biāo)的有效波高模擬結(jié)果趨勢與9月14 日后的SWIFT10 浮標(biāo)結(jié)果類似,即IC0(50%)方案出現(xiàn)異常波動,而其余方案均擬合較好。
IC0(50%)的參數(shù)設(shè)置,根據(jù)其原理,是一種不連續(xù)的處理方式,會對海冰密集度大于50%區(qū)域的海浪能量產(chǎn)生劇烈衰減,而對海冰密集度小于50%的區(qū)域衰減較微弱。根據(jù)衛(wèi)星及遙感數(shù)據(jù),浮標(biāo)所在的波弗特海域在計算時間段內(nèi)海冰密集度多為40%~70%,這就導(dǎo)致在該密集度范圍下,IC0(50%)方案對密集度變化過于敏感,結(jié)果極易產(chǎn)生突變及波動。在9月1?14 日,海冰密集度大于50%,故該時段內(nèi)IC0(50%)方案模擬結(jié)果波高與周期均為0。
圖4 不同海冰方案下兩個SWIFT 浮標(biāo)有效波高對比Fig. 4 Comparison of simulated significant wave height by two SWIFT bouys of different ice source terms
IC0(25%)方案雖然也是一種不連續(xù)的處理方式,會對海冰密集度大于75%區(qū)域的海浪能量產(chǎn)生劇烈衰減,而對海冰密集度小于25%的區(qū)域衰減較微弱,但其能在兩密集度范圍之間均勻過渡。因而雖然機(jī)理相同,但并未出現(xiàn)類似IC0(50%)方案的波高突變。
為更直觀地定量比較不同機(jī)理對有效波高的模擬水平,確定該海域最佳的海冰機(jī)理,首先分別對不同機(jī)理的模擬結(jié)果和實(shí)測數(shù)據(jù)進(jìn)行顯著性水平檢驗(yàn),假定模擬結(jié)果與實(shí)測數(shù)據(jù)之間相關(guān)性為0,計算得p值均小于0.05,證明存在相關(guān)性。其次,選取均方根誤差(Root Mean Square Error,RMSE)、相關(guān)系數(shù)(Correlation Coefficient,CC)和絕對平均誤差(Mean Absolute Error,MAE)3 個參數(shù)來統(tǒng)計分析實(shí)測值和計算值之間的誤差,分別定義為
圖5 9月1?14 日環(huán)繞SWIFT10 的海冰Fig. 5 Sea ice pictures taken around SWIFT10 from September 1st to 14th
表1 和表2 分別為SWIFT10 和SWIFT11 浮標(biāo)有效波高模擬值與浮標(biāo)值的相關(guān)系數(shù)、均方根誤差及絕對平均誤差。
由表1 和表2 可以看出,各方案下SWIFT11 浮標(biāo)處的有效波高與實(shí)測數(shù)據(jù)之間的均方根誤差和絕對平均誤差普遍小于SWIFT10,而相關(guān)系數(shù)高于SWIFT10,證明WAVEWATCHⅢ模式對有效波高的模擬效果在遠(yuǎn)離冰緣、水域相對開闊時比在被海冰控制下更好。此外,從SWIFT10 浮標(biāo)處的有效波高結(jié)果可以看出,在離海冰較近,受其影響較大的水域,IC1 方案相較其他方案而言,模擬結(jié)果誤差更小,相關(guān)系數(shù)更高。
表1 SWIFT10 有效波高模擬值與浮標(biāo)值的均方根誤差(RMSE)、相關(guān)系數(shù)(CC)和絕對平均誤差(MAE)Table 1 RMSE, CC and MAE of observed and simulated significant wave heights in SWIFT10
表2 SWIFT11 有效波高模擬值與浮標(biāo)值的的均方根誤差(RMSE)、相關(guān)系數(shù)(CC)和絕對平均誤差(MAE)Table 2 RMSE, CC and MAE of observed and simulated significant wave heights in SWIFT11
盡管已有相關(guān)研究[25]顯示,在有較為完備的海冰屬性參數(shù)時,IC3 源項(xiàng)表現(xiàn)出的對海浪的能量衰減特征更加符合實(shí)際情況,但在本次研究中,由于研究區(qū)域范圍較大,并無能力獲得除海冰密集度、厚度以外的其他海冰數(shù)據(jù),如剪切模量、黏性系數(shù)等,而這些時空均變化的海冰數(shù)據(jù)恰恰是IC3 源項(xiàng)必不可少的輸入?yún)⒘?,均取為缺省值顯然不合理。相反,IC1 源項(xiàng)對海冰數(shù)據(jù)要求較低,因此在應(yīng)用于大范圍水域、資料匱乏的前提下,計算結(jié)果反而與實(shí)測數(shù)據(jù)更為接近。
本文所采用參數(shù)的模擬結(jié)果并非最優(yōu)解,其他海冰?波浪源項(xiàng)在海冰數(shù)據(jù)充足、參數(shù)調(diào)試合理的情況下可能會得到更優(yōu)的結(jié)果。但是在以下限定條件下:(1)應(yīng)用于北極波佛特海這一特定海域的有效波高模擬,(2)缺乏海冰剪切模量、黏性系數(shù)等實(shí)測參數(shù),(3)模型參數(shù)均采用默認(rèn)值,使用IC1 海冰源項(xiàng)來計量海冰對海浪的衰減作用最為合理,能夠?yàn)楹罄m(xù)淺水地區(qū)的進(jìn)一步海浪嵌套模擬提供更加精準(zhǔn)的邊界條件。
建立典型北極風(fēng)暴下的波弗特海多重嵌套WAVEWATCH Ⅲ海浪模型,利用WAVEWATCH Ⅲ模式中的不同海冰損耗源項(xiàng)設(shè)置不同海冰方案,并利用浮標(biāo)實(shí)測數(shù)據(jù)對不同海冰源項(xiàng)在秋季波弗特海的海浪模擬效果進(jìn)行對比研究。結(jié)論如下:
(1)在受海冰影響較小的開闊水域,各海冰機(jī)理的模擬結(jié)果區(qū)別較小,且均驗(yàn)證良好;
(2)當(dāng)海冰實(shí)際存在時,若在海浪模擬中不考慮海冰的能量衰減作用,將海冰密集度較高的地方直接視為陸地,則當(dāng)海冰密集度大于某一閾值時,有效波高模擬結(jié)果會出現(xiàn)大量的向0 突變,結(jié)果不合理。因此在該種情況下,使用計入海冰對海浪的能量衰減作用的IC1、IC2、IC3 源項(xiàng)更為合理;
(3)在缺乏大范圍海冰實(shí)測數(shù)據(jù)、模型其他參數(shù)采用默認(rèn)值的前提下對波弗特海至馬更些河河口這一特定水域進(jìn)行海浪模擬時,IC1 海冰源項(xiàng)能夠更為準(zhǔn)確地表達(dá)該海域特定的冰情、冰況對海浪的能量衰減特征,能夠?yàn)楹罄m(xù)淺水地區(qū)的進(jìn)一步海浪嵌套模擬提供更加精準(zhǔn)的邊界條件。