羅志發(fā)，黃本勝，譚 超，邱 靜，黃廣靈

(1.廣東省水利水電科學(xué)研究院，廣東 廣州 510635；2.廣東省水動(dòng)力學(xué)應(yīng)用研究重點(diǎn)實(shí)驗(yàn)室，廣東 廣州 510635；3.廣東省流域水環(huán)境治理與水生態(tài)修復(fù)重點(diǎn)實(shí)驗(yàn)室，廣東 廣州 510635)

1 概述

珠江河口位于中國大陸南端，瀕臨南海，海岸線狹長，沿岸地勢較低。特殊的地理位置及氣候條件使得珠江河口成為國內(nèi)臺(tái)風(fēng)登陸最頻繁、風(fēng)暴潮災(zāi)害最嚴(yán)重的區(qū)域之一。近年來，在全球氣候變化背景下，風(fēng)暴潮災(zāi)害次數(shù)和強(qiáng)度呈增加的趨勢，登陸珠江河口的超強(qiáng)臺(tái)風(fēng)引發(fā)的風(fēng)暴潮災(zāi)害，屢次刷新沿岸潮位站歷史最高潮位和最大增水記錄。珠江河口沿岸人口稠密、經(jīng)濟(jì)發(fā)達(dá)，每年由于風(fēng)暴潮造成的經(jīng)濟(jì)損失相當(dāng)巨大，已經(jīng)成為影響人民生活質(zhì)量、制約國民經(jīng)濟(jì)高質(zhì)量發(fā)展的重要因素。

為了防御風(fēng)暴潮災(zāi)害，珠江河口沿海陸續(xù)興建了江海堤圍。盡管如此，強(qiáng)臺(tái)風(fēng)往往伴隨著強(qiáng)臺(tái)風(fēng)浪，在強(qiáng)浪沖擊作用下，江海堤圍面臨著堤身結(jié)構(gòu)損毀的風(fēng)險(xiǎn)，可能導(dǎo)致海水漫過海堤入侵內(nèi)陸，引發(fā)城市內(nèi)澇。因此，準(zhǔn)確模擬風(fēng)暴潮及臺(tái)風(fēng)浪過程，了解正面登陸臺(tái)風(fēng)影響下珠江河口近岸海域風(fēng)暴潮增水及臺(tái)風(fēng)浪變化過程，對準(zhǔn)確評估風(fēng)暴潮的災(zāi)害風(fēng)險(xiǎn)以及做好風(fēng)暴潮災(zāi)害防御與應(yīng)對工作、減少風(fēng)暴潮災(zāi)害損失有著重要現(xiàn)實(shí)的意義。

珠江河口是多種動(dòng)力因子協(xié)同作用的復(fù)雜系統(tǒng)，呈“三江匯流，八口入海”的形勢。河網(wǎng)區(qū)水網(wǎng)密布，橫向支汊發(fā)育，其間多種動(dòng)力因子相互耦合，動(dòng)力復(fù)雜多變。珠江河口風(fēng)暴潮、波浪數(shù)值模擬已開展較多的研究，如二維風(fēng)暴潮模型的建立[1]，風(fēng)浪數(shù)值模擬[2-3]，路徑、風(fēng)速對風(fēng)暴潮增水的影響[4]，地形對局部增水的影響[5]等。波浪在近岸傳播過程中受到淺水變形的影響發(fā)生破碎，產(chǎn)生的輻射應(yīng)力驅(qū)動(dòng)水面抬升形成增水，在臺(tái)風(fēng)期間波浪增水相對較大[6]。波浪與風(fēng)暴潮耦合計(jì)算已開展了相關(guān)的研究，在不同的海域，由于岸線地形及水動(dòng)力條件的不同，波浪對風(fēng)暴潮增水的影響有所差異[7-9]。針對珠江河口復(fù)雜的動(dòng)力系統(tǒng)，本文基于無結(jié)構(gòu)三角網(wǎng)格的數(shù)值模式，優(yōu)化擬合復(fù)雜岸形，以珠江河網(wǎng)—河口作為整體，綜合考慮風(fēng)、徑流、潮流、波浪等動(dòng)力因子的耦合作用，構(gòu)建珠江河口波浪—風(fēng)暴潮耦合數(shù)值模型，可有效提高風(fēng)暴潮數(shù)值模擬的準(zhǔn)確度，為進(jìn)一步探討風(fēng)暴潮—波浪耦合作用提高有效的科學(xué)手段。以1822號“山竹”臺(tái)風(fēng)為例，對臺(tái)風(fēng)過程的風(fēng)暴潮增水及波浪進(jìn)行數(shù)值模擬，定量分析考慮波浪—風(fēng)暴潮耦合作用對風(fēng)暴潮數(shù)值模擬結(jié)果的影響。

2 模型的建立與驗(yàn)證

2.1 模型簡介

本文選用海洋環(huán)流模式SELFE建立南海—珠江河口雙重嵌套模型(模型網(wǎng)格見圖1)，該模式基于非結(jié)構(gòu)化網(wǎng)格，可精細(xì)化擬合復(fù)雜岸線和地形，采用半隱式的歐拉—拉格朗日有限元算法求解N-S方程組，可有效提高模型計(jì)算效率。南海模型計(jì)算范圍為98°E～126°E，0°N～30°N，涵蓋整個(gè)南海及西北太平洋海域，南邊界至卡里馬塔海峽、北邊界至浙江省沿岸海域、東邊界至48小時(shí)警戒線。網(wǎng)格分辨率從近岸的 1 km 逐漸過渡到外海20 km，水深數(shù)據(jù)采用ETOP01全球1′×1′分辨率的地形資料。珠江河口模型上邊界為西江高要，北江石角，東江博羅，流溪河及潭江上游，外海下邊界取約100 m等深線處。模型采用無結(jié)構(gòu)網(wǎng)格，擬合復(fù)雜河岸邊界，對局部進(jìn)行加密，提高網(wǎng)格分辨率。模型網(wǎng)格共有101 752個(gè)節(jié)點(diǎn)，173 045個(gè)網(wǎng)格單元，網(wǎng)格大小為從河網(wǎng)區(qū)10 m逐漸過渡到外海的20 km；在垂向上采用Sigma坐標(biāo)，均勻分為10層；珠江三角洲網(wǎng)河區(qū)采用2005—2008年的大范圍實(shí)測地形，河口區(qū)及近岸海區(qū)采用2000—2008年海圖地形。南海模型由風(fēng)場、氣壓場及8個(gè)主要分潮(M2、S2、N2、K2、K1、O1、P1、Q1)[10]進(jìn)行驅(qū)動(dòng)，計(jì)算得到余水位，將余水位以及風(fēng)場、氣壓場、8個(gè)分潮作為珠江口模型的驅(qū)動(dòng)條件，對珠江河口風(fēng)暴潮增水進(jìn)行數(shù)值模擬，模型詳細(xì)介紹見文獻(xiàn)[7]。

圖1 模型嵌套計(jì)算網(wǎng)格示意

在風(fēng)暴潮模型的基礎(chǔ)上耦合波浪計(jì)算模塊，建立波浪—風(fēng)暴潮耦合模型。波浪模塊運(yùn)用Wind Wave Model-III(WWM-III)模式[11]，該模式主要用于模擬大洋—近岸波浪傳播過程，模擬的物理過程包括：波浪折射作用，淺水變形作用，風(fēng)成浪作用，白浪的耗散作用，水深引起的破碎作用，海底摩擦作用等?？刂品匠虨槟芰孔V平衡方程：

(1)

式中N為波作用密度譜;Cgx和Cgy分別為x和y方向的波群傳播速度;u和v分別為x和y方向的水流速度;σ為相對波頻;θ為波向;Cσ和Cθ分別為σ和θ方向的波浪傳播速度;Stot為譜密度表示的源項(xiàng)，是各種物理現(xiàn)象的源函數(shù)的疊加形式。

風(fēng)暴潮模型為波浪模型提供水位、流速、流向數(shù)據(jù)計(jì)算變化水動(dòng)力條件下的波浪要素。波浪模型輸出有效波高、平均波向、有效波周期、波長和底部波軌流速，輸出結(jié)果作為邊界波浪條件輸入風(fēng)暴潮模型，計(jì)算波流共同作用下的底切應(yīng)力及垂向混合系數(shù)。

2.2 模型驗(yàn)證

文獻(xiàn)[12]以“山竹”臺(tái)風(fēng)為算例，對模型計(jì)算的風(fēng)速、風(fēng)向、天文潮位、風(fēng)暴潮位等進(jìn)行了率定驗(yàn)證，驗(yàn)證結(jié)果良好，本文在此基礎(chǔ)上對波浪計(jì)算結(jié)果進(jìn)行進(jìn)一步驗(yàn)證。珠江口波浪模型外海波浪開邊界條件通過南海模型提供，海表面臺(tái)風(fēng)場及氣壓場通過臺(tái)風(fēng)模型及氣壓模型給出[12]，底摩擦系數(shù)設(shè)置為0.067[13]，波浪破波系數(shù)為0.78[11]。波浪計(jì)算時(shí)間步長設(shè)置為1 800 s，每個(gè)時(shí)間步長計(jì)算結(jié)束后與風(fēng)暴潮模型進(jìn)行數(shù)據(jù)反饋。模擬時(shí)間為2018年9月10日—9月20日，前5 d用于計(jì)算的穩(wěn)定，后5 d對風(fēng)暴潮及波浪過程進(jìn)行模擬。

圖2為外海浮標(biāo)站(113.9947°E，21.4948°N)有效波高驗(yàn)證示意，由圖2可知，計(jì)算結(jié)果有效波高變化過程與實(shí)測數(shù)據(jù)有效波高變化趨勢一致，相關(guān)系數(shù)為0.862，平均絕對誤差為0.53 m，二者吻合較好。計(jì)算有效波高的最大值與實(shí)測結(jié)果分別為9.4 m、9.9 m，相對誤差僅5.3%，表明模型計(jì)算結(jié)果能較好地反映有效波高的極值情況。由圖3可知，考慮波浪作用之后，潮位數(shù)值模擬結(jié)果更接近觀測值。橫門站考慮波浪作用后，最高潮位的絕對誤差由0.24 m減小為0.06 m；南沙站的模擬精度同樣有所提高。由此可見，考慮波浪—風(fēng)暴潮耦合作用后可提高風(fēng)暴潮數(shù)值預(yù)報(bào)的精度。綜上所述，模型計(jì)算結(jié)果能夠較好地反映“山竹”臺(tái)風(fēng)期間波浪及風(fēng)暴潮增水的變化過程，可為進(jìn)一步分析波浪與風(fēng)暴潮耦合提供基礎(chǔ)。

圖2 外海浮標(biāo)站有效波高驗(yàn)證示意

圖3 橫門(a)、南沙(b)站考慮和不考慮波浪作用潮位模擬結(jié)果示意

3 結(jié)果與討論

本節(jié)利用上文已建立的波浪—風(fēng)暴潮耦合數(shù)值模型，模擬“山竹”臺(tái)風(fēng)風(fēng)暴潮增水及波浪過程，討論波浪—風(fēng)暴潮耦合對風(fēng)暴潮增水?dāng)?shù)值模擬的影響。

3.1 波浪及風(fēng)暴潮增水計(jì)算結(jié)果分析

“山竹”臺(tái)風(fēng)過程對珠江口近岸波浪要素分布影響顯著，2018年9月16日13：00，臺(tái)風(fēng)中心位于珠江口南側(cè)(見圖4a)，臺(tái)風(fēng)中心風(fēng)力較小，有效波高較低，約為6～7 m；臺(tái)風(fēng)中心右側(cè)風(fēng)力較大，有效波高超過9 m。此時(shí)，伶仃洋海域?yàn)闁|北偏東風(fēng)，珠江口有效波高相對較小，約為2～4 m。16日17:00，臺(tái)風(fēng)于廣東臺(tái)山海宴鎮(zhèn)登陸，珠江口盛行東南偏南風(fēng)，有效波高有顯著增大趨勢，主要高值區(qū)更偏向近岸，有效波高最大值介于4～5 m(見圖4b)。選取燈籠山、南沙、萬頃沙、橫門站位(站位位置見圖4a)分析臺(tái)風(fēng)期間有效波高變化過程。由圖5可知，臺(tái)風(fēng)過境期間，各站位有效波高呈現(xiàn)先增大后減小的過程。磨刀門往外延伸瀕臨南海，因此燈籠山站受波浪影響較大，最大有效波高達(dá)3.2 m。南沙、萬頃沙、橫門等站位均位于伶仃洋灣內(nèi)，外海波浪傳入伶仃洋經(jīng)底摩擦、折射、繞射等作用，波能消耗，有效波高均小于2 m。

圖5 各站位臺(tái)風(fēng)期間有效波高變化

臺(tái)風(fēng)登陸前4 h，伶仃洋海域?yàn)闁|北偏東風(fēng)，伶仃洋東岸在離岸風(fēng)的驅(qū)動(dòng)下，產(chǎn)生減水，深圳灣區(qū)域最大減水值約為-1.5 m。伶仃洋西岸由于水體的橫向堆積，產(chǎn)生0.5～1.0 m的增水(見圖4c)。磨刀門、黃茅海區(qū)域?yàn)槠憋L(fēng)，有利于水體離岸輸運(yùn)，形成-0.5～-1.0 m的減水。臺(tái)風(fēng)登陸后(見圖4d)，珠江口海域普遍為東南偏南風(fēng)，珠江河口東南向的開口方向有利于水體向岸堆積并沿河道向上游輸運(yùn)，此時(shí)珠江口海域普遍達(dá)到增水的最大值。伶仃洋河口灣頂增水值較大，約為3.1 m，是由于伶仃洋河口灣喇叭狀的形態(tài)有利于水體向?yàn)稠斁奂?。其余口門區(qū)域最大增水值普遍為2.6～3.1 m。

圖4 “山竹”臺(tái)風(fēng)登陸前后波浪場及風(fēng)暴潮增水分布(a、c：登陸前，b、d：登陸后)

3.2 波浪對風(fēng)暴潮影響分析

利用風(fēng)暴潮與波浪耦合計(jì)算的總水位減去未考慮波浪作用的風(fēng)暴潮水位可得到波浪增水，波浪增水包含了由波浪輻射應(yīng)力產(chǎn)生的增水和波流相互作用產(chǎn)生的增水。圖6為臺(tái)風(fēng)登陸前后波浪增水分布示意。

臺(tái)風(fēng)登陸前(見圖6a)波浪增水影響范圍較廣，波浪增水由外海向近岸逐漸增大，外海波浪增水約為10 cm，近岸海域波浪增水約為15～20 cm，與相關(guān)研究結(jié)果較為一致[9]。近岸海域中，伶仃洋、香港島附近海域波浪增水較大，黃茅海等海域波浪增水相對較小，其空間分布與有效波高空間分布大體一致。臺(tái)風(fēng)登錄后(見圖6b)，波浪增水集中在近岸海域，主要分布在伶仃洋、磨刀門及黃茅海等海域近岸區(qū)域，波浪增水最大值約20 cm。

圖6 “山竹”臺(tái)風(fēng)登陸前后波浪增水分布示意

選取泗盛圍、南沙、橫門、燈籠山等站位分析風(fēng)暴潮增水及波浪增水過程，由圖7可知，各站位風(fēng)暴潮增水呈單峰型，增水快速增大至最大值而后迅速下降，最大增水值均超過3 m。波浪增水同樣呈現(xiàn)先增大后減小的過程，最大波浪增水提前于最大風(fēng)暴潮增水。各站位最大波浪增水均超過0.2 m，約占風(fēng)暴潮增水的6%。

圖7 各站位風(fēng)暴潮增水與波浪增水過程示意

波浪增水占風(fēng)暴潮增水的比例雖然不大，但對于設(shè)計(jì)潮位具有一定的影響。表1列舉了珠江口潮位站不同頻率下的設(shè)計(jì)潮位值，由表1可知，同一測站其不同頻率的設(shè)計(jì)潮位僅相差20～30 cm，如西炮臺(tái)、南沙、燈籠山200年一遇潮位與100年一遇設(shè)計(jì)潮位值分別相差0.17、0.25、0.22 cm，與最大波浪增水值相當(dāng)?？紤]波浪—風(fēng)暴潮耦合計(jì)算的方法，對推算沿岸設(shè)計(jì)潮位具有一定的參考意義。

表1 各測站潮位設(shè)計(jì)值

4 結(jié)語

1)本文在已建立風(fēng)暴潮數(shù)值模型的基礎(chǔ)上，耦合波浪模塊，建立了珠江河口波浪—風(fēng)暴潮耦合數(shù)值模型。采用實(shí)測波浪資料對“山竹”臺(tái)風(fēng)波浪過程進(jìn)行模擬驗(yàn)證，計(jì)算有效波高平均絕對誤差為0.52 m，最大有效波高相對誤差僅5.3%，模型計(jì)算結(jié)果較好地反映了臺(tái)風(fēng)期間的波浪過程，可用于珠江河口波浪—風(fēng)暴潮耦合數(shù)值模擬研究。

2)基于本文建立的波浪—風(fēng)暴潮耦合數(shù)值模型，以1822號臺(tái)風(fēng)“山竹”為例，討論了波浪對風(fēng)暴潮增水的影響，結(jié)果表明波浪對風(fēng)暴潮增水的影響僅限于近岸海域，最大值為10～20 cm，占風(fēng)暴潮增水的6.5%?？紤]波浪—風(fēng)暴潮耦合作用后，計(jì)算結(jié)果更接近觀測值，可提高風(fēng)暴潮數(shù)值預(yù)報(bào)的精度，同時(shí)對沿岸設(shè)計(jì)潮位的確定具有一定的參考意義。