河口海岸水沙動(dòng)力環(huán)境多尺度時(shí)空耦合數(shù)值模擬研究進(jìn)展

王芳宇，陳汨梨,，張啟博，謝守鵬，張繼生*

(1.河海大學(xué) 港口海岸與近海工程學(xué)院，南京 210024；2.中國(guó)港灣工程有限責(zé)任公司，北京 100027)

河口海岸是陸海相互作用的集中地帶，各種物理過程耦合多變，演變機(jī)制復(fù)雜，生態(tài)環(huán)境敏感脆弱。同時(shí)，該地帶經(jīng)濟(jì)較為發(fā)達(dá)、人口集居，日益增多的人類活動(dòng)會(huì)加劇河口海岸地區(qū)環(huán)境壓力。因此，河口海岸的研究、開發(fā)和保護(hù)是當(dāng)前世界沿海國(guó)家和科學(xué)家關(guān)注的熱點(diǎn)問題。

學(xué)者們通常采用物理模型試驗(yàn)、數(shù)值模擬和現(xiàn)場(chǎng)觀測(cè)等手段開展河口海岸水沙動(dòng)力環(huán)境的多尺度時(shí)空耦合作用研究。隨著計(jì)算機(jī)科學(xué)技術(shù)的快速發(fā)展，數(shù)值模擬憑借成本低、速度快、預(yù)測(cè)性強(qiáng)、結(jié)果可視化等優(yōu)勢(shì)成為當(dāng)前河口海岸水沙動(dòng)力環(huán)境研究的重要手段，并在工程實(shí)踐中得到廣泛應(yīng)用。因此，本文總結(jié)了波浪、潮流、泥沙同一物理場(chǎng)的多尺度時(shí)空耦合及不同物理場(chǎng)間的單一時(shí)空耦合模擬研究進(jìn)展，為今后河口海岸水沙動(dòng)力環(huán)境數(shù)值模擬提供研究思路及方向。

1 同一物理場(chǎng)多尺度空間耦合模擬研究

隨著河口海岸資源的持續(xù)開發(fā)與環(huán)境保護(hù)地不斷推進(jìn)，實(shí)現(xiàn)小規(guī)模沿海水沙動(dòng)力的準(zhǔn)確預(yù)測(cè)愈發(fā)重要。當(dāng)使用大范圍水沙動(dòng)力方法模擬大陸坡及小規(guī)模近海動(dòng)力環(huán)境時(shí)，由于大小范圍網(wǎng)格間的空間尺度差異，在描述垂向網(wǎng)格變化時(shí)較為困難[1]。近年來，填海工程造成的區(qū)域泥沙運(yùn)動(dòng)特征的改變以及對(duì)港口、航道和海岸生態(tài)系統(tǒng)的影響也愈來愈受到廣泛的關(guān)注，并且泥沙淤積帶來的問題尤其復(fù)雜。雖然目前同一物理場(chǎng)大范圍模擬的方法相對(duì)成熟，能準(zhǔn)確刻畫大范圍水沙動(dòng)力特性，但在局部水沙特性的刻畫上仍不夠精細(xì)。因此，如何考慮同一物理場(chǎng)不同空間尺度的相互作用，已成為準(zhǔn)確模擬河口海岸水沙動(dòng)力特性的關(guān)注重點(diǎn)。

1.1 單向嵌套

近年來，學(xué)者們發(fā)現(xiàn)雖計(jì)算流體力學(xué)方法能求解完整的Navier-Stokes方程，捕捉不同尺度的水流流動(dòng)現(xiàn)象。但是小尺度模擬對(duì)計(jì)算資源需求過大，仍不適用于工程實(shí)踐。另外，由于空間小尺度現(xiàn)象和大尺度現(xiàn)象間的相互影響不可忽視，迫切需要同步實(shí)現(xiàn)大范圍水動(dòng)力場(chǎng)的預(yù)測(cè)及局部動(dòng)力場(chǎng)的精確模擬。

基于此，學(xué)者開始著眼于同一物理場(chǎng)的多尺度空間嵌套模擬，以大范圍模型的模擬結(jié)果作為邊界條件驅(qū)動(dòng)小范圍模型。吳文挺等[2]通過嵌套模擬的方法搭建了高分辨率的模型網(wǎng)格，將臨近海域的海岸動(dòng)力過程降尺度傳遞到小尺度河口地區(qū)。隨后M?hring等[3]發(fā)展了流體動(dòng)力平面定位系統(tǒng)的多尺度建模方式，基于有限元軟件Ansys，使用Flotran求解器實(shí)現(xiàn)模型嵌套，證明了嵌套模型相較于簡(jiǎn)化模型的優(yōu)越性。Nash和Hartnett[4]利用MSN_Flood平臺(tái)開發(fā)了多尺度嵌套模型，該模型在精度和計(jì)算量間取得良好平衡且在復(fù)雜嵌套邊界上有獨(dú)特的處理方法(如圖1)。圖2展示了Nash等采用四級(jí)聯(lián)嵌套網(wǎng)格(90 m、30 m、6 m和2 m四個(gè)不同空間尺度)對(duì)較粗網(wǎng)格PG數(shù)據(jù)進(jìn)行空間內(nèi)插，以填充下一級(jí)較細(xì)子網(wǎng)格CG邊界數(shù)組(P為父網(wǎng)格，C為子網(wǎng)格)。嵌入PG90中的第一層子網(wǎng)格CG30以3:1(△x_P/△x_C)的嵌套比在重點(diǎn)區(qū)域縮小，在5:1的嵌套比下CG30為CG06提供東部邊界條件，最高分辨率的CG02則以3:1嵌套比嵌入CG06中[5]。由于小尺度水流常存在復(fù)雜邊界及源、匯項(xiàng)等問題，且模型精度受差分網(wǎng)格幾何形態(tài)限制[6]，因此在實(shí)現(xiàn)多尺度空間水流過程的高效模擬方面仍有欠缺[7]。

圖1 嵌套過程示意圖

泥沙數(shù)學(xué)模型作為研究和解決近海水域泥沙沖淤問題的有效工具已得到較為普遍的運(yùn)用[8]。一維模型計(jì)算省時(shí)但無法給出各物理量的平面范圍分布，雖二維模型能較準(zhǔn)確模擬但計(jì)算量較大經(jīng)濟(jì)性欠缺[9]。兼顧一、二維模型的優(yōu)缺點(diǎn)，Zhang[10]應(yīng)用一、二維嵌套技術(shù)，將河海連接處視為節(jié)點(diǎn)，由此成功地模擬了黃河口的演變。李玉婷等[11]考慮河口懸沙運(yùn)動(dòng)沖淤的變化，采用大小尺度雙重嵌套的方法建立了河口海岸泥沙二維模型，高峰等[12]在分析鴨綠江河口實(shí)測(cè)水文泥沙資料的基礎(chǔ)上，采用MIKE21 FM模塊建立了基于多重嵌套網(wǎng)格模式的泥沙場(chǎng)數(shù)學(xué)模型，驗(yàn)證表明該嵌套模型能較好地反映泥沙的運(yùn)動(dòng)特征，可以為今后開展進(jìn)一步的相關(guān)工程規(guī)劃和設(shè)計(jì)提供有效工具。張修忠等[13]也建立了一種河口海岸一、二維嵌套的泥沙數(shù)學(xué)模型，通過交界面實(shí)現(xiàn)數(shù)據(jù)在每一迭代步進(jìn)行耦合計(jì)算，能用較少的機(jī)時(shí)復(fù)演、預(yù)測(cè)重點(diǎn)部分的細(xì)部變形情況，驗(yàn)證結(jié)果表明在河口海岸處建立的一、二維嵌套泥沙數(shù)學(xué)模型是一種解決某些實(shí)際工程問題的可靠和高效的工具。

1.2 雙向耦合

跨空間尺度完成小范圍與大范圍網(wǎng)格之間的數(shù)據(jù)交互，準(zhǔn)確模擬多尺度空間耦合現(xiàn)象，對(duì)研究各物理場(chǎng)的本質(zhì)至關(guān)重要。因此，為更精細(xì)刻畫小尺度復(fù)雜水流過程，實(shí)現(xiàn)不同空間尺度下各流場(chǎng)信息的雙向互饋，曾季才等[14]建立多尺度空間數(shù)值離散的三維水流運(yùn)動(dòng)耦合模型，用以獲取全區(qū)域大尺度解及多子域小尺度高精度解。該模型采用了三維Darcy水頭插值算法[15]及尺度分離策略，在獲取區(qū)域大尺度粗網(wǎng)格解的同時(shí)精細(xì)刻畫局部小尺度水流運(yùn)動(dòng)過程。將不同尺度模型的空間信息關(guān)聯(lián)傳遞，以實(shí)現(xiàn)多尺度空間耦合模擬。相比傳統(tǒng)局部嵌套模型，該模型能更精確刻畫局部尺度的復(fù)雜水流運(yùn)動(dòng)過程。

在大規(guī)模沿海水流運(yùn)動(dòng)模擬中，Agnieszka等[16]首次提出POM-MSN_Flood建模耦合系統(tǒng)，并將其用于研究復(fù)雜潮汐引起的河口海岸水動(dòng)力。驗(yàn)證結(jié)果表明該模型能有效解決與水流特征尺度相對(duì)應(yīng)的流體運(yùn)動(dòng)，包括東北大西洋的大規(guī)模流動(dòng)和沿海水域精細(xì)環(huán)流。Zhang等[17]在建立近岸波流耦合的非結(jié)構(gòu)化四叉樹網(wǎng)格數(shù)值模型時(shí)，以樹狀數(shù)據(jù)結(jié)構(gòu)實(shí)現(xiàn)網(wǎng)格連通性和數(shù)據(jù)傳遞性。Wu等[18]為準(zhǔn)確模擬小范圍沿海水流的水動(dòng)力特性，將局部計(jì)算流體力學(xué)模型[19]與大范圍非結(jié)構(gòu)化網(wǎng)格海洋模型耦合，前者用于捕獲小尺度流，后者用于預(yù)測(cè)大尺度流。在局部模型中將小尺度流場(chǎng)劃分為相互重疊的雙向耦合子域，基于Schwartz交替迭代[20]和研發(fā)的質(zhì)量守恒算法[21]，實(shí)現(xiàn)子域間解的銜接過渡。如圖3所示，大范圍非結(jié)構(gòu)化網(wǎng)格海洋模型采用內(nèi)外模式求解離散化流動(dòng)，將局部模型與大范圍模型的三維內(nèi)模式相耦合，在網(wǎng)格交界面上交換速度分量，完成數(shù)據(jù)交互。

3-a 模型示意 3-b 網(wǎng)格示意

隨后Tang等[22]基于該耦合模型成功模擬了污染流從直徑為cm級(jí)的出口排放到寬度為km級(jí)的河口中并流向海域，并將模擬結(jié)果與僅采用局部模型方式計(jì)算的結(jié)果比較，證明該耦合方法的優(yōu)越性。與此同時(shí)，Tang等[23]還提出一種預(yù)測(cè)沿海海流的多尺度方法，將基于HM(Hybrid method)和DDM(Domain decomposition method)的計(jì)算流體力學(xué)和地球物理流體力學(xué)模型與Chimera交疊網(wǎng)格[24]相結(jié)合，實(shí)現(xiàn)了三維不可壓縮計(jì)算流體力學(xué)模型與大范圍非結(jié)構(gòu)化網(wǎng)格海洋模型的雙向耦合，并提供了數(shù)值示例證明其可行性，指出多尺度空間耦合模擬有很大的發(fā)展空間和應(yīng)用前景。

河口海岸泥沙輸運(yùn)是陸海物質(zhì)循環(huán)的重要過程，對(duì)物質(zhì)入海通量和全球物質(zhì)循環(huán)具有重要研究意義。劉紀(jì)根等[25]認(rèn)為當(dāng)前泥沙輸移過程隨流域尺度變化的研究缺乏應(yīng)有的廣度和深度，總結(jié)了泥沙沖刷輸運(yùn)的空間尺度相關(guān)問題，探討了不同尺度下的泥沙輸運(yùn)現(xiàn)象與規(guī)律。朱志夏等[26]針對(duì)三維泥沙數(shù)值模擬難題，采用FVCOM模型及非結(jié)構(gòu)化網(wǎng)格二維、三維嵌套方法，建立了三維泥沙數(shù)學(xué)模型。Bai等[27]對(duì)河口海岸地區(qū)泥沙輸移現(xiàn)象展開研究，采用有限元法FEM和有限差分法FDM相結(jié)合的方法提出海河河口模型，并在三維指標(biāo)系統(tǒng)中引入多層模型，最后將模型應(yīng)用于渾濁的華北海河河口。通過對(duì)比數(shù)值結(jié)果與實(shí)測(cè)資料，證明該多空間尺度泥沙模型與美國(guó)CH3D-SED模型具有可比性。

本節(jié)回顧了學(xué)者針對(duì)不同物理場(chǎng)展開的多尺度空間耦合模擬，由于小尺度和大尺度水、沙相互作用，當(dāng)進(jìn)行多尺度空間耦合模擬時(shí)，因大小網(wǎng)格的空間尺度差異性，通常會(huì)出現(xiàn)大范圍水動(dòng)力的一個(gè)網(wǎng)格涵蓋了十至百個(gè)精細(xì)網(wǎng)格的現(xiàn)象。隨著數(shù)值模擬手段的不斷深入，各物理量的空間尺度跨度是逐漸增大的，例如，泥沙的大規(guī)模輸運(yùn)、海工建筑物的局部沖淤、泥沙顆粒運(yùn)動(dòng)規(guī)律等物理現(xiàn)象的模擬，均從宏觀到微觀不斷深入，并且在該過程中空間尺度量級(jí)差是不可避免的。因此，如何在數(shù)值模擬中有效解決空間尺度量級(jí)差帶來的問題，已成為真實(shí)模擬各物理場(chǎng)的研究重點(diǎn)及難點(diǎn)。

2 同一物理場(chǎng)多尺度時(shí)間耦合模擬研究

將不同時(shí)間尺度的物理場(chǎng)進(jìn)行耦合模擬，預(yù)測(cè)未來時(shí)間內(nèi)該物理場(chǎng)的真實(shí)變化，可在某種程度上滿足當(dāng)前對(duì)長(zhǎng)時(shí)間變化物理現(xiàn)象的現(xiàn)實(shí)需求。雖前人曾采用統(tǒng)計(jì)學(xué)模型的方式對(duì)河口海岸演變進(jìn)行預(yù)測(cè)并得到?jīng)_淤平衡點(diǎn)[27-30]，但研究性質(zhì)僅限于規(guī)律，無法從機(jī)理上對(duì)長(zhǎng)期演變進(jìn)行很好地解釋，并且也無法對(duì)海岸變化的具體位置進(jìn)行預(yù)測(cè)。因此多時(shí)間尺度模擬預(yù)測(cè)未來長(zhǎng)時(shí)間變化是有必要的，如泥沙輸移的多尺度時(shí)間數(shù)值模擬被普遍認(rèn)為是研究河口海岸地貌演變的重要工具?；诖耍瑢W(xué)者們逐漸關(guān)注各物理場(chǎng)的多尺度時(shí)間耦合模擬，期盼能以當(dāng)前時(shí)間尺度預(yù)測(cè)未來時(shí)間內(nèi)物理場(chǎng)變化，對(duì)厘清河口海岸水沙動(dòng)力環(huán)境及相互作用機(jī)制具有重要意義[31]。

2.1 加速因子

當(dāng)模擬中長(zhǎng)期河口海岸動(dòng)力地貌演變時(shí)，通常對(duì)地貌演變的時(shí)間進(jìn)行加速以快速模擬動(dòng)力地貌的中長(zhǎng)期變化?；诖擞^念，Lesser[32]提出地貌加速因子法，這是一種用來解決地貌發(fā)展過程中出現(xiàn)時(shí)間尺度差異的有效方法。這一技術(shù)與Latteux[33]提出的“延長(zhǎng)潮汐”概念相類似，旨在用短時(shí)間泥沙運(yùn)動(dòng)模擬實(shí)現(xiàn)長(zhǎng)時(shí)間地貌演變模擬。

地貌加速因子法通過引入一個(gè)地貌加速因子(morphological acceleration factor)，有效擴(kuò)展地貌時(shí)間步長(zhǎng)以完成地貌變化的快速計(jì)算。通過乘上一個(gè)常量因子增加地形變化率，經(jīng)過1個(gè)潮周期的計(jì)算，模型實(shí)際上已對(duì)n(n=fMOR)個(gè)周期的地貌變化進(jìn)行了模擬，如式(1)

△tmorphology=fMOR△thyfordynamic

(1)

式中：fMOR為地貌加速因子；△t為時(shí)間步長(zhǎng)。

由于地貌演變作用的復(fù)雜性與時(shí)空尺度多樣性，當(dāng)前對(duì)河口海岸動(dòng)力地貌演變的理解與預(yù)測(cè)仍非常有限。Ganju等[34]首先利用ROMS模型展開河口地區(qū)多尺度時(shí)間嵌套的地貌模擬，概化了入海輸沙量并引入地貌加速因子。Ganju等[35-36]也對(duì)美國(guó)Suisun灣1867年—1887年地貌沖淤開展了模擬，與實(shí)際地形變化對(duì)比預(yù)測(cè)了Suisun灣2010年—2030年際地貌演變[37]，并指出模擬年代際尺度上的地形沖淤變化難度較大。

荷蘭學(xué)者們[38]利用Delft3D模型分別模擬研究了Scheldt河口西部的長(zhǎng)時(shí)間地貌演變，實(shí)現(xiàn)多時(shí)間尺度耦合模擬研究。董程等[39]從黃河潮灘地貌演化角度出發(fā)，認(rèn)為采用地貌加速因子后一個(gè)潮周期的地貌變化等價(jià)于10個(gè)潮周期的地貌，表明以多尺度時(shí)間方式模擬地貌變化是非常必要的。蔣超[40]系統(tǒng)分析了黃河入海水沙通量變化特征，揭示了黃河口不同時(shí)間尺度地貌演變的時(shí)空變化特征，為不同類型動(dòng)力地貌的多尺度時(shí)間耦合模擬建立重要紐帶。Luan等[41-42]闡述了主導(dǎo)長(zhǎng)江口年代際地貌沖淤演變的影響因素，基于Delft3D模型建立長(zhǎng)江口年代際地貌演變數(shù)值模型，預(yù)測(cè)了未來年代際演變的動(dòng)態(tài)過程。在珠江口方面，Xing和Wu[43]基于動(dòng)力地貌數(shù)值模型模擬了珠江口長(zhǎng)時(shí)間沖淤變化過程，指明河口地形的長(zhǎng)期演變與復(fù)雜的河網(wǎng)格局密切相關(guān)。吳超羽等[44]綜合考慮影響珠江口動(dòng)力地貌演變的驅(qū)動(dòng)因子，建立長(zhǎng)時(shí)間動(dòng)力地貌演變數(shù)值耦合模型，模擬了珠江口幾千年地貌、河網(wǎng)演變過程。考慮工程實(shí)際情況，竇希萍等[45]認(rèn)為在工程初期地形沖淤變化往往較劇烈，此時(shí)不能選擇大的加速因子。隨著水沙動(dòng)力與河床不斷調(diào)整，地形沖淤變化趨緩，可選擇較大加速因子以減小地貌長(zhǎng)周期模擬計(jì)算時(shí)間。當(dāng)模擬工程影響下的地貌演變時(shí)，采用隨時(shí)間變化的加速因子更能提高模擬精度和計(jì)算效率。但在中長(zhǎng)期地貌演變模擬過程中使用加速因子，模型需要設(shè)置對(duì)應(yīng)尺度的上游和外海邊界條件。如果仍采用實(shí)測(cè)邊界條件，那么水動(dòng)力過程在時(shí)間尺度上會(huì)被人為“拉長(zhǎng)”，一個(gè)潮周期對(duì)應(yīng)的地貌時(shí)間可能為數(shù)天甚至數(shù)月，顯然不符合實(shí)際物理意義，因此在選擇相應(yīng)邊界條件需要考量得當(dāng)。

2.2 長(zhǎng)期運(yùn)算

由于計(jì)算資源有限，以往的數(shù)值模擬多為短時(shí)間模擬，并且大多數(shù)基于物理機(jī)制的海岸演化模型為短時(shí)間尺度下“真實(shí)性演變過程”的一個(gè)描述，難以刻畫中長(zhǎng)時(shí)間尺度的海岸演化過程和機(jī)制[46-48]。隨著計(jì)算技術(shù)進(jìn)步，長(zhǎng)期運(yùn)算模擬手段逐漸完善。如前所述，為預(yù)測(cè)未來長(zhǎng)周期河口海岸水沙動(dòng)力變化，常用手段是使用加速因子加快時(shí)間進(jìn)程，從而縮短模擬時(shí)長(zhǎng)。但是，Roelvink等[49]指出使用加速因子的前提是水動(dòng)力計(jì)算時(shí)間步長(zhǎng)內(nèi)泥沙凈輸運(yùn)量和地形沖淤變化很小，乘以加速因子后的地形沖淤變化不會(huì)對(duì)未來年際內(nèi)地貌格局產(chǎn)生較大影響。

基于上述觀點(diǎn)，學(xué)者們嘗試采用長(zhǎng)期運(yùn)算手段模擬河口海岸水沙動(dòng)力的長(zhǎng)期演變，以期預(yù)測(cè)未來更長(zhǎng)時(shí)間尺度內(nèi)各物理量的變化過程。Roelvink等[50]基于二維(水平)模型模擬復(fù)雜海岸演變，認(rèn)為該模型可數(shù)值再現(xiàn)多年海岸形態(tài)，但計(jì)算成本高且需要豐富專業(yè)知識(shí)。為克服現(xiàn)有海岸線模型局限性，他們開發(fā)了一種全新的海岸線模擬模型(Shorelines)，旨在預(yù)測(cè)數(shù)年到數(shù)世紀(jì)的海岸線演變。研究結(jié)果證明該模型不僅能模擬漸進(jìn)的海岸變化，并且能預(yù)測(cè)未來時(shí)間尺度內(nèi)海岸的根本性轉(zhuǎn)變。該模型適用于局部到區(qū)域甚至全球領(lǐng)域，以合理的計(jì)算成本在完整生命周期中穩(wěn)健跟蹤海岸特征，捕捉海岸演化的復(fù)雜性，具有顯著的優(yōu)越性。

發(fā)展中長(zhǎng)時(shí)間尺度的海岸演化預(yù)測(cè)方法與構(gòu)建機(jī)理模型是目前國(guó)際海岸演化研究領(lǐng)域關(guān)注的焦點(diǎn)和極具挑戰(zhàn)性的主題。Hanson等[47]對(duì)目前已有的基于短時(shí)間尺度發(fā)展的20個(gè)不同類型年代際尺度海岸演化模型進(jìn)行對(duì)比分析和討論，結(jié)果表明這類模型大多不能進(jìn)行長(zhǎng)時(shí)間尺度海岸演化過程的準(zhǔn)確模擬。基于此，李國(guó)勝等[51]提出將一線模型(one-line model)與三維海洋動(dòng)力學(xué)數(shù)值模型耦合，經(jīng)過參數(shù)化后形成真正意義上的機(jī)理模型，能在長(zhǎng)時(shí)間尺度上對(duì)河口海岸演化的驅(qū)動(dòng)機(jī)制和影響機(jī)理進(jìn)行深入地研究。錢繼春等[52]基于潮流泥沙數(shù)學(xué)模型和灘涂動(dòng)態(tài)平衡的經(jīng)驗(yàn)關(guān)系，建立了長(zhǎng)歷時(shí)灘涂動(dòng)態(tài)演變的經(jīng)驗(yàn)動(dòng)力數(shù)值模型，結(jié)果表明數(shù)值結(jié)果與實(shí)測(cè)基本吻合，能真實(shí)反映灘涂的動(dòng)態(tài)演變過程。王寧舸[53]隨后對(duì)江蘇中部流域的河口潮灘開展1 a現(xiàn)場(chǎng)觀測(cè)，建立平面二維泥沙輸運(yùn)及潮灘演變概化數(shù)學(xué)模型，模擬2.5 a內(nèi)江蘇中部潮灘變化，預(yù)測(cè)了未來年際內(nèi)潮灘演變規(guī)律。程梁秋等[54]基于多年海圖數(shù)據(jù)建立了相對(duì)應(yīng)的數(shù)字高程模型，對(duì)伶仃洋西槽地貌展開長(zhǎng)達(dá)15 a的數(shù)值模擬，定量分析其長(zhǎng)期演變特征。龔政等[55]以淤長(zhǎng)型潮灘為研究對(duì)象，建立了潮灘地貌的演變概化數(shù)學(xué)模型，預(yù)測(cè)了未來16 a內(nèi)潮灘剖面的發(fā)育和演變過程，并與經(jīng)典潮灘一維數(shù)值結(jié)果對(duì)比，證明建立的長(zhǎng)期運(yùn)算模型具有較高可靠性。鄭釗[56]考慮黃河沖淤演變建立了數(shù)學(xué)模型以此預(yù)測(cè)未來50 a黃河的沖淤演變，并且驗(yàn)證了預(yù)測(cè)成果的可靠性，該預(yù)測(cè)成果對(duì)未來黃河河口治理具有重要的實(shí)際價(jià)值。

石小翠[57]從預(yù)防海洋災(zāi)害角度出發(fā)，運(yùn)用SWAN模型嵌套模擬計(jì)算了中國(guó)東部海區(qū)22 a波要素，以此構(gòu)成中國(guó)東部海區(qū)的波浪數(shù)據(jù)庫。Liang等[58]也采用SWAN模型以當(dāng)前時(shí)間尺度模擬了渤海、黃海和東海21 a的波浪要素。Langodan等[59]采用ERA-Interim全球再分析初始化WRF模式生成高分辨率區(qū)域驅(qū)動(dòng)WAVEWATCH III模式，利用30 a波浪要素追算了紅海波浪變化趨勢(shì)。徐佳麗等[60]采用TOMAWAC模型模擬中國(guó)近海40 a的波浪要素，建立相對(duì)完整的數(shù)據(jù)庫以全面描述中國(guó)海區(qū)波浪情況，以此為依據(jù)探究河口海岸地區(qū)動(dòng)力演變過程。

3 多物理場(chǎng)單一尺度時(shí)空耦合模擬研究

河口海岸水沙動(dòng)力環(huán)境的演變是個(gè)復(fù)雜的多時(shí)空尺度過程，不同時(shí)空尺度的演化過程及動(dòng)力機(jī)制存在顯著差異，相關(guān)研究的理論基礎(chǔ)及技術(shù)方法也各不相同[61]。自20世紀(jì)以來，學(xué)者致力于數(shù)值模擬各物理現(xiàn)象，但已有研究主要聚焦于各種相對(duì)孤立的物理效應(yīng)。在通常情況下，聚焦主要物理場(chǎng)時(shí)會(huì)采用概化方法模擬其他物理場(chǎng)[62]。對(duì)比學(xué)者們的模擬手段發(fā)現(xiàn)，聚焦主要物理場(chǎng)并概化次要物理場(chǎng)是常見且實(shí)用的研究手段。然而，現(xiàn)實(shí)中的物理現(xiàn)象并不是孤立發(fā)生的，只重點(diǎn)關(guān)注單一物理場(chǎng)會(huì)導(dǎo)致模擬結(jié)果與真實(shí)情況有出入。

隨著計(jì)算能力不斷提升，為更全面真實(shí)模擬河口海岸水沙動(dòng)力環(huán)境，學(xué)者們開始在數(shù)學(xué)模型中考慮多物理場(chǎng)的相互耦合作用。

3.1 單向嵌套

在河口海岸地區(qū)波浪和潮流是主要驅(qū)動(dòng)力，其過程復(fù)雜且呈非線性[63]。近岸區(qū)波浪在傳播時(shí)會(huì)受到水流的影響，波浪自身也會(huì)產(chǎn)生近岸水流驅(qū)動(dòng)力，兩者存在耦合作用，并且波浪與潮流在時(shí)間、空間尺度上差異顯著，因此研究河口海岸地區(qū)的波流相互作用具有一定意義[64]。當(dāng)前波流相互作用的研究主要采用兩種概化方法：一是把周期變化的潮流概化為具有某一特征流疊加到波浪運(yùn)動(dòng)方程中，用以模擬短時(shí)間內(nèi)波流運(yùn)動(dòng)的變化；另一種方法是把波浪作用過程概化為潮周期中具有平均意義的波浪要素疊加到潮流運(yùn)動(dòng)方程中，用以描述波對(duì)流的影響。

波流模型經(jīng)歷了從一維到二維再到三維、由單向嵌套再到雙向耦合等由簡(jiǎn)單到復(fù)雜的發(fā)展歷程[65]，Longue-Higgins和Stewart[66]首先基于輻射應(yīng)力概念闡明了波流相互作用機(jī)理以及波流的能量轉(zhuǎn)換關(guān)系；Janssen[67]隨后提出波的阻力關(guān)系，進(jìn)一步豐富波流相互作用理論。在過去幾十年，特別是近20 a，學(xué)者們從不同角度對(duì)波流相互作用展開細(xì)致模擬研究。Teles等[68]基于RANS方程將水流模塊嵌套至波浪模塊，采用任意拉格朗日-歐拉(ALE)方法模擬時(shí)變自由表面動(dòng)力學(xué)，在波浪時(shí)空尺度上同時(shí)求解波流運(yùn)動(dòng)數(shù)值解，分析自由表面流動(dòng)中波流的聯(lián)合效應(yīng)。許忠厚等[69]基于Delft3D建立了三維波流耦合模型，模擬了波浪作用下出水離岸堤附近的波浪與波生流間的關(guān)系，Zhang等[17]通過求解二維深度平均非線性淺水方程，在沿海水域的水動(dòng)力模塊加入波浪誘導(dǎo)的沖浪區(qū)近水平流，波流相互作用模擬結(jié)果與實(shí)測(cè)數(shù)據(jù)吻合良好。Sanchez等[70]采用二維深度平均方程求解波流作用下潮流從十幾公里寬的Colombia河口流向幾十萬平方公里海域的地形變化，實(shí)現(xiàn)了基于波流相互作用的地貌演變模擬。

在水沙動(dòng)力環(huán)境下，泥沙常隨水流一同運(yùn)動(dòng)。由于水流是周期性的往復(fù)運(yùn)動(dòng)流，其作用下的泥沙運(yùn)動(dòng)也是往復(fù)運(yùn)動(dòng)，運(yùn)動(dòng)機(jī)理要比河流中單向恒定作用下的泥沙運(yùn)動(dòng)機(jī)理復(fù)雜得多[71]。因此，研究近岸區(qū)泥沙輸移及岸灘演變具有十分重要的意義[72]。受限于計(jì)算流體力學(xué)理論不成熟，在過去的幾十年中水沙模型多為非耦合模型。近20 a來得益于計(jì)算效率的提高，吳偉明等[73]通過處理連接交界面讓一維和二維水沙數(shù)學(xué)模型滿足水流、泥沙運(yùn)動(dòng)連接的邊界條件，從而構(gòu)建了水沙運(yùn)動(dòng)一維、二維嵌套模型。張華慶等[74]采用正交網(wǎng)格建立了珠江河口平面二維潮流泥沙數(shù)學(xué)模型，將珠江三角洲一維河網(wǎng)模型嵌套至珠江河口二維模型，在不同空間尺度上構(gòu)建了珠江河口水沙數(shù)值模擬系統(tǒng)的核心計(jì)算模塊。張世奇[75]在河口采用一維模式并在海岸運(yùn)用二維模式，將入?？谧鳛楹雍＿B接點(diǎn)，既作為一維模擬的下邊界，又作為二維模擬的上邊界。根據(jù)水流連續(xù)原理交替使用一維、二維尺度手段求解水位流量，采用擬合的方法使其趨于一致，由此建立一維、二維連接的沖淤嵌套模型，結(jié)果表明該模型能較準(zhǔn)確模擬黃河口沖淤地形和輸沙過程。

3.2 雙向耦合

隨著科學(xué)進(jìn)步以及嵌套技術(shù)的成熟，信息的雙向傳遞已經(jīng)成為河口海岸水沙動(dòng)力環(huán)境數(shù)值模擬研究的關(guān)鍵，多物理場(chǎng)的數(shù)值模擬方式也從單向嵌套逐步轉(zhuǎn)到雙向耦合。Dietrich等[76]采用區(qū)域分解并行仿真技術(shù)，在同一非結(jié)構(gòu)化網(wǎng)格上交換網(wǎng)格頂點(diǎn)的解以實(shí)現(xiàn)SWAN和ADCIRC耦合，有效解決了波流相互作用的模擬難題。該模型具有高度可擴(kuò)展性和通信性，在較大空間梯度上提高了局部分辨率，且能準(zhǔn)確地模擬波流相互作用。王平等[77]基于非結(jié)構(gòu)化網(wǎng)格同時(shí)模擬潮流和波浪，并通過參數(shù)的同步傳遞實(shí)現(xiàn)波流的耦合計(jì)算，結(jié)果表明模型對(duì)模擬近岸波流的耦合作用有著很好的精度和適用性。王文鼎等[78]采用大尺度和小尺度兩種尺度的網(wǎng)格建立波流耦合數(shù)值模型，并將研究結(jié)論應(yīng)用于部分工程實(shí)踐，驗(yàn)證了模型的符合性、有效性和適用性。Roland等[79]基于三維水動(dòng)力模型SELFE和波譜模型WWM-II，提出一種基于非結(jié)構(gòu)網(wǎng)格的波流全耦合模型。SELFE和WWM-II模型間的耦合是在源代碼級(jí)別進(jìn)行的，在共享相同子域的基礎(chǔ)上并行MPI和域分解，從根本上保證并行效率和避免插值。Roland等[79]隨后將模擬結(jié)果和現(xiàn)場(chǎng)資料比較，結(jié)果表明即使在復(fù)雜的實(shí)際環(huán)境下，耦合模型也能捕捉到不同尺度下波流相互作用過程。

在水流問題的研究上，王智勇等[80]提出耦合邊界水位預(yù)測(cè)校正法，由此建立一維和二維耦合水動(dòng)力模型。在任一時(shí)間步計(jì)算時(shí)先在耦合連接處預(yù)測(cè)一維和二維水位邊界條件，再根據(jù)耦合邊界處的流量計(jì)算結(jié)果及當(dāng)?shù)厮?shù)對(duì)水位邊界條件反復(fù)校正，直至耦合連接條件滿足計(jì)算容差。顧杰等[81]基于MIKE-FLOOD模塊建立秦皇島入海河流與近岸海域間的一維、二維耦合河流-海岸水動(dòng)力模型，系統(tǒng)研究了秦皇島海域入海河流的水動(dòng)力特性。來志剛[82]在珠江口實(shí)現(xiàn)了一維河網(wǎng)、三維河口的水動(dòng)力耦合模擬，成功將一維、二維水動(dòng)力及泥沙連接模型基本條件和計(jì)算方法擴(kuò)展至一維、三維水動(dòng)力模型的耦合模擬。Zhou等[83]基于質(zhì)量守恒準(zhǔn)則，耦合了一維河網(wǎng)模型與三維河口模型，為珠江三角洲的水環(huán)境監(jiān)控提供可靠依據(jù)。同時(shí)，Lai等[84]著眼于長(zhǎng)江中游復(fù)雜河湖大水系，提出一種新的耦合水動(dòng)力模型CHAM，為河口海岸地區(qū)的復(fù)雜流態(tài)提供了有效模擬工具。

準(zhǔn)確描述河口海岸的水流運(yùn)動(dòng)是合理描述泥沙運(yùn)動(dòng)的基礎(chǔ)，目前相關(guān)學(xué)者在泥沙數(shù)值模擬方面已開展了相關(guān)研究并取得一定成果[85]，張修忠等[86]指出對(duì)于解決河口海岸這類具有多空間尺度的水沙問題，采用單一模型不僅困難、費(fèi)時(shí)，且不易掌握不同區(qū)域矛盾的主要特征?；诖藦埿拗业萚86]利用有限元的非結(jié)構(gòu)化網(wǎng)格特點(diǎn)和子結(jié)構(gòu)疊加的概念，建立了一種虛擬過度區(qū)間的一、二維耦合算法。蘇東升等[87]基于CFD-DEM 方法建立水沙運(yùn)動(dòng)耦合模型，張華慶等[88]考慮水沙相互關(guān)系建立了一維、二維耦合的水沙模型，將河口地區(qū)視為重點(diǎn)區(qū)域進(jìn)行計(jì)算。張鵬飛和陸建剛[89]研究了通江湖泊與外部江河的水沙交換過程，基于二維非穩(wěn)態(tài)水流和泥沙耦合數(shù)學(xué)模型，模擬了鄱陽湖與外部江河的水沙交換過程。Li等[90]基于三維水動(dòng)力模型與泥沙運(yùn)動(dòng)模型，耦合模擬了北部灣海岸地區(qū)沙波的遷移過程，模擬結(jié)果與實(shí)測(cè)基本吻合，表明多場(chǎng)耦合模型可以模擬現(xiàn)實(shí)中的沙波遷移規(guī)律。在驗(yàn)證耦合模型能較好模擬水沙運(yùn)動(dòng)后，Shen等[91]考慮河口和沿海地區(qū)水域懸浮沉積物濃度的高變異性特征，建立了泥沙模型與開源TELEMAC系統(tǒng)耦合的平衡模型。葉濤焱[92]為探究杭州灣人類活動(dòng)影響下的懸沙濃度，建立了三維水沙數(shù)值耦合模型ESed，分析了杭州灣灣口懸沙動(dòng)力與潮灘變化的互饋機(jī)理。

由此可見，在多物理場(chǎng)單一尺度時(shí)空耦合數(shù)值模擬中，通常在求解控制方程時(shí)考慮多因子的相互影響。本節(jié)總結(jié)了學(xué)者在河口海岸地區(qū)波浪、潮流以及泥沙方面的研究成果，說明了考慮不同物理場(chǎng)的相互作用能更真實(shí)反應(yīng)河口海岸的水沙動(dòng)力環(huán)境的變化。但目前數(shù)值模擬多為單一尺度時(shí)空耦合，在準(zhǔn)確刻畫不同時(shí)空尺度下各物理場(chǎng)間的相互作用上仍有所欠缺。

4 模擬難題與建議

在河口海岸水沙動(dòng)力環(huán)境的數(shù)值模擬中，模擬手段從單獨(dú)模擬到單向整合再發(fā)展至雙向耦合。以往的模型大多僅實(shí)現(xiàn)信息單向傳遞，即給定初始、邊界條件時(shí)A模型先計(jì)算，再將A模擬結(jié)果傳至B模型，此時(shí)B模型結(jié)果并未反饋給A模型。但是，較為理想的是結(jié)果相互反饋：A將結(jié)果傳遞至B，B同時(shí)反饋至A，循環(huán)迭代直到各物理場(chǎng)的計(jì)算達(dá)到穩(wěn)定。在整個(gè)模擬過程中，數(shù)據(jù)在每一時(shí)間步都進(jìn)行交換，據(jù)此實(shí)現(xiàn)兩個(gè)模型間的實(shí)時(shí)雙向耦合。

實(shí)際河口海岸水沙動(dòng)力環(huán)境是多物理場(chǎng)相互作用的結(jié)果，早在1993年Orton和Reading[93]就提出河口系統(tǒng)每一個(gè)地貌單元并不是獨(dú)立的個(gè)體，而存在頻繁的物質(zhì)交換和相互作用過程。考慮不同物理場(chǎng)相互作用時(shí)，數(shù)據(jù)相互反饋使耦合模擬結(jié)果更精確。隨著社會(huì)發(fā)展及行業(yè)需求，研究重點(diǎn)、熱點(diǎn)也更加聚焦于多物理場(chǎng)的多尺度時(shí)空耦合模擬，河口海岸地區(qū)多物理場(chǎng)間的相互作用問題也越來越廣泛受到人們的重視，例如：海冰、波浪、水流的相互作用[94]，波浪、水流、泥沙的相互作用[95]，風(fēng)、波浪、水流的相互作用[96]，大氣、海浪、海流的相互作用等[97]。與此同時(shí)，多物理場(chǎng)相互作用現(xiàn)象越發(fā)復(fù)雜，水沙動(dòng)力與結(jié)構(gòu)物響應(yīng)、波浪-海床-結(jié)構(gòu)物[98]相互作用等問題也逐漸突出，多尺度計(jì)算已經(jīng)發(fā)展成為一個(gè)具有自身特征的學(xué)科，是構(gòu)成人們探索和理解物質(zhì)世界多尺度問題的基礎(chǔ)[99]。因此，在求解多物理場(chǎng)相互作用時(shí)，如何客觀描述控制方程是提高數(shù)值模擬精度、實(shí)現(xiàn)多尺度時(shí)空耦合的研究重點(diǎn)，也是亟待解決的關(guān)鍵問題。

受限于計(jì)算技術(shù)及耦合算法的難度，多物理場(chǎng)的全耦合模擬尚未完全實(shí)現(xiàn)。結(jié)合現(xiàn)今技術(shù)手段，采用以下三種方法實(shí)現(xiàn)模型間耦合較可行。第一種方法是將兩種模型代碼直接合并，基于內(nèi)存實(shí)現(xiàn)模型間的數(shù)據(jù)傳遞，但通常兩個(gè)模型是獨(dú)立開發(fā)，因此數(shù)據(jù)的交互具有一定難度。第二種方法是利用現(xiàn)有通信工具或協(xié)議如MPI，但該方法需模型開發(fā)人員具備專業(yè)知識(shí)。第三種方法是使用模型耦合工具提供自動(dòng)傳輸數(shù)據(jù)接口，這種方法能極大減少模型開發(fā)期間修改代碼的工作量，并且非專業(yè)技術(shù)人員也能實(shí)現(xiàn)模型接口配置，相比下具有更優(yōu)性。雖第三種方法更具優(yōu)越性，但關(guān)于耦合模型接口的搭建，實(shí)現(xiàn)數(shù)據(jù)的自動(dòng)傳輸難度很大。因此如何開發(fā)模型耦合工具及使用好該工具，實(shí)現(xiàn)信息交互、耦合的自動(dòng)化以及模型結(jié)果的相互反饋是學(xué)者們非常期待的。

5 結(jié)語

河口海岸是地球四大圈層交匯、能量流和物質(zhì)流的重要聚散地帶，既是流域物質(zhì)的歸宿，又是海洋的開始，生態(tài)環(huán)境極其敏感、脆弱。高精度數(shù)值模擬手段是未來河口海岸水沙動(dòng)力耦合作用機(jī)理及其環(huán)境影響效應(yīng)評(píng)估等研究的主要手段和發(fā)展趨勢(shì)。本文回顧了河口海岸地區(qū)潮流、波浪及泥沙等多物理場(chǎng)間的耦合作用及同一物理場(chǎng)下多尺度空間、時(shí)間耦合作用，探討了多尺度時(shí)空耦合過程中周期量級(jí)差、尺度量級(jí)差、各物理場(chǎng)間的相互耦合、多尺度時(shí)空耦合技術(shù)以及多物理場(chǎng)多尺度時(shí)空耦合等關(guān)鍵難題，指出現(xiàn)階段實(shí)現(xiàn)河口海岸多物理場(chǎng)多尺度時(shí)空數(shù)值耦合的模擬仍面臨許多挑戰(zhàn)，急需展開深入系統(tǒng)研究和工程實(shí)踐。