材料第一性原理操作和试验方法.ppt

第一性原理计算的软件实现及其在材料设计中应用报告人：曹勇 博士生导 师 ：朱景川 教授2011.03.09n n1 概述n n2 计算软件CASTEP的简介及工作环境n n3 模型搭建n n4 几何优化n n5 计算结果分析n n6 应用实例1 概述n n第一原理计算指的是只依赖于基本的物理常数和原子质量，第一原理计算指的是只依赖于基本的物理常数和原子质量，通过量子力学理论求解运动方程来确定电子行为，进而预通过量子力学理论求解运动方程来确定电子行为，进而预测物质性质由于不借助任何经验、半经验理论或试验结测物质性质由于不借助任何经验、半经验理论或试验结果，因此也被称为从头算方法（果，因此也被称为从头算方法（abab initio method) initio method)n n仅需仅需5 5个基本物理常数，即电子的静止质量个基本物理常数，即电子的静止质量mm、电子电量、电子电量e e、、普朗克常数普朗克常数h h、光速、光速c c和波尔玆曼常数和波尔玆曼常数k kB B n n第一原理计算方法是从最基本和最简单的物理事实出发，第一原理计算方法是从最基本和最简单的物理事实出发，计算结果往往能够具有很高的精确度和可信度，成为独立计算结果往往能够具有很高的精确度和可信度，成为独立与试验之外研究物质性质的又一利器。

与试验之外研究物质性质的又一利器 1第一原理的可靠性Total Energy, 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 o9898o 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 PropertiesLatticeconstantMagnetism: 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.72/atomCalculated:Calculated:m mB B=1.63/atom=1.63/atomMagnetism: 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/atomMagnetic PropertiesMetalMetalCrystalCrystalstructurestructureatomic 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.65n n由于在处理多电子体系相互作用时遇到的困难由于在处理多电子体系相互作用时遇到的困难, , 第一原理第一原理在实际计算中也需要做一定程度的近似，不过这样的近似在实际计算中也需要做一定程度的近似，不过这样的近似也都是基于物理上的简单假设和数学形式上的合理简化而也都是基于物理上的简单假设和数学形式上的合理简化而得到的。

得到的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) of their constituent atoms?CSiGeNB•Structural properties?•Electronic properties?•Elastic properties?•optical properties ?•Phonon dispersion?•Magnetic properties?•Superconductivity?2第一原理所能解决的问题Hexagoal Mo2CMo2C 是合金钢中的二次硬化相，它均匀、细小、弥散的分布于基体中，具有极高的硬度和良好的高温稳定性①合金钢中合金元素种类繁多，是否存在其他的碳化物相，如Cr2C、V2C、Ni2C来取代其作用。

② Mo2C 是否具有其他更加稳定的结构③ 空位或其他原子置换对其稳定性和性能有何影响2 第一原理所能解决的问题Orthorhombic and Hexagonal Fe3C晶体结构的稳定性判定虚拟相、亚稳相、稳定相研究晶格常数和原子所处位置对晶体结构稳定性的影响对于组成元素确定的体系，可能存在不同的晶体结构分别对不同的晶体结构的能量极小值进行计算，比较不同晶体结构的能量极小值，便可确定稳定的晶体结构2 第一原理所能解决的问题Material Properties from First-PrinciplesFrom first principles!Predict new behaviors/properties of existing materialsDesign materials with desired propertiesUnderstand and explain materials properties√√Becoming realityn n目前已有一些利用第一原理对具有周期性结构的材料进行目前已有一些利用第一原理对具有周期性结构的材料进行能带计算的商业软件，如能带计算的商业软件，如VASPVASP，，ABINITABINIT，，Wien2KWien2K，，CASTEPCASTEP。

2 计算软件CASTEP的简介及运行环境n n计算软件计算软件CASTEP (Cambridge Sequential Total Energy CASTEP (Cambridge Sequential Total Energy Package) Package) 是基于密度泛函理论的第一性原理计算量子力是基于密度泛函理论的第一性原理计算量子力学程序其使用的是总能量平面波赝势的方法其使用的是总能量平面波赝势的方法n n密度泛函理论密度泛函理论 (Density Function Theory) (Density Function Theory) 巧妙的将电子巧妙的将电子之间的交换势表示为密度泛函，从而使得薛定谔方程在考之间的交换势表示为密度泛函，从而使得薛定谔方程在考虑了电子之间的复杂相互作用之后，依然可以用自洽的方虑了电子之间的复杂相互作用之后，依然可以用自洽的方法求解关于关于CASTAPCASTAP CASTAPCASTAP是特别为固体材料学而设计的一个现代的量子力学基本程序，是特别为固体材料学而设计的一个现代的量子力学基本程序，其使用了密度泛函其使用了密度泛函(DFT(DFT））平面波赝势方法，进行第一原理量子力学计算，平面波赝势方法，进行第一原理量子力学计算，以探索如半导体，陶瓷，金属，矿物和沸石等材料的晶体和表面性质。

以探索如半导体，陶瓷，金属，矿物和沸石等材料的晶体和表面性质 典型的应用包括表面化学，键结构，态密度和光学性质等研究，典型的应用包括表面化学，键结构，态密度和光学性质等研究， CASTAPCASTAP也可用于研究体系的电荷密度和波函数的也可用于研究体系的电荷密度和波函数的3D3D形式此外，形式此外， CASTAPCASTAP可可用于有效研究点缺陷（空位，间隙和置换杂质）和扩展缺陷（如晶界和位用于有效研究点缺陷（空位，间隙和置换杂质）和扩展缺陷（如晶界和位错）的性质错）的性质 Material StudioMaterial Studio使用组件对话框中的使用组件对话框中的CASTAPCASTAP选项允许准备，启动，选项允许准备，启动，分析和监测分析和监测CASTAPCASTAP服役工作服役工作计算：允许选择计算选项（如基集，交换关联势和收敛判据），作业控制计算：允许选择计算选项（如基集，交换关联势和收敛判据），作业控制和文档控制和文档控制分析：允许处理和演示分析：允许处理和演示CASTAPCASTAP计算结果这一工具提供加速整体直观化以计算结果这一工具提供加速整体直观化以及键结构图，态密度图形和光学性质图形。

及键结构图，态密度图形和光学性质图形CASTAPCASTAP计算计算 CASTAPCASTAP计算是要进行的三个任务中的一个，即单个点的能量计算，计算是要进行的三个任务中的一个，即单个点的能量计算，几何优化或分子动力学可提供这些计算中的每一个以便产生特定的几何优化或分子动力学可提供这些计算中的每一个以便产生特定的物理性能性质为一种附加的任务，允许重新开始已完成的计算以便物理性能性质为一种附加的任务，允许重新开始已完成的计算以便产生最初没有提出的额外性能产生最初没有提出的额外性能在在CASTAPCASTAP计算中有很多运行步骤，可分为如下几组：计算中有很多运行步骤，可分为如下几组：* * 结构定义：必须规定包含所感兴趣结构的周期性的结构定义：必须规定包含所感兴趣结构的周期性的3D3D模型文件，有模型文件，有大量方法规定一种结构：可使用构建晶体（大量方法规定一种结构：可使用构建晶体（Build Crystal)Build Crystal)或构建真或构建真空板空板(Build Vacuum Stab)(Build Vacuum Stab)来构建，也可从已经存在的的结构文档中引来构建，也可从已经存在的的结构文档中引入，还可修正已存在的结构。

入，还可修正已存在的结构提示：提示： CASTAPCASTAP计算所需时间随原子数平方的增加而增加因此，建议计算所需时间随原子数平方的增加而增加因此，建议是用最小的初晶胞来描述体系，可使用是用最小的初晶胞来描述体系，可使用Build\Symmetry\Primitive Build\Symmetry\Primitive CellCell菜单选项来转换成初晶胞菜单选项来转换成初晶胞 CASTAPCASTAP中中选择一项任务选择一项任务1 1 从模块面板（从模块面板（Module Explorer)Module Explorer)选择选择CASTAP\CalculationCASTAP\Calculation2 2 选择设置表选择设置表3 3 从任务列表中选择所要求的任务从任务列表中选择所要求的任务 •计算设置：合适的计算设置：合适的3D3D模型文件一旦确定，必须选择计算类型模型文件一旦确定，必须选择计算类型和相关参数，例如，对于动力学计算必须确定系综和参数，包和相关参数，例如，对于动力学计算必须确定系综和参数，包括温度，时间步长和步数选择运行计算的磁盘并开始括温度，时间步长和步数。

选择运行计算的磁盘并开始CASTAPCASTAP作业 * 结果分析：计算完成后，相关于结果分析：计算完成后，相关于CASTAPCASTAP作业的文档返回用作业的文档返回用户，在项目面板适当位置显示这些文档的一些进一步处理要户，在项目面板适当位置显示这些文档的一些进一步处理要求获得可观察量如电子结构图、光学性质求获得可观察量如电子结构图、光学性质CASTAPCASTAP能量任务能量任务CASTAPCASTAP能量任务允许计算特定体系的总能量以及物理性质能量任务允许计算特定体系的总能量以及物理性质除了总能量之外，在计算之后还可报告作用于原子上的力；也除了总能量之外，在计算之后还可报告作用于原子上的力；也能创建电荷密度文件；利用材料观测仪（能创建电荷密度文件；利用材料观测仪（Material Material VisualizerVisualizer) )允许目测电荷密度的立体分布；还能报告计算中使允许目测电荷密度的立体分布；还能报告计算中使用的用的Monkhorst-Park的的k点的电子能量，因此在点的电子能量，因此在CASTAPCASTAP分析中分析中可生成态密度图可生成态密度图。

对于能够得到可靠结构信息的体系的电子性质的研究，能量任对于能够得到可靠结构信息的体系的电子性质的研究，能量任务是有用的只要给定应力性质，也可用于计算没有内部自由务是有用的只要给定应力性质，也可用于计算没有内部自由度的高对称性体系的状态方程（即压力度的高对称性体系的状态方程（即压力- -体积，能量体积，能量- -体积关系）体积关系）CASTAPCASTAP中能量的缺损单位是电子伏特中能量的缺损单位是电子伏特( (eVeV) )，各种能量单位的换，各种能量单位的换算关系见算关系见Mohr.P.J(2000).Mohr.P.J(2000).1 eV=0.036749308 Ha=23.0605 kcal/mole=96.4853 kJ/moleCASTAPCASTAP几何优化任务几何优化任务 CASTAPCASTAP几何优化任务允许改善结构的几何，获得稳定结构几何优化任务允许改善结构的几何，获得稳定结构或多晶型物通过一个迭代过程来完成这项任务，迭代过程中调或多晶型物通过一个迭代过程来完成这项任务，迭代过程中调整原子坐标和晶胞参数使结构的总能量最小化整原子坐标和晶胞参数使结构的总能量最小化。

CASTAPCASTAP几何优化是基于减小计算力和应力的数量级，直到小几何优化是基于减小计算力和应力的数量级，直到小于规定的收敛误差也可能给定外部应力张量来对拉应力、压应于规定的收敛误差也可能给定外部应力张量来对拉应力、压应力和切应力等作用下的体系行为模型化在这些情况下力和切应力等作用下的体系行为模型化在这些情况下反反复迭复迭代内部应力张量直到代内部应力张量直到与所施加的外部应力相等与所施加的外部应力相等 几何优化处理产生的模型几何优化处理产生的模型结构与真实结构紧密相似结构与真实结构紧密相似利用利用CASTAPCASTAP计算的晶格参计算的晶格参数精度列于右图数精度列于右图CASTAPCASTAP性质任务性质任务 CASTAPCASTAP性质任务允许在完成能量，几何优化或动力学运行之后求出电性质任务允许在完成能量，几何优化或动力学运行之后求出电子和结构性质可以产生的性质如下：子和结构性质可以产生的性质如下：* * 态密度（态密度（DOSDOS）：利用原始模拟中产生的电荷密度和势能，非自恰计算价）：利用原始模拟中产生的电荷密度和势能，非自恰计算价带和导带的精细带和导带的精细MonkhorstMonkhorst-Pack -Pack 网格上的电子本征值。

网格上的电子本征值 * 带结构：利用原始模拟中产生的电荷密度和势能，非自恰计算价带和导带结构：利用原始模拟中产生的电荷密度和势能，非自恰计算价带和导带的布里渊区高对称性方向电子本征值带的布里渊区高对称性方向电子本征值 * 光学性质：计算电子能带间转变的矩阵元素光学性质：计算电子能带间转变的矩阵元素CASTAPCASTAP分析对话可用于生分析对话可用于生成包含可以测得的光学性质的网格和图形文件成包含可以测得的光学性质的网格和图形文件 * 布局数分析：进行布局数分析：进行MullikenMulliken 分析计算决定原子电荷的键总数和角动量分析计算决定原子电荷的键总数和角动量（以及自旋极化计算所需的磁矩）任旋地，可产生态密度微分计算所要（以及自旋极化计算所需的磁矩）任旋地，可产生态密度微分计算所要求的分量求的分量 * 应力：计算应力张量，并写入应力：计算应力张量，并写入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)； Si的原子位置为(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) cell(unit) celln n惯用晶胞＝惯用晶胞＝Conventional (unit) cellConventional (unit) celln n由一组基矢所决定的平行六面体所围起来的最小重复单元由一组基矢所决定的平行六面体所围起来的最小重复单元就叫就叫原胞原胞原胞原胞（或初基单胞（或初基单胞Primitive CellPrimitive Cell (固体物理常用固体物理常用) )n n原胞代表的是几何形状、体积的最小重复单元，不具有物原胞代表的是几何形状、体积的最小重复单元，不具有物理内涵n n原胞往往不能反映晶体的对称性，因而，习惯上常选择能原胞往往不能反映晶体的对称性，因而，习惯上常选择能反映晶体对称性的重复单元，这种重复单元就叫惯用反映晶体对称性的重复单元，这种重复单元就叫惯用晶胞晶胞晶胞晶胞（（conventional cellconventional cell）或非）或非初基单胞初基单胞初基单胞初基单胞（（nonprimitivenonprimitive cell) cell)。

n n晶胞一般不是最小的重复单元其体积（面积）可以是原晶胞一般不是最小的重复单元其体积（面积）可以是原胞的数倍胞的数倍 寻找已经存在的常见晶体模型寻找已经存在的常见晶体模型寻找已经存在的常见晶体模型寻找已经存在的常见晶体模型n n建立所研究材料的结构模型建立所研究材料的结构模型 ①①空间群空间群 ②②晶格参数晶格参数 ③③晶体的原子内坐标晶体的原子内坐标n n以以VCVC为例，为例， ①①其空间群为其空间群为FM-3MFM-3M；； ②②晶格常数为晶格常数为a=4.160Åa=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.8.4. 几何优化能带计算的过程晶体是一个具有周期性结构的体系，输入时只能给出一个体积有限的晶体结构模型，需要利用周期性边界条件，才能得到整个晶体的能带结构在选择了计算任务之后，需要对计算参数进行设置，如自洽场计算的精度、基组大小、k的取值等在能带计算中，基组是指用于线性组合成单电子波函数的基函数集合，在平面波赝势方法中，给出能量的截断值Ecut由于计算量的的原因计算时k值只取简约布里渊区中的有限值，以了l1× l2 ×l3的形式表示所取的k在b1、b2和b3方向上取值间隔分别为b1 / l1 , b2 / l2 , b3 / l3 , 在其他条件不变的情况下，k的取值数目增加，得到的En(k) 函数的精确度增加，但计算量显著增加自洽场计算的收敛用电子数密度或晶体总能量的收敛来标志，自洽场计算的精度是指收敛的标准能带计算的过程 利用自洽场方法求解Kohn-Sham方程，得到所设结构的晶体总能量在计算晶体的物理性质之前必须对所设的晶体结构模型进行几何优化，根据关于能量、力、应力、位移的判据来判断晶体结构是否为稳定结构（总能量最小）。

如果晶体结构不是稳定结构，重新设置晶格参数进行计算，直至得到稳定的晶体结构然后对结构优化后的晶体进行物理性质计算，输出计算结果能带计算的过程几何优化几何优化是通过调节结构模型的几何参数来获得稳定结构的过程，其结果是使模型结构尽可能的接近真实结构进行几何优化的判据可以根据研究的需要而定，一般是几个判据的组合使用常用的判据有以下几个：①自洽场收敛判据对给定的结构模型进行自洽场计算时，相继两次自洽计算得到的晶体总能量之差足够小，即相继两次自洽计算的晶体总能量之差小于设定值②力判据每个原子所受的晶体内作用力足够小，即单个原子受力小于设定最大值③应力判据每个结构模型单元中的应力足够小，即应力小于设定的最大值④位移判据相继两次结构参数变化引起的原子位移的分量足够小，即原子位移分量小于设定的最大值设置参数设置参数设置参数设置参数设置电子结构参数几何优化几何优化(a)(b)MoMoCCFig. 1 The crystal structures of Mo2C with 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成功地利用第一性原理成功地利用第一性原理从头算法研究了从头算法研究了SiSi和和GeGe的晶体结构，指出常压的晶体结构，指出常压下金刚石结构在能量上最稳定，算出的晶格常下金刚石结构在能量上最稳定，算出的晶格常数与实验值吻合的非常好数与实验值吻合的非常好n n他们还从单点能的体积变化率预言了压力导致他们还从单点能的体积变化率预言了压力导致的相变，即当压力增加体积缩小时，其他相的的相变，即当压力增加体积缩小时，其他相的单点能可能比常压下稳定相的单点能更低，从单点能可能比常压下稳定相的单点能更低，从而发生相变而发生相变n n计算结果表明，虽然六角金刚石是个次稳定的计算结果表明，虽然六角金刚石是个次稳定的结构，其单点能极小值仅次于金刚石结构，但结构，其单点能极小值仅次于金刚石结构，但在稳定结构附近的单点能体积变化率比金刚石在稳定结构附近的单点能体积变化率比金刚石结构大，而结构大，而ββ- -锡结构的单点能体积变化与金锡结构的单点能体积变化与金刚石结构一致，因此，他们预言，金刚石结构刚石结构一致，因此，他们预言，金刚石结构在静压力下将发生向在静压力下将发生向ββ- -锡结构而不是六角金锡结构而不是六角金刚石结构相变，这一预言已被试验证实。

刚石结构相变，这一预言已被试验证实 晶体的稳定性与其结合能密切相关，结合能定义为晶体的稳定性与其结合能密切相关，结合能定义为 晶体分晶体分解为单个原子所需要的能量因此结合能越大，晶体结构解为单个原子所需要的能量因此结合能越大，晶体结构越稳定其中，x，y分别为晶胞中A，B原子的个数；Etot为晶胞的总能量；EatomA,EatomB分别为自由状态下原子A，B的能量结合能FeCrMoVNbCEsolid （eV）−865.73−2468.93−1936.86−1979.07−1551.63−155.83EsolidA,EsolidB分别为固态单质单个原子A，B的能量形成能Hexagoal and orthorhombic Mo2C形成能设定任务5.2 电子结构分析n n能带结构n n能态密度n n布局分析n n电荷密度差异性能参数能带结构•单电子方程的本征值En(k)对每个n是一个对k准连续的、可区分（非简并情况）的函数，称为能带所有的能带称为能带结构•人们常常关心由原子的价电子形成的能量相对较高的几个能带，其中在绝对零度(0K)时，被价电子填满的能带称为价带，而未被填满或全空的能带称为导带。

导带底与价带顶之间的能量区间称为禁带，导带底与价带顶的能量之差称为禁带宽度能带结构•下图是CASTEP计算的闪锌矿结构的BN能带结构图，总能量的计算精度是1.0×10-5，k的取值为8×8×8，能量的截断值为250.0eV•图中虚线表示费米能级(EF=0),费米能级以下的能带是价带，费米能级以上的能带是导带导带底与价带顶位于相同点的能带结构，称为直接跃迁型能带结构导带底与价带顶位于不同点的能带结构，称为间接跃迁型能带结构•从图中可以看出，导带底与价带顶位于布里渊区的原点（点，图中以G点表示），这表明BN具有直接型能带结构能带结构•由于三维晶体的波矢k也是三维的，En(k)需要四维空间，因此，一般使波矢k沿选定的方向取值，画出二维的En(k)图所选定的直线方向一般是晶体倒易点阵的高对称方向，如立方晶体倒易点阵的⊿轴（＜100＞方向）、∑轴（＜110＞方向）和⊿轴（＜111＞方向）•全部En(k)函数给出能带结构，实际上人们最关心的是费米能级附近的一系列能带，因此，一般只计算有限个能带的En(k)函数从能带结构图上可以直观地看到在指定方向上各能带En(k)函数随k的变化、导带底与价带顶的位置、禁带宽度及禁带能隙随k的变化。

能态密度在原子中电子的本征态形成一系列分立的能级，可以具体标明各能级的能量，说明它们的分布情况然而，在固体中，电子能级非常密集，形成准连续分布，标明其中每个能级的能量是没有意义的为了描述固体能带中电子能级的分布，引入“能态密度”的概念若能量在E(E+ΔE)范围内的电子能态数目为ΔZ，则N(E)称为能态密度（或态密度）能态密度能态密度描述了电子态的能量分布，由态密度的定义式可推出，第n个能带的态密度Nn(E)为式中，为原胞体积，N为晶体中原胞数目，积分在布里渊区内进行总态密度N(E)为所有能带的态密度之和，即积分在等能面SE上进行能态密度总态密度N(E)是各能带的态密度之和，总电子数N等于N(E)从负无穷到费米能级EF的积分，即•态密度是一个十分有用的概念，使用态密度可以用对电子能量E的积分代替在布里渊区内对k的积分•态密度经常用于电子结构的快速可视分析，价带宽度、能隙及电子态密度N(E)的主要特征处强度和数目等特性有助于定性解释实验得到的光谱数据态密度分析还有助于理解电子结构的变化，如外压引起的电子结构的变化能态密度•能态密度（DOS）可分为总态密度、分波态密度(PDOS)和局域态密度(LDOS)。

局域态密度和分波态密度是电子结构分析十分有用的半定量工具•LDOS显示系统中原子的电子态对能态密度谱的每个部分的贡献PDOS根据电子态密度角动量来进一步分辨这些贡献，确定DOS的主要峰是否具有s, p或d电子的特征•LDOS和PDOS分析可对体系中电子杂化的本质和体系的XPS谱、光谱中的主要特征的来源提供定性解释•PDOS计算基于Mulliken布局分析，这种布局分析可将原子对每个能带的贡献归属于指定的原子轨道，对所有能带中某些原子轨道的贡献求和，可得到加权的态密度可选择不同的加权方式，例如，将指定原子的所有原子轨道的对各能带的贡献加起来便得到LDOS能态密度结果输出能态密度•在自旋极化体系中，α电子和β电子具有不同的空间波函数，即占据不同的能态可分别计算α电子的态密度N(E)↑和β电子的态密度N(E) ↓，它们的和给出总态密度，它们的差N(E)↑- N(E) ↓被称为自旋态密度（SDOS）•材料的磁性质与自旋态密度有关Fe是典型的自旋极化体系，用CASTEP计算α-Fe的能态密度，总能量的计算精度是1.0×10-5，k的取值为10×10×10，能量的截断值为300.0eV，分别得到α电子和β电子的能态密度，如图所示。

布局分析•对电子电荷在各组分原子之间的分布情况进行计算，称为布局分析•有多种布局分析方法，其中被广泛采用的布局分析方法是Mulliken布局分析•布局分析可以给出原子上、原子轨道上、两原子间的电子电荷分布，依次称为原子布局、轨道布局和键布局布局分析•布局分析为原子间成键提供了一个客观判据，并且两原子间的重叠布局还可用于评价一个键的共价性或离子性•键布局的值高表明键是共价的，键布局的值低表明键是一种离子相互作用还可以用有效离子价来进一步评价键的离子性，有效离子价定义为阴离子物种上原来的离子电荷与Mulliken电荷之差，若这个值为零，则表明该键是完全离子键，若这个值大于零，则表明该键的共价成分增加电荷密度差异n n利用能带理论除了计算利用能带理论除了计算E En n(k(k) )函数、能态密度等描述能带函数、能态密度等描述能带结构的物理量，还可以对固体的一些物理性质，如弹性常结构的物理量，还可以对固体的一些物理性质，如弹性常数、光学性质、热力学函数等，进行计算数、光学性质、热力学函数等，进行计算n n计算固体的物理性质时，应先进行几何优化计算固体的物理性质时，应先进行几何优化5.3 弹性性质输出与分析弹性性质•材料的弹性常数描述了它对所加应力的相应。

或者反过来说，弹性常数描述了为维持一个给定形变所需的应力应力 和应变 均为二阶对称张量，可分别用 (i=1, 2,…6)和 (j=1, 2,…6)来表示，则线弹性常数 可以表示为一个6×6的对称矩阵对于小的应力应变，有知道晶体总能量ET可算出弹性常数 利用计算得到的弹性常数 ，还可计算体弹模量、泊松比的性质待优化结构完成之后VC弹性常数弹性常数5.45.4 其他物理性质分析其他物理性质分析热力学性质•对体系热力学性质的描述基于声子，声子是晶格震动的能量子•声子的角频率与波矢q的函数关系ω(q)称为声子谱或色散关系•利用第一原理计算声子谱ω(q)的方法有两种：超胞法和线性响应法由声子谱ω(q)可计算体系的焓H、熵S自由能F和晶格热容CV，它们都是温度的函数•下图是用CASTEP计算的立方ZrO2的焓(H)、功函数(F)、温度与熵的乘积(TS)和等容热容(CV)随温度的变化光学性质光学性质6 应用实例①晶体结构稳定性对于组成元素确定的体系，可能存在不同的晶体结构分别对不同的晶体结构的能量极小值进行计算，比较不同晶体结构的能量极小值，便可确定稳定的晶体结构。

Orthorhombic and Hexagonal Fe3C(a)(b)MoMoCCFig. 1 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)②晶体导电性及结合键的强弱直接跃迁型能带结构间接跃迁型能带结构能隙的宽度结合键的共价性和离子性结合键的强度材料的磁性质与自旋态密度有关③计算晶体结构的物理性质，计算力学稳定性及材料的韧脆性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 4677.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.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, 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.55429.550.2630.2630.5210.5211.7731.773V V3 3CrCCrC4 4220.84220.84535.59535.590.2130.2130.6960.6961.4071.407V V3 3MoCMoC4 4233.73233.73558.73558.730.1950.1950.7110.7111.3071.307V V3 3NbCNbC4 4186.70186.70464.17464.170.2430.2430.5870.5871.6131.613NbNb4 4C C4 4231.04231.04552.76552.760.1960.1960.6020.6021.3131.313NbNb3 3FeCFeC4 4174.87174.87438.64438.640.2540.2540.5380.5381.7011.701NbNb3 3CrCCrC4 4204.21204.21499.21499.210.2220.2220.6480.6481.4671.467NbNb3 3MoCMoC4 4225.03225.03546.47546.470.2140.2140.5990.5991.4161.416NbNb3 3VCVC4 4220.33220.33531.18531.180.2050.2050.6030.6031.3641.364表7 金属原子M (M=Fe, Cr, Mo)置换后VC与NbC碳化物的剪切模量、杨氏模量、各向异性因子及B/G比根据参考文献，B/G大于1.75，材料表现为韧性；反之则表现为脆性。

由以上可知除V3FeC4，其余碳化物均表现为脆性④作为中间输出量，为其它计算提供参数计算的体系能量，可为热力学计算提供参数；弹性常数可作为更高尺度材料计算模拟的输入参数建立经验公式，计算硬度n n硬度 过去二十年，人们通过可从第一性原理直接计算得出的弹性常数来寻找超硬材料，事实上，硬度与弹性常数并无一一对应的关系 但是，建立基于第一性原理输出值而无需试验输入的硬度计算公式来寻找超硬材料一直是人们努力的方向 建立经验公式，计算硬度n n硬度 张发明等提出共价化合物的硬度与两个因素有关 ①键的强度②单位面积内键的个数 其中，A是常数，取740；P为布局数；Vb为键的体积Table 1 The overlap population and the hardness of VC and V4C3, where Hvcalc, Hvother and Hvexp are the calculated hardness, the calculated hardness by other person [5] and the experimental Vickers hardness[5], respectively.BondsBondsPopulationPopulation( Å( Å3 3) )Bond number in Bond number in one cellone cellBond Bond length (Å)length (Å)H Hvcalcvcalc( (GPaGPa) )H Hvothervother( (GPaGPa) )H Hvexpvexp( (GPaGPa) )VCVCV-CV-C0.740.745.765.7612122.0522.05229.629.627.227.22929TiCTiCTi-CTi-C0.790.796.516.5112122.1372.13725.725.719.119.124.724.7NbCNbCNbNb-C-C0.730.737.247.2412122.2142.21419.919.918.618.61818(a/2) [110] Edge Dislocation in SiXTEM (M. Chisholm)ball-n-stick model-[110] projection⑤模拟位错成核Model System:a/2 <110> edge dislocation[110] projection of diamond Si lattice-Dislocation Core Formation:iif局限性n n采用第一性原理方法往往又受到计算能力的限制不能处理宏观体系。

采用第一性原理方法往往又受到计算能力的限制不能处理宏观体系 “ “由于合金极为复杂，最好走理论与实验相结合的经验理论道路，单由于合金极为复杂，最好走理论与实验相结合的经验理论道路，单纯从第一原理出发可能难以解决问题纯从第一原理出发可能难以解决问题 ” ( (谢佑卿谢佑卿. . 金属材料系统科金属材料系统科学学. . 长沙：中南工业大学出版社，长沙：中南工业大学出版社，19981998））n n 从第一性原理计算的角度看，仅有几个有限的力学量，如晶体的理想从第一性原理计算的角度看，仅有几个有限的力学量，如晶体的理想强度、弹性常数等，可用第一性原理直接计算对于其他的一些宏观强度、弹性常数等，可用第一性原理直接计算对于其他的一些宏观力学性能，如韧性、断裂强度等，虽然本质上仍然由电子结构决定，力学性能，如韧性、断裂强度等，虽然本质上仍然由电子结构决定，但其与电子结构的关联并不明显，很难用第一性原理方法直接计算但其与电子结构的关联并不明显，很难用第一性原理方法直接计算除了从材料的成健类型大致判断力学性能外（如共价键材料的韧性比除了从材料的成健类型大致判断力学性能外（如共价键材料的韧性比金属健材料差、强度高），第一性原理研究材料的力学性能需要采取金属健材料差、强度高），第一性原理研究材料的力学性能需要采取一些间接的方法。

一些间接的方法非常高兴与大家分享非常高兴与大家分享谢谢大家关注！。