Nudging同化方案改進(jìn)及其在西北太平洋海域M2分潮模擬中的應(yīng)用＊

魏澤勛，王永剛，高秀敏，方國洪

（1.國家海洋局 第一海洋研究所，山東 青島266061；2.青島海洋科學(xué)與技術(shù)國家實驗室 區(qū)域海洋動力學(xué)與數(shù)值模擬功能實驗室，山東 青島266071）

潮汐是海洋最基本的運動形式之一。它與人類活動關(guān)系密切，準(zhǔn)確掌握潮汐的運動規(guī)律對漁業(yè)生產(chǎn)和生態(tài)環(huán)境都具有重要意義。此外，航運交通、海港工程、海岸防護(hù)、能源利用、環(huán)境保護(hù)以及軍事活動等，都受潮汐現(xiàn)象的影響。潮汐學(xué)的研究還對其他方面如海洋污染物擴(kuò)散、泥沙運動、水體交換、風(fēng)暴潮等問題的研究有重要作用。因此，針對潮汐進(jìn)行研究具有重要的理論意義和實際應(yīng)用價值。

對于中國近海，早期的潮汐研究主要是基于觀測資料［1－2］，這些觀測資料大多集中于沿岸和島嶼附近。隨著計算技術(shù)的發(fā)展，數(shù)值模擬方法逐漸被普遍應(yīng)用于潮汐研究［3－16］，這些研究大多集中于中國近海的某一海區(qū)，或者將渤、黃、東海作為一個整體。

近年來，觀測資料尤其是衛(wèi)星高度計資料的大量涌現(xiàn)，使潮汐學(xué)研究有了更多、更全面的資料［17－22］。如何更有效率地利用這些資料，是擺在人們面前的重要問題。多年來，人們嘗試將觀測資料與數(shù)值模型相結(jié)合，從而使數(shù)值模擬結(jié)果與觀測盡可能接近，這種將觀測資料與數(shù)值模型相結(jié)合的方法稱為資料同化方法。早在上世紀(jì)70年代，方國洪［23］就成功建立了潮波變分同化模式，并用于模擬黃海M2分潮。隨后將同化方法用于潮波模擬的研究多集中在2000年以后。韓桂軍等［24－25］采用基于最優(yōu)控制理論的伴隨法建立伴隨模型，利用驗潮站的水位資料以及TOPEX／Poseidon（以下簡稱T／P）衛(wèi)星測高數(shù)據(jù)在黃海、東海進(jìn)行變分?jǐn)?shù)據(jù)同化試驗，優(yōu)化了開邊界條件，提高了數(shù)值模擬的精度。呂咸青和方國洪［26－27］采用伴隨法，由渤海沿岸19個驗潮站的潮汐調(diào)和常數(shù)來反演渤海海域的開邊界條件，數(shù)值模擬結(jié)果較好地體現(xiàn)了渤海M1和M2潮波的基本特征。吳自庫等［28］利用伴隨同化方法，把T／P衛(wèi)星高度計資料提取得到的沿軌分潮調(diào)和常數(shù)同化到二維非線性潮汐數(shù)值模式中去，通過優(yōu)化模型中的開邊界條件和底摩擦系數(shù)，較好地模擬了南海M1和M2分潮的潮汐。張繼才和呂咸青［29－30］探討了獨立底摩擦系數(shù)如何選取，才能使得渤、黃、東海M2分潮伴隨同化數(shù)值模擬的精度提高。丘仲鋒等［31］將大約10a的T／P高度計資料，利用伴隨同化方法，同化到二維非線性潮汐數(shù)值模式中，模擬了黃海、渤海區(qū)域M2，S2，O1和K1分潮。Zu等［32］通過將T／P衛(wèi)星高度計資料同化到正壓潮波模型中，模擬了南海8個主要分潮（M2，S2，K1，O1，N2，K2，P1和Q1）。以上都是變分同化方法在潮波模式中的應(yīng)用，趨近法同化也常用于潮波數(shù)值模式中。王永剛等［33］采用趨近法，建立了渤、黃、東海潮汐同化數(shù)值模式。丘仲鋒和何宜軍［34］利用趨近插值法將高度計資料同化到二維非線性潮汐數(shù)值模式中，模擬了中國近海M2分潮的分布。朱學(xué)明等［35］利用FVCOM海洋數(shù)值模式，采用趨近法同化84個沿岸驗潮站的觀測資料，模擬了西北太平洋海域M2，S2，K1和O1四個主要分潮。

本文基于POM模式，采用Nudging同化技術(shù)，建立了西北太平洋海域潮波數(shù)值模式，并利用數(shù)值實驗討論了Nudging松弛項的差分方案及Nudging松弛系數(shù)對潮波模擬結(jié)果的影響。

1 數(shù)值模式

POM（Princeton Ocean Model）模式是美國普林斯頓大學(xué)Blumberg和Mellor自1977年建立起來的，經(jīng)過多次修改，已成為最常用的海洋數(shù)值模式。該模式水平方向采用正交曲線網(wǎng)格，網(wǎng)格配置采用Arakawa C網(wǎng)格，垂向上采用σ坐標(biāo)。

本文采用Nudging同化方法，同化方案如下：

模式原方程為

式中，ζ為水位；n為時刻；Δt為時間步長；L為連續(xù)方程中的其他項。

1）Nudging松弛項顯式差分方案

在模式原方程基礎(chǔ)上增加Nudging項，則：

式中，ζo，m為第m個觀測水位；α為松弛系數(shù)，該值愈大表明趨近觀測值的約束作用愈強(qiáng)；ωm為其權(quán)重，若dl≤DM，ωm＝exp（1－dl／D），若dl＞DM，ωm＝0（其中；D取1／6°；DM取1°；λ和φ分別為經(jīng)度和緯度）。

進(jìn)一步得：

2）Nudging松弛項隱式差分方案

在模式原方程基礎(chǔ)上增加Nudging項：

最終得：

本文以西北太平洋海域為研究區(qū)域，范圍為105°30′～141°E、14°～41°N。模式的水平分辨率為1／12°×1／12°（圖1），垂向分10層。水深采用ETOPO2（National Geophysical Data Center，NOAA）并利用當(dāng)?shù)睾D資料進(jìn)行修正。開邊界處采用水位開邊界條件進(jìn)行驅(qū)動，調(diào)和常數(shù)由TPXO 6.0［36］資料集插值得到。本研究采用正壓潮模擬，溫度和鹽度均取常數(shù)，模式從靜止?fàn)顟B(tài)開始計算，對M2分潮模擬40個周期，存儲最后一個周期結(jié)果，并對其進(jìn)行調(diào)和分析，用來研究西北太平洋海域潮汐分布及潮波傳播特征。

本文用到的觀測值采用218個沿岸和島嶼驗潮站資料和218個衛(wèi)星高度計提取得到的調(diào)和常數(shù)資料，其站位分布見圖2，其中153個沿岸和島嶼驗潮站資料和185個衛(wèi)星高度計提取得到的調(diào)和常數(shù)資料用在模式中，剩余的用于模擬結(jié)果的檢驗。

圖1 計算網(wǎng)格Fig.1 Model grid

圖2 觀測值站位分布圖Fig.2 Locations of observation stations

2 數(shù)值實驗及結(jié)果

針對Nudging松弛項顯式和隱式差分方案及Nudging松弛系數(shù)對潮波模擬結(jié)果的影響，設(shè)計了一系列的數(shù)值實驗（表1），實驗1不加Nudging松弛項，實驗2～7加Nudging松弛項且采用顯式差分方案，實驗8～15加Nudging松弛項且采用隱式差分方案。為了分析各試驗?zāi)M結(jié)果的優(yōu)劣，本文用振幅絕均差和遲角絕均差來表示計算值和觀測值之間的偏差，這種方法比較直觀。其計算公式為：

式中，H為振幅；g為遲角；下標(biāo)sim和obs分別代表計算值和觀測值；N為觀測值個數(shù)。

表1 實驗方案設(shè)置Table 1 Experiments on nudging schemes

但是我們知道，當(dāng)振幅較小時，計算和觀測的遲角都不太穩(wěn)定，同時振幅較小時的誤差對潮高計算誤差影響也較小，故式（7）中不考慮振幅的差別而將所有遲角差進(jìn)行平均并不十分合理。另一種方法不是很直觀，但更合理一些，這種方法取均方根偏差：

式中，a＝Hcosg；b＝Hsing；σ為計算值和觀測值之間的距離，它表征了計算值和觀測值的偏離程度，而它與觀測值變化性的相對偏離程度為相對偏差：

式中，r為線性回歸中的相關(guān)系數(shù)，r2表示計算值和觀測值的擬合程度。

表2給出了數(shù)值實驗?zāi)M結(jié)果與98個觀測值的偏差對比，從中可以看出：

1）潮波模式中加了Nudging松弛項（實驗2～15）的模擬結(jié)果要明顯優(yōu)于未加Nudging松弛項（實驗1）的模擬結(jié)果；

2）在松弛系數(shù)取值相等情況下，Nudging松弛項采用隱式差分方案的結(jié)果均優(yōu)于采用顯式差分方案的結(jié)果；

3）Nudging松弛項差分方案不論為顯式還是隱式，模擬結(jié)果的偏差都隨著松弛系數(shù)的增加而逐漸減小，當(dāng)松弛系數(shù)增加到最優(yōu)值（0.002）后，模擬結(jié)果的偏差都隨著松弛系數(shù)的增加而增加；

4）Nudging松弛項差分方案為顯式時，當(dāng)松弛系數(shù)太大時（＞0.005），會導(dǎo)致模式溢出；Nudging松弛項差分方案改為隱式后，能夠有效提高松弛系數(shù)的閾值。

表2 各實驗?zāi)M結(jié)果與觀測值的偏差對比Table 2 Comparison between numerical simulations and observations

圖3給出了實驗10模擬得到的M2分潮的同潮圖，與Fang等［21］相比，除了在南黃海330°和60°等遲角線走向略有不同外，傳播趨勢基本相同，較好地體現(xiàn)了西北太平洋海域的潮波傳播特征。M2分潮都有4個無潮點，其中2個位于渤海，分別位于黃河口附近和秦皇島外海；另外2個位于黃海，分別位于成山頭附近和海州灣外。最大振幅位于南黃海東部江華灣，接近3m。本文在渤海黃河口附近得到了一個完整的旋轉(zhuǎn)潮波系統(tǒng)，與王永剛等［37］給出的無潮點位置基本一致，而丘仲鋒和何宜軍［34］以及朱學(xué)明等［35］模擬得到的黃河口附近的無潮點基本都已退化到陸地上。

圖3 實驗10得到的M2分潮同潮圖Fig.3 Co－tidal chart of M2from Experiment 10

3 結(jié)語

本文基于POM模式，針對Nudging松弛項的差分方案及Nudging松弛系數(shù)進(jìn)行了數(shù)值試驗，結(jié)果表明：

1）引入Nudging松弛項的潮波模擬結(jié)果明顯優(yōu)于未加Nudging松弛項的模擬結(jié)果；

2）在松弛系數(shù)取值相等情況下，Nudging松弛項采用隱式差分方案的結(jié)果要優(yōu)于采用顯式差分方案的結(jié)果；

3）Nudging松弛項差分方案不論是顯式還是隱式，模擬結(jié)果的偏差都隨著松弛系數(shù)的增加而逐漸減小，當(dāng)松弛系數(shù)增加到最優(yōu)值后，模擬結(jié)果的偏差將隨著松弛系數(shù)的增加而加大；

4）Nudging松弛項差分方案為顯式時，當(dāng)松弛系數(shù)太大時，將導(dǎo)致模式溢出；而Nudging松弛項差分方案為隱式時，能夠有效增大松弛系數(shù)的閾值。

（References）：

［1］WYRTKI K.Physical oceanography of the Southeast Asian waters，Scientific Results of Marine Investigations of the South China Sea and the Gulf of Thailand 1959－1961［J］.NAGA Report，1961，2：195.

［2］FANG G H.Tide and tidal current charts for the marginal seas adjacent to China［J］.Chinese Journal of Oceanology and Limnology，1986，4（1）：1－16.

［3］SHEN Y J.Numerical computation of tides in the East China Sea［J］.Journal of Shandong College of Oceanology，1980，10（3）：26－35.沈育疆.東中國海潮汐數(shù)值計算［J］.山東海洋學(xué)院學(xué)報，1980，10（3）：26－35.

［4］YE A L，ROBINSON I S.Tidal dynamics in the South China Sea［J］.Geophysical Journal of the Royal Astronomical Society，1983，72：691－707.

［5］SHEN Y J，HU D M，MEI L M，et al.Numerical computation of the tides in the South China Sea［J］.Transactions of Oceanology and Limnology，1985，1：1－11.沈育疆，胡定明，梅麗明，等.南海潮汐數(shù)值計算［J］.海洋湖沼通報，1985，1：1－11.

［6］FANG G H，YANG J F.A two－dimensional numerical model of the tidal motions in the Bohai Sea［J］.Chinese Journal of Oceanology and Limnology，1985，3（2）：135－152.

［7］DING W L.Distribution of tides and tidal currents in the South China Sea［J］.Oceanologia et Limnologia Sinica，1986，17（6）：468－480.丁文蘭.南海潮汐和潮流的分布特征［J］.海洋與湖沼，1986，17（6）：468－480.

［8］DOU Z X，YANG L W，OZER J.Three－dimensional simulation of tidal currents in the Bohai Sea［J］.Acta Oceanologica Sinica，1993，15（5）：1－15.竇振興，楊連武，OZER J.渤海三維潮流數(shù)值模擬［J］.海洋學(xué)報，1993，15（5）：1－15.

［9］FANG G H，CAO D M，HUANG Q Z.Numerical simulation on tides and tidal currents in the South China Sea［J］.Acta Oceanologica Sinica，1994，16（4）：1－12.方國洪，曹德明，黃企洲.南海潮汐潮流的數(shù)值模擬［J］.海洋學(xué)報，1994，16（4）：1－12.

［10］ZHAO B R，F(xiàn)ANG G H，CAO D M.Numerical simulation on tides and tidal currents in the Bohai，Yellow and East China Seas［J］.Acta Oceanologica Sinica，1994，16（5）：1－10.趙保仁，方國洪，曹德明.渤、黃、東海潮汐潮流的數(shù)值模擬［J］.海洋學(xué)報，1994，16（5）：1－10.

［11］ZHANG Z H，WU H D.Numerical computation of tides and tidal currents in the Bohai Sea［J］.Marine Forecasts，1994，11（1）：48－54.張占海，吳輝碇.渤海潮汐和潮流數(shù)值計算［J］.海洋預(yù)報，1994，11（1）：48－54.

［12］YE A L，MEI L M.Numerical modelling of tidal waves in the Bohai Sea，the Huanghai Sea and the East China Sea［J］.Oceanologia et Limnologia Sinica，1995，26（1）：63－70.葉安樂，梅麗明.渤黃東海潮波數(shù)值模擬［J］.海洋與湖沼，1995，26（1）：63－70.

［13］WAN Z W，QIAO F L，YUAN Y L.Three－dimensional numerical modeling of tidal waves in the Bohai，Yellow and East China Seas［J］.Oceanologia et Limnologia Sinica，1998，29（6）：611－616.萬振文，喬方利，袁業(yè)立.渤、黃、東海三維潮波運動數(shù)值模擬［J］.海洋與湖沼，1998，29（6）：611－616.

［14］WANG K，F(xiàn)ANG G H，F(xiàn)ENG S Z.A 3－D numerical simulation of M2tides and tidal currents in the Bohai Sea，the Huanghai Sea and the East China Sea［J］.Acta Oceanologica Sinica，1999，21（4）：1－13.王凱，方國洪，馮士筰.渤海、黃海、東海 M2潮汐潮流的三維數(shù)值模擬［J］.海洋學(xué)報，1999，21（4）：1－13.

［15］FANG G H，KWOK Y K，YU K J，et al.Numerical simulation of principal tidal constituents in the South China Sea，Gulf of Tonkin and Gulf of Thailand［J］.Continental Shelf Research，1999，19：845－869.

［16］ZHU X M，BAO X W，SONG D H，et al.Numerical study on the tides and tidal currents in Bohai Sea，Yellow Sea and East China Sea［J］.Oceanologia et Limnologia Sinica，2012，43（6）：1103－1113.朱學(xué)明，鮑獻(xiàn)文，宋德海，等.渤、黃、東海潮汐、潮流的數(shù)值模擬與研究［J］.海洋與湖沼，2012，43（6）：1103－1113.

［17］MAZZEGA P，BERGE M.Ocean tides in the Asian semienclosed seas from TOPEX／POSEIDON［J］.Journal of Geophysical Research，1994，99（C12）：24867－24881.

［18］YANAGI T，TAKAO T，MORIMOTO A.Co－tidal and co－range charts in the South China Sea derived from satellite altimetry data［J］.La Mer，1997，35（3）：85－93.

［19］HU J Y，KAWAMURA H，HONG H S，et al.Tidal features in the China Seas and their adjacent sea areas as derived from TOPEX／POSEIDON altimeter data［J］.Chinese Journal of Oceanology and Limnology，2001，19（4）：293－305.

［20］LIU K X，MA J R，HAN G J，et al.Tidal harmonic analysis of TOPEX／POSEIDON data in North west Pacific by introducing difference－ratio relations［J］.Acta Oceanologica Sinica，2002，24（4）：1－10.劉克修，馬繼瑞，韓桂軍，等.引入差比關(guān)系法分析西北太平洋TOPEX／POSEIDON衛(wèi)星高度計測高數(shù)據(jù)［J］.海洋學(xué)報，2002，24（4）：1－10.

［21］FANG G H，WANG Y G，WEI Z X，et al.Empirical cotidal charts of the Bohai，Yellow，and East China Seas from 10years of TOPEX／Poseidon altimetry［J］.Journal of Geophysical Research，2004，109（C11006），doi：10.1029／2004JC002484.

［22］LI P L，ZUO J C，LI L，et al.Orthogonalized convolution method for analysis of South China Sea tidal data from TOPEX／POSEIDON［J］.Oceanologia et Limnologia Sinica，2002，33（3）：287－295.李培良，左軍成，李磊，等.南海 TOPEX／POSEIDON高度計資料的正交響應(yīng)法潮汐分析［J］.海洋與湖沼，2002，33（3）：287－295.

［23］FANG G H.Finite difference of tidal equations——least square method and its application of simulation on M2tide in the Yellow Sea［J］.Science in China：Series B，1985，28（4）：356－364.方國洪.潮波方程的有限差分——最小二乘方法及其對模擬黃海 M2潮的應(yīng)用［J］.中國科學(xué)：B輯，1985，28（4）：356－364.

［24］HAN G J，HE B R，MA J R，et al.Optimizing open boundary conditions of nonlinear tidal model using adjoint methodⅠ：The establishment of adjoint model and twin－experiment［J］.Acta Oceanologica Sinica，2000，22（6）：27－33.韓桂軍，何柏榮，馬繼瑞，等.利用伴隨法優(yōu)化非線性潮汐模型的開邊界條件Ⅰ：伴隨方程的建立及“孿生”數(shù)值試驗［J］.海洋學(xué)報，2000，22（6）：27－33.

［25］HAN G J，F(xiàn)ANG G H，MA J R，et al.Optimizing open boundary conditions of nonlinear tidal model using adjoint methodⅡ：Assimilation experiment for tide in the Huanghai Sea and the East China Sea［J］.Acta Oceanologica Sinica，2001，23（2）：25－31.韓桂軍，方國洪，馬繼瑞，等.利用伴隨法優(yōu)化非線性潮汐模型的開邊界條件Ⅱ：黃海、東海潮汐資料的同化試驗［J］.海洋學(xué)報，2001，23（2）：25－31.

［26］Lü X Q，F(xiàn)ANG G H.Numerical experiments of the adjoint model for M2tide in the Bohai Sea［J］.Acta Oceanologica Sinica，2002，24（1）：17－24.呂咸青，方國洪.渤海 M2分潮的伴隨模式數(shù)值試驗［J］.海洋學(xué)報，2002，24（1）：17－24.

［27］LüX Q，F(xiàn)ANG G H.Inversion of the tides on the open boundary of the Bohai Sea by adjoint method［J］.Oceanologia et Limnologia Sinica，2002，33（2）：113－120.呂咸青，方國洪.渤海開邊界潮汐的伴隨法反演［J］.海洋與湖沼，2002，33（2）：113－120.

［28］WU Z K，TIAN J W，LüX Q，et al.A numerical model of tides in the South China Sea by adjoint method［J］.Oceanologia et Limnologia Sinica，2003，34（1）：101－108.吳自庫，田紀(jì)偉，呂咸青，等.南海潮汐的伴隨同化數(shù)值模擬［J］.海洋與湖沼，2003，34（1）：101－108.

［29］ZHANG J C，LüX Q.Inversion of the spacial varying bottom friction coefficient with the adjoint method［J］.Progress in Natural Science，2005，15（9）：1086－1093.張繼才，呂咸青.空間分布底摩擦系數(shù)的伴隨法反演研究［J］.自然科學(xué)進(jìn)展，2005，15（9）：1086－1093.

［30］ZHANG J C，LüX Q.Inversion of the bottom friction coefficient in 2Dtidal model of the Bohai，the Yellow and East China Seas［J］.Chinese Journal of Computational Mechanics，2007，24（4）：430－435.張繼才，呂咸青.渤、黃、東海二維潮汐模式底摩擦系數(shù)的反演［J］.計算力學(xué)學(xué)報，2007，24（4）：430－435.

［31］QIU Z F，HE Y J，LüX Q.Tidal adjoint assimilation with the TOPEX／Poseidon altimeter data in the Huanghai and Bohai Seas［J］.Acta Oceanologica Sinica，2005，27（4）：10－18.丘仲鋒，何宜軍，呂咸青.黃海、渤海TOPEX／Poseidon高度計資料潮汐伴隨同化［J］.海洋學(xué)報，2005，27（4）：10－18.

［32］ZU T T，GAN J P，EROFEEVA S Y.Numerical study of the tide and tidal dynamics in the South China Sea［J］.Deep－Sea Research I，2008，55：137－154.

［33］WANG Y G，F(xiàn)ANG G H，CAO D M，et al.Tides of the Bohai，Yellow and East China Seas by assimilating gauging station data into a hydrodynamic model［J］.Advances in Marine Science，2004，22（3）：253－274.王永剛，方國洪，曹德明，等.渤、黃、東海潮汐的一種驗潮站資料同化數(shù)值模式［J］.海洋科學(xué)進(jìn)展，2004，22（3）：253－274.

［34］QIU Z F，HE Y J.Assimilation of the M2tides in China Seas with TOPEX／POSEIDON data［J］.Oceanologia et Limnologia Sinica，2005，36（4）：376－384.丘仲鋒，何宜軍.TOPEX／POSEIDON高度計資料用于中國近海 M2分潮的同化反演［J］.海洋與湖沼，2005，36（4）：376－384.

［35］ZHU X M，BAO X W，HUANG B G，et al.An tidal numerical assimilation modeling in northwest Pacific［J］.Acta Oceanologica Sinica，2012，34（2）：1－12.朱學(xué)明，鮑獻(xiàn)文，黃必桂，等.西北太平洋的一種潮汐數(shù)值同化模型［J］.海洋學(xué)報，2012，34（2）：1－12.

［36］EGBERT G D，BENNETT A F，F(xiàn)OREMAN M G G.TOPEX／POSEIDON tides estimated using aglobal inverse model［J］.Journal of Geophysical Research，1994，99（C12）：24821－24852.

［37］WANG Y G，WEI Z X，F(xiàn)ANG G H，et al.A numerical study on the effect of changes in water depth and coastline on M2tidal component near the Yellow River Estuary［J］.Advances in Marine Science，2014，32（2）：141－147.王永剛，魏澤勛，方國洪，等.黃河口及其鄰近海域水深和岸線變化對 M2分潮影響的數(shù)值研究［J］.海洋科學(xué)進(jìn)展，2014，32（2）：141－147.