本文關(guān)鍵詞：多沙河流的非恒定一維水沙數(shù)學模型及其應(yīng)用，由筆耕文化傳播整理發(fā)布。

第 !; 卷 第 8 期 4""$ 年 !! 月

水 科 學 進 展 7.*7EF=G HE I7%=J GFH=EF=

*&)5!;， E&58 E&KL ，4""$

多沙河流的非恒定一維水沙數(shù)學模型及其應(yīng)用

鐘德鈺，彭 楊，張紅武

!"""#$）

（清華大學水利水電工程系，北京

摘要：介紹多沙河流上非恒定一維水沙數(shù)學模型采用網(wǎng)格上流量與水位交錯布置的方案，應(yīng)用 %&’() *(+,(’,&- .,/,-0 （%*.） 格式求解基本方程， 模型具有守恒性好、無虛假振蕩、算法簡單的優(yōu)點，克服了斷面急劇變化時計算容 ,12,-3 易失穩(wěn)的困難。同時，由于采用了適合于多沙河流的泥沙運動基本方程及相關(guān)參數(shù)，因此，該一維模型能夠適應(yīng)多 沙河流高濃度輸沙及大幅度沖淤變化的特征，成功模擬了黃河干、支流上的洪水傳播。 關(guān) 鍵 詞：多沙河流；一維數(shù)學模型；交錯網(wǎng)格；%*. 格式 文獻標識碼：7 文章編號：!""!089:! （4""$） "80"9"80";

中圖分類號：%*!$456；%*!6654

一維水沙數(shù)學模型常用來模擬天然河道時間和空間尺度較大的水流、泥沙運動和河床沖淤變化。在數(shù)值格 式的研究中，近期具有代表性的研究成果是借鑒計算空氣動力學中的優(yōu)秀算法，獲得的精度高、守恒性好的多

［!］ 種計算格式 。此外，對于與泥沙運動有關(guān)的基本方程的研究和改進上取得了較大進展，特別是多沙河流懸移 ［4］ ［6， ［;］ $］ 、水流的挾沙力 、動床阻力 等方面取得了重要成果。 質(zhì)的不平衡輸沙方程

如黃河及其支流等多沙河流，以洪水暴漲陡落且洪水過程中水流挾沙濃度高，河道斷面寬深比大、河床沖 淤幅度大而著名。這就要求應(yīng)用于這類河道的數(shù)學模型應(yīng)具備：! 計算穩(wěn)定，效率高；" 守恒性好；# 能夠 模擬多沙河流特殊的高含沙輸沙現(xiàn)象和河床沖淤。由于前兩項是對模型數(shù)值方法的基本要求，應(yīng)用近期發(fā)展起 來的高性能格式是十分合適的。但筆者的研究發(fā)現(xiàn)，對于河道斷面面積、河寬急劇變化的多沙河流，簡單套用 這些方法仍然會導致較大的計算誤差，甚至計算完全失敗。對比一維水流和氣流運動方程可見，應(yīng)用守恒型方 程求解一維水流運動時，水流運動方程右端存在反映河槽坡度所導致的加速度項，而這一項在氣流運動方程中 是不存在的。在多沙河流中，斷面上高程變化劇烈，而相鄰斷面間變化更大，，局部深槽的發(fā)生和消失非常頻 繁，使兩斷面間深泓連線并不代表這一河段內(nèi)真實河床坡度，當采用兩斷面間深泓線計算河槽坡度時，往往由

［8］ 于局部深槽的出現(xiàn)而過分夸大重力影響。特別是洪水漫灘后若用深泓線計算坡降帶來更大的誤差。譚維炎等

曾建議將壓強項與重力項及河床反力合并為水面坡度，采用 %*. 格式進行計算，并對源項中的過水面積進行 特殊處理。但我們最近的計算發(fā)現(xiàn)，當斷面寬度變化劇烈或者沖淤十分強烈的情況下，仍然出現(xiàn)數(shù)值不穩(wěn)定的 現(xiàn)象。本文所介紹的多沙河流一維非恒定水沙數(shù)學模型就是基于前人的研究成果和經(jīng)驗，引入適合于多沙河流 的泥沙運動相關(guān)方程，對數(shù)值方法中存在的問題進行改進而得到。

!

基本方程及數(shù)值方法

在寬淺河道的計算中，常將河槽加以概化后進行計算，但斷面的概化多少具有任意性，因此本文介紹的模

型中并不對斷面進行概化，描述一維水流運動基本方程仍然為

收稿日期：4""60":0!4；修訂日期：4""$0"!0!" 基金項目：國家自然科學創(chuàng)新研究群體科學基金資助項目（;"44!:"6） 作者簡介：鐘德鈺 （!:9" < ） ，男，甘肅金昌人，清華大學水利系副研究員，博士，主要從事泥沙運動力學、河流數(shù)學模型研

究。 =0/(,)： >2&-3?@A ’1,-32B(5 C?B5 D-

第/期

鐘德鈺等：多沙河流的非恒定一維水沙數(shù)學模型及其應(yīng)用

4$4

!! # !$ & ’ （!） !" !% " !$ # ! ( # $ & # ! （ ) $ * )+ ）# , （"） ! !" !% 式中 ! 為過水面積； $ 為流量； ’ 為區(qū)間入流或分流流量； ( 為斷面靜水壓強； , 為河床對水流的作用力在 % 方向的投影。 ) $ % &!-. /!% 為河床比降，其中 -. 為河床高程； )+ % 0 " $ ’ $ ’ / !" 1( / ) ，其中 1 為水力半徑； 0 為 *+,,-,# 糙率系數(shù)。沖積河流的動床阻力與水流強度、床面形態(tài)有關(guān)。本文介紹的模型中，河槽部分值采用 ［.］ 根據(jù)大量室內(nèi)外觀測資料確定的動床阻力公式計算；灘地部分 0 值根據(jù)植被、地物情況而定。 王士強 如前文所指出的，多沙沖積河流兩斷面間深泓線往往不能反映真實的河床比降，同時又限于斷面測量間距 大，確定準確的河床比降是非常困難的。為避免這一困難，可將壓強項移到方程右端與河床比降、及河床反力 ［/］ 合并為水面比降，即可將方程改寫為 " !$ # ! $ & * # ! !2 # ) （)） + !" !% ! !% 式中 2 為水位。 ［"］ 對于多沙河流的模擬，對懸移質(zhì)輸沙及其導致的河床沖淤計算是關(guān)鍵。本文介紹的模型采用張紅武等 根 據(jù)黃河等多沙河流觀測資料及室內(nèi)試驗總結(jié)出的半經(jīng)驗的輸沙及河床變形方程，即 !!) # !$) & 3 （ ) " * +! ) ） （(） !""4 ! !" !% !-. （ )" * + ! ) ） （.） 56 & * 3! ! " " # !" ［)］ 式中 ) 、 )" 分別為懸移質(zhì)含沙量和水流的懸移質(zhì)挾沙力， )" 可采用適用于多沙河流的張紅武公式 或本文第 ［(］ 一作者獲得的輸沙率型多沙河流挾沙力公式 計算； 4 為河寬；" 為沉速，對于多沙河流需考慮濃度修正； 5 6 為淤積物密度。 3! 、 + ! 為經(jīng)驗系數(shù)，分別為 #

(

)

( )

(

)

3! &

!0. ! (0/. 7" $ $0. "0/. 8 "

) " ) （4） +! & ) " 式中 $ 為混水卡門系數(shù)； 7 " 為摩阻流速； 8 為流速。 由于多沙河流中細顆粒泥沙，例如黃河中下游粒徑小于 $0$". 55 的泥沙，常隨水流強度、上游來流中細 顆粒濃度的變化而在沖瀉質(zhì)與床沙質(zhì)間相互轉(zhuǎn)化，且這類細顆粒往往占懸移質(zhì)的主體，因此區(qū)分造床的床沙質(zhì) ［6］ 與不造床的沖瀉質(zhì)是十分重要的。在本文介紹的模型中采用筆者的理論公式 劃分兩者。由于該公式不僅反映 了水流強度對細顆粒在沖瀉質(zhì)與床沙質(zhì)間相互轉(zhuǎn)化的影響，而且考慮了來流中細顆粒濃度的影響，與實際情況 更為接近，應(yīng)用表明采用該公式后，細顆粒沖淤計算精度有較大幅度提高，具體對比計算結(jié)果另文總結(jié)。 隨著河床沖淤，床沙組成也相應(yīng)發(fā)生調(diào)整。對于床沙級配變化的模擬，本文采用文獻［4］介紹的方法計 算，同時對床沙分層記錄以反映床沙級配的變化歷史。 7 若令 9 &［ ! ， $， !) ， -. ］

( )

[ $0! / +123# ( ) ) ]

(

)

!0!(

（/）

, &

2 7 ［ . & ’ ，* # ! ! # )+ ， （ ) " * +! ) ） ，* 3! （ ) " * +! ) ） 3! /# 56 ］ !""4 !"" !% 則上述方程可寫為向量形式： !9 # !, & . !" !% 對式 （6） ， 可采用一階顯格式有限體積法離散， 即 　 （6）

(

)

， $) ， $] [ $ ，$ !

"

7

G7F

水

科

學

進

展

第 !H 卷

" " ’ " #! & ) * & ! # +$ （#） !" % ! $ & ! ) *# ! $ " " ( ! 通常情況下流量 , 、水位 - 、含沙量 . 等均布置在河道斷面上。但數(shù) 值計算表明，當河道斷面沿程變化劇烈，水位變化幅度較大，河床沖

(

)

淤幅度大的時候，利用動量方程計算流量時常出現(xiàn)鋸齒狀振蕩現(xiàn)象， 即便采用譚維炎等建議的計算源項中過水斷面面積的方法、水面坡降

［!， 6］ 采用迎風格式亦不能避免 。

經(jīng)分析發(fā)現(xiàn)，這種現(xiàn)象是由于流量與水位均定義于單元中點，水

圖! 變量在網(wǎng)格上的布置

面比降不能準確反映子河段壓力、河道反力對動量輸移的影響造成

［6］ $%&’! ())*+&,-,+. /0 1*)%*23,4 %+ 4.*&&,),5 &)%5 的。有鑒于此，本文模型對譚維炎等的方法 進行了改進，即將流量

與水位、含沙量交錯布置，如圖 ! 所示。 此時，動量守恒方程中水面坡度一項自然離散為 -0 ! & -0 " （!7） / # &/! % & !( !( 顯然這樣布置計算變量后，準確地反映了河段水面比降對動量輸移的影響，避免了動量計算中水位微小振 蕩所導致的流量鋸齒狀振蕩的現(xiàn)象。同時河床高程、水位、含沙量定義在斷面上，也使得河床沖淤變化易于計 算。式 （!7） （ /0 8 /0 8 ! ） 中面積" / 可用 1 " 計算。 ! （ )* # )* # ! & 2!!* ） （!!） " " 式中 2 為數(shù)值粘性系數(shù)。由于網(wǎng)格是交錯布置的，質(zhì)量守恒方程和動量守恒方程以及泥沙運動方程應(yīng)在各自 的網(wǎng)格內(nèi)確定界面通量。式 （#） 中 !’ 需滿足 <’ $’ = 條件。

［!］ 數(shù)值通量采用 9:; 格式計算 ，即

)* # ! %

對于邊界的處理，在本文的模型中邊界處采用特征格式計算，即 53 （ . 7 & .4 ） （!"） % &/ 5’ 式中 3 為 >%,-*++ 不變量。這里需要注意的是，邊界處河床比降仍需要小心處理，否則邊界處計算的不穩(wěn)定 會影響整個計算域。根據(jù)作者的計算經(jīng)驗，可采用 ?*++%+& 公式對邊界處由式 （!"） 計算得到的水位或流速進行 校正， 實際計算表明這樣的處理是可行的，而且誤差也很小。 為了說明本文格式對斷面急劇變化時的適應(yīng)性，計算了斷面急劇收縮然后放大情況下的水流運動。圖 " 所 示分別為 , 、 - 交錯布置與非交錯布置，明渠過水斷面突然縮小時兩種方法計算結(jié)果的對比。計算中渠道底 坡為零，斷面為矩形明渠，除 ( @ 6 777 - 處渠寬為 !7 - 外，其余各斷面寬度為 !77 -。對比兩圖可見交錯網(wǎng)格 有效地改善了計算精度，流量未出現(xiàn)振蕩現(xiàn)象。

（ *）交錯網(wǎng)格 9:; 格式計算結(jié)果

（2）非交錯網(wǎng)格 9:; 格式計算結(jié)果

圖"

, 、 - 非交錯布置與交錯布置計算結(jié)果對比

$%&’"

  本文關(guān)鍵詞：多沙河流的非恒定一維水沙數(shù)學模型及其應(yīng)用，由筆耕文化傳播整理發(fā)布。