CASTEP计算理论总结+实例分析

1、.CASTEP计算理论总结XBAPRSCASTEP特点是适合于计算周期性结构，对于非周期性结构一般要将特定的部分作为周期性结构，建立单位晶胞后方可进行计算。 CASTEP计算步骤可以概括为三步：首先建立周期性的目标物质的晶体；其次对建立的结构进行优化，这包括体系电子能量的最小化和几何结构稳定化。最后是计算要求的性质，如电子密度分布(Electrondensitydistribution)，能带结构 (Band structure)、状态密度分布(Densityof states)、声子能谱(Phonon spectrum)、声子状态密度分布(DOS of phonon)，轨道群分布(Orbit

2、alpopulations)以及光学性质(Optical properties)等。本文主要将就各个步骤中的计算原理进行阐述，并结合作者对计算实践经验，在文章最后给出了几个计算事例，以备参考。CASTEP计算总体上是基于DFT，但实现运算具体理论有：离子实与价电子之间相互作用采用赝势来表示；超晶胞的周期性边界条件；平面波基组描述体系电子波函数；广泛采用快速 fast Fourier transform (FFT) 对体系哈密顿量进行数值化计算；体系电子自恰能量最小化采用迭带计算的方式；采用最普遍使用的交换- 相关泛函实现DFT的计算，泛函含概了精确形式和屏蔽形式。一，CASTEP中周期性结构计

3、算优点与 MS中其他计算包不同，非周期性结构在 CASTEP中不能进行计算。将晶面或非周期性结构置于一个有限长度空间方盒中，按照周期性结构来处理，周期性空间方盒形状没有限制。之所以采用周期性结构原因在于： 依据 Bloch 定理，周期性结构中每个电子波函数可以表示为一个波函数与晶体周期部分乘积的形式。他们可以用以晶体倒易点阵矢量为波矢一系列分离平面波函数来展开。这样每个电子波函数就是平面波和，但最主要的是可以极大简化 Kohn-Sham 方程。这样动能是对角化的，与各种势函数可以表示为相应 Fourier 形式。2 k GV(GG)V(G G)V (GG )CCGGionHxci,k Gi i

4、,k G采用周期性结构的另一个优点是可以方便计算出原子位移引起的整体能量的变化，在CASTEP中引入外力或压强进行计算是很方便的，可以有效实施几何结构优化和分子动力学的模拟。平面波基组可以直接达到有效的收敛。计算采用超晶胞结构的一个缺点是对于某些有单点限缺陷结构建立模型时，体系中的单个缺陷将以无限缺陷阵列形式出现，因此在建立人为缺陷时，它们之间的相互距离应该足够的远，避免缺陷之间相互作用影响计算结果。在计算表面结构时，切片模型应当足够的薄，减小切片间的人为相互作用。CASTEP中采用的交换 - 相关泛函有局域密度近似（ LDA）（ LDA）、广义梯度近似（ GGA）和非定域交换 - 相关泛函。

5、 CASTEP中提供的唯一定域泛函是 CA-PZ， Perdew and Zunger 将 Ceperley and Alder 数值化结果进行了参数拟和。交换 -相关泛函的定域表示形式是目前较为准确的一种描述。NameDescriptionReferencePW91Perdew-Wang generalized-gradient approximation, PW91Perdew and WangPBEPerdew-Burke-Ernzerhof functional, PBEPerdew et al.RPBERevised Perdew-Burke-Ernzerhof functional

6、, RPBEHammer et al.采用梯度校正的非定域或广义梯度近似泛函与电子密度梯度d和电子密度都有关，这样可以dr同时提高能量和结构预测的准确性，但计算耗时。 CASTEP 中提供的非定域泛函有三种：PBE 泛函与 PW91泛函计算在本质上实际是相同的，但在电子密度变化迅速体系中PBE 泛函实用性更好； RPBE 是特别用来提高 DFT 描述金属表面吸附分子能量的泛函，White and Bird 描述了各种梯度校正泛函计算方法，利用广义梯度近似计算总能量使用平面波基组与定域泛函相比并不直接。包含梯度近似的交换-相关泛函计算时对电子密度数据的精度要求较高，对计算机内存占用会增大。通过采

7、用与平面波基组总能量计算中分裂交换 -相关能量采用一系列空间网格相一致的方法来定义交换-相关势。平面波基组（Plane wave basis set ）Bloch 理论表明每个k 点处电子波函数都可以展开成离散的平面波基组形式，理论上讲这种展开形式包含的平面波数量是无限多的。然而相对于动能较大的情况，动能| k+G| 2 很小时平面波系数Ck+G更重要。调节平面波基组，其中包含的平面波动能小于某个设定的截止能量，如图所示（球体半径与截止能量平方根成比例）：总能量计算会因为平面波特定能量截止而产生误差， 通过增加体系能量截止数值就可以减小误差幅度。理论上截止能量必须提高到总能量计算结果达到设定的

8、精确度为止，如果你在进行关于相稳定性的研究，而需要对比每个相能量的绝对值时，这是一种推荐计算方法。不过， 同一个结构在低的截止能量下收敛引起的差别要小于总体能量本身。 因此可以选用合适的平面波基组对几何结构进行优化或进行分子动力学研究。以上的方法对 Brillouin区取样收敛测试同样成立。有限平面波基组的校正采用平面波基组的一个问题是截止能量与基组数量的变化是间断的，一般而言在点基组（k-pointset ）中不同点对应不同能量截止（cutoffs ）时就会产生这种不连续性。此外，在截止能量不变时，晶胞形状和尺寸的变化都会引起平面波基组的间断。通过采用更加致密的点基组就可以解决这个问题，与特

9、定平面波基组相关的加权性也会消除。然而即使在基组取样很致密的情况下，这个问题依然存在，对其近似的解决方法就是引入一个校正因子（ correction factor ），利用某个状态基组计算使用了无限数量的点与实际采用的数量之间的差别来确定。晶体结构在进行几何优化时如果基组不能真正的达到绝对的收敛，有限基组的纠正就很重要。比如硅的规范保守赝势很“软”，在平面波基组截止能量是时就已经可以得到准确的计算结果了。但如果计算状态方程时使用上述截止能量（比如体积与总能量和压强都有关系），能量最小时对应的体积与体系内压为零时对应体积是不同的。在提高截止能量和增加点取样基础上重复对状态方程的计算，这两个体积之

10、间的差别会越来越小。此外截止能量.低时计算得到的E-V 曲线呈现锯齿状，提高截止能量计算的曲线连续而平滑。E-V 曲线中出现锯齿状的原因在于平面波基组在相同的截止能量时由于晶体点阵常数不同引起的平面波基组数量的间断。对总能量进行有限基组的校正，使得我们可以在一个恒定数量基组状态下进行计算，即使采用了恒定的截止能量这个更强制条件也可以纠正计算结果。Milman 等详细的讨论了这种计算方法的细节。进行这种校正所需要的唯一的参数就是dEtot /d lnEcut ，Etot是体系总能量，Ecut 是截止能量。dEtot /d lnEcut 的值很好的表示了能量截止和点取样计算收敛性质。当它的数值（每

11、个原子）小于0.01 eV/atom时，计算就达到了良好的收敛精度，对于大多数计算0.1 eV/atom就足够了。非定域交换 -相关泛函基于 LDA 或 GGA 的泛函的 Kohn-Sham 方程在计算能带带隙上存在低估。这对晶体或分子相关性质以及能量的描述是没有影响的。然而要理解半导体和绝缘体性质，就必须得到关于电子能带结构的准确的描述。 DFT 能带带隙计算误差可以通过引入经验“剪刀”校正，相对于价带而言导带产生了一个刚性的变化。当实验提供的能带带隙准确时，光学性质计算得到了较为准确的结果。电子结构实验数据缺乏时采用“剪刀”工具进行预测性研究或对能带带隙调整是不可靠的。关于DFT 计算中能

12、带带隙问题已经发展许多技术，但这些技术大多复杂而且很耗时，实际计算中最常用的是屏蔽交换（Sx-LDA ），建立在广义 Kohn-Sham 方法基础上。 广义 Kohn-Sham 泛函允许我们将总能量交换分布泛函分离为非定域、定域以及屏蔽密度组元。在CASTEP 计算中采用的广义Kohn-Sham 方法有：? HF: exact exchange, no correlation? HF-LDA: exact exchange, LDA correlation? sX: screened exchange, no correlation? sX-LDA: screened exchange, LD

13、A correlation与 LDA 和 GGA 相比 No local functionals 也有一些缺陷。在屏蔽交换泛函中不存在已知形式应力张量表达方式， 因此没有完全的非定域势可以用于单位晶胞结构优化或进行NPT/NPH 动力学。 这样利用这些泛函计算的光学性质很有可能是不准确的。在哈密顿量中引入一个完全非定域组元就可以解决这个问题，这个额外的矩阵元破坏了光学矩阵元素由位置算符转换为动量算符常用表达形式，使得哈密顿量对易很复杂。规范保守赝势和超软赝势赝势是利用平面波基组计算体系总能量中关键的一个概念，价电子与离子实之间强烈的库仑势用全势表示时由于力的长程作用很难准确的用少量的Fouri

14、er 变换组元表示。解决这个问题的另一种方法从体系电子的波函数入手，我们将固体看作价电子和离子实的集合体。离子实部分由原子核和紧密结合的芯电子组成。价电子波函数与离子实波函数满足正交化条件，全电子DFT 理论处理价电子和芯电子时采取等同对待，而在赝势中离子芯电子是被冻结的，因此采用赝势计算固体或分子性质时认为芯电子是不参与化学成键的，在体系结构进行调整时也不涉及到离子的芯电子。为了满足正交化条件全电子波函数中的价电子波函数在芯区剧烈的振荡，这样的波函数很难采用一个合适的波矢来表达。在赝势近似中芯电子和强烈库仑势替代为一个较弱的赝势作用于一系列赝波函数。赝势可以用少量的Fourier 变换系数来

15、表示。理想的赝势在芯电子区域是没有驻点的，因此需要平面波矢数量很少。众所周知的是现在将赝势与平面波矢相结合对描述化学键是很有用的。全离子势的散射性质可以通过构筑赝势得到重现，价电子波函数相位变化与芯电子角动量成分有关，因此赝势的散射性质就与轨道角动量是相关的。赝势最普遍表达方式是：.NL=|lml Vwhere | lm are the spherical harmonics and Vl is the Pseudopotential for angular momentuml .在不同角动量通道均采用同一个赝势值称为定域赝势（Local Pseudopotential），定域赝势计算效率更高

16、，一些元素采用定域赝势就可以达到准确描述。赝势的硬度（hardness ）在赝势的应用中是一个重要的概念，当一个赝势可以用很少的Fourier变换组元就可以准确描述时称为“软赝势”，硬赝势与此相反。早期发展的准确规范保守赝势很快就发现在过渡元素和第一周期元素（C、 N、O，等）中的描述十分“硬”，提高规范保守赝势收敛性质的各种方法都已经被提出，在CASTEP中采用了由 Lin 等提出的动能优化而来的规范保守赝势。Vanderbilt提出了另一种更基本的方法，放宽规范保守赝势的要求，从而生成更软的赝势。在超软赝势方法中，芯电子区的赝平面波函数可以尽可能的“软”，这样截止能量就可以大幅度的减少。超

17、软赝势与规范保守赝势相比除了“更软”以外还有其它的优点，在一系列预先设定的能量范围内遗传算法确保了良好的散射性质，从而使赝势获得更好变换性和准确性。超软赝势通常将外部芯区按照价层处理，每个角动量通道中的占据态都包含了复合矢。这样就增加了赝势的变换性和准确性，但同时是以消耗计算效率为代价的。可转移性是赝势的主要优点。赝势是通过孤立的原子或离子特定的电子排部状态下构建的，因此可以准确的描述原子在那些特定排部下芯区的散射性质。在相应条件下产生的赝势可以用于各种原子电子排部状态以及各种各样的固体中，同样也确保了在不同的能量范围内具有正确的散射状态。Milman 给出了不同化学环境和一系列结构中采用赝势

18、描述准确性事例。非定域赝势即使在最有效离散表示情况下，体系能量赝势计算依然占用了大量计算时间。此外，在倒易点阵空间采用非定域赝势会因原子数目增多而耗时以原子立方数增大，因此对于大体系是很适用的。赝势非定域性是指只有在超过原子芯区时它才会扩展，由于芯区是很小的，特别是当体系包含有许多的真空腔体时，在实空间采用赝势来计算就有很大的优势。这时计算量随体系中原子数目平方增长，因此是很适合大体系计算的。将电子划分为芯电子和价电子在处理交换- 相关相互作用时会产生新问题，在原子芯区两个亚体系叠加在赝势产生过程中很难完全去屏蔽。在赝势能量算符中与电子密度存在非线形关系的项就是交换 - 相关能。 Louie

19、等采用了一种简单的方法来处理芯电子和价电子密度之间非线性的交换-相关能。这种方法在很大程度上提高了赝势的可变换性，特别是自旋极化的计算更为准确。当准芯区电子不能简单处理为价电子时非线性核校正就很重要。另一方面将他们简单地包含在价层亚体系中从本质上可以避免 NLCC处理的必要性。规范保守赝势：采用赝势计算关键在于可以有效的对化学键的价电子进行可再现的近似，赝势与全势在超过离子实半径以后具有完全相同的函数形式。Figure 1. Schematic representation of the all-electron and pseudized wave functions and potenti

20、als两个函数平方幅度的积分数值应该是相同的，这等同于要求赝势波函数具有规范 - 保守性，比如每个赝波函数只能描述一个电子的行为。这样的条件就确保了赝势可以再现正确的散射(Scattering Properties)性质。生成赝势的典型方法如下所述：选择某个特定的电子排部状态（不一定就.是基态）全部电子计算在一个孤立的原子中进行。从而得到原子价电子能量本征值和价电子波函数。选择一个离子赝势或赝波函数参数形式，通过对参数的调节，使得赝原子计算和全电子原子赝势计算采用相同的交换 - 相关势，在超过截止半径后与价电子波函数形式相同，赝势的本征值等于价电子的本征值。如果电子波函数和赝势波函数满足正交归

21、一，两者在截止半径以外的匹配性决定了规范- 保守条件自动成立。离子赝势的截止半径是实际物理芯区的二到三倍。截止半径越小， 赝势越“硬”而适用性（ transferability ）好。计算精度和效率决定了实际中采用的截止半径的大小。规范 - 保守赝势优化在固体计算中依据能量的截止存在一系列优化赝势的方法，Lin 基于 Rappe 早期工作提出了下列赝势产生方法：在截止半径（cutoff radius ）内，赝势波函数可以表达为：PS4j( q R )li cl ( Rc )l( r )i 1 ai jl ( qi r ), jc ( qi Rc )l ( Rc )j l ( qi r ) 是球

22、形 Bessel 函数，在 r=0 和 r=Rc 之间有（ i-1）个零点。为保证赝势的实用性，截止半径越大越好。超过截止矢量qc 对动能最小化可得到系数i。Ek ( qi , qc )drPS,( r )2PS( r )qc2PS2lldqql( q )00在第一个方程中让 qc 等于 q4。其他的三个限制条件使得赝波函数在进行Lagrange连乘（ Lagrangemultipliers ）时保持正交化（ normalization ），并且使赝波函数在Rc 处的第一个二介偏微分是连续的。半径相关 Kohn-Sham 方程反转标准步骤产生的一个具有理想收敛性质的平滑赝势函数。Lee 提出了

23、进一步改进的方法，在 CASTEP数据库中固体规范保守赝势就是采用他的思想设计的。这种通用的方法消除了在特定的截止半径处赝波函数的二介偏微分必须是连续的条件，因为它是自动满足这个条件的。这样对于特定截止半径 Rc 允许我们通过调节qc 提高赝势的精度和计算效率。超软赝势（ ULTRASOFT PSEDUPOTENTIAL）为了能够使平面波基组计算中所采用的截止能量尽可能的小，Vanderbilt提出了超软赝势方法。众所周知规范 - 保守赝势在收敛优化中存在本身缺陷，所以就设计了另一种方法。超软赝势基础是在大多数情况下只有当紧密结合原子价轨道加权性分数大部分在芯区时，利用平面波基组计算才要求较高

24、的截止能量。在这种情况下，减少平面波基组的唯一方法就是解除（violate ）规范 - 保守赝势成立条件，将这些轨道中的电子从芯区移去。芯区的赝势就可以尽可能的“软”，从而使截止能量降低达到要求。从技术上讲，通过引入一个广义的正交归一化条件就可以完成。为了覆盖全部电子电荷，在芯区对由电子波函数模平方产生的电子密度进行适度放大（augmented）。电子密度划分成两部分：扩展在整个晶体中“软”部分和定域在芯区的“硬”部分。.固体中超软赝势公式超软赝势中总能量与采用其他赝势平面波方法时相同，非定域势VNL表达如下：VNl( 0)IInm , ID nmnm投影算符和系数D（ 0） 分别表征赝势和原

25、子种类的差别，指数 I对应于一个原子位置 。总能量用电子密度可以表示为：n ( r )i( r ) 2QI( r )IIinmi nm inm,I是波函数， Q(r ) 是严格位于芯区的附加函数 (Augment function) 。超软赝势完全由定域部分， Vloc ion ( r ) 和系数 D(0) , Q and 确定，这些变量计算方法在下文中将做介绍。引入一个广义正交归一条件来解除规范- 保守赝势的限制条件：i S jijS 是哈密顿重叠算符（Hermitian overlap operator ）.S1IIqnm nmnm,I系数 q 是通过对Q（ r ）积分得到，超软赝势的Ko

26、hn-Sham方程可以写为：HiiSiH 代表了动能和定域势能之和，如下所示：.IIIHTVeffnm,I Dnm nmion在 Veff 中包含离子定域势 Vloc ( r ) ，Hartree 势和交换 - 相关势等项。 通过定义一些新参数就可以将因附加（ augmented）电子密度而产生所有项全部包含在赝势的非定域部分。I(0)ID nmD nmdrVeff ( r ) Qnm ( r )与规范 - 保守赝势对比，不同之处在于在超软赝势中存在重叠算符S，波函数与D 有关而且事实上投影算符函数（ projector function ）数量要比规范 - 保守赝势中大两倍多。与附加（ au

27、gmented）电荷相关的一系列计算可以在实空间 (real space) 中进行，这与函数中定域势的性质有关。多余的步骤不会对计算效率产生较大的影响。在 Laasonen 文献中提供了超软赝势计算的详细方法以及总能量微分表达式。赝势生成：与规范- 保守赝势情况一样，在自由原子上对所有的电子进行计算，得到屏蔽原子势VAE( r )（ screened atomic potential）。每个角动量选择一系列的参考能量l ，一般两个能量参考点就足够了。这些能量参考范围必须包含良好散射性质，在每个参考能量处求解与半径相关的Kohn-Sham 方程，得到规则初始点。 选择截止半径， 对上面产生的每个

28、全电子波函数构筑一个赝势，唯一的限制条件是它必须在Rcl处与平滑相交。定义一个比所有芯区半径稍大的辅助半径R。最后就形成了定域轨道（超过 R 时就消失）：n(T V )nnloc.以及它们矩阵内积（inner products）：B nmnm这样就可以定义用于固体计算的变量(V loc ion ( r ), D (0) , Q and):Q (r )( r )( r )(r )(r ); qdrQ (r )nmnmnmnmnmn( B 1)nm mm采用去屏蔽（ descreening procedure）方法计算locion( r ), D(0)系数 :VionV loc( r )V loc

29、( r )V H( r )V xc( r )( 0)B nmm q nmdrV loc ( r ) n ( r )D nm在去屏蔽方法中可以引入非线性核校正方法（ The nonlinear core correction (NLCC) ），这与规范 - 保守赝势中所采用的方法完全一致。在以下情况下超软赝势是很适用的：赝本征值与所有电子本征值相同，在芯半径截止区以外赝轨道波函数与价电子波函数匹配一致；对于每个参考能量散射性质都是正确的，这样通过增加参考能量点数目就可以系统的提高赝势的适用性；在参考电子排部情况下，赝势价电子密度与全价电子密度相同。关于非线性核校正Louie 等人第一次提出了非线

30、性芯校正，使得赝势对磁系统的描述更准确。然而，对于非自旋极化体系中准芯区电子， NLCC也具有同样的作用。 DFT总能量准确表达需要 NLCC，如下：EtotTEionEeeExc在赝势的计算式中，电子密度分别来自于芯区电子和价电子。将芯区能量假设为一常数并切不计入计算。用一个价电子密度和由赝势计算得到的离子定域势Eion 来代替总电子密度，这样芯区电子与价电子之间所有的相互作用全部转移到赝势上。由此可以推断电子密度线性化只是对动能和简单非线性交换- 相关能的.一个近似，很明显当芯区电子和价电子在空间很好分离时是一个良好的近似。但如果两个区域电子密度的叠加密切时，计算体系本身就会产生错误，进一

31、步减弱赝势实用性。解决NLCC问题的方法就是调节赝势生成方法以及在固体中计算方法。在产生赝势时每个角动量通道对应一个屏蔽势，并且满足一定的条件，比如规范 - 保守，赝波函数本征值与全电子波函数本征值相同等。这些屏蔽势 （screened potentials）对应的原子赝波函数（atomic pseudowavefunctions）仅表示价电子。从这些波函数可以得到价层赝电子密度（ pseudo charge density），通过对势的屏蔽得到“光秃”离子势（bare ）：Vionl (r )V l (r )Vee(r )Vxc(r )由于交换 - 相关势泛函是电子密度的非线性函数，对自旋极

32、化体系采用这种方法产生的离子势与价电子排列有关。 Louie 等提出了将上面方程替换为如下表达：V l (r ) V l ( r ) V ( r ) Vxc(r )c( r )ionee在屏蔽原子势中减去总交换- 相关势。此外，在计算交换- 相关势时芯区电荷必须加到价电子中去，这个额外原子状态信息传递给CASTEP，在所有计算中芯区电荷认为（deemed）是相同的，这种做法的一个缺点是在利用赝势计算时芯区电荷很难准确的用Fourier网格表示。而且通常芯区电子密度比价电子密度大，这很容易将与价电子密度有关的影响掩盖掉。以下部分将对部分芯区校正方程建立做介绍，该方法充分的认识到价电子与芯电子密度

33、重叠的区域才是我们感兴趣的。靠近原子核的芯电子密度不会产生物理结果，虽然有如上所述的一些问题。部分NLCC采用一个在特定半径以外与c 一致的函数替代全芯电子密度，在原子核周围这个函数起伏是平滑的。在CASTEP中对一些特定元素在赝势中采用的部分芯区校正使用了数值化的芯区电子密度。在规范保守赝势中虽然有相关的内容，但在计算中并没有采用这个方法。A Introduction to DFT第一性原理（ The first principle ）计算也称为从头算起（ab-initial calculation ），由于固体的许多基本的物理性质是由其微观的电子结构决定的，因此通过求解多粒子系统的Scho

34、dinger 方程，来获取固体全部的微观信息从而预测宏观的性质。利用这个思想建立的能量的哈密顿量非相对论形式可以表示如下：h22me i2i(r ,r , r ,.r, r r.r,)k e1 e2 e3ei,n1,n2,njh221(r ,r ,r ,.r, r r .r,)2 j mj k e1 e2 e3ei, n1,n2,njnj(2ze2z z e2ejj1j2 ) (r ,r ,r ,.r,r r .r),i ireirej i,j rei rnjj ,j rk e1 e2 e3 ein,1,n2, njr1,21 2 nj1nj2i1 i2j1 j2E(r ,r ,r ,.r,r

35、 r.r),k k e1 e2 e3ein,1,n2,nj考虑到原子核与核外电子质量差别以及电子驰豫时间比原子核驰豫时间要小三个数量级，因此利用.Born-Oppenheimer近似将原子核运动和电子的运动分离，从而将体系波函数划分为电子波函数和原子核波函数两个部分，分别用和表示：(r, r , r,.r, rr.r, )(r, r,r ,.r,) (rr.r, )k e1 e2e3ei, n1,n2,nje1 e2 e3ein1,n2,nj能量的哈密顿量可以分解为如下的两个方程：2222ehezje( 2miiiirri, jr)(r, r, r,.,rE(r, r, r,.,rere1e2

36、e3ei)kee1e2e3ei)1,2eiejeinji1 i2(h212zje2zj1zj 2e2)(rr.,r) En(rr.,r)2jmnjji , jreirnjj1, j2rrn1, n2,njknn1, n2,njj1j2nj1nj 2第一性原理严格求解仅在氢分子中实现了，对于多粒子体系的计算几乎是不可能的。目前均采用不同的近似方法来实现计算，主要方法有 Hartree-Fock 近似和 DFT 近似。在 Hartree-Fock 近似中体系的哈密顿量表示如下：ETotaliHF( J ijK ij)ij iHF 为第 i 个电子的 Hartree-Fock 的轨道能，Jij 是库

37、仑积分，表示电子静电互斥能，K为交换积分。iji交换积分所代表的交换能指电子由于自旋平行而引起的电子轨道库仑能量减少的部分。密度泛函理论 (Density Functional Theory) 建立了将多电子体系化为单电子方程的理论基础，并且给出了有效势计算方法，是目前研究多粒子体系性质的一种普遍使用的重要方法。该理论认为对于处于外势场 V （ r）中相互作用的多电子系统，电子密度分布函数 (r) 是决定该系统基态物理性质的基本变量。密度泛函理论中体系的能量泛函表示如下：Et ( )T( )U ( )Exc()T () ： Kinetic energy; U ( ) :classical el

38、ectrostatic energy; E xc ( ) :exchange and correlation energy由上表达式可见体系能量是电子密度的泛函，因此可以进一步将上式表达为：.Et ( r )e2( r )( r )E xc ( r )V ( r ) ( r ) dr T ( r )2rdrdrr在上式中第一项为电子在外场中的势能，第二项为电子的动能，第三项为电子相互之间的库仑能，第四项为交换相关能，最后一项形式是未知的。系统的电子密度分布可以表达如下：( r )ii ( r ) 2利用上式可以将动能项表示为：T ( r )ni2ii2U( )表达为：U (r)n Nr(r )

39、z(r )1(r ) (r ) 1i(r(r )NzaziaR r i2 iji 1 j 2 r r1) j2aR Ra1 2a aNza11Nzaz(r1)R r2(r1) (r2)r raR Raa 11 2a aVe r(r )VN(r )VNN2E xc ( ) 形式确定有两种方法： 局域密度近似 （ LDA,Local Density Approximation ）和广义梯度近似 （ GGA, General Gradient Approximation ）。在局域密度近似（ LDA ）中采用了均匀电子气的分布函数推倒出了非均匀电子气的交换-相关能泛函，从而得到E xc ( ) 的具

40、体形式。从近期计算结果相关报道来看采用局域密度近似（ LDA ）计算在绝缘体中会产生较大的误差，而且对带隙宽的半导体等得到不正确的结果。采用局域密度近似（ LDA ）主要的缺陷现归纳如下：1 对光学跃迁带隙预测很差（一般是过低估计带隙宽度） 。这虽然对基态性质如电荷密度，总能量以及力影响不大，但在导带状态计算中却是个大问题，如关于光学性质，运输性质等的计算。在诸如光伏.装置等领域的研究中，带隙就是个很重要的问题。采用“剪刀”(Scissors)工具在固体带隙计算中很有用，但对我们未获得实验结果的物质，是不能采用这个方法的。2 对类似于二氧化硅这样的电子气分布极不均匀体系，基本假设中关于电子密度

41、分布在空间是缓慢变化的条件是不满足的，这样的体系采用LDA 处理就存在难题。3 LDA 简单的认为计算体系是顺磁性(Paramagnetic)的，对于包含未配对(Unpaired) 自旋体系采用局域自旋密度近似（LSDA ）（对自旋向上 (spin up) 和向下 (spin down) 的电子分别采用密度泛函计算）是很有用的，比如费米能级(Fermi level) 处半填充的系统。4 最后一个很少关注的领域就是玻璃陶瓷工业，LDA对弱的结合键（如偶极涨落）很难描述，氢键(Hydrogen bond) 在 LDA 中也无法获得准确的计算结果。GGA 近似则改进了L(S)DA, 将相关交换能确定

42、为电子密度极其梯度的函数，在GGA 学派中以Perdew 等人认为交换相关能的泛函形式应该以一定的物理规律为基础，构造了著名的PBE 泛函。将电子密度分布函数带入体系能量电子密度泛函中，对泛函变分求极小值，可以得到Kohn-Sham 方程：Fi ( r )Eii ( r )2FVKS ( r )V ( r )V ( r )( r )d r iE xcKSr i( r )r交换 -相关能可以按照下式计算：E xc ( r ) E xc ( r ) dr:number of particles;E xc ( r ) :exchange-correlation energy per particle

43、s in an uniform electron gas ;( r ):distribution function of electron density.Exc( r )Excxc( r )( r ) 称为交换 -相关势和，表示为：在 Castep 计算中采用了周期性边界条件，单电子的轨道波函数满足Bloch 定理，采用平面波展开式有：i (r )eiK ?Ri (r )周期性边界条件下的波函数扩展为一系列分离的平面波波矢，这些波矢与晶体的倒易点阵矢量相联系。i (r )Ci ,G eiG ?R2.2 晶体光学性质的计算基于以下原理：电磁波在真空以及某种材料介质中传播时差别可以用一个复数式的折射指数来表示：Nnik在真空中N 为实数，而且其大小为1；在其他介质中时若材料对于光是透明的则是一个纯实数，虚部对应材料的吸收系数（Adsorption Coefficient ）。它们之间的关系方程2 所示：2 kc2吸收系数表示的是电磁波通过单位厚度的材料时能量的衰减分数，通常可以用材料焦耳热的产生来衡量。反射系数（ Reflection Coefficient ）可以简单通过将垂直光束照射材料的表面引起2k21 N(n 1)2.3R(n 1)2k21 N在计算光学性质时一般先计算虚部的介电常数，其他的性质与介电常数之间建立关系。虚部