基于FLAC 3D的溶洞顶板稳定性分析

马中伟

(贵州有色地质工程勘察公司，贵阳 550002)

1 概 述

在我国西南地区，地表、地下岩溶发育，地下岩溶发育形成各种奇峰异洞，这些地下洞穴千奇百怪、美轮美奂，在一定程度上给旅游资源带来一定的有利条件，而且人民的日常生活和经济建设与岩溶又有着密切的联系[1]。另一方面，岩溶具有不利于人类工程活动的灾害性特点，在岩溶地区溶洞顶板稳定性这一困扰工程建设的大问题一直以来对工程地基稳定性有着极大影响，因此研究岩溶的发育特征及规律，探索卓有成效的治理方法对工程具有重大意义[2]。本文以贵州省凯里市台江县水投工司所委托的项目——污水处理厂区内典型溶洞为例，通过FLAC 3D 软件对溶洞进行数值模拟，对研究区溶洞顶板稳定性进行详细的分析与评价。

2 工程概况

研究区位于贵州黔东南州台江县，台江县位于贵州省东南部、黔东南苗族侗族自治州中部，E108°06′-108°29′，N26°23′-N26°52′之间。属中亚热带湿润季风气候区，场研究区整体地处云贵高原向湘西丘陵过渡的大斜坡上，地貌为深切割中心类型,原始地形为山地，较为平坦，最高处位于场地南侧，最低处位于场地北侧，高程在596.16～608.19 m之间。研究区地表地下岩溶发育。本文选取场地具有代表性的钻孔ZK-226和ZK-227钻遇两个大小不同的溶蚀洞体作为研究对象。

3 模型的建立

3.1 模型假设条件

溶洞顶板稳定性的影响因素繁多，在建模分析过程中对各种因素做减法，简化计算模型。比较常见的方法是针对具体的问题具体分析，最终得到普遍性的成果，在工程问题上，进行合理的假定具有可行性。本次计算模型的假设条件如下：

1) 溶洞岩体为位于地下的各向同性的、均匀的空间半无限体，层面、节理裂隙的影响折减至计算参数。

2) 外加荷载是导致溶洞顶板失稳的直接原因，因此溶洞在不受外加工程荷载的条件下是稳定的，且洞体围岩岩体视为同向均质弹塑性体。

3) 在一般情况下，溶洞充填物不承受应力，这是由于充填物的松胀破碎等特征导致与岩石相比具有非常低的力学特性，因此按无充填物的空溶洞进行计算。但最终所得的计算结果进行对比分析时，溶洞充填物对溶洞稳定有利。

4) 由于研究区地下水埋藏深度在100 m以下，因此可以忽略地下水对溶洞的影响。溶洞形状按椭球体考虑。见图1、图2。

3.2 计算模型的建立及参数取值

模型计算需要的岩土体参数主要有剪切模量(G)、体积模量(K)、密度(ρ)、黏聚力(C)、摩擦角(φ)、剪胀角(Ф)及抗压强度(σ)。其中，剪切、体积模量利用弹性模量、泊松比按下式换算：

式中：E为变形模量；μ为泊松比。

依据室内、原位试验数据，计算参数选取见表1。

图1 研究溶洞钻孔平面布置图

图2 研究溶洞钻孔剖面图

表1 各岩土层物理力学参数

地质概化后的模型分为上下两层：

ZK-226上覆土体黏土层，厚度为2.9 m，溶洞的跨度为17.25 m，高11.5 m，深度为12 m，顶板厚度为9.4 m。计算模型见图3。

ZK-227上覆土体黏土层，厚度为3.5 m，溶洞的跨度为3.75 m，高3.6 m，深度为7 m，顶板厚度为13.5 m，见图4。

图3 ZK-226计算模型

图4 ZK-227计算模型

4 模拟结果分析

4.1 最大不平衡力分析

模型在进行数值计算迭代计算过程中网格节点系统的内外力之差称为最大不平衡力，最大不平衡力表示在单元网格节点上所施加力的代数和，当最大不平衡力接近为零或等于零时，模型处于稳定状态；反之，当最大不平衡力从最初的变化跳跃，最终趋近于一个非零的恒定值时，模型进入塑性流动状态并趋于破坏。本文以最大不平衡力作为评价模型稳定性的一个重要指标，选取场地范围内ZK-226和ZK-227大小不同的地下溶洞为研究对象，逐级施加13 000 kN荷载作用下的最大不平衡力，见图5、图6。

图5 ZK-227钻孔地下溶洞初始地应力及荷载作用下的最大不平衡力监测图

图6 ZK-226钻孔地下溶洞初始地应力及荷载作用下的最大不平衡力监测图

对图5、图6分析可得出以下结论：

1) 由图5可知，ZK-227溶洞在初始应力作用(自重应力)及荷载作用下，最终均处于稳定状态。整个计算过程中，随时间变化，最大不平衡力呈现收敛状态，迭代效果好，系统达到平衡，节点最大不平衡力为6.635×104N。

2) 由图6可知，ZK-226地下溶洞在初始应力状态下节点只有微小的跳跃变化，并最终处于稳定状态。但施加荷载以后，计算的迭代过程相对缓慢，且节点力代数和最终是非零的值，该值不恒定但跳跃弧度小，说明在施加外部荷载以后模型被破坏。在整个计算过程中，最大不平衡力为9.516×104N。

4.2 位移变化分析

见图7-图10。

图7 ZK-227在初始应力下竖直方向位移云图

图8 ZK-227 在外荷载作用下竖直方向位移云图

图9 ZK-226在初始应力下竖直方向位移云图

图10 ZK-226在外荷载作用下竖直方向位移云图

1) 位移变化分布特征是描述模型是否发生塑性变形的又一重要物理量。对比分析图7、图8、图9和图10，荷载作用下位移分布特征较之初始应力作用下位移分布有明显变化，初始应力下的位移沉降小于外荷载作用下位移沉降。

2) 对图7、图8分析可知，在初始应力条件下及荷载作用下的稳定洞体，从宏观上分析可以看出，ZK-227遇钻溶洞位移量自上而下逐渐增加，这是由于处于稳定状态的溶洞，主要是受到地应力的影响较大所致。但从微观上分析来看，其最小位移出现在溶洞底部浅层范围，这是由于该部位受到洞体自身减压拱的作用。

3) 对图9及图10对比分析可以得出，虽在初始应力作用下溶洞处于稳定状态，但由于扁平溶洞的岩体临空方向不能提供围压而易于向内发生弯折破坏，因此从图9可以看出其最大位移出现在溶洞洞顶；在外荷载作用下的溶洞ZK-226出现顶板失稳，最大位移出现在上覆土层处，最大位移量达到81 mm，而在施加荷载处土体左右两侧出现位移上升，这是由于这两处受到外荷载及周边岩土体共同挤压作用所致，使其出现凸起。因此，外荷载的作用对于溶洞顶板的稳定性有不利的影响。

4.3 应力变形分析

见图11-图14。

图11 ZK-227在初始应力下竖直方向应力云图

图12 ZK-227在外荷载作用下竖直方向应力云图

图13 ZK-226 在初始应力下竖直方向应力云图

图14 ZK-226在外荷载作用下竖直方向应力云图

从图11至图14分析可以得出，无论是在初始应力条件下还是荷载作用下，最大主应力均以压应力为主，主要集中在溶洞洞壁左右两侧或顶板中部，这是由于受到围压作用的影响，在荷载作用下，溶洞周围岩体的竖向应力增大，ZK-227遇钻洞体最大压应力量值为-105.4 kPa，最小主应力以拉应力为主，最大拉应力量值为0.383 kPa。ZK-226遇钻洞体荷载作用下最大压应力量值为-6 291 kPa，最小主应力以拉应力为主，最大拉应力量值为101.45 kPa。因为此处围岩局部已经开始破坏，说明拉压应力是造成溶洞破坏的主要原因，且应力呈对称分布的特征。

4.4 塑性区分析

在FLAC 3D中，塑性区域分布有两种表示方法，分别是now和past。其中，正处于屈服面上的模型单元表示为now，过去曾经处于屈服面上的单元表示为past，past的区域表示此刻已经离开屈服面，处于弹性范围。模型中塑性状态为now的单元，是正处于塑性流动状态的那些区域。因此，只有now状态的单元针对实际工程问题分析塑性区域时才对模型的破坏起作用[3]见图15-图18。

图15 ZK-227在初始地应力下的塑性区分布图

图16 ZK-227在外荷载作用下的塑性区分布图

图17 ZK-226在初始地应力下的塑性区分布图

图18 ZK-226在外荷载作用下的塑性区分布图

对图15、图16、图17、和图18的分析可得出以下结论：

1) 在自重应力场条件下的ZK-227遇钻溶洞，其洞壁出现极少部分塑性破坏区，但该破坏范围不足以影响顶板稳定性；在外荷载条件下，ZK-227遇钻溶洞上覆土体浅层部位及顶板中部发生塑性破坏，但顶板中部塑性破坏区不连续，因此在荷载作用下亦不会对顶板稳定性造成影响。

2) ZK-226遇钻溶洞，初始状态稳定，而受到外荷载失稳，整个洞体处于失稳状态，且破坏形式以剪切破坏为主。

5 不同填方高度对溶洞顶板稳定性的影响

研究区填方区域面积占整个场区的70%，因此研究填方对溶洞顶板稳定性的影响具有重要意义。

见图19-图22。

图19 2 m填方高度塑性分布

图21 6 m填方高度塑性分布

图20 4 m填方高度塑性分布

图22 10 m填方高度塑性分布

图19至图22是溶洞顶板在不同填方高度下的塑性分布图。从图19-图22中可以看出，随着填方高度的增加，塑性破坏区增大，塑性破坏区集中分布在下伏土层处，2～6 m的填方高度均对溶度顶板无影响，在10 m填方高度下，顶板上部以及溶洞左右两侧出现了部分不连续的塑性破坏区，可以判断在10 m以内的填方高度下，该类型溶洞顶板较稳定。

图23至图26是溶洞顶板在不同填方高度下的竖向位移发展过程曲线。1-Zdis代表溶洞顶板下部关键点竖向位移，2-Zdis代表溶洞顶板中部关键点竖向位移，3-Zdis代表溶洞顶板上部关键点竖向位移。

图23 2 m填方竖向位移发展过程

图25 6 m填方竖向位移发展过程

图24 4 m填方竖向位移发展过程

图26 10 m填方竖向位移发展过程

从图23-图26中可以看出，位移发展过程曲线收敛，若以顶板关键点位移最大点作为位移变化的判别标准，可得到不同填方高度下关键点的最大位移值。当填方高度分别为2、4、6和10 m时，1-Zdis此关键点位移分别为-0.01、-0.18、-2.58和-5 mm。可见，虽然溶洞顶板的位移量会随填方高度的增大而增大，但减幅较小，表明当溶洞顶板较厚且埋藏较深时，其稳定性受填方高度影响很小。

6 结 语

1) 荷载作用下位移分布特征较之初始应力作用下位移分布有明显变化，初始应力下的位移沉降始终小于外荷载作用下位移沉降。

2) 拉压应力是造成溶洞破坏的主要原因，且应力呈对称分布的特征。

3) 溶洞顶板的位移量会随填方高度的增大而增大，但减幅较小。因此当溶洞顶板较厚且埋藏较深时，其稳定性受填方高度影响很小。