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

...wd......wd......wd...第二章計算方法及其根本原理介紹化學(xué)反響的本質(zhì)是舊鍵的斷裂和新建的形成，參與成鍵原子的電子殼層重新組合是導(dǎo)致生成穩(wěn)定多原子化學(xué)鍵的明顯特征。因此闡述化學(xué)鍵的理論應(yīng)當(dāng)描寫電子殼層的相互作用與重排，借助求解滿足適當(dāng)?shù)腟chrodinger方程的波函數(shù)描寫分子中電子分布的量子力學(xué)，為解決這一問題提供了一般的方法，然而，對于一些實際的體系，不引入一些近似，就不可能求解其Schrodinger方程。這些近似使一般量子力學(xué)方程簡化為現(xiàn)代電子計算機可以求解的方程。這些近似和關(guān)于分子波函數(shù)的方程形成計算量子化學(xué)的數(shù)學(xué)根基。圖2-1分子體系的坐標2.1SCF-MO方法的根本原理分子軌道的自洽場計算方法(SCF-MO)是各種計算方法的理論根基和核心局部，因此在介紹本文計算工作所用方法之前，有必要對其關(guān)鍵的局部作一簡要闡述。2.1.1Schrodinger方程及一些根本近似為了后面介紹各種具體在自洽場分子軌道(SCFMO)方法方便，這里將主要說明用于本文量子化學(xué)計算的一些重要的根本近似，給出SCFMO方法的一些根本方程，并對這些方程作簡單說明，因為在大量的文獻和教材中對這些方程已有系統(tǒng)的推導(dǎo)和闡述[1-5]。確定任何一個分子的可能穩(wěn)定狀態(tài)的電子構(gòu)造和性質(zhì)，在非相對論近似下，須求解定態(tài)Schrodinger方程〔2.1〕其中分子波函數(shù)依賴于電子和原子核的坐標，Hamilton算符包含了電子p的動能和電子p與q的靜電排斥算符，(2.2)以及原子核的動能(2.3)和電子與核的相互作用及核排斥能(2.4)式中ZA和MA是原子核A的電荷和質(zhì)量，rpq=|rp-rq|，rpA=|rp-RA|和RAB=|RA-RB|分別是電子p和q、核A和電子p及核A和B間的距離〔均以原子單位表示之〕。上述分子坐標系如圖2.1所示?？梢杂肰(R,r)代表(2.2)－(2.4)式中所有位能項之和(2.5)原子單位上述的Schrodinger方程和Hamilton算符是以原子單位表示的，這樣表示的優(yōu)點在于簡化書寫型式和防止不必要的常數(shù)重復(fù)計算。在原子單位的表示中，長度的原子單位是Bohr半徑能量是以Hartree為單位，它定義為相距1Bohr的兩個電子間的庫侖排斥作用能質(zhì)量則以電子制單位表示之，即定義me=1。Born-Oppenheimer近似可以把分子的Schrodinger方程(2.1)改寫為如下形式(2.6)由于組成分子的原子核質(zhì)量比電子質(zhì)量大103倍～105倍，因而分子中電子運動的速度比原子核快得多，核運動平均速度比電子小千倍，從而在求解電子運動問題時允許把電子運動獨立于核運動，即認為原子核的運動不影響電子狀態(tài)。這就是求解(2.1)式的第一個近似，被稱作Born-Oppenheimer近似或絕熱近似。假定分子的波函數(shù)Ψ′可以確定為電子運動和核運動波函數(shù)的乘積(2.7)其中Ф(R)只與核坐標有關(guān)，代入方程(2.2)有對于通常的分子，依據(jù)Born-Oppenheimer原理有：AΨ和A2Ψ都很小，同時MA≈103～105，從而上述方程中的第二項和第三項可以略去，于是易知也即該方程可以別離變量而成為兩個方程(2.8)(2.9)方程(2.8)為在某種固定核位置時電子體系運動方程，而方程(2.8)時核的運動方程。E(R)固定核時體系的電子能量，但在核運動方程中它又是核運動的位能。此時分子總能量用ET代表。因此，在Born-Oppenheimer近似下，分子體系波函數(shù)為兩個波函數(shù)的乘積(2.7)式。分子中電子運動波函數(shù)Ф(R)分別由(2.8)和(2.9)式確定。電子能量E(R)為分子的核坐標的函數(shù)，從(2.9)式看出它又是核運動的位能。在空間畫出E(R)隨R的變化關(guān)系稱為位能圖。單電子近似體系的電子與核運動別離后，計算分子的電子波函數(shù)Ψ歸結(jié)為求解下面的方程(2.10)(2.10)式是量子化學(xué)的根本方程，目前已有多種求解這個方程的方法。這些方法的區(qū)別首先是構(gòu)成Ψ的方式及其相應(yīng)的近似。最常用的是Hartree建議的單電子近似［6］。在多電子體系中，所有電子勢相互作用的，其中任意電子運動依賴于其它電子的運動。Hartree建議把所有電子對于每個個別電子運動的影響代換成某中有效場的作用。于是每個電子在核電荷及其余電子有效場產(chǎn)生的勢場中運動僅依賴于電子坐標。從而，電子運動分開了，對于多電子體系中每個電子可以引入單電子波函數(shù)，這種單電子波函數(shù)是(2.10)式單電子Schrodinger方程的解，其中含有算符1/rpq項，用只依賴于所研究電子坐標的有效場代替。整個多電子體系波函數(shù)等于所有電子的單電子波函數(shù)〔軌道〕乘積。電子還具有自旋角動量s，其分量sx,sy和sz滿足普通角動量算符的對易關(guān)系。算符s2和sz完全給定了電子的自旋，電子自旋波函數(shù)()滿足方程(2.11)其中是自旋坐標，通常把對應(yīng)于自旋1/2的波函數(shù)記為()，而把自旋ms=-1/2波函數(shù)記作()。在非相對論近似下和不存在外磁場時，電子的自旋和空間坐標無關(guān)，因此，因此電子的自旋軌道可取成(2.12)考慮到自旋變量的多電子波函數(shù)由自旋軌道組成，他應(yīng)當(dāng)是體系總自旋S2及其Sz的本征函數(shù)(2.13a)(2.13b)構(gòu)成體系多電子波函數(shù)Ψ時，必須考慮Ψ相對于任一對電子交換的反對稱性要求，此所謂Pauli原理［7］。因此，一般不求出Hartree方法的簡單乘積型波函數(shù)Ψ，而是求出對應(yīng)于按自旋軌道電子的所有可能置換方式的Slater行列式波函數(shù)，此為Hartree-Fock方法。對于置于n=N/2軌道的Ψ上的N電子體系，單電子近似下波函數(shù)Ψ寫為(2.14)該式的Slater行列式是保證反對稱性要求的唯一這類函數(shù)。引入單電子近似便確定了波函數(shù)Ψ的形式，用它可以求解方程(2.10)。顯然在一般的情況下，Ψ應(yīng)當(dāng)包含(2.14)型行列式的線性組合，同時滿足(2.13)式的限制。假設(shè)(2.12)式中自旋局部是單電子自旋投影算符Sz的本征值，則(2.13b)式就滿足。當(dāng)分子的n個軌道每個均為自旋反平行電子對占據(jù)時〔閉電子殼層〕，一個行列式波函數(shù)(2.14)就已滿足(2.13a)和(2.13b)。對于含有未配對的電子體系，這是做不到的，此時體系波函數(shù)是對應(yīng)于各種軌道填充方式〔不同組態(tài)〕的Slater行列式Ψl的的線性組合(2.15)當(dāng)適中選擇行列式前系數(shù)al時，條件(2.13a)和波函數(shù)的反對稱性要求均可以滿足。由于存在著電子運動的相關(guān)，不明顯處理(2.10)式中1/rpq項的單電子近似，完全忽略了這種相關(guān)效應(yīng)，所以，Hartree-Fock單電子近似使波函數(shù)的計算產(chǎn)生了誤差。變分原理上述單電子近似只是給出了所求解體系多電子波函數(shù)的一種形式，變分法提供了求解方程(2.10)的一種方法。Schrodinger方程(2.10)的解對應(yīng)于穩(wěn)定態(tài)能量。因此假設(shè)波函數(shù)Ψ是(2.10).的解，那么對于任意微小變化Ψ，取能量平均值(2.16)的變分應(yīng)等于零，即(2.17)(2.16)式中積分是對Ψ的所有變量進展的，并且已假定Ψ是歸一化的，即(2.18)由于我們尋找對應(yīng)于體系基態(tài)的波函數(shù)，總能量應(yīng)當(dāng)是極小值。因此，對單電子軌道施行變分就給出這種型式波函數(shù)，能量是極小值并滿足(2.17)式。從而求得的波函數(shù)Ψ就是多電子體系基態(tài)Schrodinger方程所欲求的解。顯然，為了施行變分，波函數(shù)Ψ的型式應(yīng)當(dāng)充分好。兩種途徑可以保證這一點：①取展開式(2.15)是從充分多項，且固定軌道Ψ只對系數(shù)al變分；②局限于盡可能少的行列式Ψl，假設(shè)有可能做到就取一個，但此時把每個ψ表成可能的簡單形式。鑒于這種選擇，區(qū)分出兩類廣泛應(yīng)用的量子化學(xué)方法，價鍵(VB)法和分子軌道法(MO).在價鍵法中，用孤立原子的原子軌道(AO)作為單電子波函數(shù)ψ去構(gòu)成Slater行列式Ψl。原子軌道的不同選擇對應(yīng)于不同的行列式Ψl。對于(2.15)施行的變分，可得到確定系數(shù)al的方程。為了充分靠近體系的能量，必須在(2.15)式中選用足夠的多項，即用多行列式波函數(shù)進展運算。用原子軌道線性組合分子軌道(LCAOMO)法提供了另外一種選擇相應(yīng)于體系能量極小的多電子波函數(shù)方法。此時，對應(yīng)于分子中單電子態(tài)的分子軌道ψi寫成原子軌道φμ〔基函數(shù)AO〕的線性組合(2.19)實際上，這種展開有完全合理的根基。因為靠近某個原子的電子所受的作用根本上是由該原子產(chǎn)生的場引起的，所以該區(qū)域中電子波函數(shù)應(yīng)當(dāng)近于原子軌道。展開該式對求解變分問題的優(yōu)點是明顯的。如果(2.15)式中選用極大數(shù)目的項，那么VB法和MO法就都給出同樣的能量E和波函數(shù)ψ，當(dāng)然表達式不完全一樣。這種唯一性的原因很簡單，因為使用LCAOMO的的每個行列式均可以展開為AO組成的一些行列式。在一般情況下，每個MO組成的行列式應(yīng)展開成AO組成的所有行列式。因為波函數(shù)ψ應(yīng)通過AO組成的行列式完全集合表達，從而，當(dāng)使用完全集合時，MO法與VB法所描述的ψ就等價。當(dāng)然，不用完全及表達時，兩種方法的等價性就破壞了。在極端性況下，某種方法中可以取一個行列式，此時可以直接看到MO法的優(yōu)越性。對于MO法，允許采用單行列式表達ψ〔至少對于閉殼層體系〕，進而，通常由一些正交分子軌道組成行列式(2.20)其中δij是Kronecker符號。從而是計算大為簡化，并能比VB法更簡單地確定(2.19)式的方程系數(shù)。同時，MO法的根本方程能很好的適應(yīng)現(xiàn)代電子計算機的能力。由于這個原因，現(xiàn)代的MO方法已經(jīng)成為最常用的計算多電子分子的電子構(gòu)造的根本方法。閉殼層體系的Hartree-Fock-Roothaan方程在分子軌道范圍內(nèi)，對閉殼層體系，在單電子近似下，用兩個自旋反平行電子填充每個分子軌道ψ，可以構(gòu)成一個Slater行列式(2.14)型波函數(shù)，選擇軌道(2.12)的自旋局部滿足(2.11)式，則保證了(1.13b)條件。根據(jù)變分原理，假設(shè)軌道ψ使得分子能量(2.16)取極小值，就求出了所研究多電子體系方程(2.10)的解。將波函數(shù)(2.14)代入(2.16)式，并進展一些推導(dǎo)［見引文1－4］，可得閉殼層分子的電子能量表達式(2.21)此處Hii是對應(yīng)于分子軌道ψ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軌道上占據(jù)兩個電子，所以乘以2。雙電子庫侖積分Jij表示ψi和ψj軌道上兩個電子間平均排斥作用。由于波函數(shù)的反對稱性要求出現(xiàn)了交換積分Kij〔在Hartree方法中不考慮它〕，減小了不同軌道ψi和ψj上平行自旋電子間相互作用，正是它們描寫了一樣自旋電子運動的交換相關(guān)。然而，在Hartree-Fock方法中，還是沒有考慮反平行自旋電子間庫侖排斥引起的電子相關(guān)效應(yīng)。為了求解Ψ的最優(yōu)近似，必須選擇一定形式的分子軌道ψi是總能量最小。這些分子軌道相互正交，在LCAO近似下表成原子軌道ψμ的展開(2.19)式。Roothaan最先解決了這一問題［8］。關(guān)于ψi對AO基展開系數(shù)的方程成為Hartree-Fock-Roothaan方程〔間記為HFR方程〕。下面簡要推導(dǎo)這個方程。首先從寫ψi的LCAO展開式分子軌道正交歸一化條件給出對MO系數(shù)的附加限制(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)變?yōu)?2.32)引入原子軌道基的電子密度矩陣元(2.33)上式化為(2.34)或者令(2.35)則可以用用矩陣形式寫為(2.36)其中電子相互作用矩陣G定義為(2.37)矩陣J(R)描寫庫侖相互作用，而K(R)描寫電子交換相互作用，其矩陣元分別為(2.38)(2.39)作用于兩個矩陣的運算Sp是指把此兩個矩陣所響應(yīng)矩陣元的乘積求和。在LCAO分子軌道法中，通常假定基原子軌道是固定的，而型式不變，因此對分子軌道變分歸結(jié)為對展開系數(shù)cμν的變分(2.40)應(yīng)當(dāng)指出，(2.40)式中的ψμ不變性在原則上是不必要的，已有把ψμ看成可變的一類計算方法，然而是極其困難和復(fù)雜的，常用的MO法一般不做這類計算。當(dāng)對分子軌道變分時，借助Lagrange乘法可以是能量極小化；這是極小化泛函即變化系數(shù)cμi，使穩(wěn)定點δG＝0。其中E由方程(2.34)定義。取G的變分，則有(2.41)由于δcμi*的任意性，所以必須有(2.42)因為分子軌道在酉變換下確實定性，可以自由選擇非對角Lagrange乘數(shù)為零。定義Fock矩陣元(2.43a)或記為矩陣形式(2.43b)其中H為Hartree-Fock矩陣的單電子局部，而G為雙電子局部。于是方程(2.42)可以寫為(2.44a)或記為矩陣形式(2.44b)此為HFR方程，其中E是Hartree-Fock算符本征值組成的對角矩陣，S是基原子軌道的重疊積分矩陣。解(2.44a)或(2.44b)，就給出按原子軌道展開的分子軌道系數(shù)和單電子分子軌道能量εi。然而，上述HFR方程是數(shù)學(xué)上廣義特征值問題。為了求解方便，先借助于下述變換將它化成標準特征值問題。由于S是Hermite矩陣，可以用酉變換θ使其對角化，即(2.45)從而可由其對角元平方根倒數(shù)構(gòu)成矩陣S-1/2，同時S1/2S-1/2和S-1/2SS-1/2=1成立。此時可做變換(2.46)(2.47)則方程(2.44)化成亦即(2.48)這是標準特征值問題，可以采用數(shù)學(xué)上標準的對角化方法處理。顯然這些方程都是系數(shù)cμi三次聯(lián)立方程組，必須用迭代法求解：當(dāng)F(0)=H時，解上述方程(2.44a)或(2.44b)就得到MO系數(shù)的零級近似C(0)。用C(0)計算儲F(1)；再代入(2.44)方程確定新的系數(shù)C(1)，然后再計算F(2)等等。最后，當(dāng)前后兩次迭代所得系數(shù)C(n)與C(n-1)符合收斂精度時，迭代過程收斂。應(yīng)用現(xiàn)代快速計算機易于實現(xiàn)這樣的計算過程。開殼層體系的非限制性Hartree-Fock方法對于含有奇數(shù)電子的分子，不可能把所有電子均成對地排布于相應(yīng)的分子軌道之中，體系中將有未配對的電子。此時，電子體系處于開殼層狀態(tài)。顯然，當(dāng)電子從閉殼層分子基態(tài)占據(jù)分子軌道躍遷到基態(tài)未占據(jù)的空軌道上去，也會產(chǎn)生類似狀態(tài)。圖2-2處理開殼層的兩種方法當(dāng)描寫開殼層體系時，波函數(shù)ψ一般應(yīng)滿足條件(2.13a)的一些Slater行列式線性組合(2.15)構(gòu)成。從而，計算方案將更加復(fù)雜。然而，對于開殼層體系對應(yīng)于極大多重度的狀態(tài)來說，可以保持波函數(shù)的單行列式表示。非限制Hartree-Fock(UHF)方法[4,9]是描寫這類體系的可能方法之一。UHF方法的根本假定是，α自旋電子所處的分子軌道不同于β電子。從而，與閉殼層體系不同，在UHF法中引入兩組分子軌道：p個α自旋電子置于分子軌道集合ψiα中，而q個β自旋電子置于分子軌道集合ψiβ中。電子體系的波函數(shù)為(2.49)由于α自旋和β自旋電子占據(jù)不同的空間軌道，所以應(yīng)用ψUHF就在每種程度上考慮了不同自旋電子的相關(guān)效應(yīng)。常用的還有另一種處理開殼層體系的限制型Hartree-Fock(RHF)方法。兩種不同方案的圖像表示如圖2.2.所示。把(2.49)代入(2.16)式中，可以求出ψUHF所描寫的體系的能量表達式(2.50)其中交換積分Kijα和Kijβ分別用分子軌道ψiα和ψiβ計算。當(dāng)ψα＝ψβ和p=q時，(2,50)式自動復(fù)原為(2.21)式，而ψUHF也就變?yōu)殚]殼層體系的波函數(shù)(2.14)式。如前所述，對于分子軌道ψiα和ψiβ分別向原子軌道φμ做LCAO展開有(2.51)此時，展開系數(shù)cα不同于cβ。ψiα和ψiβ各自滿足正交歸一化條件，因為它們不同的自旋因子保證了(2.49)式中的ψiα與ψiβ的正交性?？梢苑謩e引入α自旋和β自旋電子的原子軌道密度矩陣元(2.52)顯然，總電子密度矩陣元等于兩者的和，即(2.53)而兩者之差定義了自旋密度矩陣元(2.54)將ψiα和ψiβ的展開式(2.51)代入(2.50)式中積分表達式，可得到用原子基AO表達的能量公式，并考慮到(2.52)式，E可以寫為(2.55)相對于cμiα和cμiβ各自使用變分法極小化能量(2.55)，就導(dǎo)出聯(lián)立的兩套方程組，可計算分軌道ψiα和ψiβ的能量εiα和εiβ及系數(shù)cμiα和cμiβ。(2.56)上式中α自旋和β自旋電子的Hartree-Fock算符矩陣元為(2.57)或者寫成以下矩陣公式(2.58)正如閉殼層情況一樣，方程(2.56)也是系數(shù)cα和cβ的三次聯(lián)立方程組，只能用迭代法求解。UHF方法的缺點是，波函數(shù)ψUHF滿足條件(2.13b)，Ms=(p-q)/2，但不是S2的本征函數(shù)，即不對應(yīng)于任一個總自旋值。在一般情況下，可以把ψUHF表成具有不同自旋多重度波函數(shù)的線性組合(2.59)就數(shù)值來說，展開式(2.59)中系數(shù)Cs+m很快變小，所以別離出第一個混合態(tài)S'=S+1常常足夠了，可以借助于消滅算符［10］(2.60)用原始波函數(shù)ψUHF表達與波函數(shù)AS+1ψUHF有關(guān)的密度矩陣(2.52)、總密度矩陣(2.53)和自旋密度矩陣(2.54)式，它們是UHF方法中極重要的一些物理量，詳細的表達式見引文［10］。開殼層體系的限制性Hartree-Fock方法假定ψiα≡ψiβ，就可以從ψUHF得到滿足(2.13a)的分子電子波函數(shù)，這就導(dǎo)出一個Slater波函數(shù)，它的n1個軌道均為兩個自旋反平行電子占據(jù)〔閉殼層〕，而n2個軌道為自旋一樣電子單占據(jù)〔開殼層〕。(2.61)其圖像如圖2.2所示，波函數(shù)ψUHF是算符S2和Sz的本征函數(shù)，同時描寫極大多重度狀態(tài)S=MS=n2/2。所謂限制性Hartree-Fock法就是用(2.61)型波函數(shù)作運算的。把(2.61)代入(2.16)式，作一些變換后，可得到RHF法的電子能量表達式(2.62)其中k和l表示閉殼層局部的分子軌道，而m和n表示開殼層的分子軌道。量f表示開殼層的占據(jù)程度；a和b是Roothaan常數(shù)[11]，它們?nèi)Q于所研究電子體系的具體特征。例如，對半充滿的開殼層體系，f=1/2，a=1，b=2。公式(2.62)中其余量的含義同2.2節(jié)。在LCAO近似下，此時閉殼層和開殼層分子軌道均可按原子軌道基集合φμ展開成(2.19)式。運用關(guān)系(2.28)-(2.30)就可用AO基改寫(2.62)式。當(dāng)使用矩陣表示時，可以寫成(2.63)其中ν1=2和ν2=2f是閉殼層和開殼層的填充數(shù)；R1和R2分別是相應(yīng)于閉殼層和開殼層局部的原子基的密度矩陣［見(2.33)式］。閉殼層和開殼層局部的電子相互作用矩陣分別是G1和G2，它們是(2.64)矩陣G和G'可以用庫侖矩陣[(2.38)式]和交換作用矩陣[(2.39)]表示為(2.65)(2.66)矩陣(2.65)描述閉殼層電子相互作用，而矩陣(2.66)描述開殼層電子相互作用。用極小化(2.63)式，可以得到確定分子軌道展開成AO基的系數(shù)方程，當(dāng)然這必須在分子軌道相互正交的條件下進展，分別對閉殼層和開殼層軌道的LCAO系數(shù)變分。通常這導(dǎo)出兩個(2.43)型的Hartree-Fock矩陣，對閉殼層和開殼層分別有(2.67)還可以由不同殼層分子軌道正交條件導(dǎo)出兩個(2.44)型的方程組。然而，正如Roothaan指出的［11］，對于所研究的函數(shù)，可以將兩個(2.44)型的方程組統(tǒng)一成一個求本征值和本征向量的方程(2.68)當(dāng)然，這需要按下述公式確定上式中的Hartree-Fock矩陣(2.69)矩陣Ri'=1-Ri(i=1,2)，R3'=1-R1-R2，此處1為單位對角矩陣。可以把(2.69)式寫成與(2.63)式相應(yīng)的形式(2.70)這里B=2aJ2-bK2；R是總密度矩陣R=R1+fR2。庫侖積分和交換積分下角標1和2分別對應(yīng)于閉殼層和開殼層，包含有(2.70)式定義的Hartree-Fock算符方程(2.68)的求解歸結(jié)為逐次迭代法。矩陣R1由n1個閉殼層軌道計算，而R2由n2個開殼層軌道計算。波函數(shù)(2.61)式以及公式(2.68)-(2.70)不僅可用于描寫極大自旋多重度狀態(tài)，而且也可用于描寫開殼層軌道局部填充的一些狀態(tài)。此時，只須改變(2.70)式中a，b，f數(shù)值，也并不使問題復(fù)雜化。Roothaan論文中已給出用RHF方法描述分子狀態(tài)的這類系數(shù)［11］。Pitaer等人推廣了Roothaan的研究，除了給出原子和線型分子能量系數(shù)外，還給出了Td和Oh分子的能量系數(shù)，使該方程更一般化了［12］。SCF計算結(jié)果的分析通過SCF計算得到了收斂波函數(shù)|ψ>，就可以應(yīng)用此波函數(shù)|ψ>分析計算常用的物理量，現(xiàn)將主要的物理量計算介紹如下?？偰芰矿w系的電子總能量E是Hamilton算符的平均值<ψ|H|ψ>，即(2.71)而體系的總能量ET是電子能量加上原子核－原子核排斥能量(2.72)這是最有意義的量，因為預(yù)期穩(wěn)定分子的平衡幾何構(gòu)型在ET的極小點處。Koopmans定律電離勢的理論計算式是提供說明分子光譜的有用工具。與半經(jīng)歷計算不同，從頭算處理了全部電子，因此不必區(qū)分價電子和內(nèi)殼層電子電離勢。一般不去區(qū)分光電子光譜〔PES〕和學(xué)分析電子光譜(ESCA)。兩者在理論上是統(tǒng)一處理的。Koopmans定律提供了電離勢的理論上最簡單近似[13]，它使親態(tài)閉殼層體系的Hartree-Fock軌道能量εh的負值等于第h個電離勢(2.73)所以，按Koopmans定律，假設(shè)理論上有假設(shè)干個電子占據(jù)分子軌道，則光電子光譜應(yīng)當(dāng)存在有相應(yīng)個譜帶，粗略地說，在多數(shù)情況下實驗觀察確實如此[14]。然而，定量地講，方程并非很滿意地列數(shù)了觀測譜。這是由于Koopmans定律有三個內(nèi)在缺陷：①Koopmans定律不管從分子到離子時軌道的自洽重排，亦即忽略了與這種重排有關(guān)的能量ΔR，通常叫做重新組合能或弛豫能?？傊甂oopmans定律認為從分子到離子所有占據(jù)的分子軌道保持不變，這是一種近似。②沒有考慮從分子到離子的相關(guān)能變化ΔC。③假定分子和離子的的相對論能量一樣。首先，相對論能量隨著原子系數(shù)的增加而增大，并且對于內(nèi)殼層來說是很大的數(shù)值。例如對內(nèi)殼層電離勢來說，CH4的相對論修正為0.1eV，對于Ne的修正為0.8eV，而對于Ar對修正則是14eV。由于從頭計算不直接計算相對論修正。其次，當(dāng)忽略相對論效應(yīng)時，對于內(nèi)殼層從頭計算給出電離勢的估計值數(shù)值偏低。從而不考慮相對論效應(yīng)時，電離勢應(yīng)為(2.74)差值IPhKRC-IPhK是所謂的Koopmans缺欠，即應(yīng)當(dāng)考慮電子弛豫和相關(guān)能修正。當(dāng)采用ΔSCF方法時，可以考慮電子重排能ΔR(2.75)其中ESCF是電離前分子的SCF總能量，而EhSCF是空穴態(tài)離子的SCF能量。這說明，電離勢是絕熱SCF計算的分子和離子的總能量之差，有時稱之為豎直電離能。單電子算符的平均值諸如分子的偶極矩、四極矩、核位置的場梯度、反磁化率等可以采用單粒子算符的和描寫，一般形式為(2.76)其中g(shù)(k)是任意一種只依賴于單電子坐標的算符。單粒子算符的平均值是(2.77)從而，除密度矩陣外，只須計算單電子積分<ν|g|μ>。下面以偶極矩為例加以說明。偶極矩的經(jīng)典定義為(2.78)向量ri除的電荷為qi。偶極矩相應(yīng)的量子力學(xué)定義為(2.79)其中第一項為哪一項電子〔電荷為-1〕對偶極矩的量子力學(xué)奉獻，第二項是原子核的核電荷ZA的經(jīng)典奉獻項。電子偶極矩算符是單粒子算符。所以應(yīng)當(dāng)有(2.80)這是一個向量表達式，其分量公式可以寫為(2.81)對于μy和μz也有類似的公式。為了計算偶極矩，除了密度矩陣外Pμν外只須計算偶極矩積分(2.82)至于與電荷密度有關(guān)的問題下面加以說明。集居數(shù)分析在分子軌道計算中，人們通常把分子中電荷密度定義為(2.83)它表示空間任意點發(fā)現(xiàn)電子的幾率，常常畫出通過分子中各種截面的電荷密度輪廓圖[16]。但是，由于分子給定原子或原子核中的電子數(shù)不是唯一確定的，有時可以給出不同的集居數(shù)分析。因為(2.84)每個MO占有兩個電子。把ψi(r)的展開式代入上式就有(2.85)其中(PS)μμ為聯(lián)系于ψμ軌道的電子數(shù)，這稱為Mulliken集居數(shù)分析[17]。對于給定原子A在所有軌道上的電子數(shù)為，從而A原子的凈電荷數(shù)為(2.86)其中ZA是A原子核的核電荷數(shù)，求和取遍所有屬于A原子的軌道ψμ。但是，上述定義(2.85)不是唯一確定的，因為TrAB=TrBA，就有(2.87)對于任意a成立。當(dāng)a=1/2時，得到(2.88)其中P'μμ是對稱正交化基組表象中的密度矩陣(2.89)(2.90)此時，P'的對角元構(gòu)成Lowdin集居數(shù)分析[18](2.91)上述兩種集居數(shù)分析都是計算和文獻中常用的。參考文獻唐敖慶等.量子化學(xué),科學(xué)出版社,1981.徐光憲等.量子化學(xué)—根本原理和從頭計算(中冊),1985.SzaboA,UstlundNS.ModernQuantumChemistry：AnIntruodutiontoAdvancedElectronicStructureTheory,MacmillanPublishingCo,NewYork,1981.PopleA,BeveridgeDL.ApproximateMoleculerOrbitalTheory,McGraw-HillBookCompany，NewYork,1970.王志中,李向東.半經(jīng)歷分子軌道理論與實踐,科學(xué)出版社,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.沈爾中,王志中,江元生.分子科學(xué)與化學(xué)研究，1984,5,427.MullikenRS.J.Chem.Phys.1955,23,1833.LowdinPO.J.Chem.Phys.1953,21,374.2.2從頭算方法量子化學(xué)從頭計算〔Abinitio〕方法60年代在國際上開場流行起來，隨著計算機技術(shù)的開展，到了70年代得到了廣泛的研究和應(yīng)用。所謂“從頭計算法〞是立足于以下三條根本近似：①采用非相對論量子理論，即從Schrodinger方程出發(fā)；②使用Born_Oppenheimer近似；③軌道近似，即單粒子函數(shù)近似。除了這些近似以外，不再借助于任何經(jīng)歷參數(shù)，所以稱之為從頭計算。從頭計算實質(zhì)上是以前面所述分子軌道理論為根基的，從HFR或UHF方程出發(fā)，適當(dāng)?shù)剡x取原子軌道的線性組合法去模擬分子軌道計算時，一旦選定了基函數(shù)后，對于HFR或UHF方程所涉及的全部單電子積分和雙電子積分均嚴格地計算。一般來說，選取從頭計算基函數(shù)有兩種方法：Slater型軌道〔STO〕，用它描述電子云分布優(yōu)于其它函數(shù)，但計算積分時涉及一些無窮級數(shù)，十分復(fù)雜，所以只是早期對一些小分子實現(xiàn)了這類計算；選用Gauss型函數(shù)〔GTF〕作為基函數(shù)，使多中心積分容易計算，且有簡明的解析表達式，便于在計算機上實現(xiàn)，所以得到了廣泛的應(yīng)用。使用Gauss函數(shù)作為基函數(shù)求解Schrodinger方程，最初是Boys在1950年提出的，到了70年代初才得以迅速開展。STO與分子積分STO定義量子化學(xué)中廣泛應(yīng)用的一類基函數(shù)是Slater軌道〔STO〕：(2.92)其中Slm(θφ)是歸一化實球諧函數(shù)(2.93)式中Pl|m|是連屬Legend多項式，當(dāng)m=0時δm=0，而m≠0時δm=1，ζ稱為軌道指數(shù)。(2.92)式存在多種確定軌道指數(shù)ζ的方法。Clementi-Remond曾按極小化孤立原子能量定出了原子的ζ值［19］。Burns以最正確重復(fù)Hartree-Fock原子軌道方式定出了另一套指數(shù)[20]Gauterman-Zener用Hartree-Fock軌道給出的金屬-配位體重疊積分為標準定出一組ζ值[21]。但是最常用的還是Slater在1930年建議的計算ζ值的經(jīng)歷方法[22]，該方法是基于很好再現(xiàn)原子的前幾個電離勢。分子積分[23-35]使用STO為基函數(shù)，除了函數(shù)的表達式簡單外，它們還有疊加性。此時所需要計算的分子積分有以下四類重疊積分(2.94)動能積分(2.95)核吸引能積分(2.96)電子相互作用能積分(2.97)上述式子中，A、B等為原子核的位置向量，而r1和r2是電子的位置向量。(2.94)-(2.97)是所謂的多中心積分，其中最難計算的是(2.96)和(2.97)式中的三中心和四中心積分。實際上，到目前為止，尚無有效地解析式，允許人們用Slater軌道快速計算多原子分子問題所需要的多中心積分。大局部的單中心和雙中心積分是易于計算的。雙中心積分JAACC是庫侖積分，JABAB是交換積分，而JAAAB稱為雜化積分。詳細的積分公式和推導(dǎo)可參考文獻[23-37]。2.2.2GTF及Gauss函數(shù)積分由于STO分子積分的復(fù)雜性，限制了直接使用其進展從頭計算的應(yīng)用。1950年Boys建議用Gauss函數(shù)為基函數(shù)計算分子積分[38]。GTF的定義及其性質(zhì)徑向Gauss函數(shù)定義為(2.98)顯然，它與STO的主要差異在于e-αr2和e-αr，兩者角度函數(shù)一樣。另一種常用的笛卡爾GTF表達式為(2.99)其中歸一化常數(shù)N為(2.100)GTF的重要特點是，兩個中心的GTF相乘可以表成另一個中心GTF函數(shù)的線性組合，使得多中心積分易于計算，這便是使用GTF的優(yōu)越性。(2.101)(2.101)式詳細的證明和推導(dǎo)可參考文獻[2]。GTF的積分公式應(yīng)用(2.98)式表達研究分子的徑向GTF積分，最簡便的是所謂無節(jié)點Gauss函數(shù)〔即n=1,l=m=0,1s函數(shù)〕。這是Boys[38]和Shavitt[39]所完成的。實際計算中使用笛卡爾GTF積分公式，它是由Taketa-Huzinaga-Ohata在1964年推導(dǎo)出來的[40]。唐敖慶等著《量子化學(xué)》[1]、王志中編著的《現(xiàn)代量子化學(xué)計算方法》[41]等書中也有較詳細的推導(dǎo)。重疊積分：(2.102)(2.103)動能積分：(2.104)電子核引積分：(2.105)(2.106)雙電子相互作用積分：(2.107)(2.108)其中已定義(2.109)求和指標限制為(2.110)從所給出的公式可以看出，用笛卡爾Gauss函數(shù)計算分子積分是方便的，同時易于計算機程序化。2.2.3解HFR方程的自洽場〔SCF〕過程HFR方程的求解，實際上要求解Fock算符的本征矢{cμi}—分子軌道組合系數(shù)，以及對應(yīng)的本征值εi—分子軌道能量。但Fock算符本身又是分子軌道組合系數(shù)的二次函數(shù)，故只能用迭代方法求解，其SCF過程如下：(1)給出分子的各個核坐標，以及核電荷、電子總數(shù)，并選定基組{φμ}。(2)計算分子積分Sμν，hμν以及<μν|λσ>。(3)求出X。(4)初始猜想：即給出一組分子軌道，因基組是固定的，實際是給出展開系數(shù)，從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是否與上次一樣〔在一定限度內(nèi)〕，假設(shè)不收斂，則回到第(5)步繼續(xù)迭代。(12)如果收斂，用所得的結(jié)果c，P，F(xiàn)等，計算一些有用的量.。能量梯度和平衡幾何構(gòu)型能量梯度的解析表達式勢能對核坐標的一階導(dǎo)數(shù)稱之為勢能在該方向的梯度，勢能面上的駐點〔極大點、極小點和鞍點〕都滿足所有核在各個方向上所受的力為0，即這些點的勢能梯度向量為0向量。得到勢能的梯度值后，輔以一定的數(shù)學(xué)手段，就能夠優(yōu)化出勢能面上駐點的幾何構(gòu)型，并求得分子波函數(shù)和能量。求解能量梯度的解析方法是優(yōu)化平衡幾何構(gòu)型和過渡態(tài)構(gòu)型的根基。下面給出能量梯度的解析表達式，為與分子軌道的標號i區(qū)別，用表示某原子核在某坐標的梯度〔這里的體系是閉殼層的N電子體系〕(2.111)其中，核排斥項用梯度方法優(yōu)化平衡幾何構(gòu)型分子總能量的絕對極小值所對應(yīng)的幾何構(gòu)型稱為分子的平衡幾何構(gòu)型。有如下的三種方法(1)可變尺度法〔VariableMetricMethod〕(2.112)X為坐標向量，可以是直角坐標或內(nèi)坐標；an為尺度因子，A為嚴格推導(dǎo)的能量二階導(dǎo)數(shù)矩陣的矩陣。調(diào)整尺度因子使g到達所要求的收斂限，以求得平衡幾何構(gòu)型。(2)極端方案a.最陡下降法：An=1，為單位陣。b.Newdon-Raphson法：?。?，An為嚴格的力常數(shù)矩陣。(3)Muttagh-Sargent法：介于a和b的折衷方法，An用某種經(jīng)歷公式得到。參考文獻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.王志中.現(xiàn)代量子化學(xué)計算方法,吉林大學(xué)出版社，1998，43.2.3半經(jīng)歷SCF-MO方法以Roothaan方程及其求解為核心的SCF分子軌道方法，是從純計算的觀點〔從頭算〕討論了分子軌道理論，近似的分子軌道方法是以該理論為骨架，從兩個方面簡化計算：①忽略一些難于計算的積分，如大量的雙電子積分；②對剩下的積分，根據(jù)電子構(gòu)造的實驗資料估計一些積分值。尤其是第二方面的努力，不能僅僅看作是為了簡化計算，而是不再試圖根據(jù)量子力學(xué)原理直接推導(dǎo)分子的性質(zhì)，主要致力于解釋實驗數(shù)據(jù)〔現(xiàn)象〕內(nèi)部的相互關(guān)系，它在許多方面取得了驚人的效果。這種基于實驗結(jié)果的相互關(guān)系與嚴格量子理論的結(jié)合是半經(jīng)歷的。2.3.1零微分重疊近似〔ZDO〕[42,43]LCAO自洽場分子軌道計算的最大困難是對雙電子積分的處理，這包括兩個方面的內(nèi)容：一是積分本身計算困難，特別是多中心的雙電子積分在數(shù)學(xué)上的處理較困難；二是雙電子積分數(shù)目隨基函數(shù)數(shù)目的增加呈指數(shù)形式增長。而大量的實踐結(jié)果說明，許多這類電子排斥的數(shù)值接近于0，尤其是那些包含重疊分布，且的電子排斥積分。因此，在開展近似自洽場分子軌道方案中，一個有用的近似是系統(tǒng)地忽略數(shù)值一致的小的電子排斥積分。這一點通過零微分近似實現(xiàn)，零微分近似寫為：〔當(dāng)時〕。于是，假定包含重疊積分的電子排斥積分數(shù)值小到可以忽略，在零微分近似下其次，在分子軌道歸一化中忽略相應(yīng)的重疊積分，也就是但不忽略包含重疊積分的實積分。所討論的各種等級的近似自洽場方法的最主要區(qū)別是在處理電子排斥積分中應(yīng)用零微分重疊近似的程度。假設(shè)將零微分近似應(yīng)用于全部原子軌道，則閉殼層的Roothaan方程簡化為(2.113)這實際上就是簡化的矩陣形式Fc=cE。其中Fock矩陣元為其對角元()非對角元()這些近似大大簡化了計算，主要是忽略了許多困難的雙電子積分，尤其是全部三中心和四中心積分都變?yōu)?。2.3.2全略微分重疊〔CNDO〕[5，44-45]CNDO明確地只處理價電子，內(nèi)殼層被看成剛性實，于是原子軌道基集合是一個價基集合。其根本近似是對所有不同原子軌道積運用零微分重疊近似，即使用上面推出的簡化方程Fc=cE和Fock矩陣元公式；此外對剩下的局部積分用實驗參數(shù)代替。概括起來，CNDO有以下五步近似：(1)將Roothaan方程中重疊矩陣用單位陣代替，即在分子軌道歸一化中忽略重疊積分；(2)忽略所有雙電子積分中的微分重疊，即；(3)將剩下的庫侖型積分均簡化為每一對原子一個數(shù)值〔參數(shù)〕：，屬于A原子，屬于B原子；(4)在其它包含原子實的積分中忽略單原子微分重疊：，這里是B原子核及內(nèi)層電子〔剛性實〕產(chǎn)生的勢能，、屬于A原子，也用參數(shù)化處理；(5)取雙原子非對角實矩陣元正比于相應(yīng)的重疊積分：，屬于A原子，屬于B原子，也是參數(shù)化處理。使用所有這些近似，使計算更簡單，且參數(shù)化有不同的方案，這里不再詳細介紹。2.3.3間略微分重疊〔INDO〕[46,47]CNDO對所有不同的原子軌道運用零微分重疊近似，固然對雙電子排斥積分的處理是簡單的，但它對兩個具有平行自旋的電子間實際存在的不同相互作用未給予適當(dāng)考慮，尤其是兩個電子同屬一個原子的情況。波函數(shù)的反對稱性要求自旋平行的電子不占有空間的同一小區(qū)域，因此同一原子上不同原子軌道中的兩個電子如果自旋平行，將具有較小的排斥能，這種差異在數(shù)學(xué)上表現(xiàn)于雙電子交換積分：而在CNDO理論中，這種積分是被忽略的。這樣，原子的電子能量與電子自旋取向無關(guān)，自旋多重度無法區(qū)分，導(dǎo)致CNDO不能描述對應(yīng)于同一組態(tài)的電子能量分裂。為抑制上述缺點，最簡單的方法就是對交換積分作某些考慮，即只保存單中心積分中的單原子微分重疊。所謂單中心積分，是積分只涉及同一個原子；而單原子微分重疊，指。例如，雙電子積分中的<SAPA|SAPA>要予以保存，而<SAPA|SAPA>雖屬單原子微分重疊，但非單中心積分而予以忽略。這種微分近似方法稱為間略微分重疊〔INDO〕。2.3.4忽略雙原子微分重疊(NDDO)[48]下一等級近似的自洽場理論，是只忽略屬于不同原子的原子軌道的微分重疊，即在化簡HFR方程時，對于和屬于不同原子的原子軌道才使用零微分重疊近似。這一等級近似的主要特點是保存了偶極－偶極相互作用，因為此時保存了型積分。NDDO法雖然也忽略了全部三中心和四中心積分，但是所保存的單中心和雙中心積分數(shù)目仍很多，過渡到半經(jīng)歷水平〔參量化〕比擬困難。2.3.5NDDO的改良[49-59]歷史上，NDDO法對雙中心積分采用Slater軌道直接計算，程序很繁瑣，只能計算一些不大的分子。Dewar等人于1977年建議了改良的NDDO法，即MNDO法[54]，主要特點是用半經(jīng)歷模型處理雙中心電子排斥積分。1985年，Dewar等人又對MNDO進一步修正，主要針對MNDO高估了原子間的排斥作用，對MNDO中的核實排斥函數(shù)進展了修正，加進了一些高斯函數(shù)項，同時進展重新參數(shù)化，即AM1方法(AustinModel1)[57]。這種方法目前仍被廣泛應(yīng)用。1989年，Stewart在AM1的理論框架下，發(fā)表了新的參數(shù)化方法并進展了再參數(shù)化工作，稱為PM3(ParamiterMethod3)[59]。到目前為止，MNDO、AM1、PM3是使用最為廣泛的半經(jīng)歷方法。MNDO和AM1方法改良的重點是修改原子實間的排斥項。從靜電的角度看，原子實-原子實排斥積分為〔eV〕(2.114)然而，在范德華半徑以外，電子-電子和電子-原子核實的作用并不按照1/RAB的規(guī)律衰減，如果使用(2.114)的形式，將造成兩個中性原子或分子之間的凈排斥。為了矯正這一點，減小過高的排斥能，作如下近似MNDO對核實-核實相互作用進展了修正，采用如下形式(2.115)但是，對于O-H和N-H間的相互作用采用下式單獨處理(2.116)AM1在此根基上作了進一步的修改(2.117)其中比MNDO多出的球面Gaussian函數(shù)項的作用是，在成鍵的距離之外進一步減小多余的核實-核實排斥作用，它可以理解為是vanderWaals吸引項。(2.117)式中的a,b,c是可調(diào)參數(shù)。參考文獻PariserR,ParrRG.J.Chem.Phys.1953,21,767.PopleJA.Trans.Farady.Soc.1953,49,1375.波普爾J.A.等著,江元生譯,近似分子軌道理論方法,科學(xué)出版社,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〕在處理分子和原子的電子構(gòu)造問題變得非常熱門。DFT的主要思想是考慮總電子能量和總電子密度之間的關(guān)系。這一思想的近似模型〔TheThomas-Fermimodel〕早在上個世紀20年代就已提出。但是真正的突破性工作來自Hohengerg和Kohn于1964發(fā)表的論文，以及Slater和Sham等人的改良，他們指出體系基態(tài)的能量和其他性質(zhì)是由電子密度唯一確定的。DFT方法的主要觀點是通過電子密度函數(shù)來描述相互作用的費米子體系，而不是多體波函數(shù)。在HF理論中，多電子體系的波函數(shù)是由Slater行列式來表達的，它由N個單電子波函數(shù)〔N是分子中電子的數(shù)目〕構(gòu)成，通常含有N個電子的體系包含3N個變量。DFT也考慮單電子函數(shù)。然而，HF理論是全部計算N個電子的波函數(shù)，DFT卻只計算總的電子能量和所有電子密度的布居。無論體系多大，對于電子密度來講只有三個變量x、y、z，這樣可以大大簡化計算。DFT現(xiàn)在已經(jīng)廣泛地應(yīng)用于化學(xué)熱力學(xué)、分子光譜等領(lǐng)域，可以用來描述金屬、半導(dǎo)體、絕緣體的基態(tài)性質(zhì)，和分子軌道理論一起，成為人們探討分子體系的有力理論工具。它的成功之處在于可以成功的解釋標準的大體積宏觀物質(zhì)，而且它也能解釋諸如蛋白質(zhì)、碳納米管等復(fù)雜的物質(zhì)。DFT理論的優(yōu)點還在于它提供了許多化學(xué)上重要的概念，例如電負性、硬度、Fukui函數(shù)、相應(yīng)函數(shù)等，這些概念可以很方便的用來解釋化學(xué)性質(zhì)和分子的變化，它在預(yù)言原子、分子和固體基態(tài)的能量、電子自旋密度、鍵長和鍵角等方面已經(jīng)取得了初步的成功。2.4.1Thomas-Fermi模型1927年的Thomas-Fermi模型是最早的嘗試用電子密度函數(shù)來代替體系的波函數(shù)描述和確定體系的性質(zhì)。它的根本設(shè)想是把空間看成許多邊長為1的小立方體，每個立方體都是一個三維勢箱，則有如下公式從上述公式中可以看出，以坐標系原點為球心、R為半徑的球面上的()點為標志的量子態(tài)〔此時，都是整數(shù)〕具有一樣的能量，稱為簡并態(tài)。電子按能量沿著球面由里向外、由低到高依次填充，但只有都為正的球體中的點才代表有電子占據(jù)的能級，由此可知能量小于〔或等于〕ε的量子態(tài)數(shù)目為能量在ε和ε＋?ε間的能級數(shù)目為其中，是能級ε處的態(tài)密度。在0K時，能量在費米能級〔離子的最高占據(jù)能級〕以下的粒子組成的體系總能量為(2.118)假定每個()量子態(tài)有兩個粒子占據(jù)〔自旋方向相反〕，則離子總數(shù)為即可得到因為體系V＝0，故ΔE是體系的動能，令〔電子密度〕，對空間所有小立方體求和〔積分〕，使用原子單位，則得到，(2.119)這就是動能的電子密度的泛函表達式。對于原子，只考慮電子間的庫侖作用，有，滿足用拉格朗日乘子法，求和，令為拉格朗日乘子，由此得Euler-Lagrange方程其中(2.120)上式與條件聯(lián)立，可求解得到和。實際計算說明，分子的能量總大于其組成原子能量之和，即得不到穩(wěn)定的分子體系。因此，雖然化學(xué)家用電子密度分布概念解釋了很多現(xiàn)象，但實際結(jié)果似乎說明電子密度變化的理論缺乏物理學(xué)依據(jù)。2.4.2Hohenberg-Kohn定理現(xiàn)代密度泛函理論開場于Hohenberg-Kohn〔1964年〕定理，它指出電子體系的基態(tài)電荷密度分布唯一的決定體系總能量，可以用公式〔是體系基態(tài)的密度〕，包括第一定理和第二定理。第一定理：體系的基態(tài)的電子密度與體系所處勢場有一一對應(yīng)的關(guān)系，從而完全確定體系的所有性質(zhì)。而如果兩個體系具有一樣的基態(tài)電荷密度分布，那么它們的勢場只相差一個常數(shù)確定。該定理初步論證了用電子密度函數(shù)代替波函數(shù)來全面描述一個體系是一個可行的方法?？梢杂梅醋C法進展證明，假設(shè)兩個體系分別處于勢場和，并且它們的差值不是常數(shù)，它們的哈密頓算符可以如下表示上式中是電子體系的哈密頓算符。設(shè)和的非簡并態(tài)波函數(shù)分別為Ψ1和Ψ2，相應(yīng)的電子密度為和，則且小于依據(jù)變分原理，有即同理可得如果，則兩個不等式相加就有，這是不可能的，因此。也就是當(dāng)常數(shù)時，，這就證明了電子密度和體系所處勢場之間有一一對應(yīng)關(guān)系。當(dāng)體系電荷密度分布確定時，勢場就確定了，電子數(shù)目也確定了，所以體系的哈密頓也確定了，進而可知體系的波函數(shù)和所有其它性質(zhì)都被確定了。也就是說，體系總能量、動能、電子間相互作用能都是的泛函，記作、、，并有其中，與外勢場沒有明顯的關(guān)系。第二定理：對于任意的函數(shù)，如果滿足，，則，為體系基態(tài)能量。這一定理說明了體系基態(tài)總能量是總能量密度泛函的極小值，即，可以簡單地表示為〔是體系基態(tài)的密度，是體系處在非基態(tài)時的電子密度〕。該定理證明如下：對應(yīng)一個、和Ψ′，設(shè)Ψ′為試用變分函數(shù)，則。按拉格朗日乘子法變分，μ為拉格朗日乘子，得到Euler-Lagrange方程如果知道的具體表達式，結(jié)合的限制條件，可由上式求得及。在這里我們采用Levy提出的一種求的搜索方法即在所有給出一樣的的中，使取極小值者給出的值。于是用變分法求基態(tài)能量的過程可表示為(2.121)2.4.3Kohn－Sham方程Kohn－Sham方程主要是通過把動能泛函的一局部別離出來而得到的。假設(shè)不考慮自旋，而且粒子是相互獨立的，那么動能可表示為表示為其中，是庫侖作用能，是交換相關(guān)能，包含T和TS的差，是一個量值很小的泛函?？杀硎緸閷ψ兎智竽芰繕O值得到Kohn－Sham方程其中，，假設(shè)，求解Kohn－Sham方程，可得和從Kohn－Sham方程的形式來看，它與Hatree-Fock方程非常類似，而且物理意義上也類似，所以雖然不具有分子軌道能級的含義，通常還是把看作近似的分子軌道。2.4.4Kohn－Sham方程的求解方法DFT分子計算的過程實際上就是用變分法求解分子的Kohn-Sham方程。由Kohn-Sham方程和Hatree-Fock-Roothaan方程類似的特點可知，方程的求解也要用自洽迭代方法。但DFT的計算量較小，要用數(shù)值積分方法計算，不能采用解析法計算。盡管現(xiàn)在比擬流行的數(shù)值積分方法已經(jīng)有求積方法[77]、分區(qū)求積方法[78]和重迭區(qū)間空間分割法[79]三種方法，但是由于矩陣元的數(shù)值計算占用密度泛函方法計算的大局部資源，因此尋找高效高精度的數(shù)值積分方法仍然是一個很受重視的課題。2.4.5交換能和相關(guān)能校正均勻電子云的交換能可表示為，式中第一項為交換能，第二項為相關(guān)能，參數(shù)，γ代表自旋，表示在上自旋和下自旋密度分別為和的電子云中單個電子所具有的相關(guān)能。的準確解析形式至今還不知道，Stoll等人認為LDA的主要誤差之一是均勻電子云模型過高估計了自旋一樣的電子之間的相關(guān)能。真實體系的電荷分布是不均勻的，定域密度近似忽略了電荷不均勻性的影響，因此計算的精度也就受到了限制。本文的計算大多采用的BeckeIII的交換能校正公式：(2.122)其中，，，。一般來說，LDA過高估計相關(guān)能，與交換能的處理類似，人們引入密度梯度校正來減小LDA相關(guān)能的誤差。各種近似公式的計算結(jié)果差異不是很大，本文采用的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勢能面、平衡幾何構(gòu)型與內(nèi)稟反響坐標研究一個化學(xué)反響的反響途徑，需研究這個反響體系的勢能面，由于對整個反響勢能面的研究不易做到，因而一般只研究勢能面上駐點的性質(zhì)。駐點是指反響途徑上的能量的極小點和鞍點，它們是根據(jù)其能量的一階導(dǎo)數(shù)和二階導(dǎo)數(shù)定義的。對于有幾個變量的函數(shù)，對應(yīng)于各個變量的一階導(dǎo)數(shù)形成一個梯度的矢量，二階導(dǎo)數(shù)形成Hessian矩陣。在經(jīng)典力學(xué)中，一個粒子勢能的一階導(dǎo)數(shù)對應(yīng)于作用在粒子上的力的負值，而二階導(dǎo)數(shù)則對應(yīng)于力常數(shù)。因而，梯度的負值對應(yīng)于分子中原子或核上所受的力，二階導(dǎo)數(shù)矩陣〔Hessian矩陣〕也就稱為力常數(shù)矩陣。對于只含一個變量的函數(shù)，一階導(dǎo)數(shù)為零，對應(yīng)于一個極小點或極大點。在多維空間，駐點對應(yīng)于所有的一階導(dǎo)數(shù)必須是零，梯度矢量才為零長度。同樣，分子中原子上所有的力必須為零，這樣的點才稱為駐點。駐點的性質(zhì)可由二階導(dǎo)數(shù)決定，在多維空間中，駐點的性質(zhì)是由Hessian矩陣決定的。假設(shè)Hessian矩陣所有的本征值為正，駐點為極小，對應(yīng)一個平衡幾何構(gòu)型；假設(shè)Hessian矩陣有且只有一個負的本征值，即在一個方向為極大，其他所有方向為極小，這樣的駐點對應(yīng)于反響的過渡態(tài)。勢能面上的駐點是通過優(yōu)化方法，計算Hessian矩陣并檢驗它的本征值而找到的，其中較好的是能量梯度優(yōu)化法。本文采用Berny梯度優(yōu)化法來尋找反響的所有駐點。梯度優(yōu)化方法一般是用坐標Xk、計算的能量Ek、計算的梯度Gk和近似Hessian矩陣Hk來定義第k步的勢能面的二項式(2.123)Hessian矩陣的初始猜想為H0，隨著優(yōu)化過程而加以改良以接近真