隨機(jī)庫(kù)水位作用下病險(xiǎn)土石壩的失效概率研究

李 澤，彭 普，張小艷，譚志華

（1.昆明理工大學(xué)建筑工程學(xué)院，昆明 650500；2.昆明理工大學(xué)電力工程學(xué)院，昆明 650500；3.云南省水利水電勘測(cè)設(shè)計(jì)研究院，昆明 650021）

0 引 言

土石壩是指由當(dāng)?shù)赝潦辖?jīng)過(guò)拋填、輾壓等方法堆筑而成的擋水建筑物。目前我國(guó)在運(yùn)行的8萬(wàn)多座土石壩水庫(kù)中有接近3 萬(wàn)座為病險(xiǎn)水庫(kù)，這些病險(xiǎn)土石壩大部分都是在20 世紀(jì)50-70年代修建的中、小型土石壩。由于受當(dāng)時(shí)筑壩技術(shù)和資金投入等方面的限制，很多病險(xiǎn)土石壩工程的建設(shè)標(biāo)準(zhǔn)偏低，在長(zhǎng)期服役以后潰壩災(zāi)害時(shí)有發(fā)生，給人民財(cái)產(chǎn)安全造成重大損失［1］。病險(xiǎn)土石壩的穩(wěn)定安全問(wèn)題是關(guān)乎國(guó)計(jì)民生的重大問(wèn)題。國(guó)家相關(guān)部門已經(jīng)開(kāi)展了病險(xiǎn)土石壩的全面調(diào)查和初步治理。

中小型病險(xiǎn)土石壩具有以下特點(diǎn)：①由于特定的歷史條件，壩體的填筑材料不良、抗?jié)B性能不滿足要求，長(zhǎng)期服役以后壩體材料的物理力學(xué)性能嚴(yán)重衰減，與初始設(shè)計(jì)值相差較大；②長(zhǎng)期服役后壩體材料的物理力學(xué)性能在時(shí)間、空間上都具有不確定性和隨機(jī)性；③中小型土石壩水庫(kù)的庫(kù)容較小，一般為季調(diào)節(jié)水庫(kù)，庫(kù)水位受降雨和用水需求的影響在雨季和旱季呈現(xiàn)隨機(jī)變化；④絕大部分中小型土石壩為均質(zhì)壩體、無(wú)防滲墻，壩體的穩(wěn)定性受庫(kù)水位的影響較大。在當(dāng)前關(guān)于土石壩穩(wěn)定性的研究中還沒(méi)有能夠全面反映這些因素的研究成果［2］。

目前，在中國(guó)土石壩傳統(tǒng)工程實(shí)踐中普遍采用確定性分析方法進(jìn)行土石壩的穩(wěn)定性計(jì)算。我國(guó)土石壩規(guī)范推薦使用極限平衡分析法（LEM）對(duì)土石壩進(jìn)行穩(wěn)定性分析［3］，其通過(guò)事先假定滑裂面并對(duì)滑裂面上的土體進(jìn)行條分，并對(duì)條間力的方向進(jìn)行假定將超靜定問(wèn)題轉(zhuǎn)化為靜定問(wèn)題，LEM 法概念清晰、計(jì)算簡(jiǎn)便且工程應(yīng)用非常廣泛；但該方法未能考慮土石壩的應(yīng)力-應(yīng)變關(guān)系，因此不能給出壩體的應(yīng)力場(chǎng)、位移場(chǎng)［4，5］。同時(shí)，以有限單元法（FEM）和有限差分法（FDM）為代表的數(shù)值分析方法在土石壩穩(wěn)定性分析方面得到廣泛應(yīng)用，其能使土體同時(shí)滿足靜力平衡條件、應(yīng)變相容條件和本構(gòu)關(guān)系，可以處理任意的土層、幾何形狀、荷載和邊界，不僅可獲得壩體詳細(xì)的應(yīng)力-應(yīng)變分布而且可以獲得比LEM 方法更合理的穩(wěn)定性分析結(jié)果［6～9］；但其存在初始應(yīng)力難以確定、巖土本構(gòu)關(guān)系復(fù)雜等困難。

病險(xiǎn)土石壩的穩(wěn)定性與壩體材料的物理力學(xué)參數(shù)、庫(kù)水位隨機(jī)變化等因素密切相關(guān)。因此，土石壩的穩(wěn)定性問(wèn)題是一個(gè)與材料抗剪強(qiáng)度參數(shù)和庫(kù)水位的隨機(jī)性相關(guān)的不確定性問(wèn)題，對(duì)于不確定性問(wèn)題需采用可靠度分析方法進(jìn)行研究。我國(guó)水利水電行業(yè)規(guī)范規(guī)定土石壩的設(shè)計(jì)在單一安全系數(shù)法的基礎(chǔ)上需采用以概率理論為基礎(chǔ)、以分項(xiàng)系數(shù)表達(dá)的極限狀態(tài)設(shè)計(jì)方法進(jìn)行穩(wěn)定性驗(yàn)算［10］。同時(shí)，眾多學(xué)者基于概率理論開(kāi)展了土石壩的可靠度研究，可靠度方法在土石壩工程中得到充分重視和應(yīng)用［11-17］。

在現(xiàn)有的土石壩可靠度分析中使用的穩(wěn)定性計(jì)算方法主要以剛體極限平衡法（LEM）、有限單元法（FEM）為主［18-21］，同時(shí)塑性極限分析方法（UBM）在土石壩的可靠度分析方面得到深入研究［22-24］。相對(duì)于LEM 和FEM，塑性極限分析方法在尋求土石壩的極限狀態(tài)方面存在理論上的優(yōu)勢(shì)，根據(jù)極值定理其可以高效、準(zhǔn)確地獲取土石壩破壞時(shí)的極限狀態(tài)（應(yīng)力場(chǎng)和速度場(chǎng)）。鑒于此，本文同時(shí)考慮土石壩土體的抗剪強(qiáng)度參數(shù)的變異性以及庫(kù)水位的隨機(jī)性，將塑性極限分析上限法理論、有限元離散思想、數(shù)學(xué)規(guī)劃理論、蒙特卡洛方法結(jié)合起來(lái)對(duì)病險(xiǎn)土石壩的可靠度進(jìn)行研究。

1 土石壩的極限狀態(tài)功能函數(shù)

本文同時(shí)考慮土體抗剪強(qiáng)度參數(shù)和上游庫(kù)水位的隨機(jī)性，則土石壩的穩(wěn)定性與凝聚力、摩擦角以及上游庫(kù)水位相關(guān)；使用文獻(xiàn)［25］的方法，本文采用容重超載法使土石壩達(dá)到破壞的極限狀態(tài)。在進(jìn)行可靠度分析時(shí)，病險(xiǎn)土石壩的極限狀態(tài)功能函數(shù)可表示為：

式中：cr和φr分布是土石壩土體的凝聚力隨機(jī)量和摩擦角隨機(jī)量；Hrw是上游庫(kù)水位隨機(jī)量；km(cr，φr，Hrw)是與cr、φr和Hrw相關(guān)的安全系數(shù)隨機(jī)量；當(dāng)Z＞0時(shí)表示土石壩整體穩(wěn)定，當(dāng)Z=0時(shí)表示土石壩處于臨界狀態(tài)，當(dāng)Z＜0 時(shí)表示土石壩發(fā)生失穩(wěn)（失效）。

土石壩的上游庫(kù)水位在最低水位與最高水位之間隨機(jī)變化，庫(kù)水位一般在死水位至校核洪水位范圍內(nèi)隨機(jī)變化，并近似服從正態(tài)分布，本文采用截尾正態(tài)分布模擬庫(kù)水位的變化。庫(kù)水位隨機(jī)數(shù)按下式進(jìn)行生成：

式中：tw=(1，…，nw)，nw是土石壩上游庫(kù)水位蒙特卡洛隨機(jī)數(shù)的數(shù)量；Hrw(tw)是庫(kù)水位的第tw個(gè)隨機(jī)數(shù)；μw是庫(kù)水位的均值；σw是庫(kù)水位的標(biāo)準(zhǔn)差；Random是正態(tài)分布隨機(jī)數(shù)生成函數(shù)；Normal表示庫(kù)水位隨機(jī)數(shù)符合正態(tài)分布；Hlb是庫(kù)水位的下界，取最低水位；Hub是庫(kù)水位的上界，取最高水位。

假設(shè)土石壩土體材料的凝聚力和摩擦角均符合對(duì)數(shù)正態(tài)分布，其隨機(jī)數(shù)按下式生成：

式中：tm=(1，…，nm)，nm是土體材料凝聚力和摩擦角蒙特卡洛隨機(jī)數(shù)的數(shù)量；cr(tm)是土體材料凝聚力的第tm個(gè)隨機(jī)數(shù)；φr(tm)是土體材料摩擦角的第tm個(gè)隨機(jī)數(shù)；μc是土體材料凝聚力的均值；μφ是土體材料摩擦角的均值；σc是土體材料凝聚力的標(biāo)準(zhǔn)差；σφ是土體材料摩擦角的標(biāo)準(zhǔn)差；Random是正態(tài)分布隨機(jī)數(shù)生成函數(shù)；lognormal表示隨機(jī)數(shù)符合對(duì)數(shù)正態(tài)分布。

在同時(shí)考慮庫(kù)水位、抗剪強(qiáng)度參數(shù)的隨機(jī)性時(shí)，定義土石壩的安全系數(shù)如下：

式中：tw=(1，…，nw)，tm=(1，…，nm)，c′r(tm)是強(qiáng)度折減以后土體材料凝聚力的第tm個(gè)隨機(jī)數(shù)；φ′r(tm)是強(qiáng)度折減以后土體材料摩擦角的第tm個(gè)隨機(jī)數(shù)。

2 土石壩上限法數(shù)學(xué)規(guī)劃模型

塑性極限分析的上限定理是求解巖土體承載力的高效工具，上限法的原理是：通過(guò)構(gòu)建土石壩的機(jī)動(dòng)許可速度場(chǎng)，機(jī)動(dòng)許可速度場(chǎng)需同時(shí)滿足：三角形單元塑性流動(dòng)約束條件、公共邊速度不連續(xù)約束條件、三角形單元速度邊界約束條件；每個(gè)機(jī)動(dòng)許可速度場(chǎng)對(duì)應(yīng)有唯一的外荷載，根據(jù)上限定理可知最小外荷載必然最接近其真實(shí)的極限荷載。所以，上限法本質(zhì)上可歸結(jié)為一個(gè)求極小化的數(shù)學(xué)規(guī)劃問(wèn)題。本文在前人的研究基礎(chǔ)上［26，27］，基于上限定理建立土石壩可靠度分析的上限法數(shù)學(xué)模型。

2.1 數(shù)值離散

本文使用非共結(jié)點(diǎn)的三角形單元離散土石壩［如圖1（a）所示］［23，26］，三角形單元的每個(gè)結(jié)點(diǎn)具有沿x方向的速度uexi、y方向的速度ueyi，其中i=(1，2，3)；同時(shí)，受隨機(jī)庫(kù)水位作用單元的每個(gè)結(jié)點(diǎn)具有孔隙水壓力隨機(jī)量prei，其中i=(1，2，3)。為了構(gòu)建土石壩的機(jī)動(dòng)許可速度場(chǎng)，相鄰有限單元之間公共邊需存在速度間斷，相鄰有限單元之間公共邊的速度間斷如圖1（b）所示，一條公共邊共有標(biāo)號(hào)①②③④的4 個(gè)結(jié)點(diǎn)，其中結(jié)點(diǎn)①③屬于單元a、②④屬于單元b，公共邊的每個(gè)結(jié)點(diǎn)具有沿x方向的速度udxi和y方向的速度udyi，i=(1，2，3，4)。

2.2 土石壩的穩(wěn)定滲流計(jì)算

為了獲得土石壩中每一點(diǎn)的孔隙水壓力，本文采用穩(wěn)定滲流理論對(duì)隨機(jī)庫(kù)水位作用下的土石壩進(jìn)行滲流場(chǎng)分析［28］。二維穩(wěn)定滲流公式如下：

式中：kx是土體材料x(chóng)方向的滲透系數(shù)，ky是土體材料y方向的滲透系數(shù)，Hr是土石壩內(nèi)各點(diǎn)的隨機(jī)水頭函數(shù)。

2.3 上限法數(shù)學(xué)模型的約束條件

（1）三角形單元塑性流動(dòng)約束條件。假設(shè)土石壩土體為理想剛塑性體且符合相關(guān)聯(lián)流動(dòng)性法則，則土體由變形協(xié)調(diào)條件得到的塑性應(yīng)變分量等于由關(guān)聯(lián)流動(dòng)法則和屈服條件所得塑性應(yīng)變率分量。土石壩三角形單元塑性流動(dòng)約束條件可表示為：

式中：e=(1，…，Ne)，Ne是土石壩中所有有限單元的數(shù)量；Ae是有限單元e的面積。

由于篇幅有限，系數(shù)be1，be2，be3，ce1，ce2，ce3和Ak，Bk，Ck，k=(1，…，m)等的涵義詳見(jiàn)文獻(xiàn)［26］。

（2）公共邊速度不連續(xù)約束條件。塑性極限分析上限法規(guī)定在土體相連單元公共邊上速度可間斷但必須符合相關(guān)聯(lián)流動(dòng)法則，則相連有限單元公共邊速度不連續(xù)約束條件：

式中：d=(1，…，Nd)，Nd是土石壩中所有有限單元公共邊的數(shù)量；分別是單元公共邊d的4 個(gè)非負(fù)塑性乘子；θd是公共邊的傾角，逆時(shí)針為正；φr(tm)是土石壩土體材料摩擦角的第tm個(gè)隨機(jī)數(shù)。

（3）三角形單元速度邊界約束條件。土石壩邊界上的三角形單元速度場(chǎng)必須滿足已知的邊界條件。在速度為0的邊界上土石壩壩體三角形單元速度邊界條件為：

式中：Ab是邊界上的有限單元b的坐標(biāo)轉(zhuǎn)換矩陣；ub是邊界b的速度向量；b=(1，…，Nb)，Nb是土石壩中邊界上速度等于0 的有限單元的數(shù)量。

2.4 上限法數(shù)學(xué)模型的目標(biāo)函數(shù)

根據(jù)功能原理，在整個(gè)土石壩中內(nèi)力做的內(nèi)功功率應(yīng)該等于外力做的外功功率。當(dāng)土石壩內(nèi)有孔隙水壓力作用，并考慮容重超載時(shí)，土石壩的內(nèi)功功率與外功功率相等條件為：

式中：We是有限單元的內(nèi)功功率；Wd是有限單元公共邊的內(nèi)功功率；WG是自重在有限單元結(jié)點(diǎn)速度上所做的外功功率；W pe是孔隙水壓力在有限單元連續(xù)體內(nèi)所作的外功功率；W pd是孔隙水壓力在有限單元公共邊上所作的外功功率；kγ(tw，tm)是與第tw個(gè)隨機(jī)庫(kù)水位和第tm個(gè)抗剪強(qiáng)度參數(shù)對(duì)應(yīng)的容重超載系數(shù)。

土石壩內(nèi)有限單元的內(nèi)功功率為：

土石壩內(nèi)有限單元公共邊的內(nèi)功功率為：

式中：ld是公共邊的長(zhǎng)度。

土石壩中土體的自重在有限單元結(jié)點(diǎn)速度上所做的外功功率為：

式中：γe是單元的土體容重。

土石壩內(nèi)部的孔隙水壓力在有限單元連續(xù)體內(nèi)所作的外功功率：

土石壩內(nèi)部的孔隙水壓力在有限單元公共邊上所作的外功功率為：

式（9）是一個(gè)非線性表達(dá)式，為了避免求解非線性規(guī)劃問(wèn)題，本文將其進(jìn)行線性化處理。根據(jù)上限定理可知，機(jī)動(dòng)許可速度場(chǎng)中結(jié)點(diǎn)的速度表示單元相對(duì)于固定邊界發(fā)生塑性流動(dòng)的相對(duì)大小，其與物理學(xué)上的速度的概念不同，其僅僅表示單元發(fā)生塑性流動(dòng)的趨勢(shì)；因此，本文增加一個(gè)附加約束，設(shè)WG=1.0，將式（9）等效變換為兩個(gè)數(shù)學(xué)線性式：

根據(jù)上限定理：與機(jī)動(dòng)許可速度場(chǎng)對(duì)應(yīng)的所有外荷載中，最小的外荷載與真實(shí)荷載最接近，即：上限法是求一個(gè)外荷載的極小值的數(shù)學(xué)規(guī)劃問(wèn)題。因此，考慮土石壩土體容重超載時(shí)，可設(shè)容重超載系數(shù)為目標(biāo)函數(shù)，求解土石壩穩(wěn)定性的上限法目標(biāo)函數(shù)為：

式中：Minimize表示求“最小”。

2.5 土石壩可靠度分析上限法數(shù)學(xué)規(guī)劃模型

根據(jù)土石壩壩坡極限狀態(tài)函數(shù)式（1），目標(biāo)函數(shù)式（17），約束條件式（6）、（7）、（8）、（16），可以得到隨機(jī)庫(kù)水位作用下的土石壩可靠度分析上限法數(shù)學(xué)規(guī)劃模型：

式（18）是一個(gè)大型的線性規(guī)劃數(shù)學(xué)模型，其中的A2e、A2d、We、Wd、W pe、W pd都與土體的抗剪強(qiáng)度參數(shù)、庫(kù)水位隨機(jī)變量相關(guān)。當(dāng)使用式（2）、式（3）分別將抗剪強(qiáng)度參數(shù)隨機(jī)變量（cr，φr）、孔隙水壓力隨機(jī)變量（Hrw）離散成若干個(gè)隨機(jī)數(shù)［cr(tm)，φr(tm)，Hrw(tw)，其中tw=(1，…，nw)，tm=(1，…，nm)］以后，式（18）相應(yīng)地變成nw×nm個(gè)線性數(shù)學(xué)規(guī)劃模型，通過(guò)數(shù)學(xué)規(guī)劃算法求解可以獲得nw×nm個(gè)安全系數(shù)以及對(duì)應(yīng)的速度場(chǎng)。

3 土石壩的失效概率和可靠度指標(biāo)

在當(dāng)前，土石壩失效分析主要使用整體失效概率法，即根據(jù)土石壩的安全系數(shù)來(lái)求解土石壩的整體失效概率，其原理是：當(dāng)土石壩的安全系數(shù)≥1 時(shí)土石壩是穩(wěn)定的（整體安全）、當(dāng)土石壩的安全系數(shù)＜1 時(shí)土石壩發(fā)生失穩(wěn)（整體失效）。土石壩的失效功能函數(shù)如下：

式中：tw=(1，…，nw)，tm=(1，…，nm)，I(tw，tm)是第tw個(gè)庫(kù)水位作用下對(duì)應(yīng)著第tm個(gè)抗剪強(qiáng)度參數(shù)隨機(jī)數(shù)的土石壩失效功能函數(shù)。

土石壩在第tw個(gè)庫(kù)水位作用下的整體失效概率按下式計(jì)算：

式中：tw=(1，…，nw)，Pzf(tw)是土石壩在第tw個(gè)庫(kù)水位作用下的整體失效概率。

在所有的nw個(gè)庫(kù)水位作用下土石壩的失效概率按下式計(jì)算：

式中：PzF是土石壩在所有可能發(fā)生的地下水位作用下的整體失效概率。

土石壩第tw個(gè)庫(kù)水位作用下安全系數(shù)的均值：

式中：tw=(1，…，nw)，μk(tw)是第tw個(gè)庫(kù)水位作用下土石壩nm個(gè)安全系數(shù)的平均值。

土石壩第tw個(gè)庫(kù)水位作用下安全系數(shù)的標(biāo)準(zhǔn)差：

式中：tw=(1，…，nw)，σk(tw)是第tw個(gè)庫(kù)水位作用下土石壩nm個(gè)安全系數(shù)的標(biāo)準(zhǔn)差。

在所有的nw個(gè)庫(kù)水位作用下土石壩安全系數(shù)的均值：

式中：μk是nw×nm個(gè)安全系數(shù)的平均值。

在所有的nw個(gè)庫(kù)水位作用下土石壩安全系數(shù)的標(biāo)準(zhǔn)差：

式中：σk是nw×nm個(gè)土石壩安全系數(shù)的標(biāo)準(zhǔn)差。

4 土石壩可靠度分析上限法數(shù)學(xué)規(guī)劃模型的數(shù)值求解策略

圖2是隨機(jī)庫(kù)水位作用下病險(xiǎn)土石壩可靠度分析上限法數(shù)學(xué)規(guī)劃模型的數(shù)值求解流程，計(jì)算主要包括3個(gè)部分：前處理模塊、計(jì)算模塊和后處理模塊，具體步驟如下：

（1）確定土石壩壩體幾何參數(shù)、土石壩材料的物理力學(xué)參數(shù)、材料的抗?jié)B參數(shù)和庫(kù)水位信息；

（2）建立土石壩可靠度計(jì)算的極限狀態(tài)函數(shù)；

（3）生成水庫(kù)上游庫(kù)水位水頭的隨機(jī)數(shù)；

（4）生成土石壩壩體材料的抗剪強(qiáng)度參數(shù)的隨機(jī)數(shù)；

（5）使用三角形有限單元離散土石壩，根據(jù)步驟（3）生成的庫(kù)水位水頭信息，進(jìn)行土石壩的穩(wěn)定滲流計(jì)算，獲得有限單元節(jié)點(diǎn)的孔隙水壓力；

（6）根據(jù)式（18）構(gòu)建滿足目標(biāo)函數(shù)、約束條件的病險(xiǎn)土石壩穩(wěn)定性分析上限法線性規(guī)劃模型；

（7）使用“對(duì)偶單純形法”求解病險(xiǎn)土石壩穩(wěn)定性分析上限法線性規(guī)劃模型得到對(duì)應(yīng)的安全系數(shù)；

（8）計(jì)算土石壩安全系數(shù)的均值、標(biāo)準(zhǔn)差以及土石壩的失效概率。

5 工程實(shí)例分析

根據(jù)本文提出的方法，編制了相應(yīng)的上限法程序，并對(duì)云南省某病險(xiǎn)土石壩進(jìn)行了計(jì)算分析，將計(jì)算結(jié)果與LEM 法的計(jì)算結(jié)果進(jìn)行比較，驗(yàn)證了本文計(jì)算方法的正確性。

5.1 云南省某病險(xiǎn)土石壩概況

某病險(xiǎn)土石壩位于云南省中部，大壩為均質(zhì)壩（如圖3所示），壩頂高程為1 873.00 m，最大壩高為22.5 m，壩頂長(zhǎng)度為110 m、寬度為3.5 m，上游壩坡坡比為1∶2.75、下游壩坡坡比為1∶2.60。該壩在初期建造時(shí)采用當(dāng)?shù)夭牧线M(jìn)行筑壩，受歷史原因、投入資金的限制，該土石壩設(shè)計(jì)水平和施工質(zhì)量偏低；經(jīng)過(guò)54年的長(zhǎng)期服役，壩體材料的物理力學(xué)性質(zhì)已經(jīng)嚴(yán)重劣化，壩體的穩(wěn)定性已不滿足現(xiàn)行規(guī)范要求；經(jīng)過(guò)鑒定，該土石壩屬于B類病險(xiǎn)土石壩，需要進(jìn)行綜合治理。為詳細(xì)了解該病險(xiǎn)土石壩的實(shí)際情況，勘察單位在大壩壩體區(qū)域布置了4 個(gè)地質(zhì)鉆孔（ZK1～ZK4）對(duì)土石壩進(jìn)行了地質(zhì)勘察，勘察結(jié)果顯示：土石壩壩體黏土的力學(xué)性質(zhì)較差；壩基地層分為兩層，表層為礫石土，下層為石英質(zhì)粉砂巖，力學(xué)性質(zhì)相對(duì)較好；庫(kù)區(qū)水文工程地質(zhì)條件好，庫(kù)水不存在向庫(kù)外鄰谷滲漏問(wèn)題。表1為土石壩在設(shè)計(jì)初期和長(zhǎng)期服役以后的材料物理力學(xué)參數(shù)統(tǒng)計(jì)表。

表1 土石壩材料物理力學(xué)參數(shù)的設(shè)計(jì)值及當(dāng)前值Tab.1 Design values and current values of material parameters of earth-rock dam

本文采用三角形單元對(duì)土石壩進(jìn)行離散，得到1 235 個(gè)三角形有限單元，3 705 個(gè)不共節(jié)點(diǎn)單元，1 788 條公共邊，壩體網(wǎng)格劃分如圖4所示；建立xoy坐標(biāo)系如圖所示，水平向右為x軸正向，豎直向上為y軸正向，坐標(biāo)原點(diǎn)O的高程為1 834.59 m；圖中P1、P2、P3、P4為孔隙水壓監(jiān)測(cè)關(guān)鍵點(diǎn)。

5.2 隨機(jī)庫(kù)水位和穩(wěn)定滲流場(chǎng)

該水庫(kù)為?。ǘ┬退畮?kù)，校核洪水位1 871.39 m，正常蓄水位1 871.10 m，死水位1 857.85 m，水庫(kù)總庫(kù)容57.2 萬(wàn)m3。該水庫(kù)為季調(diào)節(jié)水庫(kù)，庫(kù)水位受降雨和用水需求的影響在雨季和旱季呈現(xiàn)隨機(jī)變化；根據(jù)該水庫(kù)歷年的水文資料統(tǒng)計(jì)分析，上游庫(kù)水位Hw在死水位到校核洪水位范圍內(nèi)隨機(jī)變化，近似服從正態(tài)分布，本文假設(shè)Hw服從截尾正態(tài)分布規(guī)律。根據(jù)式（2）離散獲得Hw的隨機(jī)數(shù)，取水庫(kù)上游死水位為水位下界，即Hlb=27.56 m（絕對(duì)高程1 862.15 m），取水庫(kù)校核洪水位為水位上界，即Hub= 40.43 m（絕對(duì)高程1 871.39 m），庫(kù)水位均值Huw=37.00 m，變異系數(shù)為0.3，取隨機(jī)數(shù)量nw= 50，上游庫(kù)水位隨機(jī)數(shù)分布情況如圖5所示。

根據(jù)生成的庫(kù)水位隨機(jī)數(shù)，計(jì)算得到該病險(xiǎn)土石壩的50個(gè)穩(wěn)定滲流場(chǎng)。圖6為在低水位、中水位、高水位作用下土石壩的孔隙水壓力等值線分布圖。隨著水位的升高，土石壩的浸潤(rùn)線隨之升高，壩體內(nèi)的孔隙水壓力值也隨之增大。由于該土石壩的壩體未設(shè)置防滲體，浸潤(rùn)線變化平緩，在高水位作用時(shí)，大部分壩體處于飽和狀態(tài)，孔隙水壓力對(duì)壩體的穩(wěn)定性有較大的不利影響。

圖7為土石壩中4 個(gè)關(guān)鍵點(diǎn)的孔隙水壓力隨水位變化關(guān)系圖。隨著水位的增加，關(guān)鍵點(diǎn)的孔隙水壓力值隨之增加。關(guān)鍵點(diǎn)的孔隙水壓力表現(xiàn)出如下規(guī)律：①在相同水位作用下，隨著關(guān)鍵點(diǎn)向下游移動(dòng)，孔隙水壓力值不斷減??；②越靠近上游面的關(guān)鍵點(diǎn)，其孔隙水壓力值變化越激烈，如：最靠近上游的P1關(guān)鍵點(diǎn)的孔隙水壓力值從-0.14 kPa升高至-114.07 kPa，而靠近下游面的P4 關(guān)鍵點(diǎn)的孔隙水壓力值從-0.02 kPa 升高至-30.54 kPa。

5.3 上限法與LEM法的對(duì)比分析

根據(jù)表1可知，該土石壩經(jīng)過(guò)50年的服役，壩體材料的抗剪強(qiáng)度參數(shù)已大幅降低（降幅在22%～25%之間），本文選用上游庫(kù)水位隨機(jī)數(shù)的部分樣本，采用本文上限法和LEM 法對(duì)土石壩在設(shè)計(jì)初期和長(zhǎng)期服役后的可靠度進(jìn)行對(duì)比分析計(jì)算。計(jì)算結(jié)果如表2和圖8、9所示。

表2 可靠度指標(biāo)的計(jì)算結(jié)果Tab.2 Calculation results of reliability index

圖8為Hrw(1) = 27.56 m的低水位作用下土石壩安全系數(shù)累積概率密度曲線，由圖8可知本文上限法與LEM 法在兩種工況下的曲線都比較接近，誤差較小，上限法曲線位于LEM 法曲線的上方，表2為土石壩在Hrw(1)、Hrw(25)、Hrw(50)時(shí)上限法和剛體極限平衡法計(jì)算所得可靠度指標(biāo)。

從圖8可以直觀地看出，不管是設(shè)計(jì)值曲線還是服役后的曲線，UBM 和LEM 法的計(jì)算結(jié)果都非常接近。服役前的均值誤差為0.48%，服役后的均值誤差為5.70%，驗(yàn)證了本文計(jì)算方法的正確性。采用兩種計(jì)算方法所得設(shè)計(jì)值累積概率密度曲線均位于服役后的曲線左側(cè)，即經(jīng)過(guò)長(zhǎng)期的服役后土石壩安全系數(shù)大幅降低，與土石壩的實(shí)際情況相符。

從計(jì)算的結(jié)果可以得出：在服役前后土石壩安全系數(shù)的均值、標(biāo)準(zhǔn)差的上限解均小于LEM 法的計(jì)算結(jié)果，但誤差較小，服役前后采用兩種計(jì)算方法均值的最大誤差為8.81%，標(biāo)準(zhǔn)差的最大誤差為11.11%；此外，服役前后采用本文計(jì)算方法所得安全系數(shù)的均值在低水位時(shí)降低24.04%、中水位時(shí)降低24.35%、高水位時(shí)降低23.95%，服役前后采用LEM 法所得安全系數(shù)的均值在低水位時(shí)降低28.02%、中水位時(shí)降低22.73%、高水位時(shí)降低23.13%。服役前后采用本文計(jì)算方法所得失效概率在低水位時(shí)增加了0.52%、中水位時(shí)增加了2.16%、高水位時(shí)增加了11.04%，服役前后采用LEM 法所得失效概率在低水位時(shí)增加了0.76%、中水位時(shí)增加了4.16%、高水位時(shí)增加了14.32%。服役前，土石壩在低、中、高水位時(shí)失效概率均較低，安全性較高；服役后，土石壩在高水位時(shí)失效概率較高，安全性較低。

圖9為土石壩失穩(wěn)時(shí)上限法速度場(chǎng)，可以看出，采用上限法所得速度矢量圍成的區(qū)域和LEM 法所得臨界滑裂面基本一致，但部分區(qū)域存在差距。這是因?yàn)椴捎肔EM 法進(jìn)行土石壩穩(wěn)定性計(jì)算時(shí)，需要事先假定滑裂面的位置，然后對(duì)假定滑裂面上土體進(jìn)行條分，通過(guò)對(duì)各條塊受力特征進(jìn)行分析，從而得到臨界滑裂面位置，從理論上講，LEM 法所得滑裂面是假定滑裂面，與真實(shí)滑裂面存在差距；上限法通過(guò)有限單元構(gòu)建土石壩可靠度上限法數(shù)學(xué)規(guī)劃模型，采用數(shù)學(xué)規(guī)劃的方法搜索得到壩體的失穩(wěn)區(qū)域，其結(jié)果更加接近于真實(shí)情況。

5.4 長(zhǎng)期服役時(shí)的壩體可靠度分析

該土石壩經(jīng)過(guò)長(zhǎng)期的服役，土石壩材料抗剪強(qiáng)度參數(shù)下降，勘察單位對(duì)該土石壩材料的抗剪強(qiáng)度參數(shù)進(jìn)行試驗(yàn)測(cè)定，根據(jù)勘察發(fā)現(xiàn)壩體土和排水棱體的黏聚力和內(nèi)摩擦角相較于設(shè)計(jì)值均有所降低，其中壩體土和排水棱體的黏聚力降低了22.86%，內(nèi)摩擦角降低了25%。本文假設(shè)黏聚力和摩擦角均服從對(duì)數(shù)正態(tài)分布，樣本數(shù)量取為nm=5 000，根據(jù)以上信息，生成壩體材料抗剪強(qiáng)度參數(shù)隨機(jī)數(shù)如圖10所示。

圖11、12分別土石壩安全系數(shù)概率密度曲線和累積概率密度曲線，圖13、14 為計(jì)算所得土石壩安全系數(shù)的均值和標(biāo)準(zhǔn)差隨水位的變化關(guān)系。根據(jù)計(jì)算結(jié)果可知，安全系數(shù)概率密度曲線基本服從正態(tài)分布，且隨著庫(kù)水位的升高，概率密度曲線和累積概率密度曲線均逐漸左移，土石壩的安全性降低，失效概率逐漸增大，可靠度指標(biāo)逐漸降低。此外，隨著庫(kù)水位的升高，安全系數(shù)的標(biāo)準(zhǔn)差也在逐步減小，概率密度曲線的分布范圍逐步變窄，累積概率密度曲線逐步變陡。

圖15、16為計(jì)算得到的安全系數(shù)均值和標(biāo)準(zhǔn)差直方圖。安全系數(shù)均值在1.55～1.57 區(qū)間內(nèi)出現(xiàn)的頻數(shù)最多，最高頻次為8次；安全系數(shù)標(biāo)準(zhǔn)差在0.267～0.270 和0.283～0.286 區(qū)間內(nèi)出現(xiàn)的頻數(shù)最多，最高頻次為7次。

本文統(tǒng)計(jì)了50個(gè)隨機(jī)庫(kù)水位作用下土石壩的失效概率（如圖17 所示）。從圖中可以看出隨著庫(kù)水位的升高，失效概率從0.52%增加到11.04%，安全性逐步降低。庫(kù)水位在27.56～34.00 m 范圍內(nèi)，失效概率隨水位變化較?。粠?kù)水位在34～38 m 范圍內(nèi)，失效概率隨水位變化開(kāi)始增大；庫(kù)水位在38.00～40.43 m 范圍內(nèi)，失效概率隨水位變化急劇增大。土石壩的失效概率隨庫(kù)水位的升高整體上呈現(xiàn)先平穩(wěn)后上升、再急劇上升的特征。庫(kù)水位主要通過(guò)影響土石壩內(nèi)部穩(wěn)定滲流場(chǎng)從而影響壩體穩(wěn)定性，在水位較低時(shí)，只有壩基和少部分壩體處于飽和狀態(tài)，隨著水位的升高，壩基和壩體處于飽和狀態(tài)的面積增大明顯，這就使得土石壩在水位較高時(shí)，安全性面臨一定的挑戰(zhàn)。對(duì)此，在土石壩實(shí)際運(yùn)行過(guò)程中，應(yīng)該采取以下措施，從而保證土石壩的安全：

（1）在土石壩運(yùn)行過(guò)程中應(yīng)該避免高水位的出現(xiàn)，在洪水期，當(dāng)水位超出38.00 m 時(shí)，應(yīng)該及時(shí)采取相應(yīng)措施，如開(kāi)閘泄水的方式來(lái)降低庫(kù)水位，從而保證土石壩的正常運(yùn)行。此外，由于土石壩在多年的運(yùn)行過(guò)程中，材料的抗剪強(qiáng)度參數(shù)會(huì)降低，土石壩安全性能降低，通過(guò)固結(jié)灌漿的方式可提高材料的強(qiáng)度，提高土石壩安全性能。

（2）由于土石壩在長(zhǎng)期的浸泡和沖刷過(guò)程中，滲透系數(shù)會(huì)顯著增大，壩體內(nèi)部孔隙水壓力值也會(huì)明顯增加，對(duì)病險(xiǎn)土石壩做防滲處理極為必要，如：通過(guò)建立防滲墻的方式，降低壩體內(nèi)部的孔隙水壓力，從而增加壩體穩(wěn)定性；將貼坡排水體改為褥墊排水從而降低壩體浸潤(rùn)線，減少壩體內(nèi)部處于飽和狀態(tài)的面積，從而增加壩體穩(wěn)定性。

7 結(jié) 論

（1）本文將上限定理、有限元離散技術(shù)、數(shù)學(xué)規(guī)劃理論和蒙特卡洛方法結(jié)合起來(lái)，同時(shí)將病險(xiǎn)土石壩壩體材料的抗剪強(qiáng)度參數(shù)和上游庫(kù)水位作為隨機(jī)變量，建立了病險(xiǎn)土石壩可靠度分析的上限法隨機(jī)規(guī)劃模型，并編制了相應(yīng)的計(jì)算程序，獲得了同時(shí)考慮土石壩壩體材料的抗剪強(qiáng)度參數(shù)和上游庫(kù)水位隨機(jī)性的土石壩安全系數(shù)、失效概率分布規(guī)律，為病險(xiǎn)土石壩可靠度分析提供了一種新方法。

（2）本文獲得了病險(xiǎn)水庫(kù)的庫(kù)水位對(duì)土石壩的安全系數(shù)分布的影響規(guī)律，在相同土體的抗剪強(qiáng)度參數(shù)樣本下，隨著庫(kù)水位的不斷升高，安全系數(shù)的均值、標(biāo)準(zhǔn)差均呈現(xiàn)下降趨勢(shì)。同時(shí)，庫(kù)水位對(duì)土石壩的整體失效概率有較大影響，失效概率隨著庫(kù)水位的升高不斷增大，且呈現(xiàn)先平穩(wěn)后上升再急劇上升的特征，對(duì)此，在實(shí)際工程中，針對(duì)病險(xiǎn)水庫(kù)，庫(kù)水位過(guò)高應(yīng)該及時(shí)通過(guò)開(kāi)閘泄水的方式使庫(kù)水位降低，也可以通過(guò)對(duì)病險(xiǎn)土石壩做防滲處理的方式，提高大壩安全性。

（3）本文是將塑性極限分析引入病險(xiǎn)土石壩可靠度分析的有益探索。在實(shí)際工程中，地震和庫(kù)水位是影響土石壩穩(wěn)定性的兩個(gè)重要因素，由于計(jì)算工作復(fù)雜，難度較大，本文研究?jī)?nèi)容僅考慮了庫(kù)水位對(duì)病險(xiǎn)土石壩的影響，未考慮地震的作用。在今后的研究工作中，應(yīng)著重對(duì)地震作用下病險(xiǎn)土石壩可靠度研究工作?！?/p>