基于随机有限元法的溶蚀坝基系统建坝适应性评价

张社荣，严 磊，王 超，孙 博，撒文奇

（天津大学 建筑工程学院，天津 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.