在上覆岩层压力及开挖力学效应影响下，地下采场及硐室顶板是整个工程中相对薄弱的部分，当开挖空间、应力分布等发生变化时，极有可能诱发其失稳破坏，对正常生产及生命安全形成极大威胁^[1-3].据统计数据分析，地下工程作业过程中，顶板事故是各种安全事故中最为突出的一类.因此，顶板稳定性的识别预测对地下工程安全作业有着极为重要的意义.

目前对顶板稳定性测试分析的方法较多，监测方法主要包括变形观测法、声发射法、位移监测分析法等^[4-6]，数据分析方法多采用经验、解析和半定量的方法^[7-10].但顶板变形破坏是一个复杂的非线性过程，从承载到失稳整个过程中所释放的信息较多，仅通过某种监测技术和分析难以全面反映顶板的应力、应变分布和破坏情况，故对顶板稳定性评价也存在较大的局限性和片面性.为此，针对顶板岩体变形过程所释放的最为直接的多方面信息，利用有效的数学方法分别对其进行量化提取并组合分析，以实现地下工程顶板稳定性的准确评价与预测.

1 顶板变形位移识别

地下工程顶板一旦被揭露，即发生相应的变形，该时段的变形称为开挖变形，随着开挖变形结束，对不同的岩体而言，顶板变形增量逐渐减小或消失，其变形处于稳定阶段，当岩体处于非稳定状态时，其变形增量又逐步增大，该时段称为加速变形阶段，直至顶板破坏失稳.因此，顶板岩体在开挖完成后整个承载过程中可能出现的位移变化阶段主要包括无变形或稳定变形阶段、加速变形阶段.

1.1 无变形阶段

该阶段说明开挖变形完成后不产生任何时效变形，如果位移值发生细微变化说明是由于测量仪器的误差而引起的，因此，实际监测的位移值的变化就是仪器测量误差值.即：

$ {d_i} = {\delta _i} $

(1)

其中d_i为每次监测的位移值，δ_i为测量误差，测量误差δ_i由仪器的精密程度引起，属于随机误差，通常是由多种相互独立的因素造成的许多微小误差的总和.其变化服从正态分布N（μ，σ）.共监测n组数据后，测量值的标准差和平均值为：

$ \hat \sigma = \sqrt {\frac{1}{{n-1}}\sum\limits_{i = 1}^n {{u_i}^2} } = \sqrt {\frac{1}{{n-1}}\sum\limits_{i = 1}^n {{{({d_i}-\bar d)}^2}} } $

(2)

$ \bar d = \frac{1}{n}({d_1} + {d_2} + \cdots + {d_n}) = \bar \delta $

(3)

其中$ {\hat \sigma } $为测量值的标准差、d为测量值的平均值，δ为误差的平均值，由式（1）和式（3）可以推出：

$ {d_i}-\bar d = {\delta _i}-\bar \delta $

(4)

进一步推得：

$ {{\hat \sigma }^2} = \frac{1}{{n-1}}\sum\limits_{i = 1}^n {{{({d_i}-\bar d)}^2}} = \frac{1}{{n-1}}\sum\limits_{i = 1}^n {{{({\delta _i} - \bar \delta )}^2}} $

(5)

式（5）中：

$ \begin{array}{l} \sum\limits_{i = 1}^n {{{({\delta _i}-\bar \delta )}^2}} = {\left( {{\delta _1}-\bar \delta } \right)^2} + {\left( {{\delta _2}-\bar \delta } \right)^2} + \cdots + {\left( {{\delta _n} - \bar \delta } \right)^2}\\ = ({\delta _1}^2 + {\delta _2}^2 + \cdots + {\delta _n}^2) - 2\bar \delta ({\delta _1} + {\delta _2} + \cdots + {\delta _n}) + n{{\bar \delta }^2} \end{array} $

(6)

此处δ_i服从正态分布N（μ，σ），所以$ {\delta _1} + {\delta _2} + \cdots + {\delta _n} = n\mu, \bar \delta = \mu $，由此，对式（6）进行整理：

$ \sum\limits_{i = 1}^n {{{({\delta _i}-\bar \delta )}^2}} = {\delta _1}^2 + {\delta _2}^2 + \cdots + {\delta _n}^2-n{\mu ^2} $

(7)

由于δ_i服从正态分布N（μ，σ），对其进行标准化变换，得到$ {T_i} = \frac{{{\delta _i}-\mu }}{\sigma } $服从标准正态分布N（0，1），进行推导可得：

$ {\delta _i} = {T_i}\sigma + \mu $

(8)

将式（8）代入式（7）得：

$ \begin{array}{l} \sum\limits_{i = 1}^n {{{({\delta _i}-\bar \delta )}^2}} = {\left( {{T_1}\sigma + \mu } \right)^2} + {({T_2}\sigma + \mu )^2} + \ldots + {({T_n}\sigma + \mu )^2}-n{\mu ^2}\\ = {\sigma ^2}({T_1}^2 + {T_2}^2 + \ldots + {T_n}^2) + 2\mu \sigma ({T_1} + {T_2} + \ldots {T_n}) \end{array} $

(9)

由于T_i服从标准正态分布N（0，1），故${T_1} + {T_2} + \ldots {T_i} = 0 $，而$ {T_1}^2 + {T_2}^2 + \ldots + {T_i}^2 $服从自由度为n的卡方分布，即：

$ {\chi ^2} = {T_1}^2 + {T_2}^2 + \ldots + {T_i}^2 $

(10)

其均值E（χ²）=n，所以$ {T_1}^2 + {T_2}^2 + \ldots + {T_i}^2 = {n^2} $，将以上结果代入式（2），整理后可得测量值的标准差：

$ \hat \sigma = \sqrt {\frac{{{n^2}{\sigma ^2}}}{{n-1}}} = \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} $

(11)

式（11）中σ为测量误差的标准差，即测量仪器的精度，n为测量次数，可见随着测量次数增多，测量值的标准差逐渐增大，所以应该尽量提高测量精度，以降低测量中的随机误差.同时，由于式（11）推导过程中测量值的标准差计算中只考虑了测量仪器误差造成的变化，而将该阶段顶板的实际变形值定为零.因此，如果顶板处于无变形阶段，则：

$ \tilde \sigma < \hat \sigma = \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} $

(12)

式（12）中，$ {\tilde \sigma } $为实际顶板监测位移值的标准差.如果$ \tilde \sigma > \hat \sigma $，则顶板不再处于无变形阶段，而过渡到其他变形阶段.

1.2 稳定变形阶段

该阶段表明顶板产生了变形，因此该阶段实际监测的位移值包含2个部分，一部分为顶板变形产生的位移，另一部分为测量仪器引起的误差值.即：

$ {d_i} = {\delta _i} + {w_i} $

(13)

式（13）中w_i为顶板实际产生的位移值，有关测量仪器引起的误差值δ_i前面已经论述，此处仅讨论顶板实际产生的位移值.

其均值为：

$ \bar w = \frac{1}{n}({w_1} + {w_2} + \ldots + {w_n}) $

(14)

其标准差为：

$ \widehat {\hat \sigma } = \sqrt {\frac{1}{{n-1}}\sum\limits_{i = 1}^n {{{({w_i}-\bar w)}^2}} } $

(15)

式（15）中$ \widehat {\hat \sigma } $为实际产生位移值的标准差，对公式中位移值与均值差的平方和进行整理，得：

$ \begin{array}{l} \sum\limits_{i = 1}^n {{{({w_i}-\bar w)}^2} = {{({w_1}-\bar w)}^2} + {{({w_2}-\bar w)}^2} + \ldots + {{({w_n} - \bar w)}^2}} \\ = (w_1^2 + w_2^2 + \cdots + w_n^2) - n{{\bar w}^2} \end{array} $

(16)

如果顶板岩体处于稳定变形阶段，则可以近似认为在等时间间隔内顶板实际产生的位移值的变化率相等，即等时间间隔内顶板位移值的增量一致.

$ ({w_1}-{w_0}) = ({w_2}-{w_1}) = \cdots = ({w_n}-{w_{n - 1}}) = t $

(17)

w₀为开始监测之前的位移值，取值为0，由此可知：

$ \begin{array}{l} w_1^2 + w_2^2 + \cdots + w_n^2 = {t^2} + {(2t)^2} + {(3t)^2} + \cdots + {(nt)^2}\\ = \frac{{n(n + 1)(2n + 1)}}{6}{t^2} \end{array} $

(18)

$ \bar w = \frac{1}{n}(t + 2t + 3t + \cdots + nt) = \frac{{n + 1}}{2}t $

(19)

将式（17）、式（18）、式（19）代入式（16）整理得：

$ \sum\limits_{i = 1}^n {{{({w_i}-\bar w)}^2} = } \frac{{n({n^2}-1)}}{{12}}{t^2} $

(20)

将式（20）代入式（15）得顶板实际位移值的标准差：

$ \widehat {\hat \sigma } = \sqrt {\frac{1}{{n-1}}\sum\limits_{i = 1}^n {{{({w_i}-\bar w)}^2}} } = t\sqrt {\frac{{n(n + 1)}}{{12}}} $

(21)

由于t为稳定变形阶段位移的增量，监测过程很难准确获得其实际值，故将其转化为最终监测的位移值，即：S=nt，S为顶板实际变化位移值的总和.故式（21）转化为：

$ \widehat {\hat \sigma } = t\sqrt {\frac{{n(n + 1)}}{{12}}} = S\sqrt {\frac{{n + 1}}{{12n}}} $

(22)

根据式（13）将测量仪器引起的误差值考虑在内，得到稳定变形阶段实际监测顶板位移值的标准差为：

$ \widehat {\hat \sigma } = S\sqrt {\frac{{n + 1}}{{12n}}} + \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} $

(23)

由此可以得到，在稳定变形阶段：

$ \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} \le \tilde \sigma \le S\sqrt {\frac{{n + 1}}{{12n}}} + \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} $

(24)

1.3 加速变形阶段

当顶板岩体处于加速变形阶段时，可近似认为在等时间间隔内顶板实际产生的位移值的变化率逐渐增大，即在相等时间间隔内，顶板位移增量逐渐变大.

$ ({w_n}-{w_{n-1}}) \ge ({w_{n-1}} - {w_{n - 2}}) \cdots \ge ({w_1} - {w_0}) $

(25)

根据1.2节推导过程，可最终得到在加速变形阶段：

$ \tilde \sigma > \widehat {\hat \sigma } = S\sqrt {\frac{{n + 1}}{{12n}}} + \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} $

(26)

通过对顶板位移监测数据进行统计分析，可以判别顶板变形阶段，但整个地下工程施工过程中，受复杂因素的影响，顶板应力状态不断变化，顶板3种变形阶段也可以相互转变，换言之，即使进入加速变形阶段，顶板也不一定失稳.仅凭位移监测难以分析顶板的稳定性状态.

2 顶板失稳声发射指标识别

岩体变形破坏过程伴随裂纹产生和扩展，裂纹产生与扩展又伴随着应变能的释放，即岩体破坏声发射现象.针对岩体破坏的声发射特征，国内外学者进行了大量的研究，通过室内试验及现场测试，揭示了岩石破坏过程与声发射特征参数的内在关系.

2.1 岩石承载过程声发射关联分维数变化规律

笔者及相关研究发现岩石破坏过程声发射特征参数具有明显分形特征^[11]，针对这一研究结果，笔者在声发射测量地应力的试验中发现^[12]：不同岩性的多组岩石试样在单轴加载全过程中，其声发射关联分维数总体在开始阶段和中间阶段，声发射关联分维数总体呈上升趋势，但在上升过程中也出现了不同程度的多次起伏，即在加载过程中多次出现关联分维数的上升和下降.当加载至岩石峰值应力前夕（即将破坏），声发射关联分维数呈现降低的特征，图 1是其试验结果曲线图.

2.2 岩石承载过程声发射长程相关指数变化规律

许福乐等^[13]提出了用声发射序列的长程波动来描述岩石破坏过程，通过不同煤岩试样单轴加载全过程声发射试验，分析了整个过程声发射强度序列的长程相关性，研究了长程相关指数H的变化规律，试验结果表明，在加载初期，H值增长缓慢，随着加载继续进行，H值不断增大，并达到最大值，随后又出现降低.但在整个过程中多次出现升高-降低的模式，直至加载至岩石达到加速非线性变形阶段（濒临破坏前夕），声发射强度长程相关指数H值出现了“最大-减小”模式，由此得出结论可以利用岩石破坏声发射强度长程相关指数H值的“最大-减小”模式判别岩体的失稳，见图 2.

2.3 岩体破坏声发射相对强弱指标变化规律

王宁等^[14]研究了评价采场岩体破坏的声发射相对强弱量化指标，对现场监测的声发射能率、总事件、大事件进行标准化处理，得到声发射事件率的相对强弱指标CR和声发射能率相对强弱指标ER，从图 3可以看出，在整个岩体承载声发射监测过程中，CR和ER总体呈现波动趋势，但在失稳前兆期，相对强弱指标出现迅速回落.

从笔者^[15]及相关人员对岩体承载声发射特性参数研究结果来看^[16]，都说明了同样一点，岩体加速破坏时，声发射监测指标频度与量度激增，岩体在失稳前夕，声发射监测指标频度与量度出现了迅速回落，对应于声发射b值的升降模式，声发射长程相关指数H值的最大-减小模式，以及声发射相对强弱指标CR、ER的起伏模式.因此岩体失稳前，其声发射监测指标频度（量度）出现激增-回落现象.但反之并不成立，从多次研究结果的变化曲线上来分析，在整个监测过程中，这种激增-回落现象在变化曲线中多次出现，也就是说这种激增-回落现象并不能说明岩体一定失稳，所以单纯利用声发射监测指标的这一性质并不能准确判断其稳定性状态.

3 位移与声发射组合模式

位移与声发射是岩体稳定性监测中最有效、最直观的伴生信息，上述分析可知单一方法无法准确判定顶板的稳定性状态，对其进行组合分析可以增加分析结果的准确性.对顶板不同的变形阶段和声发射信息指标进行组合分析，可以确定以下几种组合模式.

1）无变形稳定模式.该模式下监测周期内，顶板位移监测结果分析满足式（12），顶板处于无变形阶段，因此顶板不会破坏，也不会产生相应大量激增的声发射现象，即使有也是局部岩体受结构控制产生的小规模破裂.据此，该位移与声发射的组合模式可以分解为以下2种，见式（27）.均说明顶板处于稳定状态.

$ \left\{ \begin{array}{l} \tilde \sigma < \hat \sigma \quad 无变形阶段\\ 无声发射现象产生 \end{array} \right.或\left\{ \begin{array}{l} \tilde \sigma < \hat \sigma \quad 无变形阶段\\ 少量声发射现象产生 \end{array} \right. $

(27)

2）平稳变形稳定模式.该模式下整个监测周期内顶板位移监测结果分析满足式（24），说明顶板处于平稳变形阶段，该阶段内由于岩体产生了实际变形，变形结果可能形成2种情况，其一为岩体只发生变形并未产生内部结构的破裂，因此声发射相对“平静”；其二为岩体变形引发内部结构破坏，产生一定量的声发射现象.据此，该位移与声发射的组合模式也可分为2种情况，见式（28）.该组合分析模式也表明顶板处于稳定状态.

$ \begin{array}{l} \left\{ \begin{array}{l} \hat \sigma \le \tilde \sigma \le \widehat {\hat \sigma }\quad 平衡变形阶段\\ 无声发射现象产生 \end{array} \right.\\ 或\left\{ \begin{array}{l} \hat \sigma \le \tilde \sigma \le \widehat {\hat \sigma }\quad 平衡变形阶段\\ 一定量声发射现象产生 \end{array} \right. \end{array} $

(28)

3）加速变形稳定模式.该模式下整个监测周期内顶板位移监测结果分析满足式（26），说明顶板实际已由平稳变形阶段过渡到加速变形阶段，由于顶板岩体的加速变形，导致顶板岩体内部产生一定数量的新裂纹，并随着变形进一步增加，岩体内部裂纹逐步扩展，声发射实际监测结果也体现为其能量和频度等特征参量激增，但此时只是加速变形阶段岩体内部结构在破坏，而顶板依然可以继续承载，尚能保持稳定，组合模式见式（29）.

$ \left\{ \begin{array}{l} \tilde \sigma > \hat \sigma \quad 平衡变形阶段\\ 声发射特征参量分析结果：无激增-回落现象 \end{array} \right. $

(29)

4）加速变形失稳模式.该模式下整个监测周期内顶板位移监测结果分析满足式（26），顶板处于加速变形阶段，同时声发射特征参量分析结果出现激增-回落现象.说明此时顶板岩体大部分已经破坏，内部裂纹已经出现贯通融合，致使声发射特征参量频度（量度）回落.该组合分析模式表明顶板岩体已进入整体失稳阶段，见式（30），顶板处于失稳前兆期.

$ \left\{ \begin{array}{l} \tilde \sigma > \hat \sigma \quad 加速变形阶段\\ 声发射特征参量分析结果：激增-回落 \end{array} \right. $

(30)

4 应用实例

某铜矿对阶段水平矿柱（顶底柱）实施回采^[15]，水平矿柱上方为上部中段松散块石充填体，回采过程主要采用分层进路式上采方式，在上部散体压力作用下，进路顶板的稳定性是安全回采的关键.为此，对进路顶板实施位移与声发射组合监测，柱17^#采场7^#进路回采过后在进路中布置顶板沉降仪监测顶板位移，同时利用水平矿柱上部中段坑道在顶板上布置钻孔，安装声发射连续监测仪，测点布置见图 4.

从仪器安装至进路充填共持续监测41 d，取得位移监测数据15组，声发射监测处理数据40组.见表 1、表 2.

利用第1节中的处理方法对位移监测值进行分析，由于D₂未能取得完整监测结果，分别计算D₁、D₃测点位移监测值的标准差.

$ \tilde \sigma ({D_1}) = \sqrt {\frac{1}{{n-1}}\sum\limits_{i = 1}^n {{{({w_i}-\bar w)}^2}} } = 4.6639 $

(31)

$ \tilde \sigma ({D_3}) = \sqrt {\frac{1}{{n-1}}\sum\limits_{i = 1}^n {{{({w_i}-\bar w)}^2}} } = 0.3459 $

(32)

分别计算出D₁、D₃测点无变形阶段和稳定变形阶段的标准差上限值，测量仪器的精度σ为0.01.

$ \begin{array}{l} \hat \sigma ({D_1}) = \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} = 0.038\\ \hat \sigma ({D_3}) = \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} = 0.038 \end{array} $

(33)

$ \begin{array}{l} \widehat {\hat \sigma }({D_1}) = S\sqrt {\frac{{n + 1}}{{12n}}} + \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} = 4.639\\ \widehat {\hat \sigma }({D_3}) = S\sqrt {\frac{{n + 1}}{{12n}}} + \sigma \sqrt {\frac{{{n^2}}}{{n-1}}} = 0.582 \end{array} $

(34)

根据计算结果可知：

$ \tilde \sigma ({D_1}) > \widehat {\hat \sigma }({D_1})\;\;\;\hat \sigma ({D_3}) < \tilde \sigma ({D_3}) < \widehat {\hat \sigma }({D_3}) $

(35)

由此可知D₁测点所在区域顶板已经处于加速变形阶段，而D₃测点所在区域顶板尚处于稳定变形阶段.顶板是否出现失稳前兆应配合声发射测点的监测结果进行组合分析.

对声发射能率监测结果绘制变化曲线，见图 5.

从图 5可知，在整个监测周期内，声发射能率先后出现了3次较大的起伏，其中尤以第3次出现了声发射的激增-回落，配合该监测期内顶板位移监测结果，发现截止至12月23日，顶板D₁测点区域已经处于加速变形阶段，符合位移与声发射组合分析的加速变形不稳定模式，由此判断D₁测点附近区域顶板已经处于失稳前兆期，遂撤出该进路的设备及监测仪器，12月27日早班发现该进路端头与上盘岩体接触带附近顶板岩体小规模垮塌，27日小夜班对该进路实施充填.

5 结论

1）顶板位移监测数据数理统计分析结果表明：地下工程开挖后，顶板位移主要经历3种变形过程，即无变形过程、稳定变形过程与加速变形过程.在整个地下工程的复杂施工过程中，3个过程可以相互过渡转换.

2）岩石承载过程多次声发射特性试验结果说明，岩体失稳前，其声发射监测指标频度（量度）出现激增-回落现象.但试验全过程显示该现象出现并不能预测岩体必然失稳.

3）仅凭位移或声发射单一监测技术并不能对顶板岩体稳定性做出准确预判，对2种测试结果进行组合分析得到4种组合模式，其中加速变形失稳模式是顶板岩体失稳的真正前兆.