王 文 才，杜 薇，范 中 亞，姜 龍，許 益 新，李 一 平

(1.河海大學 淺水湖泊綜合治理與資源開發(fā)教育部重點實驗室，江蘇 南京 210098； 2.生態(tài)環(huán)境部 華南環(huán)境科學研究所，廣東 廣州 510655； 3.生態(tài)環(huán)境部 南京環(huán)境科學研究所，江蘇 南京 210042)

水流輸運的快慢及水體交換速率對于物質輸移、水生生物生長、水體自凈能力均有重要的意義[1-2]。目前研究水體輸運時間尺度的內容主要包括沖刷時間、滯留時間、轉化時間、通過時間和水齡等[3];研究手段主要包括水平衡[4]、示蹤劑、拉格朗日粒子漂浮、放射性同位素[5]及數(shù)值模擬等方法[6-7]。其中水平衡方法雖然簡單，但是不能很好地體現(xiàn)出空間上的水體動力情況；通過向水體投加示蹤劑或者通過測定水體中放射性同位素的變化來得到水體滯留時間，也存在一定的空間和時間局限性。因此目前更多的是采用數(shù)值模擬的方式來研究河流[7]、湖泊水庫[8-11]和河口海灣等各種特定地表水體的水流輸運及水力滯留[6，12]。長江作為世界第三大河流，全長6 300多km，沿江經(jīng)濟相對發(fā)達[13]。下游長江三角洲更是中國經(jīng)濟最發(fā)達的地區(qū)之一，人口密度和工業(yè)化程度也相對較高，沿江大量生活污水和生產(chǎn)廢水排放進入長江水體，造成水體中營養(yǎng)物質含量偏高，長江入?？诩敖Ｋw富營養(yǎng)化程度較高[14]，同時，沿江的化工廠也存在風險事故，存在危險化學品泄漏進入長江的可能。江蘇地區(qū)作為下游相對發(fā)達的地區(qū)，大部分的污廢水也都排放進入了長江。因此，研究長江江蘇段的水流輸運特點對于明確污染物輸運及進行河段管理都具有重要的意義。本文基于EFDC模型中水齡和拉格朗日粒子追蹤模塊，采用設計流量和下邊界潮位條件，研究長江江蘇段的水流輸運和水力滯留特點，在一定程度上闡明了長江江蘇段的水流輸運規(guī)律。

1 研究區(qū)域與研究方法

1.1 研究區(qū)域

本文以長江江蘇段(主要是南京至徐六涇江段)為研究對象(圖1)。長江下游大通以下為感潮江段，水位流速等水動力均受到長江入?？诔辈ㄏ蛏蟼鞑サ挠绊慬15]，水動力過程較復雜。大通流量通常作為下游及入海徑流量的依據(jù)[16-17]，長江南京至徐六涇江段位于長江下游，分布著南京、鎮(zhèn)江、揚州、常州、無錫、蘇州這6座城市，沿江包括南京、鎮(zhèn)江、江陰、天生港和徐六徑潮位站。

圖1 研究區(qū)域 Fig.1 Study area

1.2 研究方法

1.2.1 模型建立

本研究選取可以用于模擬包括河流、湖泊、水庫、河口、濕地、近岸等地表水體三維水流、物質輸運和生物化學過程的模型系統(tǒng)EFDC。在用其構建長江南京至徐六涇江段的水動力模型的基礎上，計算設計徑流量和下游潮位下的江段水齡變化及拉格朗日粒子輸移過程。

1.2.2 水齡計算方法

水齡的概念最早由Zimmerman提出[18]，定義為可溶性物質從入口傳輸?shù)街付ㄎ恢盟枰臅r間(通常將入口處的水齡設置為0)，可以反映水體的交換狀況，水齡越大水體交換程度越弱。EFDC模型中水齡根據(jù)示蹤劑來計算，水齡方程如下：

(1)

(2)

式中，c是示蹤劑濃度，α是初始水齡濃度，u是時空分布的流速，K是擴散張量，t為時間。水齡α的計算如下

(3)

1.2.3 拉格朗日粒子追蹤計算方法

拉格朗日粒子示蹤劑模型(LPT)，可以模擬粒子在水體中的輸移路徑。該模型考慮了質量輸運的對流擴散，三維曲線正交坐標系的質量運輸對流擴散方程是：

(4)

顆粒運動的微分方程與方程(4)一致，如下所示：

(5)

(6)

(7)

式中，dt是時間單元，p是一個由分布相同均值為0.5的隨機變量產(chǎn)生的隨機值。當轉換使用2p-1時，這個隨機部分的平均值為0并在-1～1范圍內。變化的隨機值允許擴散期限到更多顆粒+/-在平流運輸方位。

為了判定顆粒的運動軌跡，方程都被納入EFDC模型中。如上所述，數(shù)值求解方法分別被分為平流運輸和隨機成分。這個方法允許使用者開啟和關閉水平或垂直方向的隨機游走功能。微分方程可采用顯式歐拉、預測校正歐拉和四階龍格庫塔方法來解。

1.3 研究方案

長江江蘇段已開展過較多模型模擬應用，能夠較好地模擬江段水動力狀況。鑒于篇幅有限，此處不再贅述模型率定、驗證結果。由于三峽工程的建設，下游大通水文站的極值流量發(fā)生了改變，洪峰流量降低，枯水流量變大?？紤]三峽工程建設后2003～2013年最小流量和最大流量的范圍，本研究江段模型上邊界考慮10 000，20 000，30 000，40 000，50 000 m3/s和60 000 m3/s共6種設計流量；下邊界徐六涇的潮位過程根據(jù)水文年鑒中逐日高、低潮潮位數(shù)據(jù)，采用cubic插值得到逐時的潮位，之后通過調和分析得到調和常數(shù)，最后得到設計潮位過程。采用EFDC模型進行模擬過程中，拉格朗日粒子釋放時間設定在模擬開始后的0.5 d，在接近入口處隨機釋放20個粒子，通過模型計算得到各種設計情景下的江段水齡變化和拉格朗日粒子輸移過程。

2 結果分析

2.1 不同徑流量下長江江蘇段水齡模擬

統(tǒng)一截取水齡達計算基本穩(wěn)定、計算時間為第20 d的計算結果。不同徑流量下模擬結果在空間上基本呈現(xiàn)相同的特征，因此此處只給出徑流量為40 000 m3/s的結果(圖2)，水齡值從上游向下游逐漸增大，入口處南京站所在位置的水齡為0。流量越大，水齡值的變化范圍就越小。上游水齡在空間上分布相對均勻，下游接近出口徐六涇處水齡在空間上存在明顯的變化，河道斷面中間水齡小，靠近兩岸的水齡大；另外，下游出口附近，受到下邊界水流反向進入的影響，水齡的結果呈現(xiàn)出隨時間波動較大的特征，甚至會存在水齡接近于零的區(qū)域。

圖2 流量為40 000 m3/s情況下江段水齡空間分布 Fig.2 Distribution of water age under the discharge of 40 000 m3/s

根據(jù)水齡計算結果，匯總在不同設計徑流量下南京、鎮(zhèn)江、江陰、天生港、徐六涇這5個代表性站點的水齡計算值見表1。

表1 設計流量下各站點水齡計算值

Tab.1 Water age of each station under designed discharge d

站點入口流量 10000 m3/s20000 m3/s30000 m3/s40000 m3/s50000 m3/s60000 m3/s南京000000鎮(zhèn)江4.062.031.431.120.940.82江陰9.634.673.202.442.091.81天生港12.146.184.263.232.702.32徐六涇9.408.165.114.063.302.87

由于徐六涇的水齡值不穩(wěn)定，以天生港為最下游分析，可以看到流量從10 000 m3/s增加到60 000 m3/s的情況下，天生港的水齡值從12.14 d減小到2.32 d，下降約80%，其它站點水齡值也存在大幅度下降趨勢，這也就說明，豐水期和枯水期江段的水流輸運時間存在較大的變化。

2.2 拉格朗日粒子輸移結果

由于粒子是隨機釋放的，受河道斷面上流速不均勻的影響，拉格朗日粒子的輸移速率存在差異。本文中定義10%的粒子到達或通過某點為水流到達，90%的粒子通過或到達某點定義為水流完全通過。結合拉格朗日粒子的計算結果，以粒子到達和完全通過南京站的時間為第0 d，分別提取粒子到達和通過鎮(zhèn)江、江陰、天生港及徐六涇的時間(表2)。從結果可以看出，粒子從到達各站點到完全通過需要一定的時間，越往下游所需要的時間就越長。流量從10 000 m3/s增大到60 000 m3/s，粒子從南京輸移到鎮(zhèn)江、江陰、天生港、徐六涇的時間范圍分別為3.97～0.63，8.63～1.54，11.54～2.08 d和14.29～2.83 d。流量越小，輸移所需要的時間越長；越往下游，流量的影響越大。粒子完全通過的時間也呈現(xiàn)出相同的特征。與此類似的是，粒子從到達到完全通過某個站點，所需的時間也會隨流量的減小和離下游邊界距離減小而變大。

表2 不同流量下拉格朗日粒子從釋放至到達相應站點的輸移時間Tab.2 Transportation time of particle from relasee to arrival at corresponding site

3 討 論

為分析水齡的沿程變化，以上游設計流量為40 000 m3/s為例，提取結果相對穩(wěn)定的南京、鎮(zhèn)江、江陰、天生港4個站點的水齡值和相應站點與南京的距離進行二次擬合(圖3)，具有較高的擬合度。擬合函數(shù)的二次項系數(shù)為正，且從圖3中也可以看出來，越往下游，相同距離所需要的輸運時間也越長。河流中水流和粒子的輸移總體上是從上游向下游，因此水齡總是沿河流向下游不斷增大的。但是由于潮汐的原因，在長江下游潮流界內水流存在反向的特征，靠近下邊界的部分存在開邊界水向河段內倒灌的情況。因此越向下游的位置，水的輸運速度越慢，相同距離輸運所需要的時間也就越長。另外有關于長江口的水流輸運時間的研究也表明在多年平均流量條件下，水流從徐六涇輸運至河口大約需要24 d，說明下游段受到潮流的影響，水流輸運會越來越慢[19]。其它有關感潮河流、河口的研究也都呈現(xiàn)類似的變化規(guī)律[6-7]。

圖3 40 000 m3/s徑流條件下水齡的沿程變化Fig.3 Water age along the river under discharge of 40 000 m3/s

結合不同流量下各站點的水齡計算結果，采用冪函數(shù)對鎮(zhèn)江、江陰、天生港3個站點穩(wěn)定后的水齡和流量進行回歸分析(圖4)，R2均大于0.99，具有較好的擬合效果。隨著流量的增大，水齡逐漸減小，但減小的幅度在下降。這與相關研究中采用指數(shù)函數(shù)對流量和滯留時間進行回歸具有類似的效果[6,20]。在流量接近無窮大時，水齡接近0；而在流量接近0時，水齡接近無窮大，顯然采用冪函數(shù)進行回歸相比指數(shù)函數(shù)更符合實際。由于徐六涇位于模擬計算區(qū)域的下邊界，在徑流量較大情況下，隨著徑流量的增大，水齡的趨勢與其他幾個站點的趨勢基本一致；但在徑流量較小(Q=10 000m3/s)情況下，受下邊界影響較大，與其徑流量較大時的趨勢并不一致。

圖4 各站點水齡與徑流量的關系Fig.4 The relationship between water age and discharge at different stations

另外，拉格朗日粒子輸移計算的結果呈現(xiàn)出與水齡類似的分布特征。采用冪函數(shù)來對粒子到達或通過各站點的時間與流量進行回歸，同樣能夠表現(xiàn)出較好的擬合效果(見圖5)。拉格朗日粒子示蹤相比于水齡計算，在靠近下游開邊界處更加穩(wěn)定，不容易受到開邊界水倒灌的影響。其中粒子到達站點的時間與水齡計算結果接近，粒子完全通過站點的時間則要相對較大，可以反映出污染物由于河流斷面流速不均引起的粒子輸移動態(tài)變化過程。以粒子完全通過和到達的時間為界，可以反映出類似的污染物經(jīng)過相關站點所需要的時間。

圖5 設計流量條件下拉格朗日粒子到達和完全通過各站點的時間Fig.5 Time duration of particle arrived and passed the specific station under the designed discharge

4 結 論

本研究基于EFDC中的水齡和拉格朗日粒子示蹤模塊，以長江江蘇段(南京至徐六涇)為研究對象，模擬了不同設計流量及下邊界潮位條件下的江段水齡分布及拉格朗日粒子輸移時間的規(guī)律，研究結果表明：

(1) 受下游潮汐的影響，越靠近下游水體輸運越慢，水流輸運相同的距離所需要的時間也相對越長。

(2) 河段整體的水流輸運時間與徑流量呈現(xiàn)非線性關系，在流量相對較小時水流輸運受流量的影響較大。

(3) 拉格朗日粒子示蹤的方法相較于水齡計算的方法，在河段下游的水流輸運結果的穩(wěn)定性受開邊界水位影響更小。

(4) 水流沿程的輸運時間與徑流量可以采用二項式進行擬合描述。水流從入口輸運到達、完全通過特定站點以及該站點的水齡均可以采用冪函數(shù)進行回歸計算，均有較好的效果。

本研究可以為突發(fā)性污染事故的預警以及河段管理提供支撐。