EFDC模型在三白蕩二維水流數(shù)值模擬中的應(yīng)用

孫麗娜， 李一平， 張其成， 陳 星

(1.南京市水資源管理中心， 江蘇 南京 210008；2.河海大學(xué) 環(huán)境學(xué)院， 江蘇 南京 210098；3.河海大學(xué) 水文水資源學(xué)院， 江蘇 南京 210098)

在城市湖泊的設(shè)計(jì)中，岸線的設(shè)計(jì)不僅影響到外界環(huán)境包括人類活動(dòng)，污染源對(duì)湖泊的擾動(dòng)程度，而且直接影響到湖體的水動(dòng)力過程，進(jìn)而影響污染物質(zhì)在湖體內(nèi)的遷移和擴(kuò)散。故了解不同岸線形狀下的湖泊水動(dòng)力過程是湖泊規(guī)劃設(shè)計(jì)的基礎(chǔ)和前提。而湖泊水動(dòng)力數(shù)值模擬則提供了一個(gè)便利的工具，本研究以吳江市三白蕩為具體的研究實(shí)例，利用EFDC模型對(duì)三白蕩現(xiàn)狀岸線條件下的水流流場進(jìn)行模擬，探尋合理的湖泊岸線設(shè)計(jì)方案，以有效減緩湖泊與陸地環(huán)境過于劇烈的物質(zhì)和能量交換，為岸線的生態(tài)設(shè)計(jì)提供技術(shù)參考。

1 現(xiàn)狀概況

三白蕩位于太湖流域下游吳江市汾湖經(jīng)濟(jì)開發(fā)區(qū)的中心城區(qū)，是典型的城市型湖泊?，F(xiàn)狀湖泊面積6.68 km2，岸線總長19 km，湖底地形較為平坦，湖底平均高程0.07 m(吳淞高程)，平均水深僅2.81 m，屬淺水型湖泊, 全年主導(dǎo)風(fēng)向?yàn)闁|南風(fēng)，夏季主導(dǎo)風(fēng)向?yàn)闁|南風(fēng)，冬季主導(dǎo)風(fēng)向?yàn)槲鞅憋L(fēng)，多年平均風(fēng)速為3.9 m/s，水面陣風(fēng)最大達(dá)到8～9級(jí)，風(fēng)力較大，對(duì)河道湖泊岸線侵蝕作用明顯。長期以來，受城市建設(shè)及歷史圍墾影響，岸線多被蠶食[1]，形態(tài)較為曲折，對(duì)外界環(huán)境的作用響應(yīng)相對(duì)強(qiáng)烈，與外界物質(zhì)與能量交換頻繁[2]。

2 EFDC模型簡介

EFDC (environmental fluid dynamic code)是在美國國家環(huán)保署資助下由維吉尼亞海洋研究所根據(jù)多個(gè)數(shù)學(xué)模型集成開發(fā)研制的綜合模型，被用于模擬水系統(tǒng)一維、二維和三維流場、物質(zhì)輸運(yùn)(包括溫度、鹽度和泥沙的輸運(yùn))、生態(tài)過程以及淡水入流等[3-4]。EFDC 模型垂向上采用σ坐標(biāo)變換，能較好地?cái)M和近岸復(fù)雜的岸線和地形；采用Gelperin 等[5]修正的 Mellor-Yamada 2. 5 階湍封閉模式較客觀地提供垂向混合系數(shù)，避免其人為選取造成的誤差。

本研究利用EFDC建立三白蕩的二維水動(dòng)力模型，對(duì)現(xiàn)狀岸線下湖泊的流場進(jìn)行模擬，研究結(jié)果可以為湖泊岸線的設(shè)計(jì)、改造提供基礎(chǔ)支撐，并為同類湖泊提供參考。

模型中使用的動(dòng)量方程、連續(xù)方程及狀態(tài)方程[6-8]如下：

?t(mxmyHu)+?x(myHuu)+?y(mxHvu)+

?z(mxmywu)-femxmyHv=-myH?x(p+patm+φ)+

(1)

?t(mxmyHv)+?x(myHuv)+?y(mxHvv)+

?z(mxmywv)+femxmyHu=-mxH?y(p+patm+φ)+

(2)

mxmyfe=mxmyf-u?ymx+v?xmy

(3)

(τxz,τyz)=AvH-1?z(u,v)

(4)

?t(mxmyH)+?x(myHu)+?y(mxHv)+

?z(mxmyw)=QH+δ(0)(QSS+QSW)

(5)

?zp=-gH(ρ-ρ0)ρ0-1=-gHb,ρ=ρ(P,S,T)

(6)

式中：u、v、w分別為邊界擬合正交曲線坐標(biāo)x、y、z方向上的速度分量；mx、my分別為度量張量對(duì)角元素的平方根，m=mxmy是度量張量行列式的平方根；Av為垂向紊動(dòng)黏滯系數(shù)；f為科里奧利系數(shù)；p為壓力；ρ為密度；S為鹽度；T為溫度；Qu、Qv為動(dòng)量源匯項(xiàng)。本研究中鹽度S為0，并假設(shè)水為不可壓縮流體，密度ρ和水溫T為常量。

動(dòng)量方程的垂向邊界層考慮了水表面的風(fēng)拖曳力和湖底的摩擦力。即河床剪切應(yīng)力τxz和τyz取決于速度分量，可以根據(jù)二次阻力公式求得，具體方程為：

(7)

式中，Uw、Vw分別為x、y方向在水表面10 m高處的分速。風(fēng)拖曳系數(shù)求法如下：

(8)

式中，ρa(bǔ)、ρw分別為空氣和水的密度。湖底摩擦力的計(jì)算方法為：

(9)

其中，下標(biāo)1指底部的對(duì)應(yīng)流速。底摩擦系數(shù)求取方程為：

(10)

式中：k為卡門常數(shù)；Δ1底層的無量綱厚度；zo為無量綱的糙率高度，一般取值在 0.002～0.01之間。

式(1)～(10)的數(shù)值解法采用有限差分格式，空間上采用二階精度的中心差分格式，時(shí)間上采用二階三次精度的內(nèi)外模分離模式有限差分格式。使用質(zhì)量守恒格式的干濕網(wǎng)格處理方法，更好地模擬湖流的漫灘過程，保證淺水區(qū)域計(jì)算結(jié)果的準(zhǔn)確性。首先設(shè)置干濕臨界水深，當(dāng)水深小于臨界值，則認(rèn)為處于干點(diǎn)；同時(shí)設(shè)置判定開關(guān)，用于定時(shí)檢查干濕狀態(tài)。對(duì)于干點(diǎn)，則人為將單元通量設(shè)為0，濕點(diǎn)則正常運(yùn)算。

3 模型構(gòu)建

現(xiàn)狀三白蕩面積為6.68 km2，其中島嶼面積0.9 km2，岸線長度為19 km。首先對(duì)該湖泊的邊界及出入湖河道進(jìn)行概化，然后利用湖泊岸線及島嶼邊界生成網(wǎng)格，共有網(wǎng)格11 302個(gè)。本研究水平網(wǎng)格采用矩形正交坐標(biāo)，網(wǎng)格距為22.8 m，垂直方向網(wǎng)格距為23.1 m。研究區(qū)域網(wǎng)格湖底高程及主要環(huán)湖河道(作為流量邊界條件)示意圖見圖1。

圖1 三白蕩網(wǎng)格湖底高程

4 模型參數(shù)設(shè)計(jì)

EFDC具有很好的通用性，在水動(dòng)力領(lǐng)域已具備相當(dāng)強(qiáng)的模擬精度。一般情況下，EFDC 模型中的許多參數(shù)不需要修改。譬如Mellor-Yamada湍封閉參數(shù)在各個(gè)模型中是相同的。以下為主要調(diào)整的參數(shù)。

(1)湖底粗糙度Z0

湖底粗糙度Z0是EFDC模型中常需調(diào)整的參數(shù)，其默認(rèn)設(shè)置為0. 02 m。在本研究的湖流模型中，設(shè)置Z0> 0. 03 m時(shí)，模型無法運(yùn)行；設(shè)置Z0在0. 03～0. 02 m時(shí)，模擬結(jié)果相對(duì)誤差在15%左右； 設(shè)置Z0在0.01～0. 02 m時(shí)，模型結(jié)果相對(duì)誤差在6%左右；Z0為0.015 m時(shí)，模擬結(jié)果和實(shí)測(cè)值符合較好。故本模型中Z0取為0.015 m。

(2)動(dòng)邊界干濕水深設(shè)定

固定邊界模型的計(jì)算域邊界隨時(shí)間不發(fā)生變化，而動(dòng)邊界模型的計(jì)算域邊界隨水位漲落而變動(dòng)，可以模擬湖泊水位的變化過程?？梢妱?dòng)邊界模型能詳細(xì)的模擬三白蕩水位的漲落過程，對(duì)后續(xù)研究近岸湖體污染物輸運(yùn)具有重要的應(yīng)用價(jià)值。

根據(jù)比較，選擇0.1 m作為干濕網(wǎng)格的臨界水深。即當(dāng)某網(wǎng)格水深大于0. 1 m時(shí)，當(dāng)作濕網(wǎng)格處理，進(jìn)行正常的模擬計(jì)算；當(dāng)水深小于0. 1 m時(shí)，此網(wǎng)格作為干網(wǎng)格處理，不參與計(jì)算。

(3)其他參數(shù)

其他參數(shù)如時(shí)間步長、水平黏性系數(shù)、水平擴(kuò)散系數(shù)、最小黏性系數(shù)等見表1。

表1 三白蕩水動(dòng)力模型參數(shù)取值表

5 模型應(yīng)用

(1)預(yù)測(cè)方案

模型預(yù)測(cè)條件為三白蕩在規(guī)劃岸線下在不同風(fēng)向及風(fēng)速下的流場模式，具體預(yù)測(cè)方案分為6種：

方案一：東南風(fēng)，風(fēng)速3.9 m/s，現(xiàn)狀出入湖流量；方案二：西北風(fēng)，風(fēng)速3.9 m/s，現(xiàn)狀出入湖流量；方案三：東南風(fēng)，風(fēng)速10 m/s，現(xiàn)狀出入湖流量；方案四：西北風(fēng)，風(fēng)速10 m/s，現(xiàn)狀出入湖流量；方案五：無風(fēng)條件，現(xiàn)狀出入湖流量；方案六：東南風(fēng)，風(fēng)速10 m/s，無出入湖流量。

現(xiàn)狀主要入湖河道及流量：牛長涇： 31.5 m3/s；楊蕩港： 39.6 m3/s；莘塔港： 20.2 m3/s?，F(xiàn)狀主要出湖河道及流量：南莊蕩： 61.2 m3/s；北窯港： 31.5 m3/s。

(2)模型計(jì)算結(jié)果

通過對(duì)不同計(jì)算結(jié)果的分析比較，可知模型在8 d后能夠計(jì)算穩(wěn)定。圖2為三白蕩現(xiàn)狀岸線條件下6種方案的流場模擬圖。

(3)模擬成果分析

吳江市常年主導(dǎo)風(fēng)向?yàn)闁|南風(fēng)和西北風(fēng)，多年平均風(fēng)速為3.9 m/s，湖面陣風(fēng)最大達(dá)到8～9級(jí)。由方案1與方案2可以看出，在3.9 m/s的東南風(fēng)和西北風(fēng)作用下，在西北部和西南處形成微小環(huán)流，但與湖泊外界無水體交換現(xiàn)象。由方案3與方案4可以看出，在10 m/s的東南風(fēng)和西北風(fēng)作用下，湖泊西北處、西南處受風(fēng)浪的沖刷強(qiáng)烈，形成明顯的環(huán)流區(qū)，雖在局部水體交換強(qiáng)烈，但與其他水體聯(lián)系較弱，污染物進(jìn)入該區(qū)域，易發(fā)生富營養(yǎng)化現(xiàn)象。由方案5與方案6可以看出，較出入湖水體而言，風(fēng)對(duì)三白蕩的流場起主導(dǎo)作用，且方案5中在無風(fēng)情況下水流流場較弱，在西北區(qū)基本無水體流動(dòng)，水體聯(lián)系較弱，為死水區(qū)。

圖2 三白蕩現(xiàn)狀岸線條件下6種方案的流場模擬

6 湖泊形態(tài)調(diào)整后流場模擬及結(jié)果分析

根據(jù)現(xiàn)狀流場模擬的結(jié)果分析后，為改善三白蕩岸線沖刷及內(nèi)外水體聯(lián)系，規(guī)劃調(diào)整三白蕩的湖泊岸線形態(tài)；由于三白蕩現(xiàn)狀平均淤積深度為0.14 m，為保障三白蕩及其周邊的防汛安全，增加湖泊調(diào)蓄能力，規(guī)劃對(duì)三白蕩進(jìn)行生態(tài)清淤，調(diào)整湖盆形態(tài)，局部地區(qū)疏浚深度0.5 m，見圖3。模型參數(shù)及方案設(shè)計(jì)同現(xiàn)狀形態(tài)，規(guī)劃岸線流場模擬見圖4，分析結(jié)果如下：

圖3 三白蕩現(xiàn)狀與規(guī)劃岸線形態(tài)示意圖

圖4 三白蕩規(guī)劃形態(tài)條件下6種方案的流場模擬

(1)由方案1及方案2可看出，在平均風(fēng)速為3.9 m/s，有出入湖流量情況下，湖泊的流場較為穩(wěn)定和均勻。最為明顯的是岸線調(diào)整過后，三白蕩西北部形成的環(huán)流與外部的水體產(chǎn)生水流聯(lián)系。

(2)在風(fēng)速加大后，如方案3和方案4所示，三白蕩的流場強(qiáng)弱明顯。在三白蕩的入湖、出湖河道的河口處，流速均較大。由于湖泊內(nèi)湖水流聯(lián)系得到整體改善，在風(fēng)速較大的情況下，對(duì)岸線容易形成沖刷，規(guī)劃在生態(tài)護(hù)岸布設(shè)過程中，需要對(duì)這些岸線進(jìn)行加固和保護(hù)。

(3)由方案3、4、6可以看出，湖泊在較大的風(fēng)速影響下，湖內(nèi)的整體水體水流聯(lián)系較強(qiáng)。但是在岸線較為曲折的地區(qū)，分布著一些緩流區(qū)，這些緩流區(qū)在風(fēng)速較小的情況下易形成死水區(qū)，規(guī)劃中可通過生態(tài)工程的布設(shè)來改善本區(qū)內(nèi)的水生態(tài)環(huán)境，加強(qiáng)對(duì)湖體重污染物的稀釋和降解。

7 結(jié) 語

湖泊形態(tài)的設(shè)計(jì)必須從整體出發(fā)，通過全面的分析比較，最終選取最好的規(guī)劃設(shè)計(jì)方案。 一方面，三白蕩在經(jīng)過水流數(shù)值模擬后，可清晰的看到湖泊中的死水區(qū)，為研究三白蕩及太湖流域中類似淺水湖泊的水動(dòng)力過程及增強(qiáng)湖泊環(huán)境容量提供較強(qiáng)的理論支撐。另一方面，模擬結(jié)果顯示，風(fēng)生流對(duì)三白蕩的流場起主導(dǎo)作用，根據(jù)流場圖顯示的岸線受侵蝕嚴(yán)重區(qū)域的分布情況，可通過加強(qiáng)生態(tài)護(hù)岸等措施加強(qiáng)對(duì)湖泊岸線的保護(hù)。此外，水流數(shù)值模擬不僅能夠較好地模擬湖內(nèi)水流聯(lián)系、污染物遷移轉(zhuǎn)化情況，而且還可以較好地模擬湖泊的出入湖水量分配情況，對(duì)于開展湖泊水資源保護(hù)工作具有較強(qiáng)的指導(dǎo)意義。但是實(shí)際的應(yīng)用過程中，往往由于缺少入湖的徑流資料、水位變動(dòng)資料，水流模擬的仿真性值得進(jìn)一步研究，所以在實(shí)際岸線湖盆形態(tài)規(guī)劃設(shè)計(jì)中，需要水文、氣象、城市規(guī)劃等多部門間的合作。