动态核事故农业后果评估模型开发及应用研究

2021-05-25崔慧玲黄莎吕明华赵多新张俊芳李云鹏

农业环境科学学报 2021年4期

崔慧玲，黄莎，吕明华，赵多新，张俊芳，李云鹏

（中国辐射防护研究院，太原030006）

为了阻止和减小核设施对农业环境的潜在核事故后果影响，建立旨在保护食物和农产品的核事故应急反应系统非常必要。应急反应系统包含大气扩散模型、农业后果评估模型和对应的农业措施。1986年切尔诺贝利核事故由于缺乏应急反应系统，导致事故后食物里不合理的131I活度水平和之后高的甲状腺剂量及甲状腺癌记录[1]。2011 年日本福岛核事故使邻近土壤遭受了严重威胁，但是由于吸取了切尔诺贝利事故的经验，通过应急反应系统及时采取了农业应对措施，有效降低了人的吸收剂量[2]。国际上关于影响食物和农业的核应急系统开发始于2013 年，在FAO/IAEA 组织框架下进行[3]（FAO：Food and Agricul⁃ture Organization，IAEA：International Atomic Energy Agency），该计划同时开展的还有国际间关于农产品的采样分析及数据采集工作。农业的核事故应急反应，既要采取及时、恰当的农业措施，同时又不能造成大量的污染废物，这就需要农业评估模型提供可靠的评估结果。我国现有的核事故应急反应系统只有简单的食物链模型，单独针对农业环境的核事故后果评估模型开发研究较少。20世纪80年代食物链模型随着核事故后果评估系统的开发而开发，计算了核素从生态环境进入粮食作物，再通过食入途径形成食入剂量的过程，计算分为动态模型和平衡模型。平衡模型一般假设生态系统处于平衡静止的状态，核素迁移过程不存在时空差异。事实上放射性核素在生态环境中的迁移是一个动态过程，生态环境变化和生物体自身新陈代谢都会影响其最终的活度[4-6]。动态模型计算结构复杂，需要很多难以获得的动态生态参数和较高的计算性能，难以在应急情况下应用。因而计算简洁快速的平衡模型成为核事故评估系统的主要计算方法。例如英国和德国联合开发的COSYMA[7]和美国的MACCS[8]中的食物链模型都是平衡模型，只有欧洲JRODOS[9]系统中的食物链模型采用了动态模型。现在比较成熟、完全的动态食物链模型是ECOSYS 模型[10]，它是一个针对德国中部农业环境开发的动态食物链模型。由于模型应用的复杂性，ECOSYS 在其他国家的应用研究较少，主要研究工作集中在北欧地区。Liland 等[11]用包含ECOSYS 的ARGOS 和其他模式组成模式链，用于挪威塞拉菲尔德后处理厂的核事故应急评估，计算了事故条件下的食入剂量。Slavik等[12]对ECOSYS 在斯洛伐克共和国的本地化及应用做了研究，在分析收集数据的基础上，建立了新物种“苜蓿”的参数，计算了切尔诺贝利事故后各食物的活度并与检测数据进行对比，结果显示模型结果与观测数据有较好的一致性。

随着放射性核素生态环境迁移研究的深入和计算机性能的大幅提高，能更真实反映事故后核素迁移过程的动态模型逐渐成为事故后果评估的主要方法。本文基于动态模型理论和食物链模型ECOSYS 的计算方法，建立了基于我国农业生态环境的动态核事故农业后果评估模型AgrRAM，为核事故应急反应系统提供基础数据。AgrRAM 是应急反应系统的一个模块，它以大气扩散模型的结果为输入，模型从沉积计算开始完整模拟了核事故后放射性尘埃在农业生态环境中的沉积、迁移及最终形成农产品活度的过程，考虑了生态环境变化及生物体自身的生长过程，计算了随时间变化的典型农产品活度的时间、空间分布。用AgrRAM 模型对处于我国重要夏粮产地（黄淮海区域）的田湾核电厂进行了假想事故后果计算并对结果进行了讨论，说明新建模型AgrRAM 具有动态模拟计算的功能，真实反映了核事故后农产品的活度分布，定量计算了核事故后农业辐射效应。由于缺乏核事故后农作物的实际监测数据，AgrRAM 模型通过与保守的平衡模型对比进行了模型可靠性验证。

1 AgrRAM的模型结构及主要计算过程

AgrRAM 模型计算主要分两个部分：农业生态环境沉积计算和核素迁移计算。模型的结构及计算流程如图1 所示，计算始于大气扩散模型的结果也可用环境监测数据作为输入。模型的输入数据有网格化的时间积分空气活度浓度和降雨量等，输出有地面沉积活度分布和农产品比活度空间分布及时间序列。

1.1 沉积计算

放射性尘埃在农业生态环境中的沉积分两部分，一部分是植物叶面的直接沉积，另一部分是土壤的沉积。植物叶面干沉积计算考虑了大气阻尼和植物冠层阻尼，在湿沉积计算中考虑了植物的截留作用。植物沉积计算主要过程如下：

式中：Ai是对于i 植物的总沉积量，Bq·m-2；Adi是对于i植物的干沉积，Bq·m-2；Aw是对于i 植物的湿沉积量，Bq·m-2；fwi是湿沉积截留因子；Cair是空气积分浓度，Bq·m-3·s-1；Vi是对于植物i 的沉积速度，mm·s-1；Λ 是冲洗系数。沉积速度等参数的计算如下：

式中：Ra、Rc分别是大气和植物冠层阻尼系数；U 是风速，m·s-1；Zu、Z0分别是测风高度和地表粗糙度，m；Vplant，max是植物最大沉积速度（模型内部参数），m·s-1；LAIi、LAIi,max分别是计算时刻叶面指数和最大叶面指数。

对于土壤沉积的计算如下：

式中：As是土壤表面的总沉积量，Bq·m-2；Ads是土壤表面干沉积，Bq·m-2；Asoil，max是核素在农田中的最大干沉积速度，m·s-1。

1.2 农作物活度主要计算方法

核事故后农作物的活度来自叶吸收和根吸收。在叶吸收计算中，植物活度除了来自初始沉积，还需要考虑生长稀释、天气过程、核素衰减、核素易位等过程，这些过程都需通过参数化计算。根据过程不同叶吸收分为牧草、绿叶菜类、果实根茎类蔬菜和水果类、谷物类4 个基本计算模型。计算结果为各农产品的比活度（鲜质量），主要计算公式如下：

式中：Cil是来自叶吸收的活度，Bq·kg-1（FW）；Cir是来自根吸收的活度，Bq·kg-1（FW）；Yi是植物i 的单位产量，kg·m-2（FW）；λw、λr分别为天气衰减因子和核素衰减率；TFi是土壤到植物i 的转移因子；Cs是根部土壤的活度浓度，Bq·kg-1；As是土壤的总沉积，Bq·m-2；L是土壤深度，m；δ是土壤密度，g·m-3。

1.3 动物产品活度主要计算方法

农作物如牧草、玉米等作为饲料通过食入进入家禽、家畜体内，最终形成肉、蛋、奶中的活度。主要计算公式如下：

式中：Aanimal，i是动物i 食入的活度，Bq；Ck是饲料k 的活度，Bq·kg-1（FW）；Ik，i是饲料k 的食入量，kg·d-1（FW）；Ci是动物产品i的活度，Bq·kg-1（FW）；TFi是核素的转移因子；λb，ij是核素j在农产品i中的生物衰减率；ai，j和bi，j是核素在农产品中的生物半衰期对应参数。

1.4 建立模型内部参数

由上述沉积和迁移计算公式可见，计算涉及大量生态参数，而生态参数有明显的地域差异。AgrRAM模型建立后，对处于主要夏粮产地的田湾核电厂周边200 km 进行了生态参数调研和文献研究[13-15]，建立了一套基于该农业生态环境的生态参数体系。农业生态参数分为：农业区域特定参数和一般参数（与核素性质相关）。区域特定参数反映该地区农业生产的一般特征，具体如表1所示。

一般参数有：

（1）天气衰减因子（由于天气原因如风、雨、雾等导致的放射性活度的衰减），根据当地气候及降雨情况进行设置。

表1 区域特定参数Table 1 Specific ecological parameters in a region

（2）核素在土壤中的浸出率、固定率、吸收率、再悬浮率。

（3）核素从土壤到植物的转移因子（与土壤类型相关）。

（4）核素从饲料到肉蛋奶的转移因子。

（5）核素在各农作物中的易位因子。

（6）核素的生物半衰期。

除了天气因子，其他参数均选择了IAEA 技术报告[16]中的参数及JRODOS中的参数[17]。

2 田湾核电厂事故后果案例计算

2.1 核电厂址周边农业生态条件

田湾核电厂位于江苏省最北部，处于山东、江苏交界。该区域属于黄淮海地区，本区光热水分与地形条件耦合较佳，是我国传统农业区。以田湾核电厂为中心（34.68°N，119.46°E），用AgrRAM 模型进行200 km×200 km 范围假想核事故农业后果计算。该区域所属行政区有：山东省临沂市、日照市和江苏省连云港市，依据文献和统计年鉴收集特定生态参数。根据收集的环境气象数据及地理数据，对计算区域进行二级划分[18]，分为两个农业区如图2 所示，农业及气象特征如表2所示。

计算区域主要作物有：冬小麦、水稻、玉米、大豆、花生、甘薯、水果、蔬菜等，冬小麦和夏播玉米是该区的主要耕作方式。各农作物的播种收获时间如表3所示，典型农作物的叶面指数如表4 所示，单位产量来自文献和各地统计年鉴[19-21]。该区主要饲养家畜有牛、羊、猪和鸡，养殖业特点是精饲料充足，但是粗饲料特别是青饲料不足。奶牛的青饲占比为50%，产肉家畜的青饲占比为17%。建立的饲料配方见表5。

2.2 模型计算条件

不同事故发生时间对农业生产影响不同，因此选择2017 年4 月10 日、6 月1 日和10 月20 日分别进行事故后果计算。收集田湾核电厂气象铁塔及周边气象站点的观测数据，用扩散模型进行事故后的扩散计算，以扩散模型的空气积分浓度作为AgrRAM 模型的输入，进行农业事故后果计算。具体计算条件如表6所示。

表2 两个农业区的环境和农业特征Table 2 Environmental parameters and agriculture features in two sub-region

表3 两个放射生态区的农作物生长期Table 3 The growth and harvest period in two sub-region

表4 主要农作物不同生长期的叶面指数[月-日（叶面指数）]Table 4 The LAI of different plants in grown period

表5 两个区域的饲料清单（kg·d-1，FW）Table 5 Feed list of both region 001 and region 002（kg·d-1，FW）

2.3 结果讨论

2.3.1 农产品的活度分布

不同事故发生时间对于农业生态环境的影响不同，如果事故发生在作物生长期，活度来自叶面直接沉积的叶吸收和土壤的根吸收，如果事故发生在播种之前或收割之后，则没有直接沉积叶吸收，只有根吸收。为了展示活度随时间的变化，选择134Cs（ncs4）和半衰期较短的131I（nii1）进行不同事故时间活度对比（134Cs 半衰期为2 a，131I 半衰期为8 d）。图3～图5 是3个事故发生时间农产品在计算区域的比活度分布。事故发生在4 月10 日时的农产品最大比活度分布如图3 所示，青草和绿叶菜已经处于生长期，在青草和绿叶菜上发生了直接沉积。绿叶菜作为全部食用产品，叶吸收和根吸收导致绿叶菜较高的活度，而青草作为饲料通过迁移进入奶牛和牛奶形成了牛奶中的活度，由于迁移过程对活度的消减作用，最终牛奶的活度小于绿叶菜的活度。同时绿叶菜收割后很快被全部食用，核素衰减时间短，而牛奶需考虑从青草到奶牛及牛奶的迁移时间，核素衰减也导致牛奶活度比绿叶菜低。对于冬小麦，由于播种时间为9月底10月初，因此没有发生叶面直接沉积，活度来自土壤沉积后的根吸收，所以活度较低，131I对冬小麦未形成辐射危害（≤10-9Bq·kg-1FW）。在各个农产品中，134Cs活度均高于131I，这主要是由于核素半衰期的差异，从事故发生到收获农产品，131I衰减较快，导致其活度降低。

表6 AgrRAM田湾核电案例计算条件Table 6 AgrRAM simulation conditions of Tianwan NPP

图4 是事故发生在6月1日时各农产品的比活度分布。青草和绿叶菜在6 月进入成熟期，叶面指数达到最大值，当此时发生事故，放射性尘埃在作物上的直接沉积将会最大。如图所示，对于牛奶和绿叶菜事故后活度均处于较高分布，但是两种作物从4 月到6月叶面指数增长幅度较小，相同距离处活度略有增大。与图3 相比，差别较大的是131I 在冬小麦中的活度分布。这是由于事故时间6 月距离冬小麦播种时间比4 月更接近，所以对于短半衰期的131I，冬小麦活度虽比4 月增大，但也未形成辐射危害（≤10-8Bq·kg-1FW）。

事故发生在10月20日时各农产品最大比活度分布如图5 所示，此时变化最大的仍然是131I 在冬小麦中的比活度分布。由于10 月初冬小麦已经播种，发生了在幼苗上的直接沉积，最大比活度增至1×10-3Bq·kg-1FW，但仍未形成辐射危害。绿叶菜和青草在10 月仍然处于生长期，活度略有降低但变化不明显。对于其他作物例如水稻、玉米等已经收割的作物，事故发生在10 月与4、6 月相比变化较大，这些作物的活度来自下一个生长季的根吸收，对于短半衰期的131I，当下一季作物收割时，131I活度水平降至很低。

从上述计算结果可见，对于动态模型，由于农业生态环境和农产品自身变化，不同事故时间的事故后果有很大差异，同一事故不同农产品的污染水平也不相同。AgrRAM 模型展示了动态的计算过程，更真实地反映了实际事故后果。对于农作物影响最严重的事故是发生在植物生长期，放射性尘埃直接沉积在叶面进入作物，如果事故发生在临近收割之前，事故后果最严重。为了对农产品进行定量活度对比，采用日本福岛事故数据进行计算：事故为137Cs 的沉积，沉积量为100 kBq·m-2（日本福岛核事故30 km 范围内137Cs典型的沉积水平），计算3 个上述事故时间各农产品的最高活度，结果如图6所示。当事故发生在4月时，只有牛奶、羊肉、牛肉、鸡肉的比活度超出了一般食品的限值（1 000 Bq·kg-1FW）[22]，这些动物产品的活度来自饲料如青草等。当事故发生在6 月时，只有大葱、鸡肉、鸡蛋的比活度未超限值，其他产品都超过了限值。当事故发生在10 月时，青草、绿叶菜、牛奶、大葱、牛肉、羊肉的比活度超过限值。牛、羊肉中高的活度是由于饲料中的青草及储存干草的高活度。10 月一些作物的低活度是由于沉积发生在收割之后，活度来自下一个季节土壤向植物的转移。

2.3.2 典型农作物的活度时间序列

不同核素半衰期不同，图7 是事故发生在6 月1日时131I和134Cs在不同农产品中的比活度时间序列对比。比活度结果输出从各农产品的收割日开始，横坐标为时间步长（图中所示是从收割日到第60 个时间步的计算结果）。由于半衰期不同，131I 衰减速度很快，而134Cs 活度变化比较平缓。在青草中，由于生长期和收割期较长，所以核素直接沉积在青草上并形成较高的活度。但是随着时间的推移，131I 沉积发生后的175 d 活度已降至700 Bq·kg-1FW，而134Cs 活度仍然高达9.58×1010Bq·kg-1FW。对于该地区的中稻，由于沉积发生在播种之前，活度主要来自土壤向根的迁移，从沉积到水稻收割的190 d，131I 活度大幅衰减，收割时只有60 Bq·kg-1FW，而134Cs 活度达6×109Bq·kg-1FW，134Cs 在第40 a 活度才会降至几百Bq·kg-1FW。青草中高的活度浓度导致牛奶中的高活度（假设牛饲料中的青草均来自本地区），在沉积发生后第235 d，牛奶中131I 活度就降至100 Bq·kg-1FW，满足了GB/T 18871—2002 标准，但134Cs 活度一直在1×1012Bq·kg-1FW 数量级，同样大约40 a 后活度才会降至几百Bq·kg-1FW。绿叶菜与青草类似，有较长的生长期及收割期，虽然初始活度较高，但是在沉积发生后的210 d 时131I 的活度就降至几百Bq·kg-1FW，但是134Cs 活度一直较高，40～50 a 后活度将降至几百Bq·kg-1FW。

3 AgrRAM模型的验证分析

AgrRAM 模型计算流程简洁清晰，计算结果的精确性基于真实可靠的农业生态参数及输入的大气扩散模型的结果。AgrRAM 验证分析由两部分组成：第一部分是对于模型本身参数灵敏度的测试分析；第二部分是与平衡模型计算结果对比，验证模型的可靠性。

3.1 模型参数灵敏度测试

3.1.1 叶面指数

放射性尘埃沉积在植物或土壤上，沉积量的大小决定了农业环境的活度水平，而核素的沉积速度和气象条件影响了沉积量的大小。AgrRAM 中假设事故形成的是放射性气溶胶，农作物上的沉积速度计算见公式（4）。由公式可见沉积速度与作物的最大沉积速度和叶面指数相关，作物的最大沉积速度由模型内部参数给出，所以选择不同的叶面指数进行模型参数灵敏度测试。

3.1.2 转移因子

放射性核素在农业环境中的迁移包括土壤到作物的迁移（转移因子）和核素在作物内部的易位（易位因子），其中转移因子的作用更重要。IAEA EMRAS项目（2003—2011）生物工作组对11 种来源的转移因子进行了分析研究，并出版了IAEA TRS NO.472[23]等系列报告，研究指出不同来源的转移因子导致最终活度在1～2个数量级范围的变化。

选择代表性农作物和核素两组不同的转移因子进行参数灵敏度测试，改变叶面指数和转移因子，其他参数不变进行计算，计算条件如表6。图8 所示是事故发生在4月10日时的活度对比。

由图8a和图8b可见，当叶面指数增大10%时，绿叶菜的活度变化比较明显，其他作物差别不明显，这是因为农作物的叶面指数变化范围很小，进入成熟期各作物的叶面指数差别不大，所以叶面指数对活度变化的影响较小。131I和137Cs 在绿叶菜中的活度变化大于90Sr。由图8c和图8d可见，放射性核素从土壤向农作物的转移因子与作物最后的活度水平直接相关。由于计算的事故时间为大部分作物播种之前，此时核素由土壤到植物的转移因子与最终的活度水平是线性关系，当转移因子增大为模型设定数值的10 倍，计算结果也为原来的10 倍。绿叶菜的活度计算与其他作物的计算方式不同，而且活度除了根吸收还有叶面直接沉积。

3.2 验证分析

核事故后农作物的实际测量数据难以获得，而且事故发生地的生态参数、地理和气候数据也难以收集，因此无法用实测数据进行AgrRAM 模型的验证。为了检验AgrRAM 结果的合理性，采用与传统的平衡模型对比进行模型可靠性验证。平衡模型计算方法是由环境到作物直接迁移，不考虑生态环境的变化，长期被用来作为食物链模型的计算方法，为后果评估系统提供数据支持。为了评估AgrRAM 模型的可靠性，用两个模型进行平行计算，计算单位沉积量的最大比活度，结果见表7～表9。

由AgrRAM 与平衡模型的结果对比可见，两个模型对于农作物的比活度偏差在3 倍左右，对于动物产品的比活度偏差在5 倍左右，说明了AgrRAM 模型结果的可靠性。正如很多研究指出，在所有环境迁移模型中，生态环境的迁移模型不确定性最大，这是因为缺乏基础的生态参数和用于模型验证的监测数据。随着目前国际间的合作及大量基础研究的进行，模型相关的生态参数会更精确，模型的验证工作也会更有效。