隧道围岩空间变异性是围岩的固有特性，隧道开挖卸载时扰动破坏影响程度随径向距离不同而不同，更加剧了这种围岩的空间变异性．隧道位移的传统设计计算方法通常是把埋深相近和围岩等级相同的地段作为一个设计单元(或设计区间)，采用力学参数的平均值或统计计算确定值作为该单元段典型断面的力学参数进行计算，不能很好地反映空间变异性，而且围岩的这种空间变异性通常很难真实确定，具有不确定性，因此传统设计计算位移值的误差很大，围岩力学特性的空间变异性、不确定性往往成为隧道破坏和变形预测不准的主要因素，也是研究的难点．围岩的这种空间变异性通常与地质形成过程有关，有一定随机性规律，采用数值模拟仿真技术可以全面分析这种空间变异性对位移场的影响规律，因而逐渐成为研究热点．由于土质围岩的离散性相对较小，因此针对土质围岩的研究相对较少[1-4]．针对岩质围岩的空间变异性研究相对较多，文献[5-6]采用等效参数法分析了围岩空间变异性对变形的影响规律．文献[7]结合可靠度理论，分别研究了岩体参数地质强度指数(GSI)、单轴抗压强度σi和完整岩石性质mi的平均值和变异系数、几何参数和空间变异性参数对隧道的可靠度指标和破坏概率的影响规律．文献[8]通过收集与整理203篇隧道工程科研论文和工程案例，对隧道围岩力学参数、荷载参数和支护构件参数的变异系数和概率分布类型及功能函数设计表达式进行统计分析．文献[9]分析了黏聚力c、内摩擦角φ、变形模量E、衬砌厚度t、围岩稳定的极限应变率εlim的变异性对失稳概率的敏感程度．文献[10]提出了一种截断的高斯函数模拟岩石特性的随机场，分析了岩石空间变异性对隧道变形、弹塑性界面、地面反应曲线和破坏机制的影响规律．文献[11]基于蒙特卡罗随机有限差分法，将围岩密度、弹性模量和内摩擦角等参数视为三维正态随机场，分析了围岩的竖向与水平向相关距离的变化对隧道结构可靠度的影响．总体上，目前基于数值模拟研究隧道变形的较多，考虑空间变异性的理论计算方法研究成果较少[12-13]．综上所述，借助于数值仿真技术，采用随机理论进行隧道位移可靠度的设计成为有效的方法，因此本研究根据实际围岩弹性模量的空间变异性特征，兼顾其不确定性，分别提出了围岩空间变异性随机有限元模拟方法和理论简化计算方法，揭示了围岩空间随机分布特征下的位移变化规律，提出了根据围岩弹性模量空间变异性不同标准差(或变异系数)的4分区隧道变形量简化计算方法，为实际隧道工程预留变形量可靠度设计提供了理论依据．围岩力学参数满足其他概率分布函数情况的研究思路和方法相似．1 隧道位移分区计算法1.1　围岩力学参数分布特征根据实验统计，地质年代相同、岩体埋深相近、围岩种类及等级相同的岩体，其力学参数通常服从一定的概率分布形式，根据现场实测值统计拟合[14]和文献[10]等的研究，不考虑层状等特殊结构性岩层的围岩服从截断式正态分布的较多，截断式正态分布的概率密度函数为f(Mi;μMi,σMi,a,b)=1σMiφMi-μMiσMi∙Φb-μMiσMi-Φa-μMiσMi-1，(1)式中：Mi,μMi,σMi,a,b分别为满足截断式正态分布的断面上任意的点i的岩土体力学参数M设计值、断面上各点i的岩土体力学参数M平均值、标准差、下限和上限；φ(⋅)为标准正态分布的概率密度函数，φ(ξ)=(1/2π)exp(-ξ2/2)；Φ(⋅)为标准正态分布的累积分布函数，Φ(ξ)=12π∫-∞ξexp(-t2/2)dt．1.2　随机有限元模型构建根据岩土体力学参数(如弹性模量、强度指标等)空间随机分布特征，划分单元大小(围岩扰动区对稳定性影响最大，因此本研究主要针对扰动区进行空间随机规律研究)，将围岩岩体离散成各个单元(单元大小根据围岩空间变异性特征确定，理论上越小越接近真实情况)，将相似力学性质的岩体单元归为一类，采用同一力学参数，因此可把所有单元分为N类(理论上N越大越精确)．针对实际设计阶段，考虑到设计单元段内围岩的空间变异性并不能完全确定，还具有不确定性特征，因此通过编制子程序将每个单元的力学参数进行随机赋值，但所有单元中N类单元的数量总数(即单元总体积，单元划分时使每个单元体积相等)始终满足围岩空间变异性统计函数，因此可获得一个围岩力学参数空间变异性满足统计随机分布函数的有限元模型．1.2.1　随机有限元模型构建及本构模型选取圆形隧道开挖半径r0=2.5 m，左右、上下边界距圆心65 m，前后上下左右六面均为滑移边界；纵向长30 m，不考虑掌子面的空间效应(一次性开挖，即进行2.5D模拟)，主要扰动区围岩空间随机分布单元模型示意图见图1(a)．10.13245/j.hust.211220.F001图1围岩空间随机分布模型及位移统计路径图本研究采用Extended Drucker-Prager本构模型(EDP模型)，选择线性屈服准则、线性流动法则，参数的选择采用与均匀地层的摩尔库伦理论计算隧道位移值相等的等效原则确定，具体换算参数取值为：压力敏感参数[15]α=1.2，材料屈服应力σy=1.04 MPa，流动势敏感性系数α¯=0.6．为了更好地探究规律性，采用弹性模量单影响分析，其他参数不变：围岩泊松比μ=0.25，黏聚力c=0.5 MPa，内摩擦角ϕ=30°，容重γ=25 kN/m3．1.2.2　位移可靠度设计思路每次赋值均可获得一种围岩空间随机分布的有限元模拟模型，在该岩土体力学参数M空间随机分布工况下，每种单元数量保持不变，只是不同类型单元的空间分布不同．假设进行P次计算，当P趋向无穷大时，计算的工况将涵盖所有的围岩空间随机分布可能存在的情况．若对某区段洞壁处所有点的位移全部统计，则可获得该区段可靠度典型设计断面的位移概率密度函数和累积概率密度函数，从而获得不同概率下的可靠度典型设计断面的设计位移值(或称为典型设计断面的可靠度设计位移值)，本研究建议采用该值作为设计值．因此本研究把埋深相近和围岩等级相同的地段作为一个设计单元段，统计该断面上洞壁全部点的位移值(本研究为了简化计算量，选取8条路径上每隔0.2 m的点的位移值进行统计分析，如图1(b)所示的L1~L8共8条路径)．1.3　位移理论计算法根据深埋隧道围岩受力特点，提出基于环向分层总和法的位移简化计算方法思路如下．a. 针对圆形隧道设计单元段内的空间不均匀分布情况，把径向距离相同的所有围岩的某一力学参数沿环向进行统计，获得每一环的弹性模量等力学参数随机分布函数．b. 简化为平面隧道围岩位移计算模型，此时每环的弹性模量均为随机分布函数．采用随机抽样的方法，每次抽样可获得一个抽象的每环弹性模量都不相同的隧道典型断面，实际操作时，围岩的取样或波速测试沿径向总是有间距，因此引入分层总和法思想，理论上分层越多越接近真实值，获得抽象的典型断面后即可计算该断面的位移．c. 进行N次随机抽样计算(N越大越收敛于真实情况)，对N次随机计算的位移值进行统计分析，即可获得考虑围岩空间变异性的隧道位移概率密度曲线和累积概率曲线．当然这种简化为参数的随机计算忽略了围岩的空间分布性，理论计算值比三维随机有限元的误差更大些，但该计算方法简单，易为普通设计者应用，因此可以修正后作为设计参考．本研究结合修正芬纳公式进行考虑围岩弹性模量空间随机分布的理论计算方法研究，修正芬纳公式通常把围岩分为弹性区和塑性区(分层总和法的极限状态，即只有塑性区一环)，因此对设计单元段内围岩力学参数M(本研究仅分析弹性模量，其他参数可类似分析)进行总的统计分析．假设其统计规律满足相同的截断式正态分布函数(见式(1))，建立基于修正芬纳公式的位移简化计算方法如下．把统计的弹性模量随机分布函数代入修正芬纳公式计算洞壁位移，其中位移ΔR与弹性模量Ei存在对应的函数关系，ΔR=f'(Ei)，函数f'即为修正芬纳公式．由于修正芬纳公式中没有考虑塑性的非线性，通过计算发现计算的位移yi也近似服从截断式正态分布函数，位移yi的累积分布函数为F(yi)=∫-∞yif″(t,μt,σt,c,d)dt. (2)根据N次随机抽样计算的位移累积概率曲线即可进行位移可靠度设计．2 隧道位移计算案例分析2.1　有限元法和理论简化计算方法的合理性验证针对理论简化计算方法，以拟定围岩岩体弹性模量的平均值μEi=200 MPa为例，将在不同标准差σEi下电脑程序随机生成的弹性模量进行频数分布统计和正态分布拟合，理论计算法的分布弹性模量概率密度值f1(σEi)结果见图2，弹性模量N次抽样的统计规律符合实际的截断式正态分布规律，平均值误差Δ均在0.1%以内，故生成数据精度可控．10.13245/j.hust.211220.F002图2理论计算法在不同标准差下围岩弹性模量概率密度曲线针对三维随机有限元法得到的不同类型(本研究算例将围岩单元类型根据弹性模量分为9类)单元数量(由于单元大小相等，本质上即为不同类型单元的总体积空间分布规律)标准化后的统计分析结果显示：不同标准差时不同类型单元的占比均满足截断式正态分布规律，与拟合曲线的误差均小于0.1%，说明理论计算程序和有限元程序生成的模型是基本可靠的．2.2　有限元位移计算结果分析本研究算例通过不同类型单元的数量变化，即单元数量的拟合曲线标准差变化，模拟围岩弹性模量的空间变异性分布规律变化．对每种标准差控制下生成500次(通过试算，100次以上8条路径上各点的位移统计值基本趋于收敛)的计算位移结果进行统计分析．图3为某一次随机生成的围岩9种类型单元空间随机分布图．10.13245/j.hust.211220.F003图3随机生成的围岩单元空间随机分布示意图图4和图5为位移计算值Dw的概率密度曲线和10.13245/j.hust.211220.F004图4随机有限元法在不同标准差下围岩空间随机分布计算位移值的概率密度曲线10.13245/j.hust.211220.F005图5随机有限元法在不同标准差下围岩空间随机分布计算位移值的累积概率曲线累积概率曲线图，图中：f2(σEi)为标准化的围岩空间随机分布的洞壁位移概率密度值；F2(σEi)为围岩空间随机分布的洞壁位移累积概率值(随机有限元法)；σ'为高斯拟合的弹性模量平均值．随着标准差的增大，位移的概率分布曲线逐渐平缓，且分布区间逐渐增大．计算位移值的分布区间(45.5～69.0 mm)位于根据围岩均匀分布的弹性模量上限a和下限b计算的位移值29.286 7～208.843 0 mm之间．图6为不同标准差下理论计算法的位移计算值Dw的累积概率曲线，图中F1(σEi)为围岩空间随机分布的洞壁位移累积概率值(理论计算法)，计算位移值的分布区间32～202 mm，规律与三维随机有限元法计算结果基本相似． 10.13245/j.hust.211220.F006图6理论计算法在不同标准差下围岩空间随机分布计算位移值的累积概率曲线2.3　两种方法特征值对位移分布规律的影响分析图7和图8分别为随机有限元法标准差-洞壁位移值函数曲线及理论计算法标准差-洞壁位移值函数曲线，图中：UCir为随机有限元法数值模拟的洞壁位移值；UCiF为基于修正的芬纳公式的理论计算方法的洞壁位移值．10.13245/j.hust.211220.F007图7随机有限元法标准差-洞壁位移值函数曲线10.13245/j.hust.211220.F008图8理论计算法标准差-洞壁位移值函数曲线2.3.1　随机有限元法洞壁位移计算值的分布规律图7显示了不同标准差条件下位移统计特征值分布情况，大致可分为四个阶段：均匀化线性离散阶段，随机波动性离散阶段，非线性离散阶段，单一化线性离散阶段．当弹性模量正态分布拟合标准差μEi25 MPa(变异系数为12.5%)时，围岩洞壁位移统计值处于均匀化线性离散阶段，可靠度0.95对应的位移UC0.95与位移统计平均值UCu差最大占UCu的1.70%，UCu与假设围岩均匀分布时采用平均弹性模量200 MPa计算的围岩位移值UCm的差值占UCu的1.06%．可见此时考虑空间变异性的位移计算值离散度很小，位移设计值可以采用围岩弹性模量平均值计算后适当进行线性系数的修正即可，或者采用极限值方法，即以变异系数临界值12.5%的计算结果进行设计．随着弹性模量标准差的增大，位移的离散程度逐渐增大，当标准差达到112.6 MPa时，位移的离散程度又逐渐平稳，虽保持在较高水平，但此时弹性模量标准差的增加对计算位移值的影响已经不大，同样可以采取极限思想，以图7中曲线收敛值作为设计值，或通过计算均匀化处理的围岩岩体得到围岩位移值后进行修正．本案例计算结果显示：可靠度0.95所对应的位移在该阶段为假设围岩均匀性工况位移的1.1692倍，因此修正系数可取1.2．当处于随机波动离散阶段和非线性离散阶段时，有必要考虑围岩空间分布变异性对洞壁位移的影响，不能简化计算，否则误差较大．对于随机波动性离散阶段，并无明显相关性表示方式，建议采用外包络线法结合插值法确定位移设计值；对于非线性离散阶段应单独进行随机有限元的分析，当精度要求不高时也可采用拟合插值法确定位移设计值．2.3.2　理论计算法位移计算值的分布规律分析图8显示采用简化的理论计算法的位移计算特征值大致可分为三个阶段，即均匀化线性离散阶段、线性离散阶段和单一化线性离散阶段，未出现随机波动离散阶段．当标准差28.61 MPa(即变异系数为14.3%)时处于均匀化线性离散阶段，可靠度0.95对应的位移UC0.95与位移统计平均值UCu差最大占UCu的19.22%．随着弹性模量分布函数的标准差增大，计算位移特征值的离散程度逐渐增大，当标准差达到102.3 MPa时，位移的离散程度趋缓，虽保持在较高水平，但此时弹性模量分布函数的标准差的变化对计算位移特征值的影响大大减小．基本的规律与三维随机有限元法计算结果基本相似，当处于线性离散阶段时，有必要考虑围岩空间分布变异性对洞壁位移的影响，应进行单独计算，当精度要求不高时也可以采用拟合后线性插值快速计算．图9显示当假设围岩均匀分布时(即不考虑围岩的空间随机性)，随机有限元法计算值与理论计算值位移基本相等；当考虑围岩的空间随机性时，随着变异系数的增大，随机有限元计算结果逐渐小于基于修正芬纳公式的随机理论计算结果；曲线形态总体与图7相似．在弹性模量拟合变异系数为12.5%，23.0%，56.0%(标准差分别为25，46，112 MPa)处，曲线走向出现了明显的转折；当变异系数cv12.5%时，两种方法的计算结果相近，误差约为10%；当变异系数56.0%时，两种方法计算结果差距很大，比值趋于一个稳定的值，此阶段两种方法的误差接近50%．当cv介于12.5%～23.0%时，计算结果的比值在高速下降的同时呈现出波动趋势；当cv介于23%～56%时，函数呈现出非线性下降的趋势．在此阶段，两种计算结果之间的差距较大，约为30%～50%．理论上围岩空间随机分布的随机有限元法计算的结果更接近实际，简化为参数随机计算的方法无法反映围岩沿轴线的空间分布特性，由于两条曲线规律相似，理论计算方法简单快捷，因此可以采用理论计算方法结果乘以拟合经验系数的方法进行修正，作为快速设计值的参考．10.13245/j.hust.211220.F009图9两种方法的位移统计特征值比较3 结论a. 本研究考虑围岩空间变异性统计规律及其统计样本的有限性，引入随机理论，提出了考虑围岩力学参数空间变异性的隧道位移三维随机有限元模拟分析方法和理论计算方法，两者计算结果的规律基本相似，为实现隧道预留变形量的可靠度设计计算提供了理论依据．b. 本研究结合案例分析揭示了围岩弹性模量空间变异性对洞壁位移变化的影响规律，结果显示围岩弹性模量空间随机性分布函数的标准差对围岩位移的可靠度设计值有很大影响，随着标准差(或变异系数)的增加，可以把位移可靠度设计值分为四个阶段，即均匀化线性离散阶段、随机波动性离散阶段、非线性离散阶段和单一化线性离散阶段，具有不同阶段特征的围岩可以采用不同的简化计算方法．c. 考虑到理论简化计算方法便于操作和推广应用，本研究通过对比分析两种方法的计算结果，提出了结合分区计算方法，即采用理论计算值乘以分区影响系数的简化方法．总体上，采用考虑空间变异性的设计方法更科学合理．