考慮級(jí)配空間隨機(jī)特性的堆石壩變形應(yīng)力分析

朱 晟，盧知是

(1.河海大學(xué)水文水資源與水利水電工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210098；2.河海大學(xué)水工結(jié)構(gòu)研究所，江蘇 南京 210098 )

近年來，西部地區(qū)一批300 m級(jí)高心墻堆石壩即將開工建設(shè)，高堆石壩的變形控制問題成為設(shè)計(jì)研究的重點(diǎn)。已有研究表明，堆石料的級(jí)配分布顯著影響其物理性質(zhì)和力學(xué)行為，且隨荷載增大，相關(guān)性明顯提高。陳镠芬等[1]通過三軸剪切試驗(yàn)發(fā)現(xiàn)堆石料級(jí)配的分形維數(shù)越大，各圍壓下峰值強(qiáng)度越大；趙婷婷等[2]基于顆粒流方法進(jìn)行數(shù)值試驗(yàn)，認(rèn)為分形維數(shù)可用于量化分析堆石料力學(xué)特性；加力別克等[3]通過三軸壓縮試驗(yàn)，得出分形維數(shù)與其抗剪強(qiáng)度、破壞比、切線模量等鄧肯E-B模型參數(shù)之間存在二次函數(shù)關(guān)系；朱晟等[4]分析了級(jí)配與抗剪強(qiáng)度、剪脹性、壓縮性及顆粒破碎之間的規(guī)律，建立了廣義塑性模型參數(shù)與分形維數(shù)之間的函數(shù)關(guān)系。

目前高壩的變形與應(yīng)力研究成果，絕大部分以各分區(qū)壩料采用單一級(jí)配的假定為基礎(chǔ)。然而堆石料級(jí)配存在分布上的空間變異性，且受自然條件、巖石性狀、填筑過程等因素影響，表現(xiàn)出結(jié)構(gòu)上的不確定性，按單一級(jí)配計(jì)算不能完全反映工程實(shí)際。Vanmarcke[5]用齊次正態(tài)隨機(jī)場(chǎng)模擬土性剖面，建立了巖土參數(shù)隨機(jī)場(chǎng)模型。我國(guó)專家學(xué)者如張征等[6]將地質(zhì)學(xué)方法應(yīng)用于巖土參數(shù)分析，拓展了隨機(jī)場(chǎng)理論及應(yīng)用；蔣水華等[7]提出了考慮空間變異性的非均質(zhì)邊坡可靠度分析方法；劉鑫等[8]采用隨機(jī)極限平衡法發(fā)現(xiàn)邊坡破壞的演化與土體參數(shù)的空間分布密切相關(guān)；楊鴿等[9]在對(duì)筑壩料隨機(jī)場(chǎng)進(jìn)行地震響應(yīng)分析時(shí)，得出忽略材料不確定性的方法會(huì)較大程度低估壩坡危險(xiǎn)性。

已有的隨機(jī)有限元研究主要針對(duì)邊坡穩(wěn)定問題，而對(duì)于級(jí)配參數(shù)空間隨機(jī)分布的研究成果較為稀缺。為此，本文引入地質(zhì)統(tǒng)計(jì)學(xué)方法，結(jié)合施工填筑檢測(cè)資料，模擬級(jí)配參數(shù)的空間隨機(jī)分布，進(jìn)行考慮級(jí)配效應(yīng)的變形與應(yīng)力研究，為高心墻堆石壩的變形控制問題提供新的思路。

1 壩殼料級(jí)配分形分布模型

已有研究成果[10-13]表明，粗粒土在不同的粒度范圍內(nèi)均表現(xiàn)出良好的分形性狀，分形維數(shù)D已成為表征巖土材料物理力學(xué)性質(zhì)的重要指標(biāo)。對(duì)于堆石料，采用顆粒質(zhì)量的分形分布模型，即

(1)

式中：P(di)——小于di的顆粒質(zhì)量占比；di——顆粒直徑；dmax——最大粒徑。

結(jié)合建設(shè)中的兩河口心墻堆石壩檢測(cè)級(jí)配，采用分形分布擬合壩料的顆粒級(jí)配曲線如圖1所示，相關(guān)系數(shù)均大于0.989，基本滿足分形分布。

圖1 兩河口心墻堆石壩填筑級(jí)配曲線Fig.1 Grading curve of Lianghekou rockfill dam

對(duì)堆石料、過渡料檢測(cè)級(jí)配的分形維數(shù)值分別進(jìn)行K-S正態(tài)分布檢驗(yàn)，雙側(cè)相伴概率值p分別為0.395、0.194，均大于0.05，因而各分區(qū)分形維數(shù)D滿足正態(tài)分布。由于心墻、墊層、反濾料分形維數(shù)的變異系數(shù)均小于1%，級(jí)配差異相對(duì)較小，故本文著重研究堆石料、過渡料等粗粒料的空間隨機(jī)分布，假定其他分區(qū)填料為單一級(jí)配，計(jì)算時(shí)參數(shù)作均一化處理。

表1 分形維數(shù)D統(tǒng)計(jì)特性

分形分布模型能夠定性、定量描述已知級(jí)配檢測(cè)點(diǎn)的信息，然而難以掌握空間其他位置點(diǎn)的信息，大壩整體結(jié)構(gòu)仍存在很大的未確知性，需進(jìn)一步研究分形維數(shù)D在壩體內(nèi)的空間分布。

2 級(jí)配空間變異結(jié)構(gòu)分析

地質(zhì)學(xué)中，以空間坐標(biāo)x1、x2、x3為自變量的隨機(jī)函數(shù)Z(x1,x2,x3)=Z(x)稱為區(qū)域化變量[6]，其方差為

Var[Z(x)-Z(x+h)]=E[Z(x)-Z(x+h)]2=2γ(h)

(2)

式中：h——檢測(cè)點(diǎn)間距；E——期望值；γ(h)——變異函數(shù)，用來描述空間點(diǎn)之間的差異性和相關(guān)性。

散布于壩體內(nèi)、與空間位置相關(guān)的級(jí)配參數(shù)即可看作區(qū)域化變量。選取兩河口心墻堆石壩距左岸壩肩335.94 ～ 347.79 m(最大斷面壩段)的檢測(cè)級(jí)配共54組，空間位置如圖2所示。由于兩河口工程尚處于在建階段，故本文研究只應(yīng)用1/3壩高以下的檢測(cè)級(jí)配資料，通過計(jì)算較低高程處分布的壩料測(cè)點(diǎn)的空間相關(guān)性，近似模擬整個(gè)斷面級(jí)配的空間分布。

圖2 大壩填筑檢測(cè)級(jí)配分布(單位：m)Fig.2 Detecting point distribution of dam material gradation (unit: m)

為保證隨機(jī)模擬精度，需試算合理的變異函數(shù)理論模型。利用GS+、Arcgis等地質(zhì)統(tǒng)計(jì)學(xué)軟件計(jì)算變異函數(shù)，根據(jù)樣本值計(jì)算γ(h)的估計(jì)量γ*(h)，用理論變異函數(shù)模型對(duì)γ*(h)做回歸擬合，并選定誤差最小模型作為理論模型。計(jì)算時(shí)需設(shè)置最大步長(zhǎng)和Lag步長(zhǎng)等參數(shù)，最大步長(zhǎng)為分離距離的最大值hmax，通常為L(zhǎng)/2(L為區(qū)域內(nèi)沿某方向的最大尺度)[14]。該堆石壩散布的檢測(cè)點(diǎn)間距最大值為894 m，故選擇447 m作為最大步長(zhǎng)值。由于計(jì)算數(shù)據(jù)量要求對(duì)應(yīng)每個(gè)樣本間距值的點(diǎn)對(duì)至少為30對(duì)，經(jīng)試算選定合理步長(zhǎng)值為65.6 m。以堆石料為例，變異函數(shù)理論模型及模型參數(shù)模擬結(jié)果見表2。

表2 變異函數(shù)理論模型及參數(shù)模擬結(jié)果

區(qū)域化變量的相關(guān)程度有對(duì)應(yīng)的分級(jí)標(biāo)準(zhǔn)[14]，當(dāng)塊金系數(shù)C0/(C0+C)≤25%、=25%～75%、>75%時(shí)，分別對(duì)應(yīng)強(qiáng)烈、中等、微弱的3種空間相關(guān)程度等級(jí)。根據(jù)表2，球狀、指數(shù)、高斯模型的塊金系數(shù)分別為16.15%、8.69%和18.88%，表明分形維數(shù)具有強(qiáng)烈的空間自相關(guān)性。其中，指數(shù)模型離差平方和更小、R2更大、塊金系數(shù)更小，故選定指數(shù)模型作為堆石料分形維數(shù)D的變異函數(shù)理論模型，如圖3所示，橫軸為步長(zhǎng)值限定下的所有點(diǎn)對(duì)距離的平均值，縱軸為對(duì)應(yīng)于離散點(diǎn)的變異函數(shù)值。

圖3 變異函數(shù)指數(shù)模型擬合曲線Fig.3 Fitting curve of exponential model for variogram

計(jì)算得到堆石料級(jí)配參數(shù)的變異函數(shù)指數(shù)模型擬合公式為

γ(h)=0.000 229+0.000 118 9(1-e-h/63.9)

(3)

重復(fù)上述步驟對(duì)過渡料進(jìn)行變異函數(shù)計(jì)算，選擇高斯模型為變異函數(shù)理論模型，公式為

γ(h)=0.001 13+0.011 69[1-e-(h/236.4)2]

(4)

將空間距離h作為主要參數(shù)，借助變異函數(shù)理論模型明確了空間各位置點(diǎn)參數(shù)間的相關(guān)性和依賴性，在此基礎(chǔ)上對(duì)分形維數(shù)D的空間分布進(jìn)行隨機(jī)模擬。

3 級(jí)配空間隨機(jī)分布模擬

3.1 隨機(jī)有限元離散方法

采用中心點(diǎn)法對(duì)壩殼料級(jí)配參數(shù)的空間分布作離散化處理。將最大斷面壩段視為空間幾何域V，并劃分為n個(gè)子單元V1、V2、…、Vn，用Z(x)在Vi(i=1,2,…,n)中點(diǎn)xi的值Z(xi)來表征單元體Vi的特征，即可離散得n個(gè)隨機(jī)變量Z(x1)、Z(x2)、…、Z(xn)。為便于計(jì)算分析，依據(jù)壩體有限元剖分網(wǎng)格對(duì)分形維數(shù)D進(jìn)行空間離散，采用空間8節(jié)點(diǎn)等參單元，共811個(gè)單元、1 625個(gè)節(jié)點(diǎn)，如圖4所示，沿垂直于壩軸線從上游到下游為x軸正向，沿壩體高程方向?yàn)閥軸正向。

圖4 大壩有限元網(wǎng)格剖分示意圖Fig.4 Finite element meshing of the dam

3.2 克里金估計(jì)方法

克里金估計(jì)方法通過對(duì)已有數(shù)據(jù)加權(quán)線性組合，采用拉格朗日乘數(shù)法求解可得分形維數(shù)D在無(wú)偏性條件和估計(jì)方差最小原則下的線性估計(jì)量[15-16]，是較優(yōu)秀的線性無(wú)偏估計(jì)方法。估計(jì)采用克里金估計(jì)方法對(duì)壩體內(nèi)未檢測(cè)的級(jí)配不確知點(diǎn)進(jìn)行隨機(jī)估計(jì)，將估得的Z*(xi)值賦給Vi對(duì)應(yīng)的隨機(jī)變量Z(xi)。

由于克里金估計(jì)方法對(duì)正態(tài)分布的數(shù)據(jù)預(yù)測(cè)精度最高，且分形維數(shù)D服從正態(tài)分布，因此精度能夠得到保證[17]。對(duì)堆石料、過渡料級(jí)配空間離散單元進(jìn)行隨機(jī)內(nèi)插，將某次估計(jì)結(jié)果繪制成空間隨機(jī)分布圖如圖5所示，受空間不同位置點(diǎn)的相關(guān)性影響，級(jí)配參數(shù)表現(xiàn)出局部隨機(jī)性和整體結(jié)構(gòu)性，這是空間變異性的主要特征。由于級(jí)配測(cè)點(diǎn)集中于壩體較低高程，較高高程處的離散點(diǎn)間距h較大，空間依賴性降低，一定程度上影響了較高高程處的模擬精度。

圖5 壩料級(jí)配參數(shù)的空間隨機(jī)分布Fig.5 Spacial random distribution of gradation parameter of dam material

3.3 隨機(jī)模擬精度驗(yàn)證

用交叉證實(shí)法驗(yàn)證空間隨機(jī)模擬的精度，其原理為：將各已知點(diǎn)Z(x)暫時(shí)從數(shù)據(jù)系列中除去，根據(jù)其他已知點(diǎn)測(cè)值、克里金估計(jì)方法和變異函數(shù)理論模型計(jì)算該位置處的估計(jì)值Z*(x)，并將Z*(x)放回?cái)?shù)據(jù)系列，對(duì)各點(diǎn)依次遍歷可得Z*(x1)、Z*(x2)、…、Z*(xn)的新數(shù)據(jù)系列，并將其與原始數(shù)據(jù)系列Z(x1)、Z(x2)、…、Z(xn)進(jìn)行統(tǒng)計(jì)比較。對(duì)圖5模擬結(jié)果進(jìn)行交叉驗(yàn)證，結(jié)果如表3所示，堆石料、過渡料的標(biāo)準(zhǔn)差分別為0.085和0.076，數(shù)值較小，相關(guān)系數(shù)分別為0.785和0.712，認(rèn)為變異參數(shù)模型合理，插值精度較高。

表3 交叉驗(yàn)證精度評(píng)價(jià)

4 心墻堆石壩隨機(jī)有限元靜力分析

4.1 靜力計(jì)算

考慮到單次模擬的隨機(jī)性，應(yīng)重復(fù)以上方法對(duì)分形維數(shù)D的空間分布進(jìn)行多次隨機(jī)模擬，并將模擬結(jié)果應(yīng)用于有限元計(jì)算，直至主要計(jì)算值統(tǒng)計(jì)結(jié)果趨于穩(wěn)定。材料的本構(gòu)模型選取廣義塑性模型。朱晟等[4,18]以國(guó)內(nèi)外典型高壩的試驗(yàn)資料為依據(jù)，建立了廣義塑性模型參數(shù)與分形維數(shù)之間的函數(shù)關(guān)系式，可用于考慮填筑級(jí)配影響的土工計(jì)算。

基于文獻(xiàn)[4,18]中一一對(duì)應(yīng)關(guān)系，可得考慮隨機(jī)特性的廣義塑性模型參數(shù)，并針對(duì)各次模擬進(jìn)行靜力計(jì)算分析。借鑒Griffiths等[17]對(duì)于隨機(jī)模擬穩(wěn)定的判定方法，繪制主要計(jì)算結(jié)果的滑動(dòng)平均值(running average)隨模擬次數(shù)n的波動(dòng)曲線，設(shè)第i次模擬得到的計(jì)算結(jié)果為ri，則滑動(dòng)平均值定義為

(5)

圖6 下游側(cè)位移的滑動(dòng)平均值隨模擬次數(shù)的波動(dòng)示意圖Fig.6 Fluctuation diagram of moving average with simulation number for displacement in downstream

表4 材料隨機(jī)參數(shù)統(tǒng)計(jì)特性

4.2 計(jì)算結(jié)果統(tǒng)計(jì)分析

4.2.1 變形和應(yīng)力

圖7為第50次隨機(jī)模擬得到的位移和應(yīng)力分布。對(duì)500組考慮隨機(jī)特性的靜力計(jì)算結(jié)果做K-S正態(tài)分布檢驗(yàn)，得下游側(cè)位移、上游側(cè)位移、最大沉降量、第一主應(yīng)力、第三主應(yīng)力假設(shè)檢驗(yàn)的顯著性概率分別為0.932、0.474、0.864、0.580和0.973，均高于0.05，滿足正態(tài)分布。

圖7 第50組隨機(jī)有限元方法計(jì)算結(jié)果Fig.7 Results of 50th group of stochastic finite element calculation

確定性有限元方法的計(jì)算結(jié)果與隨機(jī)模擬結(jié)果的統(tǒng)計(jì)均值對(duì)比如表5所示，隨機(jī)方法計(jì)算得到最大沉降量均值、下游側(cè)位移和上游側(cè)位移分別為363.95 cm、40.119 cm和39.072 cm，均大于確定性有限元計(jì)算結(jié)果，相對(duì)誤差為10%左右。由于考慮隨機(jī)性得到的各組計(jì)算值均服從正態(tài)分布，且正態(tài)分布具有對(duì)稱特性，即p(xi，μ)=50%，若忽略堆石壩結(jié)構(gòu)的不確定性，僅按單一參數(shù)進(jìn)行有限元計(jì)算，則有超過50%的可能性會(huì)低估壩體變形和應(yīng)力，不利于探究更為安全經(jīng)濟(jì)的高堆石壩設(shè)計(jì)方案。因而認(rèn)為考慮隨機(jī)特性進(jìn)行穩(wěn)定分析和結(jié)構(gòu)設(shè)計(jì)時(shí)，結(jié)果偏于安全和可靠。

表5 隨機(jī)性有限元與確定性有限元計(jì)算結(jié)果對(duì)比

4.2.2 心墻的拱效應(yīng)

由于拱效應(yīng)的作用將使心墻豎向應(yīng)力降低，從而導(dǎo)致心墻內(nèi)產(chǎn)生裂縫或水力劈裂，合理評(píng)估心墻拱效應(yīng)是高壩變形控制中十分重要的問題。采用拱效應(yīng)系數(shù)來表征拱效應(yīng)強(qiáng)弱程度[19]，定義為

(6)

式中：R——拱效應(yīng)系數(shù)；σy——豎向土壓力；γ——土體容重；hs——上覆土層高度。

繪制確定性方法和隨機(jī)性方法[20-23]的拱效應(yīng)系數(shù)R等值線圖如圖8所示，可見隨機(jī)有限元方法計(jì)算得到的拱效應(yīng)系數(shù)R與常規(guī)有限元方法的等值線分布較為接近，能定性表征實(shí)際工程中的拱效應(yīng)問題，驗(yàn)證了隨機(jī)模擬方法的合理性；兩種方法計(jì)算結(jié)果均呈現(xiàn)出兩側(cè)大、中間小的分布規(guī)律，即位于心墻中部的土體拱效應(yīng)更為顯著，相對(duì)于兩側(cè)更易于產(chǎn)生裂縫或發(fā)生水力劈裂；從高程分布上看，位于心墻1/2高程處拱效應(yīng)系數(shù)R值較低，拱效應(yīng)作用更為顯著；拱效應(yīng)系數(shù)R越小，表示拱效應(yīng)越強(qiáng)。隨機(jī)性方法計(jì)算得到的拱效應(yīng)系數(shù)最小值為0.664 5，較傳統(tǒng)方法計(jì)算得到的最小值0.688 0偏小，表明應(yīng)用隨機(jī)性方法進(jìn)行設(shè)計(jì)計(jì)算，更有利于壩體結(jié)構(gòu)安全與穩(wěn)定。

圖8 拱效應(yīng)系數(shù)R分布等值線Fig.8 Contour map of arch effect coefficient

5 結(jié) 論

a.分形分布可較好地表征級(jí)配空間分布，選定指數(shù)型、高斯型變異函數(shù)理論模型分別描述堆石料、過渡料級(jí)配分形維數(shù)D的空間變異結(jié)構(gòu)，且塊金系數(shù)均小于或等于25%，表現(xiàn)出強(qiáng)烈的空間相關(guān)性和依賴性。

b.采用中點(diǎn)映射法進(jìn)行隨機(jī)場(chǎng)模擬與離散，結(jié)果表明，隨機(jī)模擬達(dá)500次時(shí)，計(jì)算結(jié)果趨于穩(wěn)定。

c.通過工程實(shí)例驗(yàn)證，級(jí)配分布的隨機(jī)特性對(duì)高堆石壩應(yīng)力、變形影響明顯。相較于材料參數(shù)均一的確定性有限元方法，考慮級(jí)配參數(shù)空間隨機(jī)特性時(shí)，沉降更大，拱效應(yīng)更明顯，有利于探討更為安全經(jīng)濟(jì)的高堆石壩設(shè)計(jì)方案。