-
In order to predict the release behavior of fission gas in large grain UO2 fuel and provide support for the development of accident tolerant fuel, a phase-field model is used to simulate the release behavior of fission gas in the microstructure of UO2 polycrystalline in this work. This model adopts a set of coupled Cahn-Hilliard equations and Allen-Cahn equations, using conserved field variables to represent the distribution of fission gas and vacancies, and distinguishing bubble phase from matrix phase by using order parameters. This model focuses on investigating the effects of different grain sizes, temperature conditions, and diffusion coefficients on the release behavior of fission gas, demonstrating the nucleation, growth, and fusion behavior of bubbles. Simulation results are obtained for fuel porosity, bubble coverage on grain boundaries, and average bubble radius at a certain degree of burnup. The results show that temperature and diffusion coefficient have a significant influence on porosity and bubble coverage on grain boundaries. When the diffusion coefficient is high, grain size also has a significant influence on fission gas release behavior. And when the diffusion coefficient is low, the influence of grain size is not significant. In addition, the distribution of fission gas bubbles under high burnup obtained through this model is also in good agreement with experimental result. The model can predict the behavior of fission gas release in large grain UO2 fuel.
-
Keywords:
- phase field model /
- large grain UO2 fuel /
- release of fission gas
1. 引 言
在反应堆的运行过程中, 随着裂变的进行会持续产生多种裂变产物, 例如惰性气体氙(Xe)和氪(Kr), 这两种气体每次裂变的产额大约0.25—0.3个原子[1], Xe的产生速率是Kr的将近10倍[2]. 由于这两种气体在UO2中的溶解度很低, 因此大多数的Xe和Kr通过扩散在晶内或晶间形成气泡[3,4]. 随着气泡的增加和生长, 晶内气泡会导致燃料芯块的肿胀[3], 晶间气泡连接形成通道, 裂变气体沿通道扩散到自由空间, 造成裂变气体的释放[4].
裂变气体的产生与释放是提升反应堆燃耗水平的关键制约因素之一. 裂变气体在反应堆燃料棒内的积累, 特别是反应堆燃耗达到一定水平时, 会导致反应堆燃料棒内压升高, 应力增大, 从而导致燃料棒包壳形变, 并最终导致包壳破裂, 造成放射性物质的泄漏[5]. 因此, 研究裂变气体的释放规律能够辅助改进燃料设计, 改善反应堆安全性并提升反应堆整体经济性[6].
根据研究[7,8], 提升UO2燃料的初始晶粒尺寸可以有效减少裂变气体的释放. 随着晶粒尺寸的增大, 裂变气体扩散到晶界的平均距离增大, 从而更易被晶粒滞留[5,9]. 并且由于大晶粒相对较小的塑性蠕变与热蠕变能减小燃料的肿胀, 减缓芯块-包壳相互作用[5,10], 提升反应堆的安全性. 因此, 大晶粒UO2燃料作为事故容错燃料的候选之一, 近年来受到了持续关注.
大晶粒UO2燃料的制造工艺可通过向UO2中掺杂某些氧化物实现, 如Cr2O3, TiO2, Al2O3, MgO等[11]. 以Cr2O3为例, 当UO2燃料中Cr的质量百分比达到0.16%时(接近UO2中Cr的溶解度极限), 可使晶粒平均尺寸增长到50—60 μm, 同时不过多改变燃料的材料特性[12]. 针对掺杂Cr2O3的大晶粒UO2燃料的相关研究已经展开, Killeen[7]通过实验对比了Cr2O3掺杂质量分数为0.5%的UO2燃料与未掺杂的燃料在1500 ℃、不同燃耗深度下, 燃料的肿胀与裂变气体的释放; Kashibe和Une[13]测试了不同掺杂组分, 如掺杂质量分数为0.065% Cr2O3, 掺杂质量分数为0.076% Al2O3等, 经由1100—1600 ℃的的辐照后退火实验, 对裂变气体Xe扩散系数的影响; Che等[14]使用BISON燃料性能分析代码对Cr2O3掺杂的大晶粒UO2燃料棒的性能进行了计算. Cooper等[8]在Che等[14]的基础上, 进一步研究了掺杂对UO2燃料扩散系数根源性的影响. 目前围绕大晶粒UO2开展的相关研究多是基于实验或宏观系统模型, 缺少基于介观尺度对大晶粒UO2燃料中裂变气体释放行为的研究, 通过介观尺度研究裂变气体释放时燃料微观结构的演化, 能更好理解大晶粒UO2对裂变气体释放的迟滞作用.
相场法是实现材料学和物理学中不同驱动力(如成分梯度、温度、应力应变、电场、磁场)下微观结构演化模拟的一种重要的介观尺度研究工具. 其以Ginzburg-Landau理论为基础, 通过偏微分方程建立起一种能够准确描述系统随时间演化的数学模型[15], 用以获取体系在时间和空间上的瞬时状态. 相较于蒙特卡罗法、元胞自动机法, 相场法通过将不同的界面描述纳入其公式, 来避免复杂界面的追踪问题, 且易与其他物理场(如噪声场、速度场、成分场等)耦合[16]. 使用相场法模拟常规UO2燃料裂变气体释放行为已有不少相关研究. Hu和Henager [17]采用Cahn-Hilliard方程构建了一个包含空位浓度场与间隙原子浓度场的相场模型, 该模型模拟了辐照条件下, 空洞的演化规律, 得出空洞临界半径、数量密度等一系列参数. Millett等[18]构建了一个包含Allen-Cahn方程的相场模型, 该方程通过控制序参量来区分不同相界, 同时还引入裂变气体浓度场, 研究了辐照条件下, 燃料基质中裂变气体气泡的演化规律. Aagesen等[9]采用抛物线近似的自由能来构建相场模型, 研究裂变气体气泡沿三岔晶界的生长规律, 该模型对相关参数的简化可提高计算效率, 允许模拟更广泛的区域.
采用相场法模拟裂变气体释放的相关研究已比较详尽, 但通过相场法研究晶粒尺寸、温度与扩散系数等因素对裂变气体释放行为的影响则不够充分. 本研究使用相场模型, 参考上述研究, 采用抛物线近似的自由能泛函、真实的扩散系数值与微米级的晶粒尺寸, 对辐照下不同晶粒尺寸、不同温度与扩散系数条件下UO2燃料中的裂变气体释放行为, 包括裂变气体向晶界扩散, 气泡的形核与生长等, 进行了模拟. 第2节将对本模型进行阐述, 包括本模型采用的驱动方程、自由能函数的成分、变量的定义、参数的选取等相关内容. 第3节为模拟结果, 包括常规晶粒尺寸UO2燃料的模拟结果验证、温度对裂变气体释放的影响、晶粒尺寸对裂变气体释放的影响、以及同温度下采用不同扩散系数(裂变气体在未掺杂的UO2燃料中的扩散系数与在掺杂Cr2O3的燃料中的扩散系数)对裂变气体释放的影响4个部分. 第4节对本文进行概括总结.
2. 相场模型
2.1 自变量定义
本模型考察辐照效应产生的点缺陷空位与裂变气体在多晶材料中的浓度分布, 引入两个守恒场变量: cv(r, t), cg(r, t), 分别表示空位摩尔分数和裂变气体摩尔分数, 该摩尔分数c可通过c = ρVa计算, 其中ρ为空位或裂变气体的数密度, Va为UO2晶胞中U原子体积, 其值为0.0409 nm³[9].
如引言所述, 燃料基质中的裂变气体通过扩散形成气泡, 存在两个稳定相: 燃料固体基质相和裂变气体气泡相, 使用序参量η区分该两相: 当η = 0时为基质相, η = 1时为气泡相, 如图1(a)所示. 基质相与气泡相中存在空位与气体的平衡浓度, 对于基质相, 空位与气体的平衡浓度用下式表示:
$$ c_{\text{v}}^{{\text{m,eq}}} = \exp \left( { - \frac{{E_{\text{v}}^{\text{f}}}}{{{k_{\text{B}}}T}}} \right){,} $$ (1) $$ c_{\text{g}}^{{\text{m,eq}}} = \exp \left( { - \frac{{E_{\text{g}}^{\text{f}}}}{{{k_{\text{B}}}T}}} \right){,} $$ (2) 其中, $ E_{\text{v}}^{\text{f}} $, $ E_{\text{g}}^{\text{f}} $分别为空位与裂变气体原子缺陷的形成能, 取$ E_{\text{v}}^{\text{f}} $= 5.1 eV[19], Xe原子在空位中的溶解能$ E_{\text{g}}^{\text{s}} $= 5.21 eV[20], 取裂变气体原子的形成能为空位形成能与原子在空位中的溶解能之和, $ E_{\text{g}}^{\text{f}} $=$ E_{\text{v}}^{\text{f}} + E_{\text{g}}^{\text{s}} $= 10.31 eV; kB为玻尔兹曼常量, T为温度. 对于气泡相, 取气体平衡浓度$ c_{\text{g}}^{{\text{b, eq}}} $= 0.454[9], 当cg = 0.454时, 气泡相的范德瓦耳斯自由能最小. 由于模拟过程中气体浓度cg存在涨落, 为保持气泡相中cv + cg ≤ 1, 取空位平衡浓度$ c_{\text{v}}^{{\text{b, eq}}} $= 0.9 – cg, 随着cg的增大, $ c_{\text{v}}^{{\text{b, eq}}} $将逐渐减小.
本相场模型关注裂变气体在多晶燃料中的释放行为, 使用序参量ϕi (i = 1, 2, ···, P)区分P种不同取向的晶粒, 在第i个取向的晶粒内部ϕi = 1, 如图1(b)所示, 在其他取向的晶粒相中, 该ϕi = 0. 在不同取向晶粒的接触面附近, ϕi近似等于0.5, 通过考察Φ =$ \displaystyle\sum\nolimits_{i = 1}^P {\phi _i^2} $可区分晶体内部与晶界, 当Φ = 1为晶体内, Φ ≈ 0.5为晶界, Φ = 0的区域与η = 1的区域重合, 为气泡相. 在本模型中, 不同取向的晶粒不对裂变气体与空位扩散产生相关的应力作用, 仅用于划分不同晶粒区域, 产生晶界.
2.2 自由能方程
相场法采用各相和界面自由能所组成的总自由能泛函F, 与Cahn-Hilliard定义的非均匀系统自由能一致[21]. 本研究采用的自由能泛函F形式如下:
$$ \begin{split} \;&F = \int_V {h(\eta ){f^{{\text{solid}}}}({c_{\text{v}}},{c_{\text{g}}}) + j(\eta ){f^{{\text{bubble}}}}({c_{\text{v}}},{c_{\text{g}}})} \\ & ~~+ {f^{{\text{poly}}}}(\eta ,{\phi _{1 \to P}}) + {f^{{\text{grad}}}}({c_{\text{v}}},{c_{\text{g}}},\eta ,{\phi _{1 \to P}}){\text{d}}V, \end{split} $$ (3) 式中, f solid为基质相自由能密度函数, f bubble为气泡相自由能密度函数, f poly为多晶体自由能密度函数, f grad为梯度贡献项. h(η) = (1 – η)2, j(η) = η2为权重函数, 当η = 0时, f solid对总自由能的贡献达到最大, f bubble的贡献为0, η = 1时则相反.
假设基质相的化学自由能近似为理想解, 根据文献[22], 基质相自由能密度函数写为
$$ \begin{split} & {f^{{\text{solid,ideal}}}} =\\ \;& \frac{1}{{{V_{\text{m}}}}} \big\{ RT[ {{c_{\text{v}}}\ln {c_{\text{v}}} + ( {1 - {c_{\text{v}}}})\ln ({1 - {c_{\text{v}}}})} ] + {N_{\text{A}}}E_{\text{v}}^{\text{f}}{c_{\text{v}}}\\ &+ RT [ {{c_{\text{g}}}\ln {c_{\text{g}}} + ( {1 - {c_{\text{g}}}} )\ln ( {1 - {c_{\text{g}}}} )} ] + {N_{\text{A}}}E_{\text{g}}^{\text{f}}{c_{\text{g}}} \big\}, \end{split} $$ (4) 其中Vm为单位晶格内U原子的摩尔体积, Vm = VaNA, R为理想气体常数, NA为阿伏伽德罗常数.
气泡相采用范德瓦耳斯气体处理, 其自由能密度函数由范德瓦耳斯气体的亥姆霍兹自由能给出[9,23]:
$$ {f^{{\text{bubble,vdw}}}} = {\rho _{\text{g}}}k{}_{\text{B}}T\left[ {\ln \left( {\frac{1}{{{n_{\text{Q}}}\left( {{1}/{{{\rho _{\text{g}}}}} - b} \right)}}} \right) - 1} \right] + {f_0}, $$ (5) 其中ρg为气体原子数密度; nQ = $ {\left( {\dfrac{{m{k_{\text{B}}}T}}{{2{\text{π }}{\hbar ^2}}}} \right)^{3/2}} $是Xe的量子浓度, m为Xe原子的质量, $ \hbar $为约化普朗克常数; 参数b与Xe原子体积相关, 其值为0.082 nm3/原子[9]; f0表示基质相与气泡相自由能之间的偏移, 通过f0的修正使得基质相与气泡相自由能的最小值相等, 其值为2.22 × 109 J/m3 [9].
为简化控制方程的数值计算, 对基质相自由能密度函数和气泡相自由能密度函数进行抛物线近似, 与Li等[22]采用的自由能密度函数相似, 但进行简单修正以使基质相晶界处的自由能小于晶粒内的自由能, 近似后的基质相自由能密度函数f solid的表达式如下:
$$ \begin{split} &{f^{{\text{solid}}}}({c_{\text{v}}},{c_{\text{g}}}) = {f^{{\text{solid,v}}}}({c_{\text{v}}}) + {f^{{\text{solid,g}}}}({c_{\text{g}}}) = \\ & \frac{1}{2}{f^{{\text{m,fix}}}}k_{\text{v}}^{\text{m}}{({c_{\text{v}}} - c_{\text{v}}^{{\text{m,eq}}})^2} + \frac{1}{2}{f^{{\text{m,fix}}}}k_{\text{g}}^{\text{m}}{({c_{\text{g}}} - c_{\text{g}}^{{\text{m,eq}}})^2}, \end{split} $$ (6) 其中修正项为$ {f^{{\text{m,fix}}}} = {4}/{5} + {1}/{5}\varPhi $, 假设基质相晶界处空位与气体的自由能密度是晶粒内的约十分之九, 这么做使得空位与裂变气体产生向晶界扩散的趋势.
自由能密度函数表达式中的$ k_{\text{v}}^{\text{m}} $, $ k_{\text{g}}^{\text{m}} $控制抛物线的曲率, 形式如下:
$$ k_{\text{v}}^{\text{m}} = \frac{1}{{(c_{\text{v}}^{0} - c_{\text{v}}^{{\text{m,eq}}})}}\left[ {\frac{{RT}}{{{V_{\text{m}}}}}\left[ {\ln c_{\text{v}}^{0} - \ln (1 - c_{\text{v}}^{0})} \right] + \frac{{{N_{\text{A}}}E_{\text{v}}^{\text{f}}}}{{{V_{\text{m}}}}}} \right], $$ (7) $$ k_{\text{g}}^{\text{m}} = \frac{1}{{(c_{\text{g}}^{0} - c_{\text{g}}^{{\text{m,eq}}})}}\left[ {\frac{{RT}}{{{V_{\text{m}}}}}\left[ {\ln c_{\text{g}}^{0} - \ln (1 - c_{\text{g}}^{0})} \right] + \frac{{{N_{\text{A}}}E_{\text{g}}^{\text{f}}}}{{{V_{\text{m}}}}}} \right], $$ (8) 其中假设$ c_{\text{v}}^{0} = c_{\text{g}}^{0} = 0.01 $[22]. T = 1276 K条件下, 计算得$ k_{\text{v}}^{\text{b}}= 1.7999 \times 10^{12}$ J/m³, $ k_{\text{g}}^{\text{b}}= 3.8408 \times 10^{12}$ J/m³.
类似地, 近似后的气泡相自由能密度函数f bubble的表达式如下:
$$ \begin{split} &{f^{{\text{bubble}}}}({c_{\text{v}}},{c_{\text{g}}}) = {f^{{\text{bubble,v}}}}({c_{\text{v}}}) + {f^{{\text{bubble,g}}}}({c_{\text{g}}}) \\ =\;& \frac{1}{2}f_{\text{v}}^{{\text{b,fix}}}k_{\text{v}}^{\text{b}}{({c_{\text{v}}} - c_{\text{v}}^{{\text{b,eq}}})^2} + \frac{1}{2}f_{\text{g}}^{{\text{b,fix}}}k_{\text{g}}^{\text{b}}{({c_{\text{g}}} - c_{\text{g}}^{{\text{b,eq}}})^2}. \end{split} $$ (9) 修正项:
$$ f_{\text{v}}^{{\text{b,fix}}} = \left\{ \begin{aligned} &1,&{c_{\text{v}}} \leqslant c_{\text{v}}^{{\text{b,eq}}}, \\ &{\left( {1 + {c_{\text{v}}} - c_{\text{v}}^{{\text{b,eq}}}} \right)^{10}}{,}&{c_{\text{v}}} > c_{\text{v}}^{{\text{b,eq}}}, \end{aligned} \right. $$ (10) $$ f_{\text{g}}^{{\text{b,fix}}} = \left\{ \begin{aligned} &1,&{c_{\text{g}}} \leqslant c_{\text{g}}^{{\text{b,eq}}}, \\ &{\left( {1 + {c_{\text{g}}} - c_{\text{g}}^{{\text{b,eq}}}} \right)^{10}}{,}&{c_{\text{g}}} > c_{\text{g}}^{{\text{b,eq}}}, \end{aligned} \right. $$ (11) 该修正项是根据气泡相亥姆霍兹自由能随气体浓度的增大所表现的特征而添加的, 当气体浓度cg趋近平衡浓度0.454时, 亥姆霍兹自由能将急剧增长[9]. 同时假设空位浓度cv趋近1时自由能也有类似的特征. 取f solid的曲率$ k_{\text{g}}^{\text{b}} $= 2.25 × 109 J/m3, 通过该曲率控制的近似自由能随浓度的变化与亥姆霍兹自由能在0 < cg < 0.454范围上较匹配. 另一曲率$ k_{\text{v}}^{\text{b}} $= 2 × 1010 J/m³, 这是假设气泡形核时, 基质相中空位浓度cv不超过0.02所设的值. 基质相自由能密度与气泡相自由能密度随浓度的变化如图2所示.
多晶体自由能密度函数f poly的表达式如下:
$$\begin{split} & {f^{{\text{poly}}}}(\eta ,{\phi _{1 \to P}}) \\ =\;& m\Bigg[ \sum\limits_{i = 1}^P {\left( {\frac{{\phi _i^4}}{4} - \frac{{\phi _i^2}}{2}} \right)} + \left( {\frac{{{\eta ^4}}}{4} - \frac{{{\eta ^2}}}{2}} \right)\\ &+ {a_{{\text{GB}}}}\sum\limits_{i = 1}^{P - 1} {\sum\limits_{j > i}^P {\phi _i^2\phi _j^2} } + {a_{\text{s}}}\sum\limits_{i = 1}^P {\phi _i^2{\eta ^2}} \Bigg],\\[-15pt]\end{split} $$ (12) 式中, m是自由能势垒系数, aGB和as为扩散界面系数, aGB取1.2, as取0.8.
梯度贡献项f grad的表达式如下:
$$\begin{split} {f^{{\text{grad}}}}({c_{\text{v}}},{c_{\text{g}}},\eta ,{\phi _i}) =\;&\frac{{{\kappa _{\text{v}}}}}{2}|\nabla {c_{\text{v}}}{|^2} + \frac{{{\kappa _{\text{g}}}}}{2}|\nabla {c_{\text{g}}}{|^2} + \frac{{{\kappa _\eta }}}{2}|\nabla \eta {|^2} \\ & + \frac{{{\kappa _\phi }}}{2}\sum\limits_{i = 1}^P {|\nabla {\phi _i}{|^2}} ,\\[-1pt]\end{split} $$ (13) 式中κv, κg, κη, κϕ分别是空位、气体、两种序参量的梯度项系数.
m和κϕ与界面能γint和界面宽度lint相关[24,25]:
$$ m = \frac{{6{\gamma _{{\text{int}}}}}}{{{l_{{\text{int}}}}}}, ~~$$ (14) $$ {\kappa _\phi } = \frac{3}{4}{\gamma _{{\text{int}}}}{l_{{\text{int}}}}, $$ (15) γint = 1.5 J/m2 [9], lint取0.3 μm, 该值大于参考文献[9]给出的值, 计算得m = 3.0 × 107 J/m3, κϕ = 3.38 × 10–7 J/m3. κv, κg, κη取1.69 × 10–6 J/m3.
2.3 控制方程
本相场模型的控制方程采用Cahn-Hilliard方程[21,26]控制守恒场变量cv, cg, 采用Allen-Cahn方程[27,28]控制序参量η, ϕ. cv, cg的控制方程如下:
$$ \frac{{\partial {c_{\text{v}}}}}{{\partial t}} = \nabla \cdot \left( {{M_{\text{v}}}\nabla \frac{{{\text{δ}} F}}{{{\text{δ}} {c_{\text{v}}}}}} \right) + {P_{\text{v}}}({\boldsymbol{r}},t) + \xi {}_{\text{v}}({\boldsymbol{r}},t), $$ (16) $$ \frac{{\partial {c_{\text{g}}}}}{{\partial t}} = \nabla \cdot \left( {{M_{\text{g}}}\nabla \frac{{{\text{δ}} F}}{{{\text{δ}} {c_{\text{g}}}}}} \right) + {P_{\text{g}}}({\boldsymbol{r}},t) + \xi {}_{\text{g}}({\boldsymbol{r}},t). $$ (17) 这里, Mv, Mg分别为空位与气体的迁移率,
$$ {M_{\text{v}}} = \frac{{{V_{\text{a}}}{D_{\text{v}}}{c_{\text{v}}}}}{{{k_{\text{B}}}T}}, $$ (18) $$ {M_{\text{g}}} = \frac{{{V_{\text{a}}}{D_{\text{g}}}{c_{\text{g}}}}}{{{k_{\text{B}}}T}}, $$ (19) 其中, Dv, Dg分别为空位与气体的扩散系数, 假设Dv = Dg. 未掺杂Cr2O3的UO2燃料中的裂变气体扩散系数由Turnbull的模型给出[29,30]:
$$ D_{1}^{{\text{undoped}}} = 7.6 \times {10^{ - 10}} \exp \left( \dfrac{ - 4.68 \times {{10}^{ - 19}}({\text{J}})}{{{k_{\text{B}}}T}} \right), $$ (20) $$ \begin{split} D_{2}^{{\text{undoped}}} = \;&5.64 \times {10^{ - 25}}\sqrt {\dot F} \exp \left( \dfrac{-1.91 \times {{10}^{ - 19}}({\text{J}})}{ {{k_{\text{B}}}T} } \right),\end{split} $$ (21) $$ D_3^{{\text{undoped}}} = 8 \times {10^{ - 40}}\dot F. $$ (22) 对应3种温度条件下裂变气体的扩散[29], 其中, D1为高温条件(T大于1650 K)下的本征扩散, D2为中间温度条件(T在1650—1350 K之间)下的辐射增强扩散, D3为低温条件(T小于1350 K)下无热辐射驱动扩散.
掺杂Cr2O3的UO2燃料的裂变气体扩散系数引自文献[8]:
$$ D_1^{{\text{doped}}} = \exp \left[ { - \frac{{\Delta {H_1}}}{{{k_{\text{B}}}}}\left( {\frac{1}{T} - \frac{1}{{{T_1}}}} \right)} \right]D_{1}^{{\text{undoped}}}, $$ (23) $$ D_{2}^{{\text{doped}}} = \exp \left[ { - \frac{{\Delta {H_2}}}{{{k_{\text{B}}}}}\left( {\frac{1}{T} - \frac{1}{{{T_2}}}} \right)} \right]D_2^{{\text{undoped}}}, $$ (24) $$ D_3^{{\text{doped}}} = D_3^{{\text{undoped}}}, $$ (25) 其中, ΔH1 = 0.3198 eV, ΔH2 = –0.3345 eV, T1 = T2 = 1773 K. 在模拟中, 根据不同温度、晶粒尺寸采用不同的扩散系数.
该控制方程包含两个源项, Pv, Pg分别表示空位与气体的产生速率, 裂变气体Pg产生速率由Pg =$ {V_{\text{a}}}\dot FY $计算, $ \dot F $为裂变率密度, 根据典型压水堆的运行数据, $ \dot F $取1.09 × 1013 次裂变/(cm³·s)[9], Y为U-235每次裂变裂变气体的产额, 取Y = 0.27. 计算得Pg = 1.2037 × 1010 s–1, 并假设Pv = 20Pg = 2.4074 × 1011 s–1. ξv, ξg是用简单随机函数分别表示的空位与气体的热涨落, 以使浓度非均匀分布, ξv的具体形式为$ {\xi _{\text{v}}} = {c_{\text{v}}}{R_{\text{v}}}(r, t) $, Rv为随机函数, 范围从–0.01—0.01, ξg与之类似.
序参量η, ϕ的控制方程如下:
$$ \frac{{\partial \eta }}{{\partial t}} = - {L_\eta }\frac{{{\text{δ}} F}}{{{\text{δ}} \eta }} , $$ (26) $$ \frac{{\partial {\phi _i}}}{{\partial t}} = - {L_\phi }\frac{{{\text{δ}} F}}{{{\text{δ}} {\phi _i}}}{\text{ }}(i = 1 \to P) , $$ (27) 其中, Lη, Lϕ为迁移率, Lη = Lϕ = 1.56 × 1011 m3/(J·s) [9].
模拟中使用的部分参数如表1所列. 本研究使用COMSOL Multiphysics® [31]进行建模, 该软件采用有限元法求解偏微分方程, 时间步进为默认配置, 网格尺寸0.1 μm. 对模型进行尺度变换但不做无量纲处理, 时间变换t* = 1 × 10–5 , τ 为模拟时间, 对应实际时间t = τ/t*, 自由能密度变换f * = 0.1. 模拟中采用的扩散系数$ \bar D $= D/t*, 迁移率$ \bar L $= L/t*, 自由能密度$ \bar f $= f × f *.
参数 符号 值 玻尔兹曼常量/(J·K–1) kB 1.3806 × 10–23 理想气体常数/(J·mol–1·K–1) R 8.3145 UO2晶胞U原子体积/ nm–³ Va 0.0409 空位形成能/eV $ E_{\text{v}}^{\text{f}} $ 5.1 气体原子缺陷形成能/eV $ E_{\text{g}}^{\text{f}} $ 10.31 扩散界面系数 aGB 1.2 as 0.8 自由能势垒系数/(J·m–3) m 3.0 × 107 梯度项系数/(J·m–3) κϕ 3.38 × 10–7 κv, κg, κη 1.69 × 10–6 迁移率/(m3·J–1·s–1) Lη, Lϕ 1.56 × 10–11 裂变率密度/(次裂变·m3·s–1) $ \dot F $ 1.09 × 1019 Xe产额 Y 0.27 3. 结果与分析
3.1 常规晶粒模拟结果与验证
本节讨论在5 μm的固定晶粒尺寸, T = 1276 K, D = 8.72 × 10–21 m2/s条件下, 裂变气体的释放行为, 并与相关实验与模拟结果进行对比和分析. 模型根据不同晶粒尺寸选取不同模拟区域大小, 针对直径5 μm的晶粒选取的模拟区域为10 μm × 10 μm, 周期性边界条件. 在模拟区域中构建4个正六边形晶粒, 空位与裂变气体Xe在整个区域内的初始浓度为基质相平衡浓度$ c_{\text{v}}^{{\text{m, eq}}} $, $ c_{\text{g}}^{{\text{m, eq}}} $.
T = 1276 K, 0 < τ < 550 s的晶粒与气泡分布随时间的演化如图3所示, 该分布采用Φ在模拟区域的分布, 其中红色区域为燃料的基质相, 介于红色区域间的黄色线条为晶界, 蓝色区域为气泡相. 图3(a)为初始状态, 由4个完整的六边形燃料晶粒构成. 随着裂变的进行, 燃料中不断产生裂变气体与空位, 两者浓度不断升高, 同时空位与裂变气体Xe原子向晶界扩散, 晶界处的缺陷浓度高于晶粒内, 如图4(a), (b)所示. 在模拟时间τ < 140 s内, 没有观察到气泡形核, 分析认为这一时间段为气泡演化的孕育阶段, 同Millett在文献[32]中的模拟一样, 在这一阶段, 缺陷浓度上升但不形成气泡. 随着缺陷浓度的上升, 系统总自由能逐渐增加, 在τ =160 s时, 观察到气泡在晶界处形核. 在 τ = 200 s时, 气泡在晶粒内大量形核, 当气泡数密度达到一定水平时, 形核停止, 在之后的时间里没有新的气泡产生, 现有气泡不断长大, 同时距离较近的气泡会互相融合形成更大的气泡, 气泡数密度降低, 平均直径增加, 这一过程被称为奥斯特瓦尔德熟化[33]. τ = 550 s时的空位与裂变气体的浓度分布如图4(c), (d)所示, 气泡中空位的平均浓度$ {\bar c_{\text{v}}} $≈ 0.95, 裂变气体的平均浓度$ {\bar c_{\text{g}}} $≈ 0.05, 晶界上气泡平均直径0.37 μm, 稍大于平均直径0.30 μm的晶粒内部气泡. 在晶界附近观察到明显的“无气泡区”, 与Bullough和Nelson[34]在实验中观察到的和Millett等[18]模拟得到的“空洞剥蚀区”类似, 表明晶界对缺陷和气体原子有吸收作用. 在本模型中, 由于空位与裂变气体的源项为净生成项, 且为固定正值, 无湮灭项、复合项、重融项等会减少点缺陷浓度的项, 气泡的尺寸会随着时间无限生长下去, 故选取600 s为模拟的结束时间.
经过换算, 该条件下模拟中τ = 550 s的累计裂变时间对应实际燃耗25.03 GWd/tU (UO2密度取10 t/m3), 在这一时刻, 燃料中裂变气体气泡分布(图3(f))与Zacharie等[35]实验中得到的在堆芯温度不超过1100 ℃的反应堆中燃烧至25 GWd/tU并退火的UO2燃料样品断裂表面的扫描电子显微镜照片较为相似.
整个模拟区域内平均自由能密度随时间的变化如图5(a)所示. 平均自由能密度在τ = 180 s左右存在峰值, 到达峰值后急剧下降, 随后缓慢攀升. 自由能抵达峰值到急剧下降这一时间段, 与气泡的形核期时间重合(对比图9(a)红色曲线), 表明促使裂变气体扩散形成气泡的驱动因素为整体自由能趋于最低. 从图2中自由能随浓度的变化曲线也可以观察到, 当空位与裂变气体的浓度较低时, 基质相自由能密度小于气泡相, 空位与裂变气体弥散在基质中不会形成气泡, 当浓度达到一定水平, 基质相自由能密度开始大于气泡相时, 空位与裂变气体开始形成气泡, 以满足整体自由能最小状态.
τ = 550 s时整个模拟区域总自由能密度分布如图5(b)所示, 忽略f grad的贡献. 可以观察到气泡边缘的总自由能密度较高而内部较低. 随机抽取单个气泡径向自由能密度如图5(c)所示, 气泡边缘存在总自由能密度峰值, 其值在2.3 × 108—2.7 × 108 J/m3之间, 比气泡内部高出约1/3, 气泡内部总自由能密度在1.6 × 108 —1.8 × 108 J/m3之间, 气泡外总自由能密度趋近于零. 在总自由能密度的各项贡献中, 气泡相自由能密度贡献了最高份额, 在气泡内部, 基质相自由能密度贡献几乎为0而多晶项的贡献为负数, 导致总自由能密度小于气泡相自由能密度. 气泡边缘峰值除了气泡相的贡献, 也明显受到基质相自由能密度贡献峰值的影响.
3.2 温度的影响
本节讨论在10 μm的固定晶粒尺寸下, 不同温度和扩散系数对裂变气体扩散行为的影响, 针对直径为10 μm的晶粒选取的模拟区域为20 μm × 20 μm. 选取3种温度条件: 1276, 1476, 1676 K, 分别对应第2节所述Turnbull模型中3种类型的扩散. 由(20)—(25)式计算得出的温度与扩散系数的对应关系如表2所列. 空位和Xe的初始浓度与3.1节相同.
温度T/K /(m2·s–1)$ {D^{{\text{undoped}}}} $ /(m2·s–1)$ {D^{{\text{doped}}}} $ 1276 8.72 × 10–21 8.72 × 10–21 1476 1.583 × 10–19 2.4593 × 10–19 1676 5.7461 × 10–19 5.0906 × 10–19 T = 1276, 1476, 1676 K条件下晶粒与气泡在τ = 180, 500 s的分布图如图6所示. 从图6可以观察出, 在τ = 180 s时, T = 1276 K条件下气泡尚未形核, 而T = 1476, 1676 K条件下气泡轮廓十分清晰. τ = 500 s时, T = 1276 K条件下晶界上气泡平均直径为0.35 μm, 数密度为1.38个/μm(晶界上气泡数/晶界长度), 晶粒内气泡平均直径为0.32 μm, 数密度为1.96个/μm2 (晶粒内气泡数/晶粒面积); T = 1476 K条件下晶界上气泡平均直径为0.67 μm, 数密度为0.85个/μm, 晶粒内气泡平均直径为0.47 μm, 数密度为17.87个/μm2; T = 1676 K条件下晶界上气泡平均直径为1.34 μm, 数密度为0.23个/μm, 晶粒内气泡平均直径为1.4 μm, 数密度为0.02个/μm2. 相同时刻高温条件下形成的气泡的平均半径明显大于低温形成的. 以上现象表明温度和扩散系数对裂变气体气泡的形核时间、平均半径、数量密度的影响较为显著, 温度越高, 扩散系数越大, 形核越早, 气泡平均直径越大、数量密度越少. 此外, “无气泡区”也随着温度的升高而加宽, 并且T = 1676 K下形成的“无气泡区”的宽度相较于晶粒尺寸而言较大, 宽于有气泡的区域面积, 晶粒内气泡的平均半径和数密度也因高温受到影响.
3种温度下气泡面积占整个模拟区域面积的百分比(以下简称孔隙度)与穿过气泡的晶界占整条晶界长度的百分比(以下简称晶界气泡覆盖率)随时间的演化如图7所示, 在本研究所采用的统计中, 将η ≥ 0.9 的格点均视为气泡相. 从图7(a)可以看出, 孔隙度随时间的演化分为3个阶段: 第1阶段, 维持在零水平; 第2阶段, 快速增长阶段; 第3阶段, 低速平稳增长阶段. 这3个阶段与前文所述的气泡生成的3个阶段: 孕育期, 形核期与长大期相对应, 在图7(a)中以T = 1676 K条件下的曲线为例标记了此3个阶段. 在T = 1676 K条件下, τ = 115 s时气泡开始形核, 模型孔隙度开始增长, 而在 T = 1476, 1276 K条件下则分别在τ = 125 , 160 s开始增长, 温度较高的模型, 孔隙度明显更早进入增长阶段. 由于温度越高扩散系数越大, 空位与裂变气体向晶界扩散得越快, 相同时间内晶界上聚集了更多的空位与裂变气体, 从而更早形核, 孔隙度更早进入增长阶段. 而在平稳增长阶段, 将增长趋势外推至τ = 600 s, T = 1676 K条件下的孔隙度为0.14, T = 1476 K条件下为0.137, T = 1276 K条件下为0.125, 温度较高的模型孔隙度稍高. 如果将温度与扩散系数对孔隙度的影响延伸至对肿胀率的影响, 孔隙度越大的燃料肿胀率越大, 那么该结果意味高温条件下消耗的燃料, 其肿胀率比低温条件下的稍大.
图7(b)表明, 当气泡开始形核时, 晶界气泡覆盖率快速增长, 随后增长逐渐放缓. τ = 600 s, T = 1676 K条件下的晶界气泡覆盖率为0.651, T = 1476 K条件下为0.584, T = 1276 K条件下为0.523. 相较于孔隙度, 高温条件下晶界气泡覆盖率明显大于低温条件下, 高温条件下空位与裂变气体更多地聚集在晶界处, 晶界被气泡覆盖的面积较大. 而晶界气泡覆盖率越大, 意味着裂变气体更容易从燃料芯块沿晶界排放到外界, 即发生裂变气体的释放, 假设当晶界气泡覆盖率达到0.4时裂变气体开始释放, 那么对于T = 1276, 1476, 1676 K条件下, 裂变气体释放的模拟时间分别为τ = 400, 215, 195 s, 对应实际燃耗深度分别为为18.1, 9.7, 8.8 GWd/tU.
在统计晶界上气泡覆盖率时, 采用的方法是在晶界处画一条直线, 统计直线穿过气泡的长度占整条直线长度的百分比, 视为穿过气泡的晶界占整条晶界长度的百分比, 但在演化过程中, 由于存在钉扎效应, 晶界会因气泡的存在而扭曲变形, 造成了统计中波动.
3.3 晶粒尺寸效应
本节讨论在相同的3种温度以及扩散系数条件下, 不同晶粒尺寸对裂变气体扩散的影响. 针对直径为15 μm, 20 μm的晶粒分别选取的模拟区域为30 μm × 30 μm, 40 μm × 40 μm. 选取3个温度参数: 1276, 1476, 1676 K, 对应未掺杂Cr2O3的裂变气体扩散系数如表2所列, 晶粒构造与初始浓度均与3.2节相同.
τ = 250 s, 采用未掺杂$ {D^{{\text{undoped}}}} $, 直径为5, 10, 15, 20 μm的晶粒在1476 K温度下, 气泡的分布如图8所示. 可以看到直径5 μm的晶粒内部没有形成气泡, 而其他晶粒尺寸内部形成气泡, 且分布较为类似.
相同温度下不同晶粒尺寸孔隙度与晶界气泡覆盖率随时间演化见图9. 通过对比可以看到, 如果不改变温度和扩散系数, 仅改变晶粒尺寸, 对燃料孔隙度几乎没有影响. 而晶界气泡覆盖率在某些条件下会受晶粒尺寸的影响, 具体来讲, 从图9(d)可以看出, 在模拟进行一段时间后, 直径5 μm晶粒的晶界气泡覆盖率高于其他晶粒尺寸, 且增长趋势也高于其他晶粒尺寸, 分析认为这是由于直径5 μm晶粒的尺寸较小, 在1476 K温度条件下, 较大的扩散系数使得晶粒内裂变气体全部扩散至晶界, 晶粒内没有发生气泡的形核, 而其他尺寸晶粒则不同(如图8所示), 晶粒内没有气泡捕捉裂变气体原子, 又使得源源不断产生的气体原子扩散到晶界, 被晶界上的气泡捕捉, 使得这些气泡相较于其他晶粒尺寸的晶界气泡, 吸收了更多气体原子, 半径更大, 在τ = 500 s时, 直径5 μm晶粒晶界上的气泡平均直径为0.89 μm, 数密度为0.65个/μm, 而10 μm晶粒晶界上的气泡平均直径为0.64 μm, 数密度为0.85个/μm, 5 μm晶粒晶界上的气泡平均直径与数密度的乘积大于10 μm晶粒, 意味着其晶界气泡覆盖率大于10 μm晶粒, 形成图9(d)所示的趋势. 1476 K温度下的扩散系数是1276 K下的近30倍, 在1276 K下, 5 μm内部同样产生大量气泡(如图3所示), 没有产生类似前文所述的情况, 故在1276 K条件下不同尺寸晶粒的晶界气泡覆盖率差别不大(如图9(b)所示).
此外从图9(f)看出, 1676 K温度下, 5 μm晶粒的晶界气泡覆盖率较小, 不仅小于同温度下其他尺寸晶粒的, 还小于1476 K温度条件下同尺寸晶粒的, 与3.2节的分析结论相悖, 分析后将这一现象同样归咎于5 μm晶粒的尺寸较小. 根据前文的分析, 1476 K温度下5 μm晶粒内部的裂变气体已全部被晶界上的气泡吸收, 1676 K下更大的扩散系数并不能使这些气泡吸收到更多的气体原子, 但能促进形成半径更大的气泡, 同时通过对比图7(a)、图9(c)、图9(e)可以了解到5 μm晶粒在1476 K和1676 K条件下的孔隙度十分近似, 即两者气泡面积十分近似. 晶界长度一致且穿过所有气泡, 气泡面积近似相等且形状近似为圆形, 1676 K下的气泡半径大于1476 K, 通过几何论证可以证明, 1676 K条件下穿过气泡的晶界占整条晶界长度的百分比, 即晶界气泡覆盖率, 小于1476 K条件下的晶界气泡覆盖率. 当晶粒尺寸较大, 1676 K条件下可以比1476 K吸收到更多的裂变气体原子时, 便不会出现上述情况.
3.4 大晶粒掺杂效应
考虑到直径为15, 20 μm的晶粒尺寸介于普通UO2燃料(5—10 μm)与掺杂Cr2O3的大晶粒UO2燃料(50—60 μm)之间, 本节在相同温度下采用掺杂Cr2O3的裂变气体扩散系数$ {D^{{\text{doped}}}} $对该两种尺寸晶粒再进行一次模拟, $ {D^{{\text{doped}}}} $与温度的对应关系如表2所列, 其他条件均与3.3节相同.
采用$ {D^{{\text{undoped}}}} $与$ {D^{{\text{doped}}}} $部分条件下的晶粒与气泡分布对比如图10所示, 可以看出, 相同温度下采用不同扩散系数对裂变气体气泡的演化也会产生影响, 但由于选取的扩散系数相差不大, 所以影响较小.
对于直径15 μm与20 μm晶粒, 孔隙度与晶界气泡覆盖率随时间演化如图11所示. 1476 K与1676 K条件下采用的两种扩散系数相差1.2—1.6倍, 对孔隙度与晶界气泡覆盖率的影响不大, 不同温度下采用的扩散系数相差5—30倍, 此时对孔隙度与晶界气泡覆盖率的影响才比较显著.
4. 结 论
采用Turnbull模型中给出的裂变气体扩散系数, 包括普通UO2燃料中的裂变气体扩散系数与掺杂Cr2O3的大晶粒UO2燃料中的裂变气体扩散系数, 建立了不同晶粒尺寸下晶粒与气泡演化的相场模型, 考察了不同时间、温度、扩散系数条件下燃料模型的晶粒与气泡分布、孔隙度、晶界气泡覆盖率, 得到与实验和其他相场模型相近的结果.
验证了裂变气体气泡形成的3个阶段: 孕育阶段, 形核阶段, 生长阶段. 在孕育阶段, 没有气泡产生, 孔隙度与晶界气泡覆盖率不增大; 在形核阶段, 气泡大量形成且优先在晶界上形成, 孔隙度与晶界气泡覆盖率高速增长; 在生长阶段, 气泡停止产生, 现有气泡不断生长并与临近气泡相互融合, 孔隙度与晶界气泡覆盖率平稳增长. 此外还观察到了较明显的“无气泡区”, 以及气泡边缘存在能量峰值.
研究结果表明, 温度与扩散系数对裂变气体的释放行为产生较大影响, 温度越高、扩散系数越大, 裂变气体气泡形核越早, 平均半径越大且数密度越小. 同时在晶粒半径较大的模型中, 扩散系数越大, 孔隙度与晶界气泡覆盖率越高, 晶粒尺寸对这两者的影响较小, 但在晶粒半径较小的模型中则不同, 较小的晶粒尺寸会对孔隙度与晶界气泡覆盖率产生较大影响.
在进一步的工作中, 可通过孔隙度、气泡半径及密度等数据进而计算燃料的热导率、肿胀率等参数, 从而研究裂变气体释放行为对燃料的导热性能、力学性能的影响. 此外根据文献[36–38]的相关研究, 晶粒中存在应力分布和弹性能, 该应力可能会影响气泡的生长和形状, 促使气泡成长为透镜状, 而在本文中, 气泡多为圆形, 故在之后的研究中, 需要考虑弹塑性能对气泡演化的影响.
[1] Rest J 2010 J. Nucl. Mater. 402 179Google Scholar
[2] Trinkaus H, Singh B N 2003 J. Nucl. Mater. 323 229Google Scholar
[3] Rest J, Hofman G L 1999 Nucl. Technol. 126 88Google Scholar
[4] Pastore G, Luzzi L, Di Marcello V, van Uffelen P 2013 Nucl. Eng. Des. 256 75Google Scholar
[5] Piro M H, Sunderland D, Livingstone S, Sercombe J, Revie R W, Quastel A, Terrani K A, Judge C 2020 Comprehensive Nuclear Materials (2nd Ed.) (Amsterdam: Elsevier) p248
[6] 何文, 伍晓勇, 吴璐, 温榜, 朱伟, 张伟, 潘荣剑, 王桢, 黄伟杰 2017 核动力工程 38 170Google Scholar
He W, Wu X Y, Wu L, Wen B, Zhu W, Zhang W, Pan R J, Wang Z, Huang W J 2017 Nucl. Power Eng. 38 170Google Scholar
[7] Killeen J C 1980 J. Nucl. Mater. 88 177Google Scholar
[8] Cooper M W D, Pastore G, Che Y, Matthews C, Forslund A, Stanek C R, Shirvan K, Tverberg T, Gamble K A, Mays B, Andersson D A 2021 J. Nucl. Mater. 545 152590Google Scholar
[9] Aagesen L K, Schwen D, Tonks M R, Zhang Y 2019 Comput. Mater. Sci. 161 35Google Scholar
[10] Yuda R, Harada H, Hirai M, Hosokawa T, Une K, Kashibe S, Shimizu S, Kubo T 1997 J. Nucl. Mater. 248 262Google Scholar
[11] 庞华, 辛勇, 岳慧芳, 彭航, 蒲曾坪, 邱玺, 孙志鹏, 刘仕超 2022 材料导报 36 5Google Scholar
Pang H, Xin Y, Yue H F, Peng H, Pu Z P, Qiu X, Sun Z P, Liu S C 2022 Mater. Rev. 36 5Google Scholar
[12] Delafoy C, Dewes P, Miles T 2007 Proceedings of the 2007 LWR Fuel Performance Meeting/TopFuel 2007 San Francisco, CA, United States, September 30–October 3, 2007 p1
[13] Kashibe S, Une K 1998 J. Nucl. Mater. 254 234Google Scholar
[14] Che Y, Pastore G, Hales J, Shirvan K 2018 Nucl. Eng. Des. 337 271Google Scholar
[15] Moelans N, Blanpain B, Wollants P 2008 Comput. Coupling Phase Diagrams Thermochem 32 268Google Scholar
[16] 郭灿, 康晨瑞, 高莹, 张一弛, 邓英远, 马超, 徐春杰, 梁淑华 2022 物理学报 71 096401Google Scholar
Guo C, Kang C R, Gao Y, Zhang Y C, Deng Y Y, Ma C, Xu C J, Liang S H 2022 Acta Phys. Sin. 71 096401Google Scholar
[17] Hu S, Henager Jr C H 2009 J. Nucl. Mater. 394 155Google Scholar
[18] Millett P C, El-Azab A, Rokkam S, Tonks M, Wolf D 2011 Comput. Mater. Sci. 50 949Google Scholar
[19] Jiang Y B, Liu W B, Li W J, Sun Z Y, Xin Y, Chen P H, Yun D 2021 Comput. Mater. Sci. 188 110176.Google Scholar
[20] Zhao J J, Sun D, Xi L, Chen P, Zhao J J, Wang Y Y 2023 Phys. Chem. Chem. Phys. 25 14928Google Scholar
[21] Cahn J W, Hilliard J E 1958 J. Chem. Phys. 28 258Google Scholar
[22] Li Y, Hu S, Montgomery R, Gao F, Sun X 2013 Nucl. Instrum. Methods Phys. Res. , Sect. B 303 62Google Scholar
[23] Kittel C, Kroemer H 1980 Thermal Physics (New York: WH Freeman and Company) pp287–306
[24] Moelans N 2011 Acta Mater. 59 1077Google Scholar
[25] Moelans N, Blanpain B, Wollants P 2008 Phys. Rev. B 78 024113Google Scholar
[26] Cahn J W 1961 Acta Metall. 9 795Google Scholar
[27] Allen S M, Cahn J W 1972 Acta Metall. 20 423Google Scholar
[28] Allen S M, Cahn J W 1973 Scr. Metall. 7 1261Google Scholar
[29] Turnbull J A, Friskney C A, Findlay J R, Johnson F A, Walter A J 1982 J. Nucl. Mater. 107 168Google Scholar
[30] Turnbull J A, White R J, Wise C 1989 Technical Committee on Water Reactor Fuel Element Computer Modelling in Steady State, Transient and Accident Conditions Preston, United Kingdom, September 18–22, 1988 p174
[31] INTRODUCTION TO COMSOL Multiphysics, COMSOL Co Ltd. https://cdn.comsol.com/doc/5.2/IntroductionToCOMSOLMultiphysics.zh_CN.pdf [2023-11-1]
[32] Millett P C, El-Azab A, Wolf D 2011 Comput. Mater. Sci. 50 960Google Scholar
[33] Sagui C, Grant M 1999 Phys. Rev. E 59 4175.Google Scholar
[34] Bullough R, Nelson R S 1974 Phys. Technol. 5 29Google Scholar
[35] Zacharie I, Lansiart S, Combette P, Trotabas M, Coster M, Groos M 1998 J. Nucl. Mater. 255 85Google Scholar
[36] Sheng J, Wang Y C, Liu Y, Wu S, Xu K, Chen Z H, Bo S, Liu H F, Song H F 2022 Comput. Mater. Sci. 213 111663Google Scholar
[37] Wu S, Sheng J, Yang C, Shi X, Huang H, Liu Y, Song H 2022 Front. Mater. 9 916593Google Scholar
[38] 姜彦博, 柳文波, 孙志鹏, 喇永孝, 恽迪 2022 物理学报 71 026103Google Scholar
Jiang Y B, Liu W B, Sun Z P, La Y X, Yun D 2022 Acta Phys. Sin. 71 026103Google Scholar
-
图 2 自由能密度随浓度(cv, cg)分布的变化, 曲线最低点的横坐标为平衡浓度, 随着cg增加, $ c_{\text{v}}^{{\text{b, eq}}} $逐渐减小, f bubble,v的图像(蓝色虚线)将向左平移
Figure 2. Variation of free energy density with concentration (cv, cg) distribution, the abscissa of the lowest point of the curve is the equilibrium concentration. As cg increases, $ c_{\text{v}}^{{\text{b, eq}}} $ gradually decreases, and the image of f bubble,v (blue dashed line) will shift to the left.
图 3 T = 1276 K条件下直径5 μm晶粒与气泡随时间演化分布 (a) τ = 0 s; (b) τ = 140 s; (c) τ = 160 s; (d) τ = 200 s; (e) τ = 300 s; (f) τ = 550 s
Figure 3. Distribution of grain with a diameter of 5 μm and bubble evolution over time at T = 1276 K: (a) τ = 0 s; (b) τ = 140 s; (c) τ = 160 s; (d) τ = 200 s; (e) τ = 300 s; (f) τ = 550 s.
图 4 空位(cv, (a), (c))与裂变气体(cg, (b), (d))在τ = 140 s与τ = 550 s的浓度分布(不同的颜色代表浓度的取值) (a), (b) τ = 140 s; (c), (d) τ = 550 s
Figure 4. Concentration distribution of vacancies (cv, (a) and (c)) and fission gases (cg, (b) and (d)) at τ = 140 s and τ = 550 s, different color represent the value of concentration: (a), (b) τ = 140 s; (c), (d) τ = 550 s.
图 5 (a) 整个模拟区域内的平均自由能密度随时间的变化; (b) 总自由能密度在整个模拟区域分布, 颜色栏为取值范围; (c)某一气泡径向自由能密度分布
Figure 5. (a) Variation of average free energy density over time in simulation area; (b) distribution of the total free energy density in simulation area, the color bar represents the range of values; (c) radial free energy density distribution of a certain bubble.
图 9 相同温度不同晶粒尺寸, (a), (c), (e)孔隙度与(b), (d), (f)晶界气泡覆盖率随时间演化 (a), (b) T = 1276 K; (c), (d) T = 1476 K; (e), (f) T = 1676 K
Figure 9. Evolution of (a), (c), (e) porosity and (b), (d), (f) bubble coverage on GB over time for the same temperature but different grain sizes: (a), (b) T = 1276 K; (c), (d) T = 1476 K; (e), (f) T = 1676 K.
图 6 T = 1276, 1476, 1676 K条件下直径10 μm晶粒与气泡分布图 (a) T = 1276 K, τ = 180 s; (b) T = 1476 K, τ = 180 s; (c) T = 1676 K, τ = 180 s; (d) T = 1276 K, τ = 500 s; (e) T = 1476 K, τ = 500 s; (f) T = 1676 K, τ = 500 s
Figure 6. Distribution of grains with a diameter of 10 μm and bubbles under T = 1276, 1476, 1676 K: (a) T = 1276 K, τ = 180 s; (b) T = 1476 K, τ = 180 s; (c) T = 1676 K, τ = 180 s; (d) T = 1276 K, τ = 500 s; (e) T = 1476 K, τ = 500 s; (f) T = 1676 K, τ = 500 s.
图 10 采用$ {D^{{\text{undoped}}}} $与$ {D^{{\text{doped}}}} $晶粒与气泡分布对比 (a) T = 1676 K, 15 μm, $ {D^{{\text{undoped}}}} $; (b) T = 1676 K, 15 μm, $ {D^{{\text{doped}}}} $; (c) T = 1476 K, 20 μm, $ {D^{{\text{undoped}}}} $; (d) T = 1476 K, 20 μm, $ {D^{{\text{doped}}}} $
Figure 10. Comparison of the distribution of grains and bubbles using $ {D^{{\text{undoped}}}} $ and $ {D^{{\text{doped}}}} $: (a) T = 1676 K, 15 μm, $ {D^{{\text{undoped}}}} $; (b) T = 1676 K, 15 μm, $ {D^{{\text{doped}}}} $; (c) T = 1476 K, 20 μm, $ {D^{{\text{undoped}}}} $; (d) T = 1476 K, 20 μm, $ {D^{{\text{doped}}}} $.
表 1 模拟采用的部分参数
Table 1. Parameters used in simulation.
参数 符号 值 玻尔兹曼常量/(J·K–1) kB 1.3806 × 10–23 理想气体常数/(J·mol–1·K–1) R 8.3145 UO2晶胞U原子体积/ nm–³ Va 0.0409 空位形成能/eV $ E_{\text{v}}^{\text{f}} $ 5.1 气体原子缺陷形成能/eV $ E_{\text{g}}^{\text{f}} $ 10.31 扩散界面系数 aGB 1.2 as 0.8 自由能势垒系数/(J·m–3) m 3.0 × 107 梯度项系数/(J·m–3) κϕ 3.38 × 10–7 κv, κg, κη 1.69 × 10–6 迁移率/(m3·J–1·s–1) Lη, Lϕ 1.56 × 10–11 裂变率密度/(次裂变·m3·s–1) $ \dot F $ 1.09 × 1019 Xe产额 Y 0.27 表 2 不同温度下采用的扩散系数
Table 2. Diffusion coefficients used at different temperatures.
温度T/K $ {D^{{\text{undoped}}}} $/(m2·s–1) $ {D^{{\text{doped}}}} $/(m2·s–1) 1276 8.72 × 10–21 8.72 × 10–21 1476 1.583 × 10–19 2.4593 × 10–19 1676 5.7461 × 10–19 5.0906 × 10–19 -
[1] Rest J 2010 J. Nucl. Mater. 402 179Google Scholar
[2] Trinkaus H, Singh B N 2003 J. Nucl. Mater. 323 229Google Scholar
[3] Rest J, Hofman G L 1999 Nucl. Technol. 126 88Google Scholar
[4] Pastore G, Luzzi L, Di Marcello V, van Uffelen P 2013 Nucl. Eng. Des. 256 75Google Scholar
[5] Piro M H, Sunderland D, Livingstone S, Sercombe J, Revie R W, Quastel A, Terrani K A, Judge C 2020 Comprehensive Nuclear Materials (2nd Ed.) (Amsterdam: Elsevier) p248
[6] 何文, 伍晓勇, 吴璐, 温榜, 朱伟, 张伟, 潘荣剑, 王桢, 黄伟杰 2017 核动力工程 38 170Google Scholar
He W, Wu X Y, Wu L, Wen B, Zhu W, Zhang W, Pan R J, Wang Z, Huang W J 2017 Nucl. Power Eng. 38 170Google Scholar
[7] Killeen J C 1980 J. Nucl. Mater. 88 177Google Scholar
[8] Cooper M W D, Pastore G, Che Y, Matthews C, Forslund A, Stanek C R, Shirvan K, Tverberg T, Gamble K A, Mays B, Andersson D A 2021 J. Nucl. Mater. 545 152590Google Scholar
[9] Aagesen L K, Schwen D, Tonks M R, Zhang Y 2019 Comput. Mater. Sci. 161 35Google Scholar
[10] Yuda R, Harada H, Hirai M, Hosokawa T, Une K, Kashibe S, Shimizu S, Kubo T 1997 J. Nucl. Mater. 248 262Google Scholar
[11] 庞华, 辛勇, 岳慧芳, 彭航, 蒲曾坪, 邱玺, 孙志鹏, 刘仕超 2022 材料导报 36 5Google Scholar
Pang H, Xin Y, Yue H F, Peng H, Pu Z P, Qiu X, Sun Z P, Liu S C 2022 Mater. Rev. 36 5Google Scholar
[12] Delafoy C, Dewes P, Miles T 2007 Proceedings of the 2007 LWR Fuel Performance Meeting/TopFuel 2007 San Francisco, CA, United States, September 30–October 3, 2007 p1
[13] Kashibe S, Une K 1998 J. Nucl. Mater. 254 234Google Scholar
[14] Che Y, Pastore G, Hales J, Shirvan K 2018 Nucl. Eng. Des. 337 271Google Scholar
[15] Moelans N, Blanpain B, Wollants P 2008 Comput. Coupling Phase Diagrams Thermochem 32 268Google Scholar
[16] 郭灿, 康晨瑞, 高莹, 张一弛, 邓英远, 马超, 徐春杰, 梁淑华 2022 物理学报 71 096401Google Scholar
Guo C, Kang C R, Gao Y, Zhang Y C, Deng Y Y, Ma C, Xu C J, Liang S H 2022 Acta Phys. Sin. 71 096401Google Scholar
[17] Hu S, Henager Jr C H 2009 J. Nucl. Mater. 394 155Google Scholar
[18] Millett P C, El-Azab A, Rokkam S, Tonks M, Wolf D 2011 Comput. Mater. Sci. 50 949Google Scholar
[19] Jiang Y B, Liu W B, Li W J, Sun Z Y, Xin Y, Chen P H, Yun D 2021 Comput. Mater. Sci. 188 110176.Google Scholar
[20] Zhao J J, Sun D, Xi L, Chen P, Zhao J J, Wang Y Y 2023 Phys. Chem. Chem. Phys. 25 14928Google Scholar
[21] Cahn J W, Hilliard J E 1958 J. Chem. Phys. 28 258Google Scholar
[22] Li Y, Hu S, Montgomery R, Gao F, Sun X 2013 Nucl. Instrum. Methods Phys. Res. , Sect. B 303 62Google Scholar
[23] Kittel C, Kroemer H 1980 Thermal Physics (New York: WH Freeman and Company) pp287–306
[24] Moelans N 2011 Acta Mater. 59 1077Google Scholar
[25] Moelans N, Blanpain B, Wollants P 2008 Phys. Rev. B 78 024113Google Scholar
[26] Cahn J W 1961 Acta Metall. 9 795Google Scholar
[27] Allen S M, Cahn J W 1972 Acta Metall. 20 423Google Scholar
[28] Allen S M, Cahn J W 1973 Scr. Metall. 7 1261Google Scholar
[29] Turnbull J A, Friskney C A, Findlay J R, Johnson F A, Walter A J 1982 J. Nucl. Mater. 107 168Google Scholar
[30] Turnbull J A, White R J, Wise C 1989 Technical Committee on Water Reactor Fuel Element Computer Modelling in Steady State, Transient and Accident Conditions Preston, United Kingdom, September 18–22, 1988 p174
[31] INTRODUCTION TO COMSOL Multiphysics, COMSOL Co Ltd. https://cdn.comsol.com/doc/5.2/IntroductionToCOMSOLMultiphysics.zh_CN.pdf [2023-11-1]
[32] Millett P C, El-Azab A, Wolf D 2011 Comput. Mater. Sci. 50 960Google Scholar
[33] Sagui C, Grant M 1999 Phys. Rev. E 59 4175.Google Scholar
[34] Bullough R, Nelson R S 1974 Phys. Technol. 5 29Google Scholar
[35] Zacharie I, Lansiart S, Combette P, Trotabas M, Coster M, Groos M 1998 J. Nucl. Mater. 255 85Google Scholar
[36] Sheng J, Wang Y C, Liu Y, Wu S, Xu K, Chen Z H, Bo S, Liu H F, Song H F 2022 Comput. Mater. Sci. 213 111663Google Scholar
[37] Wu S, Sheng J, Yang C, Shi X, Huang H, Liu Y, Song H 2022 Front. Mater. 9 916593Google Scholar
[38] 姜彦博, 柳文波, 孙志鹏, 喇永孝, 恽迪 2022 物理学报 71 026103Google Scholar
Jiang Y B, Liu W B, Sun Z P, La Y X, Yun D 2022 Acta Phys. Sin. 71 026103Google Scholar
Catalog
Metrics
- Abstract views: 2473
- PDF Downloads: 70
- Cited By: 0