余明輝,王睿璞,,陳小齊,劉長(zhǎng)杰,劉畫眉
(1.武漢大學(xué)水資源與水電工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室,湖北 武漢 430072;2.廣東省水利水電科學(xué)研究院,廣東 廣州 510635; 3.深圳市城市公共安全技術(shù)研究院,廣東 深圳 518000;4.哈爾濱工業(yè)大學(xué)(深圳),廣東 深圳 518000)
水系連通性反映河道干支流、湖泊及其他濕地等水流的連續(xù)性,具有自然屬性和社會(huì)屬性[1-2]。水系連通研究和連通性定量評(píng)價(jià)可為防洪、水生態(tài)保護(hù)、水環(huán)境改善及水資源開發(fā)與分配等提供基礎(chǔ)支撐。近年來,國(guó)內(nèi)外學(xué)者對(duì)連通性概念機(jī)制、定量評(píng)價(jià)及其在城市水系規(guī)劃等方面應(yīng)用開展了大量研究[3-6]。關(guān)于河流水系連通性的研究成果層出不窮,孟慧芳等[7]對(duì)比分析現(xiàn)有連通性評(píng)價(jià)方法,提出基于水流阻力和水文連通函數(shù)的評(píng)價(jià)方法可較好地評(píng)價(jià)平原河網(wǎng)區(qū)河流連通性;崔廣柏等[8]認(rèn)為水系結(jié)構(gòu)連通性是水力連通性的基礎(chǔ),水系連通條件的提高可改善水環(huán)境;王延貴等[9]構(gòu)建了包括河流邊界、水流、泥沙及生態(tài)環(huán)境等在內(nèi)的河道(湖庫(kù))連通評(píng)價(jià)體系。在筑閘建壩、河道采砂等人為干預(yù)下,部分地區(qū)的河網(wǎng)水系結(jié)構(gòu)變化和功能退化,出現(xiàn)河流連通性衰退的局面[10]。此外,河流水系的水位和流量等存在明顯的季節(jié)性差異和年際變化,河流水系連通性也會(huì)動(dòng)態(tài)變化[2]。在河口三角洲水文連通性的研究中,Tejedor等[11]分析了三角洲河網(wǎng)脆弱性,量化了流量在河網(wǎng)區(qū)河段之間的分配過程與動(dòng)態(tài)變化;Hiatt和Passalacqua[12]較為全面地考慮三角洲河道和島嶼,研究發(fā)現(xiàn)島嶼灘地在三角洲水文連通性中具有重要作用。
20世紀(jì)80年代以來,東江干流及三角洲采砂活動(dòng)大規(guī)模興起,導(dǎo)致河道演變特性由緩慢淤積為主轉(zhuǎn)變?yōu)閯×移毡橄虑?使采砂區(qū)河段河槽容積大幅度擴(kuò)大,水深明顯增加,河道橫斷面愈向窄深發(fā)展,而采砂區(qū)以上河段河床則表現(xiàn)出明顯的溯源沖刷[13-14]。東江干流及三角洲地形劇烈變化引起河道水動(dòng)力條件的改變:一方面,20世紀(jì)60年代前東江三角洲南、北水道枯水期分流比為45%和55%,但90年代后分流比逆轉(zhuǎn)為60%和40%,且上游來流越小南支流分配越大[15];另一方面,河床下切導(dǎo)致東江干流及東江三角洲同流量級(jí)下水位均有大幅度的下降,但1996年后樊屋站、石龍站、大盛站受上下游河道地形不均勻下切及河床縱比降減小的綜合影響,3站水位不同程度的上升[16-17]。現(xiàn)有研究較多分析東江干流和東江三角洲河床演變以及水動(dòng)力特征的改變,但聚焦到該區(qū)域水系形態(tài)和連通性研究的成果較少。
本文梳理1999—2009年間東江干流博羅—石龍河段及東江三角洲河網(wǎng)區(qū)地形變化,基于圖論方法描述河流結(jié)構(gòu)連通性和河流輸運(yùn)連通性,借助數(shù)值模擬手段研究徑潮條件改變和地形變化對(duì)東江干流及三角洲河流連通性的綜合影響,研究成果以期為改善東江水系河流連通性及東江防洪、水利工程規(guī)劃和水生態(tài)環(huán)境的研究提供支撐。
東江發(fā)源于江西省尋烏縣,干流博羅至石龍河段長(zhǎng)約56 km,在廣東省東莞市石龍匯入珠江三角洲,石龍以下河網(wǎng)水系稱為東江三角洲(圖1)。東江三角洲集水面積約1 380 km2,在石龍分出東江北干流和東江南支流,水系縱橫交錯(cuò),經(jīng)大盛、麻涌、漳澎、泗盛注入獅子洋經(jīng)虎門入海。
圖1 東江干流及東江三角洲示意Fig.1 Map of the lower reaches of Dongjiang River and Dongjiang River Delta (DRD)
1972年以前東江干流及東江三角洲的采砂活動(dòng)較少;20世紀(jì)80年代初改革開放以來,東江干流及三角洲地區(qū)需砂量大增,主要采砂區(qū)為東江干流、東江北干流和東莞水道等;20世紀(jì)90年代以后東江干流采砂區(qū)逐年上移到達(dá)合竹洲,這一階段的采砂長(zhǎng)期處于無序、無度狀態(tài),1997年?yáng)|江流域主管機(jī)構(gòu)介入嚴(yán)格控制采砂,實(shí)施河砂分段招標(biāo)措施,其中東江干流博羅附近河段仍劃為招標(biāo)采砂的河段,加上偷采,采砂量仍較大。2004—2007年博羅以下河段采砂活動(dòng)進(jìn)一步嚴(yán)控,2009年以后采砂規(guī)模有序減小,直至停止[18](1972—2002年間主要采砂區(qū)及采砂量變化如表1)。
表1 東江干流及東江三角洲采砂統(tǒng)計(jì)[18-19]Table 1 Statistics of sand excavation in the lower reaches of Dongjiang River and DRD
1997—2008年間,東江典型河床斷面博羅斷面6.50 m高程以下河床平均下切3.62 m,年均下切0.33 m,深泓高程降低6.69 m,10 m水位下過水?dāng)嗝婷娣e增加約59%。根據(jù)賈良文等[18]、Luo等[20]和胡曉張等[21]統(tǒng)計(jì)調(diào)查結(jié)論,1997—2009年間東江干流河床平均高程下切4.03 m,年均下切0.34 m;東江北干流石龍—新塘河段河床平均高程下切2.72 m,年均下切0.23 m;新塘—大盛河段河床平均高程下切0.07 m,但距大盛9 km的河段有沖有淤,河床平均高程抬高0.20 m;1997—2003間中堂水道河床平均高程增大1.35 m、厚街水道河床平均高程減少0.35 m。調(diào)查表明采砂主要采挖河床沙,多在河槽中進(jìn)行,而近岸淺灘挖沙較少,且主要采砂區(qū)位于東江干流、東江北干流、東莞水道等河道[22]。以“1999年地形”為基礎(chǔ),由河床平均高程年均變化累加得到1999—2009年間主槽河床高程變化,結(jié)果如表2所示,另外,潢涌水道、大汾河、厚街水道河底平均高程變化按參考值選取,入??谔幐浇佣魏哟财骄俜e0.20 m,以此推算得到2009年?yáng)|江博羅以下河道斷面地形(簡(jiǎn)稱“2009年地形”),河道深泓和典型斷面變化見圖2、圖3。
表2 1999—2009年?yáng)|江河床平均高程變化推算值 mTable 2 Calculated values of average channel elevation change in the lower reaches of Dongjiang River and DRD from 1999 to 2009
圖2 1999年和2009年?yáng)|江河道深泓高程沿程變化Fig.2 Changes of thalweg elevation in the lower reaches of Dongjiang River and DRD from 1999 to 2009
圖3 1999年和2009年?yáng)|江干流及東江三角洲典型河床斷面Fig.3 Typical cross sections in the lower reaches of Dongjiang River and DRD in 1999 and 2009
建立河網(wǎng)區(qū)一維河網(wǎng)水動(dòng)力數(shù)學(xué)模型,模型的上游邊界博羅站采用實(shí)測(cè)流量過程,增江麒麟咀站處的徑流過程在新家鋪?zhàn)鼽c(diǎn)源入?yún)R處理,下游邊界大虎站采用實(shí)測(cè)潮位過程(圖4)。高程基準(zhǔn)面采用珠江基面,建模地形采用1999年實(shí)測(cè)地形。
圖4 模型邊界及觀察斷面布置Fig.4 Model boundary and the location of observation sections
采用1998年6月21日至6月23日水文觀測(cè)資料進(jìn)行率定,得到東江河道糙率為0.023~0.031。采用同年6月25日至6月28日水文觀測(cè)資料進(jìn)行驗(yàn)證,圖5為石龍站水位驗(yàn)證結(jié)果,均方根誤差為0.09 m,未超過0.10 m[23],表明該模型可用來模擬東江干流及三角洲徑潮動(dòng)力過程。
圖5 石龍站水位模擬值與觀測(cè)值對(duì)比Fig.5 Comparison of simulated and measured water level at Shilong station
相比于枯季,洪季河流徑流較大,連通性較為顯著,選取洪季作為本文研究時(shí)期。由于缺乏與建模采用的1999年地形年份和變化后2009年地形相匹配的洪水實(shí)測(cè)資料,選取2場(chǎng)與地形年份接近的實(shí)測(cè)洪水(發(fā)生于1998年6月、2008年6月,簡(jiǎn)稱為“98·6”洪水、“08·6”洪水),洪潮過程如圖6所示。表3為2場(chǎng)洪水來流及下游潮位特性,“08·6”洪水兩站總來流量比“98·6”洪水增加約46%,但外海最大潮差減小0.04 m、平均潮差減小0.27 m。為探究不同徑潮條件和地形條件下東江干流及三角洲河網(wǎng)區(qū)結(jié)構(gòu)連通性和輸運(yùn)連通性的變化,計(jì)算過程中根據(jù)不同年份地形和不同徑潮條件,組合成4種工況:工況1為“1999年地形”+“98·6”徑潮過程實(shí)測(cè)資料;工況2為“1999年地形”+“08·6”徑潮過程實(shí)測(cè)資料;工況3為“2009年地形”+“08·6”徑潮過程實(shí)測(cè)資料;工況4為“2009年地形”+“98·6”徑潮過程實(shí)測(cè)資料。
圖6 洪水流量過程及大虎潮位Fig.6 Flood process and tidal water level in Dahu station
表3 東江流域“98·6”洪水、“08·6”洪水水文特性Table 3 Characteristics of flood “98·6”and “08·6”in the Dongjiang River basin
選用余水位、余水位比降及潮差描述東江三角洲徑潮動(dòng)力特性。余水位指1個(gè)太陰日25 h內(nèi)的平均水位,即采用傅里葉變換方法對(duì)逐時(shí)水位序列進(jìn)行以25 h為周期濾波后的低頻成分;余水位比降為單位距離下余水位的增量,反映潮汐傳播時(shí)遭受的非線性摩擦阻力,余水位比降越大對(duì)潮動(dòng)力的抑制作用越顯著;潮差指1個(gè)潮周期內(nèi),去除余水位影響后(即濾波后的高頻成分)高潮位與低潮位的差值[24-25]。
(1)結(jié)構(gòu)連通性。圖論方法描述水系的形態(tài)格局時(shí)把河流的分汊匯流點(diǎn)概化為節(jié)點(diǎn),河道概化為連接,水系概化為節(jié)點(diǎn)和連接組成的網(wǎng)狀結(jié)構(gòu)。
可選路徑數(shù)(Number of Alternative Paths,Nap)描述河道平面形態(tài)上的成環(huán)程度。子河網(wǎng)為從上游至下游出口的沿程河道組成的輸水通道。Nap為子河網(wǎng)中從上游進(jìn)口至下游出口的路徑數(shù)目,其值愈小,反映子河網(wǎng)中分汊、入?yún)R節(jié)點(diǎn)較少,成環(huán)發(fā)育程度較差[11]。
連接共享系數(shù)(Linking Sharing Index,ILS)描述連接在子河網(wǎng)間的共享程度,計(jì)算式為
(1)
式中:k為子河網(wǎng)的序號(hào);連接(ij)表示水流從節(jié)點(diǎn)i向節(jié)點(diǎn)j;Vk為子河網(wǎng)Sk的節(jié)點(diǎn)數(shù)目;bij為子河網(wǎng)Sk中連接(ij)在所有子河網(wǎng)中的重復(fù)次數(shù)[11]。
節(jié)點(diǎn)連接率(β)描述河網(wǎng)中每個(gè)節(jié)點(diǎn)與其他節(jié)點(diǎn)連接難易程度[10,26],計(jì)算式為
(2)
式中:L為河道數(shù),0 (2)輸運(yùn)連通性。輸運(yùn)連通性包含水體的流動(dòng)性和河流的連續(xù)性2個(gè)基本要素[1]。河道流量越大反映流動(dòng)性越好、連續(xù)性越優(yōu)、水量傳輸能力越強(qiáng)[27]。連通度因子計(jì)算式為 ωij=1-exp(-qij) (3) 式中:ωij為河道輸運(yùn)連通度,值愈小代表節(jié)點(diǎn)i和節(jié)點(diǎn)j之間連通程度愈差;qij為河道時(shí)均流量,m3/s。將“98·6”洪水、“08·6”洪水進(jìn)行同倍比縮小,縮小比取1960—2009年博羅站洪季多年平均流量1 053 m3/s。式(3)改寫為 ωij=1-exp(-qij/1 053) (4) 將連通度因子ωij代入鄰接矩陣A=(aij)n×n得到加權(quán)鄰接矩陣W,則判斷矩陣 (5) (6) (7) 式中:D為河網(wǎng)連通矩陣;dij為節(jié)點(diǎn)i和節(jié)點(diǎn)j之間的平均連通度;nij為節(jié)點(diǎn)i和節(jié)點(diǎn)j之間各不同長(zhǎng)度的連接的數(shù)目,nij≥1,nij∈N+;IDC為河網(wǎng)輸運(yùn)連通度。 徑潮動(dòng)力特性的變化是輸運(yùn)連通性變化的內(nèi)因。4種計(jì)算工況觀察斷面處(圖4)的余水位和潮差指標(biāo)如表4和圖7所示。分別對(duì)比工況1和工況2、工況3和工況4,均可以發(fā)現(xiàn):在同一地形條件下,當(dāng)上游來流增大46%且大虎潮差減小0.27 m時(shí),東江干流及三角洲河道余水位均增加而潮差均減小,且徑流控制作用較強(qiáng)的三角洲上游河道余水位增幅大于潮流控制作用較強(qiáng)的三角洲下游河道,而下游河道潮差的降幅則大于上游河道。 表4 計(jì)算工況余水位和潮差變化Table 4 Changes of residual water level and tidal range of model cases 分別對(duì)比工況1與工況4、工況2與工況3,當(dāng)上下游遭遇的徑潮過程相同時(shí),河道地形變化后,東江干流及三角洲河道余水位均降低,變化幅度沿程分布表現(xiàn)為位于徑流控制區(qū)愈上游的河道余水位增幅愈大,而潮差變化與河道地形變化相關(guān),大盛潮差減小0.22 m、泗盛潮差僅增大0.07 m,這與當(dāng)?shù)刂鞑酆哟哺叱唐骄俜e0.2 m的統(tǒng)計(jì)結(jié)果一致,除此之外河道潮差整體上增大,且潮流控制區(qū)愈下游的河道潮差增幅愈大。 工況1與工況3對(duì)比可探究地形變化和上下游不同徑潮遭遇對(duì)東江干流及三角洲水動(dòng)力的綜合影響,由表4和圖7可知,東江干流徑潮特性受地形變化影響較大,而東江三角洲河網(wǎng)區(qū)受徑潮過程的影響較大,表現(xiàn)為東江干流河道沿程余水位降低和潮差增大,東江北干流、東江南支流河道沿程余水位增大和潮差減小。 圖7 余水位和潮差變化Fig.7 Change of residual water level and tidal range 東江干流及東江三角洲可劃分為4個(gè)子河網(wǎng),分別用符號(hào)S1、S2、S3、S4表示(圖8),頂點(diǎn)均為博羅,出口分別為大盛、麻涌、漳澎、泗盛。基于圖論將東江干流及三角洲河網(wǎng)概化為有向圖(見圖9),其中節(jié)點(diǎn)23個(gè),連接31條。東江干流及三角洲各河道連接重復(fù)次數(shù)bij值如圖10所示。表5為結(jié)構(gòu)連通性描述指標(biāo)Nap、ILS和β計(jì)算值,子河網(wǎng)S4的Nap值最大,表明河網(wǎng)內(nèi)有諸多河道分汊、入?yún)R點(diǎn),河網(wǎng)環(huán)形輸水通道較多、成環(huán)發(fā)育程度較好,因而增加了水流入海的路徑;子河網(wǎng)S1、S2的Nap均為1,表明其流路唯一,不存在環(huán)形輸水通道,說明該子河網(wǎng)成環(huán)發(fā)育程度不高。由圖10可見,bij=4的河道僅3條,占總河道數(shù)的9.4%;不存在bij=3的河道;bij=2的河道共10條,占總河道數(shù)的31.3%;bij=1的河道共19條,占59.3%,多為子河網(wǎng)S4的輸水通道。由式(1)計(jì)算子河網(wǎng)連接共享系數(shù)值為0.28~0.55,其中,子河網(wǎng)S4的ILS值最小,說明S4子河網(wǎng)結(jié)構(gòu)上共享性較低;子河網(wǎng)S1、S2、S3的ILS均高于0.50,表明共享性較高。由式(2)計(jì)算子河網(wǎng)節(jié)點(diǎn)連接率值變化范圍為0.83~1.42,東江干流及三角洲的平均節(jié)點(diǎn)連接率β=1.39,表示平均每個(gè)節(jié)點(diǎn)有1~2個(gè)連接線,說明河網(wǎng)有較好的連接途徑。計(jì)算結(jié)果表明,整體上東江干流及三角洲河網(wǎng)成環(huán)發(fā)育程度不高,但河道共享性較高,有較好的連通途徑,總體結(jié)構(gòu)連通性較好;不同徑潮條件和地形條件下,結(jié)構(gòu)連通性未受影響。 圖8 東江干流及三角洲子河網(wǎng)示意Fig.8 Schematic of sub-network of the lower reaches of Dongjiang River and DRD 圖9 東江干流及東江三角洲概化有向圖Fig.9 Directed map of the lower reaches of Dongjiang River and DRD圖10 連接重復(fù)次數(shù)Fig.10 Value of linking sharing frequency bij 表5 東江干流及三角洲結(jié)構(gòu)連通性指標(biāo)Table 5 Topologic connectivity metrics of the lower reaches of Dongjiang River and DRD 4種工況下輸運(yùn)連通度計(jì)算結(jié)果如圖11和表6所示。分別對(duì)比工況1和工況2、工況3和工況4,其結(jié)果表明:在同一地形條件下,相比于上下游遭遇“98·6”徑潮過程,在“08·6”徑潮過程下河道輸運(yùn)連通度和河網(wǎng)輸運(yùn)連通度均增大,表明來流增大及外海潮動(dòng)力減弱使東江干流博羅—石龍河段及東江三角洲河網(wǎng)內(nèi)河道流量普遍增大,余水位比降增大反映潮波上溯的阻力增大,綜合作用下東江干流及三角洲河道潮動(dòng)力減弱,因此,來流較大、外海潮動(dòng)力較弱的“08·6”徑潮過程條件下河流輸運(yùn)連通性增強(qiáng)。分別對(duì)比工況1和工況4、工況2和工況3,其結(jié)果表明:在東江三角洲遭遇同一徑潮過程時(shí),地形變化對(duì)子河網(wǎng)輸運(yùn)連通度的影響無顯著特征,但河網(wǎng)整體輸運(yùn)連通度降低,這表明東江干流博羅—石龍河段及三角洲地形不均勻下切引起河道流量重新分配,且余水位比降減小,反映潮波上溯的阻力減小,河道內(nèi)潮動(dòng)力增強(qiáng)、水量傳輸能力下降,說明地形下切使輸運(yùn)連通度減小、輸運(yùn)連通性減弱。對(duì)比工況1和工況3,工況3河道輸運(yùn)連通度變化范圍為-0.231 3~0.272 9,其中輸運(yùn)連通度減小的河道共4條,占全部河道數(shù)的12.9%,子河網(wǎng)輸運(yùn)連通度增大范圍為0.014 3~0.391 7,整體上河網(wǎng)東江干流及三角洲輸運(yùn)連通度增大0.017 6;且相比于地形下切使河網(wǎng)輸運(yùn)連通度減小了0.004 2~0.008 1,來流及外海潮動(dòng)力帶來的增大幅度在0.021 8~0.025 7,說明洪季徑潮條件對(duì)輸運(yùn)連通性的影響程度較大,洪水流量越大則輸運(yùn)連通性越優(yōu)。 圖11 東江干流及三角洲河道輸運(yùn)連通度示意Fig.11 Path dynamic connectivity in the lower reaches of Dongjiang River and DRD 表6 東江干流及三角洲輸運(yùn)連通性指標(biāo)Table 6 Dynamic connectivity of the lower reaches of Dongjiang River and DRD 考慮采砂等人類活動(dòng)干預(yù)后東江干流及三角洲河網(wǎng)地形變化,開展人為擾動(dòng)下地形變化和徑潮條件變化對(duì)洪季河流結(jié)構(gòu)連通性和輸運(yùn)連通性兩方面的影響研究,主要結(jié)論如下: (1)結(jié)構(gòu)連通性未受地形和徑潮條件變化的影響。東江干流及三角洲成環(huán)發(fā)育程度以子河網(wǎng)博羅—泗盛為最高,河道共享性最差;4個(gè)子河網(wǎng)均有較好的連通途徑,平均節(jié)點(diǎn)連接率為1.39,總體結(jié)構(gòu)連通性較好。 (2)洪季輸運(yùn)連通性受河道地形和徑潮條件變化的影響。河道地形下切使徑流動(dòng)力減弱,潮動(dòng)力增強(qiáng),水流下泄阻力增大,水量傳輸能力降低,輸運(yùn)連通度減小0.004 2~0.008 1,河網(wǎng)輸運(yùn)連通性減弱;而徑流增大及外海潮動(dòng)力減小導(dǎo)致輸運(yùn)連通度增大0.021 8~0.025 7,河網(wǎng)輸運(yùn)連通性增強(qiáng)。3 結(jié)果及討論
3.1 徑潮動(dòng)力特性變化分析
3.2 結(jié)構(gòu)連通性對(duì)地形和徑潮條件變化的響應(yīng)分析
3.3 輸運(yùn)連通性對(duì)地形和徑潮條件變化的響應(yīng)分析
4 結(jié) 論