CN105550441A - 一种基于连续介质的工程岩体破裂劣化数值模拟方法 - Google Patents

Abstract

本发明公开了一种基于连续介质的工程岩体破裂劣化数值模拟方法，它包括建立工程对象数值模型步骤，选择工程对象数值模型中岩体的本构模型，并赋予该本构模型相应力学参数步骤，设置工程对象数值模型边界条件步骤，对工程对象数值模型进行开挖并执行工程岩体破裂劣化迭代计算，获取模拟结果，根据研究对象进行研究分析。本发明在保留原始数值迭代计算基础上，加入了围岩中应力重新分布过程中裂隙受开挖扰动的动态发展及导致围岩杨氏模量劣化，由现场实地观测或理论经验估算所得数据反演得到劣化后残余杨氏模量，能够使基于连续介质的数值模拟过程更加合理严谨，进而使数值模拟结果更加真实可靠。

Description

一种基于连续介质的工程岩体破裂劣化数值模拟方法

技术领域

本发明涉及岩石力学与工程、采矿工程领域，特别是关于裂隙岩体围岩稳定性分析与控制设计优化数值模拟技术领域。

背景技术

在岩土工程的隧道、矿业、水电、核电等领域，尤其是地下工程中，人为开挖扰动打破了围岩的原始应力状态，引起围岩中应力的重新分布。这一过程使得岩石中原生的宏观和微观裂隙发育、扩展甚至彼此贯通，使围岩的力学性质劣化并发生不同程度的变形和破坏，从而影响围岩的稳定性。围岩的稳定性控制是岩土工程领域的重要研究课题，对保证安全、高效的施工和生产有着重要意义。

近年来在分析围岩失稳机理、研究围岩控制对策、设计支护方案等领域，数值模拟研究方法被越来越多的学者所采用。当下岩土工程中广泛采用的数值模拟方法主要有如下几种：

1、有限元法

有限元法可以灵活的划分不规则的网格单元以模拟不同形状的仿真对象，具有广泛的应用性和较高的计算精度，但其计算过程复杂、计算量庞大，使其在进行较大范围和规模、介质复杂工程的模拟中存在局限。

2、有限差分法

有限差分法是目前岩土工程领域应用较广的计算发法，用均匀的网格划分求解域，用有限个离散点的差分代替微分，具有较高的计算效率，适用于较大范围和规模、介质复杂的工程模拟。

3、离散元法

离散元法是解决不连续介质问题的数值模拟方法。该方法将模拟对象视为由离散的岩体单元和单元间的节理面所组成，允许单元的平移、转动和变形，而节理面可被压缩、分离或滑动。离散元法适用于模拟已知节理结构的对象，但节理的参数很难合理估算且解决三维问题时会带来巨大的存储量和计算量。

上述方法中有限元法和有限差分法均采用连续介质力学。离散元法中节理的力学参数需要用户自行定义，但这部分力学参数很难准确估算，节理参数的选取不当会极大影响离散元法的计算结果，无法对工程实际进行准确的模拟和仿真。在不专门针对裂隙发育、扩展的围岩稳定性数值模拟中，基于连续介质的有限元发和有限差分法凭借其较高的计算效率和准确的计算结果在国内外得到更广泛的使用。

大量国内外研究指出岩石中裂隙的产状和发育程度与岩体的力学参数(尤其是杨氏模量)有直接关系，因此在围岩稳定性分析中，特别是软弱、破碎围岩环境下，岩石中裂隙受开挖扰动的动态发育、扩展及导致围岩杨氏模量劣化对围岩失稳机理、控制对策研究、支护方案设计等具有重要影响，是不可忽略的。因此，基于连续介质数值模拟方法，保证计算准确性和高效性的前提下，将围岩中裂隙发展引起的岩体杨氏模量劣化考虑在内，对提高数值模拟过程合理性和结果可靠性均有重要意义。

发明内容

为了更好分析围岩失稳机理，研究围岩控制对策和设计支护方案，本发明克服现有数值模拟技术的不足，提出了一种基于连续介质的工程岩体破裂劣化数值模拟方法，能够动态的实现岩体单元破坏状态识别和杨氏模量随裂隙发育而衰减劣化。

本发明提出的基于连续介质的工程岩体破裂劣化数值模拟方法包括以下步骤：

第一步：建立工程对象数值模型

通过有限元数值模拟软件FLAC3D建立工程对象的数值模型并划分计算网格，在重点研究位置加密网格，得到体积更小、分布更密的岩体单元格和节点，从而提高模拟精度；

第二步：选择工程对象数值模型中岩体的本构模型，并赋予该本构模型相应力学参数

具体包括以下步骤：

第2.1步：将工程对象数值模型中各类岩体的力学模型设置为FLAC3D中内置的应变软化模型，按照应变软化模型所要求的材料参数对各类岩体赋予相应的材料参数，包括：杨氏模量、泊松比、体积力、内聚力、内摩擦角、抗拉强度、残余内聚力和塑性应变量等；

第2.2步：基于地质强度指标体系GSI，使用如下公式建立岩体裂隙发育程度与残余杨氏模量的量化关系：

$E_r=\sqrt{\frac{{\mathrm{\σ}}_m}{100}}\·{10}^{\left(\frac{{\mathrm{GSI}}_t-10}{40}\right)}---\left(1\right)$

式中，E_r为岩体拉伸破坏后裂隙发育作用下的残余杨氏模量，σ_m为岩体抗压强度，由岩石物理力学性质试验和岩体强度估算得到，GSI_t为岩体受拉伸破坏产生的裂隙发育程度；

第2.3步：由工程现场实地观测或理论经验估算获得的岩体裂隙发育程度GSI_t，通过公式②得到劣化后的残余杨氏模量；

第三步：设置工程对象数值模型边界条件

根据工程对象和研究目标设置工程对象数值模型边界条件，包括施加重力、应力限定边界和位移限定边界，模拟工程对象数值模型初始应力场；

第四步：工程对象数值模型开挖并执行工程岩体破裂劣化迭代计算

具体计算包括以下步骤：

第4.1步：将工程对象数值模型需开挖部分定义为空(null)模型；

第4.2步：进行FLAC3D所采用的显式有限差分计算方法进行迭代计算，计算达某一时步数(可由使用者根据情况自由设定)暂停计算；

第4.3步：检查模型最大不平衡力是否低于默认标准值(1e-5)，若低于则模型达到平衡状态，模拟运算结束；若模型尚未达到平衡状态，则遍历工程对象数值模型中全部岩体单元并检测每个单元的破坏状态，识别并标记当前运算状态下模型中破坏状态为拉伸破坏的岩体单元；

第4.4步：对模型中全部标记为拉伸破坏的岩体单元赋予由裂隙发育导致的残余杨氏模量，其数值通过第二步得到；

第4.5步：继续进行显式有限差分迭代计算，计算达第4.2步设定的时步数后重复步骤第4.3至第4.5步，直至工程对象数值模型达到平衡状态。

第五步：模拟运算结束得到数值模拟结果并分析研究

获取模拟结果，根据研究对象进行研究分析，从而为岩土工程、采矿工程围岩失稳机理分析、稳定性控制对策研究、支护方案设计等提供有效可靠的指导依据。

本发明提供的数值模拟方法对比现有技术，具有如下有益效果：

1、采用本发明提供的数值模拟方法进行岩石工程与采矿工程相关领域，特别是裂隙岩体环境下的数值模拟计算，在保留原始数值迭代计算基础上，加入了围岩中应力重新分布过程中裂隙受开挖扰动的动态发展及导致围岩杨氏模量劣化，由现场实地观测或理论经验估算所得数据反演得到劣化后残余杨氏模量，能够使基于连续介质的数值模拟过程更加合理严谨，进而使数值模拟结果更加真实可靠。

2、采用本发明获得的数值模拟结果可以为岩土工程、采矿工程围岩失稳机理分析、稳定性控制对策研究、支护方案设计等提供有效可靠的指导依据。

附图说明

图1本发明方法流程图；

图2实施例矿井三维模型图；

图3巷道支护设计示意图；

图4岩体破坏识别算法结果图；

图5巷道掘进与工作面回采期间围岩变形规律，其中图5a代表巷道掘进期间围岩变形规律，图5b代表工作面回采期间巷道围岩变形规律；

图6巷道掘进与工作面回采期间锚杆、锚索受力规律，其中图6a代表巷道掘进期间锚杆、锚索受力规律，图6b代表工作面回采期间锚杆、锚索受力规律。

具体实施方式

下面结合附图和具体实施例对本发明进行进一步描述。具体实施例和应用中所提供的描述信息仅为示例。这里所描述的实施例的各种延伸和组合对于本领域的技术人员是显而易见的，在不脱离本发明的实质和范围的情况下，本发明定义的一般原则可以应用到其他实施例和应用中。因此，本发明不只限于所示的实施例，本发明涵盖与本文所示原理和特征相一致的最大范围。

本发明的具体实施例采用连续介质有限差分软件FLAC3D，以某地下煤矿回采巷道掘进和工作面回采过程对巷道围岩稳定性的影响为例，具体按照图1所示的步骤实施：

第一步、建立矿井三维数值模型

根据某地下煤矿某工作面实际工程地质情况，通过FLAC3D建立矿井三维数值模型，如图2，模型走向长140m，其中巷道和工作面走向长度为60m，前后各留40m边界，模型倾向宽155m，其中工作面倾向宽90m，模型高100m。

第二步：选择矿井三维数值模型中岩体的本构模型，并赋予该本构模型相应力学参数

设置矿井三维数值模型中岩体本构模型为应变软化模型，具体力学参数见表2。

巷道支护中采用FLAC3D中内置的cable结构单元模拟与现场实际相符的锚杆—锚索支护，模型巷道支护断面图如图3所示，支护体结构与力学参数见表1。

表1

表2

第三步：设置矿井三维数值模型边界条件，模拟初始地应力场

根据该工作面实际工程地质情况设置模型边界条件，在模型顶部施加15MPa的垂直应力模拟实际埋深-650m的上覆岩层载荷，X、Y方向的水平应力分别为垂直应力的0.8和1.2倍，模型四周和底部采用位移限定边界；

采用显式有限差分计算法进行迭代计算，模型的最大不平衡力低于默认标准值(1e-5)后模型达到平衡状态并停止迭代计算，即模型的初始地应力场模拟完成；

第四步：矿井三维数值模型开挖并执行工程岩体破裂劣化迭代计算

为了完整地监测记录回采巷道掘进和回采影响过程中的巷道稳定性状态，设置图2中所示位置为监测断面，记录该断面内的围岩变形量和支护体受力状态。该断面距离工作面开切眼40m，可以较完整地记录该断面在巷道掘进过程中和回采工作面推进采动影响下的围岩响应。回采巷道开挖模拟以掘进5m为一个运算循环，支护紧随巷道开挖，设置支护体并赋予力学参数后运算至模型平衡后，继续下一个掘进—支护—运算循环。巷道掘进完成后，以回采工作面后退式推进5m为一个运算循环，运算至模型平衡后继续下一个循环直到工作面推进至监测断面。

在上述模拟过程中进行迭代运算时，同时执行本发明提出的基于连续介质的工程岩体破裂劣化数值模拟方法，具体包括以下实施步骤：

(1)将矿井三维数值模型需开挖部分定义为空(null)模型。

(2)进行FLAC3D所采用的显式有限差分计算方法进行迭代计算，计算达到100时步后暂停计算。

(3)检查模型最大不平衡力是否低于默认标准值(1e-5)，若低于则模型达到平衡状态，当前模拟运算结束，执行下一步开挖；若模型尚未达到平衡状态，则遍历矿井三维数值模型中全部岩体单元并检测每个单元的破坏状态，识别并标记当前运算状态下模型中破坏状态为拉伸破坏的岩体单元，标记得到的拉伸破坏单元如图4所示。

(4)对矿井三维数值模型中全部标记为拉伸破坏的岩体单元赋予由裂隙发育导致的残余杨氏模量，该实施例中根据现场观测裂隙发育程度，设定煤层GSI_t为10，顶板岩层GSI_t为20，底板岩层GSI_t为25。

(5)继续进行显式有限差分迭代计算，计算达100时步后重复步骤第(3)至第(5)步，待矿井三维数值模型达到平衡状态后继续下一步开挖，直至完成巷道掘进和工作面回采的整个模拟过程。

第五步“提取监测数据并分析模拟结果

提取监测断面位置围岩变形量和支护体受力状态在巷道掘进和回采影响过程中的监测数据，如图5、图6所示。根据以上数据可以对回采巷道完整服务期限(始于巷道开挖，止于回采工作面推进后顶板垮落)内围岩稳定性、支护可靠性等进行评估分析，指导二次补强支护等。通过监测数据可以看出，该巷道属于围岩裂隙发育的大变形巷道，底板和两帮由软弱煤体组成，服务期限内变形剧烈，实际工程中需及时观测巷道变形量，在围岩松动破碎、稳定性差和大变形的位置实施二次补强支护。模拟中顶板锚索在服务期限末段轴向受力超过其拉伸强度，在实际工程中表现为锚索破断从而使支护失效，需要在工作面回采期间超前工作面进行假设单体液压支柱等补强支护，维持回采巷道服务期限内的正常使用。

本发明提出的基于连续介质的裂隙岩体数值模拟方法使连续介质的数值模拟过程更加合理严谨，模拟结果更加真实可靠。可为工程岩体的失稳机理分析、稳定性控制对策研究、支护设计优化等提供了更可靠的仿真模拟方法，特别是对易受拉破裂的工程环境具有普遍的应用价值。

以上所述实施例仅为本发明的其中一个应用案例，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改和等同替换等，均应包含在本发明的保护范围之内。

Claims (1)

1.一种基于连续介质的工程岩体破裂劣化数值模拟方法，其特征在于，包括以下步骤：

第一步：建立工程对象数值模型

通过有限元数值模拟软件FLAC3D建立工程对象的数值模型并划分计算网格，在重点研究位置加密网格，得到体积更小、分布更密的岩体单元格和节点，从而提高模拟精度；

第二步：选择工程对象数值模型中岩体的本构模型，并赋予该本构模型相应力学参数

具体包括以下步骤：

第2.1步：将工程对象数值模型中各类岩体的力学模型设置为FLAC3D中内置的应变软化模型，按照应变软化模型所要求的材料参数对各类岩体赋予相应的材料参数，包括：杨氏模量、泊松比、体积力、内聚力、内摩擦角、抗拉强度、残余内聚力和塑性应变量；

第2.2步：基于地质强度指标体系GSI，使用如下公式建立岩体裂隙发育程度与残余杨氏模量的量化关系：

$E_r=\sqrt{\frac{{\mathrm{\σ}}_m}{100}}\·{10}^{\left(\frac{{\mathrm{GSI}}_t-10}{40}\right)}---\left(1\right)$

式中，E_r为岩体拉伸破坏后裂隙发育作用下的残余杨氏模量，σ_m为岩体抗压强度，由岩石物理力学性质试验和岩体强度估算得到，GSI_t为岩体受拉伸破坏产生的裂隙发育程度；

第2.3步：由工程现场实地观测或理论经验估算获得的岩体裂隙发育程度GSI_t，通过公式②得到劣化后的残余杨氏模量；

第三步：设置工程对象数值模型边界条件

根据工程对象和研究目标设置工程对象数值模型边界条件，包括施加重力、应力限定边界和位移限定边界，模拟工程对象数值模型初始应力场；

第四步：工程对象数值模型开挖并执行工程岩体破裂劣化迭代计算

具体计算包括以下步骤：

第4.1步：将工程对象数值模型需开挖部分定义为空模型；

第4.2步：进行FLAC3D所采用的显式有限差分计算方法进行迭代计算，计算达某一时步数(可由使用者根据情况自由设定)暂停计算；

第4.3步：检查模型最大不平衡力是否低于默认标准值(1e-5)，若低于则模型达到平衡状态，模拟运算结束；若模型尚未达到平衡状态，则遍历工程对象数值模型中全部岩体单元并检测每个单元的破坏状态，识别并标记当前运算状态下模型中破坏状态为拉伸破坏的岩体单元；

第4.4步：对模型中全部标记为拉伸破坏的岩体单元赋予由裂隙发育导致的残余杨氏模量，其数值通过第二步得到；

第4.5步：继续进行显式有限差分迭代计算，计算达第4.2步设定的时步数后重复步骤第4.3至第4.5步，直至工程对象数值模型达到平衡状态；

第五步：模拟运算结束得到数值模拟结果并分析研究

获取模拟结果，根据研究对象进行研究分析，从而为岩土工程、采矿工程围岩失稳机理分析、稳定性控制对策研究和支护方案设计提供有效可靠的指导依据。