王 碧,杜興武,胡 成,羅 垚
(1.中國地質(zhì)大學(武漢)環(huán)境學院,湖北 武漢 430074;2.湖南省水利水電勘測設計研究總院,湖南 長沙 410001)
水庫浸沒是指水庫蓄水后,庫區(qū)周圍地下水受庫水頂托作用,排泄受阻,導致地下水水位壅高,使其可能接近或者高于地表而導致該地區(qū)土壤次生鹽堿化、沼澤化以及建筑物地基條件惡化等危害[1]。一般來說,在平原地區(qū)修建水庫,正常蓄水位一般限制在一級階地附近,也有不少水庫利用兩岸堤防抬高水頭。由于平原地區(qū)庫岸地形平坦開闊,地面高程與正常蓄水位相差不大,地下水埋深較淺,易形成浸沒問題,且浸沒面積一般較大,因此對平原地區(qū)水庫浸沒范圍進行預測很有必要[2]。浸沒問題是平原地區(qū)水庫勘察設計和建設中最重要的水文地質(zhì)問題,對于常見的沖積成因二元沉積結(jié)構的平原地區(qū),我國現(xiàn)行的《水利水電工程地質(zhì)勘察規(guī)范》(GB 50287—99)[3](以下簡稱《規(guī)范》)中采用的卡明斯基法在實際工程應用中存在諸多問題。為此,本文分析了利用卡明斯基法對二元結(jié)構平原地區(qū)水庫浸沒進行預測所存在的不足,并從黏性土的滲透性規(guī)律入手,在“啟動水力梯度”的基礎上,提出了適用于承壓含水層模型的折減系數(shù)法的計算公式,最后通過實例應用,對比分析了卡明斯基法、折減系數(shù)法和數(shù)值法3種預測方法對二元結(jié)構平原地區(qū)水庫浸沒范圍預測結(jié)果的適用性和準確性,以為類似工程水庫浸沒問題研究提供參考依據(jù)。
對于雙層結(jié)構地層,《規(guī)范》中采用卡明斯基法進行水庫浸沒范圍預測。
其計算公式如下:
(1)
根據(jù)上式計算得出的地下水水頭壅高值H,結(jié)合水動力學法可推導出地下水壅高水頭線方程[4]為
(2)
上式中:M為下部含水層(砂礫石)的厚度(m);K1為下部含水層(砂礫石)的滲透系數(shù)(m/d);K2為上部含水層(黏性土)的滲透系數(shù)(m/d);h1、H為起始斷面蓄水前、后地下水水位至含水層頂板的厚度(m);h2、y為計算斷面壅水前、后地下水水位至含水層頂板的厚度(m);Yx為壅水后任意斷面的地下水水位值(m)。
卡明斯基法以上層黏性土底板為計算基準面,水庫浸沒區(qū)計算示意圖見圖1。
圖1 水庫浸沒區(qū)計算示意圖Fig.1 Calculation diagram of the reservoir immersion zone
利用該方法對二元結(jié)構平原地區(qū)水庫浸沒進行預測存在以下不足:①卡明斯基法計算模式是地下水補給庫水,水庫蓄水后導致地下水水位整體抬高,但地下水水位仍大于水庫水位,整體的補給關系不變。而平原地區(qū)地形平坦、地面高程較小,平原地區(qū)修建的水庫水位一般在河流一級階地附近,甚至利用兩岸的堤壩抬高水頭,導致庫水位高于兩岸一級階地。水庫蓄水后庫水位往往高于地下水水位,補給關系發(fā)生了改變,由庫水補給地下水,導致地下水水位壅高;②分析卡明斯基法計算公式得出的地下水壅高水頭線(見圖1)可知,離水庫越遠地下水水位越高,其浸沒程度越大,在堤腳處浸沒程度最小。而工程實際情況是,在地形相差不大的情況下堤腳處浸沒程度最大,遠離水庫其浸沒程度逐漸減??;③平原地區(qū)二元結(jié)構地層中上層為滲透系數(shù)較小的黏土層,其厚度較大,而下層為滲透系數(shù)較大的砂礫層,兩者級差大于102,黏土層可視為弱透水層或相對隔水層,其阻水作用不能忽略,采用承壓含水層模型計算會更加符合實際。而基于潛水含水層模型的卡明斯基法計算公式顯然已不適用于這種情況,其忽略了上層黏性土的阻水性[5]。
以上三點原因是造成該規(guī)范提供的方法不適用于平原地區(qū)二元結(jié)構地層的主要原因,將會導致預測結(jié)果出現(xiàn)偏移,由于上層黏性土的特殊性,應當采取新的方法來預測水庫蓄水后的地下水水位壅高值,使得計算結(jié)果更加接近實際值。
1.2.1 考慮起始水力梯度(I0)
張忠胤[6]研究指出黏土中存在大量的結(jié)合水,其對水的流動起到阻礙作用,并提出了“起始水力梯度“的概念,即黏性土中當水力梯度很小時,不存在顯著的滲流,只有當水力梯度達到一定值后才會出現(xiàn)流動的水,這個值就是“起始水力梯度”(I0)。水在黏性土中的滲流運移具有較復雜的形態(tài),其滲透性規(guī)律可用下式表示:
v=K(I-I0)
(3)
式中:v為滲透流速;K為黏性土的滲透系數(shù);I為水力比降,即水力梯度;I0為起始水力梯度,即為v-I直線段延長在I軸線的截距,見圖2。
圖2 水在黏性土中的滲透規(guī)律曲線Fig.2 Curve of the water permeability in the cohesive soil
由公式(3)可知,當水力梯度I小于I0時,水在黏土中不產(chǎn)生顯著滲流現(xiàn)象;而當水力梯度I大于I0時,v-I曲線呈線性關系,因此I0是表征黏性土滲流特性的重要指標之一。
由于起始水力梯度的存在,承壓含水層中的地下水在頂托作用下逐漸向上覆的弱透水層滲透,在豎直方向上,弱透水層中的水頭損失與穿透厚度成正比例關系,該滲透規(guī)律已在部分工程實踐中得到了驗證[7]。如圖3所示,當鉆孔達到a點時,剛好能觀測到穩(wěn)定的地下水水位,隨著鉆孔加深至b、c點,地下水水位分別上升了ΔH1和ΔH2,當鉆孔揭穿黏土層到達下部含水層時,地下水水位上升至距黏土層底板以上處H0,顯然:
圖3 弱透水層中的含水帶示意圖Fig.3 Diagram of the water-bearing zone in the aquitard
ΔH1/L1=ΔH2/L2=ΔH3/L3
=(ΔH1+ΔH2+ΔH3)/(L1+L2+L3)
=(H0-T)/T=I
(4)
由于各點地下水水位是穩(wěn)定、靜止的,故v=0,由此得出I=I0。故可推導出任意斷面的黏土層中含水帶厚度T的表達式為
T=H0/(I0+1)
(5)
式中:H0為下層含水層的測壓水位高度;I0為黏土中的起始水力梯度;T為弱透水層中的含水帶厚度,其上界面即為水庫蓄水后地下水水位壅高位置。
1.2.2 承壓含水層模型
平原地區(qū)的二元結(jié)構地層中由于上層黏性土滲透系數(shù)很小且厚度較大,黏性土可視為弱透水層,下層可視為承壓含水層。因此,采用承壓含水層模型計算地下水水位壅高值更加符合實際情況[8]。
由于下層承壓的水頭不是一個定值,隨著距堤壩距離的增加,承壓的水頭呈線性降低,見圖4。由公式(5)可知,黏土層中的含水帶厚度與下層承壓水頭成正比例關系,因此該含水帶厚度也是一個變化值。為了計算出黏土層中含水帶厚度,只需求出下層承壓水頭的變化規(guī)律即可。
圖4 二元結(jié)構承壓含水層模型示意圖Fig.4 Diagram of the confined aquifer model with dual-formation structure
水庫蓄水前,在堤腳處通過鉆孔鉆穿上層黏土層,獲得砂礫層的承壓水頭H1;然后在計算剖面的最遠端距堤壩Lm處,同樣通過鉆孔獲得H2;利用水動力學中承壓水一維穩(wěn)定運動模型,可得出單寬流量的計算公式為
(6)
(7)
以上層黏性土底板為計算基準面。將公式(7)代入公式(5)中,可得到黏土層中含水帶厚度的計算公式為
(8)
通過上式即可計算出平原地區(qū)二元結(jié)構地層中水庫蓄水后上層黏性土中地下水水位的壅高值,從而得出水庫浸沒范圍。該黏土層中含水帶厚度T是個變化值,隨著距堤壩距離x的增加,T減小,相應的水庫浸沒程度降低,符合實際情況。
在預測二元結(jié)構平原地區(qū)水庫浸沒范圍時,因研究區(qū)域內(nèi)地下水水文地質(zhì)條件較復雜,地下水往往呈空間三維運動狀態(tài),滲流場的時空分布難以用解析方法計算。為了預測在庫水位抬升的影響下庫周地下水水位雍高及其所造成的水庫浸沒范圍,需要在分析研究區(qū)水文地質(zhì)條件的基礎上,建立水文地質(zhì)概念模型,并采用數(shù)值法求解。
利用數(shù)值法求解主要有以下4個步驟:一是通過對研究區(qū)的水文地質(zhì)條件進行合理概化后,建立符合研究區(qū)實際情況的水文地質(zhì)概念模型;二是將該水文地質(zhì)概念模型進行轉(zhuǎn)換,建立與之對應的地下水滲流數(shù)學模型;三是將該數(shù)學模型在數(shù)值模擬軟件中轉(zhuǎn)為可視化的地下水滲流數(shù)值模型。經(jīng)過上述轉(zhuǎn)換過程建立的地下水滲流數(shù)值模型是否能全面、客觀地表征研究區(qū)實際的水文地質(zhì)條件和特征,還需要對該模型進行識別與驗證,也稱之為“反演”,即數(shù)學運算中的解逆問題,就是利用水頭函數(shù)求解地下水均衡方程。由于水頭函數(shù)受地下水均衡場內(nèi)多個水文地質(zhì)條件控制,因而它是一個多元函數(shù),模型驗證時需對研究區(qū)內(nèi)各種水文地質(zhì)參數(shù)進行調(diào)整,具體做法是根據(jù)觀測孔的資料,反求水文地質(zhì)參數(shù);四是將“反演”算出的最優(yōu)參數(shù)值代入地下水滲流數(shù)值模型中,計算得到預測結(jié)果。
湘江長沙綜合樞紐工程位于湘江干流下游的長沙市望城縣境內(nèi),壩址選定在望城縣境內(nèi)的湘江蔡家洲,位于長沙市城區(qū)下游約20 km。壩址至上游長沙市區(qū)瀏陽河庫段均屬于垸內(nèi)型的庫區(qū),堤壩內(nèi)的平原區(qū)地勢低緩,地形平坦開闊,是水庫浸沒易發(fā)地帶,見圖5。
圖5 蘇托垸地區(qū)平面簡圖Fig.5 Plane diagram of the Sutuoyuan area
本文選擇水庫東岸蘇托垸地區(qū)為研究區(qū),區(qū)內(nèi)為典型的平原地區(qū)二元結(jié)構地層,上層為黏土層,下層為砂礫層,主要為第四系地層,其特點如下:
(4) 殘坡積地層(Qedl):巖性為黃褐色、紅褐色黏土、黏土夾碎石,可硬塑狀,主要分布在低山、丘陵等坡麓地帶。
由鉆孔資料得出黏土層底板高程為25 m,天然情況下,研究區(qū)段撈刀河天然水位為25.45 m,達到正常蓄水位29.7 m高度時,按照撈刀河河面的水力坡降值得出河水位為29.76 m。當撈刀河河水位為27.5 m時,堤壩處承壓水位為27.6 m,與河水位基本一致,在遠離堤壩L=2 000 m處的觀測孔內(nèi)承壓水位為27 m。由于承壓含水層與河水之間水力聯(lián)系非常緊密,且相距很近,故兩者水位近似相等。據(jù)現(xiàn)場壓水試驗結(jié)果,上層弱透水層的滲透系數(shù)K2為2.0×10-5cm/s,厚度為5 m,下層強透水層的滲透系數(shù)K1為1.5×10-3cm/s,厚度為7 m。
在研究區(qū)水庫浸沒評價時,本文根據(jù)蘇托垸地區(qū)的地形特點,并結(jié)合觀測孔布置情況,在遠離堤壩L=2 000 m處選取具有代表性的剖面(見圖5),利用卡明斯基法、折減系數(shù)法和數(shù)值法3種方法對水庫浸沒范圍進行了預測。
本文利用卡明斯基法計算公式(1),可計算得出距離河流L=2 000 m處的地下水壅高水頭值H,通過計算得出Lm處的地下水水位值為31.29 m。
將該參數(shù)代入公式(2),化簡得到水庫蓄水后地下水壅高水頭線方程為
(9)
將利用卡明斯基法計算得出的地下水水位線繪制在研究區(qū)地形剖面上,并結(jié)合前述的水庫浸沒評價標準,繪制臨界地下水水位線,見圖6。
由圖6可見,利用卡明斯基法計算公式預測得到的水庫蓄水后蘇托垸地區(qū)地下水水位線接近地表,甚至在平坦及低洼地段溢出地表,而臨界地下水水位幾乎都高于地表,也就是說,水庫蓄水后整個蘇托垸地區(qū)均為浸沒區(qū)。
圖6 利用卡明斯基法預測得到的蘇托垸地區(qū)地下水水位剖面圖Fig.6 Cross section of groundwater level using Kamenski method in Sutuoyuan area
本文通過現(xiàn)場鉆孔試驗及長觀孔數(shù)據(jù),得出起始水力梯度I0=0.36,將所有參數(shù)代入折減系數(shù)法計算公式(8),可計算得到任意斷面的黏土層中含水帶厚度,從而推導出水庫蓄水后上層黏土層中地下水壅高水位線方程:
(10)
因此水庫蓄水后地下水壅高水頭
Yx=T+25=28.5-2.2×10-4x
(11)
將利用折減系數(shù)法計算得出的地下水水位線繪制在研究區(qū)地形剖面上,并結(jié)合前述的水庫浸沒評價標準,繪制臨界地下水水位線,見圖7。
圖7 利用折減系數(shù)法預測得到的蘇托垸地區(qū)地下水水位剖面圖Fig.7 Cross section of groundwater level using the reduction factor method in Sutuoyuan area
由圖7可見,距離撈刀河越近,地下水水位壅高值越大,遠離堤壩,地下水水位降低,符合實際情況。此外,距離堤壩越近,水庫浸沒程度也越大,符合常理。
本文采用有限元數(shù)值模擬軟件FEFLOW對蘇托垸地區(qū)水庫浸沒區(qū)進行數(shù)值模擬預測。根據(jù)蘇托垸研究區(qū)實際條件,可將該地區(qū)地下水流概化為非均質(zhì)各向同性3D水流模型,邊界條件設置如下:東南西三個方向均以撈刀河為界,設為一類邊界類型,邊界上各位置的河水位依據(jù)水位站的實測數(shù)據(jù)以及河水位的平均水力坡降進行估算;北部邊界設為二類流量邊界;上邊界取為自由潛水浸潤面;底邊界取為隔水邊界。大氣降雨和潛水蒸發(fā)作為源匯項處理。模擬區(qū)在平面上共劃分為44 328個單元、29 175個節(jié)點(見圖8)。水文地質(zhì)參數(shù)按照巖性分層賦值。模型運行期為10 a。為了節(jié)省篇幅,詳細模擬過程在此省略,將另文闡述。
圖8 蘇托垸地區(qū)水庫浸沒區(qū)3D空間離散化示意圖Fig.8 3D space discretization diagram of the simulated reservior immersion zone in Sutuoyuan area
本文利用數(shù)值法通過數(shù)值模擬得到蘇托垸地區(qū)水庫蓄水1 a后、5 a后、10 a后的浸沒分區(qū)圖(見圖9至圖11)以及水庫蓄水10 a后地下水水位剖面圖(見圖12)。
圖9 蘇托垸地區(qū)水庫蓄水1 a后的浸沒分區(qū)Fig.9 Distribution of the immersed zone in Sutuoyuan area 1 a after the impoundment of the reservoir
圖10 蘇托垸地區(qū)水庫蓄水5 a后的浸沒分區(qū)Fig.10 Distribution of the immersed zone in Sutuoyuan area 5 a after the impoundment of the reservoir
圖11 蘇托垸地區(qū)水庫蓄水10 a后的浸沒分區(qū)Fig.11 Distribution of the immersed area in Sutuoyuan area 10 a after the impoundment of the reservoir
由圖9至圖12可見,水庫蓄水后,研究區(qū)內(nèi)地下水水位明顯雍高,并且部分區(qū)域地下水水位略高于撈刀河水位,具體表現(xiàn)如下:水庫蓄水初期,由于河水位的抬升,地表水補給模擬區(qū)地下水,致使區(qū)域地下水水位整體抬升,逐漸接近地表,產(chǎn)生浸沒災害;在水庫正常蓄水條件下,短期內(nèi)研究區(qū)內(nèi)出現(xiàn)大范圍的輕微浸沒災害,在近岸段浸沒程度小,而在遠岸段浸沒程度反而更大,主要原因可能是水庫蓄水前期遠岸段地下水水位本身就高,加之區(qū)內(nèi)地下水排泄受阻,因而更早出現(xiàn)浸沒現(xiàn)象;水庫蓄水5 a時,浸沒范圍大大增加,浸沒面積達70%,其中在研究區(qū)內(nèi)西南部的近岸段浸沒程度非常嚴重,東北部浸沒程度最弱;水庫蓄水10a后浸沒范圍基本穩(wěn)定,浸沒面積約占75%,與水庫蓄水5 a時浸沒情況相比,浸沒面積增加不顯著,主要集中在浸沒程度輕微區(qū)。因此,湘江長沙樞紐工程蓄水后,研究區(qū)內(nèi)近岸段浸沒程度最嚴重,遠岸段浸沒程度最輕;西南部基本都處于浸沒狀態(tài),而東北部浸沒程度最輕。
圖12 利用數(shù)值法預測得到的蘇托垸地區(qū)地下水水位剖面圖Fig.12 Cross section of groundwater level using the numerical method in Shutuoyuan area
通過對比3種方法預測得到的蘇托垸地區(qū)的地下水水位剖面圖,可以看出:卡明斯基法預測得到的地下水水位線幾乎全部高于地表面,其預測結(jié)果是水庫蓄水后蘇托垸地區(qū)均為浸沒區(qū);而折減系數(shù)法與數(shù)值法預測得到的地下水水位線相近,距離撈刀河越近,地下水水位壅高值越大,浸沒程度也越大,遠離堤壩,地下水水位降低,浸沒程度減小,這兩種方法預測得到的結(jié)果是蘇托垸部分地區(qū)為浸沒區(qū)。由此得出以下結(jié)論:
(1) 通過上述3種預測方法在具體工程實例中的應用,結(jié)果顯示卡明斯基法的預測結(jié)果偏大,計算得到的地下水水位值偏高超過約1 m。
(2) 在二元結(jié)構平原地區(qū)進行水庫浸沒預測時,需要考慮起始水力梯度I0的影響,而考慮起始水力梯度后的折減系數(shù)法預測得到的地下水水位線明顯低于卡明斯基法的預測結(jié)果,且浸沒水位變化明顯,其預測結(jié)果更加合理。
(3)I0的原始定義是“起始水力梯度”,對于這一概念目前尚存爭議,但它是折減系數(shù)法的理論基礎。為了驗證I0的存在,應當通過對已建成的水庫工程所引發(fā)的浸沒災害情況進行重訪調(diào)查,將實際浸沒情況與預測值進行對比,才能得到更加有說服力的證據(jù),并得出更加準確的計算公式或結(jié)論。
(4) 解析法(即折減系數(shù)法)主要適用于剖面狀地下水水位的計算,而數(shù)值法適用于浸沒災害的時空變化預測分析。由于二元結(jié)構平原地區(qū)地下水滲流規(guī)律較為復雜,因此采用數(shù)值法對該類地區(qū)水庫浸沒進行預測更加準確。
通訊作者:胡 成(1976—),男,博士,副教授,主要從事工程地質(zhì)、水文地質(zhì)和環(huán)境水文地質(zhì)方面的研究。E-mail:hu_cheng@cug.edu.cn