黃河下游河道準二維泥沙數學模型研究
張紅武 趙連軍 王光
摘要:本文首先簡要回顧了以往黃河河道泥沙數學模型研究成果,再以水流運動方程及經過作者修正的泥沙運動方程為基礎,同時引入與實測資料相符合的水流挾沙力、動床阻力、泥沙級配等計算公式作為補充方程,構造出黃河下游河道準二維泥沙數學模型。然后,采用1986年11月~1996年10月這10年長系列實測資料,開展了驗證計算。其結果表明,該模型不僅能計算黃河下游河道一般洪水引起的河床沖淤變形,還能成功地模擬出大沙年下游處于強烈淤積時的規律。
關鍵詞:黃河下游河道 準二維泥沙 數學模型
1 黃河河道數學模型研究的簡要回顧
黃河數學模型的研究,是與流域規劃、工程建設和管理運用等生產緊密結合的。早在1955年編制黃河綜合利用規劃技術經濟報告中,就曾在黃河三門峽水庫規劃階段用初級的一維恒定平衡輸沙模型對水庫淤積和下游河道的河床沖刷變形進行計算[1]。當時的計算結果認為,在桃花峪下游沖刷9年后河床刷深27m,沖刷量及沖刷速度顯然比實際夸大甚多。
三門峽水庫投入運用后,庫區嚴重淤積,并迅速向上游延伸,與原來的計算結果截然不同。為研究改建方案,通過實測資料建立了許多泥沙沖淤量與水文要素之間的經驗關系式。麥喬威、趙業安、潘賢娣等學者在分析下游河道沖淤變化及挾沙能力變化規律的基礎上,提出了三門峽水庫下游河床沖淤計算的方法[1]。
在1975年~1977年進行治黃規劃過程中,為研究各種規劃方案對黃河下游的減淤作用,黃河水利委員會及其有關科研、設計部門的科技人員,使用兩種河道沖淤計算的經驗模型開展了大量的計算工作,麥喬威、李保如兩位學者進行了總結。其中,在規劃方案的初步比較階段,使用的是李保如提出的框算模型,只考慮影響輸沙能力的最主要因素(來水、來沙)。對于問題比較復雜和初選的方案,則使用較詳細的經驗模型[2]。
20世紀80年代,劉月蘭、韓少發、吳知等學者在上述經驗模型的基礎上,根據不同河段汛期 、非汛期輸沙經驗公式,引入水流連續方程、動量方程及泥沙平衡方程式,建立了一套頗具特色的黃河下游河道沖淤計算模型[3]。該模型在黃河規劃中得到應用,只是在第三次加固大堤設計中,預計的下游淤積量與實際相差甚大。1987年~1990年,謝鑒衡、韋直林、魏良琰同劉月蘭、張啟衛等,開展了黃河下游洪水演進預報水動力學數學模型研究。特別是1992年之后,對國家“八五”重點科技攻關項目研究的全面啟動,使黃河數學模型研究有了進一步發展的機遇,在錢意穎、張啟衛主持下,開展了“黃河泥沙沖淤數學模型的應用”專題的研究[4],在此期間,黃科院、清華大學、武漢水利電力大學及中國水科院的有關學者,對黃河數學模型的研究都有較大貢獻。附帶指出,黃委會有關專家曾對美國最通用的泥沙模型HEC-6模型,用黃河資料進行了檢驗,認為這類模型直接用于黃河必須做很多的工作。另外,美國學者楊志達博士等1986年開發了GSTARS模型,參照我國經驗進行修正后,對黃河楊集至孫口河段進行了驗證計算。以謝鑒衡為首的專家組的驗收意見指出:“本模型要應用于包括游蕩性河段在內的黃河下游,還需進一步做更多的研究工作”。此后,世界銀行又委托沈學文對該成果進行評價,提交的書面意見認為:“模型可以用于黃河某些不平衡輸沙問題不是十分嚴重的順直河段,但橫向輸沙很麻煩且需要豐富的經驗”顯然,海外數學模型同黃河實際尚存在著很大的距離。
在黃河平面二維泥沙數學模型的研究方面,韋直林[5]采用有限元法對黃河下游夾河灘至高村河段的1982年的洪水演進過程進行了模擬:程曉陶、劉樹坤[6]以二維恒定流與不平衡輸沙理論為基礎,將顯式有限差分法與有限體積法相結合,并引入張紅武水流挾沙力公式計算泥沙沖淤變化,研制出夾河灘至高村河段灘區水沙運行數學模型。結合工程實際,張世奇、楊國錄等也對黃河京廣鐵橋以上河段的沖淤規律進行了數值計算。與此同時,張紅武、江恩惠、趙連軍、李煒、鄭邦民、楊明、鐘德鈺、張凌武等學者聯合攻關,開展了黃河下游不同河段的平面二維數學模型研究。例如,張紅武、李東風、許雨新[7]等應用黃河泥沙運動基本理論研究的最新成果,采用有限元法建立了濟南河段平面二維水流泥沙數學模型。從河床沖淤形態、流量水位過程、全域流速場等方面驗證結果可以看出,該數學模型計算結果與實測資料符合較好。
應該指出,由于黃河問題的復雜性;現有的數學模型對所存在的參數,多進行經驗處理或用實測資料率定。而黃河水沙條件及河床邊界條件如果發生很大變化,那些經驗性較強的參數的處理,往往會影響新條件下的計算結果,因而不能用來預測黃河未來沖淤趨勢。另一方面,水流挾沙力及河床糙率是黃河泥沙數學模型應考慮的兩個重要問題,不少模型中沒有解決這兩個問題。因此,即使在驗證時差別不大,但在重要的特殊洪水期的計算結果,與原型卻差異甚多,大大影響了數學模型在黃河防汛規劃中的使用價值。
再者,恢復飽和系數的確定也是黃河泥沙數學模型的關鍵和難點所在,理論值應大于1,而一般的黃河模型中的取值都遠小于1,其處理值得進一步商榷。為此謝鑒衡指出:“現階段 對黃河數學模型上的恢復飽和系數,無論理論認識或實際經驗都甚感不足,有必要在這兩方面做進一步的工作。”我們認為,該問題主要是在假定非平衡狀態下的恢復飽和
系數與平衡狀態下相等之后出現的,存在著理論上的缺陷,因此,無論對該系數如何進行改進,都難以從根本上提高泥沙數學模型的理論水平。
作者針對上述存在的問題,長期開展了黃河數學模型的研制工作,本文介紹黃河下游準二維泥沙數學模型的研究結果。
2 準二維泥沙數學模型的建立
2.1 基本方程 本模型計算除選用描述河道水流運動的基本方程外,還改正了在泥沙數學模型中占據重要位置的河床變形方程及泥沙連續方程式[8,9],即:
河床變形方程
泥沙連續方程
式(1)~(2)中:i為斷面號;j為子斷面號,河床高程最低的子斷面號j為1,最高的取j為m;A為過水面積;t為時間;x為沿流程坐標;Z為水位;ω為泥沙渾水沉速;S為水流含沙量;S為水流挾沙力;γ0為河床淤積物干容重;bij為子斷面寬度;Zbij為子斷面平均河床高程;α為平衡含沙量分布系數,由下式計算[10]:
其中:
式中:κ為卡門常數;cn為渦團參數(cn=0.375κ);C為謝才常數;u為摩阻流速;且有:
運用量綱分析及相似轉化原理,并根據動床河工模型試驗結果,可求得Ks的表達式為:
對于非飽和系數fs,通過歸納分析,再借助試驗資料,得出:
2.2 水流挾沙力及糙率計算 為保證數學模型適用于黃河下游高含沙洪水,提高多沙河流輸沙計算的精度,水流挾沙力計算采用黃河上常用的張紅武水流挾沙力公式[10]
式中:
式中:ωsk為第k組粒徑泥沙在渾水中群體沉速;ωok為第k組粒徑沉速;D50為床沙中值粒徑;d50為懸沙中值粒徑(mm);psk為第k組粒徑泥沙的重量百分數;N為非均勻沙的分組數。
經舒安平、江恩惠、朱太順及陳雪峰、陳立、李義天[11]等學者,通過黃河、長江等河流大量實測資料系統檢驗結果表明,式(8)是現有公式中最適用于黃河等天然河流的水流挾沙力公式,目前已得到廣泛應用。而且教科書中也被推薦為適用于高含沙水流的公式[12]。
本模型采用作者的動床糙率公式開展黃河洪水水位模擬計算,可較好地描述水力泥沙因子的變化對摩阻特性的影響,并能反映天然河道中各種附加糙率的影響。該糙率公式為:[13]
式中:δ為摩阻厚度,采用下游資料,得出δ與床沙中徑D50及佛汝德數Fr的關系式為:
2.3 懸移質泥沙與床沙級配的計算 從懸浮于水體中泥沙顆粒的受力分析入手,經理論推導和經驗處理可建立懸移質泥沙級配及粒徑二階圓心距的計算公式[14]
上述3式中,P(di)為小于某粒徑的泥沙所占的百分數;d50為懸沙中值粒徑;dcp為懸沙平均粒徑;m為指數;ξd為懸沙粒徑的二階圓心距,其數值可表征出泥沙組成的非均勻程度。由式(15)、(16)可知,對于變量d50、dcp、ξd,已知任意兩個,即可確定第三個變量的值,懸沙組成即可由式(14)求解出來。同時,通過黃河下游 水文站實測資料回歸分析,建立了懸沙中值粒徑d50和平均粒徑dcp之間的如下關系式:
式中:vm為渾水運動粘滯系數。
對于天然細沙河流,底沙隨時可能跳起,同懸沙的交換幾率很大。我們以躍動模式分析床沙起動時的受力入手,推導出與懸沙相類似的級配及粒徑二階圓心距公式:
式中:P(Di)為床沙中小于某粒徑的泥沙所占的百分數;Dcp為床沙平均粒徑;m為指數;ξD為床沙粒徑的二階圓心距。
從懸沙與床沙交換過程中任一粒徑組的沙量平衡入手,經理論推導建立了如下懸沙與床沙平均粒徑與二階圓心距變化規律的表達式:
式(21)~式(24)中:dL、ξL分別為側向入流泥沙的平均粒徑及粒徑的二階圓心距,P為床沙混合層厚度,dc、ξc為沖淤物平均粒徑與粒徑的二階圓心距,dc與P的具體計算方法
通過對黃河下游實測資料回歸得出如下懸沙平均粒徑沿橫向分布公式:
式中:dcpi為斷面平均懸沙平均粒徑,dcpij為斷面上任一點懸沙平均粒徑。C2為斷面形態系數,由沙量守衡可求得
2.4 懸移質含沙量橫向分布計算 通過對黃河大量實測資料分析發現,含沙量的橫向分布規律不僅與水力因子、含沙量大小有關,還與懸沙組成密切相關。懸沙粒徑越細,含沙量的橫向分布越均勻,為此除引入含沙量因子外,還引入懸浮指標ω/ku*,建立了含沙量橫向分布公式,即
式中:Vi、Vij分別為斷面平均及任意一點的流速;Hi、Hij分別為斷面平均及任意點的水深;u為斷面平均摩阻流速;C1為1左右的斷面形態系數,由沙量守恒可求得:
式中:qij為斷面任一點單寬流量,y為橫向坐標;a、b為斷面河寬兩端點起點距。
2.5 河槽在沖淤過程中河寬變化規律的模擬 本數學模型所采用的模擬河槽寬度變化過程的方法,是根據計算河段應受特定的河相關系均衡調整原理進行的。模擬河段內各橫斷面寬度修正值,除受到給定河岸條件的限制外(例如受到抗沖河岸、山嘴及河道整治工程等的限制),還要隨著河流造床過程而自動調整。本模型選用張紅武河床綜合穩定性指標[9]作為河相關系均衡調整準則計算,即:
式中:B為河寬;i為河床縱比降。作為初步探討,將河床的均衡關系式運用于河床斷面的自動調整過程,對于河寬(包括子斷面寬度)的調整以式(29)的約束關系確定其變化值。
3 驗證計算結果
該數模建立之后,主要針對黃河下游諸如“82.8”、“77.7”、“92.8”、“96.8”等不同類型洪水開展了驗證計算①。為進一步檢驗該數學模型對不同水沙系列年的適用性,特選用三門峽水庫蓄水運用及滯洪運用初期1960年9月~1965年6月的水沙系列和近期的1986年11月~1996年10月水沙系列以及相應的沖淤變形資料,對模型進行檢驗。前者見文獻[15],本文只介紹后者的驗證計算結果。
1986年~1996年,黃河下游總來水量為2979.6億m3,總沙量為85.742億t,其中汛期來水量為1404.8億m3,占總來水量的47.1%,來沙量為81.7億t,占總來沙量的95.3%。特別是1988年、1992年、1994年、1996年,這4年為豐沙年,汛期小浪底水文站平均含沙量均超過80kg/m3,汛期來沙量為48.6億t,占總來沙量的56.6%。由于水沙搭配不合理,造成了下游河道嚴重淤積,10年內下游河道共淤積25.0億t,其中鐵謝至高村河段淤積19.2億t,占總淤積量的76.8%,高村至艾山河段淤積3.6億t,占總淤積量的14.4%,艾山至利津河段淤積2.2億t,占總淤積量的8.8%。從淤積的年際年內分配上看,河道主要在汛期淤積,共淤積31.4億t,而非汛期河道基本處于沖刷狀態,共沖刷6.4億t。其中1988年、1992年、1994年、1996年四年汛期淤積最為嚴重,共淤積22.4億t,占總淤積量的89.5%。從泥沙的淤積部位看,大量泥沙主要淤積在主槽內,引起河槽高程明顯抬升,平灘流量大幅度下降。如“92.8”、“96.8”兩場洪水,花園口站洪峰流量分別為6020m3/s、7600m3/s,接近于該站洪峰的多年平均值,僅相當于中常洪水,花園口站卻先后出現有水文記載以來的最高水位,洪水漫灘范圍及傳播時間也明顯增大,洪水造成的經濟損失巨大。
采用小浪底站實測水沙過程及懸沙組成作為進口邊界條件。伊洛河、沁河入黃水沙及泥沙級配過程分別采用黑石關、小董(或武陟)站實測資料,沿黃引水引沙量采用實測資料沿程分配給出,作為已知的分流條件。初始地形采用黃河下游1986年汛后92個實測大斷面控制,并將每一個大斷面概化為12~15個子斷面。初始床沙組成根據花園口、夾河灘、高村、孫口、艾山、濼口、利津等7站的實測床沙組成按斷面間距內插或外延。出口邊界水位采用相應年份 利津站水位流量關系控制。
表1給出了鐵謝至利津不同河段不同時段計算和實測的沖淤量。計算表明,1986年11月~1996年10月,鐵謝至高村河段淤積13.588億m3,高村至艾山淤積2.445億m3,艾山至利津淤積1.654億m3,與相應的實測值13.707億m3、2.581億m3、1.571億m3接近。
從上述計算結果可以看出,無論河道處于沖刷、淤積或沖淤交替狀態,即使是遇到大沙年(包括其中的高含沙洪水)處于強烈淤積時,該模型均能很好地模擬黃河河道的輸沙特性及河床變形規律。同時,該模型對河道在不同水沙組合條件下發生的上沖下淤、上淤下沖及沿程沖刷或淤 積等規律的模擬精度也令人滿意。
4 結 語
本文首先回顧了前人關于黃河河道泥沙數學模型研究成果,再以黃河泥沙基礎理論最新成果 為基礎,構造出黃河下游河道準二維泥沙數學模型。然后采用原型實測的水沙過程及相應的河床變形資料,開展了驗證計算。結果表明:該模型不僅能模擬下游河道長系列年一般洪水引起的河床沖淤變形,還能成功復演大沙年下游處于強烈淤積時的實際情況。該數學模型已在黃河下游防洪及規劃中應用,而且還在小浪底水庫運用方式研究中發揮重要作用。目前被黃 河主管部門確定最為可靠的黃河下游泥沙數學模型。
參 考 文 獻:
[1] 麥喬威,趙業安,潘賢娣.多沙河流攔沙水庫下游河床演變計算方法[J].黃河建設.1965,(3).
[2] 李保如,華正本,樊左英,陳上群.三門峽水庫攔沙期下游河道的變化[c].見:第一次河流泥沙國際學術討論會論文集(第一卷).北京:光華出版社.1980.
[3] 劉月蘭,韓少發,吳知.黃河下游河道沖淤計算方法[J].泥沙研究.1987,(3).
[4] 錢意穎,張啟衛,李松恒.黃河泥沙數學模型的研究與應用[C].黃河泥沙.鄭州:黃河水利出版社.1996.
[5] 韋直林.黃河下游水流泥沙平面二維數學模型研究[J].武漢水利電力學院學報,1990,6:77-86.
[6] 劉樹坤,宋玉山,程曉陶,等.黃河灘區及分滯洪區風險分析和減災對策[M].鄭州:黃河水利出版社,1999.
[7] Zhang Hongwu, Li Dongfeng, Xu Yuxin, Zhang Junhua. A Two Dismesional Sediment Mathematical Model for Unsteady Flow in the Lower Yellow River [J]. International Journal of Sediment Research.1999,14(2).
[8] 張紅藝,楊明,等.高含沙水庫泥沙數學模型的研究及應用[J].水利學 報,2001,(11).
[9] 張紅武,呂昕.彎道水力學[M].北京:水利電力出版社,1993.
[10] 張紅武,江恩惠,等.黃河高含沙洪水模型的相似律[M].鄭州:河南科技出版社 ,1994.
[11] 陳雪峰,陳立,李義天.高中低濃度挾沙力公式的對比分析[J].武漢水利電力大學學報,1999,(5).
[12] 張瑞瑾,謝鑒衡,等.河流泥沙動力學[M].北京:水利電力出版社,1998.
[13] 趙連軍,張紅武.黃河下游河道水流摩阻特性的研究[J].人民黃河,1997,(9).
[14] 趙連軍,張紅武,江恩惠.沖積河流懸移質泥沙與床沙交換機理及計算方法研究[J].泥沙研究,1999,(4).
[15] 張紅武,黃玉婕,趙連軍,等.小浪底水庫運用初期下游河道沖淤數學模型預測計算[J].水力發電學報.2002,(1).