太湖風生流及污染物輸移擴散數值模擬

舒葉華，高晨晨

(1.上海市水利工程設計研究院有限公司，上海 200061； 2.上海灘涂海岸工程技術研究中心，上海 200061)

湖流從成因上可分為重力流、密度流和風生流。平原地區的淺水湖泊，重力流和密度流微乎其微，風生流是湖流的主要形式[1-3]。風生流是由湖面風場對湖泊表層水體的剪切、拖曳作用形成的流動，能顯著影響淺水湖泊底泥再懸浮與污染物輸移擴散，進而對湖泊水質產生重要影響[4-8]，因此研究淺水湖泊風生流及其驅動下的污染物輸移擴散規律具有重要的科學意義。

太湖是我國五大淡水湖之一，水域面積寬廣，但平均水深不到2 m，是典型的大型淺水湖泊。迄今為止，已有很多學者借助野外觀測或數值模擬等手段對太湖風生流開展了研究，取得了不少有價值的成果。梁瑞駒等[2]發現太湖風生流沿垂向和水平向均有較大變化，流向沿垂向和水平向也并不一致。逄勇等[9]研究了湖底地形和不規則湖區邊界對風生流場的影響。胡維平等[10]基于三維數值模擬揭示了太湖風生流存在垂直切變，表層流向與風切應力方向偏角不大，底層湖流為補償流。周杰源等[3]提出太湖夏季環流產生的主要動力因子為風，風的大小直接影響平面環流流速的大小，風向結合地形基本決定了環流的形態。王建威等[11]基于原位觀測數據，分析了太湖流場垂向分布特征，發現表層水體的流向與風向基本一致，其余各層流場與風場的對應關系并不明顯，甚至出現完全反向的情況。Li等[12]利用現場觀測的風速和湖流資料，研究了風場對太湖不同水深處流場的影響。Luo等[13]建立了太湖三維數值模型，發現在盛行風作用下湖區存在2個順時針環流和2個逆時針環流。Liu等[14]基于二維水流-污染物耦合數值模型，研究了風場驅動下的太湖環流及污染物遷移過程。Liu等[15]認為風向和風速決定了太湖平面環流的方向、強度和位置，并采用示蹤劑方法研究了定常風作用下不同位置的點源污染物泄露對太湖水質的影響。

以往對于太湖風生流及污染物輸移擴散的數值模擬還存在一些不足，如模型區域采用矩形網格進行離散，無法很好地貼合太湖復雜的岸邊界(邊界呈現出鋸齒狀)；網格精度較低，很多甚至是千米量級等，這些不足必然會影響到水流結構的計算精度。本文通過建立高精度的三維數值模型對太湖風生流及污染物輸移擴散進行精細化的數值模擬，以期為太湖水污染治理和水資源調度提供參考。

1 研究區域概況

太湖位于長江三角洲南緣(北緯30°55′40″～31°32′58″、東經119°52′32″～120°36′10″)，湖泊面積2 427.8 km2，實際水面面積2 338.1 km2，湖岸線總長405 km。西側與南側湖岸線平順、呈圓弧狀，東側與北側湖岸線曲折、灣岬相間。太湖是典型的淺水湖泊，平均水深僅1.89 m，最大水深不超過3 m。湖底地形較為平坦，湖中無深槽，也無大片的湖灘，1 m等深線靠近湖岸。環太湖現有出入湖河道200余條，其中入湖河道主要分布在湖區西部和西南部，出湖河道主要位于湖區東部。

本文結合河道流量特征、位置分布等信息，將環太湖口門概化為望亭、銅坑閘、胥口、瓜涇口、太浦閘、長兜港、長興港、大港河、城東港、浯溪橋、百瀆口、雅浦港、直湖港、梁溪河等14個，編號 1～14，具體布置見圖1。

根據太湖小雷山氣象站(圖1)2007年12月至2008年11月的逐時風資料，太湖湖面各月平均風速主要集中在4～6 m/s，年平均風速約為5 m/s。由圖2可知，太湖地區盛行東南風(SE)和西北風(NW)，頻率分別為12.5%和12.3%，其次為西北偏北風(NNW)，頻率為11.7%。西南風(SW)、西南偏西風(WSW)和西風(W)發生的頻率很低，均不超過3%。

圖1 太湖分區及環湖主要口門分布

圖2 太湖湖面風向玫瑰圖

2 模型建立與驗證

2.1 控制方程

模型采用三維不可壓縮雷諾平均Navier-Stokes方程作為控制方程，服從Boussinesq假定和靜水壓力假定。采用交替方向隱式迭代法(ADI算法)對控制方程進行差分離散求解。其連續性方程和兩個水平向動量方程為

(1)

(2)

(3)

式中：vx、vy、vz分別為x、y、z方向的速度分量，m/s；vxs、vys分別為x、y方向的點源速度分量，m/s;η為自由水面高程，m；h為全水深，m；t為時間，s；f為柯氏力參數，s-1；g為重力加速度，m/s2；ρ為水體密度，kg/m3；ρ0為水體參考密度，kg/m3；sxx、sxy、syx、syy為輻射應力張量分量，kg/(m·s2)；νt為垂向紊動黏性系數，m2/s；pa為自由表面大氣壓，Pa；Fx、Fy分別為x、y方向的水平應力；S為源匯項。

2.2 模型區域和網格劃分

模型計算范圍為整個太湖湖區，具體計算區域與網格劃分如圖3所示。由于東側與北側湖岸線較為曲折，為使網格更好地貼合岸線，采用三角形的非結構化網格對模型區域進行水平離散，其中西側與南側岸線平順處網格相對較大，最大網格邊長 1 200 m；對岬角及島嶼附近區域進行了網格加密處理，最小網格邊長約60 m。平面上共劃分了23 827個網格單元，節點總數為12 448個。考慮到太湖為淺水湖泊，為使模型垂向上也具有較高的分辨率，采用Sigma分層法將模型垂向均勻劃分為5層，即每一層的垂向網格邊長都是該位置處總水深的1/5。Sigma分層法最大的優點是垂向分層數不受水深影響，能夠保證在淺水區也具有較高分辨率，此外湖底和等Sigma面重合，可以真實地反映湖底地形。

圖3 計算區域與網格

2.3 模型參數設置

模型綜合考慮重力、風力、底床摩擦力和柯氏力的作用，其中表面風力τs、底床摩擦力τb和柯氏力f的計算公式分別為

τs=ρacd|uw|uw

(4)

(5)

f=2ωsinφ

(6)

由于環太湖各口門的空間尺度相對于整個太湖來說很小，因此模型計算時，將太湖視為一個四周封閉的水體，各主要出入湖河道的流量分別按點匯和點源處理。

2.4 模型驗證

2.4.1水位驗證

采用太湖西山站和吳溇站的實測水位資料對模型進行驗證，實測資料的時間跨度為1個月。各水位站位置見圖1，水位驗證過程如圖4所示。驗證結果表明，水位計算值與實測值基本吻合，模型精度滿足要求。

圖4 計算水位與實測水位的對比

2.4.2風生流驗證

南京地理與湖泊研究所曾基于多年的現場觀測資料，總結了太湖在北風、西風、南風及東南風4種風場作用下形成的環流方向與平面形態。現以東南風作用下的風生流為模擬對象，模型模擬時風速選擇太湖地區的平均風速5 m/s，用上述建立的三維水流模型模擬太湖風生流流場，并將模擬得到的垂向平均流場與南京地理與湖泊研究所的觀測成果做對比，結果如圖5所示。可以看出，模型模擬的垂向平均流場與南京地理與湖泊研究所實際觀測的流場在流態上較為一致，表明模型能夠較好地模擬太湖的水流運動和流場結構。

(a) 太湖風生流觀測成果

3 結果與分析

3.1 盛行風作用下的太湖風生流特征

太湖流域夏季盛行東南風，冬季盛行西北風[16-17]，因此模型中選取東南和西北兩個盛行風向的定常風進行數值試驗，風速值取太湖地區的平均風速5 m/s。圖5(b)和圖6為5 m/s東南風作用下的風生流流場模擬結果。需要指出的是，由于風向相反，5 m/s西北風作用下的垂向平均流場與各分層流場呈現跟東南風完全相反的水流形態和結構特征。

(a) 表層

在5 m/s東南風作用下達到穩定狀態的垂向平均流場中(圖5(b))，湖區西南部存在以大雷山為中心的大范圍順時針環流，該環流途經大雷山東部向東南直至太湖南岸，此后沿太湖西南側岸線折向西北，至蘭右山后再轉至大雷山北面形成閉合順時針環流。太湖湖心區出現一正(順時針)兩反(逆時針)3個環流，其中順時針環流以平臺山為中心循環流動，兩個逆時針環流分別位于順時針環流南側和東北側。此外，在東太湖、胥口灣、貢湖灣以及梅梁灣、竺山灣等湖灣區域出現幾個大小不等、方向不同的小尺度閉合環流。上述結果與王謙謙[18]和胡維平等[19]的研究結果基本一致。

圖6為5 m/s東南風作用下達到穩定狀態的分層流場。可以看出，表層由于受到湖面風場的直接剪切和拖曳作用，風生流流速較大，大致在0.005～0.100 m/s之間，其中西南沿岸以及岬角處流速較大，湖心區流速相對較小，風生流流向基本上與風向一致。底層除沿岸淺水區域以及湖灣內有局部環流存在外，其余區域流向與表層大致相反，表現為很明顯的補償流。受底部摩擦力影響，底層流速較表層流速小，在0.005～0.060 m/s之間。中層流態與底層基本類似，但流速比底層略小。這是由于從湖面到湖底依次為風應力影響控制區、壓強梯度力影響占優區以及底摩擦力影響控制區。表層水體處于風應力影響控制區，受湖面風應力直接影響，流速最大；往下流速逐漸減小，至風應力影響與壓強梯度力影響相平衡時，流速達到極小值；再往下過渡到壓強梯度力影響占優區，流速逐漸增大，至壓強梯度力與底摩擦力影響平衡區出現流速的極大值；底層水體處于底摩擦力影響控制區，流速隨水深的增加而逐漸減小[19]。

對比圖5(b)和圖6可知，垂向平均流場和分層流場存在顯著差異，因此在用平面二維模型研究淺水湖泊風生流時，應充分認識到其在揭示流場垂直結構方面的不足。

3.2 太湖風生流驅動下的污染物輸移擴散

污染物輸移擴散規律研究是太湖水質水環境研究中的一項重要內容。考慮到太湖污染物輸移擴散的復雜性，本文通過在一具有代表意義的入湖口釋放保守示蹤劑，模擬了太湖在5 m/s東南風作用下的污染物輸移擴散情況，模擬總時長為90 d。示蹤劑釋放點設置在“引江濟太”工程[20-22]的望虞河入湖口處(圖1中的概化口門1)，引水流量設置為 100 m3/s，示蹤劑釋放速率為5 kg/s，釋放時間為模擬時間90 d中的前30 d。

3.2.1污染物質量濃度的垂向分布特性

雖然在定常風作用下，太湖流場呈現出明顯的分層流動現象，但由于太湖為淺水湖泊，其垂向尺度遠小于水平尺度，在垂向擴散作用下，污染物質量濃度沿水深方向近似均勻分布。圖7為5 m/s東南風作用30 d后太湖表、中、底層的示蹤劑質量濃度分布。可以看出，對于太湖這類淺水湖泊，由于垂向擴散比較強烈，污染物質量濃度并未呈現明顯的分層現象，與楊益洪[23]對西湖污染物濃度垂向分布特性的研究結果吻合。

(a) 表層

3.2.2污染物輸移及平面分布特征

前述分析表明，風生流驅動下太湖污染物質量濃度在垂向上并未呈現明顯的分層分布，因此以下在對太湖污染物輸移及平面分布特征進行分析時示蹤劑質量濃度場均取垂向平均的質量濃度場。

5 m/s東南風作用下，示蹤劑在不同時刻的質量濃度分布情況如圖8所示。可以看出，10 d時，自望虞河入湖的示蹤劑在太湖東北部沿岸流作用下在貢湖灣北岸形成高濃度“污染帶”。此后，由于受到梅梁灣內部兩個小尺度環流的牽引影響，示蹤劑一部分進入梅梁灣，另一部分則直接越過梅梁灣口進入竺山灣。至50 d時，在以平臺山為中心的順時針環流帶動下，示蹤劑進入到湖心區，并最終在西南部順時針大環流作用下對南太湖水質產生影響。從模擬開始后的第30天直至模擬結束，太湖北部湖灣區的“污染物”質量濃度始終較高，這是由于北部的湖灣區存在多個小尺度環流，給污染物的聚集提供了水動力條件。這也從一定程度上解釋了夏季太湖藍藻水華主要集中暴發在梅梁灣、竺山灣等北部湖灣區的原因。

4 結 論

a. 太湖風生流流態較為復雜，且流場存在明顯的垂向分層。穩定后的風生流場表層流速較大，流向基本上同風向一致；底層流速較小，除沿岸淺水區域以及湖灣內有局部環流存在外，其余區域流向與表層大致相反，表現出補償流的特征。

b. 西北風作用下的太湖風生流流場呈現跟東南風作用下完全相反的水流形態和結構特征。此外，風場驅動下的太湖垂向平均流場和分層流場存在顯著差異，在用平面二維模型研究淺水湖泊風生流時，應充分認識到其在揭示流場垂直結構方面的不足。

c. 風場驅動下的太湖污染物輸移規律體現出同其風生流特征的高度一致性。太湖為大型淺水湖泊，其垂向尺度遠小于水平尺度，在垂向擴散作用下，污染物質量濃度沿水深方向近似呈均勻分布，并未呈現出垂向分層特征。

d. 太湖風生流流速較小，水動力條件較差，再加上在盛行風向的定常風作用下，湖區出現眾多大小不一的環流結構，尤其是北部的湖灣區存在多個閉合環流，具備積聚污染物的水動力條件。這也是太湖藍藻水華主要集中暴發在梅梁灣、竺山灣等北部湖灣區的一個重要原因。

(a) 10 d