1、第一性原理计算的软件实现及其在材料设计中应用报告人:曹勇 博士生导 师:朱景川 教授2011.03.09n n1 概述n n2 计算软件CASTEP的简介及工作环境n n3 模型搭建n n4 几何优化n n5 计算结果分析n n6 应用实例1 概述n n第一原理计算指的是只依赖于基本的物理常数和原子质量,第一原理计算指的是只依赖于基本的物理常数和原子质量,通过量子力学理论求解运动方程来确定电子行为,进而预通过量子力学理论求解运动方程来确定电子行为,进而预测物质性质。由于不借助任何经验、半经验理论或试验结测物质性质。由于不借助任何经验、半经验理论或试验结果,因此也被称为从头算方法(果,因此也被称
2、为从头算方法(abab initio method)initio method)。n n仅需仅需5 5个基本物理常数,即电子的静止质量个基本物理常数,即电子的静止质量mm、电子电量、电子电量e e、普朗克常数普朗克常数h h、光速、光速c c和波尔玆曼常数和波尔玆曼常数k kB B n n第一原理计算方法是从最基本和最简单的物理事实出发,第一原理计算方法是从最基本和最简单的物理事实出发,计算结果往往能够具有很高的精确度和可信度,成为独立计算结果往往能够具有很高的精确度和可信度,成为独立与试验之外研究物质性质的又一利器。与试验之外研究物质性质的又一利器。1第一原理的可靠性Total Energy
3、,Cohesive EnergyMetalMetalDFT-LDADFT-LDAExperimentExperimentNbNb0.550.550.560.56MoMo0.490.490.500.50RhRh0.450.450.420.42PdPd0.270.270.280.28Total energyof atoms(Ry.)Cohesive energy(Ry.)Structural Properties:MoleculesStructural properties of diborane a ab bd d1 1()()d d2 2()()DFT-LDADFT-LDA121121o o9
4、898o o1.321.321.201.20DFT-GGADFT-GGA122122o o9797o o1.331.331.191.19EXPEXP122122o o9797o o1.331.331.191.19Diborne(3c-2e bond)Lattice constant(a.u.)LDA predicts latticeconstants to within 2%.Bulk modulus(Mbar)LDA predicts bulk modulusto within a few percents.Structural PropertiesLatticeconstantMagnet
5、ism:BCC Fe(DFT-GGA)-2-1 0 1 2-8-6-4-2 0 2 4Energy(eV)DOSMajority spinMinority spinMeasured:Measured:m mB B=2.22/atom=2.22/atomCalculated:Calculated:m mB B=2.25/atom=2.25/atom-2-1 0 1 2-8-6-4-2 0 2 4Magnetism:HCP Co(DFT-GGA)Energy(eV)DOSMajority spinMinority spinMeasured:Measured:m mB B=1.72/atom=1.7
6、2/atomCalculated:Calculated:m mB B=1.63/atom=1.63/atomMagnetism:CCP Ni(DFT-GGA)Energy(eV)DOSMajority spinMinority spin-1.5-1-0.5 0 0.5 1 1.5-8-6-4-2 0 2 4Measured:Measured:m mB B=0.62/atom=0.62/atomCalculated:Calculated:m mB B=0.65/atom=0.65/atomMagnetic PropertiesMetalMetalCrystalCrystalstructurest
7、ructureatomic atomic configurationconfigurationMeauredMeaured (m mB B/atom)/atom)Calculated Calculated(DFT)(DFT)(m mB B)FeFeBCCBCC3d3d6 64s4s2 22.222.222.252.25CoCoHCPHCP3d3d7 74s4s2 21.721.721.631.63NiNiCCPCCP3d3d9 94s4s2 20.620.620.650.65n n由于在处理多电子体系相互作用时遇到的困难由于在处理多电子体系相互作用时遇到的困难,第一原理第一原理在实际计算中也需
8、要做一定程度的近似,不过这样的近似在实际计算中也需要做一定程度的近似,不过这样的近似也都是基于物理上的简单假设和数学形式上的合理简化而也都是基于物理上的简单假设和数学形式上的合理简化而得到的。得到的。Can the materials properties be understood knowing only the identities(atomic Can the materials properties be understood knowing only the identities(atomic numbers)of their constituent atoms?numbers)o
9、f their constituent atoms?CSiGeNBStructural properties?Electronic properties?Elastic properties?optical properties?Phonon dispersion?Magnetic properties?Superconductivity?2第一原理所能解决的问题Hexagoal Mo2CMo2C 是合金钢中的二次硬化相,它均匀、细小、弥散的分布于基体中,具有极高的硬度和良好的高温稳定性。合金钢中合金元素种类繁多,是否存在其他的碳化物相,如Cr2C、V2C、Ni2C来取代其作用。Mo2C 是否
10、具有其他更加稳定的结构 空位或其他原子置换对其稳定性和性能有何影响2 第一原理所能解决的问题Orthorhombic and Hexagonal Fe3C晶体结构的稳定性判定虚拟相、亚稳相、稳定相。研究晶格常数和原子所处位置对晶体结构稳定性的影响。对于组成元素确定的体系,可能存在不同的晶体结构。分别对不同的晶体结构的能量极小值进行计算,比较不同晶体结构的能量极小值,便可确定稳定的晶体结构。2 第一原理所能解决的问题Material Properties from First-PrinciplesFrom first principles!Predict new behaviors/proper
11、ties of existing materialsDesign materials with desired propertiesUnderstand and explain materials propertiesBecoming realityn n目前已有一些利用第一原理对具有周期性结构的材料进行目前已有一些利用第一原理对具有周期性结构的材料进行能带计算的商业软件,如能带计算的商业软件,如VASPVASP,ABINITABINIT,Wien2KWien2K,CASTEPCASTEP。2 计算软件CASTEP的简介及运行环境n n计算软件计算软件CASTEP(Cambridge Sequ
12、ential Total Energy CASTEP(Cambridge Sequential Total Energy Package)Package)是基于密度泛函理论的第一性原理计算量子力是基于密度泛函理论的第一性原理计算量子力学程序。其使用的是总能量平面波赝势的方法。学程序。其使用的是总能量平面波赝势的方法。n n密度泛函理论密度泛函理论 (Density Function Theory)(Density Function Theory)巧妙的将电子巧妙的将电子之间的交换势表示为密度泛函,从而使得薛定谔方程在考之间的交换势表示为密度泛函,从而使得薛定谔方程在考虑了电子之间的复杂相互作用
13、之后,依然可以用自洽的方虑了电子之间的复杂相互作用之后,依然可以用自洽的方法求解。法求解。关于关于CASTAPCASTAP CASTAPCASTAP是特别为固体材料学而设计的一个现代的量子力学基本程序,是特别为固体材料学而设计的一个现代的量子力学基本程序,其使用了密度泛函其使用了密度泛函(DFT(DFT)平面波赝势方法,进行第一原理量子力学计算,平面波赝势方法,进行第一原理量子力学计算,以探索如半导体,陶瓷,金属,矿物和沸石等材料的晶体和表面性质。以探索如半导体,陶瓷,金属,矿物和沸石等材料的晶体和表面性质。典型的应用包括表面化学,键结构,态密度和光学性质等研究,典型的应用包括表面化学,键结构
14、,态密度和光学性质等研究,CASTAPCASTAP也可用于研究体系的电荷密度和波函数的也可用于研究体系的电荷密度和波函数的3D3D形式。此外,形式。此外,CASTAPCASTAP可可用于有效研究点缺陷(空位,间隙和置换杂质)和扩展缺陷(如晶界和位用于有效研究点缺陷(空位,间隙和置换杂质)和扩展缺陷(如晶界和位错)的性质。错)的性质。Material StudioMaterial Studio使用组件对话框中的使用组件对话框中的CASTAPCASTAP选项允许准备,启动,选项允许准备,启动,分析和监测分析和监测CASTAPCASTAP服役工作。服役工作。计算:允许选择计算选项(如基集,交换关联势
15、和收敛判据),作业控制计算:允许选择计算选项(如基集,交换关联势和收敛判据),作业控制和文档控制。和文档控制。分析:允许处理和演示分析:允许处理和演示CASTAPCASTAP计算结果。这一工具提供加速整体直观化以计算结果。这一工具提供加速整体直观化以及键结构图,态密度图形和光学性质图形。及键结构图,态密度图形和光学性质图形。CASTAPCASTAP计算计算 CASTAPCASTAP计算是要进行的三个任务中的一个,即单个点的能量计算,计算是要进行的三个任务中的一个,即单个点的能量计算,几何优化或分子动力学。可提供这些计算中的每一个以便产生特定的几何优化或分子动力学。可提供这些计算中的每一个以便产
16、生特定的物理性能。性质为一种附加的任务,允许重新开始已完成的计算以便物理性能。性质为一种附加的任务,允许重新开始已完成的计算以便产生最初没有提出的额外性能。产生最初没有提出的额外性能。在在CASTAPCASTAP计算中有很多运行步骤,可分为如下几组:计算中有很多运行步骤,可分为如下几组:*结构定义:必须规定包含所感兴趣结构的周期性的结构定义:必须规定包含所感兴趣结构的周期性的3D3D模型文件,有模型文件,有大量方法规定一种结构:可使用构建晶体(大量方法规定一种结构:可使用构建晶体(Build Crystal)Build Crystal)或构建真或构建真空板空板(Build Vacuum Sta
17、b)(Build Vacuum Stab)来构建,也可从已经存在的的结构文档中引来构建,也可从已经存在的的结构文档中引入,还可修正已存在的结构。入,还可修正已存在的结构。提示:提示:CASTAPCASTAP计算所需时间随原子数平方的增加而增加。因此,建议计算所需时间随原子数平方的增加而增加。因此,建议是用最小的初晶胞来描述体系,可使用是用最小的初晶胞来描述体系,可使用BuildSymmetryPrimitive BuildSymmetryPrimitive CellCell菜单选项来转换成初晶胞。菜单选项来转换成初晶胞。CASTAPCASTAP中中选择一项任务选择一项任务1 1 从模块面板(从
18、模块面板(Module Explorer)Module Explorer)选择选择CASTAPCalculationCASTAPCalculation。2 2 选择设置表。选择设置表。3 3 从任务列表中选择所要求的任务。从任务列表中选择所要求的任务。计算设置:合适的计算设置:合适的3D3D模型文件一旦确定,必须选择计算类型模型文件一旦确定,必须选择计算类型和相关参数,例如,对于动力学计算必须确定系综和参数,包和相关参数,例如,对于动力学计算必须确定系综和参数,包括温度,时间步长和步数。选择运行计算的磁盘并开始括温度,时间步长和步数。选择运行计算的磁盘并开始CASTAPCASTAP作业。作业。
19、*结果分析:计算完成后,相关于结果分析:计算完成后,相关于CASTAPCASTAP作业的文档返回用作业的文档返回用户,在项目面板适当位置显示。这些文档的一些进一步处理要户,在项目面板适当位置显示。这些文档的一些进一步处理要求获得可观察量如电子结构图、光学性质。求获得可观察量如电子结构图、光学性质。CASTAPCASTAP能量任务能量任务CASTAPCASTAP能量任务允许计算特定体系的总能量以及物理性质。能量任务允许计算特定体系的总能量以及物理性质。除了总能量之外,在计算之后还可报告作用于原子上的力;也除了总能量之外,在计算之后还可报告作用于原子上的力;也能创建电荷密度文件;利用材料观测仪(能
20、创建电荷密度文件;利用材料观测仪(Material Material VisualizerVisualizer)允许目测电荷密度的立体分布;还能报告计算中使允许目测电荷密度的立体分布;还能报告计算中使用的用的Monkhorst-Park的的k点的电子能量,因此在点的电子能量,因此在CASTAPCASTAP分析中分析中可生成态密度图。可生成态密度图。对于能够得到可靠结构信息的体系的电子性质的研究,能量任对于能够得到可靠结构信息的体系的电子性质的研究,能量任务是有用的。只要给定应力性质,也可用于计算没有内部自由务是有用的。只要给定应力性质,也可用于计算没有内部自由度的高对称性体系的状态方程(即压力
21、度的高对称性体系的状态方程(即压力-体积,能量体积,能量-体积关系)体积关系)。CASTAPCASTAP中能量的缺损单位是电子伏特中能量的缺损单位是电子伏特(eVeV),各种能量单位的换,各种能量单位的换算关系见算关系见Mohr.P.J(2000).Mohr.P.J(2000).1 eV=0.036749308 Ha=23.0605 kcal/mole=96.4853 kJ/moleCASTAPCASTAP几何优化任务几何优化任务 CASTAPCASTAP几何优化任务允许改善结构的几何,获得稳定结构几何优化任务允许改善结构的几何,获得稳定结构或多晶型物。通过一个迭代过程来完成这项任务,迭代过程
22、中调或多晶型物。通过一个迭代过程来完成这项任务,迭代过程中调整原子坐标和晶胞参数使结构的总能量最小化。整原子坐标和晶胞参数使结构的总能量最小化。CASTAPCASTAP几何优化是基于减小计算力和应力的数量级,直到小几何优化是基于减小计算力和应力的数量级,直到小于规定的收敛误差。也可能给定外部应力张量来对拉应力、压应于规定的收敛误差。也可能给定外部应力张量来对拉应力、压应力和切应力等作用下的体系行为模型化。在这些情况下力和切应力等作用下的体系行为模型化。在这些情况下反反复迭复迭代内部应力张量直到代内部应力张量直到与所施加的外部应力相等。与所施加的外部应力相等。几何优化处理产生的模型几何优化处理产
23、生的模型结构与真实结构紧密相似。结构与真实结构紧密相似。利用利用CASTAPCASTAP计算的晶格参计算的晶格参数精度列于右图。数精度列于右图。CASTAPCASTAP性质任务性质任务 CASTAPCASTAP性质任务允许在完成能量,几何优化或动力学运行之后求出电性质任务允许在完成能量,几何优化或动力学运行之后求出电子和结构性质。可以产生的性质如下:子和结构性质。可以产生的性质如下:*态密度(态密度(DOSDOS):利用原始模拟中产生的电荷密度和势能,非自恰计算价):利用原始模拟中产生的电荷密度和势能,非自恰计算价带和导带的精细带和导带的精细MonkhorstMonkhorst-Pack-Pa
24、ck 网格上的电子本征值。网格上的电子本征值。*带结构:利用原始模拟中产生的电荷密度和势能,非自恰计算价带和导带结构:利用原始模拟中产生的电荷密度和势能,非自恰计算价带和导带的布里渊区高对称性方向电子本征值。带的布里渊区高对称性方向电子本征值。*光学性质:计算电子能带间转变的矩阵元素。光学性质:计算电子能带间转变的矩阵元素。CASTAPCASTAP分析对话可用于生分析对话可用于生成包含可以测得的光学性质的网格和图形文件。成包含可以测得的光学性质的网格和图形文件。*布局数分析:进行布局数分析:进行MullikenMulliken 分析。计算决定原子电荷的键总数和角动量分析。计算决定原子电荷的键总
25、数和角动量(以及自旋极化计算所需的磁矩)。任旋地,可产生态密度微分计算所要(以及自旋极化计算所需的磁矩)。任旋地,可产生态密度微分计算所要求的分量。求的分量。*应力:计算应力张量,并写入应力:计算应力张量,并写入seedname.castepseedname.castep 文档文档。启动程序CASTAP的运行环境的运行环境2.8.3 模型搭建搭建晶胞模型能带计算过程输出结果分析包括电子结构和物理性质通过计算能量极小值进行几何优化确立稳定结构MoSi2以MoSi2为例,其空间群为14/mmm;晶格常数为a=3.202,c=7.851;Mo的原子位置为(0,0,0)和(1/2,1/2,1/2);S
26、i的原子位置为(0,0,1/3),(0,0,2/3),(1/2,1/2,1/6)和(1/2,1/2,5/6)n n空间群是点阵、平移群空间群是点阵、平移群(滑移面和螺旋轴滑移面和螺旋轴)和点群的组合。和点群的组合。230230个空间群是由个空间群是由1414个个BravaisBravais点阵与点阵与3232个晶体点群系统组个晶体点群系统组合而成。合而成。n n能使无限大晶体能使无限大晶体(从微观角度来看从微观角度来看)自身重复的几自身重复的几何对称操作的集合。何对称操作的集合。原胞初级单(晶)胞原始单(晶)胞原胞初级单(晶)胞原始单(晶)胞primitive primitive(unit)c
27、ell(unit)celln n惯用晶胞惯用晶胞Conventional(unit)cellConventional(unit)celln n由一组基矢所决定的平行六面体所围起来的最小重复单元由一组基矢所决定的平行六面体所围起来的最小重复单元就叫就叫原胞原胞原胞原胞(或初基单胞(或初基单胞Primitive CellPrimitive Cell。(固体物理常用固体物理常用)n n原胞代表的是几何形状、体积的最小重复单元,不具有物原胞代表的是几何形状、体积的最小重复单元,不具有物理内涵。理内涵。n n原胞往往不能反映晶体的对称性,因而,习惯上常选择能原胞往往不能反映晶体的对称性,因而,习惯上常选
28、择能反映晶体对称性的重复单元,这种重复单元就叫惯用反映晶体对称性的重复单元,这种重复单元就叫惯用晶胞晶胞晶胞晶胞(conventional cellconventional cell)或非)或非初基单胞初基单胞初基单胞初基单胞(nonprimitivenonprimitive cell)cell)。n n晶胞一般不是最小的重复单元。其体积(面积)可以是原晶胞一般不是最小的重复单元。其体积(面积)可以是原胞的数倍。胞的数倍。寻找已经存在的常见晶体模型寻找已经存在的常见晶体模型寻找已经存在的常见晶体模型寻找已经存在的常见晶体模型n n建立所研究材料的结构模型建立所研究材料的结构模型 空间群空间群
29、晶格参数晶格参数 晶体的原子内坐标晶体的原子内坐标n n以以VCVC为例,为例,其空间群为其空间群为FM-3MFM-3M;晶格常数为晶格常数为a=4.160a=4.160;V V的原子位置为的原子位置为(0,0,0)(0,0,0);C C的原子位置的原子位置(1/2,1/2,1/2).(1/2,1/2,1/2).搭建模型选择空间群、设定晶格常数选择空间群、设定晶格常数选择空间群、设定晶格常数添加原子,输入原子内坐标建立超晶胞转换原胞设立P1空间群,降低对称性问题1 说出上述晶胞的晶胞类型,并指出后两个晶胞所含化学式VC的个数。问题2 P1空间群有何特点,什么情况下会将晶胞设置为P1空间群?2.
30、8.4.几何优化能带计算的过程晶体是一个具有周期性结构的体系,输入时只能给出一个体积有限的晶体结构模型,需要利用周期性边界条件,才能得到整个晶体的能带结构。在选择了计算任务之后,需要对计算参数进行设置,如自洽场计算的精度、基组大小、k的取值等。在能带计算中,基组是指用于线性组合成单电子波函数的基函数集合,在平面波赝势方法中,给出能量的截断值Ecut。由于计算量的的原因计算时k值只取简约布里渊区中的有限值,以了l1 l2 l3的形式表示所取的k在b1、b2和b3方向上取值间隔分别为b1/l1,b2/l2,b3/l3,。在其他条件不变的情况下,k的取值数目增加,得到的En(k)函数的精确度增加,但
31、计算量显著增加。自洽场计算的收敛用电子数密度或晶体总能量的收敛来标志,自洽场计算的精度是指收敛的标准。能带计算的过程 利用自洽场方法求解Kohn-Sham方程,得到所设结构的晶体总能量。在计算晶体的物理性质之前必须对所设的晶体结构模型进行几何优化,根据关于能量、力、应力、位移的判据来判断晶体结构是否为稳定结构(总能量最小)。如果晶体结构不是稳定结构,重新设置晶格参数进行计算,直至得到稳定的晶体结构。然后对结构优化后的晶体进行物理性质计算,输出计算结果。能带计算的过程几何优化几何优化是通过调节结构模型的几何参数来获得稳定结构的过程,其结果是使模型结构尽可能的接近真实结构。进行几何优化的判据可以根
32、据研究的需要而定,一般是几个判据的组合使用。常用的判据有以下几个:自洽场收敛判据。对给定的结构模型进行自洽场计算时,相继两次自洽计算得到的晶体总能量之差足够小,即相继两次自洽计算的晶体总能量之差小于设定值。力判据。每个原子所受的晶体内作用力足够小,即单个原子受力小于设定最大值。应力判据。每个结构模型单元中的应力足够小,即应力小于设定的最大值。位移判据。相继两次结构参数变化引起的原子位移的分量足够小,即原子位移分量小于设定的最大值。设置参数设置参数设置参数设置参数设置电子结构参数几何优化几何优化(a)(b)MoMoCCFig.1 The crystal structures of Mo2C wi
33、th 12 atoms basis in(a)Ortho-Mo2C,(b)Hexa-Mo2C.5 计算结果分析5.1 能量分析5.2 电子结构分析5.3 弹性性质分析5.4 其他物理性质分析5.1 能量分析晶体总能量晶体总能量(不包括核的动能部分)可分为两部分:一部分是原子核与内层电子组成的离子实能量,这部分能量基本上与晶体结构无关,是一个常数,赝势方法中常把总能量中这部分不变的能量设为零;另一部分是总能量与离子实能量之差,包括离子实与价电子的相互作用,离子实之间的相互作用以及价电子之间的相互作用。晶体总能量n nM.T.YinM.T.Yin和和M.L.CohenM.L.Cohen成功地利用第
34、一性原理成功地利用第一性原理从头算法研究了从头算法研究了SiSi和和GeGe的晶体结构,指出常压的晶体结构,指出常压下金刚石结构在能量上最稳定,算出的晶格常下金刚石结构在能量上最稳定,算出的晶格常数与实验值吻合的非常好。数与实验值吻合的非常好。n n他们还从单点能的体积变化率预言了压力导致他们还从单点能的体积变化率预言了压力导致的相变,即当压力增加体积缩小时,其他相的的相变,即当压力增加体积缩小时,其他相的单点能可能比常压下稳定相的单点能更低,从单点能可能比常压下稳定相的单点能更低,从而发生相变。而发生相变。n n计算结果表明,虽然六角金刚石是个次稳定的计算结果表明,虽然六角金刚石是个次稳定的
35、结构,其单点能极小值仅次于金刚石结构,但结构,其单点能极小值仅次于金刚石结构,但在稳定结构附近的单点能体积变化率比金刚石在稳定结构附近的单点能体积变化率比金刚石结构大,而结构大,而-锡结构的单点能体积变化与金锡结构的单点能体积变化与金刚石结构一致,因此,他们预言,金刚石结构刚石结构一致,因此,他们预言,金刚石结构在静压力下将发生向在静压力下将发生向-锡结构而不是六角金锡结构而不是六角金刚石结构相变,这一预言已被试验证实。刚石结构相变,这一预言已被试验证实。晶体的稳定性与其结合能密切相关,结合能定义为晶体的稳定性与其结合能密切相关,结合能定义为 晶体分晶体分解为单个原子所需要的能量。因此结合能越
36、大,晶体结构解为单个原子所需要的能量。因此结合能越大,晶体结构越稳定。越稳定。其中,x,y分别为晶胞中A,B原子的个数;Etot为晶胞的总能量;EatomA,EatomB分别为自由状态下原子A,B的能量。结合能FeCrMoVNbCEsolid(eV)865.732468.931936.861979.071551.63155.83EsolidA,EsolidB分别为固态单质单个原子A,B的能量。形成能Hexagoal and orthorhombic Mo2C形成能设定任务5.2 电子结构分析n n能带结构n n能态密度n n布局分析n n电荷密度差异性能参数能带结构单电子方程的本征值En(k)
37、对每个n是一个对k准连续的、可区分(非简并情况)的函数,称为能带。所有的能带称为能带结构。人们常常关心由原子的价电子形成的能量相对较高的几个能带,其中在绝对零度(0K)时,被价电子填满的能带称为价带,而未被填满或全空的能带称为导带。导带底与价带顶之间的能量区间称为禁带,导带底与价带顶的能量之差称为禁带宽度。能带结构下图是CASTEP计算的闪锌矿结构的BN能带结构图,总能量的计算精度是1.010-5,k的取值为888,能量的截断值为250.0eV。图中虚线表示费米能级(EF=0),费米能级以下的能带是价带,费米能级以上的能带是导带。导带底与价带顶位于相同点的能带结构,称为直接跃迁型能带结构。导带
38、底与价带顶位于不同点的能带结构,称为间接跃迁型能带结构。从图中可以看出,导带底与价带顶位于布里渊区的原点(点,图中以G点表示),这表明BN具有直接型能带结构。能带结构由于三维晶体的波矢k也是三维的,En(k)需要四维空间,因此,一般使波矢k沿选定的方向取值,画出二维的En(k)图。所选定的直线方向一般是晶体倒易点阵的高对称方向,如立方晶体倒易点阵的轴(100方向)、轴(110方向)和轴(111方向)。全部En(k)函数给出能带结构,实际上人们最关心的是费米能级附近的一系列能带,因此,一般只计算有限个能带的En(k)函数。从能带结构图上可以直观地看到在指定方向上各能带En(k)函数随k的变化、导
39、带底与价带顶的位置、禁带宽度及禁带能隙随k的变化。能态密度在原子中电子的本征态形成一系列分立的能级,可以具体标明各能级的能量,说明它们的分布情况。然而,在固体中,电子能级非常密集,形成准连续分布,标明其中每个能级的能量是没有意义的。为了描述固体能带中电子能级的分布,引入“能态密度”的概念。若能量在E(E+E)范围内的电子能态数目为Z,则N(E)称为能态密度(或态密度)。能态密度能态密度描述了电子态的能量分布,由态密度的定义式可推出,第n个能带的态密度Nn(E)为式中,为原胞体积,N为晶体中原胞数目,积分在布里渊区内进行。总态密度N(E)为所有能带的态密度之和,即积分在等能面SE上进行。能态密度
40、总态密度N(E)是各能带的态密度之和,总电子数N等于N(E)从负无穷到费米能级EF的积分,即态密度是一个十分有用的概念,使用态密度可以用对电子能量E的积分代替在布里渊区内对k的积分。态密度经常用于电子结构的快速可视分析,价带宽度、能隙及电子态密度N(E)的主要特征处强度和数目等特性有助于定性解释实验得到的光谱数据。态密度分析还有助于理解电子结构的变化,如外压引起的电子结构的变化。能态密度能态密度(DOS)可分为总态密度、分波态密度(PDOS)和局域态密度(LDOS)。局域态密度和分波态密度是电子结构分析十分有用的半定量工具。LDOS显示系统中原子的电子态对能态密度谱的每个部分的贡献。PDOS根
41、据电子态密度角动量来进一步分辨这些贡献,确定DOS的主要峰是否具有s,p或d电子的特征。LDOS和PDOS分析可对体系中电子杂化的本质和体系的XPS谱、光谱中的主要特征的来源提供定性解释。PDOS计算基于Mulliken布局分析,这种布局分析可将原子对每个能带的贡献归属于指定的原子轨道,对所有能带中某些原子轨道的贡献求和,可得到加权的态密度。可选择不同的加权方式,例如,将指定原子的所有原子轨道的对各能带的贡献加起来便得到LDOS。能态密度结果输出能态密度在自旋极化体系中,电子和电子具有不同的空间波函数,即占据不同的能态。可分别计算电子的态密度N(E)和电子的态密度N(E),它们的和给出总态密度
42、,它们的差N(E)-N(E)被称为自旋态密度(SDOS)。材料的磁性质与自旋态密度有关。Fe是典型的自旋极化体系,用CASTEP计算-Fe的能态密度,总能量的计算精度是1.010-5,k的取值为101010,能量的截断值为300.0eV,分别得到电子和电子的能态密度,如图所示。布局分析对电子电荷在各组分原子之间的分布情况进行计算,称为布局分析。有多种布局分析方法,其中被广泛采用的布局分析方法是Mulliken布局分析。布局分析可以给出原子上、原子轨道上、两原子间的电子电荷分布,依次称为原子布局、轨道布局和键布局。布局分析布局分析为原子间成键提供了一个客观判据,并且两原子间的重叠布局还可用于评价
43、一个键的共价性或离子性。键布局的值高表明键是共价的,键布局的值低表明键是一种离子相互作用。还可以用有效离子价来进一步评价键的离子性,有效离子价定义为阴离子物种上原来的离子电荷与Mulliken电荷之差,若这个值为零,则表明该键是完全离子键,若这个值大于零,则表明该键的共价成分增加。电荷密度差异n n利用能带理论除了计算利用能带理论除了计算E En n(k(k)函数、能态密度等描述能带函数、能态密度等描述能带结构的物理量,还可以对固体的一些物理性质,如弹性常结构的物理量,还可以对固体的一些物理性质,如弹性常数、光学性质、热力学函数等,进行计算。数、光学性质、热力学函数等,进行计算。n n计算固体
44、的物理性质时,应先进行几何优化。计算固体的物理性质时,应先进行几何优化。5.3 弹性性质输出与分析弹性性质材料的弹性常数描述了它对所加应力的相应。或者反过来说,弹性常数描述了为维持一个给定形变所需的应力。应力 和应变 均为二阶对称张量,可分别用 (i=1,2,6)和 (j=1,2,6)来表示,则线弹性常数 可以表示为一个66的对称矩阵对于小的应力应变,有知道晶体总能量ET可算出弹性常数 。利用计算得到的弹性常数 ,还可计算体弹模量、泊松比的性质。待优化结构完成之后VC弹性常数弹性常数5.45.4 其他物理性质分析其他物理性质分析热力学性质对体系热力学性质的描述基于声子,声子是晶格震动的能量子。
45、声子的角频率与波矢q的函数关系(q)称为声子谱或色散关系。利用第一原理计算声子谱(q)的方法有两种:超胞法和线性响应法。由声子谱(q)可计算体系的焓H、熵S自由能F和晶格热容CV,它们都是温度的函数。下图是用CASTEP计算的立方ZrO2的焓(H)、功函数(F)、温度与熵的乘积(TS)和等容热容(CV)随温度的变化。光学性质光学性质6 应用实例晶体结构稳定性对于组成元素确定的体系,可能存在不同的晶体结构。分别对不同的晶体结构的能量极小值进行计算,比较不同晶体结构的能量极小值,便可确定稳定的晶体结构。Orthorhombic and Hexagonal Fe3C(a)(b)MoMoCCFig.1
46、 The crystal structures of Mo2C with 12 atoms basis in(a)Ortho-Mo2C,(b)Hexa-Mo2C.晶体结构稳定性对于相同结构的晶体结构,研究空位或原子置换对晶体结构稳定性的影响。通过对能量极小值进行计算,比较不同处理后的能量极小值,便可确定稳定的晶体结构。B.Xiao等b还利用第一原理和试验的方法研究了W、Mo、B在Cr7C3中多元复合添加的影响,研究结果表明复合添加B和W或Mo有利于Fe16Cr12C12结构的化学稳定性。(Chemical Physics Letters,2008)晶体导电性及结合键的强弱直接跃迁型能带结构间接
47、跃迁型能带结构能隙的宽度结合键的共价性和离子性结合键的强度材料的磁性质与自旋态密度有关。计算晶体结构的物理性质,计算力学稳定性及材料的韧脆性C C1111(GPaGPa)C C12 12(GPaGPa)C C4444(GPaGPa)G(G(GPaGPa)B(B(GPaGPa)V V4 4C C4 4578.22578.22147.17147.17176.27176.27518.50518.50290.85290.85V V3 3FeCFeC4 4638.79638.79135.60135.60130.95130.95592.23592.23301.60301.60V V3 3CrCCrC4 4
48、677.78677.78128.92128.92191.09191.09633.59633.59310.62310.62V V3 3MoCMoC4 4689.40689.40115.43115.43203.90203.90653.73653.73305.57305.57V V3 3NbCNbC4 4643.42643.42129.27129.27150.84150.84599.92599.92301.13301.13NbNb4 4C C4 4721.02721.0294.3494.34188.75188.75699.39699.39303.28303.28NbNb3 3FeCFeC4 4635
49、.02635.02127.90127.90136.46136.46592.88592.88297.41297.41NbNb3 3CrCCrC4 4653.19653.19123.42123.42171.72171.72614.22614.22299.61299.61NbNb3 3MoCMoC4 4724.36724.36111.81111.81183.33183.33696.80696.80318.69318.69NbNb3 3VCVC4 4699.14699.14102.10102.10180.08180.08672.65672.65300.53300.53表6 金属原子M(M=Fe,Cr,
50、Mo)置换后VC与NbC碳化物的三个独立的弹性常数及体弹模量根据判据,以上碳化物机械稳定。发生置换后,VC的体弹模量升高。Nb的碳化物中,Nb3MoC的体弹模量最高。固体的弹性常数与材料的力学性能密切相关,可以用来表征晶体对外界应力的响应,可作为原子间相互作用力大小的判据,并且弹性常数在材料的稳定性、硬度和强度上扮演着重要的角色。G(G(GPaGPa)E(E(GPaGPa)v vA AB/GB/GV V4 4C C4 4191.04191.04470.18470.180.2310.2310.8180.8181.5231.523V V3 3FeCFeC4 4170.10170.10429.554