郭曉明，喬 明

(1.丹華水利環(huán)境技術(上海)有限公司，上海，200235；2.浙江省環(huán)境科技有限公司，杭州，310007)

在數(shù)值模擬的計算中，糙率一般都是需要經過率定的[1-3]，即在某個糙率下數(shù)模結果與實測結果或試驗結果進行比對，當相差較小時我們認為這個糙率就是模擬中可以使用的糙率。但是，在一些擬建的工程中如隧洞、渠道等，并無實測數(shù)據(jù)或試驗數(shù)據(jù)，這時僅靠邊壁的材料去確定糙率有失偏頗，特別在二維水動力學計算中會產生大的偏差。二維水動力學模擬中僅是考慮底部的影響，在模擬水深相對底部寬度比較大的隧洞和渠道時，直接用糙率系數(shù)表上的糙率進行計算的結果不準確，因此，需要增大河床單元格的綜合糙率，以近似代表邊壁糙率對流動的整體影響。

1 二維水動力學模型糙率的公式推導[4-7]

假定模擬無任何誤差，那么二維數(shù)模計算出的流速等于實際的流速(即V′=V)。利用謝才公式和曼寧公式[8]：

兩式聯(lián)立可以得出：

謝才公式是根據(jù)明渠均勻流的實測資料總結的經驗公式，在實際的工程中大多是阻力平方區(qū)的紊流，謝才系數(shù)的經驗公式是根據(jù)阻力平方區(qū)紊流大量實測資料求得，所以適合我們的計算。

由于水力坡度相同，可以得出：

矩形過水斷面情況下：

其中B為底部寬度；H為淹沒深度。

2 二維水動力學模型控制方程[9]

二維水動力模型采用的公式如下：

質量守恒方程：

X方向動量方程：

Y方向的動量方程：

以上各式中，h為基準面下水深；u、v分別為沿X、Y方向的流速；t為時間；η為床面高程；g為重力加速度；Exx、Eyy、Exy為渦動擴散系數(shù)；a為河床高程；n為曼寧糙率系數(shù)；ζ為風應力系數(shù)；Va為風速；ψ為風向；ω為地球角速度；φ為當?shù)鼐暥取?/p>

3 糙率公式驗證

為了驗證推論的合理性和正確性，我們用下莊電站引水設計的物理模型試驗結果和藕節(jié)水槽試驗結果及矩形水槽水力計算結果與二維數(shù)值模型的計算結果進行對比。

3.1 下莊電站引水設計[10]

下莊電站位于岷江支流雜谷腦河下游河段左岸，為雜谷腦河規(guī)劃十八級梯級電站中的第17級電站。由取水樞紐、泄洪閘、沖砂閘、沖砂槽、進水口、引水明渠、泄水道等建筑物組成。取水樞紐設計洪水流量916m3/s，校核洪水流量1180m3/s。為了更好地發(fā)揮電站的經濟效益和在災后重建中的作用，對電站引水系統(tǒng)進行技改，在原有雜谷腦河引水明渠旁開挖引水隧洞，引水隧洞采用馬蹄形斷面，底部寬度為6.6m，內部流動為無壓水流。

根據(jù)下莊電站的位置和所擁有的實測資料及需要計算的范圍，我們選取的模擬區(qū)域長約為1.47km，河道段長140m，引水隧洞和引水明渠長平均1.33km。采用三角形網格與四邊形網格離散計算區(qū)域，在交匯處及不規(guī)則區(qū)域主要采用三角形網格，其它則采用四邊形網格，模擬區(qū)域網格共12265個，其中三角形網格1314個，四邊形網格10951，節(jié)點38029個。

引水隧洞斷面為馬蹄形斷面，但在二維數(shù)值計算中近似取為矩形斷面。上游給以流量，下游給以水位邊界條件。在上游來水為160m3/s時，關閉所有的泄洪閘，下游控制水位為1403.5m(前池正常水位)。根據(jù)材料不同查出的初始糙率取值為：雜谷腦河率定糙率n=0.035；原引水明渠n=0.018；引水隧洞n=0.014。經過計算新建引水隧洞內平均水深約為4.3m，原引水明渠內平均水深約為3.5m。引水隧洞底部寬度為6.6m，引水明渠底部寬度為7.7m。

圖1 網格結構

再根據(jù)修正后的糙率進行模擬，得出閘前斷面CS1的水位對比和沿隧洞內前300m的水深，分別見圖2和圖3。

圖2 閘前水位對比

圖3 隧洞內水深對比

用初始的糙率進行計算的結果和實驗數(shù)據(jù)相差很大，但是用經過修正后的糙率進行計算的結果和試驗數(shù)據(jù)相差很小。圖2可以看出，未經修正糙率模擬的閘前水位與試驗數(shù)據(jù)相差60cm，而經過修正后相差不到5cm；圖3中的未經修正糙率模擬的水深相差50cm以上，但是經過修正后水深相差大幅縮小。另外，由于公式中所用水深為初始糙率進行模擬的結果，在用計算修正后糙率時使水深偏小，故可以用修正后糙率模擬出的水深重新計算糙率，會更加合適。

3.2 藕節(jié)水槽試驗

在四川大學水力學國家重點實驗室做的藕節(jié)水槽試驗，水槽長28m，坡度為0.002，有三個藕節(jié)，最寬為1.4m，最窄為0.6m，深為0.5m，材質采用抹面混凝土。模擬中取15.2m的長度，示意圖如圖4。

圖4 藕節(jié)示意

網格為四邊形網格共2000個，節(jié)點數(shù)為6421個。網格分布圖如圖5。

圖5 藕節(jié)網格劃分

進口給定流量為100L/s，出口給定水位為10.18cm，根據(jù)試驗中底部和壁面的材料，經過查詢吳持恭《水力學》糙率系數(shù)表，糙率取0.013。關于E的取值比較關鍵，在水深較淺時E應該盡可能的小，在模擬中E取25Pa/s。根據(jù)初始單一糙率0.013模擬的水深，以及不同的底部寬度可以得出修正糙率沿程的變化。

圖6 糙率沿程分布

然后根據(jù)糙率把模擬段分為8個部分，每個區(qū)域糙率取平均值，其值在0.0158～0.0165之間，相差較小。再根據(jù)修正后的糙率進行模擬，得出模擬段內的水位。

圖7 沿程水位對比

可以看出，當采用糙率0.013時其在收縮處的水位變化極大，且和實驗數(shù)據(jù)差別很大，達到4cm；當采用修正后的糙率時其差別縮小，在收縮處不到1.5cm。

3.3 理想矩形水槽

設定一個長100m，寬10m，底部坡度為1/1000的矩形水槽，糙率為0.04。假設過水流量為100m3/s、200m3/s、300m3/s，矩形水槽流動為均勻流?？梢酝ㄟ^謝才公式和曼寧公式得出：

通過試算，可以得出三個流量對應的水深H分別為6.37m、11.09m、15.62m，流速V分別為1.57m/s、1.804m/s、1.92m/s。上游進口處給定流量，下游出口給定水位(出口處地形高程為0，故水位即水深H)，對矩形槽進行模擬。計算中E取2000Pa/s。

圖8 網格及平面示意

由于水深變化較小，可以直接用下游出口水深修正糙率。

代入H=6.37m、11.09m、15.62m?？梢缘贸鋈齻€修正糙率分別為0.0692、0.0872、0.1029。根據(jù)三種不同的工況和不同的糙率進行模擬。水位、流速沿程變化圖如圖9。

(a)流量100m3/s時流速沿程變化 (b)流量100m3/s時水位沿程變化

(c)流量200m3/s時流速沿程變化 (d)流量200m3/s時水位沿程變化

(e)流量300m3/s時流速沿程變化 (f)流量300m3/s時水位沿程變化

根據(jù)均勻流的定義，水力坡度應該與底部坡度相一致且流速不發(fā)生變化。從圖9中可以看出，在流量100m3/s、200m3/s、300m3/s情況下，根據(jù)矩形槽糙率0.04所模擬出的水力坡度約為3/10000，遠遠小于底部坡度1/1000。而用修正糙率計算的結果完全符合均勻流的定義，上下游水位差0.1m，水力坡度1/1000。修正前流速沿程逐漸減小，修正后流速與均勻流公式計算出的幾乎完全一致，相差不到0.001m/s。

3.4 驗證結果分析

我們本可以通過試驗得出糙率的率定，但是不經過率定而直接通過糙率修正公式進行糙率的選取也可以得出滿意的結果，這就為如果沒有實測數(shù)據(jù)及試驗數(shù)據(jù)的河道、渠道、隧洞等計算提供了好的糙率修正準則，將大大改進數(shù)模計算的質量。

在藕節(jié)模擬中，由于斷面的收縮擴張，水流流態(tài)變化大，流動為恒定非均勻流，導致修正后的模擬數(shù)據(jù)與試驗數(shù)據(jù)存在偏差，但是通過修正后的糙率值相比直接用糙率系數(shù)表中的糙率有較大的改進。在下莊電站引水的模擬中，隧洞和明渠中的流態(tài)近似漸變流，流態(tài)比較穩(wěn)定，修正后的模擬效果更好。理想矩形水槽的數(shù)模結果與矩形水槽均勻流公式計算結果完全一致。修正公式對于均勻流和漸變流有較好的結果。

4 結論

在二維水動力模型計算中，隧洞和渠道等由于缺少實測資料和試驗資料，不能對其進行準確的率定。我們可以通過查糙率表，通過公式計算得出模型需要的糙率。公式更加適合于寬深比較大的渠系，公式為二維水動力模型計算采用糙率提供參考。即使修正后的糙率還是有一定的差別，因為修正公式使用范圍是恒定均勻流或漸變流。