本发明属于材料相变的理论计算领域，特别是涉及一种计算镧系化合物相稳定性的方法。
背景技术：
1、镧系过渡金属三元化合物lnyz是极为重要的功能材料大家族，其中ln代表镧系元素，y为过渡元素，z对应主族元素。其因在磁性、热电、拓扑绝缘体以及超导等诸多方面展现出的特性，而受到了广泛关注。由于化学元素组成的多样性，lnyz化合物能够构建出多种不同的晶体结构，包括高对称性的立方half-heusler相以及其他对称性较低的结构，例如六方相或正交相等。这些三元化合物所呈现出的丰富多样的晶体结构与性质，为材料科学的持续探索及利用提供了广阔的前景。
2、当材料因温度或压强等因素而引发相变并形成新相时，通常会伴有物理性质的改变，因此对材料的相稳定性及相变温度展开研究，对挖掘其应用潜力意义重大。在针对镧系化合物晶体结构的实验研究中，发现了比较特殊的现象，即在仅有镧系元素ln存在差异而其他组成元素相同的条件下，这一系列化合物并非仅形成单一结构，而是会伴随ln3+离子半径的变化而同时存在多种结构。例如，hartje等人对lnnisb化合物的晶体结构进行了系统研究，实验发现含有轻稀土元素(ln＝la-nd、sm)的lnnisb为六方相(空间群：p63/mmc)，而含有重稀土元素(ln＝gd-lu)的lnnisb则为立方相(空间群：)。此外，在诸如lnnio3、ln3tao7、ln2(cn2)3和ln2o3等其他镧系化合物中，也观察到材料会随ln3+半径的减小而发生结构相变。
3、借助理论手段对材料的热力学性质进行计算，有助于预测其相稳定性及相变温度，从而更好地解释相变机制。常规的热力学性质计算方法包含以下步骤：首先基于密度泛函理论通过有限位移法或密度泛函微扰理论计算材料在0k下的声子谱；其次，在声子谱收敛的前提下，基于准简谐理论(quasi-harmonic approximation,qha)考虑晶格振动的贡献后进一步计算体系不同相的热容、振动熵、自由能等热力学性质。然而，该方法对于那些在0k下声子谱存在较大虚频的体系并不适用，如图1所示。
4、由此，现在亟需一种适用范围更为广泛的热力学性质计算方法，以满足理论预测材料的相稳定性的需求。
技术实现思路
1、为解决现有技术问题，本发明旨在提供一种计算镧系化合物相稳定性的方法。通过对第一性原理分子动力学的轨迹及原子受力信息进行处理，获取各相在有限温下的力常数，从而能够预测材料不同结构的热力学性质。此方法充分考虑了晶格振动与温度效应的贡献，且在此基础上还可对材料的晶格热导率进行预测。
2、为实现上述目的，本发明采取的技术方案如下：
3、本发明为一种计算镧系化合物相稳定性的方法，包括如下步骤：
4、s1、对镧系化合物的立方相和六方相的晶体原胞进行结构优化，得到两相各自对应的平衡原胞；
5、s2、对两相的平衡原胞分别进行自洽计算，得到两相的形成焓差值和两相的振动自由能差值；
6、s3，将振动自由能差值与形成焓差值相加，获得两相的吉布斯自由能差值与温度的关系，用于预测两相的相稳定性及相变温度。
7、优选地，s1中，对镧系化合物的立方相和六方相的晶体原胞，分别采用第一性原理计算进行结构优化。
8、优选地，s1中，结构优化的参数设置为：采用投影缀加波方法，以广义梯度近似的gga-pbe作为交换相关函数，平面波截断能均设为400ev，结构优化的能量和受力的收敛精度分别设定为10-6ev和高对称k点设置为30/a+1，30/b+1和30/c+1，其中a、b和c分别对应晶体原胞的晶格参数。
9、优选地，s2包括如下子步骤：
10、s2.1，对两相的平衡原胞分别进行自洽计算，得到两相的形成焓差值
11、s2.2，对平衡原胞进行扩胞，然后进行第一性原理分子动力学计算，得到构型和受力信息，接着计算两相在有限温度下的二阶力常数，最后计算两相的声子谱；
12、s2.3，测试声子谱是否收敛，若声子谱收敛，则进入s2.4；若声子谱不收敛，则更改参数重复s2.2，直至获得稳定的声子谱再进入s2.4；
13、s2.4，基于收敛的声子谱对应的二阶力常数，得到两相的振动自由能差值。
14、优选地，s2.1中，对平衡原胞进行自洽计算时，高对称k点设置为60/a+1，60/b+1，60/c+1。
15、优选地，s2.2中，所述立方相和所述六方相超胞的扩胞倍数分别为6×6×3和6×3×3。
16、优选地，s2.2中，基于aimd的结果获得有限温力常数以计算材料的晶格振动相关的热力学性质，所述aimd计算时，选取了aimd nvt系综，温度设置为300k，模拟步数为10000步，时间步长为2fs。
17、优选地，s2.2中，利用phonopy软件处理所述s3的二阶力常数得到立方相和六方相在300k时的有限温声子谱。
18、优选地，s2.4中，基于收敛的声子谱对应的二阶力常数，使用phonopy软件计算体系在0～2000k温度区间的振动熵差值和振动焓差值得到lnnisb两相的振动自由能差值。
19、优选地，s3中，吉布斯自由能差值的计算公式为：
20、
21、其中，t为温度，表示六方相与立方相的振动焓的差值，代表两相的振动熵之差。
22、与现有技术相比，本申请的技术方案具有以下有益技术效果：
23、(1)本发明通过第一性原理分子动力学模拟得到材料的有限温力常数，可计算包含晶格振动与温度效应贡献的热力学性质，达到对材料的不同结构相稳定性的预测目的。
24、(2)本发明克服部分特定相在0k下声子谱存在虚频，导致无法计算热力学性质的困难，拓宽了适用范围。
1.一种计算镧系化合物相稳定性的方法，其特征在于，包括如下步骤：
2.根据权利要求1所述的计算镧系化合物相稳定性的方法，其特征在于，所述s1中，对镧系化合物的立方相和六方相的晶体原胞，分别采用第一性原理计算进行结构优化。
3.根据权利要求2所述的计算镧系化合物相稳定性的方法，其特征在于，所述s1中，所述结构优化的参数设置为：采用投影缀加波方法，以广义梯度近似的gga-pbe作为交换相关函数，平面波截断能均设为400ev，结构优化的能量和受力的收敛精度分别设定为10-6ev和高对称k点设置为30/a+1，30/b+1和30/c+1，其中a、b和c分别对应所述晶体原胞的晶格参数。
4.根据权利要求2所述的计算镧系化合物相稳定性的方法，其特征在于，所述s2包括如下子步骤：
5.根据权利要求4所述的计算镧系化合物相稳定性的方法，其特征在于，所述s2.1中，对所述平衡原胞进行自洽计算时，高对称k点设置为60/a+1，60/b+1，60/c+1。
6.根据权利要求5所述的计算镧系化合物相稳定性的方法，其特征在于，所述s2.2中，所述立方相和所述六方相超胞的扩胞倍数分别为6×6×3和6×3×3。
7.根据权利要求6所述的计算镧系化合物相稳定性的方法，其特征在于，所述s2.2中，基于aimd的结果获得有限温力常数以计算材料的晶格振动相关的热力学性质，所述aimd计算时，选取了aimd nvt系综，温度设置为300k，模拟步数为10000步，时间步长为2fs。
8.根据权利要求7所述的计算镧系化合物相稳定性的方法，其特征在于，所述s2.2中，利用phonopy软件处理所述s3的二阶力常数得到立方相和六方相在300k时的有限温声子谱。
9.根据权利要求8所述的计算镧系化合物相稳定性的方法，其特征在于，所述s2.4中，基于收敛的声子谱对应的二阶力常数，使用phonopy软件计算体系在0～2000k温度区间的振动熵差值和振动焓差值得到lnnisb两相的振动自由能差值。
10.根据权利要求9所述的计算镧系化合物相稳定性的方法，其特征在于，所述s3中，所述吉布斯自由能差值的计算公式为：
