基于非線性耦合的黏土心墻壩滲透系數(shù)的反分析

基于非線性耦合的黏土心墻壩滲透系數(shù)的反分析

0壩體材料參數(shù)的反演分析中國(guó)水庫(kù)超過(guò)90%是水庫(kù)。水庫(kù)是用松散的土壤、沙子和石頭填充的儲(chǔ)水建筑物。由于水庫(kù)蓄水形成的上游和下游水位差，水流將不可避免地從水庫(kù)和水庫(kù)底部滲透到下游。主要問(wèn)題是洪水破壞，這對(duì)水庫(kù)的安全非常重要。而且很多土壩蓄水運(yùn)行多年,在長(zhǎng)期運(yùn)行中,壩體結(jié)構(gòu)形態(tài)不斷進(jìn)行調(diào)整。因此有必要對(duì)壩體主要參數(shù)進(jìn)行反演分析,以評(píng)價(jià)大壩的工作性態(tài)。有的學(xué)者以土石壩滲流場(chǎng)中各測(cè)點(diǎn)水位的實(shí)測(cè)值與計(jì)算值的最優(yōu)化擬合準(zhǔn)則,建立滲透系數(shù)的反演模型,并以青山水庫(kù)主壩滲流分析為例,根據(jù)多個(gè)滲流觀測(cè)點(diǎn)數(shù)據(jù),反演出大壩各個(gè)區(qū)域的滲透系數(shù);有的學(xué)者采用遺傳算法并結(jié)合Biot固結(jié)有限元數(shù)值法對(duì)一假設(shè)路基的雙層地基的一些主要土性參數(shù)進(jìn)行反演;文獻(xiàn)基于遺傳神經(jīng)網(wǎng)絡(luò)的方法對(duì)一重力壩進(jìn)行了滲透系數(shù)的反演。以往,在分析滲流場(chǎng)時(shí)通常忽視壩體應(yīng)力應(yīng)變對(duì)其滲透水流的作用,在分析應(yīng)力應(yīng)變時(shí)則將滲透水流作為外荷載加以處理。事實(shí)上,壩體的應(yīng)力應(yīng)變將改變土體的孔隙比,從而使壩體各部位的滲透系數(shù)發(fā)生變化,影響其滲流場(chǎng)的分布;與此同時(shí),滲透水流的改變,其滲透力和浮力又影響著壩體的應(yīng)力應(yīng)變狀態(tài)。所以,在對(duì)大壩的有關(guān)材料參數(shù)進(jìn)行反演計(jì)算時(shí)必須考慮應(yīng)變場(chǎng)和滲流場(chǎng)的共同作用。本文提出的基于應(yīng)變場(chǎng)和滲流場(chǎng)耦合的滲透系數(shù)的反演方法能夠更真實(shí)的反映大壩的工作狀態(tài),為科學(xué)評(píng)價(jià)大壩安全性提供有力的依據(jù)。1方程與接力學(xué)的基本方向和解算方法的研究1.1對(duì)復(fù)合方程的研究1.1.1土體結(jié)構(gòu)及滲透系數(shù)的分析在土顆粒及水不可壓縮的條件下,平面問(wèn)題多孔介質(zhì)的滲流支配方程如下:??x[kx(εij)?h?x]+??z[kz(εij)?h?z]=St?h?t(1)??x[kx(εij)?h?x]+??z[kz(εij)?h?z]=St?h?t(1)式中:kx(εij),kz(εij))分別為土體水平、鉛直方向的滲透系數(shù),是應(yīng)變場(chǎng)的函數(shù);h為土體中任一點(diǎn)的水頭高度;St為貯水率。影響多孔巖土介質(zhì)的滲透性能的因素主要有兩個(gè)方面,一方面來(lái)源于土體中流體的性質(zhì),包括流體的密度和粘度等;另一方面是土體骨架的性質(zhì)。經(jīng)過(guò)大量的實(shí)驗(yàn)研究表明,黏性土的滲透系數(shù)的對(duì)數(shù)lgk與孔隙比e之間存在線性關(guān)系,具體的表達(dá)式為:e=a+blgk(2)a=10b(3)b=0.01Ιp+δ(4)e=a+blgk(2)a=10b(3)b=0.01Ip+δ(4)式中:a、b為常數(shù),且都與塑性指數(shù)有關(guān);Ip為土體的塑性指數(shù);δ為與土體的類型有關(guān)的常數(shù),其平均值一般可以取0.05。1.1.2土體本構(gòu)模型的雙場(chǎng)耦合在滲流作用下,以張量形式表示的土體的平衡方程為:σij,j+fi(h)+Xi=0(5)σij,j+fi(h)+Xi=0(5)式中:σij,j為σij對(duì)空間坐標(biāo)的偏導(dǎo);σij為有效應(yīng)力張量場(chǎng);fi(h)為體積力,是滲流場(chǎng)水力梯度的函數(shù),與壩體外的滲透水頭有關(guān);Xi為荷載張量。將滲流場(chǎng)數(shù)值模型、應(yīng)力場(chǎng)數(shù)值模型及滲流與應(yīng)變的經(jīng)驗(yàn)關(guān)系組合在一起,并以矩陣的形式表示,可得雙場(chǎng)耦合的數(shù)值模型為:[k]{h}+[s]{?h?t}={f}(6){dσ}=[D]{dε}=[D][B]{dw}(7)[Κ]{w}={F}+{X}(8)[k]{h}+[s]{?h?t}={f}(6){dσ}=[D]{dε}=[D][B]{dw}(7)[K]{w}={F}+{X}(8)式中:[k]為與滲透系數(shù)有關(guān)的矩陣;[S]為與貯水率St有關(guān)的貯水矩陣;{f}為已知常數(shù)項(xiàng),由已知水頭結(jié)點(diǎn)得出;{dσ}為與土體自重有關(guān)的應(yīng)力增量矩陣;[D]為彈性矩陣(或彈塑性矩陣);[B]為應(yīng)變矩陣;[K]為剛度矩陣;{w}為結(jié)點(diǎn)位移;{F}為與滲透坡降有關(guān)的滲透力矩陣;{X}為荷載矩陣。上述式(6)描述應(yīng)力場(chǎng)影響下的滲流場(chǎng);式(7)描述土體的應(yīng)力一應(yīng)變關(guān)系,可以選擇不同的土體本構(gòu)模型代入,式(8)描述滲流場(chǎng)影響下的力的平衡方程。當(dāng)采用Duncan模型時(shí),式(2)和式(6)-(8)就構(gòu)成了土體在滲流作用下的非線性彈性雙場(chǎng)耦合模型,給定邊界條件和初始條件,便可確定滲流場(chǎng)和應(yīng)力場(chǎng)的耦合效應(yīng)。1.2間接耦合法與直接耦合法在雙場(chǎng)耦合分析方面,目前大致有兩類方法,一是將兩場(chǎng)分開計(jì)算,然后,通過(guò)兩場(chǎng)的交叉迭代達(dá)到耦合的目的,即間接耦合法;另一種是通過(guò)建立以應(yīng)變場(chǎng)和滲流場(chǎng)為未知值的數(shù)學(xué)模型,通過(guò)求得解析解達(dá)到完全耦合的目的,即直接耦合法。根據(jù)本文所研究對(duì)象的復(fù)雜性和特點(diǎn),聯(lián)合求解上述式(6)-(8)有較大的困難,故采用間接耦合法。2基于合作土壤滲透系數(shù)的解決方案2.1土體初始孔隙率的確定由于滲透系數(shù)k是量級(jí)較小的數(shù)字(如10-8～10-4),為了給下面的反演計(jì)算提供方便,同時(shí)提高反演的精度,可由式(2)通過(guò)數(shù)學(xué)推導(dǎo)得到:k=c?k0?10αn0+εv1-n0-εv(9)式中:n0為土體的初始孔隙率;k0為土體的初始滲透系數(shù);εv為土體的體積應(yīng)變;α為與土體的塑性指數(shù)有關(guān)的常數(shù);c為由試驗(yàn)確定的常數(shù)。這樣,α和v即為反演的主要參數(shù)。2.2土石壩滲流量反演模型根據(jù)上面的控制方程和定解條件,并將滲流量數(shù)據(jù)和浸潤(rùn)線數(shù)據(jù)結(jié)合起來(lái),一起作為反演的控制條件,同時(shí)對(duì)控制條件進(jìn)行修正,可建立基于耦合的土石壩的滲透系數(shù)的反演模型:F(α,c)=ω11nn∑i=1|hh*-1|+ω2|QQ*-1|=min(10)式中:h是測(cè)壓管水位預(yù)測(cè)值;h*是測(cè)壓管水位觀測(cè)值;Q是滲流量預(yù)測(cè)值;Q*是滲流量觀測(cè)值;ω1和ω2是權(quán)重。2.3遺傳算法程序參數(shù)反演采用正反分析法,即把參數(shù)反分析問(wèn)題轉(zhuǎn)化為一個(gè)目標(biāo)函數(shù)的尋優(yōu)問(wèn)題。問(wèn)題的正算采用滲流場(chǎng)計(jì)算的有限單元法,反算用遺傳算法。土的本構(gòu)模型選用鄧肯-張模型,土壩壩體的應(yīng)力和滲流計(jì)算都采用有限元法。將反演程序GA與土壩的穩(wěn)定滲流程序正演程序SEEP和非線性有限元分析程序FEACD相結(jié)合,構(gòu)成基于流固耦合的土石壩土體滲透系數(shù)反分析程序。GA作為整個(gè)計(jì)算的主程序,SEEP和FEACD是其中的兩個(gè)子程序。程序采用Fortran語(yǔ)言編寫,FortranPowerstation4.0編譯,在Windows9.0以上操作環(huán)境下運(yùn)行。程序的運(yùn)行過(guò)程:耦合分析的初始應(yīng)力場(chǎng)和初始變形場(chǎng)選用土壩竣工時(shí)的應(yīng)力場(chǎng)和變形場(chǎng),土壩竣工時(shí)的滲透系數(shù)k0作為初始滲透系數(shù),由k0確定的滲流場(chǎng)即為初始滲流場(chǎng)。首先由遺傳算法程序GA隨機(jī)產(chǎn)生幾組參數(shù)值,調(diào)用基于耦合的計(jì)算的子程序即滲流計(jì)算程序SEEP和應(yīng)力應(yīng)變計(jì)算程序,用這幾組參數(shù)值進(jìn)行應(yīng)變場(chǎng)和滲流場(chǎng)耦合分析,通過(guò)迭代算出在各組滲透系數(shù)下浸潤(rùn)線上各點(diǎn)的縱坐標(biāo)值hi和滲流量Qi,然后返回GA程序中,用func=ω11nn∑i=1|hh*-1|+ω2|QQ*-1|計(jì)算各組參數(shù)的適應(yīng)度。其中h*是浸潤(rùn)線上各點(diǎn)縱坐標(biāo)的觀測(cè)值,即測(cè)壓管水位;計(jì)算浸潤(rùn)線和觀測(cè)浸潤(rùn)線上各點(diǎn)的橫坐標(biāo)對(duì)應(yīng)相等。以適應(yīng)度為控制條件,對(duì)這幾組參數(shù)進(jìn)行選擇、交換、突變等操作,得到新的幾組參數(shù),再調(diào)用SEEP和FEACD子程序。如此反復(fù)循環(huán),使待反演參數(shù)的適應(yīng)度逐漸增高,使?jié)B透系數(shù)逐漸逼近于真實(shí)滲透系數(shù)。3工程實(shí)例3.1壩體材料特性分析某水利樞紐主壩為黏土心墻壩,主壩壩頂高程為659.00m,最大壩高39.00m,壩頂長(zhǎng)458.00m。主壩基礎(chǔ)主要地質(zhì)自上而下分別是凝灰風(fēng)化頁(yè)巖、凝灰沙礫巖與凝灰?guī)r。風(fēng)化巖層深淺不一,兩岸深約10m,河床部分深約5m。施工中,已將風(fēng)化巖層挖除,用黏土夯填,并將左岸邊坡削成斜坡。1961年10月通過(guò)竣工驗(yàn)收,于1998年9月1日實(shí)行自動(dòng)化觀測(cè)。在主壩共布置了3個(gè)觀測(cè)斷面,共13支測(cè)壓管,用于觀測(cè)壩體浸潤(rùn)線,其中主壩0+041m斷面4支測(cè)壓管,編號(hào)為A1～A4;主壩0+066m斷面5支測(cè)壓管,編號(hào)為A5～A9。本文以典型斷面0+066為計(jì)算斷面。斷面及測(cè)壓管布置見圖1。典型時(shí)日水位測(cè)值和多年平均水位測(cè)值見表1。本文不考慮由于流固耦合引起的土體的力學(xué)參數(shù)的變化,假設(shè)土體的有關(guān)的力學(xué)參數(shù)是已知的常數(shù),主要反算和土體的滲透性能有關(guān)的參數(shù)。為了減少反算參數(shù)數(shù)量,提高計(jì)算精度,將對(duì)浸潤(rùn)線影響較小的壩殼滲透系數(shù)作為已知條件,主要反演心墻的滲透系數(shù)。從式(9)的推導(dǎo)過(guò)程得反算參數(shù)的取值區(qū)間為α∈[4,6],c∈[1,10],公式(10)中權(quán)重ω1和ω2都取為0.5。壩殼的滲透系數(shù)為1×10-4cm/s,黏土心墻的初始滲透系數(shù)為1×10-6cm/s。通過(guò)大型三軸試驗(yàn),得到壩體材料的計(jì)算參數(shù)如表2。選擇初始種群規(guī)模為20,進(jìn)化代數(shù)400代時(shí),目標(biāo)函數(shù)為0.0028,計(jì)算出α=5.21,a=1.56。將反算結(jié)果帶入公式(9),基于耦合進(jìn)行滲流場(chǎng)的有限元正向計(jì)算,計(jì)算結(jié)果見圖2。從圖2可以看出,用考慮耦合的反算算結(jié)果得到的正算值比較接近實(shí)測(cè)值。3.2反分析結(jié)果的統(tǒng)計(jì)檢驗(yàn)反分析完成后還需對(duì)反演結(jié)果進(jìn)行檢驗(yàn)與評(píng)價(jià),具體方法是將反演結(jié)果輸入有限元正分析模型,計(jì)算出2006年7月20日的水頭計(jì)算值與該時(shí)刻測(cè)得的實(shí)際水頭值進(jìn)行比較,然后檢驗(yàn)二者的符合程度。采用后驗(yàn)差法對(duì)反分析結(jié)果進(jìn)行統(tǒng)計(jì)檢驗(yàn)。假設(shè)實(shí)測(cè)水頭序列為hm(i),預(yù)測(cè)水頭序列為hc(i),二者之間的殘差為:ε(i)=hm(i)-hc(i)(i=1,2,?,n)(11)式中:n為特征點(diǎn)的個(gè)數(shù)。若S1為量測(cè)位移的均方差,S2為殘差的均方差,則有:S21=1nn∑i=1(hm(i)-ˉh)2(12)S22=1nn∑i=1(ε(i)-ˉε)2(13)其中:ˉh=1nn∑i=1hm(i)(14)ˉε=1nn∑i=1ε(i)(15)后驗(yàn)比值c=s1s2(16)誤差概率p=p{|ε(i)-ˉε|＜0.6475S1}(17)根據(jù)上述公式計(jì)算得到后驗(yàn)比值c為0.466,誤差概率p為0.894,說(shuō)明反分析結(jié)果是有效的、合格的。4滲透特性的參數(shù)本文通過(guò)公式推導(dǎo)得到應(yīng)變場(chǎng)和滲流場(chǎng)耦合的關(guān)系式,可以避開對(duì)量級(jí)較小的滲透系數(shù)的求解;通過(guò)比較對(duì)土壩的滲流的正向計(jì)算采用間接耦合法,基于耦合,用遺傳算法對(duì)土壩的滲透特性的有關(guān)參數(shù)進(jìn)行反演,通過(guò)實(shí)例計(jì)算得知