考慮滲透系數(shù)變異性的堤防滲透破壞概率分析方法

柯浩進(jìn)， 王 媛， 馮 迪

(河海大學(xué)巖土力學(xué)與堤壩工程教育部國家重點(diǎn)實(shí)驗(yàn)室，土木與交通學(xué)院隧道與城市軌道工程研究所, 南京 210098)

洪水是威脅中國人民生命財(cái)產(chǎn)安全的重要自然災(zāi)害之一，堤防工程作為抵御洪水的一道安全屏障，其滲流穩(wěn)定性尤為重要。然而中國的堤防流域廣、年代久，堤防本身的地質(zhì)條件存在較大的差異，且經(jīng)過歷史變遷，同一河段堤身填筑材料也不一致，以上諸多隨機(jī)因素造成了堤防土體滲透系數(shù)具有一定的變異性，對堤防滲流場有著重要影響[1-4]。由于以往確定性方法無法表征隨機(jī)因素對堤防滲流場的影響，以可靠度理論為基礎(chǔ)的概率分析方法被廣泛應(yīng)用到堤防工程中來。李錦輝等[5]基于攝動(dòng)隨機(jī)有限元法研究了水位變化、堤頂寬度等因素與堤防滲透失穩(wěn)概率的關(guān)系。王靖文等[6]將故障樹和蒙特卡洛法結(jié)合起來建立了求解堤防滲透失穩(wěn)概率的方法。王亞軍等[7]基于蒙特卡洛法和first order reliability method理論建立了堤防系統(tǒng)廣義模糊可靠度算法，分析了堤防不同失效模式下的失效概率。焦小超[8]以水頭和土體的物理力學(xué)參數(shù)作為隨機(jī)變量，計(jì)算了堤防管涌破壞極限狀態(tài)方程。張秀勇等[9]結(jié)合有限元方法與蒙特卡洛法，采用堤防單元土體有效土重力與土體底部的滲透力之差作為極限狀態(tài)方程，方程中忽略了土體黏聚力及摩擦力的影響，使得結(jié)果偏保守。邢萬波等[10]針對流土、管涌建立了不同失效模式下的極限狀態(tài)方程，但卻僅限于堤基部分，不能用于堤身的滲透破壞計(jì)算。高昂等[11]、雷鵬等[12]在有限元方法的基礎(chǔ)上，將神經(jīng)網(wǎng)絡(luò)、拉丁超立方抽樣等方法與有限元方法相結(jié)合，但卻未能將土體參數(shù)的變異性考慮到有限元模型中去，無法考慮水力梯度的不確定性。

基于以上背景，將滲透系數(shù)視為隨機(jī)變量，與有限元軟件Abaqus相結(jié)合，結(jié)合蒙特卡洛法分析了滲透系數(shù)變異性對堤防滲流場的影響?；诓煌茐膮^(qū)域，將水力梯度、土體黏聚力、摩擦角等因素作為隨機(jī)變量，建立了堤防滲流破壞極限狀態(tài)方程，并結(jié)合改進(jìn)的驗(yàn)算點(diǎn)法(JC法)實(shí)現(xiàn)堤防工程滲流破壞概率分析。

1 堤防隨機(jī)滲流場分析

1.1 堤防隨機(jī)滲流分析方法

滲透系數(shù)的不確定性對堤防滲流場內(nèi)各區(qū)域水力梯度有一定的影響，為此采取將滲透系數(shù)視為隨機(jī)變量的方法，與有限元方法相結(jié)合，分析滲透系數(shù)不同變異性下的堤防滲流場情況。利用有限元軟件Abaqus進(jìn)行建模并劃分網(wǎng)格，針對劃分好的單元網(wǎng)格，基于Python語言編寫滲透系數(shù)隨機(jī)賦值的腳本，通過腳本讀取模型的每一個(gè)網(wǎng)格單元，對應(yīng)每一個(gè)單元設(shè)置一個(gè)截面，同時(shí)隨機(jī)生成一個(gè)滿足要求的材料，并將該材料賦給截面，通過遍歷循環(huán)從而建立一個(gè)堤防隨機(jī)滲流模型。具體如圖1所示。

圖1 隨機(jī)有限元方法流程圖Fig.1 Stochastic finite element method flow chart

1.2 有限元模型建立

選取某一典型堤防斷面，模型參考Lam壩體有限元滲流分析[13]，模型如圖2所示，以堤基為基準(zhǔn)面，堤身高度為12 m，堤頂長度為4 m，兩側(cè)坡比均為1∶2，上游水位設(shè)置為10 m，下游水位為2 m，滲透系數(shù)均值設(shè)置為1×10-7cm/s，標(biāo)準(zhǔn)差為2×10-8，據(jù)上述內(nèi)容建立有限元模型，單元?jiǎng)澐纸Y(jié)果如圖2所示，共1 052個(gè)節(jié)點(diǎn)，971個(gè)單元。

圖2 堤防模型網(wǎng)格單元Fig.2 Levee model mesh units

1.3 滲透系數(shù)變異性對堤防滲流場的影響分析

將上述有限元模型滲透系數(shù)均值設(shè)定為1×10-7cm/s,變異系數(shù)分別設(shè)置為0.2、0.4、0.6來進(jìn)行計(jì)算。取1 000次計(jì)算進(jìn)行結(jié)果統(tǒng)計(jì)并使用Tecplot繪制云圖。繪制了不同變異系數(shù)條件下的水頭均值圖[圖3(a)]和標(biāo)準(zhǔn)差圖(圖4)以及水力梯度均值圖[圖3(b)]和標(biāo)準(zhǔn)差圖(圖5)，將各云圖進(jìn)行比較分析發(fā)現(xiàn)，不同變異性條件下，滲流場的水頭和水力梯度均值無明顯波動(dòng)。

圖3 均值等值線圖Fig.3 Mean isogram map

變異系數(shù)的改變對滲流場的均值無明顯影響，在變異性較小的情況下，非均質(zhì)土體滲流場與均質(zhì)土體滲流場較為接近。

圖4 不同變異系數(shù)的水頭標(biāo)準(zhǔn)差圖Fig.4 Standard deviation map of water head with different coefficient of variation

圖5 不同變異系數(shù)條件下的水力梯度標(biāo)準(zhǔn)差云圖Fig.5 Standard deviation map of hydraulic gradient under different coefficient of variation

從圖4對比可以看出，在穩(wěn)定滲流條件下，堤防水位以下的邊界位置標(biāo)準(zhǔn)差較小，這是因?yàn)榉€(wěn)定滲流時(shí)，邊界水頭是確定的，所以標(biāo)準(zhǔn)差較小。標(biāo)準(zhǔn)差較大的位置出現(xiàn)在堤防背水坡逸出點(diǎn)附近。變異系數(shù)的變化對標(biāo)準(zhǔn)差的分布形式無明顯影響，但變異系數(shù)的增大使堤防各區(qū)域水頭的標(biāo)準(zhǔn)差變大。

從圖5對比可以看出，變異系數(shù)的增大使得各區(qū)域水力梯度標(biāo)準(zhǔn)差變大，變異系數(shù)的變化對水力梯度標(biāo)準(zhǔn)差的分布形式無明顯影響。水力梯度標(biāo)準(zhǔn)差在堤防邊界處較大，堤身內(nèi)部較小，這就意味著在堤防邊界處水力梯度波動(dòng)幅度較大，而邊界處正是堤防發(fā)生破壞開始的位置，所以研究此處堤防破壞的概率是有必要的。

2 基于JC法的堤防滲流失穩(wěn)可靠性分析

2.1 JC法基本原理

工程中在進(jìn)行結(jié)構(gòu)可靠度分析時(shí)，常常會(huì)有多個(gè)因素同時(shí)影響著結(jié)構(gòu)整體的穩(wěn)定性，所以構(gòu)造出的結(jié)構(gòu)極限狀態(tài)方程往往也包含多個(gè)變量：

Z=g(X)=0

(1)

JC法實(shí)質(zhì)是改進(jìn)的驗(yàn)算點(diǎn)法，取極限狀態(tài)面上的一點(diǎn)x*，并在該點(diǎn)按泰勒級數(shù)展開保留至一次項(xiàng)：

(2)

在隨機(jī)變量空間X里，方程ZL=0 為過該點(diǎn)的極限狀態(tài)面的切平面。假設(shè)X服從正態(tài)分布，利用相互獨(dú)立正態(tài)分布的隨機(jī)變量線性組合的性質(zhì)，ZL的均值和標(biāo)準(zhǔn)差為

(3)

(4)

則有：

(5)

定義一個(gè)靈敏系數(shù)axi，令：

(6)

經(jīng)過標(biāo)準(zhǔn)化處理可得到：

(7)

(8)

(9)

JC法計(jì)算可靠度的步驟是先假定一個(gè)初始的驗(yàn)算點(diǎn)，然后檢驗(yàn)隨機(jī)變量是否服從正態(tài)分布，如果不服從，則進(jìn)行當(dāng)量正態(tài)化計(jì)算，將經(jīng)過當(dāng)量正態(tài)化轉(zhuǎn)換后的均值和標(biāo)準(zhǔn)差代入計(jì)算，得到新的驗(yàn)算點(diǎn)坐標(biāo)，并重復(fù)上述步驟進(jìn)行迭代計(jì)算，通過比較前后兩驗(yàn)算點(diǎn)坐標(biāo)之差是否小于允許值來終止計(jì)算，計(jì)算流程圖如圖6所示。

圖6 JC法計(jì)算流程圖Fig.6 JC method flow chart

2.2 堤防滲流破壞極限狀態(tài)方程

堤防的滲流破壞極限狀態(tài)方程是用于判斷堤防穩(wěn)定與否的關(guān)鍵依據(jù)。工程中對滲透破壞穩(wěn)定性的判別準(zhǔn)則是將滲流場中的最大水力梯度值Jmax和土體的臨界水力梯度Jc相比較，如果Jmax

圖7 堤防典型斷面橫剖圖Fig.7 Section plane of the levee

由于迎水坡水位較高，土體是受到水的壓力作用的，不易發(fā)生破壞，所以滲流危險(xiǎn)區(qū)域取堤基滲流出口區(qū)域①和背水坡堤身區(qū)域②。

2.2.1 堤基滲流出口區(qū)域

考慮此處土體在自下而上的滲流條件下發(fā)生頂沖破壞，單元土體的受力分析如圖8所示，得到此處的極限狀態(tài)方程為

g(x)=γ′+fφ+C-γωJY

(10)

式(10)中：γ′ 為土體的浮容重；fφ為單元土體破壞面上所產(chǎn)生的摩擦力；C為單元土體破壞面上所產(chǎn)生的黏聚力；γω為水的容重；JY為豎直方向的水利梯度?？紤]到土體在滲流作用下發(fā)生松動(dòng)，土體的摩擦阻力就不存在了，所以安全起見可以忽略摩擦力的影響，則極限狀態(tài)方程為

g(x)=γ′+C-γωJY

(11)

圖8 堤基單元土體受力分析圖Fig.8 Force analysis diagram on dyke foundation unit

2.2.2 背水坡堤身區(qū)域

背水坡土體在滲流作用下需考慮兩種破壞形式，一種是土體在滲透力作用下垂直于坡面的方向發(fā)生頂沖破壞，另一種是沿坡面方向發(fā)生下滑破壞。對于單元土體受滲流作用發(fā)生垂直于坡面方向的頂沖破壞時(shí)，其受力分析如圖9所示，此時(shí)極限狀態(tài)方程為

g(x)=γ′cosa+fφ+C-γωJYcosa-γωJxsina

(12)

式(12)中：a為坡面與水平方向的夾角；Jx為水平方向的水力梯度。

fφ=(γ′sina+γωJxcosa-γωJYsina)tanφ

(13)

式(13)中：φ為土體的內(nèi)摩擦角。

圖9 垂直坡面發(fā)生頂沖破壞受力分析圖Fig.9 Force analysis diagram of failure vertical to slope surface

對于單元土體受滲流作用沿坡面發(fā)生滑動(dòng)時(shí)，其受力分析如圖10所示，此時(shí)極限狀態(tài)方程為

g(x)=γωJYsina+fφ+C-γωJxcosa-γ′sina

(14)

fφ=(γ′cosa-γωJYcosa-γωJxsina)tanφ

(15)

圖10 沿坡面發(fā)生滑動(dòng)破壞受力分析圖Fig.10 Force analysis diagram of failure along the slope

以上為堤防在滲流作用下不同破壞形式的極限狀態(tài)方程，極限狀態(tài)方程中Jx、JY、C、φ為隨機(jī)變量，因此需要知道各隨機(jī)變量均值及標(biāo)準(zhǔn)差，結(jié)合JC法便可以求得堤防滲流破壞的可靠度。

2.3 堤防滲流穩(wěn)定可靠性分析

結(jié)合上述堤防模型及極限狀態(tài)方程，分析不同滲透系數(shù)變異性下的堤防滲流破壞概率。由于土體滲透系數(shù)的不確定性，水力梯度較大值不會(huì)只出現(xiàn)在一個(gè)位置，所以應(yīng)對危險(xiǎn)區(qū)域內(nèi)的不同位置進(jìn)行水力梯度統(tǒng)計(jì)。通過以上水力梯度云圖可以發(fā)現(xiàn)堤腳附近水力梯度較大，將該區(qū)域作為堤防滲流破壞的危險(xiǎn)區(qū)域，具體位置如圖11虛線處所示。堤防土體參數(shù)設(shè)置服從正態(tài)分布，具體如表1所示。

表1 土體材料參數(shù)Table 1 Soil materials parameters

圖11 水力梯度統(tǒng)計(jì)區(qū)域示意圖Fig.11 Statistical region sketch of hydraulic gradient

將1 000次隨機(jī)滲流分析的結(jié)果進(jìn)行統(tǒng)計(jì)分析，將計(jì)算得到的數(shù)據(jù)采用SPSS軟件作統(tǒng)計(jì)分析，經(jīng)擬合、K-S檢驗(yàn)后得到水力梯度服從正態(tài)分布。對區(qū)域處所有節(jié)點(diǎn)進(jìn)行遍歷計(jì)算，求得邊界各節(jié)點(diǎn)處的破壞概率，在堤基及堤身不同的破壞形式分別各取破壞概率最大的位置點(diǎn)進(jìn)行分析，結(jié)果如表2所示。

表2 不同破壞形式下破壞概率最大區(qū)域水力梯度統(tǒng)計(jì)結(jié)果Table 2 Statistical results of hydraulic gradient in regions with highest failure probability under different failure modes

堤基處破壞形式往往是垂直于堤基的頂沖破壞，從表2數(shù)據(jù)可看出，堤基處x方向的水力梯度均值及標(biāo)準(zhǔn)差相較y方向小很多，堤身處x向和y向的水力梯度差異則相對較小。變異系數(shù)較小時(shí)對堤防各方向水力梯度均值的影響較小，對標(biāo)準(zhǔn)差存在較大影響。不同變異條件下，堤身發(fā)生滑動(dòng)破壞與發(fā)生頂沖破壞概率最大的位置為同一點(diǎn)。

表3 堤防滲流破壞計(jì)算結(jié)果Table 3 Calculated results of seepage failure of levee

表3為堤防滲流破壞計(jì)算結(jié)果，可以看出，在以上變異性條件下，堤防的堤基處發(fā)生滲流破壞的概率較小，幾乎不會(huì)發(fā)生破壞。堤身處發(fā)生頂沖破壞的概率同樣較小，相比之下在滲流作用下，滲透系數(shù)變異性增大，堤身發(fā)生滑動(dòng)破壞的概率也會(huì)隨之增大，變異系數(shù)0.6時(shí)，堤身發(fā)生局部滑動(dòng)破壞的概率接近14%，該概率下的堤防存在發(fā)生破壞的可能。

3 結(jié)論

以一典型堤防斷面為例，構(gòu)建了堤防滲流隨機(jī)有限元計(jì)算模型，分析比較了滲透系數(shù)隨機(jī)條件下水頭和水力梯度的均值和標(biāo)準(zhǔn)差情況，根據(jù)分析結(jié)果選取堤防滲流不同破壞區(qū)域單元土體進(jìn)行受力分析并構(gòu)造了堤防滲流破壞的極限狀態(tài)方程，最后結(jié)合JC法對堤防的滲流破壞進(jìn)行可靠度分析和破壞概率計(jì)算，得到了以下結(jié)論。

(1)穩(wěn)定隨機(jī)滲流條件下，滲透系數(shù)變異性的變化對滲流場水頭和水力梯度的標(biāo)準(zhǔn)差存在一定影響，變異系數(shù)變大標(biāo)準(zhǔn)差隨之增大，變異系數(shù)增大，堤防發(fā)生局部滲流破壞的概率增大。

(2)通過不同區(qū)域的受力分析可構(gòu)建相應(yīng)的堤防滲流極限狀態(tài)方程，從而可以綜合考慮堤身與堤基區(qū)域可能發(fā)生的破壞。

(3)通過JC法理論和滲流隨機(jī)有限元方法進(jìn)行堤防滲流的可靠度分析，可得到堤防不同位置發(fā)生破壞的概率，為實(shí)際工程的穩(wěn)定分析提供可靠的依據(jù)，具有一定的應(yīng)用價(jià)值。