一种纯金属高压原子间相互作用势函数的构建方法

技术领域

本发明涉及一种势函数构建方法，特别涉及一种纯金属高压原子间相互作用势函数的构建方法，适合构造用于高压、高应变率应用环境的原子间相互作用势函数。

背景技术

金属结构材料在高速碰撞、爆炸冲击等作用下的动态响应行为广泛出现在国防和国民经济的众多领域中，例如高速飞行器关键结构部件受到砂石、飞鸟等的碰撞，卫星受到陨石、太空垃圾等的撞击，材料的快速冲压成型，车辆碰撞，常规武器领域中的防弹装甲/穿甲，战略武器中的内爆压缩过程等等。材料的动态响应在宏观上表现为动态本构关系与损伤断裂行为，而这些行为本质上是材料的微结构演化在高温、高压和高应变率等条件下的综合宏观体现。基于经验原子间相互作用势函数的原子模拟方法(如分子动力学、蒙特卡洛等)是目前理解材料在高温、高压和高应变率等极端条件下微结构演化规律最重要的一种理论途径。

由于历史原因，传统的原子间相互作用势函数几乎都是针对常压应用环境来发展的，在应用于高压环境时，经常会出现错误的预测结果，如铁中的虚假高压BCC-FCC相变和冲击塑性缺失、Nb的高压变形孪晶缺失、Pb的高压FCC-HCP相变缺失以及不准确的高压状态方程等。

传统的势函数构造思路是：在已有势模型泛函表达式的基础上，初始化模型参数，通过拟合常压下的材料物性数据，获得泛函表达式中的优化参数，然后将参数化的势模型应用于目标问题(如冲击高压环境)的模拟，检验模拟结果是否与相关实验或第一性原理计算的数据相符。如果不相符，则改变优化过程参数(如各种物性的权重等)，重新进行参数优化、验证，这个过程不断循环迭代，直到模拟结果与参考结果相符为止。但是用于优化过程的物性数据库仅包含较容易计算的物性参数，而用于检验的材料性质或行为通常很难通过简单的方法来计算得到，比如冲击高压的问题中要求势函数能准确描述熔点，而熔点的计算需要很大的计算量，放入优化的数据库中，造成优化流程的效率极其低下，很难用于复杂使用环境的势函数的构造。而冲击高压应用环境十分复杂，要求势函数不仅能准确描述材料在常温、常压下的基本物性，而且还要能准确描述高温、高压下的基本物性。通常需要考虑的典型高温、高压物性包括熔点、高压相变压力阈值、高压相结构的稳定性与状态方程等，尽管很多高压物性可以包含到参考数据库中来进行优化，但是随着参考数据库的急剧扩充，一方面传统构建方法的效率急剧降低，另一方面参考数据库中不同物性参量之间存在潜在关联(耦合)的概率在增加，尤其是高压与常压物性之间耦合，由于数据来源不可避免地存在不同程度的误差，造成这些耦合的物性可能永远都无法同时被准确拟合，最终使优化结果难以收敛。因此，常压与高压物性间的解耦方法对于高效地开发准确的高压势函数至关重要。

发明内容

本发明的目的在于克服现有技术的不足，提供一种纯金属高压原子间相互作用势函数的构建方法。

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

一种纯金属高压原子间相互作用势函数的构建方法，包括以下步骤：

S1、构建纯金属在平衡态下的势函数框架，得到平衡态势函数形式，并确定得到待定的势参数；

S2、根据步骤S1得到的势函数，结合结合能方程、单空位形成能方程、平衡方程、弹性系数方程和对势的光滑连接条件，得到势参数所满足的物理约束条件；

S3、将第一性原理计算或实验测量的结合能、晶格常数和弹性系数代入步骤S2得到的势函数的物理约束条件中，对步骤S1中得到的势参数中的部分势参数进行优化，并将剩余的势参数定义为自由参数，得到初步优化好的平衡态势函数；

S4、对步骤S3得到的初步优化好的平衡态势函数中的自由参数进行优化，并根据优化后的自由参数利用步骤S2得到的物理约束条件对非自由参数进行计算，最终完成全部势参数的确定，得到优化好的平衡态势函数；

S5、对步骤S4得到的优化好的平衡态势函数进行高压矫正，得到最终优化好的高压态势函数。

优选的，步骤S1中，所述平衡态势函数为：

其中，φ(r

优选的，所述电子密度ρ

对势截尾函数：

其中，p＝1/n

优选的，步骤S2中，所述势函数的物理约束条件为：

(1)结合能方程：

其中，ρ

(2)平衡方程

其中，n

(3)单空位形成能方程：

E

其中，N表示系统中总的原子数，ΔF＝∑

(4)弹性系数方程：

其中，

(5)对势插值节点处的光滑连接条件

φ′(r

其中，φ′(r

优选的，步骤S4具体包含以下步骤：

S41、将步骤S3得到的初步优化好的高压态势函数中的自由参数优化的目标函数定义为：

其中，Q

S42、定义始种群规模P＝M+1，种群中第j个个体A

S43、计算B

S44、重新计算新种群中的O

S45、根据优化后的自由参数利用步骤S2得到的物理约束条件对非自由参数进行计算，最终完成全部势参数的确定，得到优化好的平衡态势函数。

优选的，步骤S5具体包含以下步骤：

S51、设定r

其中，r为插值节点位置，M

S52、当r≤r

其中，D

S53、根据步骤S51和步骤S52对步骤S4得到的优化好的平衡态势函数中的与对势相关的参数χ

本发明与现有技术相比，其有益效果在于：

(1)本发明提供的纯金属高压原子间相互作用势函数的构建方法，通过有效地将常压和高压物性在势函数层面上进行了解耦，排除了潜在物性间的耦合对势函数优化过程中收敛性和精度的影响，提高了高压势函数开发的成功率；

(2)本发明通过在误差范围内微调参考数据库中的单空位形成能，很容易实现熔点的正确描述(仅需1～2次验算)，无需通过拟合-验证这种反复迭代的方式对势函数进行繁琐的矫正，极大的缩短了高压势函数的开发周期；

(3)本发明通过将常压与高压物性的拟合过程在势函数的层面上进行了有效解耦，使得采用高压物性拓展的参考数据库优化出的势函数的可移植性比原始势函数更高；

(4)与传统的势函数构建方法相比，本发明提出的势参数的优化流程采用了基于物理约束的多目标遗传进化算法，该算法可以有效地排除非物理的解出现在优化结果中，极大地提高了势函数构建的成功率。

附图说明

图1为本发明实施例中的基于多目标物理约束遗传进化算法的势参数优化流程图；

图2为采用不同数目的对势插值节点优化之后的势函数对比图；

其中，图2(a)为对势函数，图2(b)为嵌入能函数，图2(c)为原子电子密度函数，且不同势函数用额外插值节点数N

图3为势函数预测的广义结合能和对势与实验结果(Ref.)的对比图；

其中，图3(a)为势函数预测的广义结合能与实验结果(Ref.)的对比图；图3(b)为势函数预测的对势与实验结果(Ref.)的对比图；且N

图4为优化的势函数图；

其中，图4(a)为对势；图4(b)为嵌入能；图4(c)为原子电子密度函数；

图5为计算的声子散射关系与实验结果的对比图；

其中，

图6为势函数的广义层错能曲线；

图7为分子动力学结合势函数EAM-II，模拟预测的Al的熔化行为图；

图8为势函数EAM-II预测的Al的固相和液相的自由能-温度曲线；

图9为预测的能量与压力形式的高压状态方程及其与相关参考结果的对比(

具体实施方式

下面将参照附图更详细地描述本公开的示例性实施例。虽然附图中显示了本公开的示例性实施例，然而应当理解，可以以各种形式实现本公开而不应被直立阐述的实施例所限制。相反，提供这些实施例是为了能够更透彻地理解本公开，并且能够将本公开的范围完整的传达给本领域的技术人员。

本发明实施例以常用的计算效率较高的嵌入原子模型(EAM)势为例，来阐述高压势函数的构建方法及其在高压势函数构建过程中的应用。相关原理也可以用于其它类型的原子间相互作用势模型，例如对势和含有键角贡献的响应嵌入原子模型势等。

本发明实施例提供了一种纯金属高压原子间相互作用势函数的构建方法，具体包括以下步骤：

S1、构建纯金属在平衡态下的势函数框架，得到平衡态势函数形式，并确定得到待定的势参数；

平衡态势函数框架为：

其中，φ(r

电子密度ρ

其中，f(r

一般来说，势函数的构建实质就是根据参考物性来确定φ(r

其中，p＝1/n

和相应的对势截尾函数：

其中，r

其中，r

上述方程组可以唯一确定λ

原子电子密度函数一般地取为Ackland型三次样条插值函数：

其中，R

综上所述，在选定插值节点位置的条件下，上述模型中的参数包括n,F

S2、根据步骤S1得到的势函数，结合结合能方程、单空位形成能方程、平衡方程、弹性系数方程和对势的光滑连接条件，得到势函数的物理约束条件；通过物理约束条件来约束势参数，使参数空间的有效维度降低，达到提供优化效率的目的。同时，通过选取平衡态近邻原子间距作为插值节点位置(以下简称平衡节点)，可以有效地将常压与高压物性进行(部分)解耦。例如，平衡态的物性主要与平衡态晶格位置或平衡原子间距相对应的势能及其导数相关，而高压下的物性由其它位置处势能及其导数决定，因此可以在构建出的常压势函数基础上，通过在平衡节点之间增加额外节点的方式来改进高压物性描述的准确性，与此同时，势函数对常压物性的描述精度可以得到完全的保持。此外，这种方法还有一个独特的优点，即，可以在针对常压物性及特殊缺陷能量优化好的势函数基础上进一步其对改善高压状态方程的描述，这个过程不需要重新优化势参数，只需修改近程(小于第一平衡近邻原子间距)对势部分。

通过势函数的具体形式，可以分别推导出结合能、平衡条件、弹性系数以及光滑连接条件的方程，最终得到的相应物理约束条件如下：

(1)结合能方程：

其中，ρ

(2)平衡方程

其中，n

上述方程可进一步表示成关于势参数{y

(3)单空位形成能(E

E

其中，N表示系统中总的原子数，ΔF＝∑

(4)弹性系数(C

其中，

上述方程具体通过以下方法得到：

Ω

其中，Ω

上式中，

其中，将上述前两个方程表示成关于势参数{y

其中，当

当

(5)对势插值节点处的光滑连接条件

φ′(r

对φ(r)的表达式(4)进行求导，分别取r

S3、将第一性原理计算或实验测量的结合能、晶格常数和弹性系数分别代入步骤S2中的方程(10)、(12)和(14)中，然后联立这些方程以及光滑连接条件，总共N+4个方程，对步骤S1中得到的2N个势参数({y

S4、对步骤S3得到的初步优化好的平衡态势函数中的自由参数进行优化，并根据优化后的自由参数利用步骤S2得到的物理约束条件对非自由参数进行计算，完成全部势参数的确定，得到优化好的平衡态势函数。这个过程具体包含以下步骤：

S41、将步骤S3得到的初步优化好的高压态势函数中的自由参数优化的目标函数定义为：

其中，Q

S42、为了方便优化算法的阐述，定义始种群规模P＝M+1，种群中第j个个体A

优化流程的核心思想是通过不断更新换代，使种群的平均目标函数值逐代降低。为此，从第k代到第k+1代的种群演化按如下规则进行：

S43、计算B

(i)如果F

(ii)如果F

(iii)如果F

(iv)如果F

通过上述步骤，可以完成一次种群的更新换代。

S44、重新计算新种群中的O

F(ρ

其中，α为任意实数，在实践中，α的取值依赖于f(r

具体的基因变异算法为：取种群中每个个体与f(r

图1表示出了多目标遗传进化算法与物理约束条件间的交互关系，其中多目标遗传进化算法用于优化自由参数。

S45、根据优化后的自由参数利用步骤S2得到的物理约束条件对非自由参数进行计算，完成全部势参数的确定，得到优化好的平衡态势函数。

S5、对步骤S4得到的优化好的平衡态势函数进行高压矫正，得到最终优化好的高压态势函数；

由于在对势表达式(4)和(5)中仅涉及r≥r

S51、设定一个适当的r

其中，r为插值节点位置，M

S52、当r≤r

其中，D

S53、根据步骤S51和步骤S52对步骤S4得到的优化好的平衡态势函数中的与短程对势相关的参数χ

下面通过应用案例对本发明实施例提供的势函数的构建方法进行阐述

一、应用于Al的高压势函数构建与优化

以一个模型金属Al为例，采用本发明实施例的方法构建其高压势函数，阐述本发明所提出的方法如何高效地实现势函数精度的逐级改进以及解耦方法在其中所发挥的作用。

对于Al的高压势函数，需要考虑的最基本的物性包括常压下的基本物性(见下表1)、高压状态方程以及平衡熔点。首先，确定模型截断距离。这里取第三近邻模型(N＝3)，此时对势函数至少包含4个节点(其中3个对应于平衡的近邻原子间距，最后一个对应于势能截断距离)。其次，考虑额外的插值节点数目和位置。由于较多的模型参数有利于更好地拟合参考数据库，然而模型参数太多会导致优化程序效率极具降低，甚至找不到满意的最优解。以下将采用逐级优化的思路：先采用少的额外插值节点来优化势函数，然后再在前一次优化的基础上逐渐增加额外的插值节点，已到达不断改进势函数精度的目的。在这个过程中，由于每次额外节点的插入都是在上一次优化结果的基础上进行的，因此优化程序很容易收敛并找到更加优化的解。例如，依次采用2、3和4个额外插值节点(节点位置如图2a所示)来构建势函数，原子密度函数的插值节点数都为3，其位置如图2c所示。经过逐级优化之后，每个势函数的函数曲线如图2所示，可以发现随着节点数的不断插入，对势曲线的形状逐级变得复杂。表1中给出了每种势函数对常压物性参考数据库的拟合程度，不难看出，三种势函数均可较好地描述晶格性质，值得注意地是，随着插值节点的不断增加，振动性质的拟合程度越来越好。此外，从计算的广义结合能关系(图3a)中也可以看出，随着插值节点的不断增加，拟合程度在不断增加。这种特点有助于高可移植性势函数的开发。

表1优化之后势函数的计算结果与参考数据的对比表

注：Ref.表示参考数据库中待拟合物性的参考值。

此外，对于针对平衡态物性优化好的势函数，其对高压段的广义结合能关系描述可以得到进一步优化。例如，通过图3可以看出，通过简单调节χ

关于平衡熔点的问题，需要说明地是：经大量数值模拟测试发现，采用本发明所提出的构建方法，再结合至少包含表1中所列的物性的参考数据库，优化出的势函数的熔点几乎总是正确的，究其原因是平衡熔点主要依赖于单空位形成能，而势函数构造过程中已经拟合了该物理性质。

为了对这种方法的效果进行较为全面的检验，接下来在上述优化的势函数基础上，通过拓展参考数据库，使之包含稳定的层错能(γ

/>

其中，上式仅含一个模型参数θ，r

表2.优化的势模型参数

二、结果与验证

采用三种新的势函数预测的常压物性见表3。从中可以发现，EAM-I与EAM-II的结果非常接近，并且都与参考值(Ref.)接近，尤其是它们预测的单空位形成能(E

表3新构建的势函数的预测结果与其它相关结果的对比

a

b

c

d

e

f

g

h

综上所述，本发明实施例提供的纯金属高压原子间相互作用势函数的构建方法，可用于高精度高压势函数的高效构建，通过将常压与高压物性的拟合过程在势函数的层面上进行了有效解耦，使得采用高压物性拓展的参考数据库优化出的势函数的可移植性比原始势函数更高，而且本发明提出的势参数的优化流程采用了基于物理约束的多目标遗传进化算法，该算法可以有效地排除非物理的解出现在优化结果中，极大地提高了势函数构建的成功率。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。