1. RMC模拟与其初始构型.cfg文件
逆蒙特卡洛(RMC)模拟是根据一维衍射数据(X-ray, Neutron)反推出三维原子排列的有效手段,该算法的实现可以借助Orsolya Gereben等人开发的RMC_POT程序包。RMC_POT做RMC模拟所需的文件包括:主程序文件(.dat),初始构型文件(.cfg),实验数据文件(格式用户自己定义),sfactorcube文件(控制原子移动最大位移)。在这其中.cfg初始构型文件是极其重要的。以反推金属玻璃(Metallic glasses,MG)中原子排布为例,总的来说.cfg文件的创建需要满足以下几个原则:i) 构型内原子的数密度需要符合实验测定的材料密度;ii) 原子无规堆积,缺少长程有序结构; iii) 原子两两之间的最小间隔满足定义的cutoff。如何事半功倍的构建出一组.cfg初始构型文件对进行后续独立的RMC模拟实验是直观重要的。
.cfg文件的格式(以四元MG PdNiCuP为例)如下所示:
(Version 3 format configuration file) !file created by SimpleCfg::save !
MG_Pd43Ni20Cu27P10 metallic glass@RMC_POT
0 0 0 0 moves generated, tried, accepted, potential-accepted
0 configurations saved
5400 molecules of all types
4 types of molecule
1 is the largest number of atoms in a molecule
0 Euler angles are provided
F (box is cubic)
Defining vectors are:
20.815651 0.000000 0.000000
0.000000 20.815651 0.000000
0.000000 0.000000 20.8156512322 molecules of type 1
, 1 atomic sites
0.000000 0.000000 0.0000001080 molecules of type 2
, 1 atomic sites
0.000000 0.000000 0.0000001458 molecules of type 3
, 1 atomic sites
0.000000 0.000000 0.000000540 molecules of type 4
, 1 atomic sites
0.000000 0.000000 0.000000-0.868839994050115 -0.991980000410556 -0.436619997026109
-0.808699995279919 -0.878679998204117 -0.613059997554325
-0.760199993841675 -0.952720001332732 -0.487940013432776……
每个原子的初始位置以(-1,1)fraction坐标表示,构型满足周期性边界条件约束。
2. AIMD模拟和XDATCAR文件
AIMD模拟不依赖于势函数的限制,因而其适合处理多组元MG的模拟问题。但是AIMD建模MG存在理论及技术上的限制:a) 从玻璃弛豫现象上来讲,金属玻璃在低温状态处于绝对的亚稳定过程,弛豫时间近似无限长,理论上任何的计算机模拟/实验研究都不能遍历如此缓慢的动力学过程;b) AIMD模拟基于求解薛定谔方程,其运算速度慢,所需计算核心多,很难完整勾勒出MG实际形成过程的结构重组行为,因而基于melting-quenching创建的MG构型严格而言与实验得到的结构存在不小的差异。基于以上所述,基于AIMD输出的构型,利用衍射实验的数据做RMC精修,是捕捉实际MG玻璃结构的有效手段。
XDATCAR 文件是VASP进行AIMD计算输出的每模拟步的原子坐标集合。AIMD输出的每步构型中原子排列满足.cfg文件构建的三个重要原则,但是我们并不能直接将XDATCAR中的构型拿来进行RMC模拟,最关键的原因在于AIMD模拟构型非常小,一般<300原子数,而RMC模拟构型原子数在1000~10000之间,二者存在数量级差异。如何对AIMD的输出构型改造进行RMC模拟是个有趣的问题。
XDATCAR文件的格式(以四元MG PdNiCuP为例)如下所示:
Random MG from XJ Han
1.00000000000000
13.8765888450000006 0.0000000000000000 0.0000000000000000
0.0000000000000000 13.8765888450000006 0.0000000000000000
0.0000000000000000 0.0000000000000000 13.8765888450000006
Pd Ni Cu P
86 40 54 20
Direct
0.1967512433709550 0.0120249114084959 0.8450697621079324
0.2869488440860818 0.1819650879617962 0.5803966400872402
0.3597098725726727 0.0709252575579089 0.7680881170493263……
每个原子的初始位置以(0,1)fraction坐标表示,构型满足周期性边界条件约束。
3. 技术路线
<1>首先取出XDATCAR中任一一步原子构型,将其用VESTA打开,另存为.cif格式文件
<2>利用materials studio打开上述.cif文件,利用build-supercell命令将其扩充(一般选择3x3x3)
<3>将上述超胞保存为.cif文件,再次利用VESTA打开,将其另存为POSCAR格式(其与上述XDATCAR格式完全一样,但是只包含一步的原子坐标),保存时选择笛卡尔坐标
<4>利用以下的MATLAB程序段将其转变为.cfg格式
% This code was programmed in order to transform the POSCAR data into .cfg files for RMC simulation.
% created by Ge Xuan at 2020/03/09
% Copyright interpretation belongs to Dr. Ge
tic
format long
clear;
close all;
clc;
%以下循环是将模拟过程中XDATCAR构型数据依次输进程序中
fp=fopen('filepath\POSCAR_3x3x3','rt');
% c=fscanf(fp,'%c');
%以字符串的形式将带题头的全部数据输入到XDATCAR矩阵中
[b1,b2,b3]=textread('filepath\POSCAR_3x3x3','%s %s %s');
XDATCAR=[b1,b2,b3];
fclose(fp);
num_atoms=input('please give the number of atoms in youe simulation box: \n');
num_steps=size(XDATCAR,1)/(num_atoms+10) % 计算出模拟步数;
if num_steps~=fix(num_steps)
fprintf('you gave a wrong input')
return;
end
for t=1:num_steps %%以下的循环是按照RMC_cfg的格式输出
dir = strcat('E:\a. 我的论文\AlNiCo熔体结构数据\周焕乙\MG_PdNiCuP_', num2str(t), '.cfg');
fid=fopen(dir,'wt');
%在文件第一行输出固定格式字符串
fprintf(fid,' %s\n','(Version 3 format configuration file) !file created by SimpleCfg::save !');
%在文件第二行输出固定格式字符串和模拟步数
fprintf(fid,'%s\n','MG_Pd43Ni20Cu27P10 metallic glass@RMC_POT'); %this sentence should be identified every times
%第三/四行输出空白行
fprintf(fid,' \n\n');
%第五行输出固定的表达式
fprintf(fid,' %d %d %d %d %s\n',0,0,0,0,'moves generated, tried, accepted, potential-accepted');
%第六行输出固定的表达式
fprintf(fid,' %d %s\n \n',0,'configurations saved');
%第七行输出固定的表达式
fprintf(fid,' %d %s\n',num_atoms,'molecules of all types');
%第八-第十一行输出固定的表达式
fprintf(fid,' %d %s\n %d %s\n %d %s\n \n',4,'types of molecule',1,'is the largest number of atoms in a molecule',0,'Euler angles are provided');
%第十二行到第十七行输出模拟盒子的大小
fprintf(fid,' %s\n %s \n %8.6f %8.6f %8.6f\n %8.6f %8.6f %8.6f\n %8.6f %8.6f %8.6f\n\n',...
'F (box is cubic)','Defining vectors are:',...
0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+3,1)),0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+3,2)),0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+3,3)),...
0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+4,1)),0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+4,2)),0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+4,3)),...
0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+5,1)),0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+5,2)),0.5*str2double(XDATCAR((t-1)*(num_atoms+10)+5,3)));
%第十八行到第二十一行输出第一类原子的信息,这个每次都要修改一下
fprintf(fid,' %d %s\n, %d %s \n %7.6f %7.6f %7.6f\n\n',2322,'molecules of type 1',1,'atomic sites',0,0,0);
%第二十二行到第二十五行输出第二类原子的信息,这个每次都要修改一下
fprintf(fid,' %d %s\n, %d %s \n %7.6f %7.6f %7.6f\n\n',1080,'molecules of type 2',1,'atomic sites',0,0,0);
%第二十五行到第二十八行输出第三类原子的信息,这个每次都要修改一下
fprintf(fid,' %d %s\n, %d %s \n %7.6f %7.6f %7.6f\n\n',1458,'molecules of type 3',1,'atomic sites',0,0,0);
%第二十八行到第三十一行输出第四类原子的信息,这个每次都要修改一下
fprintf(fid,' %d %s\n, %d %s \n %7.6f %7.6f %7.6f\n\n',540,'molecules of type 4',1,'atomic sites',0,0,0);
%提取盒子大小
x_size = str2double(XDATCAR((t-1)*(num_atoms+10)+3,1));
y_size = str2double(XDATCAR((t-1)*(num_atoms+10)+4,2));
z_size = str2double(XDATCAR((t-1)*(num_atoms+10)+5,3));
for j=(t-1)*(num_atoms+10)+11:1:t*(num_atoms+10)
fprintf(fid,' %16.15f %16.15f %16.15f\n',str2double(XDATCAR(j,1))/x_size*2-1,str2double(XDATCAR(j,2))/y_size*2-1,str2double(XDATCAR(j,3))/z_size*2-1);
end
fclose(fid);
end
toc
上述代码中的空格和换行是按照.cfg文件格式设置的,不建议随便删改;程序运行时间在4~5s之间(原子数5400个),因而上述代码可以迅速根据XDATACR文件创造出一系列满足要求的RMC初始构型文件。
上述代码如有使用不变,请send Email to:gejade9254@163.com。谢谢