Pu在膨润土层中的反应性迁移模拟研究(1)——地下水-膨润土体系演化的地球化学模拟

刘东旭,黄流兴,赵振华,胡立堂,司高华,叶远虑

(1.西北核技术研究所,陕西 西安 710024;2.北京师范大学 水科学研究院,地下水污染控制与修复教育部工程研究中心,北京 100875;3.生态环境部 核与辐射安全中心,北京 100082)

高放废物中往往含有毒性大、半衰期长的超铀核素(如239Pu等),其安全处置需实现核素与生物圈的长期有效隔离(时间长达1万年甚至更长)[1]。高放废物处置系统安全评价涉及多相流、热传导、溶质运移和化学反应等复杂的过程,这些过程通过线性和非线性关系紧密耦合[2]。这种耦合过程的数值模拟可分为3类:流体动力学模拟、地球化学模拟和反应性迁移模拟,其中反应性迁移模拟是最具实效性的模拟技术,它通常建立在前两类模拟研究基础之上[3-4]。

高放废物处置多重屏障系统主要包括工程屏障(废物包装材料、膨润土缓冲/回填材料等)和自然屏障(花岗岩、黏土岩、凝灰岩等)。在安全处置所考虑的长时间尺度下,随着工程屏障系统(EBS)各部分与周围环境介质(如地下水、气体)之间的物理-化学相互作用,这些屏障材料的性能和地下水的化学条件可能发生变化,进而影响核素在地下水中的迁移行为和处置设施的长期性能[4-5]。膨润土/黏土材料的性能在很大程度上取决于其特定的矿物组成和化学性质,如膨润土的膨胀特性与其组成矿物的性质和丰度有关,其中的蒙脱石矿物有利于材料吸水膨胀、封闭裂缝,而伊利石的形成则会影响材料的膨胀性,同时矿物的丰度和性质也会影响材料的渗透性、扩散性和吸附性[4]。此外,这些材料的机械变形、蒙脱石等矿物的溶胀和热变形及矿物溶解-沉淀等过程,也可能导致材料孔隙度发生演变,而在反应性迁移模拟中则较多关注于化学演化条件下矿物溶解-沉淀引起的孔隙增大或堵塞及其所导致的渗透率和流动特性变化[6-7]。

近30年来,科学界持续研究了膨润土屏障性能评价相关的前沿交叉科学问题[4,8-14],如多场耦合条件下材料的化学演化、多过程反应性迁移建模分析。Liu等[15]和Zheng等[16]针对黏土岩处置库和Kunigel-V1膨润土工程屏障建立了概念模型,利用TOUGHREACT-FLAC3D构建热-水-力-化学(THMC)多场耦合模型预测处置库近场演化,模拟结果表明,该膨润土演化过程中主要发生伊利石化作用(即蒙脱石溶解和伊利石形成)。Xie等[17]利用PHREEQC和MIN3P-THCm对加拿大地质处置多重屏障系统的长期地球化学演化研究表明,相对于混凝土,MX-80膨润土的反应性较低,矿物体积分数和孔隙度变化很小,然而在不同材料的界面附近会发生大量矿物沉淀、溶解,导致孔隙度减小甚至堵塞。李娜娜等[18]利用PHREEQC和GWB软件模拟给出了高庙子钠基膨润土与甘肃北山地下水相互作用后的矿物平衡体系和pH值演化规律(pH值最终维持在7.45或7.75左右)。Sun等[19]通过实验分析研究表明,高庙子膨润土的膨胀应变受到水溶液组成及浓度的显著影响,含有钾离子的溶液(特别是高碱性溶液)可显著降低该膨润土的膨胀性能。张虎元等[20-21]基于室内批式固-液静态反应的实验研究表明,高庙子膨润土矿物变化与碱溶液pH值的相关性较大,在新鲜水泥浸出液(pH>12.6)作用下,高庙子膨润土会发生蒙脱石矿物的溶解、伊利石/蒙脱石混层次生矿物的形成,从而降低高庙子膨润土的缓冲性能。总之,已有的研究显示:1) 各种类型膨润土材料的化学演化相当复杂,后续水力性能演化的原因(温度升高或地球化学反应)仍存在争议;2) 虽然场地尺度模拟实验中膨润土的伊利石化(即蒙脱石溶解、伊利石形成)有待证实,但并不意味着膨润土保持完好无损,事实上它们可能经历各种演化,包括Na-蒙脱石转化为Ca-蒙脱石、Fe-蒙脱石溶解等;3) 温度在膨润土化学演化中的作用尚不明显,加热条件下,大多数研究中观察到渗透率增加,但膨胀性能的变化尚无定论;4) 国外在膨润土演化的机理分析、模型预测等系统研究中具有优势,而国内更多是实验测试分析和技术总结跟进,多场-多过程动态耦合数值模拟研究成果相对欠缺。鉴于此,对于特定的处置系统,仍然需要结合特定的工程信息、材料特征等数据来进行模拟预测。

我国某区含Pu高放废物处置研究中,采用柯尔碱膨润土(若非特别说明,后续简称膨润土)作为回填/缓冲材料,前期已开展了膨润土矿物/化学成分、含水率、渗透性、比表面积、扩散等实验表征及核素迁移研究[22-26],然而,尚未定量分析蒙脱石、孔隙度、pH值等演化过程对膨润土膨胀性、渗透性及Pu反应性迁移的影响。为此,本文基于地下水入渗情景和概念模型,采用地球化学模拟手段,重点分析膨润土-地下水体系的长期演化及其潜在影响,为后续Pu随地下水迁移的模拟预测提供技术依据,同时丰富工程屏障材料性能演化机理分析与动态模拟研究成果。

1 材料与方法

1.1 概念模型

结合某区高放废物多重屏障处置系统设计,建立概念模型,如图1所示。将含Pu废物置于金属罐和混凝土组成的包装体中,四周回填厚度为1 m的柯尔碱膨润土,即包装体+膨润土层作为工程屏障,实现废物与外界环境的隔离。随着时间推移,包装体会逐渐腐蚀、失效,同时,地下水可能通过围岩裂隙通道渗入,与工程屏障层相接触,并与膨润土、废物体之间发生水-岩相互作用,最终可能导致Pu随水流而扩散迁移出屏障层。

图1 工程屏障系统概念模型Fig.1 Conceptual model of engineering barrier system

实际中,工程屏障材料时空演化涉及复杂的多场耦合过程,然而受工程设计、参数获取、模型不确定度等因素所限,若同时考虑这些过程,则难以建立模型。由于水动力过程或水动力-化学过程是数值建模中的重要环节,因此本文暂不考虑废物体腐蚀、膨润土侵蚀、气体产生和迁移等潜在事件和过程,重点分析地下水-膨润土体系的地球化学演化,即随着地下水的渗入,近场环境发生再饱和过程,水-岩相互作用后地下水成分特征、膨润土矿物组分也会发生演化,同时膨润土矿物的溶解-沉淀可能导致其孔隙度、渗透性发生变化。

1.2 地球化学模拟软件

采用TOUGHREACT模拟地下水-膨润土的长期演化。TOUGHREACT是在TOUGH2模拟软件的基础上添加了地球化学反应模块,可用于模拟一维、二维、三维地质介质中多相流体运动、热运移、溶质运移和化学反应过程,它考虑的主要机制包括:不同温度、压力、离子强度、水饱和度、pH/Eh值条件下,热、物理、化学和生物的耦合过程,反应动力学、生物降解、线性吸附和衰变;在局部平衡假设下的水溶液配位反应、氧化还原、阳离子交换以及表面配位反应;局部平衡或动力学条件下的矿物溶解-沉淀反应及其对地质介质孔隙度、渗透率的影响等[27]。TOUGHREACT提供了7类流体属性模块,用于模拟多相流体渗流和地球化学迁移相关的多种问题,如废物地下处置、CO2地质封存等,该软件最新发展包括耦合地质力学效应和反应性迁移的MPI-OpenMP并行化[3-4]。

1.3 计算模型构建及数据参数

1) 计算模型建立

基于工程屏障概念模型,利用TOUGHREACT建立二维剖面模型,如图2如示。模拟区域为5 m×5 m,采用矩形网格剖分并细化关注区域,共剖分4 900个网格;设定水流从上部入渗,在膨润土层与废物体界面处A点、膨润土中部B点、膨润土层与花岗岩界面处C点分别设置3个观测点分析水饱和过程和地球化学演化,在膨润土层底部D点处设置观测点分析Pu扩散迁移的浓度变化。

图2 二维剖面模型示意图Fig.2 Schematic diagram of 2D profile model

2) 环境初始条件及数据

柯尔碱膨润土和花岗岩的矿物组成如表1所列,研究区地下水的参数如表2所列。利用表1、2数据进行平衡模拟,得到模型中初始地下水和膨润土孔隙水组成,并为后续反应性迁移模拟提供初始条件。

表1 柯尔碱膨润土和花岗岩的矿物组成Table 1 Mineral composition of Kerjian bentonite and granite

表2 地下水参数Table 2 Groundwater parameter

数值模拟中主要关注膨润土的再饱和过程、地球化学演化及Pu的迁移规律,模拟的压力条件为1个大气压、温度为25 ℃,模型中参数取值列于表3。水-气两相流模拟中,相对渗透率krl和毛细压力pcap分别采用Corey公式和van Genuchten公式描述,结合文献数据给定相关参数[28-29]。在膨润土饱和后的渗流迁移模拟中,废物体材料设定为混凝土。Jacques等[30]的研究表明,随着时间的演变,受侵蚀降解等影响,混凝土的渗透率变化较大,约在10-13～10-19m2量级,据此本文中废物体渗透率取典型值10-17m2。在膨润土低渗透性条件下,渗流和迁移以分子扩散为主导,根据文献[26,31-32]分析,水流与Pu的分子扩散系数分别取为1.0×10-9m2/s和5.0×10-10m2/s。

表3 模型参数Table 3 Model parameter

TOUGHREACT中的地球化学模型考虑了各种平衡和动力学化学反应,所用热力学数据库是基于EQ3/6数据库修改而来,EQ3/6计算机代码曾用于美国尤卡山高放废物处置系统性能评价等多项研究,使用SUPCRT92软件对热力学数据进行评估和更新。TOUGHREACT中矿物的溶解-沉淀动力学数据已得到不同学者的反复验证和应用,并在文献[16,27]中有所介绍。

2 结果与讨论

2.1 水饱和过程演化

膨润土层中3个观测点处的液相饱和度演化过程如图3所示。由图3可见,在水流由饱和花岗岩区域向非饱和膨润土层渗透过程中,C点基本处于水饱和状态,B点先于A点达到饱和状态;膨润土完全饱和约需20 a,这与国际上对Kunigel-V1膨润土的水力演化研究结果[15-16]相一致。

图3 柯尔碱膨润土层水饱和度随时间的演化Fig.3 Temporal evolution of water saturation in Kerjian bentonite layer

2.2 蒙脱石的伊利石化分析

膨胀性是各种膨润土、黏土类缓冲/回填材料的关键属性,在降低EBS材料的渗透性、提供EBS对周围区域的力学支撑以及开挖破损区自密封方面起着重要作用。这些膨胀能力与其中蒙脱石含量直接相关,蒙脱石溶解而转化为伊利石的过程非常重要,因为它会导致蒙脱石的损失,进而使得材料的膨胀能力降低[15]。

柯尔碱膨润土长期演化的模拟结果示于图4。由图4可见,膨润土层中产生的伊利石最大体积分数仅为10-10量级,表明伊利石化作用并不明显,可认为膨润土层的膨胀性几乎没有变化,这有利于保持膨润土的缓冲、阻滞、封闭等屏障性能。伊利石化作用很弱,其原因可能在于膨润土和地下水中K的含量较低、蒙脱石的溶解份额很少等[15]。

图4 柯尔碱膨润土中伊利石体积分数随时间的演化Fig.4 Temporal evolution of illite volume fraction in Kerjian bentonite

2.3 柯尔碱膨润土的孔隙度和渗透率变化

处置库关闭后,围岩地下水将重新饱和工程屏障中的膨润土,这一再饱和过程会引起矿物的溶解和沉淀。由于矿物溶解-沉淀反应导致的孔隙度变化如图5所示,采用Carman-Kozeny公式计算的渗透率变化如图6所示。由图5可见,1 000 a时孔隙度变化最大,由初始值0.4先增加至约0.403,再减小至约0.397,这一变化量很小、不足初始值的2%,反映出孔隙度随时间的变化并不显著;不同位置处孔隙度的变化有所不同,C点孔隙度变化最大,这可能与不同位置处的水-岩相互作用强度不同相关,即入流边界处较强烈的水-岩反应、界面处的边界效应等因素[16],使得靠近花岗岩界面处的孔隙度变化更明显。由图6可知,由于Carman-Kozeny公式中渗透率与孔隙度正相关,因此渗透率变化也不显著,其变化规律与孔隙度变化类似。总体上,在再饱和过程的1 000 a中,孔隙度和渗透率均没有明显变化。通常认为膨润土完全饱和后,水-岩相互作用接近平衡状态,除非发生地震等低概率事件[16],否则化学演化将非常缓慢,因此,可认为化学反应不会导致膨润土孔隙度及渗透性能发生显著改变。

图5 孔隙度随时间的演化和空间分布Fig.5 Temporal evolution and spatial profile of porosity

图6 渗透率随时间的演化和空间分布Fig.6 Temporal evolution and spatial profile of permeability

2.4 地下水pH值的变化

3个观测点处地下水pH值变化如图7所示。由图7可见,A点处pH值随时间先增加后减小,最大值约为10.3,膨润土完全饱和后pH值又逐渐减小至约9.5;B点处pH值变化与A点处类似,不同之处在于pH值从最大值10.3减小至9.5的过程较为平缓;C点处pH值先增加至约9.4,之后随着膨润土完全饱和而迅速降低至初始值,其原因可能在于膨润土经历水合作用、与花岗岩围岩的相互作用、边界效应等[16],使得膨润土层与花岗岩界面附近区域有最高的反应性,导致该点pH值变化有所不同。总之,通过水-岩相互作用,地下水的pH值呈上升趋势并维持在偏碱性氛围。由于Pu在这种条件下通常以难溶、难迁移的形式存在,因此,膨润土的这一pH值缓冲性能可提供偏碱性的水文地质环境,有利于实现放射性废物的长期安全隔离。

图7 pH随时间的演化Fig.7 Temporal evolution of pH

2.5 膨润土主要矿物的演化

膨润土中主要矿物在地下水中溶解-沉淀行为的模拟结果示于图8、9。由图8、9可看出:膨润土中蒙脱石、长石、石英的体积分数变化很小,分别约为0.03%、-0.05%、-0.001%,表明蒙脱石有少量沉淀,长石和石英少量发生溶解;方解石和白云石的体积分数变化相对较大,分别约为0.8%和-0.6%,表明在地下水-膨润土体系中,由于水中离子浓度与矿物的不平衡,导致白云石的溶解和方解石的形成;C点处矿物的体积分数变化较大,可能是由于水合作用、水-岩相互作用、边界效应等因素,导致不同材料界面处溶解-沉淀反应性较强[16]。通常,质子参与许多反应,其最重要的反应是蒙脱石(montmorillonite)、伊利石(illite)、白云石(dolomite)和方解石(calcite)等矿物的溶解-沉淀,以及具有表面吸附和离子交换过程的质子化反应,去白云石化反应(式(1))生成的方解石对伊利石化(式(2))的影响较小,但对控制近场环境的酸碱条件很重要[15-16],本文结果也证实了这一点,随着去白云石化反应的进行,H+浓度减小、pH值增大,这也解释了上文所述pH值的演化特征。蒙脱石溶解是伊利石化的重要机制,由于地下水和膨润土中Na的浓度较高,从而抑制了蒙脱石的溶解和伊利石的形成。

图8 观测点A、B、C处膨润土主要矿物的时间演化Fig.8 Temporal evolution of bentonite minerals at observation point of A, B and C

图9 膨润土主要矿物的空间分布Fig.9 Spatial distribution of primary minerals in bentonite

(1)

0.6K+=illite+0.26H2O+0.08Mg2++

0.33Na++0.5SiO2(aq)

(2)

2.6 地下水主要离子浓度的演化

图10 液相中K+和随时间的演化Fig.10 Temporal evolution of K+ and in aqueous

3 结论

1) 地球化学模拟表明,柯尔碱膨润土层大约需20 a达到完全饱和;膨润土演化过程中蒙脱石含量几乎没有溶解损失,这有利于维持膨润土的膨胀性能和处置工程的安全性;矿物溶解-沉淀反应引起的孔隙度变化并不明显,最大变化量不到初始值的2%,渗透率变化也有类似特征;地下水pH值呈上升趋势并维持在偏碱性的8.1～10.3范围内,膨润土的这一pH值缓冲性能将有利于阻滞核素向外迁移、实现Pu的长期安全隔离。

感谢北京师范大学王金生教授、吉林大学许天福教授和田海龙教授在数值模拟方面给予的支持。