放射性核素在包气带中的迁移模拟

张 雪， 王小元， 胡继华

(中国核电工程有限公司， 北京 100840)

放射性核素通过环境介质对人类及生物产生辐射累积效应，地下水环境是主要的途径之一。资料显示，由于历史原因，美国、英国、日本、西班牙等多个国家的部分核设施周边地下水中发现了放射性核素，不仅给当地地下水环境带来了较大影响，同时也给这些核设施后期退役工作带来了诸多难题，凸显了核设施的地下水环境影响评价的重要性。对于放射性核素来说，包气带对放射性核素的阻滞作用非常明显，李书绅等[1]的研究表明，在自然条件下，可吸附核素的迁移距离非常短，因此，研究放射性核素在包气带中迁移是必要的。

世界上最早的核素迁移研究工作是1958年进行的放射性废物的吸附研究[2]。20世纪70—80年代，美国、英国、日本、加拿大、中国、苏联等国家开展了以下方面的研究工作：①进行包气带水分运移现场试验研究，如美国在新墨西哥州的试验[3-4]；②在实验室进行核素在原状土柱中的迁移的模拟试验研究，如日本在沙土和壤土的原状土柱中做过60Co、85Sr和137Cs的迁移模拟实验[5]；③为了了解核素迁移的实际情况，美国、英国等国家做过野外核素迁移的示踪试验[6-9]。中国开展的比较著名的研究为1988年1月至1993年1月中国辐射防护研究院与日本原子力研究所合作进行的60Co、85Sr、134Cs和137Cs在非饱和黄土中迁移的野外试验、实验室模拟实验和非饱和水分运移规律现场观测等[1]。中国学者主要开展了：①核素在包气带中迁移相关参数测定及其影响因素的分析研究[10-11]；②核素在包气带中迁移的野外实验与室内实验研究[12-15]；③核素在包气带中迁移的数值模拟研究[16-18]。

在充分吸收中外已有研究成果的基础上，同时借鉴非放射性核素迁移模拟相关研究内容[19]，从地下水环评(包气带-饱和带耦合模拟)角度出发，通过重点分析潜水面处(包气带-饱和带交界)水流及核素迁移特征，说明包气带介质对核素迁移的迟滞作用，进而阐明地下水环评中评价核素在包气带中迁移的必要性，同时通过分析包气带厚度和泄漏时间等对核素迁移影响，为辐射防护、应急措施等的制定提供借鉴和参考。

放射性核素在包气带中的迁移主要以垂直向为主，且与地下水环评的评价范围相比，水平向的运动几乎可以忽略，因此，现主要通过野外试验获得包气带参数，利用Hydrus-1D软件对包气带中放射性核素迁移进行模拟，充分掌握放射性核素在包气带中的迁移扩散规律，以期为更加准确地开展地下水环境影响评价工作提供技术支持，同时为核设施选址、辐射防护设计等提供技术支撑。

1 研究区概况

研究区位于甘肃省酒泉市，属大陆性干旱气候，降水少，蒸发强烈。地势北高南低，北部为基岩山地，海拔1 250.0～1 440.0 m，南部为碎石戈壁，海拔一般为1 229.0～1 275.0 m。夏季干热而短促，冬季寒冷而漫长。多年平均降雨量62.9 mm，多年平均蒸发量为2 338.9 mm研究区地层从上至下依次为全新统冲洪积中砂层、上更新统冲洪积砾砂混角砾层、海西中期黑云母花岗闪长岩和黑云母花岗岩。地下水水位埋深在2.80～20.06 m，水位标高在1 205.25～1 227.8 m。

2 原位土-水特征曲线测试及包气带渗水试验

2.1 原位土-水特征曲线测试

为了获取包气带土-水特征曲线，为数值模拟提供参数，开展了原位土-水特征曲线测试工作。

2.1.1 试验原理

采用负压计法，通过注水产生地表积水条件模拟放射性液体泄漏，包气带介质达到饱和后停止注水并清除地表积水，通过监测包气带脱水过程的含水量θ和吸力值h获取土-水特征曲线。

2.1.2 试验步骤

(1)开挖长约3 m、宽约1 m、深4 m的区域作为试验地块。

(2)于深度0.2、0.4、0.7、1.0、1.5、2.0、2.5、3.0、4.0 m处布设土壤水分测定仪(TDR)传感器并在各深度采集适量包气带介质进行颗粒分析，同时于深度0.2、0.4、1.0、2.0、3.0、4.0 m处布设水势传感器，布设完成后将介质回填，并尽可能将其恢复至开挖前状态。

(3)在地表嵌入高度30 cm、直径1 m的铁环，确保土壤水分测定仪、水势传感器包含在铁环所包围的垂直柱体中。

(4)向铁环中注水使其保持地表积水状态24 h，持续监测土壤含水量和土壤水吸力值，直至土壤含水量降至注水前状态，试验照片如图1所示。

2.1.3 试验结果

(1)颗粒分析。由颗粒分析试验结果(表1)可知，各深度包气带介质粒径>2 mm颗粒质量在总质量25%～50%，粒径大于0.075 mm的颗粒质量超过总质量的50%，根据《岩土工程勘察规范》(GB50021—2001)(2009年版)，将包气带介质定名为砾砂。

表1 包气带介质颗粒分析试验结果Table 1 Particle analysis test results of vadose medium

(2)土-水特征曲线测试结果。土-水特征曲线选取脱湿过程的数据进行分析，如图2所示，由图2可知，各深度包气带介质土-水特征曲线差异不大，结合颗粒分析试验结果，可将包气带介质概化为均质的砾砂。

图2 各深度介质土-水特征曲线Fig.2 Soil-water characteristic curve of medium

2.2 现场包气带渗水试验

为了获取包气带渗透系数，开展了包气带渗水试验，包括单环法、双环法及IAM(the inverse auger method)法，结果如表2所示。

表2 包气带渗水试验结果Table 2 Permeability test results of vadose zone

3 放射性核素迁移数值模拟研究

3.1 源项条件假设

假设放射性液体在地表发生瞬时泄漏，选择不被吸附的H-3和被吸附的Sr-90作为研究对象，H-3的浓度假定为1.2×109Bq/L，Sr-90浓度为1.78×10-1Bq/L。

3.2 水流与核素迁移概念模型

研究区气候干旱，降水稀少，蒸发强烈，包气带厚度在2.80～20.06 m，试验场地包气带厚度约10 m，根据颗分试验结果，包气带岩性为砾砂，建立概念模型时，对核素在包气带中迁移做如下假设：①包气带中水流和核素为垂向一维流，核素的释放源在地表；②包气带介质对核素的吸附是可逆的，为一级线性反应；③Sr-90的衰变子核素为Y-90，半衰期为2.67 d，不考虑Sr-90衰变子核素对环境的影响；④忽略温度和气体的影响。

根据上述假设，以地表为零基准面，建立向上为正的概念模型，模拟区域为-H≤z≤0，其中z为坐标，H为包气带厚度，如图3所示。

图3 水文地质概念模型Fig.3 Hydrogeology conceptual model

3.3 水流与核素迁移数学模型

3.3.1 水流模型

(1)

式(1)中：θ为体积含水率；h为土壤水压力水头；z为垂向坐标(向上为正)；t为时间；K(h)为非饱和水力传导度；Tp(t)为表土入渗强度，包括降水和泄漏；Tε(t)为表土蒸发强度；qs为地表水分通量；hb为各时刻地表压力水头；h0为初始时刻各点压力水头。

3.3.2 核素迁移模型

(2)

3.4 定解条件

3.4.1 水流模型的初始和边界条件

(1)水流模型初始条件。本次模拟计算中，水分迁移模型的初始值由压力水头给出：

h(z,0)=h0(z)， -H≤z≤0

(3)

假设初始条件下水流处于稳定状态，设潜水面处压力水头为0，向上递减，地表处压力水头为-10 m。

(2)水流模型边界条件。上边界为变流量大气边界：

z=0,t≥0

(4)

下边界选择潜水面为变水头边界：

h(z,t)=hb，z=-H，t≥0

(5)

3.4.2 核素迁移模型的初始和边界条件

(1)核素迁移初始条件。初始时刻选在放射性液体刚刚从地表泄漏，初始浓度设为0。

(2)核素迁移边界条件。下边界为零浓度梯度边界。上边界为定浓度通量边界：

(6)

当边界流量为土壤水蒸散量/降水量时，C=0；当边界流量为泄漏量时，C取泄漏放射性液体中放射性核素浓度。

3.5 模型主要参数的确定

参数主要包括土-水特征曲线参数、分配系数Kd、衰变常数λ、弥散度αL、分子扩散系数D等。

(1)土-水特征曲线参数

饱和含水量θs、残余含水量θr、形状系数n、形状系数α通过原位土-水特征曲线，利用MATLAB软件，使用Van Genuchten公式拟合，Ks取现场包气带渗水试验平均值，即K=7.33 m/d，土-水特征拟合曲线如图4所示，取值如表3所示。

图4 土-水特征拟合曲线Fig.4 Soil-water characteristic curve

表3 土-水特征曲线参数取值Table 3 Soil-water characteristic curve parameters

(2)分配系数Kd。H-3不被吸附，Kd为0，Sr-90的分配系数Kd选用《Groundwater Monitoring Guidance for Nuclear Power Plants》[20]推荐的经验值，Kd=15 cm3/g。

(3)弥散度αL。弥散度选用多孔示踪试验值αL=0.77 m。

(4)衰变常数λ。H-3和Sr-90的半衰期分别为12.3 a和29.1 a，衰变常数分别为1.54×10-4/d和6.53×10-5/d。

3.6 模拟方案

包气带厚度为10 m，假设泄漏时长为3 d，模拟至包气带对饱和带停止补给。

3.7 模拟结果及分析

3.7.1 潜水面处水流速度随时间变化

潜水面处水流速度随时间变化如图5所示。

图5 潜水面处水流速度随时间变化曲线Fig.5 Graphs of the water flow velocity change with time at the water table

结果显示：①初始时刻(刚发生泄漏)，潜水面处的水流速度极小(趋于0)；②第0.36天(8.64 h)时水流速度急速增大，第0.37天(8.88 h)时已增大至饱和渗透系数7.33 m/d；③第3天时，泄漏停止，由于没有地表泄漏补给，水流在重力和基质吸力共同作用下向下迁移，随着包气带中水体向饱和带迁移，含水率降低，基质吸力迅速增大，潜水面处水流速度快速下降，1 d内水流速度由7.33 m/d降至0.286 m/d，随后水流速度缓慢降低，包气带缓慢脱水，第10天时降为3.83×10-2m/d，第1年时降至3.31×10-4m/d，直至第120年降为0，包气带对饱和带无补给。

3.7.2 水面处H-3的迁移规律

(1)潜水面处H-3浓度随时间变化。潜水面处H-3浓度随时间变化如图6所示。

图6 潜水面处H-3浓度随时间变化曲线Fig.6 Graphs of the H-3 concentration change with time at the water table

结果显示：①H-3在对流作用下流入饱和地下水，第0.35天(8.4 h)时，H-3随水流迁移到潜水面处，浓度极低为1.53×10-33Bq/L；②随后H-3浓度迅速增大，由于H-3不发生吸附，第1.3天时，潜水面处H-3浓度达到峰值1.2×109Bq/L，与地表泄漏液浓度相同；③随后在衰变的作用下H-3浓度缓慢降低，在包气带对饱和带停止补给时(120 a)降为1.29×106Bq/L。

(2)地表及潜水面处H-3通量随时间变化。为了解包气带对H-3迁移的延迟作用，对地表(无包气带的迟滞作用)及潜水面处H-3通量随时间变化进行对比分析，如图7所示。

图7 地表及潜水面处H-3通量随时间变化曲线Fig.7 Graphs of the H-3 flux change with time at the ground surface and water table

结果显示：①地表处H-3通量全部集中在泄漏发生的3 d内：初始时刻(t=0)，地表处即存在H-3通量，浓度值在8.8×1012Bq/(m2·d)上下波动；直到第0.37天(8.88 h)浓度通量达到稳定值8.8×1012Bq/(m2·d)；第3天时，渗漏补给停止，浓度通量为0；②潜水面处H-3通量存在滞后，且时间上更加分散：第0.35天(8.4 h)时，潜水面处开始有H-3由包气带流入饱和带，浓度通量为1.3×10-30Bq/(m2·d)；第1.3天时浓度通量达到峰值8.8×1012Bq/(m2·d)，并趋于稳定；3 d后由于水流通量快速降低，潜水面处H-3的通量随之快速降低，第4天时已降至3.46×1011Bq/(m2·d)，随后缓慢降低，直至对饱和带无H-3的补给。

(3)地表及潜水面处H-3累积通量随时间变化。

为了解包气带对H-3的滞留，对地表及潜水面处H-3累积通量随时间变化进行对比分析，如图8所示。

图8 地表及潜水面处H-3累积通量随时间变化曲线Fig.8 Graphs of the H-3 accumulated flux change with time at the ground surface and water table

结果显示：①地表处H-3的累积通量在前3天呈线性增大，第3天时达到峰值，随后保持恒定值2.64×1013Bq/m2；②潜水面处H-3的累积通量第0.35天(8.4 h)时开始呈线性增大，第3天时累积通量达到2.25×1013Bq/m2，随后缓慢增大，第120年达到2.54×1013Bq/m2；③即使H-3不发生吸附，但由于包气带介质对放射性液体存在基质吸力，仍有部分H-3被滞留在包气带中，第3天时潜水面处的H-3累积通量(2.25×1013Bq/m2)比地表(2.64×1013Bq/m2)减少了3.9×1012Bq/m2(14.8%)，第120年时潜水面处的H-3累积通量(2.54×1013Bq/m2)比地表处减少了1.0×1012Bq/m2(3.8%)，其中，3 d～120 a时间段内H-3由包气带渗入到饱和地下水(累积通量2.9×1012Bq/m2)的过程极其缓慢，地下水对其具有极强的稀释作用，使其在地下水中迁移过程中浓度大大降低。

3.7.3 潜水面处Sr-90的迁移规律

(1)潜水面处Sr-90浓度随时间变化。潜水面处Sr-90浓度随时间变化如图9所示。

图9 潜水面处Sr-90浓度随时间变化曲线Fig.9 Graphs of the Sr-90 concentration change with time at the water table

结果显示：①由于吸附作用产生的延迟现象，使Sr-90流入饱和带的时间比不发生吸附的H-3长，第0.9天(21.6 h)时潜水面处才有Sr-90出现，浓度极低为1.33×10-33Bq/L；②随后浓度逐渐增大，第3天时浓度增大至8.05×10-12Bq/L；③3 d后地表无泄漏补给，由于吸附作用产生的迟滞效应，被吸附的Sr-90缓慢释放，潜水面处Sr-90的浓度仍在缓慢增加，第10天时增至3.84×10-11Bq/L，第10年时增至9.21×10-11Bq/L，第65年时增至9.71×10-11Bq/L，随后Sr-90浓度趋于稳定。不同于H-3，由于Sr-90发生吸附，其峰值浓度(9.71×10-11Bq/L)远小于泄漏的放射性液体中Sr-90浓度(1.78×10-1Bq/L)。

(2)地表及潜水面处Sr-90通量随时间变化。地表及潜水面处Sr-90通量随时间变化如图10所示。

图10 地表及潜水面处Sr-90通量随时间变化曲线Fig.10 Graphs of the Sr-90 flux change with time at the ground surface and water table

结果显示：①初始时刻地表即存在Sr-90通量，浓度值在1.31×103Bq/m2·d上下波动；第0.37天(8.88 h)浓度通量达到稳定值1.31×103Bq/(m2·d)；第3天渗漏补给停止，通量为0；②第0.9天(21.6 h)时，潜水面处开始有Sr-90由包气带流入饱和带，浓度通量为4.86×10-30Bq/(m2·d)；随后逐渐增大，第3天时达峰值5.79×10-8Bq/(m2·d)，由于吸附作用，其值远小于地表处通量峰值[1.31×103Bq/(m2·d)]；随后由于水流通量快速降低，Sr-90通量随之降低，直至第120年通量为0。

(3)地表及潜水面处Sr-90累积通量随时间变化。地表及潜水面处Sr-90累积通量随时间变化如图11所示。

图11 地表及潜水面处Sr-90累积通量随时间变化曲线Fig.11 Graphs of the Sr-90 accumulated flux change with time at the ground surface and water table

结果显示：①地表处Sr-90的累积通量在前3天呈线性增大，第3天时达到峰值，随后保持恒定值3.91×103Bq/m2；②潜水面处Sr-90的累积通量0.9 d(21.6 h)时开始快速增大，第3天时累积通量达到8.05×10-12Bq/m2，随后缓慢增大，第120年时达到9.64×10-11Bq/m2。可知潜水面处累积通量(3 d时为8.05×10-12Bq/m2，第120年时为9.64×10-11Bq/m2)远小于地表的累积通量(3 d达到峰值3.91×103Bq/m2)，包气带对发生吸附的Sr-90具有极好的缓冲和滞留作用，几乎将所有的Sr-90滞留在包气带中，只有极少量的Sr-90入渗补给到饱和地下水环境。

3.8 包气带厚度和泄漏时长对核素在包气带中迁移的影响模拟研究

3.8.1 模拟方案

为了研究包气带厚度和泄漏时长对核素在包气带中迁移的影响，结合场地实际条件，共设计2组模拟方案，如表4所示。

表4 模拟方案设计表Table 4 Design table of simulation plan

3.8.2 模拟结果与影响因素分析

(1)包气带厚度对核素在包气带中迁移的影响。包气带厚度对核素在包气带中迁移影响模拟结果如表5所示。

表5 包气带厚度对核素在包气带中迁移影响模拟结果Table 5 Simulation results of the effect of vadose thickness on nuclides migration

由表5可知：包气带厚度越大对核素的迟滞作用越强；对于不被吸附的H-3，包气带厚度主要影响其到达潜水面的时间，对补给到饱和地下水的总量影响较小，包气带厚度分别为3、10、20 m时，H-3到达潜水面处的时间分别为0.08、0.35、0.99 d，潜水面处的累积通量差异不大，分别为2.59×1013、2.50×1013、2.33×1013Bq/m2；对于被吸附的Sr-90，包气带厚度对其到达潜水面的时间和补给饱和地下水总量的影响均较大。包气带厚度分别为3、10、20 m时，Sr-90到达潜水面的时间比不被吸附的H-3迁移时间更长，分别为0.11、0.9、2.99 d；潜水面处的累积通量呈数量级递减，分别为3、1.03×10-7、9.17×10-22Bq/m2。

(2)泄漏时长对核素在包气带中迁移的影响。泄漏时长对核素在包气带中迁移影响模拟结果如表6所示。

表6 泄漏时长对核素在包气带中迁移影响模拟结果Table 6 Simulation results of the effect of leakage duration on nuclides migration

由表6可知，泄漏持续时间越长，流入饱和地下水中的核素总量越大；对于不被吸附的H-3，泄漏时长分别为1、2、3 d时，潜水面处累积通量分别为7.79×1012、1.66×1013、2.50×1013Bq/m2，基本符合线性关系；对于被吸附的Sr-90，泄漏时长分别为1、2、3 d时，累积通量相差极大，呈数量级递增，分别为2.65×10-22、3.40×10-12、1.03×10-7Bq/m2。

根据模拟结果，包气带厚度越大，对核素迁移的迟滞作用越强，可为选址提供借鉴和参考；同时可为应急响应提供技术支持：应加强监测，一旦发生放射性液体的泄漏，应立即采取措施切断补给源来缩短泄漏持续时间，并可采取降低地下水位增大包气带厚度的应急措施来降低放射性核素对地下水的污染。

4 结论

为了能够为更加准确地开展地下水环境影响评价工作提供技术支持，分析放射性核素在包气带中的迁移特征及包气带厚度及泄漏持续时间对核素迁移的影响，选用Hydrus-1D软件，利用数值模拟方法对H-3和Sr-90两种放射性核素在包气带中的迁移行为进行了研究。得到如下结论。

包气带对核素的迁移存在明显的迟滞作用(延迟和滞留)，主要表现在核素迁移至饱和地下水的时间、浓度及通量三个方面，且对被吸附的核素的影响大于不被吸附的核素。

(1)由于水流迁移至饱和带需要一定时间(0.37 d)，因此核素在对流的作用下迁移至饱和带也存在短暂的延迟，且被吸附的Sr-90到达饱和带的时间(0.9 d)比不发生吸附的H-3(0.35 d)长；泄漏停止后，由于包气带介质对水体的基质吸力，使赋存于水体中的核素滞留于包气带介质中不能瞬时迁移至地下水，而是缓慢迁移至饱和地下水。

(2)由于H-3不被吸附，流入饱和带的H-3的峰值浓度与泄漏的H-3浓度相同，但由于包气带介质对H-3的滞留作用和H-3的衰变作用，流入饱和带的H-3浓度会发生缓慢的降低；而对于被吸附的Sr-90，除了包气带介质通过影响水流迁移导致的Sr-90迁移的迟滞，由于吸附作用的存在使其到达饱和带的浓度峰值(9.71×10-11Bq/L)远小于地表泄漏液的核素浓度值(0.178 Bq/L)。

(3)由于包气带介质对含H-3液体的基质吸力，部分H-3被滞留在包气带中，第3天时潜水面处的H-3累积通量比地表的减少了3.9×1012Bq/m2(14.8%)，第120年时潜水面处的H-3累积通量比地表处减少了1.0×1012Bq/m2(3.8%)，其中，第3天至第120年时间段内H-3由包气带渗入到饱和地下水(累积通量2.9×1012Bq/m2)的过程极其缓慢，地下水对其具有极强的稀释作用，使地下水中H-3浓度大大降低；对于发生吸附的Sr-90，包气带介质的基质吸力和吸附的共同作用下，且潜水面处Sr-90的累积通量(3 d时为8.05×10-12Bq/m2，120 a时为9.64×10-11Bq/m2)远远小于地表Sr-90的累积通量(3d达到峰值3.91×103Bq/m2)，包气带对Sr-90具有极好的缓冲和滞留作用，几乎将所有的Sr-90滞留在包气带中，只有极少量的Sr-90入渗补给到饱和地下水环境。可知以往不考虑包气带对放射性核素(特别是发生吸附的核素)迁移的迟滞作用进行环境影响评价的结果过于保守，可补充开展核素在包气带中的迁移模拟研究。

包气带厚度和泄漏时长对核素在包气带中迁移影响较大，主要包括以下几个方面。

(1)包气带厚度越大对核素的迟滞作用越强。对于不被吸附的H-3，随着包气带厚度增大，其到达潜水面的时间呈非线性增大，补给到饱和地下水的总量基本不变。对于被吸附的Sr-90，随着包气带厚度增大，其到达潜水面的时间也呈非线性增大，且比不发生吸附的H-3时间长，补给饱和地下水总量呈数量级递减。

(2)泄漏持续时间越长，流入饱和地下水中的核素总量越大。对于不被吸附的H-3，泄漏时长与流入饱和地下水中的核素总量基本符合线性关系；对于发生吸附的Sr-90，随着泄漏时长的增大，累积通量呈数量级递增。

根据模拟结果，包气带厚度越大，对核素迁移的迟滞作用越强，可以为选址提供借鉴和参考；同时可为应急响应提供技术支持：应加强监测，一旦发生放射性液体泄漏，应立即采取措施切断补给源来缩短泄漏持续时间，并可采取降低地下水位增大包气带厚度的应急措施来降低放射性核素对地下水的污染。