第一性原理计算原理和方法

PAGE45-第二章计算方法及其基本原理介绍化学反应的本质是旧键的断裂和新建的形成，参与成键原子的电子壳层重新组合是导致生成稳定多原子化学键的明显特征。因此阐述化学键的理论应当描写电子壳层的相互作用与重排，借助求解满足适当的Schrodinger方程的波函数描写分子中电子分布的量子力学，为解决这一问题提供了一般的方法，然而，对于一些实际的体系，不引入一些近似，就不可能求解其Schrodinger方程。这些近似使一般量子力学方程简化为现代电子计算机可以求解的方程。这些近似和关于分子波函数的方程形成计算量子化学的数学基础。图2-1分子体系的坐标2.1SCF-MO方法的基本原理分子轨道的自洽场计算方法(SCF-MO)是各种计算方法的理论基础和核心部分，因此在介绍本文计算工作所用方法之前，有必要对其关键的部分作一简要阐述。2.1.1Schrodinger方程及一些基本近似为了后面介绍各种具体在自洽场分子轨道(SCFMO)方法方便，这里将主要阐明用于本文量子化学计算的一些重要的基本近似，给出SCFMO方法的一些基本方程，并对这些方程作简略说明，因为在大量的文献和教材中对这些方程已有系统的推导和阐述[1-5]。确定任何一个分子的可能稳定状态的电子结构和性质，在非相对论近似下，须求解定态Schrodinger方程（2.1）其中分子波函数依赖于电子和原子核的坐标，Hamilton算符包含了电子p的动能和电子p与q的静电排斥算符，(2.2)以及原子核的动能(2.3)和电子与核的相互作用及核排斥能(2.4)式中ZA和MA是原子核A的电荷和质量，rpq=|rp-rq|，rpA=|rp-RA|和RAB=|RA-RB|分别是电子p和q、核A和电子p及核A和B间的距离（均以原子单位表示之）。上述分子坐标系如图2.1所示。可以用V(R,r)代表(2.2)－(2.4)式中所有位能项之和(2.5)原子单位上述的Schrodinger方程和Hamilton算符是以原子单位表示的，这样表示的优点在于简化书写型式和避免不必要的常数重复计算。在原子单位的表示中，长度的原子单位是Bohr半径能量是以Hartree为单位，它定义为相距1Bohr的两个电子间的库仑排斥作用能质量则以电子制单位表示之，即定义me=1。Born-Oppenheimer近似可以把分子的Schrodinger方程(2.1)改写为如下形式(2.6)由于组成分子的原子核质量比电子质量大103倍～105倍，因而分子中电子运动的速度比原子核快得多，核运动平均速度比电子小千倍，从而在求解电子运动问题时允许把电子运动独立于核运动，即认为原子核的运动不影响电子状态。这就是求解(2.1)式的第一个近似，被称作Born-Oppenheimer近似或绝热近似。假定分子的波函数Ψ′可以确定为电子运动和核运动波函数的乘积(2.7)其中Ф(R)只与核坐标有关，代入方程(2.2)有对于通常的分子，依据Born-Oppenheimer原理有：AΨ和A2Ψ都很小，同时MA≈103～105，从而上述方程中的第二项和第三项可以略去，于是易知也即该方程可以分离变量而成为两个方程(2.8)(2.9)方程(2.8)为在某种固定核位置时电子体系运动方程，而方程(2.8)时核的运动方程。E(R)固定核时体系的电子能量，但在核运动方程中它又是核运动的位能。此时分子总能量用ET代表。因此，在Born-Oppenheimer近似下，分子体系波函数为两个波函数的乘积(2.7)式。分子中电子运动波函数Ф(R)分别由(2.8)和(2.9)式确定。电子能量E(R)为分子的核坐标的函数，从(2.9)式看出它又是核运动的位能。在空间画出E(R)随R的变化关系称为位能图。单电子近似体系的电子与核运动分离后，计算分子的电子波函数Ψ归结为求解下面的方程(2.10)(2.10)式是量子化学的基本方程，目前已有多种求解这个方程的方法。这些方法的区别首先是构成Ψ的方式及其相应的近似。最常用的是Hartree建议的单电子近似［6］。在多电子体系中，所有电子势相互作用的，其中任意电子运动依赖于其它电子的运动。Hartree建议把所有电子对于每个个别电子运动的影响代换成某中有效场的作用。于是每个电子在核电荷及其余电子有效场产生的势场中运动仅依赖于电子坐标。从而，电子运动分开了，对于多电子体系中每个电子可以引入单电子波函数，这种单电子波函数是(2.10)式单电子Schrodinger方程的解，其中含有算符1/rpq项，用只依赖于所研究电子坐标的有效场代替。整个多电子体系波函数等于所有电子的单电子波函数（轨道）乘积。电子还具有自旋角动量s，其分量sx,sy和sz满足普通角动量算符的对易关系。算符s2和sz完全给定了电子的自旋，电子自旋波函数()满足方程(2.11)其中是自旋坐标，通常把对应于自旋1/2的波函数记为()，而把自旋ms=-1/2波函数记作()。在非相对论近似下和不存在外磁场时，电子的自旋和空间坐标无关，因此，因此电子的自旋轨道可取成(2.12)考虑到自旋变量的多电子波函数由自旋轨道组成，他应当是体系总自旋S2及其Sz的本征函数(2.13a)(2.13b)构成体系多电子波函数Ψ时，必须考虑Ψ相对于任一对电子交换的反对称性要求，此所谓Pauli原理［7］。因此，一般不求出Hartree方法的简单乘积型波函数Ψ，而是求出对应于按自旋轨道电子的所有可能置换方式的Slater行列式波函数，此为Hartree-Fock方法。对于置于n=N/2轨道的Ψ上的N电子体系，单电子近似下波函数Ψ写为(2.14)该式的Slater行列式是保证反对称性要求的唯一这类函数。引入单电子近似便确定了波函数Ψ的形式，用它可以求解方程(2.10)。显然在一般的情况下，Ψ应当包含(2.14)型行列式的线性组合，同时满足(2.13)式的限制。若(2.12)式中自旋部分是单电子自旋投影算符Sz的本征值，则(2.13b)式就满足。当分子的n个轨道每个均为自旋反平行电子对占据时（闭电子壳层），一个行列式波函数(2.14)就已满足(2.13a)和(2.13b)。对于含有未配对的电子体系，这是做不到的，此时体系波函数是对应于各种轨道填充方式（不同组态）的Slater行列式Ψl的的线性组合(2.15)当适当选择行列式前系数al时，条件(2.13a)和波函数的反对称性要求均可以满足。由于存在着电子运动的相关，不明显处理(2.10)式中1/rpq项的单电子近似，完全忽略了这种相关效应，所以，Hartree-Fock单电子近似使波函数的计算产生了误差。变分原理上述单电子近似只是给出了所求解体系多电子波函数的一种形式，变分法提供了求解方程(2.10)的一种方法。Schrodinger方程(2.10)的解对应于稳定态能量。因此若波函数Ψ是(2.10).的解，那么对于任意微小变化Ψ，取能量平均值(2.16)的变分应等于零，即(2.17)(2.16)式中积分是对Ψ的所有变量进行的，并且已假定Ψ是归一化的，即(2.18)由于我们寻找对应于体系基态的波函数，总能量应当是极小值。因此，对单电子轨道施行变分就给出这种型式波函数，能量是极小值并满足(2.17)式。从而求得的波函数Ψ就是多电子体系基态Schrodinger方程所欲求的解。显然，为了施行变分，波函数Ψ的型式应当充分好。两种途径可以保证这一点：①取展开式(2.15)是从充分多项，且固定轨道Ψ只对系数al变分；②局限于尽可能少的行列式Ψl，若有可能做到就取一个，但此时把每个ψ表成可能的简单形式。鉴于这种选择，区分出两类广泛应用的量子化学方法，价键(VB)法和分子轨道法(MO).在价键法中，用孤立原子的原子轨道(AO)作为单电子波函数ψ去构成Slater行列式Ψl。原子轨道的不同选择对应于不同的行列式Ψl。对于(2.15)施行的变分，可得到确定系数al的方程。为了充分靠近体系的能量，必须在(2.15)式中选用足够的多项，即用多行列式波函数进行运算。用原子轨道线性组合分子轨道(LCAOMO)法提供了另外一种选择相应于体系能量极小的多电子波函数方法。此时，对应于分子中单电子态的分子轨道ψi写成原子轨道φμ（基函数AO）的线性组合(2.19)实际上，这种展开有完全合理的基础。因为靠近某个原子的电子所受的作用基本上是由该原子产生的场引起的，所以该区域中电子波函数应当近于原子轨道。展开该式对求解变分问题的优点是明显的。如果(2.15)式中选用极大数目的项，那么VB法和MO法就都给出同样的能量E和波函数ψ，当然表达式不完全相同。这种唯一性的原因很简单，因为使用LCAOMO的的每个行列式均可以展开为AO组成的一些行列式。在一般情况下，每个MO组成的行列式应展开成AO组成的所有行列式。因为波函数ψ应通过AO组成的行列式完全集合表达，从而，当使用完全集合时，MO法与VB法所描述的ψ就等价。当然，不用完全及表达时，两种方法的等价性就破坏了。在极端性况下，某种方法中可以取一个行列式，此时可以直接看到MO法的优越性。对于MO法，允许采用单行列式表达ψ（至少对于闭壳层体系），进而，通常由一些正交分子轨道组成行列式(2.20)其中δij是Kronecker符号。从而是计算大为简化，并能比VB法更简单地确定(2.19)式的方程系数。同时，MO法的基本方程能很好的适应现代电子计算机的能力。由于这个原因，现代的MO方法已经成为最常用的计算多电子分子的电子结构的基本方法。闭壳层体系的Hartree-Fock-Roothaan方程在分子轨道范围内，对闭壳层体系，在单电子近似下，用两个自旋反平行电子填充每个分子轨道ψ，可以构成一个Slater行列式(2.14)型波函数，选择轨道(2.12)的自旋部分满足(2.11)式，则保证了(1.13b)条件。根据变分原理，若轨道ψ使得分子能量(2.16)取极小值，就求出了所研究多电子体系方程(2.10)的解。将波函数(2.14)代入(2.16)式，并进行一些推导［见引文1－4］，可得闭壳层分子的电子能量表达式(2.21)此处Hii是对应于分子轨道ψi的核实Hamilton量Hcore(1)的单电子矩阵元(2.22)而Hcore包含电子动能算符和分子中原子核对电子的吸引能算符(2.23)下面两式分别表示库仑积分Jij和交换积分Kij(2.24)(2.25)积分取遍电子1和2的全部空间坐标。从(2.22)至(2.25)可以看出(2.21)式中各项的物理意义。显然，单电子积分Hii表示在核势场中分子轨道上电子能量，由于每个ψi轨道上占据两个电子，所以乘以2。双电子库仑积分Jij表示ψi和ψj轨道上两个电子间平均排斥作用。由于波函数的反对称性要求出现了交换积分Kij（在Hartree方法中不考虑它），减小了不同轨道ψi和ψj上平行自旋电子间相互作用，正是它们描写了相同自旋电子运动的交换相关。然而，在Hartree-Fock方法中，还是没有考虑反平行自旋电子间库仑排斥引起的电子相关效应。为了求解Ψ的最优近似，必须选择一定形式的分子轨道ψi是总能量最小。这些分子轨道相互正交，在LCAO近似下表成原子轨道ψμ的展开(2.19)式。Roothaan最先解决了这一问题［8］。关于ψi对AO基展开系数的方程成为Hartree-Fock-Roothaan方程（间记为HFR方程）。下面简要推导这个方程。首先从写ψi的LCAO展开式分子轨道正交归一化条件给出对MO系数的附加限制(2.26)其中Sμν是原子轨道ψμ和ψν间的重叠积分(2.27)我们不难用AO基写出(2.23)至(2.25)式(2.28)(2.29)(2.30)其中Hμν是核实Hamilton量(2.23)相对原子轨道ψμ和ψν的矩阵元。双电子相互作用积分<μν|λσ>表达式为(2.31)它表示电子云分布ψμψν和ψλψσ间相互作用。从而能量E的表达式(2.21)变为(2.32)引入原子轨道基的电子密度矩阵元(2.33)上式化为(2.34)或者令(2.35)则可以用用矩阵形式写为(2.36)其中电子相互作用矩阵G定义为(2.37)矩阵J(R)描写库仑相互作用，而K(R)描写电子交换相互作用，其矩阵元分别为(2.38)(2.39)作用于两个矩阵的运算Sp是指把此两个矩阵所响应矩阵元的乘积求和。在LCAO分子轨道法中，通常假定基原子轨道是固定的，而型式不变，因此对分子轨道变分归结为对展开系数cμν的变分(2.40)应当指出，(2.40)式中的ψμ不变性在原则上是不必要的，已有把ψμ看成可变的一类计算方法，然而是极其困难和复杂的，常用的MO法一般不做这类计算。当对分子轨道变分时，借助Lagrange乘法可以是能量极小化；这是极小化泛函即变化系数cμi，使稳定点δG＝0。其中E由方程(2.34)定义。取G的变分，则有(2.41)由于δcμi*的任意性，所以必须有(2.42)因为分子轨道在酉变换下的确定性，可以自由选择非对角Lagrange乘数为零。定义Fock矩阵元(2.43a)或记为矩阵形式(2.43b)其中H为Hartree-Fock矩阵的单电子部分，而G为双电子部分。于是方程(2.42)可以写为(2.44a)或记为矩阵形式(2.44b)此为HFR方程，其中E是Hartree-Fock算符本征值组成的对角矩阵，S是基原子轨道的重叠积分矩阵。解(2.44a)或(2.44b)，就给出按原子轨道展开的分子轨道系数和单电子分子轨道能量εi。然而，上述HFR方程是数学上广义特征值问题。为了求解方便，先借助于下述变换将它化成标准特征值问题。由于S是Hermite矩阵，可以用酉变换θ使其对角化，即(2.45)从而可由其对角元平方根倒数构成矩阵S-1/2，同时S1/2S-1/2和S-1/2SS-1/2=1成立。此时可做变换(2.46)(2.47)则方程(2.44)化成亦即(2.48)这是标准特征值问题，可以采用数学上标准的对角化方法处理。显然这些方程都是系数cμi三次联立方程组，必须用迭代法求解：当F(0)=H时，解上述方程(2.44a)或(2.44b)就得到MO系数的零级近似C(0)。用C(0)计算储F(1)；再代入(2.44)方程确定新的系数C(1)，然后再计算F(2)等等。最后，当前后两次迭代所得系数C(n)与C(n-1)符合收敛精度时，迭代过程收敛。应用现代快速计算机易于实现这样的计算过程。开壳层体系的非限制性Hartree-Fock方法对于含有奇数电子的分子，不可能把所有电子均成对地排布于相应的分子轨道之中，体系中将有未配对的电子。此时，电子体系处于开壳层状态。显然，当电子从闭壳层分子基态占据分子轨道跃迁到基态未占据的空轨道上去，也会产生类似状态。图2-2处理开壳层的两种方法当描写开壳层体系时，波函数ψ一般应满足条件(2.13a)的一些Slater行列式线性组合(2.15)构成。从而，计算方案将更加复杂。然而，对于开壳层体系对应于极大多重度的状态来说，可以保持波函数的单行列式表示。非限制Hartree-Fock(UHF)方法[4,9]是描写这类体系的可能方法之一。UHF方法的基本假定是，α自旋电子所处的分子轨道不同于β电子。从而，与闭壳层体系不同，在UHF法中引入两组分子轨道：p个α自旋电子置于分子轨道集合ψiα中，而q个β自旋电子置于分子轨道集合ψiβ中。电子体系的波函数为(2.49)由于α自旋和β自旋电子占据不同的空间轨道，所以应用ψUHF就在每种程度上考虑了不同自旋电子的相关效应。常用的还有另一种处理开壳层体系的限制型Hartree-Fock(RHF)方法。两种不同方案的图像表示如图2.2.所示。把(2.49)代入(2.16)式中，可以求出ψUHF所描写的体系的能量表达式(2.50)其中交换积分Kijα和Kijβ分别用分子轨道ψiα和ψiβ计算。当ψα＝ψβ和p=q时，(2,50)式自动还原为(2.21)式，而ψUHF也就变为闭壳层体系的波函数(2.14)式。如前所述，对于分子轨道ψiα和ψiβ分别向原子轨道φμ做LCAO展开有(2.51)此时，展开系数cα不同于cβ。ψiα和ψiβ各自满足正交归一化条件，因为它们不同的自旋因子保证了(2.49)式中的ψiα与ψiβ的正交性。可以分别引入α自旋和β自旋电子的原子轨道密度矩阵元(2.52)显然，总电子密度矩阵元等于两者的和，即(2.53)而两者之差定义了自旋密度矩阵元(2.54)将ψiα和ψiβ的展开式(2.51)代入(2.50)式中积分表达式，可得到用原子基AO表达的能量公式，并考虑到(2.52)式，E可以写为(2.55)相对于cμiα和cμiβ各自使用变分法极小化能量(2.55)，就导出联立的两套方程组，可计算分轨道ψiα和ψiβ的能量εiα和εiβ及系数cμiα和cμiβ。(2.56)上式中α自旋和β自旋电子的Hartree-Fock算符矩阵元为(2.57)或者写成下列矩阵公式(2.58)正如闭壳层情况一样，方程(2.56)也是系数cα和cβ的三次联立方程组，只能用迭代法求解。UHF方法的缺点是，波函数ψUHF满足条件(2.13b)，Ms=(p-q)/2，但不是S2的本征函数，即不对应于任一个总自旋值。在一般情况下，可以把ψUHF表成具有不同自旋多重度波函数的线性组合(2.59)就数值来说，展开式(2.59)中系数Cs+m很快变小，所以分离出第一个混合态S'=S+1常常足够了，可以借助于消灭算符［10］(2.60)用原始波函数ψUHF表达与波函数AS+1ψUHF有关的密度矩阵(2.52)、总密度矩阵(2.53)和自旋密度矩阵(2.54)式，它们是UHF方法中极重要的一些物理量，详细的表达式见引文［10］。开壳层体系的限制性Hartree-Fock方法假定ψiα≡ψiβ，就可以从ψUHF得到满足(2.13a)的分子电子波函数，这就导出一个Slater波函数，它的n1个轨道均为两个自旋反平行电子占据（闭壳层），而n2个轨道为自旋相同电子单占据（开壳层）。(2.61)其图像如图2.2所示，波函数ψUHF是算符S2和Sz的本征函数，同时描写极大多重度状态S=MS=n2/2。所谓限制性Hartree-Fock法就是用(2.61)型波函数作运算的。把(2.61)代入(2.16)式，作一些变换后，可得到RHF法的电子能量表达式(2.62)其中k和l表示闭壳层部分的分子轨道，而m和n表示开壳层的分子轨道。量f表示开壳层的占据程度；a和b是Roothaan常数[11]，它们取决于所研究电子体系的具体特征。例如，对半充满的开壳层体系，f=1/2，a=1，b=2。公式(2.62)中其余量的含义同2.2节。在LCAO近似下，此时闭壳层和开壳层分子轨道均可按原子轨道基集合φμ展开成(2.19)式。运用关系(2.28)-(2.30)就可用AO基改写(2.62)式。当使用矩阵表示时，可以写成(2.63)其中ν1=2和ν2=2f是闭壳层和开壳层的填充数；R1和R2分别是相应于闭壳层和开壳层部分的原子基的密度矩阵［见(2.33)式］。闭壳层和开壳层部分的电子相互作用矩阵分别是G1和G2，它们是(2.64)矩阵G和G'可以用库仑矩阵[(2.38)式]和交换作用矩阵[(2.39)]表示为(2.65)(2.66)矩阵(2.65)描述闭壳层电子相互作用，而矩阵(2.66)描述开壳层电子相互作用。用极小化(2.63)式，可以得到确定分子轨道展开成AO基的系数方程，当然这必须在分子轨道相互正交的条件下进行，分别对闭壳层和开壳层轨道的LCAO系数变分。通常这导出两个(2.43)型的Hartree-Fock矩阵，对闭壳层和开壳层分别有(2.67)还可以由不同壳层分子轨道正交条件导出两个(2.44)型的方程组。然而，正如Roothaan指出的［11］，对于所研究的函数，可以将两个(2.44)型的方程组统一成一个求本征值和本征向量的方程(2.68)当然，这需要按下述公式确定上式中的Hartree-Fock矩阵(2.69)矩阵Ri'=1-Ri(i=1,2)，R3'=1-R1-R2，此处1为单位对角矩阵。可以把(2.69)式写成与(2.63)式相应的形式(2.70)这里B=2aJ2-bK2；R是总密度矩阵R=R1+fR2。库仑积分和交换积分下角标1和2分别对应于闭壳层和开壳层，包含有(2.70)式定义的Hartree-Fock算符方程(2.68)的求解归结为逐次迭代法。矩阵R1由n1个闭壳层轨道计算，而R2由n2个开壳层轨道计算。波函数(2.61)式以及公式(2.68)-(2.70)不仅可用于描写极大自旋多重度状态，而且也可用于描写开壳层轨道部分填充的一些状态。此时，只须改变(2.70)式中a，b，f数值，也并不使问题复杂化。Roothaan论文中已给出用RHF方法描述分子状态的这类系数［11］。Pitaer等人推广了Roothaan的研究，除了给出原子和线型分子能量系数外，还给出了Td和Oh分子的能量系数，使该方程更一般化了［12］。SCF计算结果的分析通过SCF计算得到了收敛波函数|ψ>，就可以应用此波函数|ψ>分析计算常用的物理量，现将主要的物理量计算介绍如下。总能量体系的电子总能量E是Hamilton算符的平均值<ψ|H|ψ>，即(2.71)而体系的总能量ET是电子能量加上原子核－原子核排斥能量(2.72)这是最有意义的量，因为预期稳定分子的平衡几何构型在ET的极小点处。Koopmans定律电离势的理论计算式是提供说明分子光谱的有用工具。与半经验计算不同，从头算处理了全部电子，因此不必区分价电子和内壳层电子电离势。一般不去区分光电子光谱（PES）和学分析电子光谱(ESCA)。两者在理论上是统一处理的。Koopmans定律提供了电离势的理论上最简单近似[13]，它使亲态闭壳层体系的Hartree-Fock轨道能量εh的负值等于第h个电离势(2.73)所以，按Koopmans定律，若理论上有若干个电子占据分子轨道，则光电子光谱应当存在有相应个谱带，粗略地说，在多数情况下实验观察的确如此[14]。然而，定量地讲，方程并非很满意地列数了观测谱。这是由于Koopmans定律有三个内在缺陷：①Koopmans定律不管从分子到离子时轨道的自洽重排，亦即忽略了与这种重排有关的能量ΔR，通常叫做重新组合能或弛豫能。总之Koopmans定律认为从分子到离子所有占据的分子轨道保持不变，这是一种近似。②没有考虑从分子到离子的相关能变化ΔC。③假定分子和离子的的相对论能量相同。首先，相对论能量随着原子系数的增加而增大，并且对于内壳层来说是很大的数值。例如对内壳层电离势来说，CH4的相对论修正为0.1eV，对于Ne的修正为0.8eV，而对于Ar对修正则是14eV。由于从头计算不直接计算相对论修正。其次，当忽略相对论效应时，对于内壳层从头计算给出电离势的估计值数值偏低。从而不考虑相对论效应时，电离势应为(2.74)差值IPhKRC-IPhK是所谓的Koopmans缺欠，即应当考虑电子弛豫和相关能修正。当采用ΔSCF方法时，可以考虑电子重排能ΔR(2.75)其中ESCF是电离前分子的SCF总能量，而EhSCF是空穴态离子的SCF能量。这表明，电离势是绝热SCF计算的分子和离子的总能量之差，有时称之为竖直电离能。单电子算符的平均值诸如分子的偶极矩、四极矩、核位置的场梯度、反磁化率等可以采用单粒子算符的和描写，一般形式为(2.76)其中g(k)是任意一种只依赖于单电子坐标的算符。单粒子算符的平均值是(2.77)从而，除密度矩阵外，只须计算单电子积分<ν|g|μ>。下面以偶极矩为例加以说明。偶极矩的经典定义为(2.78)向量ri除的电荷为qi。偶极矩相应的量子力学定义为(2.79)其中第一项是电子（电荷为-1）对偶极矩的量子力学贡献，第二项是原子核的核电荷ZA的经典贡献项。电子偶极矩算符是单粒子算符。所以应当有(2.80)这是一个向量表达式，其分量公式可以写为(2.81)对于μy和μz也有类似的公式。为了计算偶极矩，除了密度矩阵外Pμν外只须计算偶极矩积分(2.82)至于与电荷密度有关的问题下面加以说明。集居数分析在分子轨道计算中，人们通常把分子中电荷密度定义为(2.83)它表示空间任意点发现电子的几率，常常画出通过分子中各种截面的电荷密度轮廓图[16]。但是，由于分子给定原子或原子核中的电子数不是唯一确定的，有时可以给出不同的集居数分析。因为(2.84)每个MO占有两个电子。把ψi(r)的展开式代入上式就有(2.85)其中(PS)μμ为联系于ψμ轨道的电子数，这称为Mulliken集居数分析[17]。对于给定原子A在所有轨道上的电子数为，从而A原子的净电荷数为(2.86)其中ZA是A原子核的核电荷数，求和取遍所有属于A原子的轨道ψμ。但是，上述定义(2.85)不是唯一确定的，因为TrAB=TrBA，就有(2.87)对于任意a成立。当a=1/2时，得到(2.88)其中P'μμ是对称正交化基组表象中的密度矩阵(2.89)(2.90)此时，P'的对角元构成Lowdin集居数分析[18](2.91)上述两种集居数分析都是计算和文献中常用的。参考文献唐敖庆等.量子化学,科学出版社,1981.徐光宪等.量子化学—基本原理和从头计算(中册),1985.SzaboA,UstlundNS.ModernQuantumChemistry：AnIntruodutiontoAdvancedElectronicStructureTheory,MacmillanPublishingCo,NewYork,1981.PopleA,BeveridgeDL.ApproximateMoleculerOrbitalTheory,McGraw-HillBookCompany，NewYork,1970.王志中,李向东.半经验分子轨道理论与实践,科学出版社,1981.HartreeD.ClaculationalofAtomicStructure,Wiley,1957.PauliW.Z.Phys.1925,31,765.RoothaanCCJ.Rev.Mod.Phys.1951,23,231.PopleJA,NesbetRK.J.Chem.Phys.1954,22,517.AmosT.J.Chem.Phys.1964,41,1773;1965,43,3670.RoothaanCCJ.Rev.Mod.Phys.1960,32,179.PitzerRM.O.S.U.-T.C.G.ReportNo.101.KoopmansT.Phsica.1933,1,104.TurnerDW,etal.MolecularPhotoelectronSpectroscopy,WileyIntescience,London,1973.CadeE,etal.J.Chem.Phys.1966,44,1973.沈尔中,王志中,江元生.分子科学与化学研究，1984,5,427.MullikenRS.J.Chem.Phys.1955,23,1833.LowdinPO.J.Chem.Phys.1953,21,374.2.2从头算方法量子化学从头计算（Abinitio）方法60年代在国际上开始流行起来，随着计算机技术的发展，到了70年代得到了广泛的研究和应用。所谓“从头计算法”是立足于以下三条基本近似：①采用非相对论量子理论，即从Schrodinger方程出发；②使用Born_Oppenheimer近似；③轨道近似，即单粒子函数近似。除了这些近似以外，不再借助于任何经验参数，所以称之为从头计算。从头计算实质上是以前面所述分子轨道理论为基础的，从HFR或UHF方程出发，适当地选取原子轨道的线性组合法去模拟分子轨道计算时，一旦选定了基函数后，对于HFR或UHF方程所涉及的全部单电子积分和双电子积分均严格地计算。一般来说，选取从头计算基函数有两种方法：Slater型轨道（STO），用它描述电子云分布优于其它函数，但计算积分时涉及一些无穷级数，十分复杂，所以只是早期对一些小分子实现了这类计算；选用Gauss型函数（GTF）作为基函数，使多中心积分容易计算，且有简明的解析表达式，便于在计算机上实现，所以得到了广泛的应用。使用Gauss函数作为基函数求解Schrodinger方程，最初是Boys在1950年提出的，到了70年代初才得以迅速发展。STO与分子积分STO定义量子化学中广泛应用的一类基函数是Slater轨道（STO）：(2.92)其中Slm(θφ)是归一化实球谐函数(2.93)式中Pl|m|是连属Legend多项式，当m=0时δm=0，而m≠0时δm=1，ζ称为轨道指数。(2.92)式存在多种确定轨道指数ζ的方法。Clementi-Remond曾按极小化孤立原子能量定出了原子的ζ值［19］。Burns以最佳重复Hartree-Fock原子轨道方式定出了另一套指数[20]Gauterman-Zener用Hartree-Fock轨道给出的金属-配位体重叠积分为标准定出一组ζ值[21]。但是最常用的还是Slater在1930年建议的计算ζ值的经验方法[22]，该方法是基于很好再现原子的前几个电离势。分子积分[23-35]使用STO为基函数，除了函数的表达式简单外，它们还有叠加性。此时所需要计算的分子积分有下列四类重叠积分(2.94)动能积分(2.95)核吸引能积分(2.96)电子相互作用能积分(2.97)上述式子中，A、B等为原子核的位置向量，而r1和r2是电子的位置向量。(2.94)-(2.97)是所谓的多中心积分，其中最难计算的是(2.96)和(2.97)式中的三中心和四中心积分。实际上，到目前为止，尚无有效地解析式，允许人们用Slater轨道快速计算多原子分子问题所需要的多中心积分。大部分的单中心和双中心积分是易于计算的。双中心积分JAACC是库仑积分，JABAB是交换积分，而JAAAB称为杂化积分。详细的积分公式和推导可参考文献[23-37]。2.2.2GTF及Gauss函数积分由于STO分子积分的复杂性，限制了直接使用其进行从头计算的应用。1950年Boys建议用Gauss函数为基函数计算分子积分[38]。GTF的定义及其性质径向Gauss函数定义为(2.98)显然，它与STO的主要差别在于e-αr2和e-αr，两者角度函数相同。另一种常用的笛卡尔GTF表达式为(2.99)其中归一化常数N为(2.100)GTF的重要特点是，两个中心的GTF相乘可以表成另一个中心GTF函数的线性组合，使得多中心积分易于计算，这便是使用GTF的优越性。(2.101)(2.101)式详细的证明和推导可参考文献[2]。GTF的积分公式应用(2.98)式表达研究分子的径向GTF积分，最简便的是所谓无节点Gauss函数（即n=1,l=m=0,1s函数）。这是Boys[38]和Shavitt[39]所完成的。实际计算中使用笛卡尔GTF积分公式，它是由Taketa-Huzinaga-Ohata在1964年推导出来的[40]。唐敖庆等著《量子化学》[1]、王志中编著的《现代量子化学计算方法》[41]等书中也有较详细的推导。重叠积分：(2.102)(2.103)动能积分：(2.104)电子核引积分：(2.105)(2.106)双电子相互作用积分：(2.107)(2.108)其中已定义(2.109)求和指标限制为(2.110)从所给出的公式可以看出，用笛卡尔Gauss函数计算分子积分是方便的，同时易于计算机程序化。2.2.3解HFR方程的自洽场（SCF）过程HFR方程的求解，实际上要求解Fock算符的本征矢{cμi}—分子轨道组合系数，以及对应的本征值εi—分子轨道能量。但Fock算符本身又是分子轨道组合系数的二次函数，故只能用迭代方法求解，其SCF过程如下：(1)给出分子的各个核坐标，以及核电荷、电子总数，并选定基组{φμ}。(2)计算分子积分Sμν，hμν以及<μν|λσ>。(3)求出X。(4)初始猜测：即给出一组分子轨道，因基组是固定的，实际是给出展开系数，从Fock矩阵表示式来看，也就是给出密度矩阵.最简单的情况是令P≡0，即F=h。(5)计算G矩阵，即P及双电子积分<μν|λσ>的计算。(6)形成F=h+G。(7)计算F′＝X+FX。(8)对角化F′得到c′和E。(9)计算c＝Xc′。(10)由c通过形成新的P矩阵。(11)决定是否终止计算，即新P是否与上次相同（在一定限度内），若不收敛，则回到第(5)步继续迭代。(12)如果收敛，用所得的结果c，P，F等，计算一些有用的量.。能量梯度和平衡几何构型能量梯度的解析表达式势能对核坐标的一阶导数称之为势能在该方向的梯度，势能面上的驻点（极大点、极小点和鞍点）都满足所有核在各个方向上所受的力为0，即这些点的势能梯度向量为0向量。得到势能的梯度值后，辅以一定的数学手段，就能够优化出势能面上驻点的几何构型，并求得分子波函数和能量。求解能量梯度的解析方法是优化平衡几何构型和过渡态构型的基础。下面给出能量梯度的解析表达式，为与分子轨道的标号i区别，用表示某原子核在某坐标的梯度（这里的体系是闭壳层的N电子体系）(2.111)其中，核排斥项用梯度方法优化平衡几何构型分子总能量的绝对极小值所对应的几何构型称为分子的平衡几何构型。有如下的三种方法(1)可变尺度法（VariableMetricMethod）(2.112)X为坐标向量，可以是直角坐标或内坐标；an为尺度因子，A为严格推导的能量二阶导数矩阵的矩阵。调整尺度因子使g达到所要求的收敛限，以求得平衡几何构型。(2)极端方案a.最陡下降法：An=1，为单位阵。b.Newdon-Raphson法：取＝1，An为严格的力常数矩阵。(3)Muttagh-Sargent法：介于a和b的折衷方法，An用某种经验公式得到。参考文献ClememtiE,RaimondD.C.J.Chem.Phys.1963,38,2886.BurnsG.J.Chem.Phys.1964,41,1521.GautermanM,ZernerM.Theor.Chim.Acta.1966,4,44.SlaterJC.Phys.Rev.1930,36,57.RoothaanCCJ,etal.Rev.Mod.Phys.1960,32,186.GuseinovII.J.Chem.Phys.1977,67,3837.AlperJS.J.Chem.Phys.1971,55,3770;3780.MussoGF,etal.J.Phys.B:AtomandMolecularPhys.1970,4,1415.GuseinovII.J.Phys.B:AtomandMolecularPhys.1970,4,1399.SilverstoneHJ,etal.J.Chem.Phys.1968,48,4108.SilverstoneH.J.J.Chem.Phys.1968,48,4098.SilverDM,etal.J.Chem.Phys.1969,49,4301.SalmanLS,etal.Int.Quant.Chem.1972,6,353.GuseinovII.J.Chem.Phys.1978,69,4990.GuseinovII.Phys.Rev.(A)1985,31,2851.GuseinovII.Int.Quant.Chem.Symp.1986,19,149.CondonI,ShortleyG.TheoryonAtomicApectra.Combridge,1935.BarnettMP,CoulsonCA.Phil.Trans.R.Soc.,(London),SeriesA,1951,243,221.HarissonFE,MichelHH.J.Phys.Chem.,1965,43,165.BoysSF.Proc.Roy.Soc.A.(London),1950,200,542.ShvittI.MethodsinComputationalPhycis,Vol.2,AcademicPress,NewYork,1963,91.TaketaH,HuzinagaS,O-ochataK.J.Phys.Soc.(Japan)1966,21,2313.王志中.现代量子化学计算方法,吉林大学出版社，1998，43.2.3半经验SCF-MO方法以Roothaan方程及其求解为核心的SCF分子轨道方法，是从纯计算的观点（从头算）讨论了分子轨道理论，近似的分子轨道方法是以该理论为骨架，从两个方面简化计算：①忽略一些难于计算的积分，如大量的双电子积分；②对剩下的积分，根据电子结构的实验资料估计一些积分值。尤其是第二方面的努力，不能仅仅看作是为了简化计算，而是不再试图根据量子力学原理直接推导分子的性质，主要致力于解释实验数据（现象）内部的相互关系，它在许多方面取得了惊人的效果。这种基于实验结果的相互关系与严格量子理论的结合是半经验的。2.3.1零微分重叠近似（ZDO）[42,43]LCAO自洽场分子轨道计算的最大困难是对双电子积分的处理，这包括两个方面的内容：一是积分本身计算困难，特别是多中心的双电子积分在数学上的处理较困难；二是双电子积分数目随基函数数目的增加呈指数形式增长。而大量的实践结果表明，许多这类电子排斥的数值接近于0，尤其是那些包含重叠分布，且的电子排斥积分。因此，在发展近似自洽场分子轨道方案中，一个有用的近似是系统地忽略数值一致的小的电子排斥积分。这一点通过零微分近似实现，零微分近似写为：（当时）。于是，假定包含重叠积分的电子排斥积分数值小到可以忽略，在零微分近似下其次，在分子轨道归一化中忽略相应的重叠积分，也就是但不忽略包含重叠积分的实积分。所讨论的各种等级的近似自洽场方法的最主要区别是在处理电子排斥积分中应用零微分重叠近似的程度。若将零微分近似应用于全部原子轨道，则闭壳层的Roothaan方程简化为(2.113)这实际上就是简化的矩阵形式Fc=cE。其中Fock矩阵元为其对角元()非对角元()这些近似大大简化了计算，主要是忽略了许多困难的双电子积分，尤其是全部三中心和四中心积分都变为0。2.3.2全略微分重叠（CNDO）[5，44-45]CNDO明确地只处理价电子，内壳层被看成刚性实，于是原子轨道基集合是一个价基集合。其基本近似是对所有不同原子轨道积运用零微分重叠近似，即使用上面推出的简化方程Fc=cE和Fock矩阵元公式；此外对剩下的部分积分用实验参数代替。概括起来，CNDO有以下五步近似：(1)将Roothaan方程中重叠矩阵用单位阵代替，即在分子轨道归一化中忽略重叠积分；(2)忽略所有双电子积分中的微分重叠，即；(3)将剩下的库仑型积分均简化为每一对原子一个数值（参数）：，属于A原子，属于B原子；(4)在其它包含原子实的积分中忽略单原子微分重叠：，这里是B原子核及内层电子（刚性实）产生的势能，、属于A原子，也用参数化处理；(5)取双原子非对角实矩阵元正比于相应的重叠积分：，属于A原子，属于B原子，也是参数化处理。使用所有这些近似，使计算更简单，且参数化有不同的方案，这里不再详细介绍。2.3.3间略微分重叠（INDO）[46,47]CNDO对所有不同的原子轨道运用零微分重叠近似，固然对双电子排斥积分的处理是简单的，但它对两个具有平行自旋的电子间实际存在的不同相互作用未给予适当考虑，尤其是两个电子同属一个原子的情况。波函数的反对称性要求自旋平行的电子不占有空间的同一小区域，因此同一原子上不同原子轨道中的两个电子如果自旋平行，将具有较小的排斥能，这种差别在数学上表现于双电子交换积分：而在CNDO理论中，这种积分是被忽略的。这样，原子的电子能量与电子自旋取向无关，自旋多重度无法区分，导致CNDO不能描述对应于同一组态的电子能量分裂。为克服上述缺点，最简单的方法就是对交换积分作某些考虑，即只保留单中心积分中的单原子微分重叠。所谓单中心积分，是积分只涉及同一个原子；而单原子微分重叠，指。例如，双电子积分中的<SAPA|SAPA>要予以保留，而<SAPA|SAPA>虽属单原子微分重叠，但非单中心积分而予以忽略。这种微分近似方法称为间略微分重叠（INDO）。2.3.4忽略双原子微分重叠(NDDO)[48]下一等级近似的自洽场理论，是只忽略属于不同原子的原子轨道的微分重叠，即在化简HFR方程时，对于和属于不同原子的原子轨道才使用零微分重叠近似。这一等级近似的主要特点是保留了偶极－偶极相互作用，因为此时保留了型积分。NDDO法虽然也忽略了全部三中心和四中心积分，但是所保留的单中心和双中心积分数目仍很多，过渡到半经验水平（参量化）比较困难。2.3.5NDDO的改进[49-59]历史上，NDDO法对双中心积分采用Slater轨道直接计算，程序很繁琐，只能计算一些不大的分子。Dewar等人于1977年建议了改进的NDDO法，即MNDO法[54]，主要特点是用半经验模型处理双中心电子排斥积分。1985年，Dewar等人又对MNDO进一步修正，主要针对MNDO高估了原子间的排斥作用，对MNDO中的核实排斥函数进行了修正，加进了一些高斯函数项，同时进行重新参数化，即AM1方法(AustinModel1)[57]。这种方法目前仍被广泛应用。1989年，Stewart在AM1的理论框架下，发表了新的参数化方法并进行了再参数化工作，称为PM3(ParamiterMethod3)[59]。到目前为止，MNDO、AM1、PM3是使用最为广泛的半经验方法。MNDO和AM1方法改进的重点是修改原子实间的排斥项。从静电的角度看，原子实-原子实排斥积分为（eV）(2.114)然而，在范德华半径以外，电子-电子和电子-原子核实的作用并不按照1/RAB的规律衰减，如果使用(2.114)的形式，将造成两个中性原子或分子之间的净排斥。为了矫正这一点，减小过高的排斥能，作如下近似MNDO对核实-核实相互作用进行了修正，采用如下形式(2.115)但是，对于O-H和N-H间的相互作用采用下式单独处理(2.116)AM1在此基础上作了进一步的修改(2.117)其中比MNDO多出的球面Gaussian函数项的作用是，在成键的距离之外进一步减小多余的核实-核实排斥作用，它可以理解为是vanderWaals吸引项。(2.117)式中的a,b,c是可调参数。参考文献PariserR,ParrRG.J.Chem.Phys.1953,21,767.PopleJA.Trans.Farady.Soc.1953,49,1375.波普尔J.A.等著,江元生译,近似分子轨道理论方法,科学出版社,1976.ClackDW.J.Chem.Soc.(A)1970,2688;Ibid,1971,299;2783.PopleJA,etal.,J.Chem.Phys.1967,47,2026.DixonRN.Molec.Phys.1967,12,83.PopleJA,SegalGA.J.Chem.Phys.1965,43,S136.DewarMJS,etal.J.Chem.Phys.1969,50,1260;1275.DewarMJS,etal.J.Am.Chem.Soc.1969,91,352.DewarMJS,etal.J.Am.Chem.Soc.1970,92,590;3854.DewarMJS,etal.J.Am.Chem.Soc.1972,94,5296.DewarMJS,etal.J.Am.Chem.Soc.1975,97,1285-1311.DewarMJS,ThielW.J.Am.Chem.Soc.1977,99,4499.DavisLP,etal.J.Comp.Chem.1981,2,433.DewarMJS,McKeeML,RzepaHS.J.Am.Chem.Soc.1978,100,3607.DewarMJS,ZoebischEG,HealyEF."AM1:ANewGeneralPurposeQuantumMechanicalMolecularModel,"J.Am.Chem.Soc.1985,107,3902.StewartJJP.J.Comp.Chem.1989,10,209.StewartJJP.J.Comp.Chem.1989,10,209;221.2.4DFT方法理论简介[60-79]自上个世纪80年代末到90年代（Parr1983；Wimmer1997）开始，密度泛函理论（Densityfunctionaltheory,DFT）在处理分子和原子的电子结构问题变得非常热门。DFT的主要思想是考虑总电子能量和总电子密度之间的关系。这一思想的近似模型（TheThomas-Fermimodel）早在上个世纪20年代就已提出。但是真正的突破性工作来自Hohengerg和Kohn于1964发表的论文，以及Slater和Sham等人的改进，他们指出体系基态的能量和其他性质是由电子密度唯一确定的。DFT方法的主要观点是通过电子密度函数来描述相互作用的费米子体系，而不是多体波函数。在HF理论中，多电子体系的波函数是由Slater行列式来表达的，它由N个单电子波函数（N是分子中电子的数目）构成，通常含有N个电子的体系包含3N个变量。DFT也考虑单电子函数。然而，HF理论是全部计算N个电子的波函数，DFT却只计算总的电子能量和所有电子密度的布居。无论体系多大，对于电子密度来讲只有三个变量x、y、z，这样可以大大简化计算。DFT现在已经广泛地应用于化学热力学、分子光谱等领域，可以用来描述金属、半导体、绝缘体的基态性质，和分子轨道理论一起，成为人们探讨分子体系的有力理论工具。它的成功之处在于可以成功的解释标准的大体积宏观物质，而且它也能解释诸如蛋白质、碳纳米管等复杂的物质。DFT理论的优点还在于它提供了许多化学上重要的概念，例如电负性、硬度、Fukui函数、相应函数等，这些概念可以很方便的用来解释化学性质和分子的变化，它在预言原子、分子和固体基态的能量、电子自旋密度、键长和键角等方面已经取得了初步的成功。2.4.1Thomas-Fermi模型1927年的Thomas-Fermi模型是最早的尝试用电子密度函数来代替体系的波函数描述和确定体系的性质。它的基本设想是把空间看成许多边长为1的小立方体，每个立方体都是一个三维势箱，则有如下公式从上述公式中可以看出，以坐标系原点为球心、R为半径的球面上的()点为标志的量子态（此时，都是整数）具有相同的能量，称为简并态。电子按能量沿着球面由里向外、由低到高依次填充，但只有都为正的球体中的点才代表有电子占据的能级，由此可知能量小于（或等于）ε的量子态数目为能量在ε和ε＋∆ε间的能级数目为其中，是能级ε处的态密度。在0K时，能量在费米能级（离子的最高占据能级）以下的粒子组成的体系总能量为(2.118)假定每个()量子态有两个粒子占据（自旋方向相反），则离子总数为即可得到因为体系V＝0，故ΔE是体系的动能，令（电子密度），对空间所有小立方体求和（积分），使用原子单位，则得到，(2.119)这就是动能的电子密度的泛函表达式。对于原子，只考虑电子间的库仑作用，有，满足用拉格朗日乘子法，求和，令为拉格朗日乘子，由此得Euler-Lagrange方程其中(2.120)上式与条件联立，可求解得到和。实际计算表明，分子的能量总大于其组成原子能量之和，即得不到稳定的分子体系。因此，虽然化学家用电子密度分布概念解释了很多现象，但实际结果似乎表明电子密度变化的理论缺乏物理学依据。2.4.2Hohenberg-Kohn定理现代密度泛函理论开始于Hohenberg-Kohn（1964年）定理，它指出电子体系的基态电荷密度分布唯一的决定体系总能量，可以用公式（是体系基态的密度），包括第一定理和第二定理。第一定理：体系的基态的电子密度与体系所处势场有一一对应的关系，从而完全确定体系的所有性质。而如果两个体系具有相同的基态电荷密度分布，那么它们的势场只相差一个常数确定。该定理初步论证了用电子密度函数代替波函数来全面描述一个体系是一个可行的方法。可以用反证法进行证明，假设两个体系分别处于势场和，并且它们的差值不是常数，它们的哈密顿算符可以如下表示上式中是电子体系的哈密顿算符。设和的非简并态波函数分别为Ψ1和Ψ2，相应的电子密度为和，则且小于依据变分原理，有即同理可得如果，则两个不等式相加就有，这是不可能的，因此。也就是当常数时，，这就证明了电子密度和体系所处势场之间有一一对应关系。当体系电荷密度分布确定时，势场就确定了，电子数目也确定了，所以体系的哈密顿也确定了，进而可知体系的波函数和所有其它性质都被确定了。也就是说，体系总能量、动能、电子间相互作用能都是的泛函，记作、、，并有其中，与外势场没有明显的关系。第二定理：对于任意的函数，如果满足，，则，为体系基态能量。这一定理说明了体系基态总能量是总能量密度泛函的极小值，即，可以简单地表示为（是体系基态的密度，是体系处在非基态时的电子密度）。该定理证明如下：对应一个、和Ψ′，设Ψ′为试用变分函数，则。按拉格朗日乘子法变分，μ为拉格朗日乘子，得到Euler-Lagrange方程如果知道的具体表达式，结合的限制条件，可由上式求得及。在这里我们采用Levy提出的一种求的搜索方法即在所有给出相同的的中，使取极小值者给出的值。于是用变分法求基态能量的过程可表示为(2.121)2.4.3Kohn－Sham方程Kohn－Sham方程主要是通过把动能泛函的一部分分离出来而得到的。假设不考虑自旋，而且粒子是相互独立的，那么动能可表示为表示为其中，是库仑作用能，是交换相关能，包含T和TS的差，是一个量值很小的泛函。可表示为对变分求能量极值得到Kohn－Sham方程其中，，若已知，求解Kohn－Sham方程，可得和从Kohn－Sham方程的形式来看，它与Hatree-Fock方程非常类似，而且物理意义上也类似，所以虽然不具有分子轨道能级的含义，通常还是把看作近似的分子轨道。2.4.4Kohn－Sham方程的求解方法DFT分子计算的过程实际上就是用变分法求解分子的Kohn-Sham方程。由Kohn-Sham方程和Hatree-Fock-Roothaan方程类似的特点可知，方程的求解也要用自洽迭代方法。但DFT的计算量较小，要用数值积分方法计算，不能采用解析法计算。尽管现在比较流行的数值积分方法已经有求积方法[77]、分区求积方法[78]和重迭区间空间分割法[79]三种方法，但是由于矩阵元的数值计算占用密度泛函方法计算的大部分资源，因此寻找高效高精度的数值积分方法仍然是一个很受重视的课题。2.4.5交换能和相关能校正均匀电子云的交换能可表示为，式中第一项为交换能，第二项为相关能，参数，γ代表自旋，表示在上自旋和下自旋密度分别为和的电子云中单个电子所具有的相关能。的精确解析形式至今还不知道，Stoll等人认为LDA的主要误差之一是均匀电子云模型过高估计了自旋相同的电子之间的相关能。真实体系的电荷分布是不均匀的，定域密度近似忽略了电荷不均匀性的影响，因此计算的精度也就受到了限制。本文的计算大多采用的BeckeIII的交换能校正公式：(2.122)其中，，，。一般来说，LDA过高估计相关能，与交换能的处理类似，人们引入密度梯度校正来减小LDA相关能的误差。各种近似公式的计算结果差别不是很大，本文采用的LYP公式为：其中，，，a=0.04918，b=0.132，c=0.349。参考文献SalahubDR,ZernerMC.TheChallengeofdandfElectrons,ACS,Washington,D.C.,1989.ParrRG,YangW.Density-functionaltheoryofatomsandmolecules,OxfordUniv.Press,Oxford,1989.HohenbergP,KohnW.InhomogeneousElectronGas,Phys.Rev.B.136,864,1964.KohnW,ShamLJ.Self-ConsistentEquationsIncludingExchangeandCorrelationEffects,Phys.Rev.A1965,140,1133.SlaterJC.QuantumTheoryofMolecularandSolids.Vol.4:TheSelf-ConsistentFieldforMolecularandSolids,McGraw-Hill,NewYork,1974.PopleJA,GillPMW,JohnsonBG."Kohn-Shamdensity-functionaltheorywithinafinitebasisset,"ChemicalPhysicsLetters,1992,199,557.LeeC,YangW,ParrRG.DevelopmentoftheColle-Salvetticorrelation-energyformulaintoafunctionaloftheelectrondensity,PhysicalReviewB1988,37,785.BeckeAD.Phys.Rev.A1988,38,3098.BeckeAD.Density-functionalthermochemistry.III.Theroleofexactexchange,J.Chem.Phys.1993,98,5648.VoskoSH,WilkL,NusairM.Accuratespin-dependentelectronliquidcorrelationenergiesforlocalspindensitycalculations:acriticalanalysis,,CanadianJ.Phys.1980,58,1200.MiehlichB,SavinA,StollH,etal.Chem.Phys.Lett.1989,157,200.JohnsonBG,FrischMJ.Animplementationofanalyticsecondderivativesofthegradient-correcteddensityfunctionalenergy,J.Chem.Phys.1994,100,7429.JohnsonBG,FrischMJ.Analyticsecondderivativesofthegradient-correcteddensityfunctionalenergy:Effectofquadratureweightderivatives,Chem.Phys.Lett.1993,216,133.StratmannRE,BurantJC,ScuseriaGE,etal.J.Chem.Phys.1997,106,10175.PerdewJP.Density-functionalapproximationforthecorrelationenergyoftheinhomogeneouselectrongas,Phys.Rev.B1986,33,8822.BurkeK,PerdewJP,WangY.inElectronicDensityFunctionalTheory:RecentProgressandNewDirections,Ed.DobsonJF,VignaleG.andDasMP.(Plenum,1998).PerdewJP,BurkeK,WangY.Phys.Rev.B1996,54,16533.RosenA,EllisDE.J.Chem.Phys.1976,65,3629.BoerrigterPM,IdeGT,BarendsEJ.Int.J.QuantumChem.1988,33,87.BecksAD.J.Chem.Phys.1988,88,2547.2.5势能面、平衡几何构型与内禀反应坐标研究一个化学反应的反应途径，需研究这个反应体系的势能面，由于对整个反应势能面的研究不易做到，因而一般只研究势能面上驻点的性质。驻点是指反应途径上的能量的极小点和鞍点