楼晓明 黄 慎 韩雪靖 周 平

（1.福州大学紫金矿业学院，福建福州350116；2.福州大学爆炸技术研究所，福建福州350116）

矿山采空区在矿柱流变、分化及剥落等因素的影响下，空区顶板会发生大面积冒落，引发严重的矿山事故［1-2］。因此建立合理的力学模型并预测采空区顶板稳定时间，对于指导采空区进行及时有效的支护、充填等工作具有非常重要的现实意义。目前，已有很多学者对采空区稳定性进行了研究。在采空区稳定时间方面，来兴平［3］等采用BP神经网络等方法，构建出具有多个自由度的神经元，将其耦合成适用于力学问题的网络，对采空区围岩损伤进行计算，得到时间序列模式下采空区围岩环境损伤破坏的演化趋势；李爱兵［4］等以结构风险最小化为基础，采用SVM算法，对现有工程的实际数据进行学习回归，确定了采空区稳定时间；张源［5］等利用工业电子窥视仪对顶板进行摄像观察，通过统计分析其摄像数据，得出顶板裂隙演化规律及顶板破坏形式，预测了顶板的稳定时间。但由于国内矿山的数据资料存在着不连贯、不及时、不可靠等问题，所以采用基于数据库的算法确定的采空区稳定时间存在较大误差。另外，实地摄像、监测等方法存在实施过程复杂，推广性不强等缺陷。在采空区顶板—柱体系力学模型方面，常见的一类模型是将采空区顶板简化为岩梁，将矿柱视为Winkler弹性地基基础［6-7］。岩梁模型虽简单易算，但存在未考虑顶板空间效应及顶板各向异性的缺点。为解决岩梁模型存在的问题，基于弹性薄板理论，J A Wang［8］等将顶板简化为弹性薄板，矿柱视为均布弹簧建立力学模型，但其未考虑矿柱本身的流变；孙琦［9］等同样将顶板简化为弹性薄板，并将矿柱视为Burgers体建立力学模型，但在进行解析计算时，存在未完全求解的缺点；李铁［10］等在建立力学模型时，将矿柱视为广义Kelvin体，得到了采空区变形发展与时间的关系，但其模型参数较少。

为了更好反映岩石特性，得到更具有代表性的采空区流变方程，选用Hooke-Kelvin-kelvin模型建立矿柱—顶板流变力学体系，结合弹性薄板理论，分析顶板流变失稳过程，预测采空区顶板的稳定时间。

1 流变力学模型的建立与求解

1.1 模型建立

将采空区顶板简化为弹性矩形薄板（如图1），设该弹性矩形薄板长为2a，宽为2b（b≤a），厚为h。矿柱视为Hooke-Kelvin-kelvin模型（如图2），采空区矿柱—顶板体系简化模型见图3。假设矿柱等间距分布，矿柱总的数量为n，高度为h，平均面积为A。

则Hooke-Kelvin-kelvin模型的本构方程可表示为：

式中，E1，E2，E3为弹性系数；η1，η2为黏性系数。令：

则矿柱本构关系可简化为：

矿柱应变与位移关系为：

1.2 矿柱—顶板体系微分方程的建立

基于弹性薄板理论，采空区顶板控制方程为：

其中，

式中，D为薄板抗弯刚度；∇为拉普拉斯算子；ω为薄板挠度；ξ为将矿柱应力等效成均布应力系数；σ为矿柱内应力；q为顶板所受均荷载；E为弹性模量；H为薄板厚度；μ为泊松比。

求采空区的顶板变形转化为求式（2），（3），（4）解的问题。将Hooke-Kelvin-kelvin模型本构方程，矿柱应变与位移关系式，代入顶板控制方程式可得：

式（7）即为采空区顶板—矿柱体系的微分方程。其中，ω̈，ω̇分别为顶板下沉挠度ω对时间的二阶导数和一阶导数。

1.3 矿柱-顶板体系微分方程的求解

采用Galerkin法求解式（7），设挠度有如下解析解的形式：

式中，ω0（t）为采空区顶板中心点的下沉位移，亦即采空区顶板最大下沉位移。则有

则有：

式（11）为一个二阶非齐次线性微分方程，其通解为：

式中，C1，C2为积分常数；r1，r2表达式如下：

1.4 积分常数的求解

空区顶板的变形破坏过程：初始状态，空区顶板四周边界固支，此时空区顶板变形沉降导致边界破坏，随后长短边发生破坏，边界状态由固支转为简支。具体见图4。

（1）当采空区顶板未发生破坏时，其边界为固定边界，边界条件为：

设顶板挠度的解析解的形式为:

显然f（x，y）满足所有边界条件，将f（x，y）代入式（10）可得：

矿柱在初始状态受压时会产生瞬时弹性变形，此时的矿柱有效刚度为k2，其位移量［11］为：

随着时间的推移，顶板位移逐渐增加，其破坏条件［12］为：

（2）内部破裂边界条件为：

假设顶板的挠度解析解形式为：

则可解得：

此时的破坏条件［12-13］为：

式中，σxmax为顶板中点x方向最大拉应力；σymax为顶板中点y方向最大拉应力。

2 工程实例分析

福建省龙岩地区某矿于2008年2月对其矿山矿产资源运用房柱式采矿方法进行开采，采空区平面尺寸为2a×2b=162×143 m，矿柱高为6 m，矿柱总的横截面面积与采空区顶板的面积比值为35.7%。采空区顶板上覆岩层厚度134 m，容重为22 kN/m3，采空区顶板厚度为H=34 m，容重为25 kN/m3，采空区顶板岩石平均弹性模量值为E=48 GPa，泊松比为μ=0.26，采空区顶板岩石顶板平均抗拉强度值为［σT］=3.6 MPa，矿柱流变参数为E1=13.4 GPa/m，E2=87 GPa/m，E3=296 GPa/m，η1=32 GPa·h，η2=5.4×105GPa·h。

由以上参数可得，采空区顶板上覆岩层及采空区顶板自重总的均布载荷为：

采空区顶板岩石的抗弯刚度为：

初始时刻采空区顶板挠度为：

采空区顶板初始下沉速度为：

由式（15）可知，固定边界变为简支边界顶板破坏条件如下。

对于采空区顶板长边：

对于采空区顶板短边：

将C1，C2，ω0a，ω0b代入式（12）可求得t1≈34月。

当第二阶段发生破坏时，根据式（17）可求得t2≈47月。

由此可求得该采空区总的稳定时间t为81月，预测该铜矿采空区将在2014年11月发生采空区大面积顶板冒落。现场监测数据表明，在2014年9月，该铜矿采空区发生了采空区顶板局部冒落现象。该采空区的现场监测数据验证了该流体力学模型的可靠性。

3 结论

（1）和以往研究者采用的Burgers体、Kelvin体、广义Kelvin体建立模型相比，本研究建立模型采用的Hooke-Kelvin-kelvin模型能更好地反映岩石的蠕变特性，不仅考虑了岩石的初始弹性及流变特性，且参数适中，较好地兼顾了理论计算与工程实际的需要。

（2）基于Hooke-Kelvin-kelvin模型建立的采空区流变力学模型对采空区的稳定时间进行了预测，并将其与现场实际的监测数据进行对比和分析，实践结果验证了建立的采空区流变力学模型的合理性与可靠性。

（3）利用Hooke-Kelvin-kelvin模型建立空区模型进行计算时，没有考虑矿柱及顶板受节理、地下水、分化作用等外界因素的影响，如果要对这类外界因素对采空区顶板稳定性的影响进行研究，可通过对一些相关的系数进行适当的折减。