袁 科
(中鐵第四勘察設(shè)計(jì)院集團(tuán)有限公司,湖北武漢 430063)
水庫水位變化和降雨會(huì)導(dǎo)致坡體穩(wěn)定性的下降。對(duì)于在這些條件下坡體穩(wěn)定性計(jì)算方法,首先要確定浸潤(rùn)線的位置,然后根據(jù)浸潤(rùn)線來確定滲透壓力,再進(jìn)行穩(wěn)定性分析。然而目前浸潤(rùn)線的確定大多根據(jù)經(jīng)驗(yàn),這樣可能造成治理工程的不安全。因此,浸潤(rùn)線的計(jì)算有必要進(jìn)一步研究。本文對(duì)在庫水位等速下降時(shí)岸坡地下水浸潤(rùn)線的公式進(jìn)行推導(dǎo)計(jì)算。
①含水層各向同性、均質(zhì),側(cè)向無限延伸,具有水平不透水層。
②潛水流為一維流,水位下降前原始潛水面水平。
③水位以v0等速下降。
依據(jù)以上假設(shè)條件,可得到各向同性土體的一維非穩(wěn)定滲流運(yùn)動(dòng)的基本微分方程[1]
(1)
上式為二階非線性偏微分方程,求其解析解,通常是采用簡(jiǎn)化方法,對(duì)其進(jìn)行線性化,然后再求解。簡(jiǎn)化方法一般是把含水層厚度H看作為一常量,即用始、末時(shí)段潛水流的平均厚度hm來代替,可得到簡(jiǎn)化的滲流運(yùn)動(dòng)方程為
(2)
(3)
式中k——滲透系數(shù)/(m/d);
hm——含水層的平均厚度/m;
μ——給水度或貯水率;
t——水位下降持續(xù)時(shí)間/d。
浸潤(rùn)線計(jì)算如圖1所示。初始時(shí)刻時(shí),即t=0,由假設(shè)條件②可知,區(qū)內(nèi)各點(diǎn)初始水位均為h0,0。設(shè)距岸坡x處t時(shí)刻的地下水位變幅為u(x,t)=h0,0-hx,t,t=0時(shí)刻時(shí)的水位變幅為:u(x,0)=0。由假設(shè)條件③,當(dāng)水位以v0勻速下降時(shí),在x=0處,t時(shí)刻地下水變幅為:u(0,t)=v0t;在x=∞處,可以認(rèn)為水位無變化,u(∞,t)=0。
圖1 浸潤(rùn)線計(jì)算
根據(jù)上述的邊界和初始條件,把上述水位下降的半無限含水層中地下水非穩(wěn)定滲流歸結(jié)為以下模型方程
(4)
(5)
由定積分的分部積分法可知式(5)左邊
(u(x,0)=0)
式(5)右邊
因此
(6)
將方程組(4)中第二個(gè)方程兩端同乘e-st并積分得
(7)
因此由以上變換可將模型方程組(4)變換為如下方程組
(8)
(9)
對(duì)式(9)進(jìn)行拉普拉斯(Laplace)逆變換可得
(10)
最終可得水位下降時(shí)浸潤(rùn)線隨時(shí)間變化曲線方程
h(x,t)=h0,0-u(x,t)
(11)
從R(λ)的計(jì)算式可以看出,R(λ)的計(jì)算十分復(fù)雜,需要積分才能求得,不利于工程應(yīng)用。為了得到便于應(yīng)用的表達(dá)式,對(duì)圖2所示的曲線進(jìn)行多項(xiàng)式擬和,得到的擬和公式如下
圖2 λ與R(λ)關(guān)系曲線
R(λ)=
(12)
于是可得到庫水位下降時(shí)浸潤(rùn)線計(jì)算的簡(jiǎn)化公式,其表達(dá)式為
hx,t=
(13)
為了驗(yàn)證所得簡(jiǎn)化公式的適用性,運(yùn)用3D-Modflow對(duì)前述公式進(jìn)行驗(yàn)證。Visual-Modflow是由加拿大Waterloo水文地質(zhì)公司在Modflow基礎(chǔ)上,應(yīng)用現(xiàn)代可視化技術(shù)研制而成,是目前世界上最通行且被各國(guó)同行一致認(rèn)可的可視化專業(yè)軟件系統(tǒng)。該系統(tǒng)自1994年8月首次在國(guó)際上公開發(fā)行以來,在工程建設(shè)、環(huán)境保護(hù)、城鄉(xiāng)發(fā)展規(guī)劃、水資源利用等許多領(lǐng)域的科研和生產(chǎn)中得到了廣泛的應(yīng)用[3~5],后期經(jīng)過不斷完善,目前已經(jīng)成為最為普及的地下水運(yùn)動(dòng)數(shù)值模擬程序。
為了方便驗(yàn)證,采用Visual-Modflow分析軟件,對(duì)不同滲透系數(shù)、水位下降速度及給水度等因素的影響效應(yīng)進(jìn)行對(duì)比驗(yàn)證分析。
岸坡巖土體為均質(zhì)直立岸坡,各向滲透系數(shù)相同,其值為0.25 m/d,庫水位下降高度為30 m,下降速度為1 m/d,給水度為0.03,孔隙率為0.25,含水層厚度取升降前后的平均值為40 m,模型的范圍為400 m×60 m×60 m,長(zhǎng)度400 m。
利用程序與簡(jiǎn)化公式計(jì)算對(duì)模型進(jìn)行地下水變動(dòng)預(yù)測(cè),預(yù)測(cè)時(shí)間為30 d,滲透系數(shù)按0.05 m/d,0.1 m/d,0.5 m/d,1.0 m/d,2.0 m/d計(jì)算,水位變化計(jì)算結(jié)果對(duì)比曲線見圖3。
圖3 V-M程序與簡(jiǎn)化公式計(jì)算結(jié)果
經(jīng)對(duì)比分析可知,對(duì)于庫水位等速下降時(shí),滲透系數(shù)越小,二者預(yù)測(cè)浸潤(rùn)線相差越??;相反,滲透系數(shù)越大,此時(shí)在坡體前部相差不大,在坡體后部二者有較明顯差異,簡(jiǎn)化公式計(jì)算結(jié)果偏高,但總體而言二者擬合度較好。
參數(shù)設(shè)置同前,庫水位下降速度按0.25 m/d,0.5 m/d,1.0 m/d,2.0 m/d,5.0 m/d計(jì)算,計(jì)算結(jié)果見圖4。
圖4 V-M程序與簡(jiǎn)化公式計(jì)算結(jié)果
由計(jì)算結(jié)果可知,二者擬合度較高,其偏差只是在坡體后部偏差要大些。
參數(shù)設(shè)置同前,給水度分別按0.005,0.008,0.03,0.08,0.12進(jìn)行計(jì)算,計(jì)算結(jié)果見圖5。
圖5 V-M程序與簡(jiǎn)化公式計(jì)算結(jié)果
由計(jì)算結(jié)果可知,二者計(jì)算結(jié)果基本吻合。隨給水度增大,預(yù)測(cè)結(jié)果與簡(jiǎn)化公式計(jì)算結(jié)果偏差降低。
綜合可知,文中得到的簡(jiǎn)化公式可以滿足計(jì)算精度的要求,且應(yīng)用計(jì)算方便。
由上述分析,令水位的下降高度Δh=v0t,則時(shí)間t=Δh/v0,此時(shí)有
(14)
從圖2可以看出,隨λ的增大,R(λ)值越小,R(λ)為減函數(shù)。也就是說,λ越小,坡體自由水面下降越快;反之,λ越大,坡體自由水面下降越慢。當(dāng)λ=0,R(λ)=1,坡體中的自由水面與庫水位同步升降;當(dāng)λ=∞,R(λ)=0,即坡體中自由水面不受庫水位升降的影響。當(dāng)λ≥2時(shí),R(λ)已非常接近0,此時(shí)庫水位的變化對(duì)地下水浸潤(rùn)線的影響可以忽略。因此可以把該處作為庫水位升降對(duì)地下水浸潤(rùn)線產(chǎn)生影響的水平范圍,即
由前面的計(jì)算分析可知,浸潤(rùn)線的確定與坡體滲透系數(shù)、給水度、含水層平均厚度、水位升降速度及水位升降高度有密切關(guān)系,因此在計(jì)算坡體浸潤(rùn)線時(shí),這些因素的確定就顯得十分重要。
給水度是一個(gè)十分重要常用的水文地質(zhì)參數(shù),它的大小應(yīng)通過實(shí)際測(cè)試加以確定。我們通常把它定義為地下水位下降一個(gè)單位深度,從地下水位延伸到地表面的單位水平面積巖土柱體,在重力作用下釋出的水的體積。在當(dāng)實(shí)際工作中當(dāng)無實(shí)驗(yàn)資料時(shí),往往借助經(jīng)驗(yàn)公式。
別申斯基(1960)根據(jù)砂礫石土樣(最小粒經(jīng)d=0.25 cm,最小滲透系數(shù)k=1.87×10-2cm/s)的試驗(yàn)結(jié)果分析得出給水度μ的簡(jiǎn)單經(jīng)驗(yàn)公式為
式中:k為滲透系數(shù)/(m/d),此式只適合粗顆粒土,不適用于顆粒小于細(xì)砂的土質(zhì)。
給水度代表巖土體的排泄水量的能力,它的大小除了跟滲透系數(shù)密切相關(guān)外,還與土的密實(shí)度有關(guān),即與孔隙率密切相關(guān)。文獻(xiàn)[2]根據(jù)國(guó)內(nèi)外砂礫石和黏性土的試驗(yàn)資料,分析求得給水度的經(jīng)驗(yàn)公式為
μ=1.137n(0.000 117 5)0.607(6+lgk)
式中n——孔隙率;
k——滲透系數(shù)/(cm/s)。
在無試驗(yàn)資料時(shí),可用此式來確定給水度的大小。
由前面分析可知,為了求解非穩(wěn)定滲流的基本微分方程,通常采用始、末段潛水流厚度的平均值hm代替式中H,這樣的假定對(duì)于升降高度比較小的情況適用,但對(duì)于升降高度比較大時(shí)誤差就比較大,為此可采用下面的方法來確定含水層平均厚度。
(1)不透水層為水平面
根據(jù)滲流分析研究成果,坡體中的浸潤(rùn)線可簡(jiǎn)化為一條拋物線(拋物線方程y2=ax+b),這個(gè)假定在堤壩的分析研究中有應(yīng)用。為了得到浸潤(rùn)線的平均厚度,可用拋物線模型來確定潛水流的含水層平均厚度,如圖6為庫水位下降和上升兩種情況下含水層厚度計(jì)算簡(jiǎn)圖。以水位下降為例,假設(shè)初始含水層厚h0,水位下降高度為h,影響范圍為s,把點(diǎn)0(0,0),P(s,h)代入拋物線方程,則得到下降后的浸潤(rùn)線方程為
(15)
圖中浸潤(rùn)線osp的面積為
將含水層轉(zhuǎn)化為寬度為s的矩形,可得到庫水位下降時(shí)含水層平均厚度hm,即
(16)
圖6 含水層厚度計(jì)算
同理可以推導(dǎo)出庫水位上升時(shí)含水層的平均厚度為
(17)
式中:h為水位升降高度。
從得到的公式可以看出,當(dāng)水位升降時(shí)含水層的平均厚度僅與水位升降高度和初始含水層厚度有關(guān)。在這種情況下很容易確定含水層的厚度,在實(shí)際工程應(yīng)用中十分方便。
(2)不透水層為不規(guī)則面的情況
在實(shí)際工程中經(jīng)常遇到不透水層形狀不規(guī)則(圖7),此時(shí)含水層平均厚度一般可采用如下公式來進(jìn)行近似計(jì)算
(18)
由于一般地質(zhì)剖面都是用AutoCAD繪制,在AutoCAD里確定OABC區(qū)的面積是非常容易的,因此計(jì)算平均含水層厚度也就變的十分簡(jiǎn)單。
圖7 含水層厚度計(jì)算
(1)根據(jù)非穩(wěn)定滲流運(yùn)動(dòng)基本微分方程和邊界條件,推導(dǎo)庫水位等速下降條件下坡體內(nèi)浸潤(rùn)線的計(jì)算公式,并對(duì)其進(jìn)行多項(xiàng)式擬合得到簡(jiǎn)化公式,便于實(shí)際工程應(yīng)用。
(2)通過Visual-Modflow預(yù)測(cè)結(jié)果與簡(jiǎn)化公式的預(yù)測(cè)結(jié)果對(duì)比分析可知,雖然簡(jiǎn)化公式在考慮諸多因素情況下,與有限元模擬結(jié)果有一定偏差,但多數(shù)情況下偏差很小,說明該簡(jiǎn)化公式的應(yīng)用是可行的。但考慮到各因素的綜合影響效應(yīng),在應(yīng)用中應(yīng)需要對(duì)其進(jìn)行適當(dāng)?shù)男拚?/p>
[1]顧慰慈.滲流計(jì)算原理及其應(yīng)用[M].北京:中國(guó)建材工業(yè)出版社,2000:12-22
[2]毛昶熙.滲流計(jì)算分析與控制(第二版)[M].北京:中國(guó)水利水電出版社,2003:352-354
[3]魏玉虎,許 模,等.長(zhǎng)江三峽工程奉節(jié)庫岸斜坡地下水滲流場(chǎng)模擬分析[J].四川地質(zhì)學(xué)報(bào),2003,23(1):18-22
[4]李 平,盧文喜,等.Visual MODFLOW在地下水?dāng)?shù)值模擬中的應(yīng)用[J].工程勘察,2006(3):24-27
[5]周念清,等.MODFLOW在宿遷市地下水資源評(píng)價(jià)中的應(yīng)用[J].水文地質(zhì)與工程地質(zhì),2000(6):9-13
[6]王立久.尾礦壩滲流浸潤(rùn)線計(jì)算的實(shí)用方法[J].金屬礦山,1987(5)
[7]吳漢杰,林建南,侯曙光.降雨入滲條件下公路邊坡滲流數(shù)值分析[J].勘察科學(xué)技術(shù),2005(3)
[8]鄭穎人,時(shí)衛(wèi)民,孔位學(xué).庫水位下降時(shí)滲透力及地下水浸潤(rùn)線的計(jì)算[J].巖石力學(xué)與工程學(xué)報(bào),2004(18)
[9]張友誼,胡卸文.庫水位等速上升作用下岸坡地下水浸潤(rùn)線的計(jì)算[J].水文地質(zhì)工程地質(zhì),2007(5)