中国的离子型稀土(或称风化壳淋积型稀土)分布在江西、福建、湖南、广东、广西、云南和浙江等南方7省，而江西赣南不但是稀土开采技术的发祥地，也是国内离子型稀土的重要资源地与供应地，占到全国同类资源的30 %，素有“稀土王国”之美誉^[1-4].然而，赣南稀土历经40多年的开采，在追求经济效益的同时，付出的环境代价是相当巨大且惨重的^[5].

离子型稀土矿开采先后经历了池浸、堆浸和原地浸矿3种开采工艺^[6]，池浸和堆浸对地表植被破坏巨大，易造成水土流失；原地浸矿通过注液孔注入矿体母液(硫酸氨等化学成分)，将矿物中稀土离子交换浸出，然后，通过集液沟收集稀土母液^[7].然而，无论使用哪种采矿工艺，都会对地下水造成污染，尤其是现在推广的原地浸矿工艺，由于受地质条件及开采技术限制，母液进入矿层不能全部回收，部分浸矿液会随着时间经土壤入渗到地下造成不同程度的水土环境污染^[8-10].由于国外稀土矿开发少，稀土矿产资源开发对水土污染影响尚未进行深入、全面系统的研究.国内学者对稀土矿矿开采引起的水土污染做了大量研究，主要集中在稀土矿不同开采工艺对水土污染的影响，涉及到的污染因子有氮化物污染、开采后尾矿残留稀土元素的污染、重金属污染及放射性污染等^[11-13].目前，对离子型稀土原地浸矿通过数值模拟研究其对水土污染规律尚存不足.为此，文中以江西龙南足洞矿区为例，利用FEFLOW 6.1数值模拟软件^[14]，对离子型稀土矿开采在地下水中氨氮进行数值模拟，预测地下水氨氮污染空间分布与特征.

1 研究区概况

矿区位于江西省龙南县城东南相距10 km处，面积为23.306 km²，开采深度为250~573 m.根据水文地质勘查范围和足洞矿区所在区域，本次研究区的面积为90.6 km²(见图 1)，该矿研究区域为山区，地形起伏变化大，地面标高为200~450 m.

2 水文地质条件

根据当地的气象资料，区内年平均降雨量1 510.8 mm，其中每年的4~6月为丰水期，占全年降雨量的56.4 %，10月至翌年的元月为枯水期，占全年降雨量的14.2 %，而2月、7月、8月、9月等4个月为贫水期，当地年降雨量还与地貌、地形的高低有关，从平地到山地有降雨量随地势的增高而增大趋势.

研究区的地下水主要赋存于基岩裂隙与第四系松散层中，基岩裂隙含水层，在调查区内广泛分布.稀土矿选场、原地浸矿的注液孔落在该含水层之上.基岩裂隙水含水层以风化带网状裂隙水为主，其赋存在花岗岩的风化裂隙中，含水不均，富水性差，地下水水位埋深为2.85~4.26 m，单井涌水量为1.06~16.44 m³/d，水力性质为潜水，局部为微承压.

第四系冲洪积砂砾石孔隙潜水含水层，主要分布在调查区河谷两侧及山间、山前低洼处.稀土矿采、选矿尾水处理厂就座落在该含水层之上.水力性质为潜水，地下水水位埋深一般为0.5~1.0 m，局部可达2.0 m，单井涌水量一般为5~20 m³/d，水量贫乏.

3 研究区地下水渗流数值模型 3.1 水文地质概念模型

研究区西部边界为东江河，东北部为濂江河，河水与河流沉积的第四系松散地层具有直接的水力联系.为此，将该部分边界定为一类水头边界，其余边界则定为二类流量边界.

把整个矿区看做一个完整的水文系统和独立的水文地质单元，以保持系统的完整性，提高地下水数值模拟预测的精度.由于地下水埋深变化幅度较大，基岩裂隙水含水层风化裂隙往往透水而不含水，未形成统一地下水位，因此，模型结构按一层潜水含水层考虑.研究可概化为二维非均质、各向同性、单一结构的非稳定地下水流系统^[15].

3.2 水文地质数学模型

根据上述水文地质条件，可用下述偏微分方程及定解条件来描述研究区内水流运动：

$ \left\{ \begin{array}{l} \mu \frac{{\partial h}}{{\partial t}} = \frac{\partial }{{\partial x}}(K(h - b)\frac{{\partial h}}{{\partial x}}) + \frac{\partial }{{\partial y}}(K(h - b)\frac{{\partial h}}{{\partial y}}) + p + \varepsilon \;\;\;x,y \in \Omega ,\begin{array}{*{20}{c}} {}&{t > 0} \end{array}\\ {\left. {h(x,y,t)} \right|_{t = 0}} = {h_0}(x,y)\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;x,y \in \Omega \\ {\left. {H(x,y,t)} \right|_{{\Gamma _1}}} = \phi (x,y,t)\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;x,y \in {\Gamma _1},\begin{array}{*{20}{c}} {}&{t > 0} \end{array}\\ {K_n}{\left. {\frac{{\partial h}}{{\partial \bar n}}} \right|_{{\Gamma _2}}} = q(x,y,t)\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;x,y \in {\Gamma _2},\begin{array}{*{20}{c}} {}&{t > 0} \end{array} \end{array} \right. $

其中：Ω表示评价区域；h为含水层地下水的水位标高(m)；（m）为底板高程(m)；μ为含水层给水度；t为时间(d)；K为含水层渗透系数(m/d)；K_n为边界面法向方向的渗透系数(m/d)；ε为地下水开采强度(m/d)；p为为降雨入渗量及蒸发量的代数和(m/d)；h₀(x, y)为含水层的初始水位分布(m)；Γ₁为一类边界；Γ₂为二类边界；n为边界面的法线方向；φ(x, y, t)为一类边界水头(m)；H(x, y, t)为含水层边界的水头(m)；q(x, y, t)为二类边界的单位面积流量(m³/(d·m²))，流入为正，流出为负.

3.3 网格划分

运用基于有限元原理的FEFLOW 6.1平台软件对地下水流进行数值模拟.对研究区进行三角网格剖分，网格剖分时，将水位统测点、地下水开采井、河流及其主要支流放置于结点上，在地下水开采处和足洞矿区内的拟采矿区域进行适当加密.将研究区的单元剖分划分，共计剖分三角单元18337个，结点18594个，如图 2所示.

网格剖分后，依次将地面高程、初始地下水水位、含水层底板高程等基础数据输入模型.根据研究区已有的水文地质勘查资料，尤其是抽水试验、压水试验获得的含水层渗透系数，进行水文地质参数分区，并合理赋值，作为初始水文地质参数.模型要素全部输入模型后，即可运行模型，得到研究区各结点的计算水位.

4 模拟研究区域

根据探明的稀土保有储量及下达足洞稀土开采总量指标，足洞矿区的开发利用方案为：稀土矿块分8年开采，共计79块，开采面积合计300.54万m².每年开采的矿块数及面积详见表 1.足洞矿区各开采矿块的空间分布见图 3.

5 污染源强确定

为研究原地浸矿母液对周边地下水环境影响，开展了室内稀土淋溶试验和标准矿块开采现场试验，以确定开采经济指标和污染源指标.

5.1 室内稀土淋溶试验

吸附在矿物表面的稀土阳离子在遇到化学性质更活泼的阳离子时，能被更活泼的阳离子交换解吸而进入溶液，阳离子交换方程为^[16-18]：

2(高岭土)^-3 RE⁺³+3(NH₄)₂⁺¹SO₄^-2→

2(高岭土)^-36(NH₄)⁺¹+ RE₂⁺³(SO₄)₃^-2

根据试验块室内的稀土淋溶试验数据，每开采1 (mol)的稀土元素，需要1.5 (mol)的硫酸铵.因此，该试验中阳离子交换需要的硫酸铵为0.18 (mol)，计硫酸铵质量为24.0 g，占留在稀土矿岩体中硫酸铵为19%，清水淋溶出的硫酸铵占留在稀土矿岩体中硫酸铵为67%.

5.2 标准矿块现场试验

1) 渗入含水层中水量.现场试验以6000 m²标准矿块为开采方案，开采按1年期计.提取出标准块中的稀土元素需要顶水量加母液36500 m³，试验得出顶水量+母液回收率为90%.则开采后10%的顶水量+母液进入包气带和含水层中.按照最大风险考虑该部分水全部渗入含水层中.标准块开采1年，母液+顶水渗入含水层中的量为3650 m³/a，计10 m³/d.

2) 单位体积含水层污染源强.根据标准矿块的现场试验方案，需要注入硫酸铵190 t，回收的硫酸铵质量为109.36 t；开采结束后，留在稀土矿岩体中的硫酸铵质量为80.64 t；渗入含水层中污染物的量为54.02 t，其中氨氮质量为11.46 t，硫酸根为为39.29 t.

6 地下水污染预测

由室内稀土淋溶试验和标准矿块现场试验，可计算出各矿块的面积渗入含水层中的水量以及渗入含水层中氨氮的量，得出每天单位体积含水层污染源强，即可构建各开采矿块开采过程中对地下水污染的溶质运移模型.将各开采矿块视为地下水污染面源、滞留于含水层中的水量和各污染因子的污染源质量(g·m^-3/d)赋与模型.结合足洞稀土矿开采完所有矿块需要8a方案，故对地下水环境影响预测年限设定为20a，氨氮浓度空间分布的时间点分别为100 d、1 000 d、10 a和20 a.即可模拟预测各矿块开采过程中及开采后，地下水氨氮污染空间分布和特征.

6.1 地下水氨氮污染空间分布

本次地下水氨氮浓度预测中，按照地下水质量标准(GB/T14848-93)^[19]，地下水质量Ⅲ类标准氨氮为0.2 mg/L，当地下水中的浓度大于Ⅲ类标准时，即认为稀土开采对地下水环境产生了影响.以地下水氨氮的检出限Ⅰ类和Ⅱ类标准(均为0.02 mg/L)影响界限，标识最小影响范围.

经模型预测，足洞矿区的稀土开采对地下水中氨氮浓度的影响分布见图 4，预测地下水氨氮浓度空间分布的时间点分别为100 d、1000 d、10a和20 a.

由图 4可知，随着时间的推移(100 d→1 000 d→10 a→20 a)，地下水中氨氮浓度的影响范围越来越大，尤其是在矿块开采过程中(100 d→1 000 d→10 a)，其污染范围成2倍数以上增长；当矿块开采完成后(10 a→20 a)，由于受到降雨入渗的稀释作用其浓度又呈现降低的趋势.仅靠自然稀释作用，其超标面积、影响面积及矿区外超标面积均呈扩大趋势，与开采过程中对比影响范围虽放缓，但并不能在短期内有效阻止地下水氨氮的污染.

6.2 特征点浓度变化分析

在足洞矿区周围有60个村庄，100余眼地下水开采井，地下水主要作为饮用水源和生活用水.根据对这一区域的长期水质监测，地下水中氨氮污染较为突出，测出所在区域地下水中氨氮超出地下水Ⅲ类标准的近1 000倍.

为了进一步研究氨氮浓度变化，选择16个点作为特征点，特征点选择原则为：位于影响范围内的井点作为特征点；另外在开采矿块边缘且靠近足洞矿区界线设置4个特征点.开采井和特征点分布见图 5.各特征点的位置及氨氮浓度变化过程中的最大值详见表 2.浓度变化曲线见图 6.考虑特征点的属性及位置，将特征点分为4组，第1组位于矿块边缘，即BY1、BY2、BY3、BY4；第2组位于超标区范围内的水井，即L35、L109、L24、L112、L113、L114；第3组位于外围的未超标区域内的水井，即L37、L110、L117；第4组特征点为距离矿块最近的河流边缘点，即XL1、XL2、XL3.

从图 5可看出氨氮20a超标范围和未超标范围，大部分水井位于氨氮影响范围外，仅小部分水井位于影响范围内.

根据模型预测的地下水氨氮浓度结果，矿区开采对地表水5号、10号、11号溪流会产生一定程度的污染，氨氮进入河流的最大污染浓度为达到1 400.8 mg/L.矿区边缘的6个自然村的居民用水井也受到不同程度的污染，氨氮最大污染浓度为49.18 mg/L.

结合图 5、图 6和表 2可知：

1)  第1组特征点点位于矿块边缘，第4组特征点距矿块较近，直接受矿块开采的影响，矿块开采过程中地下水中氨氮浓度快速升高，开采期第8年左右浓度达到最高值，而后呈指数曲线形式快速下降，并渐趋平稳.

2) 第2组特征点距开采矿块较远，与第1组特征点相比，其氨氮浓度的升高明显滞后，且最高浓度值仅为49.18 mg/L，显著低于第1组的最高浓度值；矿区边缘的6个自然村的居民用水井在开采期基本不受影响，当矿块开采完后，残留在矿体中的氨氮随降雨入渗到居民用水井中，导致居民用水井氨氮浓度有上升趋势.

3)  与第1、2组特征点相比，第3组特征点的氨氮浓度最低.处于同一流域内第2、3组特征点中，第3组特征点浓度升高时刻滞后于第2组特征点.

4)  受降雨淋溶影响，各特征点氨氮浓度变化呈现出丰水期低、枯水期高的变化特点，凸显了降雨淋溶的稀释作用.

7 结论

1) 采用FEFLOW软件对矿区地下水氨氮污染进行数值模拟预测，软件较准确地反映了地下水氨氮浓度动态变化，以图形输出方式定量表现出稀土矿开采期间及开采完后地下水氨氮变化趋势及其影响范围，为地下水污染防治措施提供了理论依据.

2) 稀土矿区开采对地表水5号、10号、11号溪流会产生一定程度的污染，矿区边缘的6个自然村的居民用水井也受到不同程度的氨氮污染.

3) 矿区开采过程中随着时间的推移，地下水中氨氮浓度的影响范围越来越大.地下水中氨氮浓度超标范围也不断增加.在矿块开采过程中地下水中的氨氮浓度快速上升；当矿块开采完成后，降雨入渗的稀释作用使其浓度又呈现降低的变化趋势.但仅靠自然稀释作用，并不能在短期内有效阻止地下水的氨氮污染.

4) 矿块边缘直接受矿块开采的影响，矿块开采过程中地下水中氨氮浓度快速升高，而后呈指数曲线形式快速下降，并渐趋平稳.

5) 通过在稀土矿区建立完善的地下水模拟预测模型，对制定科学合理的开采方案，优化矿区地下水环境管理，推动生态环境与经济协调发展，从而实现赣南离子型稀土矿区的可持续发展具有重要的理论和现实意义.