2.3.2 静态自洽计算的输出文件#
在这里,你需要跟着教程初步了解VASP的主要输出文件及其含义。
提交计算任务之后,如果任务正常开始计算时,我们会看到当前路径下多了很多文件,用ls命令查看
(base) [zpwu@mu003 Si]$ ls
2617356.err DOSCAR KPOINTS POTCAR vasprun.xml
2617356.out EIGENVAL OSZICAR REPORT WAVECAR
CHG hostfile OUTCAR result XDATCAR
CHGCAR IBZKPT PCDAT SYMMETRY
CONTCAR INCAR POSCAR vasp641.slurm在这里,我们重点需要关注的是result、OUTCAR、OSZICAR这三个文件,除了第一个文件名字取决于脚本里的写法,其余文件名字都是VASP程序固定的命名。
(1)result文件
Result文件最上面几行会给出此次计算调用的计算资源(核数)、vasp版本等。之后会进入自洽迭代步骤(俗称,电子步)。请注意,result文件电子步前面部分可能有warning,一般可以不管,感兴趣可用自行去阅读这些warning内容。但是如果最后出现了warning,通常意味着计算可能有问题。这一点我们会在第三章详细说明。
(2)OUTCAR文件
VASP的全部计算信息都会在这个文件显示。输入cat OUTCAR并回车,然后借助鼠标滚轮上滑从第一行开始看。OUTCAR中前面依次会有vasp的版本信息,计算资源(核数),INCAR内容,POTCAR主要内容(这里不做展示)。
继续往下,看到了结构信息
Lattice vectors:
A1 5.4687280655, 0.0000000000, 0.0000000000)
A2 = 0.0000000000, 5.4687280655 0.0000000000)
A3 = 0.0000000000, 0.0000000000, 5.4687280655)
Analysis of symmetry for initial positions (statically):
Subroutine PRIcEL returns following result:
LATTYP: Found a face centered cubic cell.
ALAT 5.4687280655
Lattice vectors:
A1 = 2.7343640327, -2.7343640327, 0.0000000000)
A2 = 2.7343640327, 0.0000000000, 2.7343640327)
A3 = 0.0000000000, -2.7343640327, 2.7343640327)
4 primitive cells build up your supercell.注意,这里vasp默认对结构寻找了对称性,因此找到了一个更小单元的胞进行计算,这个小胞里只有两个原子。
Primitive cell
volume of cell : 40.8883
direct lattice vectors reciprocal lattice vectors
2.734364033-2.734364033 0.000000000 0.182857876 -0.182857876 -0.182857876
2.7343640330.000000000 2.734364033 0.1828578760.182857876 0.182857876
0.000000000 -2.734364033 2.734364033 -0.182857876-0.182857876 0.182857876
length of vectors
3.866974700 3.866974700 3.866974700 0.316719132 0.316719132 0.316719132
position of ions in fractional coordinates (direct lattice)
0.750000000 0.750000000 0.750000000
0.000000000 0.000000000 0.000000000以及K点的信息
Subroutine IBzKPT returns following result:
Found 4 irreducible k-points:
Following reciprocal coordinates:
Coordinates Weight
0.125000 0.125000 0.125000 8.000000
0.375000 0.125000 0.125000 24.000000
0.375000 0.375000 0.125000 24.000000
0.375000 0.375000 0.375000 8.000000
Following cartesian coordinates:
Coordinates Weight
0.022857 0.022857 0.022857 8.000000
0.068572 0.022857 0.022857 24.000000
0.068572 0.068572 0.022857 24.000000
0.068572 0.068572 0.068572 8.000000回想我们设置的KPOINTS文件,取的k点应该有64个。由于体系具有对称性,很多k点是等价的,实际只需要计算那些不可约的点,这里只有4个。在IBZKPT文件中也可以查看到这些不可约k点。当我们把每个k点的权重Weight相加,就得到了64个k点。在k点信息之后,还会给出完整的INCAR参数,因此实际计算中INCAR取的参数值也可以在OUTCAR这里查看。
SYSTEM = unknown system
POSCAR = Si8
Startparameter for this run:
NWRITE = 2 write-flag & timer
PREC = normal normal or accurate (medium, high low for compatibility)
ISTART job :0-new 1-cont 2-samecut
ICHARG 2 charge: 1-file 2-atom 10-const
ISPIN 1 spin polarized calculation?
LNONCOLLINEAR = F non collinear calculations
LSORBIT = F spin-orbit coupling
INIWAV: 1 electr: 0-lowe 1-rand 2-diag
LASPH F aspherical Exc in radial PAw
Electronic Relaxation 1
ENCUT = 400.0 eV 29.40 Ry 5.42 a.u. 8.92 8.92 8.92*2*pi/ulx,y,z
ENINI = 400.0 initial cutoff
ENAUG = 322.1 eV augmentation charge cutoff
NELM 60; NELMIN= 2; NELMDL= -5 # of ELM steps继续下滑,我们会看到这部分内容,左边就是4个不可约k点的坐标,最右边是每个k点对应的波函数所需要的平面波(可以理解为需要这么多平面波将该电子态的波函数做展开,这些平面波数量直接取决于INCAR中设置的ENCUT参数,参数越大,考虑的平面波越多,结果也越精确,代价就是贵(对于二维材料需要在z方向加真空层,真空层越大,所需要的平面波越多,后面我们会进一步了解)。
k-point 1 0.1250 0.1250 0.1250 plane waves: 2969
k-point 2 0.3750 0.1250 0.1250 plane waves: 2962
k-point 3: 0.3750 0.3750 0.1250 plane waves: 2970
k-point 4: 0.3750 0.3750 0.3750 plane waves: 2983
maximum and minimum number of plane-waves per node 2983 2962
maximum number of plane-waves: 2983
maximum index in each direction:
IXMAX= 8 IYMAX= 8 IZMAX= 8
IXMIN= 6- IYMIN= 6- IZMIN= 6-
serial 3DFFT for wavefunctions
parallel 3D FFT for charge:
minimum data exchange during FFTs selected (reduces bandwidth)继续下滑,我们可以看到:
Iteration 1( 1)
POTLOK: cpu time 0.0075: real time 0.0130
SETDIJ: cpu time 0.0011: real time 0.0011
EDDAV: cpu time 0.0504: real time 0.1203
DOS: cpu time 0.0007: real time 0.0007
LOOP: cpu time 0.0597: real time 0.1352这是程序开始进行自洽迭代计算了,括号里的(1)意思是第一次循环,叫做第一个电子步。每次电子步算完,都会输出能量信息
Free energy of the ion-electron system (eV)
alpha Z PSCENC 13.17272099
Ewald energy TEWEN -907.77962997
-Hartree energ DENC = -45.29355852
-exchange EXHF 一 0.00000000
-v(xc)+e(xc) XCENC -74.74068619
PAW double counting 704.99039692 -564.24414232
entropy T*S EENTRO = -0.00205010
eigenvalues EBANDS 73.94756885
atomic energy EATOM = 824.52596983
Solvation Ediel_sol 0.00000000
free energy TOTEN 24.57658948 eV
energy without entropy = 24.57863958 energy(sigma->0) : 24.57761453经历多次电子步迭代,当前后两次能量差达到我们INCAR中设置的收敛阈值EDIFF时,程序停止迭代,开始计算电子结构和原子受力等信息。
比如这里输出了费米能级E-fermi和每个k点在每个态band(能带)上的电子占据occupation情况。
E-fermi : 5.6033 XC(G=0): -9.3267 alpha+bet :-11.7317
Fermi energy: 5.6032749100
k-point 1: 0.1250 0.1250 0.1250
band No. band energies occupation
1 -5.9979 2.00000
2 -2.9678 2.00000
3 -2.9678 2.00000
4 -2.9678 2.00000
5 -1.1620 2.00000
6 -1.1620 2.00000
7 -1.1620 2.00000
8 2.2585 2.00000
9 2.2585 2.00000继续下滑,会看到计算达到收敛相关的信息输出
aborting loop because EDIFF is reached
CHARGE: cpu time 0.0048: real time 0.0048
FORLOC: cpu time 0.0003: real time 0.0003
FORNL cpu time 0.0117: real time 0.0117
STRESS: cpu time 0.0319: real time 0.0320
FORCOR: cpu time 0.0029: real time 0.0029
FORHAR: cpu time 0.0011: real time 0.0011
MIXING: cpu time 0.0003: real time 0.0003
OFIELD: cpu time 0.0000: real time 0.0000随后开始计算压强和受力
FORcE on cell =-STRESS in cart. coord. units (eV):
Direction XX ZZ XY YZ ZX
Alpha Z 13.17272 13.17272 13.17272
Ewald -302.59373 -302.59373 -302.59373 -0.00000 -0.00000 -0.00000
Hartree 20.23882 20.23882 20.23882 -0.00000 -0.00000 -0.00000
E(xc) -100.98605 -100.98605 -100.98605 -0.00000 -0.00000 0.00000
Local -118.08797 -118.08797 -118.08797 0.00000 -0.00000 0.00000
n-local 309.70644 312.31673 307.06329 -3.39890 -2.83265 -2.51354
augment -46.57837 -46.57837 -46.57837 0.00000 0.00000 0.00000
Kinetic 218.93795 228.79224 227.35762 -4.58272 -6.30846 -4.07425
Fock 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
Total -0.10981 -0.10981 -0.10981 0.00000 0.00000 -0.00000
in kB -1.07575 -1.07575 -1.07575 0.00000 0.00000 -0.00000
external pressure -1.08 kB Pullay stress = 0.00 kB其中,Total后面的六个数字是位力(Virial),单位是eV,六个数值分别对应XX, YY, ZZ, XY, YZ, ZX六个分量。请自行补“位力”相关的物理背景知识。in kB后面的六个数字也是Virial,单位是 kbar(1GPa=10kB,1B就是1bar,即一个标准大气压),数值上直接等于Total那一栏的数字除以当前这一帧结构的体积,注意单位换算。这个单位下的位力,可以简单理解为当前结构的静态压强。可以发现此时整个体系会有个负压强external pressure = -1kB(对应0.1GPa,1GPa=10kB,1B就是1bar,即一个标准大气压),说明晶格在此时有自发向内收缩的趋势;如果这个数值是正的,意味着晶格在当前压强下有膨胀趋势。对于晶体来说,0.1GPa对晶格常数、带隙等性质影响非常弱。
压强之后会输出原子受力信息:
POSITION TOTAL-FORCE (eV/Angst)
1.36718 1.36718 1.36718 0.000000 000000*0 -0.000000
0.00000 2.73436 2.73436 -0.000000 0.000000 0.000000
1.36718 4.10155 4.10155 000000*0 -0.000000 -0.000000
0.00000 0.00000 0.00000 0.000000 -0.000000 0.000000
4.10155 1.36718 4.10155 -0.000000 0.000000 -0.000000
2.73436 2.73436 0.00000 -0.000000 000000*0 0.000000
4.10155 4.10155 1.36718 0.000000 -0.000000 -0.000000
2.73436 0.00000 2.73436 -0.000000 -0.000000 0.000000
total drift: -0.000000 -0.000000 0.000000由于我们的结构POSCAR中原子位置非常高对称,因此每个原子受力TOTAL-FORCE都为0。这只是我们当前例子比较特殊导致的巧合。
最后是关键的能量信息:
FREE ENERGIE OF THE ION-ELECTRON SYSTEM (ev)
free energy TOTEN -43.39328400eV
energy without entropy= -43.39328400 energy(sigma->0) 43.39328400这之后就是整个计算所运行的时间(1.78秒)
General timing and accounting informations for this job:
Total CPU time used (sec): 1.783
User time (sec): 1.399
System time (sec): 0.384
Elapsed time (sec): 2.802
Maximum memory used (kb): 221260.
Average memory used (kb) : N/A
Minor page faults: 21065
Major page faults:
Voluntary context switches: 4565当你看到OUTCAR最后出现这个信息时,意味着计算成功运行且结束了,但并不意味着自洽计算部分收敛了。需要结合前面提到的“EDIFF is reached”信息以及OSZICAR或者result文件中的电子步。
(3)OSZICAR,输入cat OSZICAR并回车
(base) [zpwu@muoo3 Si]$ cat OSZICAR
N E dE d eps ncg rms rms(c)
DAV: 1 0.245765894816E+02 0.24577E+02 -0.11363E+04 256 0.658E+02
DAV: 2 -0.432737377084E+02 -0.67850E+02 -0.67296E+02 384 0.786E+01
DAV: 3 -0.440666097061E+02 -0.79287E+00 -0.79287E+00 384 0.112E+01
DAV: 4 -0.440672202112E+02 -0.61051E-03 -0.61051E-03 512 0.371E-01
DAV: 5 -0.440672202424E+02 -0.31149E-07 -0.31149E-07 320 0.312E-03 0.948E+00
DAV: 6 -0.436028846451E+02 0.46434E+00 -0.30374E-01 416 0.144E+00 0.580E+00
DAV: 7 -0.433887736873E+02 0.21411E+00 -0.59457E-01 384 0.215E+00 0.326E-01
DAV: 8 -0.433920258292E+02 -0.32521E-02 -0.13379E-02 320 0.432E-01 0.155E-01
DAV: 9 -0.433931954126E+02 -0.11696E-02 -0.13474E-03 480 0.138E-01 0.106E-01
DAV: 10 -0.433932757836E+02 -0.80371E-04 -0.66202E-05 320 0.313E-02 0.792E-03
DAV: 11 -0.433932822670E+02 -0.64834E-05 -0.83780E-06 480 0.103E-02 0.568E-03
DAV: 12 -0.433932839156E+02 -0.16486E-05 -0.26691E-06 320 0.646E-03 0.227E-03
DAV: 13 -0.433932839974E+02 -0.81813E-07 -0.83379E-09 320 0.422E-04
1 F=-.43393284E+02 E0= -.43393284E+02 dE =-.200187E-14说明程序运行了13个电子步达到收敛,每个电子步有一个dE,当dE小于INCAR中设置的EDIFF=1E-6时,scf结束。
并且输出了基态总能量E0,以及F还有一个d E=-0.829E-3(这个d E和另外一个dE不同,前者是以正值EDIFFG为收敛阈值的)
我们再查看一下OUTCAR中最后一个电子步给出的能量信息
Free energy of the ion-electron system (eV)
alpha Z PSCENC 13.17272099
Ewald energy TEWEN -907.77962997
-Hartree energ DENC -60.71687757
-exchange EXHF 0.00000000
-v(xc)+e(xc) XCENC -72.23214521
PAW double counting 3311.55222047 -3171.62851930
entropy T*S EENTRO -0.00000000
eigenvalues EBANDS 19.71297677
atomic energy EATOM 824.52596983
Solvation Ediel_sol 0.00000000
free energy TOTEN -43.39328400eV
energy without entropy = -43.39328400 energy(sigma->0) = -43.39328400关注最下面的三个能量:free energy,energy without entropy和energy(sigma-0)。这里的free energy就是OSZICAR中的F,energy(sigma-0)就是OSZICAR的E0. Free energy是对蓝色框中的全部能量求和,其中EATOM来源于POTCAR(可以自行查看POTCAR前面几行)。E0就是在F的基础上扣掉了电子熵对能量的贡献(entropy TS)并且对展宽sigma取极限得到的,所以基态能量取E0。由于金刚石Si在常压下是半导体,因此不需要考虑电子熵的贡献问题(费米面上没有电子态占据);但如果是金属,电子熵影响就需要考虑了。无论如何,在常规的计算中,以E0即energy(sigma-0)作为体系的能量。特殊地,某些高温下的金属体系模拟,是需要考虑电子熵的。关于能量E0的理解可以阅读:https://mp.weixin.qq.com/s/1rsCjl-f_w6Es1wZ6kw_Dw
需要强调的是,这里直接算出来的能量可以简单理解为是电子和离子相互作用的势能,整个体系被认为是完全静止的。没有考虑原子振动对能量的贡献,这部分能量需要通过计算声子获得,超出了这一节范畴。
其余输出文件除了WAVECAR之后,其余都可以直接以文本形式打开查看。在大部分计算中,用到较多的输出文件除了上面三个result,OUTCAR,OSZICAR之外,还有CONTCAR,CHGCAR,WAVECAR(不能直接打开查看),DOSCAR,EIGENVAL,IBZKPT,XDATCAR(结构驰豫或者AIMD过程中的原子轨迹)。这些文件的内容请查阅官网,直接在vaspwiki搜索栏输入文件名(大写)即可。