1、sii si14si.4.hgh #赝势的文件名#ENDINP文件的内容为:# Crystalline cubic Si #ndtset 11 #说明下面将有11组数据acell: 3*9.8112 # 晶格常数a=b=c,将从9.8112. a.u.开始增加acell+ 3*0.09 #晶格常数将以0.09 a.u.的间隔进行增加#Ground state calculationkptopt 1 #在 k点网格取样时根据对称性来取样,并由下面的 # ngkpt和kptrlatt, 或者nshiftk和shiftk来确定k点的数目iscf 5 #采用CG方法对能量进行优化,用在基态计算中。#D
2、efinition of the unit cellrprim 0.0 0.5 #下面三行定义了原胞的基矢,本例子中Si是fcc结构#Definition of the atom typesntypat 1 #定义原胞中原子的类别的数目,本例子中只有1类原子znucl 14 #定义原胞中原子的核电荷数#Definition of the atomsnatom 2 #定义原胞中原子的总个数,本例子中有2个原子typat 2*1 #定义每类原子的个数,本例子中第一类原子有2个xred #下面定义了原胞中原子的坐标 0.25 0.25 0.25#Gives the number of band, e
3、xplicitely (do not take the default)nband 16 #定义了要计算的能带的数目,最好按这样来设置: # nband 原胞中总的价电子数目/2 + 10#Exchange-correlation functionalixc 1 #定义交换关联函数,本例子中,采用的是Teter Pade 参数化的LDA形式#Definition of the planewave basis setecut 20.0 #定义了平面波的切断动能 #Definition of the k-point gridngkpt 8 8 8 #下面定义了k点网格取样的大小nshiftk 4s
4、hiftk 0.5 0.5 0.5 0.5 0.0 0.0 0.0 0.5 0.0 0.0 0.0 0.5#Definition of the SCF procedurenstep 60 #电子自洽迭代的最大步数diemac 12.0 #介电常数设置tolvrs 1.0d-20 #电子自洽收敛的标准- END-计算完后,得到OUT文件,用下面的命令:grep volume OUT得到如下的内容:Unit cell volume ucvol= 2.3610688E+02 bohr3 2.4266422E+02 bohr3 2.4934185E+02 bohr3 2.5614088E+02 boh
5、r3 2.6306239E+02 bohr3 2.7010748E+02 bohr3 2.7727725E+02 bohr3 2.8457279E+02 bohr3 2.9199518E+02 bohr3 2.9954553E+02 bohr3 3.0722493E+02 bohr3然后用下面的命令:Etotal Etotal= -7.92750029752797E+00 Etotal= -7.92997465524506E+00 Etotal= -7.93167675973445E+00 Etotal= -7.93266612552653E+00 Etotal= -7.93299797094
6、926E+00 Etotal= -7.93272412167288E+00 Etotal= -7.93189304516315E+00 Etotal= -7.93055314505404E+00 Etotal= -7.92874706916830E+00 Etotal= -7.92651686675655E+00 Etotal= -7.92390242137627E+00因此,Volume 和Etotal对应的关系为:2.3610688E+02 -7.92750029752797E+002.4266422E+02 -7.92997465524506E+002.4934185E+02 -7.93
7、167675973445E+002.5614088E+02 -7.93266612552653E+002.6306239E+02 -7.93299797094926E+002.7010748E+02 -7.93272412167288E+002.7727725E+02 -7.93189304516315E+002.8457279E+02 -7.93055314505404E+002.9199518E+02 -7.92874706916830E+002.9954553E+02 -7.92651686675655E+003.0722493E+02 -7.92390242137627E+00下面就可
8、以用 Birch-Murnaghan 3阶状态方程进行 (Birch F, Phys. Rev. 71, p809 (1947)拟合得到体弹性模量和平衡状态下的体积:V0 = 263.276940709097 a.u.3B0 95.497(GPa)第二章 ABINIT参数设定与收敛测试 ABINIT计算晶体时主要参数的确定(切断动能和k点网格)。 在采用平面波赝势法进行固体的电子结构计算时,为了确保计算的精度和计算结果的可靠性,在计算晶体的物理性质之前,要进行几个重要参数的测试,以保证这些参数的选取使得计算结果有很好的收敛性,这些参数是平面波的切断动能和k点网格取样的大小。下面以采用ABINI
9、T计算立方的ZrO2晶体为例:采用的赝势是40zr.psp_mod和8o.psp_mod(它们均是LDA的TM赝势)。输入文件,in.files的内容为:INP OUT zroi zroo zro ./40zr.psp_mod ./8o.psp_mod 在测试平面波切断动能的收敛性时,我们通过设置平面波切断动能从20 Ha开始,以2 Ha递增,直到58 Ha,其他参数不变计算ZrO2的总能。在ABINIT的输入文件中很方便的通过ndtset来设置。输入文件INP如下(紫色标示):# Crystalline ZrO2-cubic ndtset 20 #表示有20组数据 ecut: 20.0 #平
10、面波切断动能从20 Ha开始 ecut+ 2 #以2 Ha递增,也就是Ecut为20.0 + i*2.0, (i从1到20) #Definition of the unit cell acell 3*9.65285 #设置晶格常数a=b=c为9.65285 a.u. rprim 0.0 0.5 0.5 #同上面的acell确定了原胞的基矢 0.5 0.0 0.5 0.5 0.5 0.0 #Definition of the atom types ntypat 2 #设置原胞的原子种类数,这里有2类原子 znucl 40 8 #每类原子的核电荷数Z natom 3 #原胞中总的原子数目,这里原胞
11、总共有3个原子 typat 1 2*2 #第一类原子有1个,第二类原子有2个 xred #下面的按分数坐标给原胞中原子的坐标位置 0.00 0.000 0.000 0.25 0.25 0.25 0.75 0.75 0.75 #Definition of the k-point grid kptopt 1 #设置生成k点的方法,这里表明有ngkpt和nshiftk来确定k点网格的大小 ngkpt 8 8 8 #设置对布里渊区进行8x8x8网格的划分 nshiftk 4 #对划分得到的k点按下面的偏移量进行平移 shiftk 0.5 0.5 0.5 0.5 0.0 0.0 0.0 0.5 0.0
12、0.0 0.0 0.5 #Definition of the SCF procedure #ecut 40.0 iscf 5 #自洽迭代中的算法,这里采用cg方法优化。toldfe 1.0d-10 #自洽迭代迭过程总能收敛的标准 diemac 3.0 #体系宏观的介电常数 nstep 60 #最大自洽迭代步数 #Definition of the outfile information prtwf 0 #不输出波函数文件。运行下面的命令 abins OUT 对OUT采用进行分析,用下面的命令取出一些数据:ecut OUT summary Total energy(eV)得到summary类似如
13、下数据:2.00000000E+01 -2.16713704994688E+03 2.20000000E+01 -2.17244914205186E+03 2.40000000E+01 -2.17629539439455E+03 2.60000000E+01 -2.17888687489711E+03 2.80000000E+01 -2.18045824630388E+03 3.00000000E+01 -2.18139622512928E+03 3.20000000E+01 -2.18192168325676E+03 3.40000000E+01 -2.18219516090481E+03
14、3.60000000E+01 -2.18232706572898E+03 3.80000000E+01 -2.18237668482667E+03 4.00000000E+01 -2.18239923089573E+03 4.20000000E+01 -2.18241955706857E+03 4.40000000E+01 -2.18244775901834E+03 4.60000000E+01 -2.18248714583389E+03 4.80000000E+01 -2.18253519502611E+03 5.00000000E+01 -2.18258901301065E+03 5.20
15、000000E+01 -2.18264245486986E+03 5.40000000E+01 -2.18269485143783E+03 5.60000000E+01 -2.18274313138312E+03 5.80000000E+01 -2.18278585826504E+03 画图可以看到当Ecut取40 Ha时,体系的总能有很好的收敛了。下面测试对k点网格的,in.files的输入文件同上。我们这里把k点网格从6x6x6开始增加,以2x2x2递增。INP的输入文件(紫色标示):ndtset 10 ngkpt: 6 6 6 ngkpt+ 2 2 2 acell 3*9.65285 0
16、.0 0.5 0.5 ntypat 2 znucl 40 8 natom 3 typat 1 2*2 xred kptopt 1 #ngkpt 8 8 8 nshiftk 4 shiftk 0.5 0.5 0.5 ecut 40.0 iscf 5 toldfe 1.0d-10 diemac 3.0 nstep 60 prtwf 0 运行 abinis OUT&计算完后对OUT进行分析, ngkptcomment 取出得到的数据如下:28.0000 -2182.416564 60.0000 -2182.416514 110.0000 -2182.416506 182.0000 -2182.416
17、526 280.0000 -2182.416499 408.0000 -2182.416512 570.0000 -2182.416506 770.0000 -2182.416514 1012.0000 -2182.416515 1300.0000 -2182.416504 画图可以看到8x8x8的k点网格就能保证体系的总能很好的收敛了。第三章 ABINIT计算晶体的能带结构采用第一原理的电子结构计算方法来计算晶体的能带结构一般来说,要进行两个步骤。这不论是采用VASP、PWSCF还是ABINIT这些程序。步骤为:先进行自洽的电子结构迭代得到自洽计算得到的电荷密度,然后读入这个自洽得到的电荷密
18、度,进行非自洽的计算得到体系的本征值。下面以立方的ZrO2晶体为例,采用ABINIT来计算。下面采用ndset这个关键词,在输入文件中输入这两步的控制参数,一次性计算完成得到能带结构。并对相关的相关的输入参数进行解释:采用的赝势是:40zr.psp_mod和8o.psp_mod,输入文件 in.files的内容为 此in.files文件的内容设置了主要输入文件的名称和赝势文件所在的目录。主要输入文件INP的内容为(紫色标示):ndtset 2 #表示有两组控制参数:第一组控制参数用来设置自洽计算,第二组是非自洽的本征值计算 #Dataset 1 : usual self-consistent
19、calculation kptopt1 1 #自洽计算中设置k点网格取样的方法,表明采用ngkpt和shfitk来设置k点网格。nshiftk1 4 #使生成的k点进行平移。shiftk1 ngkpt1 8 8 8 #K点网格取样,网格划分的分割数 prtden1 1 #表明输出电荷密度文件。tolvrs1 1.0d-20 # 自洽计算收敛的标准 iscf1 5 #自洽迭代计算时,采用 CG方法来优化有效势 #Dataset 2 : #用来设置本征值计算时的参数 iscf2 -2 #表明非自洽计算 getden2 -1 #读入上一组数据进行自洽计算得到的电荷密度文件 kptopt2 -5 #负
20、数,表示下面计算能带计算时,有5段特殊线(由6个特殊k点来确定) ndivk2 10 12 18 8 8 #每段特殊线上分几等份 kptbounds2 #特殊k点的坐标 0.5 0.25 0.75 # W point 0.5 0.0 0.0 # L point 0.0 0.0 0.0 # Gamma point 0.5 0.0 0.5 # X point 0.0 0.0 0.0 # Gamma enunit2 0 #输出本征值时,本征值的单位,这里为0表示是以Hatree为单位给出 prteig2 1 #表明输出本征值到文件中 acell 3*9.48196 #设置晶常数 0.0 0.5 0.
21、5 #设置计算原胞的基矢(同上面的acell一起构成原胞的基矢) ntypat 2 #原胞的原子种类数目 znucl 40 8 #每类原子的核电荷数 natom 3 #原胞中总的原子数目 typat 1 2*2 #表示第一个原子是第一类的,后面两个原子是第二类的。xred #以分数坐标给出原子的位置 nband 30 #在计算中考虑多少条能带 ecut 40.0 #平面波切断动能 diemac 3.0 #体系的宏观介电常数,给一个近似值就可以了。nstep 60 #自洽迭代时的最大步数 计算得到的本征值文件zroo_DS2_EIG内容为:Eigenvalues (hartree) for nk
22、pt= 57 k points:kpt# 1, nband= 30, wtk= 1.00000, kpt= 0.5000 0.2500 0.7500 (reduced coord) -1.41786 -0.62938 -0.62938 -0.61620 -0.26378 -0.26378 0.14567 0.14721 0.23910 0.23910 0.27711 0.32275 0.49670 0.50302 0.57458 0.57458 0.61760 0.78769 0.90930 0.92430 0.92430 1.05159 1.05159 1.10749 1.20089 1.28891 1
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1