DFT电子结构计算李震宇

  电子结构理论与计算李震宇中国科学技术大学合肥微尺度物质科学国家实验室Version:20110812  1主要内容第一章电子结构理论基础BO近似，HF理论，电子关联与后HF理论，密度泛函理论（DFT），交换关联泛函，赝势，半经验电子结构理论第二章电子结构计算方法原子基组，平面波基组与周期性边界条件，自洽场（SCF）第三章电子结构分析电荷密度，分子轨道，能量第四章电子结构理论与计算的应用动力学与统计力学性质附录量子力学基本原理，常见数值算法  2前言理念：1.实用：介绍计算参数的意义和选择策略；2.简洁：注重物理概念，不过分纠缠于细节。参考资料：[1]RichardMMartin,ElectronicStructure:BasicTheoryandPracticalMethods.CambridgeUniversityPress,Cambridge,2004.[2]ChristopherJCramer,EssentialsofComputationalChemistry:TheoriesandModels.2ndEd.Wiley,England,2004.感谢所有帮助我将本讲义电子化的同学，以及所有其他选修本课程的同学。讲义中错误和疏漏之处在所难免，如有任何意见与建议，请发送到[e-mail]，十分感谢！  3第一章电子结构理论基础由于计算能力较弱，早期的科学研究以理论和实验的互动为主，其中理论研究通常使用唯象理论。随着科技的发展，现在已经有可能完全从量子力学基本原理出发来计算体系的各种性质。这就使得计算成为一种重要的独立的科学研究手段，与理论和实验形成三足鼎立之势。  图示：新的科学研究范式 三者之间有紧密的联系：计算可以帮助理解实验结果，例如扫描隧道显微镜(STM)观测到的是局域态密度，要从中获得几何结构的信息，需要计算模拟来配合；计算也可以验证理论猜想，例如通过计算机实验来检验材料设计原理。我们希望通过科学计算获得理论设计与计算表征的能力。 体系的大多数性质都决定于它的电子结构，即电子的空间与能量分布。  1.1基于波函数的量子化学方法1.1.1波恩-奥本海默(BO)近似在原子中，核的质量比电子大很多，运动也比电子慢得多。1amu(12C/12)是电子静止质量的1822.83倍。因此，通常在研究某一瞬间电子结构时可以忽略原子核的速度，假定它和原子核长期固定在某位置时的电子结构一样。这就是BO近似。其核心思想是把核的运动和电子的运动分开处理：处理电子运动时，认为核是固定不动的；处理核运动时，认为快速运动的电子建立一个平均化了的负电荷分布，核在电子的负电荷平均场中运动。 对于包含许多电子和原子核的体系，其定态薛定谔方程可写为： ,,HErRrR 其中  222222222IJIIiIiiIijIJIiiIijIJeZZeZeHmmrrr＝- 上述哈密顿量中第一项核动能项为小项。引入BO近似，对原子核和电子进行分离变量 ,;NelrRRrR elelelHERrRr  其中NR是描写原子核状态的波函数，它只和所有核的位置有关；而elHR和  4;elrR是描写电子状态的哈密顿量和波函数，以原子核位置为参数【输入参数：POSCAR】。ER和原子核的位置有关，通常被称为势能面(PES)。知道势能面以后可以从其极小点得到几何构型，从其势阱深度得到反应能，从其鞍点得到过渡态，等等。BO近似又称为“绝热近似”，它认为电子总是不断地适应核位置的变化，或者说跟随着核且始终围绕着核来运动，运动过程中不会发生从一个势能面E1至另一个势能面E2的跃迁，即吸收或放出能量发生电子能级跃迁。这是因为电子态涉及到的能量量级为电子伏特，远大于室温对应的能量，所以电子通常处于基态。  图示：石墨烯能量随C-C键长变化的势能面   1.1.2单电子近似与平均场近似由于电子之间存在相互作用，所以当电子数较多时，即使在BO近似下薛定谔方程仍然很难求解，需要作进一步的简化。昀粗糙的近似是假设电子之间没有任何相互作用。这样，体系的哈密顿算符简化为： iiHh 其中      222IiiIiiIeZhmr 于是，我们可以先只考虑体系中的某一个电子，其定态薛定谔方程为iiiih，本征函数i为单电子波函数，通常也称为分子轨道(MO)，本征值i为对应的分子能级。下面的问题是怎么从单电子波函数写出总的多体波函数呢？    5对无相互作用电子，考虑泡利不相容原理，体系的波函数及其对应的能量本征值为 12...HPn iiE 这种由分子轨道直接相乘得到的多体波函数称为Hartree积(HP)。 单电子近似并不意味着需要完全忽略电子间的相互作用。对严格的Hamiltonian采用HP多体波函数作为试探波函数能量(ˆelH在多体波函数下的期望值)变分可以得到新的单电子哈密顿算符 2222IiiiIiiIeZhvmr 其中 jijjiijvdrr 2jje 显然，该单电子哈密顿算符中的第三项可以看成体系中其他电子的库仑相互作用的一个平均场近似，该平均场不随时间而改变。因为在求平均场时其他电子的坐标都是积分变量，所以ˆih仅与第i个电子的坐标有关。对这种采用HP多体波函数的平均场近似，体系的总能量是 222||||12ijHPiijiijijeEddrrr 其中第二项为重复计数修正。 我们看到，在平均场近似中，需要求解单电子本征方程才能得到体系的解，但要求解单电子方程又必须知道各个轨道的电荷分布。这看起来自相矛盾。在实际操作中我们可以从N个初始单粒子波函数0i出发，建立相应的(0)ˆih，然后通过解单电子本征方程，得到一组新的分子轨道1i，再用这组1i来建立新的(1)ˆih，重复这个过程，直到昀后得到的ni与上一步的1ni相等或近似相等。这样一个过程称为自洽场(SCF)方法【输入参数：EDIFF】。  1.1.3Slater行列式与Hatree-Fock自洽场显然，HP多体波函数不满足电子（费米子）全同性所要求的波函数的反对称性。为简单起见，先讨论无相互作用两电子体系 12Hhh 其中 221112hvmr   6222222hvmr 假设,12ababE是两电子体系本征方程的解，则 12ˆ122112abababP 也是体系本征方程的解。这样，对应于同一个本征值，有两个线性独立的本征函数，但这两个本征函数都不满足反对称条件，为此作一线性组合，得： 112212Aabab 为体系满足波函数反对称条件的解，用行列式的形式表示为 111222abSDab 该行列式被称为Slater行列式(SD)。对于由N个电子组成的全同粒子体系，可以将上述结果进行直接推广，相应的定态波函数为 12121212121112221,,,,!NNSDNNNNNNNrrr 上式的反对称性是很明显的，因为任意两个粒子坐标的互换，相当于行列式中相应的两列元素的互换，将改变行列式的符号。同时，若有两个或两个以上单粒子波函数相等，则0SD，意味着这样的状态是不存在的，这就是Pauli不相容原理。 Fock建议将Hartree-SCF过程推广到采用Slater行列式作为多体波函数的情形。对两电子体系的电子库仑相互作用算符，我们有【作业】 12121SDSDababddJKrrr 其中 22121212121112abHPHPabJddddrrrrrr＝ 121211122abababKddrrr＝ 分别代表电子间的库仑排斥能和交换能。交换能Kab是采用Slater行列式以后比HP多出来的同自旋电子间的一种相互作用。Roothann提出可以基组表象来建立了HF‐SCF的矩阵代数。如果将分子轨道写成一些基函数的线性组合，则单电子哈密顿量在这组基组下可以表示为一个矩阵。常用的基组有原子轨道和平面波。对闭壳层体系，昀后可以得到久期方程det(F-ES)=0，其中单电子哈密顿量f通常被称为Fock算符，其矩阵元为   7221122IIeIFZPmr 其中 212121122eddrrr 2occiiiPaa 由于Fock矩阵的表达式中包含密度矩阵，我们需要自洽求解。 至此，我们通过如下步骤得到了作为量子化学基础的Hartree‐Fork理论。   1.1.4电子关联与post-HF方法HF理论采用单行列式多体波函数，它假定一个电子在由原子核和其它电子形成的平均势场中独立运动，没有考虑电子间的瞬时相关作用，也就是电子关联。但是，HF理论中采用的Slater行列式波函数满足Pauli不相容原理，所以它在一定程度上包含了同自旋电子的关联作用，这部分电子关联被单独命名为交换相互作用。所以，通常所说电子关联能是指精确的基态能量与HF能量之差：Ecorr=E-EHF。 关联能一般只占体系电子总能量的0.3~2%，但是它可以与化学反应的反应热或者活化能相当甚至更大。一般说来，化学反应途径上电子关联程度不同。所以，在考虑化学反应和电子激发等问题时电子关联十分重要。由电子瞬时作用产生的电子关联称为动态关联。此外，还有非动态关联，它指体系几个对称性相同的电子组态函数（单行列式）接近简并时导致HF方法给出很大的误差。常用的电子关联计算方法【输入参数：计算方法】有 组态相互作用（CI）：Configuration Interaction 耦合簇方法（CC）：Coupled‐cluster method 多体微扰理论（MPx）：Moller‐Plesset perturbation theory 密度泛函理论（DFT）：Density function theory 前面三种都是以HF的为基础，经常被称为post‐HF方法。而DFT基于完全不同的思路，将在下一节中专门介绍。 1.CI理论可以通过引入激发组态来考虑电子关联   800||||...viroccviroccaaababiiijijaiabijCCC 一般用S, D, T, Q等字母来表示从低到高的激发组态，其中昀简单的情况是只引入单激发组态(CIS)。CIS可用来计算激发能，但对基态能量没有修正【作业】。昀复杂的情况是考虑所有激发组态，即Full CI，它是特定基组下的精度极限，一般只能处理非常小的体系。 考虑到精度与计算量的平衡，一种有效的方法是将分子轨道分为不同的类型。对昀重要的分子轨道，如Fermi能级附近的轨道，考虑所有可能的电子组态。这些轨道组成空间称为全活性空间(CAS)。对次重要的轨道，考虑一定限制条件下的可能的电子组态(RAS)。在确定需要考虑的电子组态以后，通常在优化CI组态系数的同时优化单粒子轨道，即进行有多个组态的自洽场(MCSCF)。    在MCSCF的基础上还可以再做CI。这时，CI不以HF单行列式为基础。而是由MCSCF计算的波函数中选择组合系数较大的组态作为参考组态，从这些参考组态的占据轨道上移动电子到虚轨道上去，产生新的激发组态函数。由于有多个参考组态，这种方法称为MRCI。 2.Coupled-Cluster方法耦合簇理论是目前昀流行的包含电子关联效应的量子化学方法之一。耦合簇理论通过对波函数施加一个如下的ansatz来改进对体系的解： ˆTccHFe 其中ˆT为簇算符。如