分类:
2009-10-22 08:54:22
Abinit:采用状态方程拟合晶格常数和体弹性模量。
| |
下面举例如何做好输入文件计算金刚石结构的Si的状态方程,由此也得到Si的晶格常数和体弹性模量。 大致的步骤是在Si的晶格常数的实验值附近取11个数据点,也就是说取11个晶格常数或体积,然后计算在这些晶格常数下的总能。在计算得到总能后,采用状态方程拟合得到状态方程、平衡态时的体积(或晶格常数)和体弹性模量。 本例子中采用的是LDA-HGH赝势。赝势的文件名为:14si.4.hgh。 in.files的内容为: #####Begin INP #设置关键词的文件名为INP OUT #主要的输出文件为OUT,该文件将被写入计算最重要的结果 sii si si 14si.4.hgh #赝势的文件名 ######END INP文件的内容为: # 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 calculation kptopt 1 #在 k点网格取样时根据对称性来取样,并由下面的 # ngkpt和kptrlatt, 或者nshiftk和shiftk来确定k点的数目 iscf 5 #采用CG方法对能量进行优化,用在基态计算中。 ####################################################################### #Definition of the unit cell rprim 0.0 0.5 0.5 #下面三行定义了原胞的基矢,本例子中Si是fcc结构 0.5 0.0 0.5 0.5 0.5 0.0 #Definition of the atom types ntypat 1 #定义原胞中原子的类别的数目,本例子中只有1类原子 znucl 14 #定义原胞中原子的核电荷数 #Definition of the atoms natom 2 #定义原胞中原子的总个数,本例子中有2个原子 typat 2*1 #定义每类原子的个数,本例子中第一类原子有2个 xred #下面定义了原胞中原子的坐标 0.0 0.0 0.0 0.25 0.25 0.25 #Gives the number of band, explicitely (do not take the default) nband 16 #定义了要计算的能带的数目,最好按这样来设置: # nband > 原胞中总的价电子数目/2 + 10 #Exchange-correlation functional ixc 1 #定义交换关联函数,本例子中,采用的是Teter Pade 参数化的LDA形式 #Definition of the planewave basis set ecut 20.0 #定义了平面波的切断动能 #Definition of the k-point grid ngkpt 8 8 8 #下面定义了k点网格取样的大小 nshiftk 4
shiftk 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 procedure nstep 60 #电子自洽迭代的最大步数 diemac 12.0 #介电常数设置 tolvrs 1.0d-20 #电子自洽收敛的标准 ------------ END-------------- 计算完后,得到OUT文件,用下面的命令: grep 'volume' OUT 得到如下的内容: Unit cell volume ucvol= 2.3610688E+02 bohr^3 Unit cell volume ucvol= 2.4266422E+02 bohr^3 Unit cell volume ucvol= 2.4934185E+02 bohr^3 Unit cell volume ucvol= 2.5614088E+02 bohr^3 Unit cell volume ucvol= 2.6306239E+02 bohr^3 Unit cell volume ucvol= 2.7010748E+02 bohr^3 Unit cell volume ucvol= 2.7727725E+02 bohr^3 Unit cell volume ucvol= 2.8457279E+02 bohr^3 Unit cell volume ucvol= 2.9199518E+02 bohr^3 Unit cell volume ucvol= 2.9954553E+02 bohr^3 Unit cell volume ucvol= 3.0722493E+02 bohr^3 然后用下面的命令: grep 'Etotal' OUT 得到如下的内容: >>>>>>>>> Etotal= -7.92750029752797E+00 >>>>>>>>> Etotal= -7.92997465524506E+00 >>>>>>>>> Etotal= -7.93167675973445E+00 >>>>>>>>> Etotal= -7.93266612552653E+00 >>>>>>>>> Etotal= -7.93299797094926E+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+00 2.4266422E+02 -7.92997465524506E+00 2.4934185E+02 -7.93167675973445E+00 2.5614088E+02 -7.93266612552653E+00 2.6306239E+02 -7.93299797094926E+00 2.7010748E+02 -7.93272412167288E+00 2.7727725E+02 -7.93189304516315E+00 2.8457279E+02 -7.93055314505404E+00 2.9199518E+02 -7.92874706916830E+00 2.9954553E+02 -7.92651686675655E+00 3.0722493E+02 -7.92390242137627E+00 下面就可以用 Birch-Murnaghan 3阶状态方程进行 (Birch F, Phys. Rev. 71, p809 (1947))拟合得到体弹性模量和平衡状态下的体积: V0 = 263.276940709097 a.u.^3 B0 = 95.497(GPa) |
原文链接http://zhaoyunawuli.blog.163.com/blog/static/70213976200811194275577/