资料来自小木虫论坛;H.J.Zhao汇总
CASTEP获得电荷密度等高线的Matlab作图法
MATLAB程序用来处理MS电荷密度的等高线做法,分四步:
1.到电荷密度图的存储文件,一般是在自己所建的计算project中,隐藏的,castep模块是.charg_frm格式的,可以用txt打开后另存为txt格式。(这个地方注意把文件单独拷贝出来,以免破坏原始文件)
2.读取数据,在matlab命令里输入一下命令:
[path,fn]=uigetfile('*.txt','Open');
fp=fopen([fn,path],'r');
head=fscanf(fp,'%s',4);
data=fscanf(fp,'%f',[4,18*18*144]);
fclose(fp);
这是读取刚才保存的txt数据的,其中18*18*144是可以改的,分别对应a,b,c的重复单元。
3.作矩阵。
a=reshape(data(1,:),[18,18,144]);
b=reshape(data(2,:),[18,18,144]);
c=reshape(data(3,:),[18,18,144]);
1.到电荷密度图的存储文件,一般是在自己所建的计算project中,隐藏的,castep模块是.charg_frm格式的,可以用txt打开后另存为txt格式。(这个地方注意把文件单独拷贝出来,以免破坏原始文件)
2.读取数据,在matlab命令里输入一下命令:
[path,fn]=uigetfile('*.txt','Open');
fp=fopen([fn,path],'r');
head=fscanf(fp,'%s',4);
data=fscanf(fp,'%f',[4,18*18*144]);
fclose(fp);
这是读取刚才保存的txt数据的,其中18*18*144是可以改的,分别对应a,b,c的重复单元。
3.作矩阵。
a=reshape(data(1,:),[18,18,144]);
b=reshape(data(2,:),[18,18,144]);
c=reshape(data(3,:),[18,18,144]);
d=reshape(data(4,:),[18,18,144]);
其中18*18*144是可以改的,分别对应a,b,c的重复单元,即长宽高的值。
4.作等高线图
layer=1;
contour(reshape(a(:,layer,:),[18,144]),reshape(c(:,layer,:),[18,144]),reshape(d(:,layer,:),[18,144]),8);
layer=1(2,3,4....) 可以修改作图的层数,分别是从前面到后面。最后的8表示显示线条数目,可以自己修改成9,10或其它。
hold on 表示累积作图,hold off表示清空前面的 作图
举例;做表面层和第九层的等高线图
>> layer=9;
contour(reshape(a(:,layer,:),[18,144]),reshape(c(:,layer,:),[18,144]),reshape(d(:,layer,:),[18,144]),8);
>> hold on
>> layer=1;
其中18*18*144是可以改的,分别对应a,b,c的重复单元,即长宽高的值。
4.作等高线图
layer=1;
contour(reshape(a(:,layer,:),[18,144]),reshape(c(:,layer,:),[18,144]),reshape(d(:,layer,:),[18,144]),8);
layer=1(2,3,4....) 可以修改作图的层数,分别是从前面到后面。最后的8表示显示线条数目,可以自己修改成9,10或其它。
hold on 表示累积作图,hold off表示清空前面的 作图
举例;做表面层和第九层的等高线图
>> layer=9;
contour(reshape(a(:,layer,:),[18,144]),reshape(c(:,layer,:),[18,144]),reshape(d(:,layer,:),[18,144]),8);
>> hold on
>> layer=1;
contour(reshape(a(:,layer,:),[18,144]),reshape(c(:,layer,:),[18,144]),reshape(d(:,layer,:),[18,144]),8);
这是沿侧面做等高线(即xz平面),其它方向作图类似。只要修改其中的参数就可以了,如把a(:,layer,:),[18,144]),reshape(c(:,layer,:),[18,144]),reshape(d(:,layer,:),[18,144]),8);
改为a(layer,:,:),[18,144]),reshape(c(layer,l:,:),[18,144]),reshape(d(layer,:,:),[18,144]),8);
这是沿侧面做等高线(即xz平面),其它方向作图类似。只要修改其中的参数就可以了,如把a(:,layer,:),[18,144]),reshape(c(:,layer,:),[18,144]),reshape(d(:,layer,:),[18,144]),8);
改为a(layer,:,:),[18,144]),reshape(c(layer,l:,:),[18,144]),reshape(d(layer,:,:),[18,144]),8);
MS 4.3磁性体系LDA+U计算设置
在最新发布的MS版本中增加了LDA+U的计算,这是MS计算性能特别是CASTEP模块的一
次巨大提高,终于解决了已久的d,f轨道带隙分裂问题,也不会再被审稿人轻视认为是简单DFT计算结果了。
这里讲解一下MS.43里面element设置方面的问题,我们选择了Fe3O4体系,该体系是铁磁性的,计算该体系主要涉及到两个参数设置,一个是自旋设置(SPin parameters),其次是Fe Hubbard参数设置。
首先介绍需要用到的菜单:
打开MS之后看到主菜单,打开如下图所示的菜单:(点击此图,可放大观看)
这里讲解一下MS.43里面element设置方面的问题,我们选择了Fe3O4体系,该体系是铁磁性的,计算该体系主要涉及到两个参数设置,一个是自旋设置(SPin parameters),其次是Fe Hubbard参数设置。
首先介绍需要用到的菜单:
打开MS之后看到主菜单,打开如下图所示的菜单:(点击此图,可放大观看)
由于计算体系位Fe3O4,计算采用USPP完成,Fe的USPP价电子参数为3d64s2,O是2s22p6,结构中的磁性主要来自于Fe原子磁矩的定向排列,计算中自旋参数和Hubbard参数主要是针对Fe调整。再建造好Fe3O4结构之后(注MS结构库中包含此结构),为一个面心立方晶格,一个惯用晶胞包含8个Fe3O4单位,因此整个结构有24个Fe原子,32个O。计算可以采用FCC结构的初级晶胞结构,也可以采用惯用晶胞。我们计算采用后者。
打开Modify-electronic configuration菜单,看到Spin和Hubbard U两个子菜单
如图:
SPin主要是自旋控制参数,由于再铁磁性结构中Fe自旋是高度极化的,因此这个地方设置为High spin state,即3d6结构为4个自旋平行同向电子。自旋方向不管是Up还是Down当然不会改变Fe原子自旋极化的性质,但会改变总的磁矩(自旋电子数目,这点会再CASTEP计算参数Spin设置中体现出来,可以尝试把部分Fe原子自旋改成High,SPin up,部分为High,Downmatlab等高线填充颜,然后查看CASTEP计算控制面板的SPin数值)。
打开Modify-electronic configuration菜单,看到Spin和Hubbard U两个子菜单
如图:
SPin主要是自旋控制参数,由于再铁磁性结构中Fe自旋是高度极化的,因此这个地方设置为High spin state,即3d6结构为4个自旋平行同向电子。自旋方向不管是Up还是Down当然不会改变Fe原子自旋极化的性质,但会改变总的磁矩(自旋电子数目,这点会再CASTEP计算参数Spin设置中体现出来,可以尝试把部分Fe原子自旋改成High,SPin up,部分为High,Downmatlab等高线填充颜,然后查看CASTEP计算控制面板的SPin数值)。
下面按照Fe3O4具体情况设置,Fe自旋为High spin, 方向Up。自旋电子数目+4,O原子不予考虑。HUbbard带参数设置,主要针对d轨道,Fe的Hubbard参数再CASTEP里面默认数值是2.5eV,根据一些文献也可以选取其他参数,比如再FeO计算文献中选用了6eV,也有用4-5eV的,根据具体计算情况和结果可以调整,这可能会要求计算几次,比较结果。但为了保证结果准确性,是必要的。我们计算选Fe的d U=5eV。如下所示:
元素的电荷这里不需要考虑,实际上一般QM计算原子电荷状态是中性原子状态,再MD计算
元素的电荷这里不需要考虑,实际上一般QM计算原子电荷状态是中性原子状态,再MD计算
中,可能对这个参数有要求,或者计算原子必须要带电荷。一般默认即可!
版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。
发表评论