淺水湖泊平面二維水質模型研究

朱文謹，周 凱

（淮海工學院土木工程學院，連云港 222006）

淺水湖泊平面二維水質模型研究

朱文謹，周 凱

（淮海工學院土木工程學院，連云港 222006）

磷氮元素在水體中遷移對湖泊富營養化研究具有重要性。采用有限元伽遼金加權余量法建立了簡單的二維水質模型，并計算了在風場作用下的蘇州市澄湖豐水期和枯水期的氮磷的分布以及輸移。計算中引入了可由實測資料率定的可為正負的綜合衰減參數，綜合考慮了氮磷元素的沉降速率和釋放速率，這樣減少了水質多參數難確定的困難。計算的結果表明，流態與實際情況吻合較好，氮磷元素的濃度計算值與實測值相差不大。

磷氮元素；水質模型；綜合衰減系數；澄湖

Biography：ZHU Wen-jin（1981-），male，doctor.

近年來，隨著城市經濟的發展，城市湖泊水體富營養化不斷加劇，嚴重破壞了湖泊的生態環境。為了保護水庫湖泊的水質，科學地利用水資源、必須對水庫湖泊的富營養化狀況進行有效的監測和預測。磷氮元素在水體中遷移規律對湖泊富營養化研究具有重要意義。現擬建一簡單二維水質模型并計算澄湖（蘇州市）中氮磷元素在水體中的遷移規律，通過對于模型的驗證，證實建立模式的準確性，以期今后的研究中為淺水湖泊的富營養化預測和防治提供一定的參考［1-3］。

1 平面二維水流、水質模型

考慮到湖泊一般水深較淺，垂向流速變化較小，平面二維淺水湖泊水動力控制方程為［4-7］

式中：u，v分別為x和y方向上的速度；h為水深；g為重力加速度；zb為河床高程；ρ為水密度；ζ為風應力的經驗系數；W為風速；ψ為x正方向與風向之間的夾角；εij為渦粘性系數；C為謝才系數；ω為地球自轉角速度；φ為當地緯度。

平面二維水質基本方程為

根據伽遼金（Galerkin）有限元法則描述控制方程（1）～（3）可寫為

利用面積坐標的等參單元函數進行坐標變換，單元矩陣微分形式可表示為下面的形式（以水質方程為例）

式中：［k］為單元系數矩陣；［t］為時間矩陣；｛c｝為未知結點濃度向量；｛f｝為單元源、匯矢量。整體矩陣微分方程

采用Crank-Nicholson差分格式對該方程進行時間離散

式中：θ為隱式差分系數；Δt為時間步長；J，J+1分別表示時刻。

2 模型的應用

澄湖地處蘇州市東郊，西臨太湖，北穿吳淞江與陽澄湖息息相通。位于東經120°41′，北緯31°22′，沿湖周長20 km，跨吳中區、昆山、吳江三地，全湖呈三角形，東西長度為9 km，南北長度為10.4 km，全湖面積56.06 km2，約5萬畝，平均水深1.83 m，容積為0.73×108m3。澄湖河道大多數由西和西北方向注入湖泊，經東和東南方向排出泄入淀山湖輾轉入黃浦江。湖泊水位主要受本地的降雨和吳淞江水位的雙重影響。澄湖水位緩漲緩落，年內變幅較小，年中10月至次年2月為枯水期，水深1.70 m左右，3～5月先后入梅，水位上漲，水深1.80 m左右，6～9月為汛期，水位陡漲，水深2.50 m左右。

2.1 澄湖氮磷濃度的測定

以2005年12月下旬（枯水期）和2006年6月上旬（豐水期）所測的澄湖氮磷濃度為計算初始值，計算了枯季5個月，洪季4個月的氮磷濃度分布情況。在湖中心處布置一個采樣點，河流入湖口周邊布置了25個測點，每天采樣2次的濃度平均后作為當天的濃度。同時測定相應水體的水溫、pH值、溶解氧、總磷、總氮、亞硝酸氮、硝酸氮、氨氮、總有機碳、CODMn、葉綠素。

2.2 參數的選取及初邊值條件

根據以往的研究成果、2次監測資料選取參數如下［8-9］：

（1）水動力力參數：曼寧糙率n值為0.015～0.019，渦粘性系數E取10～15 m2/s，緯度31.3°。

（2）水質參數：縱、橫向擴散系數為1.0 m2/s，氮磷綜合衰減系數（可為正負，綜合考慮沉降速率釋放速率。懸沙濃度為0.055 g/m3，中值粒徑0.006 3 mm）采用實測資料估算和率定為：磷豐水期取8×10-4d-1，枯水期4×10-4d-1；氮豐水期取 4.7×10-3d-1，枯水期 3.6×10-3d-1。

計算中的初邊值條件設定如下：

（1）初始條件：假定初始湖水是靜止的，湖面是水平的，即：u（x，y）=v（x，y）=0。

（3）考慮了風對湖流的影響。根據實測資料，12月份為5～7級的北風，6月份為2～3級東南風。

2.3 數值計算及模型驗證

進出澄湖的河道較多，而有些河道的流量很小甚至是沒有流量，模型考慮主要幾個影響較大的進出口。本文采用的有限元方法結合三角形網格，可以較好地貼合自然邊界，提高了計算精度和計算效率。模型的網格在進出口處和部分較窄地方加密；由于湖中心附近水域的水深變化不大，流速較小，因此湖中心附近的網格較大，面積變化范圍為32.53～40 529.84 m2。模型共布置三角形網格單元6 964個，結點數14 465個，網格劃分如圖1所示，計算區域見圖2。

圖1 計算水域網格圖Fig.1 Mesh of computational domain

圖2 計算區域示意圖Fig.2 Sketch of Chenghu Lake

圖3 澄湖12月份流場圖Fig.3 Flow field of Chenghu Lake in dry season

圖4 澄湖6月份流場圖Fig.4 Flow field of Chenghu Lake in flood season

為了能夠清楚的顯示出湖流，流場中大于0.01 m/s的速度均調整為0.01 m/s。由圖3、圖4可以看出，洪、枯季澄湖河道大多數由西和西北方向注入湖泊，經東和東南方向排出，各個進出口區域流場表現為吞吐流流態，而其余廣大湖區流場表現為風生流，基本反映了澄湖的水動力特性。入湖流量小，湖中心位置的由于冬季的風速較夏季的大，湖中心的流速也較大一些。枯水期，總磷平均濃度0.533 mg/L、總氮平均濃度7.7 mg/L；豐水期，總磷平均濃度0.265 mg/L、總氮平均濃度6.26 mg/L（圖5～圖8）。

圖5 澄湖豐水期總磷（TP）分布圖Fig.5 Concentration distribution of TP in flood season

圖6 澄湖枯水期總磷（TP）分布圖Fig.6 Concentration distribution of TP in dry season

圖7 澄湖豐水期總氮（TN）分布圖Fig.7 Concentration distribution of TN in flood season

圖8 澄湖枯水期總氮（TN）分布圖Fig.8 Concentration distribution of TN in dry season

從模型計算的結果與實測值對比（表1～表3）可知，水動力模式計算的湖泊流速在清水港與大姚橋附近的精度較高，湖泊中心處水流流速的精度相對較低，這可能與湖泊中心處流速較小有關。水質模式計算的澄湖總磷和總氮的濃度值與實測值相差較小，基本能夠反映洪枯季澄湖水質的總體情況。枯季澄湖地區的風速較大，從而湖泊中心的流速和風浪也較大，這可能造成湖泊底泥污染物的再懸浮和釋放。而本文的水質模型采用單一綜合衰減系數雖然可以考慮眾多因素，但仍然不能夠充分反映底泥釋放的影響，所以造成了枯季湖泊中心處總磷濃度的計算值誤差稍微偏大，達到18.181%（表3）。

表1 模型計算結果與實測結果驗證流速表Tab.1 Verification of computation value and observation data for velocity m/s

表2 模型計算結果與實測結果驗證總氮表Tab.2 Verification of computation value and observation data for TN concentrationmg/L

表3 模型計算結果與實測結果驗證總磷表Tab.3 Verification of computation value and observation data for TP concentration mg/L

3 結論

本文采用有限元伽遼金加權余量法對二維水流水質控制方程進行求解。模型對風場作用下的蘇州市澄湖豐水期和枯水期的氮磷的分布以及輸移進行了計算模擬。計算中引入了可由實測資料率定的一個綜合衰減參數（可為正負），它綜合考慮了氮磷元素的沉降速率和釋放速率，這樣減少了水質多參數難確定的困難。結果表明，數值計算值與實測值吻合較好，從而為澄湖的富營養化研究提供了實際可行的計算模型。但泥沙對污染物的吸附與釋放是個動態過程，而模式中的綜合衰減參數如何反應這一動態過程，尚需作更多的研究和實際驗證。

［1］梅新敏，阮曉紅，張蘭芳，等.調引太湖水改善蘇州市水動力條件研究［J］.環境科學與管理，2006，36（1）：60-62.

MEI X M，RUAN X H，ZHANG L F，et al.Research of Water Diversion From Taihu Lake to Improve the Water Environmental in Suzhou［J］.Environmental Science and Management，2006，36（1）：60-62.

［2］趙艷艷.金山湖二維水流水質數學模型研究［D］.南京：河海大學，2005：23-28.

［3］張錫輝.水環境修復工程學原理與應用［M］.北京：化學工業出版社，2002：13-15.

［4］李瑞杰，嚴以新，宋志堯.太平水道懸移質輸運數學模型［J］.泥沙研究，2003（4）：5-7.

LI R J，YAN Y X，SONG Z Y．Calculation of Suspended Sediment Transport in Taiping Waterway［J］．Journal of Sediment Research，2003（4）：5-7．

［5］汪德爟.計算水力學理論與應用［M］.南京：河海大學出版社，1989：11-20.

［6］方子云.水資源保護工作手冊［M］.南京：河海大學出版社，1988：43-57.

［7］櫆文信.河流海岸環境學［M］.武漢：武漢大學出版社，2004：78-92.

［8］高永霞，蔡琳琳，趙林林，等.豐水期環太湖河流與湖區水質比較研究［J］.環境科學，2011，32（10）：2 840-2 848.

GAO Y X，CAI L L，ZHAO L L，et al．Water Quality Comparison Between Lake Taihu and Contribute River During High Waterlevel Period［J］.Environmental Science，2011，32（10）：2 840-2 848．

［9］李曉靜，婁安剛，王璟，等.青島堿廠溫排水對婁山河口附近水環境的影響預測［J］.水道港口，2011，32（2）：128-134.

LI X J，LOU A G，WANG J，et al．Prediction of environment effect near Loushan Estuary of thermal discharge by Qingdao Alkali Corporation［J］.Journal of Waterway and Harbor，2011，32（2）：128-134.

Two-dimensional model of water quality and its application in shallow lake

ZHU Wen-jin，ZHOU Kai
（School of Civil Engineering，Huaihai Institute of Technology，Lianyungang222006，China）

In order to predict the lake eutrophication，it is necessary to study the nutrients distribution such as phosphorus and nitrogen pollution.In this paper,the weighted residual method of finite element was used to establish the two dimensional depth-averaged water quality model.Then，the hydrodynamical movement and transport of nutrient material in Chenghu Lake were simulated with important functions.A comprehensive decay parameter of pollutant confirmed by field data was introduced in the calculation，which synthetically considered the settlement velocity and release speed of nitrogen and phosphorus at the same time.The results are in fair agreement with numerical data and observed data.

phosphorus and nitrogen；water quality model；comprehensive decay parameter；Chenghu Lake

X 824；O 242.1

A

1005-8443（2012）04-0353-05

2011-11-28；

2011-12-22

朱文謹（1981-），男，江蘇省江都市人，博士，主要從事港口航道專業。