陳偉倫王 偉
(1.河海大學港口海岸與近海工程學院,江蘇南京 210098;2.中交海洋投資控股有限公司,海南三亞 572000)
泥沙沖淤計算主要依賴以下4種方法:淤現場資料分析法;于物理模型試驗法;盂數值模擬法;榆經驗公式法。目前多采用數值模型法和經驗公式法。在河口治理工程對水沙動力和沖淤演變的影響方面,竇希萍等[2]研究了長江口深水航道治理工程實施后的航道回淤情況;王兆華等[3]分析了治理工程對北槽攔門沙河段的影響;王艷姣等[4]利用GIS技術對長江口北港河段的沖淤演變進行了可視化分析;陸忠民等[5]評價了青草沙水庫建設的可行性;時連強等[6]探討了沒冒沙形成及演變規(guī)律;陳界仁等[7]計算出三角洲淤積影響平均寬度和延伸的平均距離。而在河床泥沙回淤經驗公式的研究方面,劉家駒等[8]提出了用于淤泥質海岸航道開挖后泥沙回淤計算的方法;曹祖德等[9]根據水流歸槽和平衡泥沙質量濃度的概念,推導出順岸式碼頭前沿的港池淤積計算公式;樂培九等[10]提出了適用于河口海灣地區(qū)挖槽回淤估算方法;王義剛等[11]提出了計算邊灘圍墾淤積的一種簡便方法。
本文通過建立大范圍無結構三角網格的高精度數值模型,探究青草沙水庫工程建設前后周邊河口范圍水流變化和沖淤演變規(guī)律,以期為河口水庫工程的建設提供更為可靠的工程咨詢,并對只能應用于計算回淤強度的王義剛公式進行改進,使其可以同時計算河床沖淤。
依托MIKE21水動力模型軟件中HD模塊及其附屬的MT模塊來建立高精度的數值模型,其連續(xù)方程:
式中:x、y為右手Cartesian坐標系;t為時間;濁為水位;h為總水深;d為靜水深;u、v分別為流速在x、y方向上的分量;Pa為當地大氣壓;籽為水密度,籽0為參考水密度;f為 Coriolis 參量;Sxx、Sxy、Syx、Syy為輻射應力分量;Txx、Txy、Tyx、Tyy為水平黏滯應力項;子sx、子sy分別為x和y方向上的表面切應力分量;子bx、子by分別為x和y方向上的底部切應力分量;S為源匯項;us和vs分別為x和y方向源匯項水流流速;Dx、Dy分別為x、y方向上的泥沙擴散系數;Fs為泥沙沖淤函數;酌d為床沙干容重;濁b為海底床面的豎向位移(即沖淤變化量);贅為地球自轉角速率,贅=0郾729伊10-4s-1;漬為地理緯度。
在缺少實測現場干容重資料的情況下,可以通過泥沙中值粒徑公式[7]估算。對于長江口區(qū)域,干密度取值相對穩(wěn)定,酌d取值720 kg/m3。經過數值模型的反復調試,絮凝沉降的相關控制參數取值為:泥沙密度籽s取值2650 kg/m3,絮凝開始發(fā)生時的泥沙質量濃度(含沙量,下同)cfloc取值0郾01 kg/m3,阻礙沉降開始發(fā)生時的泥沙質量濃度chin取值10 kg/m3,凝膠點泥沙質量濃度cget取值50 kg/m3。由Richardson鄄Zaki公式[12]結合侯志強等[13]的觀點,沉速系數棕取為0郾0005 m/s。選取臨界淤積切應力在0郾10~0郾25Pa之間,泥沙臨界沖刷切應力在 0郾15~0郾5 Pa之間。粗糙高度選用模型推薦值0郾001 m。
本文模型涵蓋長江口和杭州灣大范圍區(qū)域,上邊界取在江陰驗潮站,南邊界至浙江象山縣,北邊界取在江蘇如東縣中部,東部外海開邊界位于123毅34憶E,從北至南約為50~60 m等深線位置。模型東西長約310 km,南北長約為311 km,整個模型面積為5郾75伊104km2。此外,在大模型內嵌套建立了僅包含長江口區(qū)域的小模型,北邊界、南邊界及東邊界為水邊界,給定潮位過程,該區(qū)域潮流受東中國海潮波系統控制,開邊界的潮位過程由全球潮波數學模型提供,上游邊界條件采用江陰潮位站實測潮位資料。
考慮到實際研究區(qū)域地形的復雜性,模型網格采用非結構三角形網格,網格生成考慮地形邊界和潮流的走向,整個計算區(qū)域包括個28232個網格節(jié)點,54028個單元,網格間距在200~5 000 m之間,研究區(qū)域水下地形圖及網格如圖1所示。本文模型以1985國家高程基準為計算基面,坐標系統為北京1954坐標系,中央經線123毅。
圖1 研究區(qū)域水下地形及計算網格
1.4.1 潮流模型驗證 選取2002年3月1—9日大、小潮全潮水文測驗的結果進行驗證。各驗潮站、測流垂線及流量斷面如圖2所示,部分潮位站大潮驗證結果如圖3所示,部分測點大潮驗證結果如圖4所示,北港、南港斷面大潮驗證結果如圖5所示。本文約定落潮流速及落潮流量為正值。
1.4.2 懸沙模型驗證
圖6反映了2002年3月1日0時至10日0時Y6和Z9測點大、小潮過程中泥沙質量濃度計算值與實測值的比較,表明所建模型基本能反映不同潮型下泥沙質量濃度隨時間變化的情況,計算結果與實測值基本接近,滿足后續(xù)工作的精度需要。
1.4.3 河床沖淤驗證
基于三角形無結構網格建立水流泥沙模型和潮流泥沙數值模擬,利用2007年12月至2009年10月實測水下地形進行沖淤演變計算驗證,見圖7。
由圖7可以看出,計算得到的沖淤區(qū)域及沖淤強度與實測結果基本吻合,表明數值模型中相關參數的取值合理,數值模擬的計算結果滿足工程研究的需要。
圖2 長江口水文測驗布置示意圖
圖3 高橋站大潮驗證結果
圖4 X3測點流速及流向大潮驗證結果
圖5 部分斷面流量大潮驗證結果
圖6 部分測點泥沙質量濃度大、小潮驗證結果
圖7 模型計算沖淤結果與實測地形對比(單位:m)
由于工程的影響,工程區(qū)域附近的流速和流向均發(fā)生了改變。根據實測資料,水庫建成后,長興島南側航道流速呈增大趨勢,小部分區(qū)域流速下降;水庫工程前沿取水口處,水流流向發(fā)生改變,方向沿著水庫圍堤向北側航道偏轉,同時流速增加;水庫工程北側航道圍堤附近流速微增,圍堤較遠處流速下降;水庫工程尾部流速下降,工程對水流的影響范圍呈扇形向外擴散,距離工程區(qū)域越近,流速和流向改變越明顯。南側航道的分流比增加,水庫工程的建設使得北側航道流速整體下降,南側航道流速整體上升。
2.2.1 數值模擬法
水庫工程的實施對周邊泥沙沖淤強度分布的影響主要集中在水庫工程的兩側區(qū)域,分別布置了15個采樣點(圖8)進行定量分析。青草沙水庫工程實施后,工程對水下地形沖淤影響如圖8和圖9所示:水庫口門A區(qū)域至B區(qū)域,由于工程圍堤前沿的挑流作用,口門前圍堤區(qū)域水流流向急轉,流速增加,在工程圍堤拐角處B區(qū)域出現沖刷,沖刷呈現由口門到工程圍堤拐角處增大的趨勢,沖刷量最大達5郾0 m/a左右。C區(qū)域到D區(qū)域以及K區(qū)域,由于工程實施導致工程北側水道的分流比和流速下降,南側水道主要呈淤積狀態(tài),最大淤積出現在C區(qū)域,最大淤積強度達到3郾2 m/a,由C至D區(qū)域,淤積強度不斷下降,D區(qū)域下降到大約1郾8 m/a。在工程的圍堤尾部E區(qū)域,由于圍堤的導流作用,流速下降,泥沙呈淤積狀態(tài),最大淤積強度達到3郾5 m/a。F區(qū)域由于處在橋墩后方,呈沖刷狀態(tài)。由于H區(qū)域位于限流潛堤的前端,底部流速下降,G區(qū)域受到限流潛堤的導流作用導致流速下降,故H、G區(qū)域呈淤積狀態(tài),最大淤積強度達到2郾0 m/a。工程南側水道由于分流比增加,流速增加,南側水道主要呈現出沖刷趨勢,L、M、N、I沿線水道沖刷強度一般在1郾3 m/a左右。如圖9所示,南側水道I區(qū)域由于靠近岸邊,工程前水深僅有1 m左右,但由于工程后流速的增加I區(qū)域受到強烈的水流沖刷,最大沖刷強度達到6 m/a。綜上所述,青草沙水庫工程的實施整體上使得工程南側航道呈現沖刷狀態(tài),北側航道呈現淤積狀態(tài),南側航道變化量小于北側航道,工程的實施對北側航道影響較大。
圖8 工程前水深、地形、潮量斷面及沖淤點布置
圖9 工程后水庫周邊年平均沖淤強度分布(單位:m)
2.2.2 經驗公式法
2.2.2.1 王義剛公式簡介
王義剛等[11]結合竇國仁[14]的潮汐水流懸沙運動微分方程、羅肇森[15]的航道回淤計算公式和劉家駒等[8]以及樂培九等[10]的研究成果,提出了計算邊灘圍墾淤積的簡便方法[16],其表達式如下:
式中:P為第1年回淤強度;琢為沉降概率;n為1年中的潮數;棕為泥沙沉速;T為潮周期;S*1、S*2分別為工程前、后對應的全潮平均泥沙質量濃度。
一般來說,河口邊灘區(qū)域水深相對較淺,單個潮周期內的泥沙沉降距離遠大于邊灘區(qū)域內最大平均水深,即棕T>H,令棕T抑k0H,k0臆1,則式(6)可簡化為
式中:H為全潮平均水深;u1、u2分別為工程前、后全潮平均流速;k0為考慮波浪及其他因素的綜合影響系數;m為計算參數,大小隨著地區(qū)不同而變化。k0、m均由實測資料率定得出。
灘面的促淤標高一般都在平均潮位以下,灘面高程越接近平均潮位,其年回淤強度就越小,所以促淤工程建后第2年的泥沙回淤量必然小于第1年的回淤量。在外界來沙條件不變的情況下,可以認為當年回淤強度與去年回淤強度的比值同去年平均水深與前年平均水深的比值成正比[16],即
式中:Pi為當年的回淤強度;Pi+1來年的回淤強度;Hi為當年的全潮平均水深;Hi-1為去年的全潮平均水深;K為比例系數,工程應用中一般取K=1郾0。
當u1逸u2時,表示工程施工后水域的流速減小,水動力減弱,相應的挾沙力下降,該區(qū)域出現淤積。當u1 2.2.2.2 王義剛公式的改進 由于缺乏泥沙起動流速的資料,無法判別泥沙是否起動,且由于計算時采用全潮平均流速,可能出現平均流速小于臨界起動流速,而最大流速大于臨界起動流速的現象,因此本文采用折中的處理方式,在原公式中加入沖刷系數茁。該系數由李瑞杰等[18]提出,是一個與沉速和泥沙質量濃度相關的系數,而沉速和泥沙質量濃度可以通過實測數據率定得出,所以沖刷系數的引入并不需要通過試驗得到起動流速。沖刷系數茁的取值應在0~1之間,考慮到工程地區(qū)底床表層為淤泥層,處于沖淤頻繁的狀態(tài),因此沖刷系數應取值較大,本文采取試算方法得到沖刷系數取值0郾85較為適宜。在綜合考慮泥沙起動和長時間尺度的影響后,提出的改進方法計算公式如下: 式中駐P指單個潮周期的回淤強度。 改進方法沖淤計算結果見圖10(b)。與改進前的計算結果(圖10(a))相比,改進后的王義剛公式能夠進行沖刷計算,且沖刷和淤積計算結果與圖10(c)中實測值和圖9中數值模擬結果基本吻合。由此可知,在改進經驗公式時,加入適當的沖刷系數茁是可行的。結合董曉偉[19]的計算結果,認為茁取值在0郾70~0郾85之間是合理的,如在其他河口海岸地區(qū),需通過當地的地形演變資料對茁合理取值后參考應用。 a.應用MIKE 21軟件,采用三角形網格,建立了長江口區(qū)域平面二維水沙數學模型,從潮位、流速、流向、泥沙質量濃度過程線、地形沖淤等角度對其進行驗證。模擬結果顯示模型精度滿足研究需要,因此可利用此模型進行水庫工程對周邊水域沖淤影響的計算研究。 圖10 河床沖淤深度計算結果與實測值對比(單位:m) b.青草沙水庫工程建成后,水庫北側航道的分流比減小,南側航道的分流比增大,水庫工程的建設使得北側航道流速整體減小,南側航道流速整體增大。 c.工程實施后,整體上使得水庫南側航道呈現沖刷狀態(tài),北側航道呈現淤積狀態(tài),南側航道變化量小于北側航道,水庫工程前沿口門處發(fā)生嚴重沖刷,而尾部區(qū)域發(fā)生嚴重淤積,但對距工程較遠的區(qū)域的水下地形影響不大。 d.在王義剛公式中添加沖刷系數茁(茁值由實測資料率定),可使改進后的王義剛公式同時計算水庫工程建設后對沖刷和淤積的影響,且滿足工程應用精度要求。 [1]孫東坡,李國慶,朱太順,等.治河及泥沙工程[M].鄭州:黃河水利出版社,1999. [2]竇希萍,李來,竇國仁.長江口全沙數學模型研究[J].水利水運科學研究,1999(2):136鄄145.(DOU Xipin,LI Lai,DOU Guoren.Numerical model study on total sediment transport in the Yangtze River Estuary[J].Journal of Nanjing Hydraulic Research Institute,1999(2):136鄄145.(in Chinese)) [3]王兆華,杜景龍.長江口深水航道一、二期工程建設以來北槽河段的沖淤演變[J].海洋通報,2006(6):55鄄62.(WANG Zhaohua,DU Jinlong.Evolvement of the north passage after executing the Deep鄄Water Channel Project of Yangtze Estuary[J].Marine Science Bulletin,2006(6):55鄄62.(in Chinese)) [4]王艷姣,張鷹,陳榮.基于GIS長江口北港河段沖淤變化的可視化分析[J].河海大學學報(自然科學版),2006,34(2):223鄄226.(WANG Yanjiao,ZHANG Ying,CHEN Rong.GIS鄄based visualized analysis of scouring and siltation in North Channel of Changjiang Estuary[J].Journal of Hohai University(Natural Sciences),2006,34(2):223鄄226.(in Chinese)) [5]陸忠民,盧永金,宋少紅,等.青草沙水庫建設與長江口綜合整治關系研究[J].給水排水,2009(1):55鄄58.(LU Zhongmin,LU Yongjing,SONG Shaohong,et al.Study of Qingcaosha Reservoir construction and comprehensive regulation of the Yangtze Estuary relationship[J].Water&Wastewater Engineering,2009(1):55鄄58.(in Chinese)) [6]時連強,李九發(fā),應銘,等.長江口沒冒沙演變過程及其對水庫工程的響應[J].長江流域資源與環(huán)境,2006(4):458鄄466.(SHI Lianqiang,LI Jiufa,YING Ming,et al.Evolutional process in Momaosa of the Yangtze River Estuary and itsresponse to reservoirproject[J].Resources and Environment in the Yangtze Basin,2006(4):458鄄466.(in Chinese)) [7]陳界仁,陳國祥.黃河河口三角洲演變的數學模擬[J].河海大學學報(自然科學版),2001,29(1):103鄄106.(CHEN Jieren,CHEN Guoxiang.Mathematical simulation of the Yellow River mouth delta processes[J].Journal of Hohai University(Natural Sciences),2001,29(1):103鄄106.(in Chinese)) [8]劉家駒,張鏡潮.淤泥質海岸航道、港池淤積計算方法及其推廣應用:兼論連云港建西大堤后的泥沙問題[J].水利水運科學研究,1993(4):301鄄320.(LIU Jiaju,ZHANG Jinchao.The calculation method and its application of silt coastal waterway and harbor basin,and the sediment problem after the construction of the west bank of lianyungang.[J].Journal of Nanjing Hydraulic Research Institute,1993(4):301鄄320.(in Chinese)) [9]曹祖德,王桂芬.波浪掀沙、潮流輸沙的數值模擬[J].海洋學報,1993,15(1):107鄄118.(CAO Zude,WANG Guifen.Numerical simulation of wave lifting sand and tidal current sediment transport[J].Acta Oceanologica Sinica,1993,15(1):107鄄118.(in Chinese)) [10]樂培九,張華慶.河口海灣地區(qū)航道挖槽回淤的估算方法[J].水道港口,1993(3):1鄄11.(YUE Peijiu,ZHANG Huaqing.An estimate method of siltation for dredged channels in estuaries[J].Journal of Waterway and Harbor,1993(3):1鄄11.(in Chinese)) [11]WANG Yigang,LI Xi,LIN Xiang.Analysis on suspended sediment deposition rate for muddy coast of reclaimed land[J].China Ocean Engineering,2001(1):147鄄152. [12]WILLIAMS D J,PEIRCE T J.Experiments on certain aspects of sedimentation of estuarine muds[C]//[S.l.]:Institution of Civil Engineers,1966:391鄄402. [13]侯志強,楊華.黃驊港外航道驟淤分析[J].水道港口,2004,25(4):213鄄215.(HOU Zhiqiang,YANG Hua.Analysis of sudden siltation of outer channel of Huanghua Harbor[J].Journal of Waterway and Harbor,2004,25(4):213鄄215.(in Chinese)) [14]竇國仁.潮汐水流中的懸沙運動及沖淤計算[J].水利學報,1963(4):13鄄23.(DOU Guoren.Suspended sediment movement and scouring and silting calculation in tidal currents[J].Journal of Hydraulic Engineering,1963(4):13鄄23.(in Chinese)) [15]羅肇森.河口航道開挖后的回淤計算[J].泥沙研究,1987(2):13鄄20.(LUO Zhaosen.Computation of siltation in dredged channels in estuaries[J].Journal of Sediment Research,1987(2):13鄄20.(in Chinese)) [16]王義剛,林祥,馮衛(wèi)兵.計算淤泥質海岸圍墾多年回淤強度的一種簡便方法[J].河海大學學報(自然科學版),2000,28(6):100鄄102.(WANG Yigang,LIN Xiang,FENG Weibin.A simple method for calculating long鄄term deposition rate after construction of warping bank in estuary[J].Journal of Hohai University(Natural Sciences),2000,28(6):100鄄102.(in Chinese)) [17]韓其為.泥沙起動規(guī)律及起動流速[J].泥沙研究,1982(2):11鄄26.(HAN Qiwei.Characteristics of incipient sediment motion and incipient velocity[J].Journal of Sediment Research,1982(2):11鄄26.(in Chinese)) [18]李瑞杰,羅鋒,周華民,等.泥沙沉降機率和沖刷系數分析[J].泥沙研究,2010(1):63鄄66.(LI Ruijie,LUO Feng,ZHOU Huamin,et al.Analysison sediment setting probability and erosion coefficient[J].Journal of Sediment Research,2010(1):63鄄66.(in Chinese)) [19]董曉偉.中尺度海岸演變動力模型與行為模型對比研究[D].南京:河海大學,2014.3 結 論