陳偉倫王 偉
(1.河海大學(xué)港口海岸與近海工程學(xué)院,江蘇南京 210098;2.中交海洋投資控股有限公司,海南三亞 572000)
泥沙沖淤計算主要依賴以下4種方法:淤現(xiàn)場資料分析法;于物理模型試驗法;盂數(shù)值模擬法;榆經(jīng)驗公式法。目前多采用數(shù)值模型法和經(jīng)驗公式法。在河口治理工程對水沙動力和沖淤演變的影響方面,竇希萍等[2]研究了長江口深水航道治理工程實施后的航道回淤情況;王兆華等[3]分析了治理工程對北槽攔門沙河段的影響;王艷姣等[4]利用GIS技術(shù)對長江口北港河段的沖淤演變進行了可視化分析;陸忠民等[5]評價了青草沙水庫建設(shè)的可行性;時連強等[6]探討了沒冒沙形成及演變規(guī)律;陳界仁等[7]計算出三角洲淤積影響平均寬度和延伸的平均距離。而在河床泥沙回淤經(jīng)驗公式的研究方面,劉家駒等[8]提出了用于淤泥質(zhì)海岸航道開挖后泥沙回淤計算的方法;曹祖德等[9]根據(jù)水流歸槽和平衡泥沙質(zhì)量濃度的概念,推導(dǎo)出順岸式碼頭前沿的港池淤積計算公式;樂培九等[10]提出了適用于河口海灣地區(qū)挖槽回淤估算方法;王義剛等[11]提出了計算邊灘圍墾淤積的一種簡便方法。
本文通過建立大范圍無結(jié)構(gòu)三角網(wǎng)格的高精度數(shù)值模型,探究青草沙水庫工程建設(shè)前后周邊河口范圍水流變化和沖淤演變規(guī)律,以期為河口水庫工程的建設(shè)提供更為可靠的工程咨詢,并對只能應(yīng)用于計算回淤強度的王義剛公式進行改進,使其可以同時計算河床沖淤。
依托MIKE21水動力模型軟件中HD模塊及其附屬的MT模塊來建立高精度的數(shù)值模型,其連續(xù)方程:
式中:x、y為右手Cartesian坐標(biāo)系;t為時間;濁為水位;h為總水深;d為靜水深;u、v分別為流速在x、y方向上的分量;Pa為當(dāng)?shù)卮髿鈮?籽為水密度,籽0為參考水密度;f為 Coriolis 參量;Sxx、Sxy、Syx、Syy為輻射應(yīng)力分量;Txx、Txy、Tyx、Tyy為水平黏滯應(yīng)力項;子sx、子sy分別為x和y方向上的表面切應(yīng)力分量;子bx、子by分別為x和y方向上的底部切應(yīng)力分量;S為源匯項;us和vs分別為x和y方向源匯項水流流速;Dx、Dy分別為x、y方向上的泥沙擴散系數(shù);Fs為泥沙沖淤函數(shù);酌d為床沙干容重;濁b為海底床面的豎向位移(即沖淤變化量);贅為地球自轉(zhuǎn)角速率,贅=0郾729伊10-4s-1;漬為地理緯度。
在缺少實測現(xiàn)場干容重資料的情況下,可以通過泥沙中值粒徑公式[7]估算。對于長江口區(qū)域,干密度取值相對穩(wěn)定,酌d取值720 kg/m3。經(jīng)過數(shù)值模型的反復(fù)調(diào)試,絮凝沉降的相關(guān)控制參數(shù)取值為:泥沙密度籽s取值2650 kg/m3,絮凝開始發(fā)生時的泥沙質(zhì)量濃度(含沙量,下同)cfloc取值0郾01 kg/m3,阻礙沉降開始發(fā)生時的泥沙質(zhì)量濃度chin取值10 kg/m3,凝膠點泥沙質(zhì)量濃度cget取值50 kg/m3。由Richardson鄄Zaki公式[12]結(jié)合侯志強等[13]的觀點,沉速系數(shù)棕取為0郾0005 m/s。選取臨界淤積切應(yīng)力在0郾10~0郾25Pa之間,泥沙臨界沖刷切應(yīng)力在 0郾15~0郾5 Pa之間。粗糙高度選用模型推薦值0郾001 m。
本文模型涵蓋長江口和杭州灣大范圍區(qū)域,上邊界取在江陰驗潮站,南邊界至浙江象山縣,北邊界取在江蘇如東縣中部,東部外海開邊界位于123毅34憶E,從北至南約為50~60 m等深線位置。模型東西長約310 km,南北長約為311 km,整個模型面積為5郾75伊104km2。此外,在大模型內(nèi)嵌套建立了僅包含長江口區(qū)域的小模型,北邊界、南邊界及東邊界為水邊界,給定潮位過程,該區(qū)域潮流受東中國海潮波系統(tǒng)控制,開邊界的潮位過程由全球潮波數(shù)學(xué)模型提供,上游邊界條件采用江陰潮位站實測潮位資料。
考慮到實際研究區(qū)域地形的復(fù)雜性,模型網(wǎng)格采用非結(jié)構(gòu)三角形網(wǎng)格,網(wǎng)格生成考慮地形邊界和潮流的走向,整個計算區(qū)域包括個28232個網(wǎng)格節(jié)點,54028個單元,網(wǎng)格間距在200~5 000 m之間,研究區(qū)域水下地形圖及網(wǎng)格如圖1所示。本文模型以1985國家高程基準(zhǔn)為計算基面,坐標(biāo)系統(tǒng)為北京1954坐標(biāo)系,中央經(jīng)線123毅。
圖1 研究區(qū)域水下地形及計算網(wǎng)格
1.4.1 潮流模型驗證 選取2002年3月1—9日大、小潮全潮水文測驗的結(jié)果進行驗證。各驗潮站、測流垂線及流量斷面如圖2所示,部分潮位站大潮驗證結(jié)果如圖3所示,部分測點大潮驗證結(jié)果如圖4所示,北港、南港斷面大潮驗證結(jié)果如圖5所示。本文約定落潮流速及落潮流量為正值。
1.4.2 懸沙模型驗證
圖6反映了2002年3月1日0時至10日0時Y6和Z9測點大、小潮過程中泥沙質(zhì)量濃度計算值與實測值的比較,表明所建模型基本能反映不同潮型下泥沙質(zhì)量濃度隨時間變化的情況,計算結(jié)果與實測值基本接近,滿足后續(xù)工作的精度需要。
1.4.3 河床沖淤驗證
基于三角形無結(jié)構(gòu)網(wǎng)格建立水流泥沙模型和潮流泥沙數(shù)值模擬,利用2007年12月至2009年10月實測水下地形進行沖淤演變計算驗證,見圖7。
由圖7可以看出,計算得到的沖淤區(qū)域及沖淤強度與實測結(jié)果基本吻合,表明數(shù)值模型中相關(guān)參數(shù)的取值合理,數(shù)值模擬的計算結(jié)果滿足工程研究的需要。
圖2 長江口水文測驗布置示意圖
圖3 高橋站大潮驗證結(jié)果
圖4 X3測點流速及流向大潮驗證結(jié)果
圖5 部分?jǐn)嗝媪髁看蟪彬炞C結(jié)果
圖6 部分測點泥沙質(zhì)量濃度大、小潮驗證結(jié)果
圖7 模型計算沖淤結(jié)果與實測地形對比(單位:m)
由于工程的影響,工程區(qū)域附近的流速和流向均發(fā)生了改變。根據(jù)實測資料,水庫建成后,長興島南側(cè)航道流速呈增大趨勢,小部分區(qū)域流速下降;水庫工程前沿取水口處,水流流向發(fā)生改變,方向沿著水庫圍堤向北側(cè)航道偏轉(zhuǎn),同時流速增加;水庫工程北側(cè)航道圍堤附近流速微增,圍堤較遠(yuǎn)處流速下降;水庫工程尾部流速下降,工程對水流的影響范圍呈扇形向外擴散,距離工程區(qū)域越近,流速和流向改變越明顯。南側(cè)航道的分流比增加,水庫工程的建設(shè)使得北側(cè)航道流速整體下降,南側(cè)航道流速整體上升。
2.2.1 數(shù)值模擬法
水庫工程的實施對周邊泥沙沖淤強度分布的影響主要集中在水庫工程的兩側(cè)區(qū)域,分別布置了15個采樣點(圖8)進行定量分析。青草沙水庫工程實施后,工程對水下地形沖淤影響如圖8和圖9所示:水庫口門A區(qū)域至B區(qū)域,由于工程圍堤前沿的挑流作用,口門前圍堤區(qū)域水流流向急轉(zhuǎn),流速增加,在工程圍堤拐角處B區(qū)域出現(xiàn)沖刷,沖刷呈現(xiàn)由口門到工程圍堤拐角處增大的趨勢,沖刷量最大達5郾0 m/a左右。C區(qū)域到D區(qū)域以及K區(qū)域,由于工程實施導(dǎo)致工程北側(cè)水道的分流比和流速下降,南側(cè)水道主要呈淤積狀態(tài),最大淤積出現(xiàn)在C區(qū)域,最大淤積強度達到3郾2 m/a,由C至D區(qū)域,淤積強度不斷下降,D區(qū)域下降到大約1郾8 m/a。在工程的圍堤尾部E區(qū)域,由于圍堤的導(dǎo)流作用,流速下降,泥沙呈淤積狀態(tài),最大淤積強度達到3郾5 m/a。F區(qū)域由于處在橋墩后方,呈沖刷狀態(tài)。由于H區(qū)域位于限流潛堤的前端,底部流速下降,G區(qū)域受到限流潛堤的導(dǎo)流作用導(dǎo)致流速下降,故H、G區(qū)域呈淤積狀態(tài),最大淤積強度達到2郾0 m/a。工程南側(cè)水道由于分流比增加,流速增加,南側(cè)水道主要呈現(xiàn)出沖刷趨勢,L、M、N、I沿線水道沖刷強度一般在1郾3 m/a左右。如圖9所示,南側(cè)水道I區(qū)域由于靠近岸邊,工程前水深僅有1 m左右,但由于工程后流速的增加I區(qū)域受到強烈的水流沖刷,最大沖刷強度達到6 m/a。綜上所述,青草沙水庫工程的實施整體上使得工程南側(cè)航道呈現(xiàn)沖刷狀態(tài),北側(cè)航道呈現(xiàn)淤積狀態(tài),南側(cè)航道變化量小于北側(cè)航道,工程的實施對北側(cè)航道影響較大。
圖8 工程前水深、地形、潮量斷面及沖淤點布置
圖9 工程后水庫周邊年平均沖淤強度分布(單位:m)
2.2.2 經(jīng)驗公式法
2.2.2.1 王義剛公式簡介
王義剛等[11]結(jié)合竇國仁[14]的潮汐水流懸沙運動微分方程、羅肇森[15]的航道回淤計算公式和劉家駒等[8]以及樂培九等[10]的研究成果,提出了計算邊灘圍墾淤積的簡便方法[16],其表達式如下:
式中:P為第1年回淤強度;琢為沉降概率;n為1年中的潮數(shù);棕為泥沙沉速;T為潮周期;S*1、S*2分別為工程前、后對應(yīng)的全潮平均泥沙質(zhì)量濃度。
一般來說,河口邊灘區(qū)域水深相對較淺,單個潮周期內(nèi)的泥沙沉降距離遠(yuǎn)大于邊灘區(qū)域內(nèi)最大平均水深,即棕T>H,令棕T抑k0H,k0臆1,則式(6)可簡化為
式中:H為全潮平均水深;u1、u2分別為工程前、后全潮平均流速;k0為考慮波浪及其他因素的綜合影響系數(shù);m為計算參數(shù),大小隨著地區(qū)不同而變化。k0、m均由實測資料率定得出。
灘面的促淤標(biāo)高一般都在平均潮位以下,灘面高程越接近平均潮位,其年回淤強度就越小,所以促淤工程建后第2年的泥沙回淤量必然小于第1年的回淤量。在外界來沙條件不變的情況下,可以認(rèn)為當(dāng)年回淤強度與去年回淤強度的比值同去年平均水深與前年平均水深的比值成正比[16],即
式中:Pi為當(dāng)年的回淤強度;Pi+1來年的回淤強度;Hi為當(dāng)年的全潮平均水深;Hi-1為去年的全潮平均水深;K為比例系數(shù),工程應(yīng)用中一般取K=1郾0。
當(dāng)u1逸u2時,表示工程施工后水域的流速減小,水動力減弱,相應(yīng)的挾沙力下降,該區(qū)域出現(xiàn)淤積。當(dāng)u1 2.2.2.2 王義剛公式的改進 由于缺乏泥沙起動流速的資料,無法判別泥沙是否起動,且由于計算時采用全潮平均流速,可能出現(xiàn)平均流速小于臨界起動流速,而最大流速大于臨界起動流速的現(xiàn)象,因此本文采用折中的處理方式,在原公式中加入沖刷系數(shù)茁。該系數(shù)由李瑞杰等[18]提出,是一個與沉速和泥沙質(zhì)量濃度相關(guān)的系數(shù),而沉速和泥沙質(zhì)量濃度可以通過實測數(shù)據(jù)率定得出,所以沖刷系數(shù)的引入并不需要通過試驗得到起動流速。沖刷系數(shù)茁的取值應(yīng)在0~1之間,考慮到工程地區(qū)底床表層為淤泥層,處于沖淤頻繁的狀態(tài),因此沖刷系數(shù)應(yīng)取值較大,本文采取試算方法得到?jīng)_刷系數(shù)取值0郾85較為適宜。在綜合考慮泥沙起動和長時間尺度的影響后,提出的改進方法計算公式如下: 式中駐P指單個潮周期的回淤強度。 改進方法沖淤計算結(jié)果見圖10(b)。與改進前的計算結(jié)果(圖10(a))相比,改進后的王義剛公式能夠進行沖刷計算,且沖刷和淤積計算結(jié)果與圖10(c)中實測值和圖9中數(shù)值模擬結(jié)果基本吻合。由此可知,在改進經(jīng)驗公式時,加入適當(dāng)?shù)臎_刷系數(shù)茁是可行的。結(jié)合董曉偉[19]的計算結(jié)果,認(rèn)為茁取值在0郾70~0郾85之間是合理的,如在其他河口海岸地區(qū),需通過當(dāng)?shù)氐牡匦窝葑冑Y料對茁合理取值后參考應(yīng)用。 a.應(yīng)用MIKE 21軟件,采用三角形網(wǎng)格,建立了長江口區(qū)域平面二維水沙數(shù)學(xué)模型,從潮位、流速、流向、泥沙質(zhì)量濃度過程線、地形沖淤等角度對其進行驗證。模擬結(jié)果顯示模型精度滿足研究需要,因此可利用此模型進行水庫工程對周邊水域沖淤影響的計算研究。 圖10 河床沖淤深度計算結(jié)果與實測值對比(單位:m) b.青草沙水庫工程建成后,水庫北側(cè)航道的分流比減小,南側(cè)航道的分流比增大,水庫工程的建設(shè)使得北側(cè)航道流速整體減小,南側(cè)航道流速整體增大。 c.工程實施后,整體上使得水庫南側(cè)航道呈現(xiàn)沖刷狀態(tài),北側(cè)航道呈現(xiàn)淤積狀態(tài),南側(cè)航道變化量小于北側(cè)航道,水庫工程前沿口門處發(fā)生嚴(yán)重沖刷,而尾部區(qū)域發(fā)生嚴(yán)重淤積,但對距工程較遠(yuǎn)的區(qū)域的水下地形影響不大。 d.在王義剛公式中添加沖刷系數(shù)茁(茁值由實測資料率定),可使改進后的王義剛公式同時計算水庫工程建設(shè)后對沖刷和淤積的影響,且滿足工程應(yīng)用精度要求。 [1]孫東坡,李國慶,朱太順,等.治河及泥沙工程[M].鄭州:黃河水利出版社,1999. [2]竇希萍,李來,竇國仁.長江口全沙數(shù)學(xué)模型研究[J].水利水運科學(xué)研究,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]王兆華,杜景龍.長江口深水航道一、二期工程建設(shè)以來北槽河段的沖淤演變[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].河海大學(xué)學(xué)報(自然科學(xué)版),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]陸忠民,盧永金,宋少紅,等.青草沙水庫建設(shè)與長江口綜合整治關(guān)系研究[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ā),應(yīng)銘,等.長江口沒冒沙演變過程及其對水庫工程的響應(yīng)[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]陳界仁,陳國祥.黃河河口三角洲演變的數(shù)學(xué)模擬[J].河海大學(xué)學(xué)報(自然科學(xué)版),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]劉家駒,張鏡潮.淤泥質(zhì)海岸航道、港池淤積計算方法及其推廣應(yīng)用:兼論連云港建西大堤后的泥沙問題[J].水利水運科學(xué)研究,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]曹祖德,王桂芬.波浪掀沙、潮流輸沙的數(shù)值模擬[J].海洋學(xué)報,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].水利學(xué)報,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)兵.計算淤泥質(zhì)海岸圍墾多年回淤強度的一種簡便方法[J].河海大學(xué)學(xué)報(自然科學(xué)版),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]李瑞杰,羅鋒,周華民,等.泥沙沉降機率和沖刷系數(shù)分析[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].南京:河海大學(xué),2014.3 結(jié) 論