复杂多相体系陶瓷型核燃料的致密化烧结过程的模拟方法

摘要

本发明公开了一种复杂多相体系陶瓷型核燃料的致密化烧结过程的模拟方法，通过第一性原理和实验的方法得到烧结过程中所需材料物性参数，考虑晶粒、颗粒和气孔三相相互作用，构建含第二相颗粒多晶多相烧结体系的自由能函数，通过相场变量平衡解的方法确定材料物性参数与自由能参数的关系得到总自由能函数并进行变分，同时考虑烧结过程中晶粒的刚体运动，建立修正的演化方程，采用显示Euler法和有限差分法求解，实现对含第二相颗粒的多晶多相烧结致密化过程、晶界与颗粒或气孔相互作用的精准模拟。通过本发明，解决了以往相场模型缺少致密化与第二相颗粒同时作用的问题；扩展了相场模型的适用范围，提出一种多晶多相体系陶瓷型核燃料烧结的新思路。

说明书

技术领域

本发明涉及计算材料学领域。具体而言，本发明涉及一种用于复杂多相体系陶瓷型核燃料的致密化烧结过程的模拟方法。

背景技术

相场模型基于扩散界面模型，可以有效地追踪烧结过程中的复杂界面的演变过程，进而高效地揭示其烧结微观组织演变，近年来逐渐被大量应用于烧结过程的分析和研究。随着烧结技术不断进步，烧结过程逐渐变的更加复杂，除了固相与气相的相互作用，还存在引入第二相颗粒作为催化剂或不可避免夹杂第二相颗粒的作用机制。然而，以往的相场模型只考虑了两相相互作用，并没有同时考虑气孔与第二相颗粒对晶界的相互作用以及对烧结过程中晶粒长大的影响。以往的相场模型也并未考虑烧结过程气孔收缩致密化的作用，导致无法保证烧结模拟过程的精准性。目前，在多晶多相体系陶瓷型核燃料烧结过程缺少适用性较好的模型。

现有技术方法，主要解决如下技术问题：

相场模型中系统自由能函数包含局部能量密度函数与界面梯度能函数，界面梯度能由特定公式积分得到，因此，局部能量密度函数的构造对于烧结相场模拟至关重要。在相场模型一直以来的实践中，存在由晶粒长大模型发展而来的气孔与晶界相互作用的局部能量密度函数，也存在第二相颗粒与晶界相互作用的局部能量密度函数，但都没有同时考虑局部能量密度函数中晶粒、颗粒、气孔三相相互作用。然而，烧结是一个复杂多相的过程，现有的模型无法对复杂多相体系烧结过程作较准确的描述。

以往烧结相场模拟的扩散机制主要以表面扩散机制开展，然而以表面扩散主导的相场模型并不能适用模拟烧结后期的致密化过程，烧结致密化的机理有待进一步揭示。

发明内容

为了克服上述现有技术存在的问题，本发明提出一种复杂多相体系陶瓷型核燃料的致密化烧结过程的模拟方法，用于解决或至少部分解决相场模型在含第二相颗粒条件下模拟气孔收缩致密化的问题。

为了解决上述技术问题，本发明采用如下技术方案：

一种用于复杂多相体系陶瓷型核燃料的致密化烧结过程的模拟方法，包括：

S1：通过第一性原理与实验的方法获取需要的材料物性参数，其中材料物性参数包括材料晶界能、表面能、晶界扩散率、表面扩散率和晶格扩散率；

S2：根据步骤S1获取的材料物性参数，考虑晶粒、颗粒和气孔三相共存相互作用，采用相场变量进行多相耦合，构建含第二相颗粒多晶多相烧结体系的总自由能函数，总自由能函数包括局部自由能密度函数与界面梯度自由能；推导一维无限半双晶体系与一维无限半单晶体系的系统自由能，通过相场变量平衡解的方法确定材料物性参数与自由能参数的关系式；再根据质量守恒，建立浓度场演化Cahn-Hilliard方程、取向场演化Allen-Cahn方程描述晶粒、第二相颗粒与气孔的相互作用；

S3：根据步骤S2构建的一组演化方程，考虑多晶多相烧结过程中晶粒平移旋转的刚体运动，引入描述粒子刚体运动的平流通量项，得到修正的演化方程，修正的演化方程与步骤S2构建的总自由能函数形成了多晶多相相场模型；烧结过程中，描述粒子刚体运动的平流通量项由平移和转动两部分组成；为了便于计算机模拟计算，对多晶多相相场模型的变量和常数进行无量纲化；

S4：将步骤S1中获取的材料物性参数代入步骤S2所确定的材料物性参数与自由能参数的关系式中得到总自由能函数表达式，将步骤S1中获取的材料物性参数代入步骤S3所修正的演化方程中得到演化方程表达式，输入模拟采用的边界条件、进行网格划分；根据多晶多相烧结体系模拟所需的微观结构随机生成紧密的单峰分布、双峰分布的初始化微观组织；通过演化方程考虑烧结过程中表面扩散、晶界扩散、晶格扩散，以及过饱和空位被晶界和自由表面吸收过程，使用显式Euler法和有限差分法对步骤S3中修正的演化方程进行求解；

S5：将步骤S4中修正的演化方程的求解结果进行可视化处理，获得多晶多相烧结过程中晶粒、气孔、颗粒微观组织的演化图像；根据可视化的微观组织结构，跟踪内部气孔收缩度的时间演变来监控烧结动力学过程；统计烧结颈的形成增长尺寸、平衡态时晶界二面角、气孔的尺寸、晶粒等效半径、界面曲率；

S6：根据步骤S5的演化图像，观察分析烧结模拟过程中气孔随晶界移动、颗粒随晶界移动的演变，探究烧结过程中晶界与气孔、晶界与第二相颗粒的相互作用规律；根据步骤S5的统计结果，计算气孔收缩度随时间的演化，探究烧结致密化气孔收缩的规律。

在一种实施方式中，S2建立多晶多相体系相场模型具体包括：

S2.1：根据步骤S1获取的材料物性参数，考虑晶粒、颗粒和气孔三相共存相互作用，采用相场变量进行多相耦合，构建含第二相颗粒多晶多相烧结体系的总自由能函数，总自由能函数包括局部自由能密度函数与界面梯度自由能：

总自由能函数表达式为：

式中，F为多晶多相烧结体系总自由能；f(ρ,η

S2.2：推导一维无限半双晶体系与一维无限半单晶体系的系统自由能，通过相场变量平衡解的方法确定材料物性参数与自由能参数的关系式：

在一维无限半双晶体系中晶界能的积分为：

式中，γ

再由多元函数的偏微分得：

对于一维无限双晶体系，边界条件为：

式(7)的边界条件代入式(6)求解得：

由式(3)与式(8)得：

将式(9)变形得：

根据双晶几何结构的对称性η

晶界宽度由相场变量变化的长度来表述，在一维无限半双晶体系中的积分为：

式中，δ

在一维无限半双晶体系中气孔表面能的积分为：

式中，γ

再由多元函数的偏微分得：

对于一维无限半单晶体系，边界条件为：

式(18)的边界条件代入式(17)求解得：

由式(14)与式(19)可得：

将式(19)变形得：

在气孔界面处，由ρ＝η

气孔表面宽度由相场变量变化的长度来表述，在一维无限半单晶体系中的积分为：

式中，δ

假定晶界与气孔表面有相同的宽度，即δ

化简得：

自由能参数与材料物性参数的关系如下：

式中，δ为扩散界面宽度；

S2.3：根据质量守恒，建立浓度场演化Cahn-Hilliard方程、取向场演化Allen-Cahn方程描述晶粒、第二相颗粒与气孔的相互作用：

浓度场演化Cahn-Hilliard方程：

式中，M为迁移率张量函数，表达式为：

式中，v

D

式中，D

式中，I表示单位张量；符号

D

式中，D

式中，n

D

D

ξ(ρ)＝ρ

式中，D

式中，L为表征晶界迁移率的系数，在各向同性假设下为常数，与材料的具体性质关系如下：

式中，M

在一种实施方式中，S3所述的方法具体包括：

S3.1：根据步骤S2构建的一组演化方程，考虑多晶多相烧结过程中晶粒平移旋转的刚体运动，引入描述粒子刚体运动的平流通量项，得到修正的演化方程：

式中，r表示空间中的位置；v

S3.2：烧结过程中，描述粒子刚体运动的平流通量项由平移和转动两部分组成：

平流通量项中的平流速度表示为：

式中，v

S3.3：为了便于计算机模拟计算，对多晶多相相场模型的变量和常数进行无量纲化：

根据浓度场演化Cahn-Hilliard方程含有的变量为t,M,F,

由此可得迁移率张量M的无量纲化为

式中，l

在一种实施方式中，在对演化方程求解时，S4具体包括：

S4.1：将步骤S1获取的晶界能、表面能代入步骤S2所确定的材料物性参数与自由能参数的关系式中得到总自由能函数表达式，将晶界迁移率、表面扩散率、晶格扩散率代入步骤S3所修正的演化方程中得到演化方程表达式，输入周期性边界条件，进行网格划分；

S4.2：根据多晶多相烧结体系模拟所需的微观结构随机生成紧密的单峰分布、双峰分布的初始化微观组织：

对于多晶多相模拟，基于正太分布随机生成晶粒半径单峰分布、双峰分布的晶粒颗粒分布，具体为：

(1)在划分的网格中按照掺杂比例与晶粒尺寸分布随机位置生成互不重叠的晶粒；

(2)将第(1)步生成的晶粒先向x轴方向平移至晶粒紧密接触，再将晶粒向y轴方向平移至晶粒紧密接触，重复操作直至晶粒位置无明显变化；

(3)在平移后的网格空白区域重复第(1)步、第(2)步至晶粒分布无明显变化；

(4)在第(3)步结束后的网格空白区域按照颗粒掺杂比例与颗粒尺寸分布随机生成互不重叠的颗粒；

(5)将第(4)步生成的颗粒先向x轴方向平移至颗粒紧密接触，再将颗粒向y轴方向平移至颗粒紧密接触，重复操作直至颗粒位置无明显变化；

(6)在平移后的网格空白区域重复第(4)步、第(5)步至颗粒分布无明显变化；

S4.3：通过演化方程考虑烧结过程中表面扩散、晶界扩散、晶格扩散，以及过饱和空位被晶界和自由表面吸收等过程，使用显式Euler法和有限差分法对步骤S3中修正的演化方程进行求解；

在一种实施方式中，S5具体包括：

S5.1：将步骤S4中修正的演化方程的求解结果进行可视化处理，获得多晶多相烧结过程中晶粒、气孔、颗粒微观组织演化的图像：

表征烧结模拟过程可视化变量：

式中，

S5.2：根据可视化的微观组织结构，跟踪内部气孔收缩度的时间演变来监控烧结动力学过程，其定义如下：

式中，A(0)和A(t)分别为在时间0和时间t时刻的气孔面积；R(0)和R(t)分别为在时间0和时间t时刻的气孔等效半径；

S5.3：统计烧结颈的形成增长尺寸，烧结过程中烧结颈的增长曲线满足规律：

式中，l和D分别表示烧结颈的长度和晶粒的直径；α和β为常数；n为生长指数；

计算平衡态时晶界二面角，平衡态时晶界二面角方程为：

界面处的局部曲率由下式给出：

式中，n是晶粒i在晶界处的法线，第i个晶粒的曲率K用序参量显式表达为：

式中，下标x和y表示导数。

在一种实施方式中，当采用Euler法求解时，施加周期性边界条件，将刚体运动部分与相场部分循环求解，直到控制指标满足误差要求为止。

本申请实施例中上述一个或多个技术方案，至少具有如下一种或多种技术效果：

通过第一性原理计算与实验获得烧结数值模拟所需要的材料物性参数，建立热力学自由能函数，将材料物性参数代入建立的多晶多相系统自由能函数中，得到总自由能函数表达式，对总自由能函数进行变分并考虑晶粒平移旋转刚体运动，采用显式Euler法求解修正的演化方程，实现对含第二相颗粒的烧结致密化过程、晶界与颗粒或气孔相互作用的精准模拟。

通过本发明，极大的拓展了相场模型的适用范围，解决了以往相场模型缺少含第二相颗粒与晶粒气孔同时相互作用的问题，补充了传统相场烧结模拟后期缺少气孔收缩致密化的空缺。

进一步地，提出热力学自由能函数耦合第二相颗粒作用的新方法，首次开发了晶粒、颗粒、气孔三相相互作用的相场模型。

进一步地，将烧结致密化过程机理加入相场模型演化方程中，进一步拓展了相场模型的适用范围。

进一步地，采用张量形式同时考虑了3种各向异性地扩散机制(表面扩散、晶界扩散和晶格扩散)，使得相场模型扩散机制更加全面。

附图说明

图1是本发明的模拟方法流程图。

图2(a)为模拟步长2×10

图2(b)为模拟步长70×10

图2(c)为模拟步长350×10

图2(d)为模拟步长400×10

图3(a)为模拟步长1×10

图3(b)为模拟步长25×10

图3(c)为模拟步长50×10

图3(d)为模拟步长100×10

图4(a)为模拟步长1×10

图4(b)为模拟步长40×10

图4(c)为模拟步长80×10

图4(d)为模拟步长120×10

图4(e)为模拟步长160×10

图4(f)为模拟步长200×10

图5(a)为模拟步长1×10

图5(b)为模拟步长10×10

图5(c)为模拟步长20×10

图5(d)为模拟步长30×10

图5(e)为模拟步长40×10

图5(f)为模拟步长50×10

图6不同颗粒体积分数下平均晶粒尺寸演变图。

具体实施方式

下面结合附图和具体实施例对本发明作进一步详细说明。

本发明的主要构思如下：

如图1所示，通过第一性原理和实验的方法得到烧结过程中所需的材料物性参数，构建含第二相颗粒多晶多相烧结体系的总自由能函数，将材料参数代入建立的多晶多相烧结体系的总自由能函数中，得到总自由能函数表达式，对总自由能函数进行变分并考虑多晶烧结过程中晶粒平移旋转刚体运动，采用显式Euler法求解修正的演化方程，实现对含第二相颗粒的多晶多相烧结致密化过程、晶界与颗粒或气孔相互作用的精准模拟。通过本发明，极大的拓展了相场模型的适用范围，解决了以往模型缺少含第二相颗粒与晶粒气孔同时相互作用的问题，补充了传统相场烧结模拟后期缺少气孔收缩致密化的空缺，提出一种多晶多相体系陶瓷型核燃料烧结的新思路。

为了使本发明的目的技术方案及优点更加清晰，以下结合附图及其实例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用解释本发明，并不局限于本发明。

(1)构建含第二相颗粒多晶多相烧结体系的总自由能函数

烧结过程中相场模型体系的总自由能泛函数F采用以下形式：

式中，F为多晶多相烧结体系总自由能；f(ρ,η

(2)建立平移旋转刚体运动平流通量项

烧结过程中粒子的刚体运动由平移和转动两部分组成，平流速度表示为：

式中，v

(3)考虑烧结过程颗粒平移旋转的刚体运动修正演化方程

为了描述烧结过程的气孔收缩和致密化，在演化方程中加入平流通量项：

式中，r表示空间中的位置；v

式中，v

(4)多晶多相相场模型无量纲化

为了便于计算机模拟计算，对多晶多相相场模型的变量和常数进行无量纲化。根据浓度场演化Cahn-Hilliard方程含有的变量为t,M,F,

由此可得迁移率张量M的无量纲化为

式中，l

(5)可视化处理并统计模拟结果

表征烧结模拟过程可视化变量：

式中，

可以跟踪内部气孔收缩度的时间演变来监控烧结动力学过程，其定义如下：

式中，A(0)和A(t)分别为在时间0和时间t时刻的气孔面积；R(0)和R(t)分别为在时间0和时间t时刻的气孔等效半径。

为了使本发明的目的、技术方案及优点更加清楚，以下结合附图及实施例，对本发明进行进一步详细说明。

本发明实施例提供了一种用于复杂多相体系陶瓷型核燃料的致密化烧结过程的模拟方法，按照具体实施方式中的流程，利用相场模型对二氧化铀陶瓷粉末在2000K温度下的烧结过程进行了数值模拟，模拟使用参数如表1所示：

表1二氧化铀的物理参数

注:D

构建含第二相颗粒多晶多相烧结体系的总自由能函数，自由能函数表达式为：

式中，F为多晶多相烧结体系总自由能；f(ρ,η

建立烧结过程描述气孔相的演化方程：

建立烧结过程描述晶粒相的演化方程：

式中，M为迁移率张量函数，包括表面扩散、晶界扩散和晶格扩散；L为表征晶界迁移率的系数，在各向同性假设下为常数；r表示空间中的位置；v

输入获取的材料物性参数、模拟的边界条件并划分网络，为了能够准确的模拟烧结过程的演化图像，网格尺寸应该足够细，网格尺寸与模拟尺寸的比值不得大于0.5；通过显式Euler法与有限差分法求解演化方程，并对数值结果使用Matlab或Paraview可视化处理，其中，图中颜色深浅代表相场变量的值；晶粒与颗粒气孔相互作用演化图像如图2(a-d)所示，晶界收缩受曲率影响，颗粒与气孔先促进后抑制晶界移动；四颗粒气孔收缩致密化过程演化图像如图3(a-d)所示，致密化气孔收缩最终形成稳定的三叉晶界；多晶烧结致密化过程演化图像如图4(a-f)所示，演化模拟结果与实验观察结果基本一致；多晶多颗粒多孔烧结过程演化图像如图5(a-f)所示，不同颗粒掺杂对烧结过程平均晶粒尺寸演变影响如图6所示，适量的颗粒掺杂量提高了晶粒的生长速度。

综上所述，本发明是一种结构严谨、实施流程清晰且适用范围较广的模拟方法，该发明可为复杂多晶多相核燃料致密化烧结过程提供理论与技术支撑，实现多晶多相陶瓷核燃料烧结致密化过程、晶界与颗粒或气孔相互作用的精准模拟。

与现有技术相比，本发明的复杂多相体系陶瓷型核燃料的致密化烧结过程的模拟方法考虑了晶粒、气孔、颗粒三相共存相互作用，创建含第二相颗粒多晶多相烧结体系的总自由能函数，解决了以往烧结模型无法考虑三相同时作用的问题；创建同时考虑表面扩散机制与晶粒平移旋转刚体运动的混合机制，解决了烧结后期致密化的问题，保证了烧结前期粘合与后期致密化整个过程的模拟演变。