一种空位团簇跳出超势能谷的模拟方法

技术领域

本发明属于核材料辐照损伤模拟技术领域，具体涉及一种空位团簇跳出超势能谷的模拟方法。

背景技术

辐照损伤是由于中子、带电粒子或电磁波等和固体材料的点阵原子发生一系列碰撞，引起材料内部出现大量原子尺度的缺陷的过程，这个过程在很短时间内就会发生。这些缺陷经过长时间的迁移、聚集和复合等形成缺陷团簇、空位等，引起材料微观组织变化，使材料的宏观力学、热学等性能退化。因此，理解辐照缺陷团簇(如空位团簇)的行为对认识辐照损伤微结构演化至关重要。

由于含有空位团簇体系往往具有多个局域亚稳态，当体系从一个基态跃迁到近邻基态时需要历经一系列亚稳态。相对于基态而言，亚稳态由于具有较高的能量，导致空位团簇激发以后很容易返回基态或陷入基态附近由多个亚稳态构成的超势能谷中，从而阻碍体系进一步演化。目前通常是采用势能谷自动填充法(ABC)或格点动力学蒙特卡洛(LKMC)方法访问体系势能面，但是当体系中含有诸如一个亚稳态到另一个亚稳态速率快得多的事件时，这些方法往往因体系“自转”而失效，不能真正地使体系跳出势能谷。

发明内容

有鉴于此，本发明有必要提供一种空位团簇跳出超势能谷的模拟方法，该模拟方法通过记录跃迁过程中的跃迁序列和禁止循环跃迁，可有效访问基态附近的势能面特征，可使得空位团簇跳出含有多个亚稳态的超势能谷，有效促使体系朝着近邻基态跃迁。

为了实现上述目的，本发明采用以下技术方案：

本发明提供了一种空位团簇跳出超势能谷的模拟方法，包括以下步骤：

S1、获取空间团簇周围可能的跃迁路径能垒，建立速率表格；

S2、获取跃迁序列seqtr，所述跃迁序列seqtr中记录有体系内连续跃迁的初末态空位编号，将所述速率表格与所述跃迁序列seqtr进行比对，若存在循环跃迁则修改所述速率表格，禁止循环跃迁；

S3、在修改后的速率表格中选择跃迁事件并执行，其中，每个跃迁事件被选中的几率与其速率成正相关；

S4、对比两次跃迁，更新跃迁序列seqtr；

S5、判断①是否跃迁至另一基态；②质心是否发生运动；

若①②均为是，则进入步骤S6，否则返回步骤S1；

S6、S1-S5迭代结束，获得体系跃迁路径、能量和结构变化。

进一步的，步骤S1中，所述可能的跃迁路径能垒的获得，具体包括如下步骤：

获取空位团簇周围可能的跃迁路径：假设空位团簇的跃迁是空位与空位周围第一近邻和第二近邻通过原子交换而进行的，对于含有n个空位的团簇V

|rvi-r|

其中，disn为第一近邻距离disn1或第二近邻距离disn2，ε的范围在

获取可能的跃迁能垒：根据获得的可能的跃迁路径计算可能的跃迁能垒。

优选的，所述的根据获得的可能的跃迁路径计算可能的跃迁能垒的步骤，具体为：调用NEB方法对每一条可能的跃迁路径计算其跃迁能垒。

进一步的，步骤S1中，所述的建立速率表格的步骤，具体为：假设获得的第i条跃迁路径的跃迁能垒为Esad，则根据式(2)获得速率表格中第i项速率rate

其中，T为温度，v

进一步的，步骤S2的具体步骤为：假设速率表格中第i项跃迁的初末态编号分别为gi、gj，比对gj是否与所述跃迁序列中的某项相同，若相同，则修改速率表格，令该项速率为0。

进一步的，步骤S3中，所述的选择跃迁事件采用蒙特卡洛方法进行选择。

进一步的，步骤S4中，所述更新跃迁序列seqtr的具体方法为：设上一次跃迁为从gi到gj，本次为从gi1到gj1，若gj＝gi1，表明体系沿着路径gi-gj-gj1连续跃迁，则将gj1添加到跃迁序列seqtr中；否则跃迁序列seqtr清空，然后将gi1、gj1添加到跃迁序列seqtr中。

进一步的，步骤S5中，所述的是否跃迁至另一基态的判断方法，具体为：

假设按照步骤S3跃迁后的团簇能量为E

|E

其中，ε的范围为0.05-0.08eV。

进一步的，步骤S5中，所述的质心是否发生运动的判断方法，具体为：

假设同时若程序运行初始团簇构型坐标为rVn0，按照步骤S3跃迁后为rVn1，每个坐标方向求和并除以n，获得质心坐标，进而得到质心位移dMC，若满足式(4)，则质心发生运动

进一步的，步骤S6中，所述体系跃迁路径、能量和结构变化通过迭代中输出的相邻两次跃迁的结构、跃迁能量曲线获得。

与现有技术相比，本发明具有以下有益效果：

本发明通过在通常格点动力学蒙特卡洛的框架中，记录含空位团簇体系跃迁过程中的跃迁序列，禁止基态附近可能重新返回基态的亚稳态之间的低能垒跃迁，避免了重复“自转”的无效跃迁，可使得体系跳出含有多个亚稳态的超势能谷。同现有的模拟方法相比，该模拟方法一方面保持了格点动力学蒙特卡洛在模拟团簇长时间跃迁行为方面的优势，另一方面在LKMC框架中抑制了循环跃迁的事件，可有效促使体系朝着近邻基态跃迁。

附图说明

图1为本发明中一些较佳实施方式中空位团簇跳出超势能谷的主要步骤流程图；

图2为图1中的算法流程图；

图3为含空位团簇体系附近势能面特点以及本发明模拟方法的原理示意图；

图4为应用实施例1中含四个空位的团簇(V

具体实施方式

为了便于理解本发明，下面将结合具体的实施例对本发明进行更全面的描述。但是，本发明可以以许多不同的形式来实现，并不限于本文所描述的实施方式。相反地，提供这些实施方式的目的是使对本发明的公开内容理解的更加透彻全面。

除非另有定义，本文所使用的所有的技术和科学术语与属于本发明的技术领域的技术人员通常理解的含义相同。本文中在本发明的说明书中所使用的术语只是为了描述具体的实施方式的目的，不是旨在于限制本发明。

本发明第一方面提供了一种空位团簇跳出超势能谷的模拟方法，在一些实施方式中，主要步骤如图1中所示的，包括以下步骤：

S1、获取空间团簇周围可能的跃迁路径能垒，建立速率表格；

S2、获取跃迁序列seqtr，所述跃迁序列seqtr中记录有体系内连续跃迁的初末态空位编号，将所述速率表格与所述跃迁序列seqtr进行比对，若存在循环跃迁则修改所述速率表格，禁止循环跃迁；

S3、在修改后的速率表格中选择跃迁事件并执行，其中，每个跃迁事件被选中的几率与其速率成正相关；

S4、对比两次跃迁，更新跃迁序列seqtr；

S5、判断①是否跃迁至另一基态；②质心是否发生运动；

若①②均为是，则进入步骤S6，否则返回步骤S1；

S6、S1-S5迭代结束，获得体系跃迁路径、能量和结构变化。

本发明通过记录含空位团簇体系跃迁过程中的跃迁序列seqtr，避免了重复“自转”的无效跃迁，可使得体系跳出含有多个亚稳态的超势能谷，从而有效促进体系朝着近邻基态跃迁。其中，本文中的所述的“势能面”指的是含空位团簇确定的体系的势能随其核位置改变的图形，“超势能谷”指的是在体系基态附近形成的多个亚稳态。请参阅图3，图3(a)示意的是连接体系两个基态S

在本发明的一些实施方式中，步骤S1中，所述可能的跃迁路径能垒的获得，具体包括如下步骤：

获取空位团簇周围可能的跃迁路径：假设空位团簇的跃迁是空位与周围第一近邻和第二近邻通过原子交换而进行的，对于含有n个空位的团簇V

|rvi-r|

其中，disn为第一近邻距离disn1或第二近邻距离disn2，ε为一个微小量，其选择的原则是使得空位与空位周围原子距离稍微比要考虑的近邻距离大，具体范围可为

获取可能的跃迁能垒：根据获得的可能的跃迁路径计算可能的跃迁能垒。

在获得跃迁路径后，根据跃迁路径计算跃迁能垒，该计算方法没有特别的限定，可以是本领域中常规的计算方法，在本发明的一些具体的实施方式中，所述的根据获得的可能的跃迁路径计算可能的跃迁能垒的步骤，具体为：调用弹性带算法NEB方法对每一条可能的跃迁路径计算其跃迁能垒，由于NEB方法计算跃迁能垒的方法为本领域中的常规方法，因此，这里不再具体阐述。

在本发明的一些实施方式中，步骤S1中，所述的建立速率表格的步骤，具体为：假设获得的第i条跃迁路径的跃迁能垒为Esad，则根据式(2)获得速率表格中第i项速率rate

其中，T为温度，v

在本发明的一些实施方式中，通过跃迁序列记录器记录跃迁序列seqtr，步骤S2的具体步骤为：假设速率表格中第i项跃迁的初末态编号分别为gi、gj，比对gj是否与所述跃迁序列seqtr中的某项相同，若相同，说明该跃迁已经循环至之前某个态了，则修改速率表格，令对应的第i项跃迁的速率为0。

在本发明的一些实施方式中，步骤S3中，选择跃迁事件采用蒙特卡洛方法进行选择，可以理解的是，本发明中，选择跃迁事件的方法还可以采用本领域中任意的常规方法，这里没有特别的限定。跃迁事件被选择后，执行跃迁事件，则团簇中空位运动至相应位置。

进一步的，步骤S4中，所述的更新跃迁序列seqtr的在一些实施方式中，其具体方法为：设上一次跃迁为从gi到gj，本次为从gi1到gj1，若gj＝gi1，表明体系沿着路径gi-gj-gj1连续跃迁，则将gj1添加到跃迁序列seqtr中；否则跃迁序列seqtr清空，然后将gi1、gj1添加到跃迁序列seqtr中。具体的说，跃迁序列seqtr中记录的是一系列跃迁初末态，因此上一次的末态就是下一次的初态，假设跃迁序列seqtr中第一次记录的是A-B，第二次记录的若是B1-C跃迁，而且B1＝B则跃迁序列seqtr中为：A-B-C，若下一次不是从B连续跃迁则A-B清空，跃迁序列seqtr记录为B1-C。

进一步的，在本发明的一些具体实施方式中，步骤S5中，所述的是否跃迁至另一基态的判断方法，具体为：

假设跃迁后的团簇能量为E

|E

其中，ε为一个微小量，其选择与计算体系能量的系统误差有关，在本发明的一些实施方式中，取为0.05-0.08eV，优选的，取为0.05eV。

在本发明的一些实施方式中，步骤S5中，所述的质心是否发生运动的判断方法，具体为：

假设同时若初始团簇构型坐标为rVn0，跃迁后为rVn1，每个坐标方向求和并除以n，获得质心坐标，进而得到质心位移dMC，若满足式(4)，则质心发生运动

进一步的，步骤S6中，所述体系跃迁路径、能量和结构变化通过迭代中输出的相邻两次跃迁的结构、跃迁能量曲线获得。

图2中示出了本发明一实施例中模拟方法的算法流程图，具体流程如下：

S101、开始，进入S102；

S102、读入晶体构型坐标，以及空位团簇初始基态坐标，进入步骤S103；

S103、计算空位团簇周围可能的跃迁路径，并根据可能的跃迁路径计算获得对应的跃迁能垒，建立速率表格，进入步骤S104；

S104、将速率表格与跃迁序列记录器中记录的跃迁序列进行比对，判断空位是否循环跃迁，若是，则进入步骤S105，若否，则进入步骤S106；

S105、已经进入循环跃迁的空位对应速率修改为0，进入步骤S106；

S106、获得更新后的速率表格，进入步骤S107；

S107、选择跃迁事件并执行跃迁事件，进入步骤S108；

S108、更新跃迁序列，进入步骤S109；

S109、判断是否找到另一质心运动的基态，若是，则转入步骤S110，若否，则转入步骤S103。

S110、迭代结束，输出体系跃迁路径、能量和结构变化。

下面结合具体的应用实施例对本发明中空位团簇跳出超势能谷的模拟方法进行更加清楚完整的说明。

应用实施例1

本应用实施例采用本发明一较佳实施方式模拟体心立方金属铁(100)表面上空位团簇从一个基态跃迁到另一个近邻基态的过程，其中，体心立方金属铁(100)的晶格常数a为

具体模拟过程如下：

S1、以满足公式rvi-r

此时，跃迁速率记录器初始化为seqtr:-1-1,gi＝-1,gj＝-1。

S2、在上述速率表格中进行跃迁事件选择，选中行号为2的项，执行后，seqtr＝311295，gi＝311,gj＝295。空位跃迁后，速率表格如下：

经比对后，发现其中行号为6的事件出现了311-295-311循环跃迁，则修改速率表格，禁止该跃迁，修改后的速率表格如下所示：

修改后，上表中0行事件被选中，执行后seqtr＝327 311,gi＝327,gj＝311；

空位跃迁后，速率表格如下：

其中7行事件被禁止，禁止后8行事件被选中执行，同时检测到体系已经跃迁到另一个基态，而且团簇初始坐标为：

最终坐标为：

质心坐标分别为：

MC0＝[11.4220 11.4220 7.1550]

MC1＝[11.4212 8.5659 7.1549]

计算得到质心为dMC＝2.8561，满足跃迁条件，程序停止。

通过以上实施例，可以看出，通过记录含空位团簇体系跃迁过程中的跃迁序列，禁止基态附近可能重新返回基态的亚稳态之间的低能垒跃迁，避免了重复“自转”的无效跃迁，可使得体系跳出含有多个亚稳态的超势能谷，可有效促使体系朝着近邻基态跃迁。

以上所述实施例的各技术特征可以进行任意的组合，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。

以上所述实施例仅表达了本发明的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。因此，本发明专利的保护范围应以所附权利要求为准。