煤炭资源是内蒙古的支柱性矿产,已探明储量居全国首位。近年来内蒙古一些大型煤炭企业积极和电力行业融合发展,大力开展煤电一体化项目[1]。其中,鄂尔多斯地区煤炭资源丰富,火力发电厂集中意味着鄂尔多斯地区的燃煤量将大大增加。但在燃煤发电的同时,电厂产生的污水和固体废弃物经淋溶作用产生的淋滤液极易对土壤造成污染,甚至通过运移扩散作用对地下水产生影响。土壤和地下水具有埋藏性和系统复杂性,污染问题不如大气和地表水污染表现得直观,因而长期受到忽视,使得污染问题更加突出,因此燃煤电厂周边土壤和地下水的污染问题应引起足够重视。
近年来,已有相关学者对燃煤电厂所在地区土壤产生的影响进行研究,宁东基地燃煤电厂周边表层土壤程轻微生态风险,深层土壤存在重金属污染风险[2];吉林省某电厂周边农田土壤受多环芳烃污染,来源于化石燃料、木材以及高分子化合物的燃烧[3];早期建造的燃煤火电厂大多已接近或超过退役年限,即将面临关闭,大量的工业构筑物和废弃物遗留在场地中,对土壤和地下水存在一定的环境污染风险[4]。然而这些研究大多并未对污染物进一步污染深层土壤甚至地下水的可能性进行研究,缺乏将电厂污染物与所在地区土层相结合的一体化研究。而且研究的污染物类型集中为重金属和多环芳烃,忽视了燃煤电厂产生的其他种类污染物对土壤和地下水的影响。近几年对电厂废水排放检测后发现,在一些废水处理设施出口以及排放口处氨氮的检出率增大。而现有的电厂废水处理系统除了生活污水处理,其他都不具备处理氨氮的能力。据调查氨氮废水的来源主要有生活污水、脱硫废水和化学精处理废水[5]。粉煤灰是燃煤电厂排放的固体废弃物,目前虽然许多大型电厂都采用了干贮灰技术,但是雨水的淋滤作用仍然会引起泄漏。灰渣淋滤液中主要污染因子包含氟,燃煤电厂在湿式除尘过程中产生大量氟浓度超标的废水[6-7],一旦发生管道泄漏、防渗层破损或雨水冲刷等情况,氨氮、氟离子则会进入非饱和带,存在地下水被污染的风险。因此,研究燃煤电厂污染物氨氮和氟离子在电厂区域土层中运移规律具有重要意义。本研究区域选择在燃煤电厂较为集中的内蒙古鄂尔多斯市,该地区某燃煤电厂位于鄂尔多斯市东部,毛乌素沙漠东南端,地处黄土高原。地表被大厚度砂质粉土所覆盖,一般厚度为20~90 m,最大厚度约为160 m,土壤孔隙较大,含少量姜结石,湿度较小。另外,区内部分地段地表有厚度不等的风积细砂覆盖,结构松散、湿度小、持水性差。
目前,溶质在非饱和带土壤中的迁移规律一般采用室内土柱淋滤实验和数值模拟进行研究,Porro等[8]在瞬时供给溶质和恒定不变供给方式下通过室内土柱试验研究了溶质运移规律。对于包气带中溶质运移的数值模拟,国内外有很多学者采用Hydrus-1D模拟宏观和微观尺度上一维饱和及非饱和介质中的水流、溶质、热、二氧化碳等的运移及反应情况,该模拟软件以确定性机理模型为主体[9]。Jellali等[10]使用Hydrus-1D软件建立一维运移模型,发现平衡模型模拟在模拟氨氮迁移规律中较非平衡模型效果更好。
因此,本文针对鄂尔多斯市某燃煤电厂污染物氨氮、氟离子在细砂、砂质粉土层中的运移规律进行研究,并运用Hydrus-1D平衡溶质运移模型模拟预测非饱和带对污染物的截留能力,以期为燃煤电厂污染物在非饱和带介质中的运移规律研究提供参考。
土样采自燃煤电厂所在地区的非饱和带土壤——细砂、砂质粉土,去除植物残根、石灰结核等杂质,自然风干、碾磨、过1 mm筛备用。采用激光粒度分析仪和X-射线衍射仪对土样进行颗粒分析和检测,发现2类土样主要含结晶态物质为石英砂,土样基本参数见表1。
试验中污染物浓度的确定参考当地多家电厂中水样污染物浓度,以用于数值模拟,进行淋滤实验时适当增大污染物浓度以提高数据精度,检测得到氨氮平均浓度为30 mg/L,氟离子浓度为3.5 mg/L。
表1 非饱和带土壤理化性质分析
Table 1 Physical and chemical properties analysis
土壤类型含水率/%pH值总有机碳/%干密度/(g·cm-3)土壤粒度分布砂粒/(%,>20 μm)细粉粒/(%,20~2 μm)黏粒/(%,<2 μm)细砂0.137.010.721.4897.062.940砂质粉土10.87.360.3781.196.1977.7216.09
通过吸附动力学试验和等温吸附试验得到土样对污染物最大吸附容量及分配系数。动态淋滤采取土柱实验的方法,其主要目的是模拟研究溶质运移过程,和相同水文地质条件下的土壤及地下水污染问题[8]。土柱设计及填装按照土柱实验标准进行,土柱材质为有机玻璃,为了避免产生侧壁优先流,实验前用砂纸对柱壁进行均匀打磨,填装土柱时采用干堆法。
对溶质在非饱和带运移的模拟实验分为模拟柱弥散实验和污染物淋滤实验,前者通过NaCl溶液在一维土柱中的穿透曲线,计算2种土样的饱和弥散系数DL与弥散度α;后者通过进行污染物定水头淋滤实验,在土柱底端口取样,分析污染物在取样口浓度随时间变化。
进行弥散实验时,首先使用超纯水对土柱进行淋滤,去除土柱本底值,淋滤结束后配制1 g/L NaCl溶液,置于马氏瓶中,进行NaCl弥散实验。细砂土柱间隔5 min取样,砂质粉土柱间隔1 h取样,使用电导率仪测量其电导率,后期取样时间可适当延长,根据标准曲线绘制Cl-穿透曲线并计算参数弥散系数DL与弥散度α。进行污染物淋滤实验时,同样使用超纯水淋滤去除本底值,配制100 mg/L的NH4Cl溶液以及15 mg/L的NaF溶液,细砂土柱间隔10 min取样,砂质粉土柱间隔30 min取样,使用紫外分光光度计和离子色谱分别检测NH3-N和F-浓度,后期取样时间可适当延长,直到达到其初始浓度为止,绘制穿透曲线。
在进行Cl-、NH3-N、F-的淋滤实验后,选择Hydrus-1D软件进行污染物在非饱和带的运移过程模拟。根据溶质和非饱和带的特点,选择VG水流运移模型和平衡溶质运移模型,设置各项参数初始值,进行参数反演、调参,建立2种特征污染物在典型非饱和带介质中的运移模型。
根据van Genuchten-Mualem水流模型,水分特征曲线所需参数由粒度分析实验(Bettersize 2000激光粒度仪)得到通过国际粒度分级标准分类得到非饱和带介质的颗粒粒径组成,骨架密度由土柱装填质量计算得到,然后通过 Hydrus-1D软件中的神经网络预测功能获取(表2)。
表2 土壤水力特征参数
Table 2 Soil hydraulic characteristic parameters
种类Qr/(cm3·cm-3)Qs/(cm3·cm-3)Alpha/m-1nlKs/(cm·min-1)细砂0.04830.35180.03462.85910.50.1209砂质粉土0.07530.47680.00591.63310.50.043882
注:Qr为土壤剩余含水率;Qs为土壤饱和含水率;Ks为饱和导水率;Alpha、n、l为VG模型公式中的经验常数。
将水力特征参数建立水流模型与Cl-弥散实验的实测穿透曲线进行拟合,调参后将校准值弥散度α运用到溶质运移模型的建立,选用不间断供水模式,设置NH3-N污染物模拟时长分别为400,1200 min,设置上边界为100 mg/L的定浓度边界;F-模拟时长分别为480,720 min,设置上边界为15 mg/L的定浓度边界。溶质运移模型中的反应参数按照静态吸附实验中获得的Langmuir等温式中吸附参数设置,具体数值见表3。
表3 溶质反应参数设置
Table 3 Parameters’ setting for solute reaction
污染物种类土样种类α/cmKd/(mg·mL-1)Nu/(mL·mg-1)氨氮细砂0.9382.9150.1砂质粉土0.754.0445.7氟化物细砂0.9385.226砂质粉土0.756.6628
进行溶质运移模型的验证和调参,使用Inverse Solution进行参数反演,运行模型。为进一步提高模拟数据与实测值的拟合度调整参数,选取由静态吸附实验获得的溶质反应参数Kd和Nu作为反演参数对象,进行参数反演并反复调参,优化模型。
2.1.1 Cl-穿透曲线
根据NaCl弥散实验结果绘制穿透曲线(图1a),采用插值法计算获得土壤弥散系数DL与弥散度α(表4),以此参数为基础进行模型建立及调参。由NaCl穿透曲线可发现:Cl-相对浓度逐渐增大,且穿透速率均先缓慢上升,达到峰值后,缓慢降低。
—细砂; —砂质粉土。
注:W=土柱底端流出污染物浓度/进样污染物浓度。
图1 动态淋滤试验穿透曲线
Fig.1 Penetration curve of dynamic leaching test
2.1.2 NH3-N穿透曲线
NH3-N淋滤实验的穿透曲线(图1b)与NaCl弥散实验穿透曲线相似。由于吸附作用的原因导致完全穿透时间较长,在细砂土柱中需280 min完全穿透,在砂质粉土柱中需1000 min完全穿透,说明NH3-N在介质中的运移远远落后于水分。因为Cl-在非饱和带中仅受到水流的对流弥散作用,而NH3-N则受到各种物理化学作用和对流弥散共同作用。NH3-N在2种土壤介质中运移产生的穿透曲线均表现为快速增加再变缓,可大致分为3部分:第1阶段相对浓度为0,在此阶段非饱和带中各部分的吸附位点充分,大部分被土壤所吸附且土壤吸附能力较强,到达底部时浓度为0;第2阶段中土柱中各部分的吸附位点逐渐被占据,淋滤液污染物浓度迅速上升,来自上层的淋滤液迅速代替土柱中原有水分,吸附速率降低;第3阶段吸附速率继续降低,淋滤液浓度增加缓慢,吸附量几乎达到最大值,溶液浓度也达到初始值[10]。其中砂质粉土土柱的氨氮穿透曲线在中后期出现剧烈波动,随着时间推移,波动幅度逐渐变小,这是由于土柱的吸附位点已基本饱和,随着污染液对土层的冲刷,土柱内部结构发生变化,形成一些细小的裂隙,使得淋滤液通过一些渗透性较强的通道,导致浓度不断波动最终趋于平稳[11]。
2.1.3 F-穿透曲线
F-淋滤的穿透曲线(图1c)与NaCl弥散试验穿透曲线形态相似,由于吸附作用导致完全穿透的时间较长。在细砂土柱和砂质粉土柱中分别需300,700 min完全穿透。F-的穿透曲线与Cl-相似,但存在时间差。F-在土柱中向下运移,除了对流弥散作用,还有吸附作用,土壤介质对F-的吸附较NH3-N小,因此F-的穿透曲线落后于NH3-N穿透曲线[12-13]。
2.2.1 参数校验
在正向模型的基础上,使用Inverse Solution将实测数据与模拟曲线进行拟合,Hydrus-1D模拟运行结果显示,相关系数分别为0.9801和0.9905,相关性较好。为了进一步提高相关性,选取纵向弥散系数α作为反演参数对象进行参数反演,优化得到模拟和实测值对比(图2),调参前后的拟合结果见表4。
表4 反演参数结果
Table 4 Parameter inversion results
土壤类型实测值α/cm反演前相关系数校准值α/cm反演后相关系数细砂0.8920.98010.9380.99676砂质粉土0.790.99050.750.99858
通过微调纵向弥散系数α进行参数反演,使得模拟穿透曲线与实测穿透曲线相关系数增大,提高拟合度,从而使得参数不断优化,以进行后期的预测模拟。由图2可以看出,在调参前后,模拟穿透曲线与实测值吻合度提高,在参数微调的前提下,相关系数均增大至0.996以上,可以采用参数纵向弥散系数0.938 cm、0.75 cm对污染物在细砂土、砂质粉土中运移规律进行预测模拟研究。
水动力弥散系数是用来描述一定水流流速下,多孔介质对某种溶质弥散能力大小的参数,一定程度上反映了多孔介质中水流过程和空隙结构特征对溶质运移过程的影响。介质的弥散度受到许多因素的影响,包括介质的含水率,介质隙的大小、数量、连通性,骨架密度以及溶质性质等,同时,尺度效应也是影响弥散系数的因素。本次模拟中发现反演的最佳弥散度与实测弥散度存在一定差异,分析是由装填土柱的压实程度、密度的微小差异等造成的[14]。
注:实线为模拟值,空心点为实测值,下同。
图2 NaCl弥散实验拟合结果
Fig.2 Fitting curves of NaCl tracer test
2.2.2 模型校准
针对NH3-N和F-淋滤实验的数值模拟进行参数反演,选择溶质反应参数——分配系数Kd和Nu作为反演参数对象,运行模型并进行多次调参。模拟图见图3—4,反演前后参数对比见表5。
图3 氨氮和氟化物淋滤试验拟合结果
Fig.3 Fitting curves of NH3-N and fluoride leaching test
由图3、4可以看出,对反演参数Kd和Nu进行微调后,在调参后的拟合图中,模拟值与实测值的吻合度显著提高,相关系数大幅增加,进行优化后的参数Kd和Nu可运用至后期预测模拟中,预测不同污染情况下污染物在细砂土、砂质粉土中的迁移规律。
经过对NH3-N、F-在非饱和带中运移模型的验证和参数修正,建立了能够预测不同工况下污染物运移规律的模型。本节将基于实际工况,结合反演调参后的相关参数,预测污染物穿透非饱和带厚度随污染物变化、持续排放污染物和突发性污染事故下污染物运移规律。
图4 穿透厚度与污染物浓度预测关系
Fig.4 Relationship between soil interception thickness and pollutants concentration
在时间步长固定的前提下,探究污染物穿透土层厚度随污染物浓度变化情况,根据检测到的电厂实际污染物排放浓度设计相应浓度梯度。细砂土层模拟时长为90 d,砂质粉土层为1年的污染物浓度与截污厚度关系见图4。
表5 反演前后参数对比
Table 5 Parameter comparison table
污染物种类非饱和带介质Kd/(mg·mL-1)Nu/(mL·mg-1)相关系数实测值校准值实测值校准值反演前反演后NH3-N细砂2.912.550.150.90.8950.924砂质粉土4.043.745.7470.8960.967氟化物细砂5.24.832628.910.9730.996砂质粉土6.665.922829.160.9440.987
注:Kd为最大分配系数;Nu为Langmiur反应参数。
污染物在细砂土层中运移快速,在运行时间为90 天前提下,氨氮浓度为43.4 mg/L时将穿透100 m厚土层,两者呈线性正相关,并得到线性方程如图4(注:D-穿透厚度,c-污染物浓度);氟化物浓度为164.6 mg/L时将穿透100 m厚土层。
砂质粉土对污染物的截留效果较好,针对其特性设定模拟时长为一年,根据模拟结果得出当氨氮浓度为50.24 mg/L时将穿透厚度为100 m的土层;氟化物浓度为134.9 mg/L时穿透厚度为100m的土层。
考虑到燃煤电厂入渗的实际条件,选定非饱和带平均厚度10 m作为模拟厚度,在模拟过程中,忽略弥散系数的尺度效应(认为弥散系数不变)。根据燃煤电厂采集水样结果发现,污染物平均浓度中NH3-N为30 mg/L,F-为3.5 mg/L。模型设置中将非饱和带厚度设为10 m,设置情况见表6,模型运行结果见图5。
表6 模拟输出时间一览表
Table 6 Output-time list by the simulation
T0/dT1/dT2/dT3/d细砂03.5710.5砂质粉土0102030
注:输出时间为T1代表细砂土模拟输出时间为第3.5天,砂质粉土模拟输出时间为第10天。
观察对比特征污染物在2种非饱和带介质剖面上的分布规律发现:在相应的模拟步长中,污染物在砂土介质中的运移速率和运移深度明显大于砂质粉土介质,这主要是由于土壤粒径分布差异所造成的,土壤颗粒越细,土壤黏粒含量越多,土壤渗透系数越小,吸附性能越大。
观察污染物在非饱和带中随时间变化的运移规律,可发现污染物下渗速率随着时间推移而降低。对于细砂而言,T1时NH3-N、F-运移深度分别为380,180 cm;T2时分别为630,210 cm;T3时分别为970,240 cm。对于砂质粉土来说,T1、T2、T3时NH3-N、F-运移深度分别为230,210 cm;430,250 cm;620,400 cm。
突发性污染事故发生时污染物的剂量较大,可迅速入渗到非饱和带,后续无污染物下渗。本次模拟与3.2节情况污染物总量相等,进行污染物一次性入渗模拟。污染物在1 d内全部释放,观察污染物在180 d 内的垂向分布特征,一次性入渗浓度NH3-N为5400 mg/L,氟化物为630 mg/L,分3次输出时间(表7)运行模型,得到NH3-N和F-在非饱和带剖面中的垂向分布图(图6)。
——T0; ——T1; ——T2; ——T3。
图5 持续排放下污染物垂向分布
Fig.5 Vertical distribution of pollutants under continuous discharge
表7 模拟输出时间一览表
Table 7 Output-time list by the simulation
T0/dT1/dT2/dT3/d细砂0510180砂质粉土0510 180
——T0; ——T1; ——T2; ——T3。
图6 突发性污染事故下污染物垂向分布
Fig.6 Vertical distribution of pollutants under a simulated sudden pollution accident
污染物一次性入渗情况下,由图6a、b可知,两种污染物在T2(10天)时均会穿透10 m厚度的细砂土层进入地下水中, 运行至T3(180天)时已完全穿透土层。由从图6c、d可知,运行至T3(180天)时,NH3-N污染物穿透10 m厚的砂质粉土层,氟化物未穿透土层。在突发性污染事故情况下,一次性大剂量的污染物下渗至土层中,砂质粉土层的截污效果优于细砂土层。
1)土柱淋滤实验表明:污染物在细砂土中运移速度较大,NH3-N、F-污染物完全穿透细砂土柱所需时间比砂质粉土柱分别少720,400 min。
2)运用Hydrus-1D建立模型,发现模拟值与实测值相关性较好,相关系数为0.950~0.996,模型校准后,获得用于预测的最佳参数值,细砂土与砂质粉土弥散度α为0.937,0.75 cm,NH3-N和F-在细砂土与砂质粉土中运移时的溶质反应参数Kd、Nu分别为2.5 mg/mL、50.9 mL/mg和4.83 mg/mL、28.91 mL/mg。
3)模拟预测NH3-N浓度为43.4 mg/L、F-为164.6 mg/L在90 d后穿透100 m厚细砂土层,浓度为50.24 mg/L NH3-N、134.9 mg/L的F-在一年后穿透100 m厚砂质粉土层;对持续排放污染物情况进行模拟可知,污染物可在短时间内穿透10 m厚细砂土层,而砂质粉土层对污染物截留能力相对较优;模拟发生突发性污染事故,高浓度污染物一次性入渗后,污染物在10天时均会穿透10m厚度的细砂土层,180天时NH3-N污染物穿透砂质粉土层,氟化物未能穿透,砂质粉土层的截污效果较好。
4)砂质粉土层截污能力优于细砂土层,判断污染物是否能够穿透非饱和带进入地下水取决于污染物浓度、土层质地与厚度、污染物排放时间等因素,本文的模拟预测研究可提供一定的参考。
本文以鄂尔多斯区域燃煤电厂为背景进行污染物迁移研究,该地区土壤质地较为松散,黏粒含量较少。下一步工作将深入研究黏土等黏粒含量较高的土层对污染物的截留能力,以选取能有效截污的土质作为燃煤电厂防渗层。
[1] 王道涵,杨亚利,任鹏,等. 内蒙古某电厂周围土壤汞分布特征及其影响因素[J].环境工程学报,2015,12(9):6135-6140.
[2] 罗成科,张佳瑜,肖国举,等. 宁东基地不同燃煤电厂周边土壤5种重金属元素污染特征及生态风险[J].生态环境学报,2018,27(7):1285-1291.
[3] 崔政武,王洋,于锐,等. 吉林省电厂周边农田土壤中多环芳烃含量特征及风险评价[J].环境污染与防治,2018,40(7):806-811.
[4] 刘虎,刘丹丹,吕倬,等. 燃煤火电厂生态恢复景观设计策略研究[J].天津大学学报(社会科学版),2017,19(4):377-379.
[5] 乐园园.火电厂氨氮废水来源及处理可行性探讨[J].电力科技与环保,2016,32(4):16-18.
[6] 傅嘉媛, 董榕.粉煤灰中主要有害元素的浸出与迁移试验研究[J].福州大学学报(自然科学版),2004,32(1):123-126.
[7] 白继红,张永波. 电厂粉煤灰场氟离子对地下水影响的试验[J].地球科学与环境学报,2008,30(4):408-411.
[8] Porro I, Wierenga P J. Transient and steady-state solute transport through a large unsaturated soil column[J]. Groundwater, 1993, 31(2): 193-200.
[9] 杨欢. 基于平衡和非平衡模型的包气带土壤中氨氮运移过程研究[D]. 北京:中国地质大学, 2015.
[10] Jellali S, Diamantopoulos E, Kallali H, et al. Dynamic sorption of ammonium by sandy soil in fixed bed columns: evaluation of equilibrium and non-equilibrium transport processes[J]. Journal of Environmental Management, 2010, 91(4): 897-905.
[11] Turkeltaub T, Kurtzman D, Dahan O. Real-time monitoring of nitrate transport in the deep vadose zone under a crop field-implications for groundwater protection[J]. Hydrology & Earth System Sciences Discussions, 2016, 20: 1-31.
[12] 陈初雨,王延亮,孙春,等. 粉煤灰-粉质粘土对氟离子的吸附作用[J]. 吉林地质, 2013(3): 90-92,95.
[13] 焦新庆. 氟、砷污染物在包气带(非饱和带)中的迁移规律的研究[J]. 环境科学与技术, 1988(2): 4-12.
[14] 曹巧红,龚元石. 应用Hydrus-1D模型模拟分析冬小麦农田水分氮素运移特征[J]. 植物营养与肥料学报, 2003, 9(2): 139-145.