基于欧拉方法的核素近海迁移计算

2020-06-21李子超司广成秦雪猛

核科学与工程 2020年2期

李子超，周涛，司广成，秦雪猛

(1.华北电力大学核科学与工程学院，北京 102206；2.东南大学能源与环境学院，江苏南京 211096；3.中国科学院海洋研究所，山东青岛 266071)

发展核电是我国能源发展的战略需求，也是解决环境污染、优化能源结构的有效途径。但是核电也有其发展的局限性，一旦发生核泄漏事故，将会造成极大的影响，日本福岛核电厂受地震海啸影响，发生严重事故，大量放射性核素释放，导致严重后果，造成不可估量的损失。我国核电厂分布在沿海地区，如果发生核泄漏事故，将严重威胁附近海域居民安全，污染海洋生态环境，有效应对海洋核污染已经成为一个重要议题。目前，全球第二台第三代核电机组AP1000[1]于2018年8月成功并网发电[2]，AP1000是先进第三代核电技术[3]，安全系数提高了一个量级。尽管如此，还不能完全排除发生核泄漏事故的可能性，尤其是该机组建于山东海阳，紧邻黄海，由于黄海是一个半封闭海域，一旦发生放射性泄漏事件，造成的后果将十分严重。福岛事故后，国外学者对核素在海洋中的输运进行了一系列研究。Inoue等人[4，5]分析了事故前后日本沿岸放射性物质的观测数据，结果显示事故发生后沿岸的核素浓度迅速升高，2011年6月因核素被带到远海而开始降低。Nakano等人[6]建立了全球放射性物质输运模式，分辨率为2.0°×2.0°，对137Cs的长期的输运路径进行了预测。Dietze等人[7]采用全球环流模式，对福岛近海进行了局部加密，模拟了137Cs的迁移扩散情况。国内学者也进行了不少数值研究。何晏春等人[8]利用海洋环流模式MICOM模拟了核素在海洋中的长距离输送，结果显示不同的排放情景和气象资料对核素在表层和次表层的迁移路径没有显著影响。赵昌等人[9]用普林斯顿海洋环流模式POM建立了全球放射性物质输运模型，对137Cs的输运进行了长期模拟，预测了三十年的迁移路径，表明四五年后核素扩散到美国西海岸，八九年后扩散到整个北太平洋。王辉等人[10]运用分辨率为1/8°×1/8°的北太平洋环流模式，对日本福岛核泄漏物质在海洋中的运输进行了十年的模拟和预测。乔清党等人[11]通过对国内外海流和放射性物质输运方式的分析，针对海洋环境问题，确立了核事故评价系统的总体方案。目前，对放射性核素在海洋中输运的研究多基于全球尺度，对近海域的核素扩散研究较少。因此，有必要基于AP1000核电机组，根据实时气象数据，建立高精度的近海核素扩散模型，分析核泄漏事故后核素的迁移路径，为核泄漏事故后的应急响应提供参考。

1 研究对象

海阳核电厂规划建设6台百万千瓦级压水堆机组，留有两台扩建余地，单机容量1 250 MW，热功率为3 415 MW。其中，一期工程建设2台美国西屋电气公司第三代核电技术AP1000百万千瓦级压水堆核电机组。核电厂位置及假想泄漏位置为东经121.3度，北纬36.7度。

2 计算模型

2.1 源项模型

为了实现快速响应，裂变产物的活度可以用简化的模型[12]进行估算，对于辐照时间远大于其半衰期的核素，产物可以达到平衡状态，堆芯积存量如公式(1)所示。

A=310YP

(1)

式中:A——活度，1012Bq；

Y——核素的裂变产额，%；

P——反应堆热功率，MW。

若核素半衰期明显长于照射时间，则活度将随时间线性增加，如公式(2)所示。

A=210YtP/T0.5

(2)

式中:T0.5——半衰期，s；

t——照射时间，s。

2.2 环流模型

ROMS是一个比较新的三维的、自由表面的随底坐标模式[13]，可以较好地表征近海地形变化，在温跃层和底层边界有更高的解析度，被广泛应用于海洋近岸和河口海洋环境预报，因此基于ROMS进行建模。在垂向静压近似和Boussinesq假定下，通过有限差分方法，近似求自由表面Reynolds平均的原始Navier—Stokes方程。模型在水平方向使用正交曲线网格，垂向采用地形拟合的可伸缩S坐标系统。

在水平笛卡尔坐标系与垂向σ坐标系下，运动方程如公式(3)至公式(5)所示，连续性方程如公式(6)所示，温度或盐度的对流扩散控制方程如公式(7)所示，海水状态方程如公式(8)所示。

(3)

(4)

(5)

(6)

(7)

ρ=f(C)

(8)

式中:u——x方向的流速分量，m/s；

v——y方向的流速分量，m/s；

Ω——垂向流速，m/s；

η——相对平均海平面的水位，m；

Hz——单元格高度，m；

f——科氏力参数；

P——压力，Pa；

ρ——海水密度kg/m3；

ρ0——参考密度，kg/m3；

g——重力加速度，m/s2；

υ——黏滞系数；

υθ——扩散系数；

C——温度、盐度或悬浮泥沙含量；

S——辐射应力项；

Csource——源汇项；

Dc——水平扩散项；

Du、Dv——水平黏性项。

2.3 衰变模型

放射性核素活度为单位时间内衰变的核子数，活度随时间的变化关系如公式(9)所示；衰变常数λ与半衰期之间的关系如公式(10)所示。

A=λN0e-λt=A0e-λt

(9)

λ=0.693/T0.5

(10)

式中：A0——原来活度，Bq；

A——经过t时间衰变后的活度，Bq；

T0.5——半衰期，s；

λ——衰变常数。

3.模型配置

3.1 核素释放量

根据2.1节中的源项模型，选择海阳核电站正在运行的一台AP1000机组，计算了137Cs和131I的堆芯积存量，如图1所示。

图1 核素堆芯积存量

从图1可以看出，由于137Cs半衰期很长，短时间内137Cs堆芯积存量基本随时间线性增加；131I半衰期较短，堆芯积存量很快能达到平衡。海阳核电厂AP1000机组热功率为3 415 MWt。设定运行时间为一年，则总的堆芯积存量131I为3.18×1018Bq，137Cs为1.48×1017Bq。根据美国反应堆安全研究报告[14]和福岛核事故放射性核素在海洋中的泄漏情况[15]，设定释放的131I为1×1015Bq，137Cs为1×1015Bq，与扩散时间相比，释放时间较短，设定为瞬时释放。

3.2 中国近海模型配置

中国近海水动力模型的地形数据采用GEBCO[16](General Bathymetric Chart of the Oceans)地形数据；初始场采用WOA13[17](World Ocean Atlas 2013)的气候态温盐场；海面边界采用COADS[18](Comprehensive Ocean-Atmosphere Data Set)风场和热通量；开边界采用SODA[19](Simple Ocean Data Assimilation)模型的水动力结果。

3.3 核电站近海水动力模型配置

核电站近海域水动力模型分辨率为1/60°×1/60°，垂向为地形跟随的坐标系统，分为12层。初始场采用中国近海水动力模型的温度和盐度结果；采用NCEP—DOE气象数据集中的风速、气压、气温，降雨等数据作为本模型的海面边界条件；中国近海水动力模型的温盐和海流结果作为开边界场，并在边界上引入潮汐。

3.4 海床基位置

收集了2014年6月海床基的实测水位数据，用来检验模型计算的结果，海床基布放位置为东经122.25度，北纬36.6度。

4 模型验证及结果分析

4.1 模型验证

根据2.2节中水动力模型及3.3中模型配置，基于实时气象数据[20]，计算了2014年6月8日至24日的半个月水动力结果，对应一个完整的潮周期，并与中国科学院海洋研究所的海床基实测水位进行了比较，海床基处模拟水位与实测水位如图2所示。

图2 模拟水位(蓝线)与实测水位(红线)结果

从图2可以看出，模拟结果与实测结果吻合较好。在6月14日至16日大潮期间，水位差基本控制在10 cm之内。

4.2 核电站近海域水动力场

黄海地区夏季多为西南风，典型的水文环境是黄海冷水团；冬季多为西北风，典型的水文环境是黄海暖水舌，因此黄海近海域水动力场会随着不同季节的风场和温盐场的变化而变化，首先分析夏季黄海近海域的水动力场。基于2.2节中水动力模型及3.2节中边界条件，模拟了黄海夏季流场和温度场。中国近海域水动力模拟结果较好地刻画了黄海冷水团等真实的海洋现象，从而也说明了模型的合理可靠性。基于2.2节中水动力模型及3.3节中边界条件，以气候态中国近海域水动力场结果作为边界，以六小时预报风场为海面边界条件，在山东半岛南部海区再建立核电站近海域水动力模型，计算了2014年6月8日至24日核电站近海域的水动力场。15日3时大潮期间涨急流场如图3(a)所示，15日5时大潮期间涨急流场如图3(b)所示，15日7时大潮期间落急流场如图3(c)所示，15日9时大潮期间落急流场如图3(d)所示。

从图3可以看出，涨潮期间，核电厂近海域水位上涨，海水向岸边流动；落潮期间，核电厂近海域水位下降，海水背离岸边流动。大潮期间涨急流和落急流较好地刻画了潮汐对海水运动和水位的影响，与真实海洋现象符合较好，核电厂近岸水动力模型的建立，将为下一步核物质扩散路径和分布的模拟奠定良好的基础。

4.3 不考虑衰变的核素扩散

核泄漏事故后放射性核素半衰期不同，137Cs半衰期为30年，而131I为8.3天。为了模拟短期核素扩散情况，实现快速应急响应，对于半衰期特别长的核素137Cs，衰变对其放射性的影响可以忽略，不考虑污染物的衰变，只考虑污染物的扩散过程。计算了2014年6月8日至24日两周内137Cs的迁移路径。137Cs在核电厂近海域扩散如图4所示。

从图4可以看出，受西南风影响，核素总体向东北方向迁移。从图4(a)可以看出，137Cs释放后第二天，核素基本在核电厂附近向远处扩散，随着与源项距离的增大，核素浓度逐渐减小。从图4(b)可以看出，污染物释放后第三天，核素向东方向迁移。从图4(c)可以看出，一周后，核素已经扩散到南黄岛村附近。从图4(d)可以看出，两周后核素已经跨过南黄岛村，向荣成海岸迁移。可以看出，海阳核电厂位置特殊，六月期间，一旦发生核泄漏事故，核素很容易向东北方向沿着海岸线迁移扩散，而这一地区多为旅游区，因此，其安全性至关重要。

(a)大潮期间涨急流场(15日3时)

(b)大潮期间涨急流场(15日5时)

(c)大潮期间落急流场(15日7时)

(d)大潮期间落急流场(15日9时)

(a)137Cs释放第二天

(b)137Cs释放第三天

(c)污染物释放一周

(d)污染物释放两周

为了分析不同位置核素浓度随时间的变化关系，设定了三个不同站位，三个站位的核素浓度变化如图5所示。

图5 不同站位核素浓度变化

从图5可以看出，核素在海流和扩散的作用下，先后迁移到站位01、站位02和站位03。开始时每个站位核素浓度基本为零，随着时间的增加，当核素扩散到相应站位时，核素浓度迅速升高，到达峰值后，缓慢降低。

4.4 考虑衰变核素扩散结果

对于半衰期特别短的核素131I，衰变对其放射性的影响不可忽略，考虑污染物的衰变和扩散过程，计算了2014年6月8日至24日两周内131I的扩散规律。131I在核电厂近海域扩散规律如图6所示。

从图6可以看出，131I在核电厂近海域的迁移方向与137Cs基本一致，都是沿着海岸线向东迁移，同时向南扩散。但是131I半衰期为8.3天，两周后放射性活度变为原来的30%左右，因此，随着扩散时间的增加，浓度不但因扩散而降低，还会因衰变而减小。对比图4可以看出，核素释放前两天，131I和137Cs的放射性活度分布没有明显区别；释放一周和两周后，131I的放射性活度明显低于137Cs的放射性活度。

4.5 应对措施

福岛核事故后，日本政府没有及时地科学预测核素扩散范围，应急区划分不明确，暴露出响应迟缓、部署混乱等诸多问题。基于福岛核事故的教训，沿海核电厂一旦发生核泄漏事故，应将高放废液进行吸附处理后再排放。如果进入海洋的放射性核素数量较大，应基于八小时预报风场，快速滚动计算并预测核素迁移路径和扩散范围，做好人员疏散等应急预案。此外，平时工作中应进行政府、厂区和公众的联合应急演练，以便发生事故时可以高效有序部署。