鄂尔多斯盆地侏罗纪地层是我国重要的含煤地层，其上覆岩层多为弱胶结地层，广泛发育富水性强或较强的基岩孔隙−裂隙含水层和第四系松散含水层。大量开采实践证明，厚煤层大跨度工作面开采后对顶板覆岩扰动强度大，弱胶结覆岩破坏严重，工作面回采产生的导水通道易与上覆含水层沟通，造成工作面顶板突水事故的发生。因而，对黄土梁峁地区砂岩−泥岩互层覆岩突水类型的有效判别，直接决定了现场采取的防治水方向与措施。

近年来，有学者先后提出了多元统计分析数量化理论^[1]、二级模糊综合评判数学模型^[2]、ESN正则化模型^[3]、Fisher判别法^[4-5]、改进的集对分析模型^[6]、混沌麻雀搜索RF模型^[7]、自适应粒子群优化极限学习机模型^[8]、贝叶斯判别理论^[9-11]、灰色关联度理论^[12]等水源判别方法。这些研究成果对于指导煤矿安全生产，防治突水事故具有重要的指导意义^[13]。但上述成果主要针对石炭−二叠纪煤层开采时坚硬–较坚硬覆岩突水水源类型的定性识别，且多选用常规水化学指标作为模型的输入变量，在一定程度上忽略了指标间的相关性。而未确知信息作为一种新的不确定性信息，不同于模糊信息、随机信息和灰色信息，是一种将诸多因素不确定信息定量化评价的数学方法^[14]，相比于其他评价方法，未确知信息满足“归一性”和“可加性”原则，判别结果更加科学可信^[15]。基于此，笔者选取了7组常规、同位素特征比因子，在主客观赋权的基础上采用测度函数判别突水水源数据，推导不同疑似水源所占比例，以期为黄土梁峁地区矿井顶板水害防治提供理论依据。

1   研究区概况

朱家峁煤矿位于陕北侏罗纪煤田榆横矿区东北部，毛乌素沙漠东南缘和陕北黄土高原接壤地带，地貌为典型的黄土梁峁区，地表被第四系松散沉积物覆盖，断层不发育，褶皱较发育。钻孔揭露的地层由老至新是三叠系上统瓦窑堡组，侏罗系下统富县组、中统延安组和直罗组，第四系(图1)。侏罗系中统延安组是主要含煤地层，目前开采延安组第三段顶部的3号煤层，是唯一可采煤层。研究区处于黄土梁峁地貌区，第四系松散层不含水或赋水性极弱，所以地表水及顶板砂岩水是构成该煤层开采的主要充水水源。

图  1  研究区工作面剖面

Figure  1.  profile of the working face in the Zhujiamao Coal Mine

下载: 全尺寸图片 幻灯片

研究区地势南高北低，地形支离破碎，沟壑纵横。区内水系发育，距离朱家峁煤矿最近的无定河，位于煤矿北边界，由西向东流过，年平均流量18.7 m³/s；煤矿3盘区中部的河流沙空沟，丰水期流量60~100 m³/h，枯水期流量10~40 m³/h，由南向北汇入无定河；煤矿西北部的大川沟水库，库容量约187万m³，库底所压3号煤层底板高程约895 m，距库底115 m，煤层上部未有较好隔水层存在。

2   材料与方法

2.1   样品采集与测试

样品采集：为了判别朱家峁煤矿3号煤层开采时的涌水来源，2020年7—11月开展了流域地表水、矿井水采样。共采集26组水样，其中顶板砂岩水样5组，沙空沟水样6组，大川沟水库水样5组，无定河水样5组，突水水样5组。

测试方法：采用便携式多参数仪测试pH、总溶解性固体(TDS)，采用滴定法测试HCO₃ ⁻，采用电感耦合等离子体原子发射光谱法测试K⁺、Na⁺、Ca²⁺、Mg²⁺，采用离子色谱法测试Cl⁻、SO₄ ²⁻、NO₃ ⁻、F⁻，采用Finnigan Delta V Advantage同位素比率质谱仪和Flash 2000 HT元素分析仪测试δ²H、δ¹⁸O。

2.2   分析方法

2.2.1   T-球型模糊聚集算子TOPSIS理论

令$T = \left( {{t_1},{t_2}, \cdots ,{t_k}} \right)$为决策专家集合，专家${t_k}$针对每一个方案$ {A_i} = \left( {{A_1},{A_2}, \cdots ,{A_q}} \right) $在属性${C_j} = \left( {{C_1},{C_2}, \cdots ,{C_n}} \right)$的评分为$ r_{ij}^o $，构成决策矩阵${{\boldsymbol{R}}^k} = \left( {{\boldsymbol{r}}_{ij}^k} \right)\left( i = 1, 2, \cdots , q;j = 1, 2, \cdots ,n \right)$。假设${P_k} = \left( {s_k},{l_k},{t_k} \right)$为T-球型模糊数集，则T-球型模糊数加权平均聚集算子为^[16-17]：

式中：$w_k $为下球型模糊数集合$P_k $的权重向量。

采用$ f_{ij}^ + $、$ f_{ij}^ - $分别表示第j个属性C_j下专家信息f_ij的正理想解和负理想解，计算n个属性的群体最优效益值$ w_j^ + $、遗憾值$ w_j^ - $和妥协值W_j^[16-17]：

妥协值W_j反映了第j个属性的相对重要性，W_j越小，则该属性越重要。

2.2.2   粗糙集理论

粗糙集理论利用上近似集合和下近似集合概念描述粗糙集。给定一个知识表达系统$S{\text{=}}\left( {U,C \displaystyle\cup D,V,f} \right)$，其中，$ U{\text{ = }}\left( {{x_1},{x_2}, \cdots ,{x_n}} \right) $是一个有限的非空对象集合；C和D分别为非空有限的条件属性集和决策属性集，交集为空；V为属性值域；f为信息函数。令X为U的任意非空子集，则存在一个不可辨识关系ind(B)和非空集合Y，使得集合X的B上近似集和下近似集分别表示为^[18-19]：

X的B边界域可用以下3种方式表示：记为：$ bnd\left( X \right) = \bar BX - \underline B X $，$ pos\left( X \right) = \underline B X $，$neg\left( X \right) = U - \bar BX$。显然可知，其上近似集是正域和边界域的并集，即$ \bar BX = pos\left( X \right) \cup bnd\left( X \right) $。

对于任意属性子集$ B \subseteq C \cup D $决定了一个二元等价关系$ind\left( B \right) = \left\{ {\left( {x,y} \right) \in U \times U\left| {\forall a \in B,f\left( {x,a} \right) = f\left( {y,a} \right)} \right.} \right\}$，则属性a相对于B和D重要度定义为：

2.2.3   D-S证据理论

D-S证据理论通过Dempster组合规则实现多个证据的融合，采用m₁和m₂表示证据融合的基本概率指派，则Dempster组合规则定义为^[20]：

式中：Q_j为T-球型模糊聚集算子TOPSIS理论求解的指标权重；R_j为粗糙集理论求解的指标权重；g为冲突系数。

2.2.4   单指标未知测度函数

单指标未确知测度函数是指指标按照未确知测度定义构造的单指标函数${\textit{z}}\left( x \right)\left( {{x_{ij}} \in {L_k}} \right)\left( i = 1,2, \cdots ,n; j = 1, 2, \cdots ,m;k = 1,2, \cdots ,t \right)$，则第i个评价因素x_i对应的第j个评价指标的观测值x_ij构成的未确知测度矩阵：

定义水源识别空间为$ L{\text{ = }}\left\{ {{L_1},{L_2}, \cdot \cdot \cdot ,{L_t}} \right\} $，d_k和d_k+1分别为识别空间中L_k量值范围的阈值，则单指标未知测度函数^[16]：

假设t元联系数µ的同一、差异、对立度的隶属度矩阵${{\boldsymbol{M}}_i} = \left( {{a_i},{b_{i1}},{b_{i2}}, \cdots ,{b_{i\left( {t - 2} \right)}},{c_i}} \right)$，系数分量${{\boldsymbol{E}}}_{i}= {\left(1\;{p}_{i1}\;{p}_{i2}\;\cdots \;p{}_{i\left(t-2\right)}\;{q}_{i}\right)}^{{\rm{T}}}$，则第i个水样隶属于各潜在突水水源的混合比例为^[16]：

3   结果与讨论

3.1   水化学特征

研究区所有水样呈弱碱性。与沙空沟、大川沟水库、无定河水样相比，顶板砂岩水样阳离子以Na⁺(77.26%，为毫克当量百分数，以下同)为主，阴离子以SO₄ ²⁻(93.08%)为主，高于地表水Na⁺(57.89%)、SO₄ ²⁻(37.68%)占比。推测地表水Na⁺、SO₄ ²⁻可能来源于盐岩的蒸发和溶解，在地表水向下渗入过程中不断溶解钠长石和硫酸盐，致使砂岩水Na⁺、SO₄ ²⁻显著增加。Piper三线图将水样点划分为Ca·Mg-SO₄·Cl和Na-SO₄·Cl，以Na-SO₄·Cl为主(图2)。在地表水水样点中，无定河水化学类型为Ca·Mg-SO₄·Cl和Na-SO₄·Cl混合，不同于沙空沟和大川沟水库的水−岩作用机理，其在形成过程溶解了更多的Ca²⁺和Mg²⁺。地表水水化学类型以Na-SO₄·Cl占优势，地下水(砂岩水样和突水水样)水化学类型均为Na-SO₄·Cl，表明地表水与地下水间存在一定的水力联系。

图  2  Piper三线图

Figure  2.  Piper diagram

下载: 全尺寸图片 幻灯片

3.2   水化学控制机制

如图3所示，沙空沟、大川沟水库水样位于以岩石风化为主的区域，无定河水样位于岩石风化与蒸发结晶交界区，砂岩水、突水水样位于蒸发结晶控制区。可知，地表水水化学成分受岩石风化和蒸发结晶共同控制，而地下水水化学成分主要受控于蒸发结晶作用，且随高程降低，地下水受蒸发结晶作用控制越明显。

图  3  Gibbs图

Figure  3.  Gibbs diagram

下载: 全尺寸图片 幻灯片

如图4所示，地表水主要离子均与TDS显著正相关，其中Ca²⁺、Mg²⁺、HCO₃ ⁻占主导，在p≤0.001水平下与TDS相关系数分别为0.96、0.88、0.80。地下水TDS与Na⁺呈显著正相关关系(p≤0.05)，与K⁺、Ca²⁺(p≤0.01)、Mg²⁺、HCO₃ ⁻(p≤0.05)呈负相关关系。表明地表水离子来源与地下水存在差异，地下水中发生了强烈的离子交替吸附作用。

图  4  水化学指标相关性

Figure  4.  Correlation among main hydrochemical indexes in surface water and groundwater

下载: 全尺寸图片 幻灯片

3.3   常规水化学离子比特征

如图5a所示，地表水和地下水水样点全部落在γ(K⁺+Na⁺)∶γCl⁻=1∶1(γ为离子毫克当量浓度)上方，但地表水水样点更接近于1∶1线，表明地表水Na⁺主要来自硅酸盐岩的风化。砂岩水样点分布在c[(HCO₃ ⁻+SO₄ ²⁻)−(Ca²⁺+Mg²⁺)]∶c(Na⁺+K⁺−Cl⁻)=−1∶1线附近，表明砂岩水发生了不同程度的阳离子交换作用(图5b)。地表水呈弱碱性，地下水呈弱酸性，弱酸性环境易于NO₃ ⁻的富集(图5c)。由图5d知，地表水c(Ca²⁺+Mg²⁺)∶c(HCO₃ ⁻)(c为浓度)平均比为0.79，c(Ca²⁺+Mg²⁺)∶c(HCO₃ ⁻+NO₃ ⁻)为0.78，c(Ca²⁺+Mg²⁺)∶c(HCO₃ ⁻+SO₄ ²⁻)为0.50，表明加入SO₄ ²⁻，地表水阴离子过剩状态得以改善，这个结果解释了衍生硝酸盐参与矿物风化的可能性。地下水c(Ca²⁺+Mg²⁺)∶c(HCO₃ ⁻)平均比为3.52，c(Ca²⁺+Mg²⁺)∶c(HCO₃ ⁻+NO₃ ⁻)为3.36，c(Ca²⁺+Mg²⁺)∶c(HCO₃ ⁻+SO₄ ²⁻)为0.14，表明加入SO₄ ²⁻，地下水由阳离子过剩状态向阴离子过剩状态发展。如图5e、图5f所示，地表水氟化物处于非饱和状态，而砂岩水中氟浓度基本达到了上限值，氟丰度受含氟矿物(萤石)和含钙矿物(方解石)共同控制。位于趋势线2右下角的部分突水水样，萤石与方解石质量比小于1∶200，表明含氟矿物中的氟通过淋溶作用进入地下水。

图  5  常规水化学特征离子比

Figure  5.  Characteristic ratios between different ions

下载: 全尺寸图片 幻灯片

3.4   同位素空间变异性

同位素空间变异函数如图6所示，模型拟合参数见表1。块金值(C₀)与基台值(C₀+C)的比值反映了变量的空间自相关性。若[C₀/(C₀+C)]≤25%，相关性显著，结构因素占主导；若25%<[C₀/(C₀+C)]<75%，相关性中等，结构因素和随机因素共同控制；若[C₀/(C₀+C)]>75%，相关性弱，随机因素占主导。

图  6  稳定同位素空间变异函数

Figure  6.  The spatial variograms of stable isotopes in surfer water samples and groundwater samples

下载: 全尺寸图片 幻灯片

表  1  模型拟合参数

Table  1.  Fitting models for stable isotopes and their optimal parameters

水样点	稳定同位素
	δ2H/‰					δ18O/‰
	拟合模型	块金值(C0)	基台值(C0+C)	[C0/(C0+C)]		拟合模型	块金值(C0)	基台值(C0+C)	[C0/(C0+C)]
沙空沟	球型	0.010	5.450	0.002		球型	0.023	0.186	0.124
大川沟水库	高斯	0.010	6.150	0.002		高斯	0.0001	0.111	0.001
无定河	线性	0.759	0.759	1.000		高斯	0.012	0.113	0.106
砂岩水	高斯	0.950	6.909	0.138		高斯	0.036	0.740	0.049
突水水样	线性	0.100	39.810	0.003		高斯	0.104	0.950	0.109

下载: 导出CSV 

| 显示表格

突水水样的同位素[C₀/(C₀+C)]值介于0~0.25，空间自相关性强，表明突水水源δ²H、δ¹⁸O的空间变异性受结构因素控制，且δ²H在高程较大区域表现出更大的负值，在高程较小区域表现出相对小的负值，表明突水水样δ²H存在较大的空间变异性，δ¹⁸O空间变异性较小。

3.5   同位素指示特征

当地雨水线性方程(LMWL)与全球雨水线性方程(GMWL)如下^[21]：

为评估朱家峁煤矿矿井水来源，了解蒸发对同位素组成的影响，基于δ¹⁸O与δ²H相关关系，构建地表水水线和地下水水线。

地表水、地下水水样点均位于LMWL和GMWL线之下，但地下水水样点比地表水水样点更接近LMWL(GMWL)，地表水水线的斜率和截距小于地下水水线，表明陕北地区地表水经历了强烈的蒸发作用，地下水经历了蒸发作用或受到了蒸发水的混合。突水水样水线的斜率与截距接近GMWL，表明突水水样的同位素值与降水同位素值相似，较好地保存了降水同位素特征(图7)。

图  7  氢氧稳定同位素拟合

Figure  7.  The correlation between δ¹⁸O and δ²H for surface water samples and groundwater samples

下载: 全尺寸图片 幻灯片

3.6   基于T-TOPSIS-RST-DSET-SIUMF的混合突水水源判别模型

将T-球型模糊聚集算子TOPSIS方法(T-TOPSIS)、粗糙集理论(RST)、D-S证据理论(DSET)与单指标未知测度函数(SIUMF)相结合，对突水水源混合比例进行判别，具体方案如下：

(1) 确定混合突水水源判别模型的特征指标和潜在突水水源类别

特征指标X={X₁，X₂，···，X₇}={γ(Na⁺+K⁺)/γCl⁻，γ(Na⁺−Cl⁻)/γ[(Ca²⁺+Mg²⁺)−(HCO₃ ⁻+SO₄ ²⁻)]，γ(NO₃ ⁻+Cl⁻)/γHCO₃ ⁻，c(Ca²⁺+Mg²⁺)/c(HCO₃ ⁻+NO₃ ⁻)，SIF/SICa，{lg[F⁻]活性}/{lg[Ca²⁺]活性}，δ²H/δ¹⁸O}，潜在突水水源类别Y={顶板砂岩水(A)，沙空沟(B)，大川沟水库(C)，无定河(D)，AB，AC，AD，BC，BD，CD}，2个字母表示2种水源混合。

(2) 采用T-球型模糊聚集算子TOPSIS方法求解特征指标主观权重

基于专家评分信息，构建T-SFWAA算子矩阵，计算正理想解f_ij⁺和负理想解f_ij⁻：

依据式(2)和式(3)，计算特征指标群体最优效益值和遗憾值，见表2。

表  2  $w^-_j $和$w^+_j$计算结果

Table  2.  The results of $w^-_j $ and $w^+_j $

指标	X1	X2	X3	X4	X5	X6	X7
wj−	0.310	0.400	0.286	0.476	1.323	1.000	0.027
wj+	0.240	0.512	0.526	0.696	1.000	1.266	0.025

下载: 导出CSV 

| 显示表格

依据式(4)，计算妥协值，归一化处理后，得到特征指标主观权重，见表3。

表  3  特征指标主、客观权重

Table  3.  Subjective and objective weights of characteristic indexes

判别指标	X1	X2	X3	X4	X5	X6	X7
主观权重(Ws)	0.171	0.133	0.107	0.123	0.173	0.134	0.158
客观权重(Wo)	0.111	0.333	0.111	0.111	0.111	0.111	0.112

下载: 导出CSV 

| 显示表格

(3) 采用粗糙集理论求解特征指标客观权重

依据特征指标空间阈值，构建各类潜在突水水源隶属于X₁, X₂, ···, X₇的属性函数：

确定论域U关于指标X₁, X₂, ···, X₇的等价划分集合和正域：

计算指标X₁, X₂, ···, X₇的相对重要度，归一化处理后，得到特征指标客观权重见表3。

(4) 基于D-S证据理论的特征指标赋权耦合

采用Dempster证据合成规则计算指标X₁, X₂, ···, X₇的组合权重：

(5) 单指标未知测度函数判别潜在突水水源混合比例

单指标未知测度模型的系统层包括：γ(Na⁺+K⁺)/γCl⁻、γ(Na⁺−Cl⁻)/γ[(Ca²⁺+Mg²⁺)−(HCO₃ ⁻+SO₄ ²⁻)]、γ(NO₃ ⁻+Cl⁻)/γHCO₃ ⁻、c(Ca²⁺+Mg²⁺)/c(HCO₃ ⁻+NO₃ ⁻)、SIF/SICa、{lg[F⁻]活性}/{lg[Ca²⁺]活性}、δ²H/δ¹⁸O，特征指标未知测度函数如图8所示。

图  8  单指标未知测度函数

Figure  8.  Single-index unascetained functions

下载: 全尺寸图片 幻灯片

将突水水样特征指标值代入相应的测度函数进行计算，在考虑同一性、差异性和对立关系的基础上，判别突水水样混合比例，结果见表4。

表  4  潜在突水水源混合比例判别结果

Table  4.  Identification results of mixing proportions for various water inrush sources

突水水源	不同水源占比/%
	顶板砂岩水	沙空沟	大川沟水库	无定河
1	28.96	50.75	15.17	5.12
2	32.30	44.34	17.84	5.52
3	30.49	49.82	15.31	4.37
4	37.72	14.24	45.71	2.32
5	24.10	20.96	47.41	7.52

下载: 导出CSV 

| 显示表格

3.7   三维高密度电法探测

由上述判别可知，朱家峁煤矿开采3号煤层，突水来源主要为顶板砂岩水和沙空沟水，为证明模型判别准确性，采用三维高密度电法技术对沙空沟进行探测(图9)。低电阻率特征表示裂隙发育，而较完整的岩体结构则表现出高电阻率特征。依据三维高密度电法探测成果显示，Zone B是斑裂线发育区(图10)，Zone A是未开采区域。发现B区存在明显的低阻异常，3302工作面采空区存在明显的低阻异常，低阻异常自地下130 m深度发育至地表，表明沙空沟水由斑裂线通过Zone B进入采空区。

图  9  沙空沟三维高密度电法探测

Figure  9.  Results of three dimensional high density electrical prospecting in Shakonggou riverbed

下载: 全尺寸图片 幻灯片

图  10  朱家峁煤矿开采形成的地表斑裂线

Figure  10.  Surface zebra linear cracks formed by mining of Zhujiamao Coal Mine

下载: 全尺寸图片 幻灯片

4   结 论

a. 地表水水化学类型以Na-SO₄·Cl占优势，地下水(砂岩水样和突水水样)水化学类型为Na-SO₄·Cl，表明地表水与地下水间存在一定的水力联系。

b. 地表水水化学成分受岩石风化和蒸发结晶共同控制，而地下水水化学成分主要受控于蒸发结晶作用，且随高程降低，地下水受蒸发结晶作用控制越来越明显。

c. 突水水样的同位素空间自相关性强，其空间变异性主要受结构因素控制，且δ²H存在较大的空间变异性，δ¹⁸O空间变异性较小。

d. T-TOPSIS-RST-DSET-SIUMF混合突水水源判别模型表明，朱家峁煤矿开采3煤层开采突水供水水源包括顶板砂岩水、沙空沟水、大川沟水库及无定河，但主要为顶板砂岩水和沙空沟水，两者之和所占比例大于50%，其中沙空沟水通过斑裂线进入采场。

e. 研究表明，该方法不仅克服了传统的Piper三线图难以识别水化学类型(Na-SO₄·Cl)相同但水源不同的缺点，而且能够揭示不同含水系统水化学演化规律和水源补给关系，以及其受控因素的差别。