岩溶区排泥库泥浆深层泄漏机理

欧孝夺，秦金喜，全守岳，江 杰*，杨杉楠

(1.广西大学土木建筑工程学院，南宁 530004；2.广西大学工程防灾与结构安全教育部重点实验室，南宁 53004；3.广西金属尾矿安全防控工程技术研究中心，南宁 530004；4.中铁二院工程集团有限责任公司南宁勘察设计研究院，南宁 530012)

随着中国经济快速增长，金属工业产品的需求逐年上升，而开采量的增加导致排泥库数量和库容大幅增加。排泥库作为金属冶炼废弃物堆存设施，是金属非金属矿山的重大危险源。特别在岩溶地区修建排泥库时，地下溶洞交错连接，常因底部溶洞未探查清楚，致使排泥库在长周期运行后或遇暴雨、地震等其他自然因素作用下发生泄漏事故。排泥库泥浆富含重金属离子，一旦发生泄漏事故，库内的废渣、废水就会污染地下水，甚至淹没下游的农田村庄，不仅破坏生态环境，同时给人民的生命财产带来巨大的安全隐患[1-3]。据宋园生[4]统计，溶洞塌陷占造成排泥库泥浆泄漏事故直接原因的48.1%，因而开展岩溶区排泥库深层泥浆的泄漏机理研究非常有必要，为后期尾矿库防治工作提供理论基础。

为研究流体的泄漏机理，国内外学者进行了大量有益的研究。Ger-manopoulos等[5]、Lin[6]、Alonso等[7]和中国城镇供水协会[8]分别建立了不同情况下的管网泄漏量计算模型，在一定条件下为不同区域的流体泄漏提供了计算方法，但并未针对岩溶区尾矿库区泥浆流体的泄漏提供相应的计算模型。王勇智[9]认为泥石流这类流体的浓度有一定的变化范围，其相对黏滞系数与固体体积浓度、颗粒大小均有关系。徐黎明等[10]研究得出影响泥石流流速的因素包括泥沙粒径、泥浆深度、泥浆密度等，并利用此参数建立BP(back propagation)神经网络。费祥俊等[11]研究了泥石流中不同颗粒大小的运动轨迹等参数，证明固液两相流的运动形式受颗粒大小影响，并确定了两相流的分界粒径。沈寿长[12]研究了应该采取何种计算模型才能准确反映泥石流固液相的相互作用，并提出了相应的本构关系方程。泥石流等流体相关性质虽同尾矿泥浆相似，但目前国内外学者并未针对尾矿泥浆深层泄漏机理进行相应研究。宋园生[4]分析了处于岩溶区的排泥库面临溶洞塌陷时采取的险情处置方式，给出了岩溶地区排泥库在建设和运行过程中的建议并得到相应的应用。但现目前，国内外针对排泥库的泥浆的探测研究鲜有报道，而针对大坝、水库等处于岩溶地区中的设施泥浆泄漏通常采用地球物理法[13]、综合示踪法[14-16]等进行检测，虽方法效果显著但较为复杂。对于流体泄漏，学者主要集中研究泥石流溃坝及水库泄漏或仅针对喀斯特地区尾矿库泄漏所存在的灾情处置手段等方面，而针对岩溶区排泥库底部泥浆泄漏规律的研究仍处于空白。

为解决岩溶区排泥库底部泄漏的快速探测问题，自主研发了一种探测装置[17]，利用该探测装置所测的试验数据运用灰色关联建立了泥浆参数计算模型，并用ANSYS进行泄漏孔直径、泥浆含水量、动力黏度对泄漏点的影响半径数值模拟研究[18]。该研究创新性提出新型排泥库泄漏探测方式，建立了泥浆对探测器拉力的多元线性回归计算模型，对快速准确查找排泥库泄漏点位置有极大帮助。通过数值模拟计算了深层泥浆的不同物理性质对影响半径的影响作用。本文对于排泥库深层泄漏点快速探测及明确泄漏区域泥浆流动状态都具有重要的理论指导意义。

1 试验设计

试验主要测量泥浆探测装置所受渗流力与泥浆流速、泥浆含水量和泥浆密度之间的关系，以期能寻找规律并建立相关理论。

1.1 试验装置

试验装置主要由测量装置、流动装置、探测装置、储存装置等组成，如图1所示。

图1 试验装置模型图

1.1.1 泥浆渗流力测量装置

泥浆渗流力采用WH-A04L型数显拉力计进行测量，测量精度为0.001 kg。数显拉力计下端连接自制的泥浆探测器。泥浆探测器采用聚氯乙烯塑料块加工制成，为增大泥浆探测器和流动泥浆的摩擦阻力，使泥浆对探测器的拉力作用更加敏感，将探测器做成扇贝形状，并在探测器表面进行粗糙条纹处理，探测器的内部进行机械掏空，可通过填充不同重量的填充物改变其密度，使之沉到不同的泥浆层中，探测器如图2所示。试验中为使泥浆流速计算误差减小，使之更接近出口处泥浆流速，探测器放置在聚氯乙烯(PVC)管底部靠近孔口处。

图2 可调重量探测器

1.1.2 泥浆流动装置

泥浆流动装置由内径150 mm，长度1 000 mm的PVC白色塑料管组成。将塑料管竖直放置，在管底使用塑料盖进行封底。为得到不同的泥浆流速，对盖子进行不同直径大小的钻孔，孔口直径分别设为10、15、20、25、30、35 mm。因排泥库泥浆随着深度的增加，含水量有很大区别，研究表明排泥库泥浆静置抽干上部清水后的表层泥浆含水量基本维持在120%，排泥库发生泄漏的泥浆一般为新排入排泥库的浑浊泥浆，由于远距离泵送功量需要，吹填入库的泥浆需保证一定的流动性，使所吹填泥浆含水量比静置的泥浆要高很多，根据此情况对泥浆进行分组，分别配置含水量为100%、110%、120%、130%、140%、150%、160%、170%、180%的泥浆。

1.2 试验步骤

测量泥浆泄漏渗流力的试验步骤如下。

(1)采用常规土工密度实验对不同含水量的泥浆进行密度测量，得到相应的数值。试验前，将管底泄漏孔口预先用东西封住，往PVC管倒满一定含水量的泥浆。

(2)打开底部泄漏孔，同时开始计时并进行上部注浆，保持泥浆水头不变，进而控制底下流速保持一致。由于流速较快，故只记录前面5 s的流量，5 s后停止泄漏，将泄漏出来的泥浆进行称重，记录5 s钟泄漏的流量和上面的拉力读数。

(3)重复上面的试验，更换不同含水量的泥浆和不同大小的泄漏孔口，得到正交数据。

2 试验结果分析

2.1 关联度分析

为了解泥浆的各物理力学参数与泥浆渗流力探测装置的关联程度，通过对多种方法进行比较分析，选择灰色关联分析法从泥浆的各物理力学参数中筛选出与探测装置所受泥浆渗流力关联度较高的影响因子作为定量因子。

通过以上室内土工实验和泄漏模型试验，获取了大量的排泥库泥浆物理力学性质数据和泥浆泄漏流动数据，共计54组原始数据。

原始数据处理具体分析步骤如下。

2.1.1 流速换算

由于试验所测数据为泥浆流量，需要将其换算为探测器位置的泥浆流速。换算公式为

(1)

式(1)中：v为泥浆流速，cm/s；m为泄漏出来的泥浆质量，g；ρ为泥浆密度，g/cm3；t为时间，此处取t=5 s；A为内径150 mm的PVC圆管截面面积。

2.1.2 确定参考数列和比较数列

能反映系统行为特征的数据序列称为参考数列X0，本文参考数列需要代表泥浆探测装置的受力指标，此处选择泥浆探测器所受渗流力大小作为参考数列X0。影响系统行为的因素组成的数据序列，称比较数列，此处选择泥浆流动参数和泥浆物理性质指标作为比较数列Xi，分别为泥浆的流速v(X1)、天然含水量ω(X2)、天然密度ρ(X3)。

2.1.3 原始数据无量纲化

由于系统中各因素的物理意义不同，导致数据的量纲不相同，不便于比较计算，难以得出正确的结论。因此在进行灰色关联度分析时，需要对数据进行无量纲化处理后才可进行关联性分析。常用的无量纲化处理方法有三种，分别为初值化、均值化和区间值化处理。经过分析，本次试验数据没有“0”，且含水量是从100%开始计算，采用初值法进行无量纲化处理，可以使计算简单、不易出错，故对排泥库泄漏试验所测原始数据使用初值法进行无量纲化处理。

2.1.4 对序列求差处理

计算参考数列与各比较序列差值的绝对值，计算公式为

Δi(k)=|X0(K)-Xi(K)|

(2)

求出差值绝对值的最大值和最小值，计算公式如式(3)、式(4)。

max[maxΔi(k)]=1.52

(3)

min[minΔi(k)]=0.00

(4)

2.1.5 求关联系数与关联度

用式(5)求参考数列与比较数列的关联系数。

(5)

因为关联系数是比较数列与参考数列在各个时刻的关联程度值，所以它的数不止一个，而信息过于分散不便于进行整体性比较。因此有必要将关联系数集中为一个值，即求其平均值，作为比较数列与参考数列间关联程度的数量表示。

通过上述数据处理手段进行试验原始数据的处理得到如图3所示的各因素与渗流力关联系数。

图3 各因素与渗流力关联度

用式(6)求解关联度ri：

(6)

式(6)中：ri为X0与Xi的关联度，其值越接近1，则关联度越好。各指标与泥浆所受渗流力关联度分别为r1=0.818 1、r2=0.705 2、r3=0.846 7，由此可知各关联度大小顺序为r3>r1>r2。

由灰色关联分析结果可知，排泥库泥浆含水量、密度和流速三者对泥浆渗流力的关联系数均较高，其中排泥库泥浆含水量与泥浆探测器所受渗流力关联度较低，而泥浆密度和泥浆探测器所受渗流力关联度较高，相互之间的关系较为密切，均可作为计算模型考虑的因素。

2.2 泥浆泄漏参数多元线性回归计算模型

通过运用多元线性回归理论对前期试验数据进行分析，计算时取系数个数p=3，样本数量取n=54，根据数值分析相关回归计算公式计算得到：

y=-18.535 5-0.017 3x1+0.013 6x2+

10.802 3x3

(7)

式(7)中：y为泥浆对探测器的渗流力大小，N；x1为泥浆流速，：cm/s；x2为泥浆含水量，%；x3为泥浆密度，g/cm3。

这就是排泥库泥浆深层泄漏时泥浆对探测器的渗流力与泥浆内部参数的计算模型。

对此计算模型进行显著性检验。运用F分布的方差分析法计算本次数据的F=740.667，查F分布上侧分位数表知F0.05(3，50)<2.84，显而易见F>F0.05(3，50)，证明此次计算结果线性回归显著，结论可靠。

3 数值模拟分析

在建立泥浆探测器渗流力与泥浆参数的计算模型后，为能快速准确查到泄漏点位置，还需确定泄漏点的影响半径。

3.1 模型建立

3.1.1 基本假定

排泥库泥浆由于在重力作用下发生沉积，其含水量随深度发生变化并出现分层现象，模拟泥浆分层较为复杂，为简化计算，缩短计算时间，故假设模拟过程泥浆均为含水量单一的介质。且运用VOF模型时需控制体积必须充满单一流体相或存在多相的联合但其有一相作为可压缩的理想气体，故本模型上部排泥库库面以上假设为可压缩的理想气相。壁面处为非滑移边界条件，热流量假设为零热流量，即壁面假设为绝热壁面、光滑壁面。

3.1.2 软件及模型选取

本次数值模拟使用ANSYS ICEM CFD作为前处理器，Fluent 15.0作为求解器和后处理软件。该模型模拟的两相是空气和泥浆，选用Fluent中的VOF模型，故只能绝对速度的压力基瞬态求解器。模拟泥浆旋转效应不是特别明显，采用常用的标准k-ε湍流模型。由于讨论的泥浆流动是瞬态的，采用更适宜瞬态计算的PISO算法，收敛标准取10-4，有效提高收敛性和缩短计算时间。

3.1.3 建模尺寸

通过对广西百色地区尾矿库区调研发现，为贴近现场实际情况，更真实模拟排泥库泥浆泄漏情况，模型模拟尺寸为50 m×50 m×30 m。整个模型模拟的是一个底部泄漏的排泥库，模型上部模拟的是排泥库库面，模型侧面和底面模拟的是排泥库壁面，圆孔模拟的是排泥库泥浆泄漏孔，几何模型示意如图1所示。

3.1.4 关键参数及边界条件

(1)关键参数

通过前期查阅相关国内外学者对于流体泄漏研究内容以及对排泥库现场调研发现，排泥库泄漏点的影响半径与泄漏孔直径、含水量、动力黏度、泄漏孔形状、泥浆均匀程度、库区是否闭库等因素有关。但目前学者研究流体泄漏影响半径主要有泄漏孔直径、含水量、动力黏度，而泥浆泄漏方面缺乏研究。因此文章借鉴流体泄漏方面研究，主要考虑泄漏孔直径、含水量、动力黏度对泥浆泄漏影响半径的影响。泥浆泄漏孔直径取1、2、3 m，含水量取120%、150%、180%，动力黏度取1、5、10 Pa·s。

(2)边界条件

空气、泥浆两相流中的气体为理想可压缩气体，气体是可压缩的，泥浆不可压缩，空气随着泥浆泄漏过程不断下降，体积变大，密度和压力均不变，故模型的上表面采用压力输入入口。模型的出口存在大气压力，设为压力出口，模型壁面设为固体壁面，壁面出为非滑移边界条件。由于Interface的连接面之间的节点可以不一致，选用Interface类型。

模型网格底部划分如图4所示。

图4 网格底部

3.2 模拟结果分析

数值模拟对泄漏点主要影响因素不同大小泄漏孔口、泥浆含水量、动力黏度的分析结果如下。

3.2.1 泥浆在不同大小泄漏孔发生泄漏时的影响半径分析

为了更真实地模拟排泥库泥浆的泄漏情况，本次模拟的排泥库尺寸取50 m×50 m×30 m，泥浆的含水量设为120%，密度数值大小采用试验所测数值，见表1，泄漏孔直径大小分别设为1、2、3 m。

表1 泥浆含水量与密度的关系

由于泥浆探测器能探测到的最小泥浆速度数值未知，泥浆本身具有一定的动力黏度，流动现象没有清水那么明显，而流速太小又不能被探测到，故设置泥浆流速0.5 m/s作为流动半径的界定标准。为同模型试验相对应，使研究规律较明显，选取泄漏10 s时的状态作为研究对象。根据观测模型试验状态和现场实际结合，推测在泄漏孔泄漏时间在20、60 s…后，由于泥浆介质相对均质单一，泄漏孔的影响半径会保持不变。而当泥浆下沉到泄漏孔的影响范围时，出现一些肉眼可见的现象：表层逐渐出现漩涡，由外而内流速逐渐变大，在泄漏孔中心流速最大，下沉最多，呈明显的漏斗状分布。图5是模拟排泥库泥浆在不同泄漏孔直径下泄漏10 s时的速度云图。

图5 不同直径泄漏孔泥浆速度云图

根据以上对不同泄漏孔直径模拟结果，测出影响半径如表2所示。

表2 模拟影响半径结果

由图5、表2可知，在泥浆高度一样，孔口大小不一的情况下，洞口处的泄漏泥浆流速大小有所差别，洞口越大，流速越大。由于泥浆存在一定的动力黏度，在排泥库底部发生泄漏时，泥浆泄漏孔存在一定的影响半径，水平影响半径平均为泄漏孔半径的4.61倍，竖直影响半径平均为泄漏孔半径的5倍，泄漏孔中心流速最大，在库深30 m，孔口大小为1～3 m情况下，最大流速平均为22.8 m/s。由此能解释排泥库泥浆泄漏时，泥浆表面没有漩涡和漏斗、表层水或泥浆流动不明显的原因。

根据模拟结果可推测泥浆深度较大的排泥库泥浆流动过程如下：

(1)泥浆表层整体下沉，表面无漩涡和漏斗。

(2)泥浆下沉到泄漏孔的影响范围时，泥浆表层逐渐出现漩涡，由外而内流速逐渐变大，在泄漏孔中心流速最大，下沉最多，呈明显的漏斗状分布。

(3)泥浆泄漏完毕，由于排泥库底部泥浆密度和动力黏度均较大，含水量较低，泥浆内部存在相互作用力，无法完全从泄漏孔流出，在排泥库底发生淤积，上部密度和动力黏度较小，含水量较大的泥浆则全部从泄漏孔流出。

3.2.2 不同含水量泥浆发生泄漏时影响半径分析

由于排泥库泥浆从排矿口排出静置一段时间后，在重力作用下发生沉积，泥浆含水量发生变化，上部泥浆含水量较小，底部泥浆含水量较大。使用软件模拟泥浆分层较为复杂，故模拟过程均是采用含水量单一的介质。为了研究不同泥浆含水量对模拟结果是否产生较大影响，取泄漏孔直径为3 m的排泥库模型，分别模拟泥浆在含水量分别为120%、150%、180%情况下的泄漏影响半径。

由数值模拟可知，在泄漏孔口直径为3 m，含水量分别在120%、150%、180%情况下，泄漏孔影响范围基本相同，三种不同含水量速度云图大致如图6所示，只有数值上有些许差别，证明不同含水量对排泥库泥浆泄漏影响范围不大。

图6 含水量泥浆速度云图

3.2.3 不同动力黏度的泥浆发生泄漏时影响半径分析

由于泥浆与普通清水参数不一样，它具有较大的动力黏度，在泥浆内部发生相对剪切运动时，会产生相互阻碍的力，使得泥浆产生流动需要克服泥浆内部的阻力。为了研究泥浆在不同的动力黏度情况下是否会对泥浆泄漏影响范围产生影响，对不同动力黏度的泥浆泄漏情况进行数值模拟。模拟过程中，取泥浆泄漏孔口大小为3 m，泥浆动力黏度μ分别设置为1、5、10 Pa·s。

由图7可知，在泄漏孔口直径3 m，泥浆动力黏度μ分别为1、5、10 Pa·s情况下，泄漏孔处泥浆最大速度依旧仅有微小变化，三者数值分析图片只有数值上有些许差别，故只放置了一张云图。从得到的结果分析泥浆动力黏度对泥浆泄漏影响范围变化不明显，说明泥浆动力黏度对泥浆泄漏影响范围影响不大。

图7 动力黏度泥浆速度云图

4 结论

以广西铝土排泥库为研究背景，通过室内模型试验和Fluent 15.0流体计算软件研究了排泥库深层泥浆发生泄漏时的流动机理，得到了以下结论。

(1)排泥库泥浆含水量、密度和流速三者对泥浆渗流力的关联系数均较高，均可作为排泥库泥浆泄漏计算模型考虑的因素。

(2)运用多元线性回归计算理论建立了泥浆探测器所受渗流力与泥浆流速、含水量和密度之间的多元线性回归模型。并运用数理统计理论对此多元线性回归模型进行验证，证明结果显著。

(3)通过Fluent 15.0流体计算软件计算了在排泥库底部发生泄漏时，泥浆泄漏孔存在一定的影响半径，水平影响半径平均为泄漏孔半径的4.61倍，竖直影响半径平均为泄漏孔半径的5倍。并可以根据模拟结果推测厚度较大的排泥库泥浆流动过程。不同含水量和动力黏度对排泥库泥浆泄漏范围影响不大。

(4)通过对泄漏孔半径的影响因素分析，影响的最关键因素是泄漏孔半径，而泥浆含水量和动黏度系数则对泄漏孔影响半径不大。随着排泥库底部泄漏孔半径的增加当泥浆发生泄漏时其影响半径最大，会出现泥浆表面下沉、产生旋涡等现象。因此在进行排泥库选址时应充分探测熟知该地区的地质构造情况，必要时对危害排泥库安全的溶洞进行施工处理，这在源头上进行防控和治理有极大的意义。

首次对排泥库泥浆泄漏机理进行研究，采用自制探测装置对排泥库泥浆深层泄漏点的快速查找，弥补了排泥库泄漏机理研究领域的空白，研究成果的实际应用将大幅减少泄漏事故的危害和损失，对排泥库泄漏的安全防控和后期治理起到极大的指导意义。