FLAC3D中锚杆支护的数值模拟研究综述

杜学领

(贵州理工学院 矿业工程学院，贵州 贵阳 550003)

锚杆作为一种岩土锚固中的重要材料，被广泛应用于边坡、坝体、基坑、隧道、巷道等工程场合，锚杆的使用距今已有100多年的历史[1]。在中国煤矿领域，从早期的砂浆锚杆，到近年来被广泛应用的树脂锚杆和玻璃钢锚杆，煤矿巷道使用锚杆作为主要支护材料的时间已超过60 a。近年来，为适应深部开采、复杂工程地质条件等，可伸长锚杆、中空注浆锚杆、锚杆桁架等被研发并成功应用于多种工程场合[2]。2018年，《煤矿巷道锚杆支护技术规范》颁布，进一步规范了煤矿锚杆的施工、检测及监测要求等[3]。此外，还有D锚杆、Garford锚杆、Yield-Lok锚杆、Roofex锚杆等几十种吸能锚杆[4]。与之相对应的，FLAC3D采用混合离散法进行显式、动态有限差分解算，加之其具有丰富的本构关系、多种结构单元、可进行开放式二次开发等特点，是岩土、矿业领域求解连续介质非线性、大变形、动态问题中使用最为频繁的数值模拟软件之一[5-10]，以FLAC3D为载体的锚杆支护问题也被大量研究。本文以FLAC3D6.0版本的命令流为例，结合国内FLAC3D中锚杆支护的研究案例，在分析模拟技术的基础上，按照FLAC3D中锚杆支护的主要应用场景总结已有研究进展，并进一步分析FLAC3D中锚杆支护数值模拟存在的问题。

1 FLAC3D中锚杆建模概述

1.1 锚杆模型的构建

锚杆锚固系统一般由杆体、锚固剂、托盘等构成，而锚索锚固系统一般由钢绞线、锚固剂、托盘等构成。锚杆、锚索在煤矿巷道使用中往往还需要配合金属网、钢带、螺母等一起使用。尽管锚杆和锚索是两种不同的支护材料，但在FLAC3D中，较多的研究中使用相同的结构单元来模拟锚杆和锚索，以下如无特别说明，所述的“锚杆”在FLAC3D模拟中代指锚杆和锚索。在FLAC3D中，较为常用的锚杆建模方式主要有以下三种：

1) 采用Cable结构单元建模。Cable和Pile都是FLAC3D中内置的结构单元。采用Cable创建锚杆时，可通过以下三种方式之一进行创建：①直线法(by-line)，通过指定锚杆在坐标系中的首尾位置来创建；②节点法(bynodeids)，通过指定锚杆首尾两端2个端点的结点ID来创建(注：一般将nodes译为结点，将gridpoint译为节点，二者以示区分)；③射线法(by-ray)，通过指定锚杆起点、方向(点)、锚杆长度三个参数来创建[10]。方法①适用于创建水平或平行于坐标轴的锚杆，方法③适用于创建与坐标轴斜交且易于确定锚杆起点和所经过点的锚杆。此外，FLAC3D中的Beam、Pile这两种结构单元也可采用上述三种方法进行建模。

采用Cable创建的锚杆较多应用于巷道掘进支护设计、边坡锚固效果分析等，该方法在FLAC3D模拟锚杆支护的研究中较为多见，以至于一些研究中直接将Cable翻译为锚杆，实际上Cable更倾向于译为锚索。Cable结构单元可赋予的属性参数有12个，其中较为常用的参数主要有锚杆杆体密度(density)、杨氏模量(young)、横截面积(cross-section-area)、抗压强度(yield-compression)、抗拉强度(yield-tension)等，锚固剂或注浆体单位长度的黏聚力(grout-cohesion)、摩擦角(grout-friction)、单位长度的刚度(grout-stiffness)、外露周长(grout-perimeter)等。Cable结构单元可受拉和受压，但不能抵抗弯矩。由于网格划分的疏密不同，对于长度为2 m左右的锚杆，一些拉拔实验模拟研究中单一Cable结构单元的长度设置为0.05～0.2 m[11-13]，即一根完整的锚杆被分为10～40个Cable结构单元。而在工程类模拟中单一Cable结构单元的长度则更大。一般，对应于FLAC3D中的一个区域(zone)，存在Cable结构单元的一个结点(node)即可，过多的Cable结构单元会影响到整个锚杆的失效表现[10]。但笔者的初步研究表明，FLAC3D帮助手册中的这条建议并不可靠，Cable数量对锚杆性能的影响已另行撰文讨论。Cable结构单元的锚固参数默认是全长锚固，对于非全长锚固的锚杆(如端头锚固)，一般对此类锚杆施加预紧力，由于该类锚杆的大量应用，下文将单独介绍其建模方法。

2) 采用Pile结构单元建模。Pile结构单元是FLAC3D手册中明确提出可用于模拟岩石锚杆(rock bolt)加固的构件，其被翻译为“桩”。除具备梁(Beam)的结构特性外，Pile结构单元还可与网格发生法向和切向的摩擦作用，具有抗弯特性，因而Pile又被视为是Beam和Cable的综合，但Pile结构单元不可直接施加预紧力[10]。目前国内煤矿巷道支护方面使用Pile结构单元构建锚杆的研究要比使用Cable的少，可能是与Pile的属性参数相对复杂有关。每个Pile结构单元具备基本参数20个，即6个材料参数、4个截面几何特性参数、10个耦合特性参数。当使用锚杆逻辑时，还有额外的7个参数需要考量。与Cable仅有12个属性参数而言，使用Pile结构单元模拟锚杆功能有27个参数，尽管这些参数并非都需要一一设置，但还是会增加Pile模拟锚杆的复杂性和难度。

Pile在基坑、桩锚等应用相对较多，与其本身“桩”的设定有一定关系。在矿山支护相关的应用中，也有部分研究用其来模拟锚杆的力学性能(如拉拔试验、剪切破断特性)及巷道支护等[14-16]。

Cable和Pile结构单元的坐标系统均采用局部坐标系，X轴正向为建立结构单元时所指定的起始结点1到终点结点2所指向的方向，所不同的是，Cable为2个自由度，而Pile与Beam相同，具有12个自由度。二者的作用机理如图1所示。

图1 FLAC3D中Cable和Pile结构单元的作用机理

3) 采用实体结构单元建模。更少的研究采用FLAC3D中的实体单元(如Cylinder、brick等)配合接触面(Interface)来实现锚杆的模拟[17-18]。该类型的建模一般用于锚杆性能表征的研究中，在支护工程中采用此方法构建锚杆的案例鲜有报道。

上述建模方法的典型研究案例如图2所示，其中Cable和Pile在外观上是相同的，其差异在于不同结构单元的属性参数，而实体单元建模的锚杆往往与本构模型、材料参数、接触面特性等有关。

图2 不同建模方法的典型案例

1.2 预应力锚杆的实现

预应力被认为是锚杆支护中的重要参数之一，并有研究认为提高预应力可显著减小巷道变形[19]。因此，预应力锚杆的模拟也是FLAC3D中的研究热点之一。FLAC3D中的Cable和Pile结构单元经常被用于模拟锚杆，已公开文献显示，FLAC3D中预应力场的形成方法主要有以下几种：

1) 单向拉伸法。对于Cable结构单元建立的锚杆，陈育民等[8]基于FLAC手册中预紧力的施加方法，提出锚杆预紧力的实现可通过以下三种途径：①删除锚杆自由段端头的原有Link，并在该位置新建与Zone刚性连接的Link；②将锚杆自由段端头部位的Cable结构单元设置为锚固剂参数极大；③删除锚杆自由段端头的原有Link后，新建Liner结构单元和新Link，并将新Link刚性连接到Liner上。实施以上三种方法后，均需对锚杆自由段部分施加单向拉伸载荷，以获得预应力场。同时认为，相同拉拔载荷下以上三种方法获得的预应力场效果是基本相同的。但在其提出的方法中，并没有指出何时取消对锚杆自由段的拉伸载荷。而根据软件开发公司的指导手册[10]，预紧力施加后会影响到原有的应力场状态，需要获得新的平衡后取消对预紧力的施加。从工程实践而言，预紧力也不可能无限地作用在锚杆和围岩之上的，也需要在施加预紧力后的一定运算步中撤除预紧力。但对于工程模拟而言，不可能每次安装锚杆之后都运算新的平衡，这样不仅耗时，而且与工程实际中开挖、支护的连续过程也不相符。从这一点而言，一些使用Cable结构单元并施加预紧力的数值模拟中，如果一直保持预紧力存在，可能在方法上是错误的。

2) 对向拉伸法。文献[14,20]针对FLAC3D中Pile结构单元本身不具备施加预紧力的问题，提出利用Pile结构单元施加预紧力的主要思路如下：①在锚杆位置建立上下两个锚杆分段，上分段进行端部锚固，下分段的端部刚性连接到Liner结构单元；②在两分段施加一对与预紧力大小相等、方向相反的结点力对锚杆上下分段进行拉拔，并固定网格；③模型运算平衡后，建立两分段之间的中间锚杆分段，移除结点力并继续计算新的平衡；④释放网格，获得锚杆作用后在围岩中形成的预应力场。利用该方法可获得类似“两压一拉”的应力集中区，但该方法中有几个细节并未披露，如上分段的锚固长度如何确定、中间分段的长度如何选取、网格固定的范围为何、添加中间分段对围岩应力分布有何影响、施加结点力的范围多大等。此外，这种方法的合理性还有待于考证。由于Cable结构单元中也包含有node，因此对向拉伸法也可用于Cable结构单元。

图3给出了单向拉伸法中方法2设置自由段端头锚固剂参数极大获得的预应力场分布，对向拉伸法形成的预应力场可参考文献[14,20]。由图3可知，通过拉拔作用，锚杆受力带动围岩运动，因此锚杆周围围岩中的应力发生重新分布。上述两种方法的共同点在于均通过拉拔锚杆的形式获得预应力场。其中单向拉伸法拉伸方向与锚杆锚固深度方向相反，实现拉伸过程的主要命令是“structure cable apply tension value”；而对向拉伸法拉拔段有两个，实现拉伸过程的主要命令是“structure node apply force”。但存在的问题是：如图3所示的预应力分布是否符合真实的场景；这种拉拔作用是否与结构单元数量、拉拔位置等是否有关；预应力场的应力量级对于工程而言是否具备数量优势；采用结点拉拔与整段拉拔又有哪些差异。诸如此类的问题还非常多，以至于虽然上述两种方法在一定程度上实现了锚杆的预紧力施加，但却未能从原理层面上阐明预应力场形成的机理、预应力场的作用范围和作用时间等问题，目前看到的一些预应力研究，缺乏基础性的数据呈现，也缺乏物理实验或工程实践一一对应的支撑。

图3 单向拉伸法形成的预应力场

需要指出的是，在一些研究中将预应力和预紧力混为一谈。在数值模拟中，预紧力一般是施加给锚杆的拉拔力，这个力一般作用于锚杆轴向，而预应力是因施加预紧力后改变了围岩内的应力分布而形成的应力重新分布形式，即预应力场。二者并不相同，而且施加的预紧力也并不能完全转化为预应力。

1.3 其他支护相关要素的构建

除锚杆外，其他支护相关要素的构建方法主要有：①地层的构建一般以实体单元为主，其中应用最多的为Brick。使用外部导入的模型时，也可能使用Wedge网格。②巷道的形状以半圆拱形、矩形、圆形、半圆形为主，其中半圆拱形、矩形应用较多，一般不考虑地层倾角的变化。通过Brick及Radial brick、Radial tunnel、Radial cylinder等网格建模后，使用null并配合range命令实现巷道的开挖及矿产的开采。③锚杆支护辅助构件的建模。相当数量的研究并没有提及金属钢、托盘等在FLAC3D中如何处理，在明确辅助构件的建模中，主要采用以下几种方法：一是采用Beam构建钢带、钢筋梯子梁、混凝土喷层等[11，21]；二是采用Liner构建蝶形托盘和W钢带[14-15]；三是改变围岩参数来实现锚固后围岩性质改变的等效[22]。

2 FLAC3D中锚杆支护的研究进展

由于FLAC3D一般作为验证性、辅助性手段使用，往往是对实验、实践内容进行合理性的证明，尽管FLAC3D中的Cable、Pile参数众多，但众多的研究往往从工程或应用视角出发证实相关设计的合理性，而较少从FLAC3D模拟软件本身出发来探究锚杆模拟参数对支护效果的影响。在FLAC3D手册中，提供了一些与锚杆支护相关的案例，其中采用Cable建模的案例主要有重力作用下梁的加固、路基的加固、全长锚固的拉拔试验、直墙边坡开挖支护、沉箱挡土墙支护、浅埋隧道支护等；采用Pile建模的案例主要有轴向加载桩、横向加载桩、全长锚固的拉拔及剪切试验、塑性硬化模型(Plastic-Hardening model)下的隧道开挖支护等[10]。现有的一些研究主要在这些案例的基础上进行拓展，如研究锚杆的拉拔或剪切特性、不同工程中的支护效果、支护参数优化等。以下结合FLAC3D锚杆支护研究中的不同侧重点展开论述。

2.1 锚杆及其辅助锚固材料的特性

树脂锚杆破坏主要发生在杆尾螺纹与垫圈、交界面(如外露段与钻孔交界面、锚固段与自由段交界面、杆体所穿越的围岩裂隙面)。其中杆尾螺纹段失效是煤矿巷道支护中最为常见的破坏形式，拉弯应力复合作用是导致该破坏的主要原因。此外，材料缺陷、加工工艺、元件的相互配合、施工技术等也是造成锚杆失效的重要原因[23]。锚杆材质强度的提高有助于提高锚杆的极限抗剪强度。围岩强度提高有助于提高锚杆的屈服载荷，但对极限抗剪强度影响较小。锚固剂及锚固方式对锚杆支护效果具有较大影响，一般锚固剂强度越高，锚固系统的初期抗剪刚度越高，但对锚固系统的后期影响则相对较小。钻孔与锚杆的间隙不宜超过10 mm，以6～8 mm为宜[24]。左旋无纵筋螺纹钢锚杆的拉拔力要高于右旋全螺纹钢锚杆，因此当前我国煤矿以左旋无纵筋锚杆应用为主[25]。对于左旋无纵筋螺纹钢锚杆而言，横肋高度增加有助于提高锚杆的拉拔力；横肋间距增大有助于提高残余锚固力，但同时锚杆安装过程的最大扭矩及安装推力会减小、锚杆拉拔力会先增大后减小；托板承载力与其拱高有关，应使拱高与支护设计相匹配；螺母承载性能应不低于锚杆杆体，球形垫圈强度应高于托板；树脂锚固剂在锚固范围、界面接触、安装位置等方面对支护效果产生影响，树脂锚固剂存在最佳工作温度，超过或低于最佳温度时其锚固力都会下降；锚杆支护系统在淋水条件下会随着淋水量的增加而性能降低，甚至是严重恶化，锚杆的偏心安装会造成锚杆受力的非对称性，并造成局部锚固剂内应力集中、塑性区扩大[23-25]。W形钢带的性能与其厚度、宽度、钢带孔几何特征、加工工艺等有关，长圆孔钢带的承载性能要显著低于无孔及圆孔钢带，钢筋托梁的承载力与其直径、力学特性等密切相关[26]。

拉拔及剪切试验一直是评价锚固系统性能的重要方法，锚杆的试验模拟在FLAC3D中也被大量研究，典型的研究模型如图4所示。进行拉拔试验时，一般对锚杆端部施加位移控制条件或力的控制条件；进行剪切试验，一般通过对围岩施加力或围岩控制条件来观察剪切模拟中锚杆的响应。

图4 拉拔及剪切模拟模型

锚杆的拉伸过程一般包括弹性、屈服、强化、颈缩四个阶段，锚杆加工的热处理工艺可使屈服阶段缩减甚至无屈服阶段，其相应的拉断延伸率也会降低[24]。江文武等[27]通过拉拔数值模拟的分析认为，拉拔过程中，自由段处轴力最大，并沿杆体向深部逐渐衰减，拉拔前端是受力的主承载区域，增强锚固剂的内摩擦角、粘结力和有效围压可增强锚固效果，而锚固剂所受的剪应力分布并不均匀，传统的采用均匀化处理的设计方式而得到的锚固段长度与实际工程并不符合。郑卫锋等[28]采用Cable结构单元研究单根锚杆的拉拔模拟，认为增加围压会提高锚杆的屈服强度，并据此提出可通过二次高压注浆的方式来提高锚杆承载力。康红普等[25]采用实体单元和Interface建立了锚杆拉拔的数值模拟模型，并对杆尾施加80 kN的拉拔载荷，研究认为圆钢锚杆和螺纹钢锚杆周围锚固剂的受力、变形有较大的区别，螺纹锚杆的峰值剪应力分布在横肋，且存在应力集中，横肋更好地发挥了锚固剂力的传递作用。但当钻孔孔径和锚固剂厚度增大后，轴向力的传递效率会降低，过小的孔径和锚固剂厚度则易出现滑移变形。锚杆安装偏心会使得应力分布更为复杂及应力峰值向深部转移。同时，其认为围岩介质属性会影响到锚固剂的破坏形式，坚硬岩石中以拉伸破坏为主，而软岩中复合了拉伸和剪切破坏。陈瑶[29]建立了类似的锚杆拉拔分析模型，认为剪应力沿杆体先增大后减小，锥形分布的塑性区分布在锚固段前端，锚固段剪应力分布存在一个范围，使得增大锚固长度并不一定能增加锚固效果；增加锚杆直径会降低峰值剪应力，而锚固剂厚度对锚杆锚固力学性能则影响有限；相同拉拔载荷下，围岩弹性模量越大，形成的剪应力峰值会越大。尽管高德军等[30]认同增大锚固剂弹性模量可以改善锚固效果的观点，但同时其指出改变锚固剂长度的影响并不理想，主要原因在于锚杆的主承载区域位于锚杆的前端。对拉拔应力场的进一步研究认为，拉拔作用下，单根锚杆的应力分布形态为“两峰夹一谷”，且中部应力较低，2根锚杆的支护应力场为单根应力场的叠加，变为“四峰夹一盆地”的特征[31]。李桂林等[32]利用粘滞边界研究了动载扰动下Cable在输入动载为三角载荷的冲击波影响，结果表明，距离动载源越近，围压越大，受扰动影响下围岩首先承受爆炸应力波较大的正应力影响，再受到较小的拉应力作用。对锚杆而言，锚端、中部受压力作用，而锚头则受拉力作用，从锚端到锚头，体现出压力降低、拉力升高的特点。言志信等[33]在不考虑界面滑移的情况下采用弹性本构模型对动载导波的特性进行了研究，结果表明，导波激发频率大于50 kHz时应力波在锚杆中随频率增大而显著增加，但均低于纵波波速，纵波在锚杆中衰减很大，而stonely波沿锚杆轴向则不会衰减。

剪切试验表明，剪切破坏的变形易发生在围岩节理面，变形范围是杆体直径的2～4倍，高强度岩体中易发生拉剪破坏，软弱岩体中易发生拉弯破坏。锚杆橫肋的存在可提高锚杆与锚固剂之间的摩擦效应，预紧力对不同类型的锚杆影响有所差异，一般提高预紧力可在一定范围内提高屈服强度、极限强度和锚杆的整体抗剪刚度，但如果预紧力过高，可能导致锚杆过早进入屈服阶段，有可能造成极限剪切载荷降低[24]。煤矿巷道支护中推荐的预紧力为锚杆屈服强度的30%～60%[3]。基于Pile结构单元建模的研究表明：Pile中无法实现锚杆破断过程中的颈缩现象，剪切位移也低于实验测定的数值；锚杆直径存在临界效应，适当增大锚杆直径有助于提高抗剪强度，但过大的直径则对轴向受力、抗剪强度提升影响较小[24]。高杰[34]在施加围压的条件下研究了锚杆直径为16～24 mm的支护效果，但其认为增加锚杆直径可提高控制围岩塑性区发育和减小围岩变形的能力。

借鉴陈昌富等[35,10]对锚固界面载荷传递的总结，锚杆性能表征的本构关系一般可分为如图5所示的5种类型：①理想弹塑性模型，该模型是FLAC3D中Cable、Pile结构单元的默认本构关系，当载荷达到峰值载荷时，载荷保持在峰值载荷不变，并继续发挥承载作用。需要指出的是，该模型虽然为默认的本构关系，但与真实条件下锚杆达到峰值载荷后破断并失去全部或部分承载能力的现实是不相符合的，因此并不适用于常规锚固剂锚固的情况。②峰后破断型，该模型可直接应用于Pile结构单元中，通过设置峰后的应变软化情况，可直接实现锚杆在峰后达到相应的应变时应力直接跌落到某一水平，实现类似锚杆破断的效果(对于破断的情况，一般将峰后承载力设为0或较小值)。对于Cable结构单元，则需要借助fish语言来实现类似的破断效果。对可破断锚杆的模拟中，一般采用此种本构关系[11-16]。③残余强度型(应变软化型)，该模型实际上是峰后破断型的一种变体，最典型的为三折线峰后软化模型，第一折线为峰前的弹性变形阶段，第二折线为峰后的应力降低阶段，第三折线为应力跌落后维持在某一较低残余强度的阶段。实际使用中，还可根据需要对峰后的软化进一步细分为多个软化阶段，实现分段软化和渐进性破坏。④应变硬化型，与残余强度型类似，在锚杆达到峰值载荷时，实现承载能力的分段提升。该类型的本构关系可用于模拟锚杆杆体的变形破坏过程，应变硬化特征与钢材破坏过程中出现的弹性阶段、屈服阶段、强化阶段等相类似。与残余强度型综合使用，还可实现强化阶段后锚杆的突然断裂。⑤理想弹性型。与锚固剂、围岩相比，锚杆杆体材料的弹性模量一般非常大，因此，一些采用实体单元建模的研究中将锚杆杆体视为理想弹性材料，其承载力随着变形的增加而无限增加。该模型的缺点在于，此时锚杆永远不会破断，弹性材料也不会达到屈服状态。为了改变这种情况，可以采用fish语言设定一定判定条件，达到判定标准时更改锚杆的本构关系和力学性能表现。

上述研究中，在相应的研究模型中往往聚焦于锚杆的某一个方面，如较多的研究对锚杆的载荷-位移曲线给予了充分关注，却没有把整个锚固系统的特性阐述清楚，如塑性区如何演化、塑性区演化过程与载荷位移之间具有怎样的关系等。由于锚固系统是由锚杆、锚固剂、围岩构成的综合体系，仅仅关注锚杆的力学响应，有可能对锚固系统的失效表现判断失误。类似的研究也出现在锚杆的物理实验当中，对拉拔载荷给予了充分关注，对锚固剂破碎形态、破碎范围、破碎粒径、锚杆断口位置及大小等其他细节的关注则较少。另一个需要指出的问题是，很多研究是基于FLAC3D中默认的全长锚固进行的，但并没有对端头锚固或加长锚固的情况进行对比分析。或即便建立了相应的分析模型，分析过程中更侧重结果数据的对比呈现，而缺少对破坏失效演化过程的分析，以至于很多研究呈现的只是数据结果的堆砌，而缺少过程性研究。

2.2 预应力锚杆的研究进展

预应力锚杆支护技术在支护工程中已被广泛采用。在技术规范方面，我国的预应力锚杆与欧美日等国家的预应力锚杆规范还存在一定差异，以至于尽管我国工程领域预应力锚杆的应用范围、应用质量已走在世界前列，但在相应的理论研究中却依然与发达国家存在一定差距[36]。以煤矿锚杆的预应力而言，在实验方面：林健等[37-38]通过对单根锚杆的室内测试表明，单根锚杆在围岩中所形成的空间应力场形似“石榴”、“酒精灯”，并在其长度范围内形成“两压一拉”的应力集中区，即使改变预紧力，预应力场的宏观形态是相似的；周逸群[39]进一步将该研究的锚杆根数扩大到2根、4根，研究不同围岩复合条件下的预应力场分布，相似试验结果表明复合围岩条件下的预应力场分布形式与林健等的研究结论是类似的，并通过试验证明组合构件的影响范围不能达到锚固深部，而仅作用于围岩表层附近；黄雪峰等[40]对桩锚支护结构的受力监测表明，锚杆未施加预应力时，在开挖初期其自由段轴力几乎为零，锚固段则先增后减，而在施加预应力后，自由段的轴力最大且分布均匀、锚固段内的轴力则逐渐衰减；张剑[41]对井下不同锚杆的预紧扭矩及预应力关系展开测试，表明二者呈正相关，且锚杆的偏斜度、安装角、围岩强度等也会影响预应力场的大小。理论方面：周辉等[42-43]认为预应力锚杆的作用体现在杆体抗拉的“轴压”效应及抗剪切的“销钉”效应，并可削弱裂隙尖端应力集中；何思明等[44]根据弹性理论，获得不考虑与考虑界面滑脱下的锚杆剪应力分布特性，不考虑滑脱时，剪应力单调下降，而考虑滑脱时，滑脱后的剪应力峰值点会向锚固深部转移；吴永等[45]在此基础上指出，通常峰值剪应力主要来源于外载荷，锚杆长度对峰值剪应力的影响不大；罗卫华等[46]利用半无限体柱状孔洞内受均布压力作用下的弹性解推导了位移和剪应力分布形式，并通过参数分析认为拉拔载荷具有增大剪应力、位移的效果，适当提高拉拔载荷有助于提高锚固效果；王洪涛等[47]认为，施加预紧力应留设一定自由段长度以利于有效锚固承载结构的形成，如果锚杆间距较大，可通过提高预紧力、合理减小锚固长度来增强预紧力的效果。

与之相对应的，数值模拟技术也被广泛应用于预应力锚杆的评价中。

早期的研究当中，曾有研究人员将预紧力产生的预应力场按照不同阶段进一步细分为初锚力、工作锚固力、残余锚固力。近年来的研究中，对预应力场并未做如此严格的区分。冯光明等[48]采用允许大变形的方式模拟端锚和全锚的效果，并施加不同吨位的预紧力，其获得的预应力场表现为应力等值线靠近锚杆两个端部稠密、中部稀疏，在托板中间均匀分布、两侧呈“X”形。与端锚相比，全锚能提供更大的剪应力，但端锚支护中提高初锚力所获得的预应力场提升效果要比全锚时更明显，该值一般在0～0.65 MPa之间变化。王同旭等[49]采用Cable、Interface模拟锚杆和节理，研究认为预紧力的主要功能在于改善围岩应力状态，预紧力同样存在临界效应，施加的力如果超过临界预紧力，则在改善围岩稳定性方面效果并不明显。金爱兵等[50]通过对锚杆张拉顺序的研究表明，围压会改变围岩中的应力分布，但并不会改变锚杆杆体的受力状态，相邻锚杆的预紧会造成首拉锚杆的轴向力损失，锚杆间距越小，首拉锚杆轴向力损失越多。按前后顺序张拉要优于间隔张拉和从中间向两边张拉。高峰[51]采用FLAC2D并通过改变围岩属性的方式来实现不同预应力水平的影响，但这种模拟从方法上而言是错误的。对于给定的围岩条件，预应力只是改变了围岩的受力状态，并没有改变围岩的物性。围岩的应力、位移变化是由开挖、支护等共同决定的，预应力是开挖后产生的，而不是开挖前产生的。所以，要想准确模拟预应力的效果，应首先在模拟流程上还原真实的施工条件。李厚恩[52]以基坑工程为基础对比分析了预应力锚杆对土钉支护的影响，其利用Brick建立土体、利用Cable建立土钉，并对非全长锚固的土钉施加预紧力，结果表明施加预紧力具有控制坡面变形的效果，但却并不能完全控制临空面的变形，预应力能对局部土体起到“紧箍”作用而提高作用位置附近土体的强度和稳定性，但预应力作用的范围是有限的。徐松山[53]进行了类似的研究，进一步的，徐松山还对数值模型进行了分级加载模拟，但其并未给出加载的详细过程。其获得的结论与李厚恩获得的结论类似，即：预应力锚杆可在有限空间内起到控制变形的作用。涂兵雄等[54]以大连胜利广场深基坑工程为背景，研究预应力锚杆柔性支护法的力学行为。其利用Cable和Beam分别建立预应力锚杆和喷射混凝土层面，并按照“开挖-加锚杆-喷射混凝土-锚杆张拉”的过程进行模拟。其认为预应力对变形的控制作用与预应力数值大小有关，预应力较小时，增大预应力可减小变形，但预应力较大时，持续增加预应力并不能获得持续降低变形的效果。预应力锚杆柔性支护法中预应力锚杆的作用体现在控制变形和塑性区发育、提高稳定性(对抗滑移)等方面。郑军[55]建立了未考虑原岩应力的巷道支护数值模拟模型，并认为锚杆预应力在一定范围内呈衰减性扩散，预应力分布在锚杆尾部的一端而非两端，预紧力增大后，围岩中的压应力也会随之增大。在群锚当中，锚杆间排距较大时，会造成所产生的压应力与单根锚杆形成的压应力场类似。当锚杆间距小到一定距离时，可形成有效的围压应力圈，即在围岩中形成稳定的加固拱。预应力的作用体现在恢复围岩三向应力状态、提高残余强度、控制变形等方面，预应力锚杆作用的发挥与锚杆的高强度有关。丁书学等[56]利用线性叠加原理研究裂隙影响下的预应力场分布，认为裂隙影响因素中由强到弱依次为裂隙的长度、倾角和位置。何富连等[57]认为，提高锚杆刚度或直径，可扩大其应力分布范围，并可降低应力峰值。李超等[58]在复合岩层模型中研究预应力的传递过程，认为岩层组合形式对预应力场的分布具有重要影响，在“软-硬”界面会出现阻隔现象而影响预应力的扩散。

在预应力的数值模拟研究中，存在最明显的问题是自由段的预应力场问题。如图6为以往物理实验获得的预应力场分布，真实的场景当中，锚杆自由段是在钻孔中临空存在的，自由段除两端外不与其他实体产生接触，因此绝大部分的自由段中是不应该存在预应力场的(杆体中存在轴向作用力，但临空的杆体周围不存在应力场)，至少在自由段杆体外边缘到钻孔直径范围以内，是不可能形成应力场的。这一问题在数值模拟研究中更为普遍，在模拟当中，由于Cable、Pile发挥作用需要依赖与实体单元建立的link和处于实体单元中的node来发挥作用，当Cable、Pile结构单元临空时，就无法发挥出其作用。但显然，目前大部分研究中所呈现的预应力场都未考虑自由段临空的情况，以至于其呈现的模拟结果中出现自由段周围也存在预应力场。关于锚杆自由段临空的问题，笔者已另行撰文讨论。另一方面，预应力场的形成与预紧力作用密不可分，但对于预紧力的作用范围、作用时间效应、作用机制等并没有阐述清楚，而仅仅是呈现了预紧力的作用结果。特别的，实际上预紧力所形成的预应力场在量级上而言是相对较小的，在地应力、采动应力作用下，预应力的作用机理究竟在于何处，预应力在整个锚固周期如何发挥作用，目前在此方面依然没有达成普遍共识。

笔者采用线上、线下相结合的混合式教学模式：首先开始上课时，我们对线上内容进行测验或课前问题答案展示，这样做可以监督学生线上自觉进行微课视频的学习，增强了同学们的成就感和满足感，激发了同学们学习编程的兴趣和自信心，同时也给学生提供了锻炼自我的机会。

图6 以往物理实验研究获得的预应力场(文献[37])

2.3 煤矿巷道支护系统视角下的锚杆支护模拟

除支护材料本身的性能影响外，支护材料之间的相互配合也会对支护效果产生重要影响，包括锚杆的间距、排距、倾角、布置形式、组合方式等。单仁亮等[59]结合马兰矿的条件综合使用应变软化和摩尔库伦模型对不同支护方案进行模拟，并通过围岩变形量、塑性区发育等对支护效果展开评价，其采用一掘一锚的方式对掘进长度为11 m的工程进行模拟，研究认为上下帮锚杆分别与水平面上下斜交15°时安装锚杆可分别改善顶板和底板的支护效果，锚杆间排距的影响要比锚杆长度、直径的影响更大，同时，锚索排距对巷道变形也有重要影响。肖福坤等[60]模拟研究认为，增加底角锚杆有助于控制巷道的底鼓，但其并未指出底角锚杆的安装角度、锚固长度情况等。郭建明等[61]则认为，锚杆长度控制围岩变形的效果不如锚杆间排距。张益东[62]通过相似模拟实验及数值模拟研究证明，增大支护密度或锚杆直径有助于控制围岩的变形，但增加锚杆长度仅在一定范围内有增强效应，超出此长度范围则锚固增强作用不会继续增加。胡志云[63]对比分析了单一锚杆支护、锚杆锚索联合支护、全断面短锚索支护的差异，认为全断面短锚索的支护性能要优异其他两种支护形式。但应注意到，在其研究中对前两种方案施加了50 kN的预紧力，而第三种方案当中则施加了120 kN的预紧力。对于这种差异，究竟是预紧力的作用，还是支护形式的作用，并不能形成有效的对比。耿亚东等[64]利用Hoek-Brown准则将岩石参数转换成岩体参数，认为锚杆超前安设的角度为15～25°。常聚才等[65]结合望峰岗矿的条件，采用应变软化模型和正交试验法分析不同因素的影响程度，按其重要程度，其认为依次为锚杆间排距、直径、长度和预紧力。陈登红[66]利用正交试验探究采煤工作面长度、采高、埋深、支护构件长度的影响，并建立长度达上百米的分析模型，认为应变软化模型与深部围岩的力学变化特性相符，深部围岩在掘进初期的应变软化特性对其自身的稳定性具有重要影响，巷道埋深是影响围岩稳定性最关键的因素。李小裕等[67]采用Cable结构单元构建锚杆，并通过正交试验证明对控制顶板变形而言，支护密度、预紧力、布置形式、锚杆长度的重要性依次减弱。唐红等[68]根据锚杆锚索协同支护理论在FLAC3D中研究锚杆锚索的协同支护效果，具体做法为锚杆支护运行300时步后再添加锚索。其研究认为，锚索预紧力、长度、锚固长度对巷道变形的影响力依次递减，20 kN锚杆预紧力与100 kN锚索预紧力的组合可发挥很好的协同支护效果。但在其研究中，锚固长度系数均不超过0.4，也就无法体现加长锚固、全长锚固的效果。白启树[69]通过数值模拟认为，预应力锚杆的作用体现在形成压应力区上，锚杆的间排距设计应能形成有效压应力区叠加才有利于锚固作用的发挥。严红[70]提出了支护强度因子的概念，该因子综合反映锚杆材质、长度、直径、密度、角度、预紧力、锚固剂、钢筋网、钢带等，但其并没有给出该因子的具体计算方法以及每种因素在该因子中所占的比重。

此外，锚杆支护在煤矿领域还体现在特殊地质条件、特殊开采条件等方面。张斌川等[71]结合新汶潘西矿的条件，采用应变软化模型模拟煤层、摩尔库伦模型模拟其他岩层，研究认为强力支护能更好控制围岩变形。孙海英[72]利用Pile、Liner分别建立锚杆和钢带，分析了长度为100 m的工作面、5 m煤柱沿空留巷工艺下的锚杆受力，表明锚杆轴力在巷道中间大、边角小，剪力分布则与轴力相反。乔伟[73]对沿空留巷的支护情况进行了数值模拟，认为锚杆轴力会随着工作面的推进而呈现不同程度的增大。冯俊伟等[74]则认为，大倾角煤层中的巷道高帮上部变形较大，应加强支护，增强锚杆长度可提高控制围岩变形的能力。张伟光等[75]建立了近距离煤层的巷道支护分析模型，研究结果则表明，锚杆长度对于巷帮水平位移和顶底板移近量的影响要比锚杆间排距更重要，锚索则表现出临界长度的效应。肖同强等[76]研究了不同侧压系数下的支护效果，研究认为构造应力增大后，会增大肩角锚杆失效几率，且肩角锚杆在煤岩交界处轴力最大，安装角越大，轴力越大，掘进初期应采用较小的安装角。刘建伟等[77]采用Cvisc蠕变模型和Cable结构单元研究深埋巷道中锚杆的蠕变响应，认为锚杆应力变化要经历瞬时增阻、静力平衡、衰减增阻、等速增阻四个阶段。此外，软岩巷道、冲击地压影响、多煤层开采、孤岛工作面、复采、断层等与支护有关的特殊条件也被学者们在FLAC3D中进行一定程度的研究。

对大量煤矿锚杆支护的模拟研究进行调研分析可以发现，当前大量的巷道支护研究选取部分巷道作为分析模型，仅仅研究局部巷道开挖后的支护问题，而没有从整个生产服务周期的视角去探究锚杆在全生产周期中的表现。另一方面，即便是局部巷道，也几乎没有从这些模型中看到类似于周期来压的效果，以至于对为什么选择这么长的巷道进行研究也无法言明。此外，面对复杂的地质条件和开采条件，关于锚杆长度、锚固长度、安装角度等问题还存在一些相矛盾的观点，是模拟的疏漏，还是由地质异常所引起的差异，还需要进一步深入研究。特别是对于大部分研究而言，由于仅建立了部分分析模型，在边界条件类似、锚杆属性相近的条件下，应该获得一些近似的结论。对于这种差异，还有必要进行对比分析。一些数值研究中的过程、观点也存在疏漏的嫌疑，如真实的巷道形状为4.2 m×3 m，却能在数值模拟中建立5.2 m×3.5 m的模型[78]；掘进期间的最大变形仅为29.4 mm，回采期间最大仅为30.9 mm[79]，这在量级上应是存在错误的；还有研究认为无采动影响可不进行支护[80]，这显然也是不符合安全生产要求的。由于在数值模拟中进行了适当的简化，出现与现实不完全一致的结论是非常正常的。但数值模拟在趋势性的推演、观测变量的量级方面应该是与现实较为接近时，才能获得较为合理的结论。

2.4 岩土工程中的锚杆支护模拟

FLAC3D锚杆支护模拟在岩土工程中主要应用于边坡加固、隧道支护、基坑开挖支护等方面，在锚杆参数的优化方面，如锚杆直径、锚杆长度、预应力的影响等方面，与前述的很多结论是近似的，因此不再赘述。此外，岩土工程中的模拟，由于支护形式、支护要求的不同，其研究的侧重点也有所不同。以下结合岩土工程中的一些代表性研究，予以说明。

吴礼舟等[81]利用Cable、Beam等建立锚杆框架梁的数值分析模型，并将锚杆与框架梁进行刚性联接，通过fish实现地下水位的模拟，结果表明对于二级开挖边坡，上下段应分开设计，锚杆间距、布设角对变形影响较大，而锚杆长度超过一定值时则影响有限。倪红梅[82]认为，对于基坑工程，锚杆长度应控制在1.5～2倍开挖深度，安装倾角不宜超过30°以免最大水平位移的影响加重，设置深度以支护结构顶部1 m左右为宜。张成龙[83]利用Liner建立衬砌结构来模拟预应力锚杆柔性支护的效果，施加初始预紧力200 kN，并进行分步开挖、求解、加锚杆，其认为，锚杆轴力在自由段与锚固段的结合处最大。胡建敏等[84]将围岩的本构关系设定为Drucker-Prager模型研究锚杆间排距、长度等的影响，研究认为间排距过大时锚杆的作用就不明显。张钦喜等[85]利用Cable、Beam分别模拟锚杆和护坡桩，在不考虑锚杆重力的情况下分析认为锚杆自由段是承受拉力的主要部分。张华等[86]利用Cable、Pile、Beam等研究锚杆、锚索、框架梁加固边坡的效果，其对预应力的实现方式处理略显不同：在模拟中，首先建立锚固段并赋参，此时自由段锚固参数为0，但对自由段施加拉力。计算一定时步后，再赋予自由段锚固参数并计算至平衡。研究认为：锚杆对边坡深度的加固效果有限，但框架梁形式的加固可显著控制顺坡向最大水平位移。同时，锚杆索的间距及预应力大小对边坡的位移和应力有重要影响，前者负相关，后者正相关。刘波等[87]采用砌体节理模型来描述破碎围岩，并利用Beam来建立拉杆，研究了锚杆支架支护顶板的效果，并认为对于节理化顶板，60°安装角更为有利，但锚杆轴向力则对围压拉应力的作用有限。吴静[88]通过Pile、Cable分别建立桩和锚杆来模拟基坑工程中超前钢管桩及预应力锚杆的支护效果，其采用了RBF神经网络位移反分析的方法确定模型参数，分层开挖及支护的研究及敏感性分析表明，敏感度较大的参数包括岩体参数(容重、内摩擦角、粘聚力、抗拉强度)、预应力锚杆参数(预应力、自由段长度、倾角)、钢管桩的长度及位置以及开挖深度和施工载荷等。高美玲[89]分别利用Pile、Cable、Shell建立桩、土钉、面板层研究单支点双排桩复合锚杆支护结构的效果，认为锚杆自由段拉力相同，锚固段则逐渐衰减。隋明昊等[90]在边坡失稳判据方面综合了以静力平衡收敛和剪应变增量两个因素，在其模拟中，采用Cable来模拟植被根系，并使用Shell来模拟植生基材。对坡度变化的研究表明，坡度增大，边坡稳定性系数普遍降低，复合型锚杆-土工网垫喷播结构有助于提升边坡稳定性。戚铧钟等[91]进一步对加筋层厚度和锚杆间距的研究指出：基质中的加筋层厚度应被尽量降低，合理选择锚杆间距以避免过大或过小的不合理间距。卢涛[92]针对生态护坡这一特殊的边坡问题，分别采用Cable、Geogrid建立锚索和土工网垫，并采用实体单元建立喷播植生层土体和边坡岩体，分析认为植生层坡向位移与其厚度、含水量和坡度正相关，但与植生密度负相关。张鹤[93]采用Cable和Beam建模，对参数的研究表明，边坡平台宽度增加，锚杆轴力和边坡最大水平位移会减小，锚杆长度的作用体现在一定范围内，且不宜过长，直立开挖的边坡，支护倾角以5～10°为宜，粘聚力和内摩擦角的增加均可导致安全系数的增大，施加预紧力可提高土体粘聚力和边坡稳定性。吴文清等[94]综合CAD、ANSYS、Fortran程序等建模，并在FLAC中利用Cable、Beam来模拟锚杆和框格梁，其中框格梁与锚杆刚性连接，认为大气对土体的影响深度为2.5 m，并相对于原状土进行了50%的强度折减。研究认为，边坡倾角降低、锚杆倾角增大、锚杆长度增加、垂直间距减小等有助于提高边坡稳定和安全系数。方高奎[95]利用Cable、Shell、Beam等建立了双排桩复合锚杆的分析模型，分析认为双排桩可与锚杆发挥协同支护作用，其中前排桩是主要受力构件。史江伟[96]利用Pile、Shell等建立双排桩的数值分析模型，并认为在双排桩中加入锚杆，在控制位移、弯矩等方面具有积极意义。同时，其研究了双排桩-锚杆支护的影响因素，认为适当增加桩径、桩长、锚杆长度可控制位移，但超过一定限度时，控制的效果并不会十分明显。排距过大或过小均不利于双排桩效果的发挥，排距应控制在3～6倍桩径内。同时，基坑深度及锚杆预应力对支护效果的发挥也有重要影响。此外，边坡研究中一般还对强度折减和安全系数格外关注[97]。图7所示为岩土工程中典型的研究案例。

图7 岩土工程中的典型研究案例[88,89,92]

与煤矿巷道支护相类似，边坡、基坑等岩土工程的数值模拟中一般建立的分析模型也比较简单，如：余科等[98]建立了分级开挖的相对精细化的模型，类似的研究则较少。岩土工程中一般锚杆的长度更长、安装角度要根据边坡的形状进行调整，而且桩、框架等结构被应用于不同工况。较多的研究在指定平台上研究支护形式、开挖坡度等，但没有从整个边坡系统来考虑不同的开挖进度、不同组合开挖与支护形式的影响，边坡上下台阶之间的变化相对比较单一、简单，没有精细化刻画出实际边坡中非均匀性、非线性的特点。未来的研究中，有必要在精细化建模、开挖工程从技术安全走向审美与技术安全协调设计、多因素耦合分析等方面进一步深入研究。

2.5 特殊锚杆的数值模拟

图8 变径锚杆模拟模型

此外，还有学者模拟了自旋锚杆、弹柔锚杆、GFRP锚杆、多锚点锚杆等[109-112]，但在这些研究中，却没有阐明如何实现不同功能锚杆的模拟、如何控制模拟计算的时步等关键问题，以至于降低了这些研究结论的可信度和参考价值。

由于锚固理论众多、锚杆支护设计的形式多样、深部锚固不断出现新的挑战和新的工程实践[113]，还有大量的研究通过FLAC3D软件对相应的观点进行了验证。但其建模方式与研究视角一般与上述研究类似，在此不再赘述。而张宇翔[114]将数值模拟与遗传算法相耦合，并以此获得最优的崖体加固方案；王茂源[115]综合fish语言、C号语言接口、BP神经网络模型等，开发了以FLAC3D为技术支撑主体的混合智能系统，这些研究进一步扩展了FLAC3D在跨学科、跨平台方面的应用。

3 FLAC3D中锚杆支护研究存在的问题

3.1 复杂工程问题的动态还原能力相对较弱

当前煤矿、岩土工程等领域运用FLAC3D的研究主要聚焦于拉拔试验、抗剪切特性、局部支护问题等，研究内容以静态、局部研究居多，无法精细化还原真实的复杂工程动态过程，而且研究中还存在对复杂工程的适当简化处理。如在拉拔试验中，真实过程包括钻孔、锚杆安装、拉拔这几个步骤，而在数值模拟中，钻孔及锚杆安装过程对围岩的影响一般是被忽略的；又如巷道支护的模拟中，大部分研究都是选取局部的巷道作为分析重点，并研究一次支护中的支护参数，忽略施工进度等动态因素，甚至托盘、金属网等辅助支护构件在很多研究中也被忽略掉。但真实的工程中，掘进、一次采动后，巷道还要受到开采及其他辅助工序的影响，FLAC3D中构建的局部巷道并未能继承该动态过程所形成的复杂应力条件，同一矿井掘进支护与采动影响、二次支护的模型往往是分段、相互独立的，与工程实际不对应。此外，受软件本身及计算机硬件的限制，目前高度还原真实工程环境的精细化建模研究相对较弱，无法完美的在FLAC3D中重现复杂断层、层理、节理、褶皱、水、气体等真实的地质赋存条件，多场耦合的研究和应用相对较少。

3.2 时步与时间不对应

FLAC3D中最为常用的是“时步-三维空间”的解算模式，而真实的工程场景通常可以简化为牛顿的“一维时间-三维空间”四维绝对时空惯性参考系，FLAC3D中的时步与真实的时间是不对应的。比如巷道掘进过程，FLAC3D中是以“设置开挖距离—运算时步—进行支护—运算时步—继续开挖—运算时步—循环上述步骤”来实现的，在这一过程中尽管有事件发生时间的先后，但并没有体现出时间的度量意义，就会造成时步的多少需要由人来控制，而不是真实的时间流逝和客观反映。现实的工程场景中，掘进支护工作可以看做是一个“连续—间隔—连续—间隔”的循环过程，而且真实工程中还会存在检修班，就会造成更长时间上的间隔。此外，使用Cable结构单元施加预紧力通常需要“施加预紧力—运算时步达到平衡—关闭预紧力设置”三个步骤来完成，如果不关闭预紧力设置，预紧力的设置就会覆盖法向力—位移准则，如何判断达到平衡、需要多少时步来达到平衡显然也是因人而已的。由于目前诸多研究中并未考虑时步与时间的匹配性问题，FLAC3D可以被用来做趋势性分析，定量化分析中的量与真实世界的量之间存在差异就在所难免了。

3.3 研究细节缺失及二次开发的非公开性问题

如前所述，Cable和Pile结构单元包含了众多的属性参数，但大量的研究却没有呈现使用FLAC3D的研究细节。一些研究中罗列了网格数，却没有给出网格尺寸划分的详情。甚至在一些研究中提出采用增大围岩参数的方法来实现对围岩强化效果的模拟，但增大多少却是未知之数。而在利用fish语言的二次开发中，一些研究公开了开发原理，却没有公开程序源代码，造成研究的可重复性较差、代码的正确与否不得而知、真实性受到怀疑。特别是在一些研究中，仅展示部分研究结果用以证明其研究的有效性，而没有公开参数和设置的细节，会让这种证明的效果大打折扣。由于与复杂工程不匹配、时间与时步不对应，加之研究细节的缺失，很多研究中得到与实验或实践完全相同的结果，反而是值得怀疑的。试问：一个简化的模型中，如何可能得到与真实复杂条件完全相同的结果呢？非公开性问题造成另一个结果是，即便采用完全相同的数据、类似的数值分析模型，有时也不能获得类似的模拟结果，使得研究不能沿着纵向深入高效开展。

3.4 对理论研究的支撑作用相对薄弱

目前的研究较多的从工程优化视角出发去探讨支护参数、支护形式以及由此带来的支护效果改变，却较少从原创性理论出发去探讨相应的支护机理、灾变演化规律。在支护领域，不论是支护理论还是新的支护工具，都在快速发展。但FLAC3D中对各种理论研究的支撑作用相对薄弱。如锚杆支护中的悬吊理论，在FLAC3D中怎样体现悬吊、怎样体现悬吊所发挥出的支护作用、怎样区分悬吊理论与其他理论的差异等，都没有得到很好的模拟。又如众多研究推导了锚杆轴向应力、剪切应力的分布形式，但在数值模拟中却未能进一步将相应的理论重现，模拟结果仍然以较为简单的“大”“小”来评价。又如锚杆支护中的松动圈理论，实际设计时应避开松动圈，但在模拟中，往往同一水平的岩层是相同的属性，并没有体现出松动圈的意义。对于锚杆支护过程中的能量的转移和能量释放等问题，也未能进行深入的研究。

4 FLAC3D中锚杆支护的研究展望

作为矿业、岩土等领域被广泛应用的软件，尽管FLAC3D目前还存在一些缺憾，但由于其庞大的用户群和产品的不断完善，可以预见，FLAC3D在未来依然拥有广阔的应用前景，如：

1) 新型支护理论和支护技术的验证。支护问题既是一个基础性问题，需要保障工程临时支护或永久支护的有效性；同时又不断面临着新的挑战，如深部开采中高应力、软岩巷道、地温影响、溶蚀性环境、复杂构造、冲击地压等。截至目前，锚杆支护已经提出悬吊、组合梁、组合拱、松动圈、中性点、锚固强化等多种理论[113]，与之相对应的就存在不同的锚杆支护设计方法。随着新型支护材料的研发、新的支护理论的提出，FLAC3D依然是验证新型支护技术的重要工具。此外，针对FLAC3D本身内置模型的修正[11-16]，也需要使用FLAC3D来进行对比分析。但在未来的应用中应注意：新技术新理论新在哪、在FLAC3D中实现这种新技术的途径以及如何体现模拟与现实的有机结合。

2) 跨平台建模及跨平台研究。长期以来，数值模拟研究一般局限于同一软件平台，而没有发挥出跨平台建模及研究的优势。目前，FLAC3D的网格导入支持ABAQUS、ANSYS格式直接导入，还可通过挤出命令建模，或通过Midas、Kubrix、Rhino等建模导入。如图9所示为笔者采用Rhino建模后利用FLAC3D进行实体单元拉拔模拟结果，可见跨平台建模可以更精细地还原实际情况。但由于跨平台目前的应用还相对较少，在研究过程中也会遇到很多新问题。FLAC3D6.0版本后，可实现与PFC的协同分析，探索区域(Zones)与颗粒(Particles)之间的相互作用。FLAC3D7.0版本中可以实现和UDEC的耦合分析。从原理上讲，一个网格单元形成后，在运算过程中就会在网格节点与区域内形成不同的作用参数，将原本在一个平台进行的模拟运算拆分成不同阶段，调用前一阶段的数据用于后一平台的分析，如此持续往复，则可以实现不同平台间的数据共享和连续分析。如锚杆螺纹的磨损可以通过ANSYS来研究，相应的结果则可以被调用到FLAC3D中，实现跨平台的数据共享和连续分析。与单一模拟平台相比，跨平台研究为精准建模、综合分析提供了可能，也可能催生更准确的定量化描述。

图9 Rhino联合建模研究结果

3) 精细化建模与多因素耦合研究。复杂精细化建模技术之所以不能快速发展，与民用计算机硬件的配置关系密切。未来，随着计算机高速计算、快速传输、海量存储、云计算、跨平台数据共享等的发展和普及，复杂精细化建模有望得到更广阔的应用，从而可以在FLAC3D中表征更丰富的信息。如图10所示为笔者结合复杂地表形貌建立的分析模型，但从锚杆支护而言，实际上这种精细化建模只是更加精细的刻画了围岩，而非锚杆，其优点在于精细化建模后有可能与实际情况更贴近，但缺点在于非线性、非均匀性等问题会更加突出，在分析结果时可能要比传统的模拟更为费时费力。目前在FLAC3D中，锚杆与锚索一般仅通过系统的属性参数予以区分，流固耦合、气固耦合等研究还相对较少，不能很好地实现水压作用、水的软化、气体的吸附解吸、高含硫腐蚀等。未来，随着跨平台建模的逐渐完善和精细化建模的发展，多因素耦合分析有可能借助FLAC3D这一平台实现新的突破。在此基础上，可以提高对复杂工程问题的动态还原，并借助高精度还原的运算结果，建立运算时步与一维时间、客观工程实际之间的关系。

图10 基于复杂地表形貌建立的分析模型

5 结 语

当前的锚杆支护中较多的采用Cable或Pile结构单元，虽然参数众多，但在多场耦合方面仍然存在不足。尽管如此，FLAC3D依然是锚杆支护数值模拟的首选工具之一。目前，Itasca咨询公司拥有3DEC、FLAC、PFC、UDEC等众多岩土分析软件，有利于跨平台数据的融通，不同模拟软件之间的交互使用在未来变得更具可行性。随着技术的发展，未来有望实现复杂工程实际的高度还原和跨平台耦合建模，扩大FLAC的应用前景。对于数值模拟中研究细节缺失问题，则需要从数据共享、期刊要求等方面予以改进，在透明、公开、知识产权保护的前提下，更好地推动科技发展。