软件
产品
在之前视频中，我们介绍了用VESTA操作显示差分电荷的内容。有朋友在评论区留言说，如果沿着一个方向去求和，怎么把电荷密度累加起来？这里给大家提供一个MATLAB的程序（讲解视频在这里）。
当我们在INCAR里边设置LCHARG=T以后，一般就会看到两个跟电荷有关的文件，CHG和CHGCAR。后者的有效数字多一些，每行是5列，我们一般是处理这个文件。以下面分子为例，我们用VESTA做完差分的（也可以根据今天的程序自己做差分）。这里重点讲解Matlab程序的导入。可以看到，文件的前半部分是结构信息，对应POSCAR部分，这里4个碳6个氢，有晶胞参数，每个原子的具体坐标。这里要注意的是196，168，140这三个整数，对应电荷密度沿不同方向的划分，下面的数据就是每一个网格点上的电荷密度啊，这里有正有负的话呢，是因为做了差分。我们需要用程序把这些数据导成一个196×168×140的矩阵，然后我们就可以沿着任意一个方向去进行累加。这里要注意我们选的晶胞是长方体，就才可以用现在的程序（附录A）。如果比较复杂的晶胞的话，有通用的程序（@五班最强路人）。
文件举例
首先，我们要把文件前面的结构信息和三个整数那一行全删掉，保证下面是一个矩阵，这里用load的指令导入数据。这里a是921984行，5列，累计有4609920个数据，就是196x168x140，这里VASP就是按一定顺序把数据写出来的。这里的for循环就是把这些数据变成1行，4609920列，然后继续转化。这里nn=5是因为CHGCAR文件每行是5列数据。
数据读入
在程序中，不同文件，每次大家修改n2,n1,n3三个整数即可。这里有三个for循环，就是 把对应的数据放到矩阵all的相应的位置，最核心的就是这里，这里主要是根据VASP的写文件顺序确定的，z方向呢，就是140个分割点。
数据转化
生成完all矩阵后，就可以取其中的一层来画图，也可以把它加起来。变量tp的初值是第一层，从第二层一直到N3，全部把它累加，那就把电荷全加起来了，是一个二维矩阵，行列分别是168×19。可以用contour做等高线，也可以用image函数做彩图。这里用了flipud函数是为了显示和分子结构的方向一致。注意，matlab不同版本的colormap是不一样的，旧的版本image在1~60区分度很高，但其实比较新的版本是1~200多才有区分度。大家可以看到这些蓝的地方就是失电荷比较多，红的地方得电荷比较明显。在二维材料、表面等体系的电荷转移，这个显示效果就会好一点。
数据显示
附录：
clear
a=load('CHGDIFF.vasp');
sa=size(a);
b=zeros(1,sa(1)*sa(2));
nn=5;
for ii=1:length(a)
b(1,nn*(ii-1)+1:nn*ii)=a(ii,:);
end
% 196 168 140
n2=196；n1=168；n3=140；
all=zeros(n1,n2,n3);
for ii=1:n3
for jj=1:n1
for kk=1:n2
all(jj,kk,ii)=b(1,n1*n2*(ii-1)+n2*(jj-1)+kk);
end
end
end
tp=all(:,:,1);
for ii=2:n3
tp=tp+all(:,:,ii);
end
contour(tp)
figure
tp=flipud(tp);
image((tp-min(min(tp)))*200/(max(max(tp))-min(min(tp))))
colormap(jet)
免责声明：本文系网络转载或改编，未找到原创作者，版权归原作者所有。如涉及版权，请联系删
技术文档
推荐好文
155-2731-8020
