一种利用数值模拟分析非均匀分布裂隙巷道稳定性的方法与流程

本发明涉及采矿技术领域，尤其是一种利用数值模拟和现场监测方法分析非均匀分布裂隙巷道稳定性的方法。

背景技术：

岩土工程的隧道、矿业、水电、核电等领域，尤其是地下工程中，岩体裂隙发育具有随机性、形态多样性和空间组合复杂性等特点，从而使得工程岩体的力学参数如弹性模量等具有较高的非均质性，岩体参数的非均质性对工程稳定产生重要影响，严重时会造成井下巷道等支护工程的意外失效。研究非均匀裂隙围岩稳定性控制是岩土工程领域的重要研究课题，对保证安全、高效的施工和生产具有重要指导意义，但是常规的数值模拟方法往往将岩体视为均匀介质，忽略岩体的非均匀性对工程的影响，因此无法对非均匀裂隙发育巷道围岩控制提供准确可靠的支持。

目前，分析围岩失稳机理、研究围岩控制对策、设计支护方案等领域，岩土工程中广泛采用的数值模拟方法主要有：(1)有限元法，有限元法可以将不规则的网格单元拼接成不同形状的研究对象，有限元法应用性强且计算精度较高，但是计算过程复杂、计算量庞大，因此在模拟复杂的工程问题时有限元法存在一定的局限。(2)有限差分法，有限差分法用均匀的网格划分求解域，用有限个离散点的差分代替微分，具有建模方便计算效率高的优点，尤其适用于范围和规模较大的工程模拟。(3)离散元法，离散元法是解决不连续介质问题的数值模拟方法。该方法将裂隙岩体视为由完整岩体单元和节理面组成，单元可以发生平移、转动和变形，节理面允许被压缩、分离或滑动。离散元法适用于模拟节理面空间分布结构与节理面参数明确的裂隙岩体。但是节理空间分布形态以及参数很难进行合理的估算，而且计算维问题时计算量与存储量巨大。有限元法和有限差分法均采用连续介质力学。离散元法中节理的力学参数需要用户自行定义，但这部分力学参数很难准确估算，节理参数的选取不当会极大影响离散元法的计算结果，无法对工程实际进行准确的模拟和仿真。在不专门针对裂隙发育、扩展的围岩稳定性数值模拟中，基于连续介质的有限元法和有限差分法凭借其较高的计算效率和准确的计算结果较好。

由于岩体中裂隙的产状和发育程度影响岩体的强度与变形参数，岩体裂隙发育状态严重影响巷道的变形规律，但是由于岩体裂隙发育的非均匀性造成不同区域巷道围岩力学参数差异性明显，围岩力学参数的非均匀分布对围岩失稳机理、支护方案设计等产生的影响是不可忽略的。因此，在连续介质数值模拟方法的基础上，并在保证计算准确性和高效性的前提下，考虑围岩裂隙的非均匀分布对岩体参数的影响，并提出反演岩体非均质性参数的方法，从而进一步提高数值模拟过程的合理性和结果可靠性。

技术实现要素：

为提高围岩裂隙的非均匀分布条件下的巷道围岩稳定性数值模拟的合理性和可靠性，本发明提供了一种利用数值模拟分析非均匀分布裂隙巷道稳定性的方法，具体技术方案如下。

一种利用数值模拟分析非均匀分布裂隙巷道稳定性的方法，其特征在于，步骤包括：

步骤一.利用有限元模拟软件建立三维数值模型，包括数值模型的尺寸设计和网格划分；

步骤二.非均匀岩体概率分布参数设定，包括确定概率分布函数，模拟岩体的非均质性；

步骤三.测量并计算巷道围岩的地质强度指标，包括钻探取芯获取岩石质量指标和节理条件，地质强度指标的概率分布拟合；

步骤四.对数值模型的材料参数进行赋值，所述材料参数包括弹性模量、泊松比、内聚力、内摩擦角、抗拉强度、残余内聚力和塑性应变量；

步骤五.设定边界条件，模型开挖计算，所述边界条件包括施加重力模拟初始应力，设置应力边界条件和位移边界条件，所述开挖计算使用有限差分计算方法；

步骤六.提取监测数据，进行巷道稳定性评估和二次支护设计。

优选的是，步骤一中的三维数值模型为长方体，在三维数值模型中部开挖巷道，巷道中部划分为顶板监测区域。

优选的是，步骤二中概率分布函数取韦布尔概率分布函数：

其中u是单个元素的力学参数值，u0是尺度参数，m为材料的均质性系数。

进一步优选的是，步骤二中使用蒙特卡洛法模拟岩体的非均质性，包括：a.在[μmin，μmax]之间随机取一个值μi，计算μi的概率密度函数f(μi)；b.在f(μmin)与f(μmax)之间随机产生一个概率密度函数f(μ)，当f(μ)<f(μi)时，导出μi，c.继续进行迭代计算，最终得到μi的数据组，μi的数据组服从韦布尔概率分布规律。

优选的是，步骤三中的地质强度指标gis与岩石质量指标rqd和节理条件jcond89的关系为：gsi＝1.5jcond89+rqd/2。

进一步优选的是，岩石质量指标rqd取钻探获取的岩芯完整度来表示，具体为长度大于10cm的岩芯占钻孔总长度的比值；所述节理条件jcond89的取值根据巷道围岩裂隙长度、宽度、粗糙度、充填类型和风化程度确定。

优选的是，步骤四中弹性模量e和地质强度指标gis之间的关系为：

其中σm为岩石单轴抗压强度。

优选的是，步骤五中模型开挖使用有限差分法进行迭代计算，当模型的最大不平衡力当低于1e-5时模拟运算结束。

优选的是，步骤六中二次支护设计包括，在巷道变形值达到120mm以上的区域，采取注浆加固和增加锚杆预紧力的支护措施。

本发明的有益效果包括：

(1)本方法结合现场监测和数值模拟对非均质岩体的巷道裂隙分布和顶板下沉规律进行模拟，模拟结果和现场监测相互验证，利用韦布尔概率分布函数模拟岩体的非均质性，提高了模拟的合理性和可靠性。

(2)根据钻探结果确定地质强度指标从而更好的还原了工程现场的巷道围岩条件，保证了模拟结果的准确性，尤其是对地质强度指标gis的统计值进行概率分布拟合，确定了顶板的模拟参数，从而能够更准确的计算巷道顶板下沉变形值。

(3)根据巷道模拟结果，提取巷道变形量监测数据，对巷道围岩的稳定性和可靠性进行评估分析，从而进一步的指导巷道二次支护保证了巷道的安全和支护的有效性。

另外本方法还具有操作简单，计算准确高效，支护可靠性高等优点。

附图说明

图1是利用数值模拟分析非均匀分布裂隙巷道稳定性的方法流程结构图；

图2是三维数值模型结构示意图；

图3是巷道顶底板岩性特征和厚度柱状图；

图4是不同均质性系数的韦布尔分布曲线图；

图5是韦布尔分布参数拟合示意图；

图6是数值模拟结果巷道顶板下沉示意图。

具体实施方式

结合图1至图6所示，本发明提供的一种利用数值模拟分析非均匀分布裂隙巷道稳定性的方法，具体实施方式如下。

实施例1

一种利用数值模拟分析非均匀分布裂隙巷道稳定性的方法，其特征在于，步骤包括：

步骤一.利用有限元模拟软件建立三维数值模型，包括数值模型的尺寸设计和网格划分；

在本步骤中的三维数值模型为长方体，在三维数值模型中部开挖巷道，巷道中部划分为顶板监测区域。

步骤二.非均匀岩体概率分布参数设定，包括确定概率分布函数，模拟岩体的非均质性；

在本步骤中概率分布函数取韦布尔概率分布函数：

其中u是单个元素的力学参数值，u0是尺度参数，m为材料的均质性系数。

其中使用蒙特卡洛法模拟岩体的非均质性，包括：a.在[μmin，μmax]之间随机取一个值μi，计算μi的概率密度函数f(μi)；b.在f(μmin)与f(μmax)之间随机产生一个概率密度函数f(μ)，当f(μ)<f(μi)时，导出μi，c.继续进行迭代计算，最终得到μi的数据组，μi的数据组服从韦布尔概率分布规律。利用韦布尔概率分布函数模拟岩体的非均质性，提高了模拟的合理性和可靠性。

步骤三.测量并计算巷道围岩的地质强度指标，包括钻探取芯获取岩石质量指标和节理条件，地质强度指标的概率分布拟合；

在本步骤中的地质强度指标gis与岩石质量指标rqd和节理条件jcond89的关系为：gsi＝1.5jcond89+rqd/2。岩石质量指标rqd取钻探获取的岩芯完整度来表示，具体为长度大于10cm的岩芯占钻孔总长度的比值。其中节理条件jcond89的取值根据巷道围岩裂隙长度、宽度、粗糙度、充填类型和风化程度确定。钻探结果确定地质强度指标从而更好的还原了工程现场的巷道围岩条件，保证了模拟结果的准确性，尤其是对地质强度指标gis的统计值进行概率分布拟合，确定了顶板的模拟参数，从而能够更准确的计算巷道顶板下沉变形值。

步骤四.对数值模型的材料参数进行赋值，材料参数包括弹性模量、泊松比、内聚力、内摩擦角、抗拉强度、残余内聚力和塑性应变量；

在本步骤中弹性模量e和地质强度指标gis之间的关系为：

其中σm为岩石单轴抗压强度。

步骤五.设定边界条件，模型开挖计算，边界条件包括施加重力模拟初始应力，设置应力边界条件和位移边界条件，所述开挖计算使用有限差分计算方法；

在本步骤中模型开挖使用有限差分法进行迭代计算，当模型的最大不平衡力当低于1e-5时，说明模型达到平衡状态，模拟运算结束。

步骤六.提取监测数据，进行巷道稳定性评估和二次支护设计。

在本步骤中的二次支护设计，在巷道变形值达到120mm以上的区域，采取注浆加固和增加锚杆预紧力的支护措施。根据巷道模拟结果，提取巷道变形量监测数据，从而对巷道围岩的稳定性和可靠性进行评估分析，从而进一步的指导巷道二次支护保证了巷道的安全和支护的有效性。

实施例2

本实施例结合某矿回采巷道的工程实际，对本发明的具体实施做进一步的说明，采用连续介质有限差分软件flac^3d进行数值模拟，根据数值模拟结果指导巷道的二次支护，该方法流程如图1所示，具体实施如下。

在某矿的回采巷道中利用数值模拟分析非均匀分布裂隙巷道稳定性的方法指导二次支护的具体步骤包括：

步骤一.利用有限元模拟软件建立三维数值模型，包括数值模型的尺寸设计和网格划分。

在本步骤中的三维数值模型为长方体，如图2所示，设置模型长度为200m，模型宽度为80m，高度为62m，在三维数值模型中部开挖巷道。巷道为长方形，巷道长度为200m，巷道两侧到边界的距离取35m巷道中部100m的区域划分为顶板监测区域，巷道顶底板岩层根据图3所示的巷道顶底板岩性特征和厚度柱状图设置。

步骤二.非均匀岩体概率分布参数设定，包括确定概率分布函数，模拟岩体的非均质性。

在本步骤中概率分布函数取韦布尔概率分布函数，韦布尔概率分布能够模拟一个区间内的随机分布岩体非均质参数，不同均质性系数m的韦布尔分布曲线如图4所示，韦布尔分布函数如下：

其中u是单个元素的力学参数值，u0是尺度参数，随着尺度参数的改变函数曲线放大或缩小，m为材料的均质性系数，随着均质性系数的改变曲线的形状随之改变。

其中使用蒙特卡洛法模拟岩体的非均质性，步骤包括：a.在[μmin，μmax]之间随机取一个值μi，计算μi的概率密度函数f(μi)；b.在f(μmin)与f(μmax)之间随机产生一个概率密度函数f(μ)，当f(μ)<f(μi)时，导出μi，c.继续进行迭代计算，最终得到μi的数据组，μi的数据组服从韦布尔概率分布规律，将该方法通过fish语言编写，在flac^3d中运行该程序。利用韦布尔概率分布函数模拟岩体的非均质性，提高了模拟的合理性和可靠性。

步骤三.测量并计算巷道围岩的地质强度指标，包括钻探取芯获取岩石质量指标和节理条件，地质强度指标的概率分布拟合。在回采巷道内选取52个测点，对测点的地质强度指标gis进行统计处理，统计处理的步骤包括：a.将地质强度指标gis统计值制作概率分布直方图，如图5所示；b.利用origin软件处理数据，进行韦布尔weibull概率分布参数拟合，得到gsiave＝40.72，均质性系数m＝1.96；c.根据拟合结果取gsimax＝10与gsimin＝80。

在本步骤中的地质强度指标gis与岩石质量指标rqd和节理条件jcond89的关系为：gsi＝1.5jcond89+rqd/2。岩石质量指标rqd取钻探获取的岩芯完整度来表示，具体为长度大于10cm的岩芯占钻孔总长度的比值。其中节理条件jcond89的取值根据巷道围岩裂隙长度、宽度、粗糙度、充填类型和风化程度确定。钻探结果确定地质强度指标从而更好的还原了工程现场的巷道围岩条件，保证了模拟结果的准确性，尤其是对地质强度指标gis的统计值进行概率分布拟合，确定了顶板的模拟参数，从而能够更准确的计算巷道顶板下沉变形值。

步骤四.对数值模型的材料参数进行赋值，材料参数包括弹性模量、泊松比、内聚力、内摩擦角、抗拉强度、残余内聚力和塑性应变量，各层的材料参数如表1所示。

表1.材料参数取值表

在本步骤中弹性模量e和地质强度指标gis之间的关系为：

其中σm为岩石单轴抗压强度，根据该公式计算顶板岩体的弹性模量。

利用flac^3d中的内置cable结构单元模拟锚杆和锚索支护结构，设置的支护结构力学参数如表2所示。

表2.支护参数取值表

步骤五.设定边界条件，模型开挖计算，边界条件包括施加重力模拟初始应力，设置应力边界条件和位移边界条件，开挖计算使用有限差分计算方法。在模型顶部施加15mpa的垂直应力，模拟巷道实际埋深650m的实际情况，在水平方向上分别施加水平应力12mpa和18mpa，模型的四周和底部设置位移边界条件。采用显式有限差分计算法进行迭代计算，模型的最大不平衡力低于1e-5时，说明模型达到平衡状态，停止计算，初始应力场模拟完成。

在本步骤中模型开挖使用有限差分法进行迭代计算，具体是：将开挖区域定义为空模型null，巷道的开挖以掘进10m为一个循环，支护紧随巷道开挖，设置支护体并赋予力学参数，运算至模型平衡后，进行下一个10m的开挖循环。当模型的最大不平衡力当低于1e-5时，说明模型达到平衡状态，模拟运算结束。

步骤六.提取监测数据，进行巷道稳定性评估和二次支护设计。具体是利用hist命令对巷道顶板的变形进行监测，当巷道开挖结束后提取巷道变形量的监测数据，输出巷道顶板下沉结果如图6所示。

在本步骤中的二次支护设计，在巷道变形值达到120mm以上的区域，采取注浆加固和增加锚杆预紧力的支护措施。根据巷道模拟结果，提取巷道变形量监测数据，从而对巷道围岩的稳定性和可靠性进行评估分析，尤其是巷道变形量差异明显的位置，对巷道变形量大的位置进行重点加强支护，从而提升支护效率和支护效果，进一步的指导巷道二次支护保证了巷道的安全和支护的有效性。

当然，上述说明并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的技术人员在本发明的实质范围内所做出的变化、改型、添加或替换，也应属于本发明的保护范围。