基于隨機有限元法的溶蝕壩基系統建壩適應性評價

張社榮，嚴 磊，王 超，孫 博，撒文奇

（天津大學 建筑工程學院，天津 300072）

1 引 言

巖溶地基在我國分布廣泛，巖溶發育區域面積約占全國土地面積的 1/3，隨著工程建設的開展，必將遇到越來越多的復雜巖溶地基。在巖溶地基上修建大壩等建筑物通常會遇到地基承載力不足、不均勻沉降以及抗滑穩定性不足等諸多問題[1]，而如何評價復雜溶蝕壩基條件下建壩適應性則是現階段水利建設的關鍵問題。目前，對于復雜地質特別是溶蝕發育條件下建壩適應性的研究主要集中在兩方面[2]：一是依靠試驗手段，研究復雜壩基巖體地質特性及溶蝕發育對巖體物理力學參數的影響；二是借助有限元方法，研究復雜溶蝕地基條件下壩體的強度、變形和穩定性，從溶蝕壩基系統作用效應方面評價建壩適應性。而利用有限元分析理論研究復雜地基條件下壩體的工作性態是當前評價復雜地質條件下建壩適應性的主流。但有限元分析多以確定性理論為主，難點在于對溶蝕體的模擬。

當前常用的溶蝕地基有限元分析模型為確定性的實體模型，即先對復雜溶蝕地基實體建模，建立溶洞單體確定性的幾何模型和力學模型，然后對不同種類的巖體所對應的單元賦予不同的力學參數，將壩基系統作為不同種類巖體的復合體[3]。但由于地質勘測的局限性，地質工程師只能通過有限的鉆孔資料分析線溶蝕率，同時根據溶蝕發育規律推斷出溶蝕發育情況，繪制的工程地質圖并不能完全準確地表達地下溶蝕的實際發育情況。此外，不同區域和巖層的溶蝕發育方式也不盡相同，溶蝕巖體表現出強烈的非均質性和空間變異性，確定性的有限元方法無法很好地反映這些特點。基于此，在溶蝕壩基的數值模擬方法上，張菊明等[4]采用三維隨機橢球洞體數學模型生成并描述了溶蝕地質體；陳祥軍等[5]利用溶蝕巖體隨機結構模型對小規模溶蝕巖體滲漏進行了評價；劉宏等[6]采用二維有限元模型，考慮了溶孔的空間分布特點，對溶蝕地基穩定性進行了研究。總體上看，壩基溶蝕的模擬方法經歷了從定性到定量的研究過程，采取了各種簡化方法，而應用數學和地質統計學方法，引入隨機場理論、模糊理論等克服確定性模型中對溶蝕體模擬的人為假定，更為客觀真實地模擬和表達出溶蝕展布特征，應是壩基溶蝕體模擬的必然趨勢。

本文基于地質統計資料，在分析溶蝕展布特征的基礎上確定壩基溶蝕帶，建立溶蝕帶材料參數隨機場模型。該模型以溶蝕帶隨機場內溶蝕單元的分布等效模擬溶蝕體的空間分布；視溶蝕體由一個或幾個“溶蝕元”組成，以溶蝕帶隨機場的相關距離近似模擬溶蝕元幾何尺寸。基于APDL語言進行二次開發，采用改進的Monte Carlo隨機有限元法進行數值模擬。最后，研究溶蝕帶隨機場模型中壩體工作性態的統計特性，進而基于該模型評價巖溶壩基建壩適應性。

2 隨機場的基本原理

隨機場的理論[7]最早是由VanMarcke提出的，用來描述自然界中分布無序系統（D.D.系統）。其核心思想是用一個方差折減系數把土體的點變異性和空間平均性質結合在一起，使土體的點性質轉化為空間平均性質。溶蝕帶巖體便可以視為這樣一種分布無序系統。

2.1 相關函數與方差折減系數

對于某個隨機場X(?)，點特性的變異性由標準差σ描述，空間平均特性X(z)的變異性由標準差σz描述。σz將隨空間平均范圍的增加而減小，可以定義一個無量綱函數Γ(z)表示由于空間平均而使該區域方差產生折減的大小，局部空間隨機場的均值等于隨機場的點均值，空間均值方差可表示為點方差乘以方差折減系數Γ2(z)，Γ2(z)即為方差折減系數。

相關函數ρ (τ)描述的是用某點的參數值來預測距離該點為τ的另一點的參數值的準確程度。方差折減函數與相關函數之間的關系為

式中：θ為隨機場相關偏度。

2.2 相關距離求解

VanMarcke指出，如果存在

那么，對于充分大的z，有近似公式：

式中：δ為相關距離，是巖土體的一種基本屬性。在相關距離范圍內，土性強烈相關，大于該距離可以認為隨機場性質基本不相關。當z＞2δ時，就有下面的近似公式：

3 隨機有限元方法

溶蝕巖體的物理力學參數受多種因素影響而具有較強的空間變異性。巖土參數的空間變異性是隨機性和結構性的統一體，可采用隨機場（radom field）理論進行描述[8－9]。將隨機場理論和有限元方法相結合即為隨機有限元法（SFEM）。隨機有限元較確定性有限元有幾個突出的優勢[10]：隨機有限元能夠考慮因變量在某一自變量整個隨機參數空間內的所有梯度值；能夠依靠程序自動實現對自變量統計特性的概率輸入；能夠同時改變所有自變量的輸入值，考慮自變量間的相關性，分析結果更加合理。

在各種隨機有限元法中，蒙特卡羅隨機有限元法具有直觀、精確、信息量大的特點，Shinozuka和Strill首先將該方法引入，并推動了其工程應用[11]。經過許多學者的努力，不斷改進了抽樣方法和縮減方差技術，改善和提高了模擬效率，使得改進后的Monte Carlo隨機有限元更趨于實用，改進比較成功的是Newmann展開隨機有限元和拉丁超立方（LHS）抽樣Monte Carlo隨機有限元[12]。本文基于ANSYS軟件APDL語言進行了二次開發，編制程序，采用改進的拉丁超立方抽樣蒙特卡羅隨機有限元法進行分析，大大提高了模擬效率。

拉丁超立方抽樣是一種多維分層抽樣方法，其原理是：①定義參與計算運行的抽樣數目；②把每一次輸入等概率的分成N列，

③對每一列僅抽取一個樣本，各列中的樣本位置是隨機的。

3.3.4 硫酸核糖霉素硫酸鹽含量與酸堿度結果相關性分析 將所測19批次樣品的pH值與硫酸鹽含量（HPLC-ELSD法）進行數據擬合，結果為y=-4.364 6x+52.579 4（r=0.94，其中 y 為硫酸鹽含量（%），x 為 pH 值，6.0≤x≤8.0）。結果表明，硫酸鹽含量與酸堿度結果呈顯著負相關，硫酸鹽含量在17.66%～26.39%范圍內和溶液酸堿度存在線性關系，此時分子式中對應的硫酸數量n值的范圍為1.0～1.7。可見酸堿度可以達到控制硫酸鹽含量的目的，《中國藥典》標準項目的設置有一定合理性。

本文具體研究思路是：（1）分析研究勘測資料和巖土體試驗數據，確定壩基溶蝕帶區域的展布特征。（2）視溶蝕帶子區域的巖體力學參數為隨機場，進行隨機場空間變異性分析和隨機場離散，建立溶蝕帶隨機場模型。（3）設計輸入隨機變量，計入其統計特性，對隨機場中的單元體賦值。（4）隨機有限元數值分析，獲得壩體響應量的樣本數據及其統計特性，通過溶蝕地基的作用效應評價其建壩適應性。基于隨機有限元的溶蝕地基建壩適應性評價的基本流程如圖1所示。

圖1 基于隨機有限元的溶蝕地基建壩適應性評價流程Fig.1 Flow chart of evaluating the dissolution foundation’s adaptability based on SFEM

4 工程實例

4.1 工程基本參數

通過對多種勘探調查資料的分析歸納，壩基溶蝕巖體分布主要呈現出兩種主要特征[2]，一是，為壩基巖體溶蝕較強烈的層位，溶蝕巖體主要為鈣質礫巖和砂巖；二是壩基巖體主要呈現順層溶蝕特征，共存在4條大的溶蝕帶，水平方向主要分布于河床區域，豎直方向主要集中在廠房壩段下部850～1050 m高程范圍內。壩基主要巖土體類型及溶蝕帶的展布特征如圖2所示。

圖2 壩基溶蝕帶分布特征圖Fig.2 Distribution of the dissolution belt in dam foundation

4.2 空間變異性分析及隨機場離散

選取壩基溶蝕最嚴重的典型壩段作為研究對象，視壩基4條溶蝕帶內巖體力學參數為隨機場，進行隨機場空間變異性分析和隨機場離散。相關距離是表征巖土體空間內兩點有無相關性的臨界距離，是對巖土體力學參數等變量在空間范圍內相關程度的度量。據文獻[13]研究成果，對于同一類巖土體，其不同的巖土參數具有基本相同的變異性和相關距離。本文在進行空間變異性分析時，認為溶蝕帶巖體的不同力學參數（彈性模量、黏聚力、內摩擦角等）的相關距離近似相同。據文獻[14]研究成果，每個隨機場單元可以包含一個或幾個有限元單元，本文隨機場單元網格與有限元單元網格的劃分一致。

計算模型的有限元離散網格如圖3所示，本次計算采用的坐標系為：X向由上游指向下游；Z軸豎直向上；Y軸由壩體右岸指向左岸，符合右手螺旋定則。材料本構模型地基材料采用理想的彈塑性模型，壩體混凝土采用線彈性模型。

圖3 計算模型有限元離散網格Fig.3 Discrete grids of computational model

4.3 隨機參數設計與賦值

若將所有力學參數均作為隨機變量處理，反而會影響模擬效率[15]，因此，文中首先進行隨機參數的設計。通過敏感性分析，只選取對壩體強度和穩定性起主要作用的巖體變形模量、內摩擦角、黏聚力作為主要的輸入隨機變量，計入輸入隨機變量樣本數據的統計特性；其余參數設計為常值。壩基各類巖體的參數見表 1，溶蝕帶的巖體參數隨機變量統計特性見表 2，其中，輸入變量的分布類型及相關參數參考《水利水電工程結構可靠度設計統一標準》，并使輸入變量值盡可能涵蓋溶蝕帶中的各類巖體材料參數。本文計算中采用最不利的確定性滲流場分析模型，簡化考慮溶蝕壩區滲流對應力場的影響。

4.4 對比方案設計

采用對比計算方案研究溶蝕模擬的隨機有限元方法與確定性有限元方法的差別。方案1根據地質統計資料確定溶蝕帶展布規律，視壩基溶蝕帶力學參數為隨機場，采用隨機有限元方法計算，此時，溶蝕帶參數按表2取值，壩體及非溶蝕巖體不考慮其材料隨機性，按表1取值；方案2根據地質工程師按照鉆孔資料及線溶蝕率繪制的地質工程圖，對溶蝕體進行確定性實體建模，采用確定性有限元方法計算，模型各參數均按表1取值。以正常蓄水位工況為典型計算工況，對兩種計算方案壩體響應量進行對比分析，計算方案對比設計如表3所示。

表1 模型參數表Table 1 Model parameters table

表2 輸入變量統計特性Table 2 Statistical properties of input variables

表3 計算方案對比設計Table 3 Contrastive design of calculation schemes

5 結果分析

5.1 溶蝕帶隨機場模型響應量統計特性

基于改進的Monte Carlo方法得到壩體特征部位應力、位移響應量的樣本數據，分析樣本的統計特性：以樣本均值描述壩體響應的集中趨勢；以樣本標準差、方差刻畫壩體響應的離散程度；以樣本偏度和峰度描述特征部位應力位移響應的分布形態。

壩體特征部位應力位移響應的統計特性見表4，其中，下游壩面底部折坡點第3主應力和壩頂順水流向位移的隨機輸出變量的概率分布及頻數直方圖如圖5所示。可以看出，溶蝕帶材料參數隨機條件下，壩體特征部位的應力、位移響應量離散程度較小，樣本數據基本在各自均值左右對稱分布。壩踵處第1主應力均值約為1.29 MPa，壩體主壓應力最大值出現在下游壩面底部折坡點區域，達到-6.25 MPa；位移特征方面，壩體順水流向和豎直向位移最大值均出現在壩頂區域，其中，順水流向位移均值為20.89 mm，豎直方向壩體最大沉降量均

圖4 壩體特征部位Fig.4 Feature zone of the dam

表4 壩體特征部位應力位移響應統計特性Table 4 Statistical properties of stress and displacement response of the dam

圖5 隨機輸出變量的概率分布及頻數直方圖Fig.5 Probability distribution and frequency histogram of output variables

同時，隨著重分析壩頂最大沉降值輸出變量的變化規律。圖6給出了基于隨機有限元方法求得的壩頂豎向位移變化曲線。圖中橫坐標表示隨機抽樣次數，縱坐標表示壩頂最大豎向沉降值。可見，方案1壩頂豎向沉降在-30.14～-31.70 mm之間波動。

基于隨機有限元分析的應力場成果，計算正常蓄水位工況下壩體沿建基面的作用效應函數S(?)與抗力函數R(?)，采用承載能力極限狀態法反算該條件下壩體沿建基面的抗滑穩定結構系數，該典型壩段沿建基面抗滑穩定結構系數樣本值在3.51左右波動，且波動較小。

圖6 基于Monte Carlo法的壩頂豎向位移變化曲線Fig.6 Change curve of vertical displacement at dam crest based on Monte Carlo method

5.2 兩種模型溶蝕效應對比分析

在方案2的確定性溶蝕實體模型條件下，壩踵第 1主應力為 0.85 MPa，壩趾第 3主應力為-1.79 MPa，壩頂的豎向沉降最大值為28.57 mm，壩體沿建基面的抗滑穩定反算的結構系數為5.19。與溶蝕確定性實體模型相比，溶蝕帶隨機場模型條件下，壩踵處的主拉應力和壩趾處主壓應力最大值均有所增大；在位移方面，隨機有限元模擬方法下壩頂沉降也有所增大，影響率在5.49%～10.96%范圍內；溶蝕隨機模型下壩體沿建基面的抗滑穩定結構系數有較大降低，減小了約32.37%。

一方面，與溶蝕體確定性實體模型相比，溶蝕帶隨機場模型對壩體的位移特性和穩定指標產生了一定的影響，影響率分別達到 5.49%～10.96%和32.37%。究其原因，溶蝕帶隨機場模型引入隨機場理論，克服了確定性模擬中對線溶蝕率的人為假定，全面考慮了工程地質圖中因不能完全、準確表達溶蝕的實際分布而造成的確定性模型溶蝕體缺失現象，更合理地揭示了溶蝕展布特征。另一方面，對于本工程而言，不管采用何種溶蝕模型，壩體變形特性和穩定指標均滿足規范要求，壩體的工作性態是安全的；同時考慮到局部軟弱巖帶對壩基荷載傳遞的影響會通過壩基開挖置換等措施得到進一步改善，因此，該工程場區的溶蝕帶不會影響壩基的建壩適應性。

6 結 論

（1）由于受到地質勘探局限性的影響，勘探人員無法對整個地基的溶蝕分布做出詳細的描述。采用隨機有限元方法，以溶蝕隨機場內巖體力學參數的變化等效模擬空間內隨機分布的溶蝕體系統，克服了確定性模擬中對線溶蝕率的人為假定，避免了工程地質圖中因不能完全準確地表達溶蝕的實際分布而造成的確定性模型中溶蝕體缺失現象，能夠在有限元模型中合理模擬和表達溶蝕展布特征。

（2）與溶蝕體確定性實體模型相比，溶蝕帶隨機場模型模擬了地質工程圖上未能表達的溶蝕體，因而對壩體的位移特性和穩定指標產生了一定的影響。針對本文工程實例，影響率分別達到 5.49%～10.96%和32.37%。采用該方法對壩體工作形態和溶蝕壩基建壩適應性的評價更準確。

（3）溶蝕帶材料參數隨機場模型，在合理模擬地基溶蝕的局部化效應基礎上，對壩體響應量的統計特性進行分析，賦予了壩體工作性態及建壩適應性評價的統計學意義。

（4）采用隨機有限元方法，反映了該工程場區區域溶蝕地質對壩基巖體變形參數和抗剪強度的影響，而且能給工程人員認識、判斷溶蝕對壩基巖體的影響程度以直觀印象。

需要指出的是，溶蝕地區建壩，壩基滲流場對大壩工作性態的影響不可忽略，限于篇幅，本文采用最不利確定性滲流場分析模型計入了溶蝕壩區滲流對應力場的影響，結果偏于保守。而建立考慮滲流應力耦合場隨機性的隨機耦合模型，客觀準確地描述巖體的滲流狀態，對于合理評價巖溶地區建壩的適宜性非常重要，這也是下一步研究的重點。

[1] 周建普, 李獻民. 巖溶地基穩定性分析評價方法[J]. 礦冶工程, 2003, 23(1): 4－8.ZHOU Jian-pu, LI Xian-min. Methods of stability analysis in karst foundation[J]. Mining and Metallurgical Engineering, 2003,23(1): 4－8.

[2] 楊燕. 觀音巖水電站壩基砂、礫巖溶蝕損傷力學特性與工程應用研究[碩士論文D]. 四川: 成都理工大學,2009.

[3] 胡小榮, 唐春安. 巖土力學參數隨機場的空間變異性分析及單元體力學參數賦值研究[J]. 巖石力學與工程學報, 2000, 19(1): 59－63.HU Xiao rong, TANG Chun’an. Spatial variation analysis of the random field of mechanical parameters for rock and soil and the parameter estimation of elements[J].Chinese Journal of Rock Mechanics and Engineering,2000, 19(1): 59－63.

[4] 張菊明, 王思敬, 曾錢幫, 等. 溶蝕巖體三維隨機洞體數學模型的設計[J]. 工程地質學報, 2004, 12(3): 237－242.ZHANG Ju-ming, WANG Si-jing, ZENG Qian-bang, et al.Mathematic model design for three-dimensional random caverns of karstified rock[J]. Journal of Engineering Geology, 2004, 12(3): 237－242.

[5] 陳祥軍, 馬鳳山, 王思敬, 等. 溶蝕巖體隨機結構模型建立及其在巖體滲漏評價中的應用[J]. 工程地質學報,2004, 12(2): 194－198.CHEN Xiang-jun, MA Feng-shan, WANG Si-jing, et al.The stochastic structure model of karstic rock and its applications to leakage evaluation[J]. Journal of Engineering Geology, 2004, 12(2): 194－198.

[6] 劉宏, 趙躍平, 鄔相國, 等. 強溶蝕帶巖溶地基穩定性研究[J]. 人民長江, 2009, 40(20): 56－58.LIU Hong, ZHAO Yue-ping, WU Xiang-guo, et al.Stability analysis of intensively soluted karst ground[J].Yangtze River, 2009, 40(20): 56－58.

[7] VANMARCKE E H. Random fields: Analysis and synthesis[M]. Cambridge: MIT Press, 1983.

[8] 龐小朝, 周小文, 溫慶博, 等. 隨機場的模擬及其在堤坡可靠性分析中的應用[J]. 長江科學院院報, 2002,19(4): 27－29.PANG Xiao-chao, ZHOU Xiao wen, WEN Qing-bo.Random fields simulation and its application to calculation of reliability of embankment[J]. Journal of Yangtze River Scientific Research Institute, 2002,19(4): 27－29.

[9] 朱登峰, 高大釗. 土性平穩隨機場的空間統計特性分析[J]. 巖土力學, 2003, 24(3): 455－462.ZHU Deng-feng, GAO Da-zhao. Space statistical characteristics analysis of smooth stochastic field of soil characteristics[J]. Rock and Soil Mechanics, 2003, 24(3):455－462.

[10] 丁浩, 蔣樹屏, 楊林德. 外水壓下隧道圍巖與襯砌的隨機有限元分析[J]. 巖土工程學報, 2009, 31(4): 643－647.DING Hao, JIANG Shu-ping, YANG Lin-de. Stochastic finite element analysis of surrounding rock and lining of tunnels under external water pressure[J]. Chinese Journal of Geotechnical Engineering, 2009, 31(4): 643－647.

[11] SHINOZUKA M, DEODATIES G. Response variability of stochastic finite element systems[J]. Journal of Engineering Mechanics, ASCE, 1988, 114(3): 499－519.

[12] ENGELUND S, RACKWITZ R. A benchmark study on importance sampling techniques in structural reliability[J].Structural Safety, 1993, 12(4): 255－276.

[13] 鐘文華. 巖土參數隨機場特性及最優估計研究[碩士論文D]. 天津: 河北工業大學, 2002.

[14] 傅旭東, 鄧建剛. 土性參數隨機場的模擬和隨機場的離散[J]. 工業建筑, 2002, 32(1): 32－36.FU Xu-dong, DENG Jian-gang. Random field modeling for soil properties and random field discretization[J].Industrial Construction, 2002, 32(1): 32－36.

[15] 胡小榮, 唐春安. 巖土力學參數隨機場的離散研究[J].巖土工程學報, 1999, 21(4): 450－455.HU Xiao-rong, TANG Chun′an. Research on discretion of the random field of rock and soil mechanical parameters[J]. Chinese Journal of Geotechnical Engineering, 1999, 21(4): 450－455.