汾河入黃口水流特性與污染物擴(kuò)散模擬

李 娟，馮民權(quán)，劉淑文

（1.西安理工大學(xué) 陜西省西北旱區(qū)生態(tài)水利工程重點(diǎn)實(shí)驗(yàn)室，陜西 西安710048；2.山西中部引黃水務(wù)開(kāi)發(fā)有限公司，山西 太原030012；3.山西水利職業(yè)技術(shù)學(xué)院 水利工程系，山西 運(yùn)城044004）

汾河是黃河的一級(jí)支流，汾河入黃口屬于干、支流交匯的河段。天然河流中干支流相匯處由于干、支流兩股水流相互頂托，交匯區(qū)內(nèi)水流紊動(dòng)摻混作用強(qiáng)烈，能量損失很大，可能形成河口淺灘或江心灘（洲）。因此，在干支流匯合口流段，無(wú)論是水流現(xiàn)象、污染物濃度分布都具有復(fù)雜多變的特點(diǎn)。馮鏡潔，李然［1］等采用k—ε雙方程模型和VOF模型模擬流場(chǎng)，自由面，研究了分離區(qū)的出現(xiàn)和幾何尺度。得出分離區(qū)出現(xiàn)在匯口下游主流右側(cè)，形狀隨著入?yún)R角的減小而變得狹長(zhǎng)，尺寸隨著入?yún)R角度和匯流比的增大而增大。李霞［2］等研究了當(dāng)入?yún)R角為60°時(shí)，匯流比對(duì)水位，流速以及泥沙淤積的影響。河口與海岸環(huán)境中的水體流動(dòng)具有明顯的空間三維特性，且三維模型不僅能反映真實(shí)的水流運(yùn)動(dòng)，而且可以從三維計(jì)算結(jié)果中提取某斷面或者某層的結(jié)果作為二維結(jié)果，因此很多學(xué)者對(duì)三維模型展開(kāi)了研究。劉明［3］用ADI三維潮流模型，模擬和預(yù)測(cè)河口處的流場(chǎng)，并建立污染物的輸運(yùn)模型，模擬污染物的現(xiàn)狀分布和預(yù)測(cè)污染物的變化。Isabel［4］等研究了不同流動(dòng)狀態(tài)對(duì)葡萄牙杜若河口污染物擴(kuò)散的影響。Ng［5］等提出了由GIS和三維水動(dòng)力，泥沙和重金屬運(yùn)輸?shù)募赡Ｐ?，?dòng)態(tài)演示水動(dòng)力變化和污染物擴(kuò)散過(guò)程，并很好地應(yīng)用于實(shí)例計(jì)算。對(duì)于水平尺度遠(yuǎn)大于垂向尺度的寬淺型河口，或者窄深型的河口可以將三維模型簡(jiǎn)化為二維模型。韓龍喜［6］采用非穩(wěn)態(tài)淺水方程和垂向平均的二維對(duì)流擴(kuò)散方程建立水動(dòng)力及水質(zhì)耦合模型；王昆［7］在詳細(xì)分析遼東灣沿岸各重點(diǎn)排污口排污現(xiàn)狀基礎(chǔ)上，建立了遼東灣海域污染物輸移的對(duì)流擴(kuò)散模型；徐明德［8］在矩形網(wǎng)格系統(tǒng)下建立二維水動(dòng)力水質(zhì)模型。萬(wàn)清華［9］對(duì)甌江河口區(qū)域建立了Delft 3維水動(dòng)力和WASPS水質(zhì)耦合模型。這些學(xué)者采用所建立的模型模擬了不同工況下的流速分布，并探討排污口尾水排放對(duì)河口污染物的遷移轉(zhuǎn)化的影響。侯建華［10］等采用矩形網(wǎng)格劃分地形，研究了50%，75%頻率條件下汾河入黃口的水流水質(zhì)特性。本研究將汾河看作源相匯入黃河，結(jié)合入黃口處的復(fù)雜地形環(huán)境，采用非結(jié)構(gòu)三角網(wǎng)格來(lái)劃分地形，更能適合形狀不規(guī)則的區(qū)域，并在侯建華［10］的基礎(chǔ)上增加90%頻率年的計(jì)算條件，更加全面地考慮不同工況下，入黃口處的水流特性和COD的分布情況。同時(shí)分析入黃口處不同匯流比對(duì)形成的河口灘面積的影響。

1 模型建立

1.1 基本方程的建立

模型的數(shù)學(xué)基礎(chǔ)是雷諾平均化的N—S方程，k—ε混合模型作為計(jì)算工具。

連續(xù)性方程：

動(dòng)量守恒方程：

在對(duì)水體水質(zhì)模擬時(shí)，在質(zhì)量及動(dòng)量方程的基礎(chǔ)上還需要加一個(gè)濃度的對(duì)流擴(kuò)散方程。形式為：

k—ε混合模型形式為：

式中：xi——坐標(biāo)系統(tǒng)；ui——xi方向的流速分量；t——時(shí)間；ρ——水的密度；p——壓力；vt——紊動(dòng)黏 性 系 數(shù)，；c——經(jīng) 驗(yàn) 常 數(shù)；δ ——uijKronecker函數(shù)；k——紊動(dòng)動(dòng)能；ε——紊動(dòng)能耗散率；gi——重力加速度在xi方向的分量；C——濃度（mg／L）；Ei——濃度混合擴(kuò)散系數(shù)（m2／s）；k1——污染物降解系數(shù)；cε1，cε2，σk，σε——經(jīng)驗(yàn)常數(shù)（表1）。

表1 特征參數(shù)值［10］

1.2 模型求解

模型采用交替方向隱式迭代法（ADI）對(duì)水流連續(xù)性及動(dòng)量守恒方程進(jìn)行積分。具體離散用半隱式，求解采用追趕法。

1.3 模型驗(yàn)證

采用2010年2月份黃河和汾河實(shí)測(cè)的流量數(shù)據(jù)，選取兩個(gè)斷面，并將模擬值和實(shí)測(cè)值比較來(lái)驗(yàn)證模型的可靠性。兩個(gè)斷面分別為入黃口上游斷面和入黃口下游斷面。兩斷面的高程、水位、流速變化的模擬結(jié)果如圖1—5所示。由圖1—5可知，入黃口上游斷面貫穿了黃河以及汾河，但此時(shí)黃河沒(méi)有江心洲，兩條河水位和流速最大的部分主要分布河道中心，局部偏左或偏右。入黃口下游斷面只貫穿黃河，由于江心洲的影響，水深和流速最大值分布在江心洲兩側(cè)。兩斷面水位模擬的平均絕對(duì)誤差分別為±0.014 6和±0.000 2m，相對(duì)誤差分別為±4.09×10－5，±8.9×10－6；兩斷面流速模擬的平均絕對(duì)誤差分別為±0.03和±0.009m，相對(duì)誤差分別為±9.4%和±3.6%。因此模擬結(jié)果與實(shí)測(cè)資料相吻合，故建立的二維水動(dòng)力模型是可靠的。

圖1 入黃口上游和下游兩斷面河底高程分布

圖2 入黃口上游斷面水位模擬結(jié)果

圖3 入黃口下游斷面水位模擬結(jié)果

圖4 入黃口上游斷面流速模擬結(jié)果

圖5 入黃口下游斷面流速模擬結(jié)果

2 模型計(jì)算條件

2.1 計(jì)算區(qū)域

計(jì)算區(qū)域的上游邊界設(shè)在黃淤斷面60處，下游邊界設(shè)在黃淤斷面60斷面與59斷面之間，在入黃口下游的4.5km處，計(jì)算區(qū)域面積約9.96×106m2。交匯口處黃河河寬與汾河河寬的比例約為3∶1，且河中心有江心洲，此時(shí)形成的平面區(qū)域不規(guī)則，采用三角網(wǎng)格能夠更好地貼合邊界。單個(gè)網(wǎng)格面積約為1 000m2，對(duì)網(wǎng)格進(jìn)行100次迭代光滑，以便達(dá)到較好的網(wǎng)格質(zhì)量。

2.2 設(shè)計(jì)條件確定

2.2.1 水動(dòng)力模型的初始條件和邊界條件確定 模型定解條件包括初始條件和邊界條件。對(duì)于初始條件，根據(jù)實(shí)際的高程值初始水位取值為365m，初始流速場(chǎng)采用“零啟動(dòng)”的方法給定，污染物濃度初始場(chǎng)采用全場(chǎng)均勻分布，并等于進(jìn)口處濃度值的方式給定。對(duì)于邊界條件，上游進(jìn)口給定流量，下游出口給定水位，將汾河設(shè)為源項(xiàng)匯入黃河。根據(jù)黃河龍門(mén)水文站和汾河下游河津水文站所提供的1988—2008年的實(shí)測(cè)資料來(lái)計(jì)算50%，75%，90%頻率年的流量、流速，并作為水文設(shè)計(jì)條件，此時(shí)黃河流量分別為503.6，378.5，303.4m3／s；汾河的流量分別為9.54，6.33，4.4m3／s。3種工況下，匯流比（（支流流量與主流流量比值）分別為0.018，0.016.0.014。

2.2.2 物質(zhì)輸運(yùn)模型的初始條件和邊界條件確定物質(zhì)運(yùn)移模型穩(wěn)定狀態(tài)下運(yùn)行，因此可不考慮初始值的設(shè)定。根據(jù)實(shí)測(cè)資料，2007年龍門(mén)斷面COD平均濃度為13.4mg／L，隨著COD的不斷衰減與降解，計(jì)算可得到入黃口處濃度約值為10mg／L，取10mg／L作為濃度背景值。按照山西省地表水環(huán)境功能區(qū)劃結(jié)果，汾河入黃處汾河水質(zhì)按Ⅴ類(lèi)水標(biāo)準(zhǔn)，即COD濃度取40mg／L。從2005—2010年月平均實(shí)測(cè)COD值可得，汾河水質(zhì)較好時(shí)，COD濃度能夠達(dá)到水功能區(qū)劃的要求。根據(jù)黃河水功能區(qū)劃，汾河入黃處的黃河的水質(zhì)目標(biāo)為Ⅲ類(lèi)，COD的濃度為20mg／L。以此來(lái)研究汾河入黃口以后混合區(qū)COD的輸移擴(kuò)散規(guī)律。通過(guò)對(duì)汾河運(yùn)城段20多個(gè)排污口進(jìn)行調(diào)查，給出3種工況下汾河排入黃河COD的量都為3 276.631t／a。

3 結(jié)果與分析

根據(jù)水文水質(zhì)資料，采用二維水動(dòng)力水質(zhì)模型分別對(duì)50%，75%，90%頻率條件下水流水質(zhì)進(jìn)行模擬?？傻闷淞鲃?dòng)特性為：支流受上游來(lái)流的擠壓和河岸固壁的限制，在交匯口下游附近形成分離區(qū)。（1）黃河與汾河交匯處入?yún)R角約30°，受汾河來(lái)水的影響，入?yún)R口水流總會(huì)趨于平緩，流速相對(duì)于主流其他河段較小，導(dǎo)致易在交匯口形成新的河口灘。（2）交匯口處所形成分離區(qū)，能夠影響泥沙沉積，分離區(qū)的尺寸與交匯形狀、入?yún)R角、干支流比和糙率等因素有關(guān)。水流模擬過(guò)程中入?yún)R角，糙率等因素不變，但是不同頻率年的干支流比卻不相同，因此，干支流比決定著河口灘面積的大小。隨著三種工況下匯流比的減小，河口灘的面積逐漸減小及50%頻率年形成河口灘大于75%頻率年所形成的河口灘，而90%水文頻率年下未形成新河口灘。（3）由于支流的頂托作用，3種工況下匯口上游主流發(fā)生壅水現(xiàn)象。

以流場(chǎng)分布特性為基礎(chǔ)分別對(duì)50%，75%，90%頻率年COD濃度分布進(jìn)行模擬，可得混合帶COD濃度沿程變化范圍和變化幅度。COD濃度分布規(guī)律大致為：（1）隨著水流的混摻運(yùn)動(dòng)，污染物質(zhì)不斷混合稀釋?zhuān)S河河中間的江心洲，明顯阻礙了COD的擴(kuò)散稀釋?zhuān)虼巳N頻率年條件只給出江心洲左側(cè)的COD分布情況。（2）不同頻率年的COD分布詳見(jiàn)表2。表2中縱向分布表示沿黃河的水流方向污染物從入黃處到達(dá)Ⅲ類(lèi)水的起止坐標(biāo)。隨著水文頻率的增大，在黃河縱向分別經(jīng)過(guò)2.35，3.00，4.20km的距離后水質(zhì)達(dá)標(biāo)。橫向分布表示污染物在黃河橫向擴(kuò)散的范圍，也由起止坐標(biāo)來(lái)表示，即隨著水文頻率的增大分別經(jīng)過(guò)約0.5，0.6，0.7km 后達(dá)標(biāo)。水文頻率變化對(duì)橫向分布影響較縱向小。

以上分析結(jié)果可知，COD遷移擴(kuò)散過(guò)程中在縱向的遷移距離遠(yuǎn)大于橫向距離，故而，COD總體分布在一個(gè)狹長(zhǎng)的區(qū)域內(nèi)。污染物在河流擴(kuò)散時(shí)，在橫向擴(kuò)散起主導(dǎo)作用，而在縱向既有對(duì)流又有擴(kuò)散，且對(duì)流作用遠(yuǎn)大于擴(kuò)散作用。模擬結(jié)果符合污染物遷移轉(zhuǎn)化規(guī)律。

表2 汾河入黃口不同頻率年COD分布情況

4 結(jié)論

（1）采用ADI法在非結(jié)構(gòu)三角網(wǎng)格系統(tǒng)求解連續(xù)方程及動(dòng)量方程，可以較好地適用不規(guī)則的天然河道，并且保證水量、動(dòng)量完全守恒。

（2）通過(guò)對(duì)二維水動(dòng)力驗(yàn)證，表明所采用的模型模擬效果良好，可以進(jìn)行流場(chǎng)模擬和污染分布計(jì)算。

（3）分別對(duì)3種頻率年進(jìn)行水流水質(zhì)模擬，能夠給出汾河水質(zhì)達(dá)標(biāo)后COD在黃河分別經(jīng)過(guò)約2.35，3.00，4.20km的距離后，水質(zhì)能夠達(dá)到黃河水功能區(qū)劃的要求。

［1］ 馮鏡潔，李然，王協(xié)康，等.河流交匯分離區(qū)特性研究［J］.水動(dòng)力學(xué)研究與進(jìn)展（A輯），2009，24（3）：320-325.

［2］ 李霞，侯極.入?yún)R角60°時(shí)入?yún)R口泥沙淤積試驗(yàn)研究［J］.水電能源科學(xué)，2011，29（8）：75-77.

［3］ 劉明.遼河口污染物擴(kuò)散數(shù)值模擬及總量控制研究［D］.遼寧 大連：大連海事大學(xué)，2006.

［4］ Isabel C A，Adriano A B，Pedro M D.Influence of river discharge patterns on the hydrodynamics and potential contaminant dispersion in the Douro estuary（Portugal）［J］.Water Research，2010，44（3）：3133-3146.

［5］ Sandy M Y N，Onyx W H W.Integration of a GIS and a complex three-dimensional hydrodynamic，sediment and heavy metal transport numerical model［J］.Advances in Engineering Software，2009，40（6）：391-401.

［6］ 韓龍喜，易路，劉軍英，等.連云港近岸海域污染物輸移規(guī)律［J］.河海大學(xué)學(xué)報(bào)：自然科學(xué)版，2011，39（3）：248-253.

［7］ 王昆，王年斌.遼東灣海域水體污染物輸運(yùn)的數(shù)值模擬［J］.水動(dòng)力學(xué)研究與進(jìn)展（A輯），2010，25（4）：493-498.

［8］ 徐明德.黃海南部近岸海域水動(dòng)力特性及污染物輸移擴(kuò)散規(guī)律研究［D］.上海：同濟(jì)大學(xué)，2006.

［9］ 萬(wàn)清華.甌江平面二維水量—水質(zhì)耦合模型研究及應(yīng)用［D］.江蘇 南京：河海大學(xué)，2007.

［10］ Xing Xiaopeng，Hou Jianhua，F(xiàn)eng Minquan，et al.Research on Flow and Pollutant Diffusion at the Fen River Entrance into Yellow River［C］∥Computer and Management（CAMAN），2011International Conference on IEEE，2011：1-4.