密度泛函理論的基本計(jì)算方法研究進(jìn)展

唐楠楠，張 姝，李強(qiáng)林

(1.成都紡織高等?？茖W(xué)校，四川成都611731;2.四川師范大學(xué)化學(xué)與材料科學(xué)學(xué)院，四川成都610068)

0 引言

量子力學(xué)是現(xiàn)代物理學(xué)的基石。Heisenberg在1925年首先提出了矩陣力學(xué)，即是量子力學(xué)的矩陣表達(dá)形式。隨后，Schr?dinger在1926年建立了波動(dòng)力學(xué)，即量子力學(xué)的波動(dòng)方程形式。二者雖然在形式上并不相同，但其實(shí)分屬量子力學(xué)的兩種表達(dá)形式，Schr?dinger通過(guò)嚴(yán)格的數(shù)學(xué)推理證明了二者的等價(jià)性。相比矩陣形式，Schr?dinger的波動(dòng)方程形式更易于被化學(xué)工作者所接收。1927年，Heitler和London將量子化學(xué)波動(dòng)形式運(yùn)用到氫分子系統(tǒng)的研究工作中，成功得到了精確描述氫分子系統(tǒng)的波函數(shù)解析解，并對(duì)氫分子系統(tǒng)的化學(xué)鍵本質(zhì)做出了量子力學(xué)解釋，這標(biāo)志著量子化學(xué)的誕生。

量子化學(xué)的發(fā)展給物理學(xué)和化學(xué)提供了嶄新的理論，對(duì)于解釋微觀世界的物質(zhì)規(guī)律和反應(yīng)本質(zhì)起到了決定性的作用。但是人們很快發(fā)現(xiàn)，對(duì)于氫分子以外的大量體系，都很難精確求得解析解。這主要是由于量子化學(xué)雖然在理論上嚴(yán)格的描述了微觀世界的運(yùn)動(dòng)規(guī)律，但是在實(shí)際求解多粒子體系時(shí)會(huì)遇到許多困難，有些求解在當(dāng)時(shí)看來(lái)甚至是不可能完成的工作。

基于此，科研人員開始研究如何簡(jiǎn)化計(jì)算，通過(guò)對(duì)模型和方程的近似處理來(lái)得到可靠精度范圍內(nèi)的數(shù)值解。主要的研究思路是如何將一個(gè)多粒子問(wèn)題轉(zhuǎn)化為多電子問(wèn)題，進(jìn)而再將多電子問(wèn)題轉(zhuǎn)化為單電子問(wèn)題，如果能夠完成這個(gè)步驟，會(huì)使多電子Schr?dinger方程轉(zhuǎn)化為單電子方程，這樣就使求解過(guò)程被大大的簡(jiǎn)化。目前被研究人員廣泛使用的是 Hartree-Fock 近似處理方法［1，2］和密度泛函理論(DFT)方法［3，4］。

這些理論和方法構(gòu)成了第一性原理計(jì)算(first-principles calculation)［5］。隨著計(jì)算機(jī)科學(xué)的迅猛發(fā)展，尤其使得基于DFT的第一性原理得到了廣泛重視，成為物理學(xué)、化學(xué)、生物學(xué)和材料學(xué)等學(xué)科的重要研究方法，被廣泛運(yùn)用于各類原子、分子、團(tuán)簇、晶體和物質(zhì)表面等的電子結(jié)構(gòu)計(jì)算和反應(yīng)過(guò)程分析中［6-8］。密度泛函理論計(jì)算是完全以量子化學(xué)從頭算(ab initio)作為理論依據(jù)的計(jì)算方法，為了區(qū)別于其它量子化學(xué)從頭算方法，有時(shí)也把基于DFT的計(jì)算稱為狹義的第一性原理計(jì)算。本文闡述了DFT理論的發(fā)展歷程和理論核心，并對(duì)其基本計(jì)算原理和發(fā)展趨勢(shì)做出重點(diǎn)分析，希望能對(duì)DFT理論的進(jìn)一步完善提供一定的理論參考。

1 DFT計(jì)算方法發(fā)展歷程和現(xiàn)狀

近20年來(lái)，DFT計(jì)算方法被廣泛地應(yīng)用于物理學(xué)、化學(xué)、材料科學(xué)和生命科學(xué)等領(lǐng)域。DFT方法因其計(jì)算量相對(duì)較小、計(jì)算結(jié)果合理可靠而被人們廣泛接受。DFT方法也是目前計(jì)算科學(xué)當(dāng)中的主流方法之一，它在量子力學(xué)原理和實(shí)際應(yīng)用之間架起一座橋梁，具有不可替代的意義。

密度泛函理論的實(shí)質(zhì)是將波函數(shù)轉(zhuǎn)換為電子密度函數(shù)作為研究的基本單元，從而簡(jiǎn)化運(yùn)算過(guò)程，并且使多電子問(wèn)題得以轉(zhuǎn)化為單電子問(wèn)題進(jìn)行求解。假設(shè)電子數(shù)目為N，則波函數(shù)中總共的變量數(shù)為3N，求解極其復(fù)雜，而采用密度泛函理論，可以將變量數(shù)縮減到3個(gè)，即只包含了三個(gè)空間變量，這樣既可以極大的簡(jiǎn)化計(jì)算過(guò)程，又可以確保計(jì)算精度。

1.1 Thomas-Fermi模型

Thomas和Fermi在1927年建立了基于均勻電子氣模型的Thomas-Fermi模型［9］。均勻電子氣模型的出發(fā)點(diǎn)是假設(shè)電子之間無(wú)相互作用且無(wú)外力干擾。在此基礎(chǔ)上得到電子運(yùn)動(dòng)的Schr?dinger方程為:

解得:

引入0K下的電子排布規(guī)律，可經(jīng)過(guò)推導(dǎo)得到電子密度和單電子總能的表達(dá)式，分別為:

則體系的動(dòng)能密度為:

至此，再引入對(duì)電子間的庫(kù)倫勢(shì)和外場(chǎng)的描述，可以推到出僅由電子密度函數(shù)決定的電子體系的總能表達(dá)式。依據(jù)變分學(xué)概念，我們可以認(rèn)為此時(shí)總能是電子密度函數(shù)的泛函。Thomas-Fermi模型第一次引入密度泛函的概念，成為后來(lái)DFT的雛形。這個(gè)模型的成立是基于理想狀態(tài)下的均勻電子氣假設(shè)，雖然大大的簡(jiǎn)化了計(jì)算形式和過(guò)程，但并沒(méi)有考慮電子間的交換作用，沒(méi)有對(duì)動(dòng)能項(xiàng)進(jìn)行精確描述，所以這個(gè)模型非常粗糙，以此算得的結(jié)果在很多體系中并不合理。但這個(gè)模型給后來(lái)的研究者提供了很好的研究思路，經(jīng)過(guò)多年的努力，研究人員在20世紀(jì)60年代基本完善了密度泛函理論的內(nèi)容，并最終建立起嚴(yán)格意義的密度泛函計(jì)算理論。對(duì)DFT而言，Hohenberg-Kohn定理和Kohn-Sham方程的提出對(duì)本理論的形成和完善起到了關(guān)鍵作用，被譽(yù)為DFT的兩大基石，現(xiàn)分別作如下介紹。

1.2 Hohenberg-Kohn 定理

Hohenberg-Kohn定理是由 Hohenberg和Kohn在系統(tǒng)的研究Thomas-Fermi模型之后提出的，主要針對(duì)了非均勻電子氣模型，其定理的核心思想是:體系中的所有物理量都可以通過(guò)只包含電子密度的變量來(lái)唯一決定，而實(shí)現(xiàn)方法是通過(guò)變分原理來(lái)求得到體系基態(tài)。Hohenberg-Kohn定理由兩條已被嚴(yán)格證明的子定理構(gòu)成［10］。

(1)處在外勢(shì)(除電子相互作用以外的勢(shì))的忽略自旋的電子體系，其外勢(shì)可通過(guò)電子密度唯一決定。

(2)對(duì)于給定的外勢(shì)，系統(tǒng)基態(tài)能量即能量泛函的最小值。

根據(jù)已上定理，假設(shè)一個(gè)體系受到外勢(shì)的影響，則此時(shí)體系的哈密度量由電子動(dòng)能，外場(chǎng)中電子的勢(shì)能和庫(kù)倫作用組成。這三部分中只有外勢(shì)是根據(jù)不同體系而變化的，當(dāng)對(duì)外勢(shì)的描述確定之后，系統(tǒng)的哈密頓量也就隨之確定，而哈密頓量通過(guò)Schr?dinger方程確定整個(gè)體系的波函數(shù)，至此，我們可以看出，電子密度的確能夠唯一的描述總波函數(shù)以及體系基態(tài)。此時(shí)可以對(duì)體系的能量泛函做出如下的描述:

式(1-6)中的四項(xiàng)分別代表外場(chǎng)中的電子勢(shì)能、動(dòng)能項(xiàng)、電子間的庫(kù)倫作用和交換關(guān)聯(lián)能?？梢钥闯?，Hohenberg-Kohn定理嚴(yán)格證明了電子密度函數(shù)的意義，建立了以電子密度函數(shù)為變量的求算形式，并指出了以變分原理為主的求解路徑，從這一點(diǎn)看是非常成功的。但是此定理并沒(méi)有給出求算形式中電子密度函數(shù)、動(dòng)能泛函和交換關(guān)聯(lián)泛函的具體表達(dá)式，這使得我們對(duì)體系的具體求算還是無(wú)法進(jìn)行。直到1965年，Kohn和沈呂九建立了Kohn-Sham方程［11］，得到了各項(xiàng)的具體描述形式，使密度泛函理論開始進(jìn)入實(shí)際應(yīng)用階段。其中Kohn因?yàn)樗贒FT計(jì)算方法上的卓越貢獻(xiàn)榮獲諾貝爾獎(jiǎng)。

1.3 Kohn-Sham 方程

依據(jù)上文的結(jié)論，考慮到體系粒子數(shù)守恒，將原方程變分表達(dá)式等價(jià)為:

Kohn和Sham針對(duì)式(1-7)中仍然沒(méi)有確定具體表達(dá)形式的動(dòng)能泛函項(xiàng)提出使用已經(jīng)具有明確表達(dá)形式的互相無(wú)影響的粒子動(dòng)能泛函來(lái)近似代替，而將二者的差異納入交換關(guān)聯(lián)泛函的未知項(xiàng)中，這樣就使得體系僅剩下一個(gè)未知量，即新的交換關(guān)聯(lián)泛函項(xiàng)Excρ(r[ ])。定義密度函數(shù):

可將動(dòng)能泛函項(xiàng)描述為:

將拉格朗日乘子替換為Ei，根據(jù)變分原理可得:

式(1-8)、(1-10)、(1-11)合稱為 Kohn-Sham方程。通過(guò)選取初始電子密度泛函表達(dá)式，對(duì)體系進(jìn)行自洽迭代運(yùn)算得到更為精確的表達(dá)式，重復(fù)此計(jì)算直到滿足精度要求，即可以求出體系的近似波函數(shù)。Kohn-Sham方程給除了交換關(guān)聯(lián)泛函以外的各項(xiàng)都賦予了明確的表達(dá)形式，并把復(fù)雜作用項(xiàng)也歸并在此項(xiàng)中。至此，求算難度得到大大簡(jiǎn)化，所有的工作都圍繞著如何描述交換關(guān)聯(lián)泛函而展開。

2 交換關(guān)聯(lián)泛函

由上文可知求解問(wèn)題的關(guān)鍵是建立精確和簡(jiǎn)潔的交換關(guān)聯(lián)泛函表達(dá)式。現(xiàn)對(duì)目前常用的方法做出介紹。

2.1 局域密度近似(LDA)

LDA方法也是由Kohn和沈呂九在1965年提出，目的是為了使未知交換關(guān)聯(lián)項(xiàng)能被近似表達(dá)，使得DFT能夠適用于實(shí)際計(jì)算。LDA的核心思想是使用已有的均勻電子氣密度函數(shù)來(lái)求算非均勻電子氣條件下的交換關(guān)聯(lián)項(xiàng)。假設(shè)體系中電子密度隨空間的變化率極小，則此時(shí)非均勻電子氣的交換關(guān)聯(lián)項(xiàng)可以表示為:

εxc[ρ(r)]表示粒子在均勻電子氣中的交換關(guān)聯(lián)能。則相應(yīng)的交換關(guān)聯(lián)勢(shì)可表示為:VLxcDA(r)

可得:

式(2-2)和(2-3)被稱為Kohn-Sham局域密度泛函方程。通過(guò)已上表達(dá)式我們不難看出:

式(2-4)說(shuō)明，交換關(guān)聯(lián)泛函又可以拆解為交換項(xiàng)和關(guān)聯(lián)項(xiàng)兩個(gè)部分來(lái)分別求算。LDA方法求解過(guò)程簡(jiǎn)單，運(yùn)算成本較低，非常適合電子密度變化率較小的系統(tǒng)，在早期求算一般金屬和半導(dǎo)體能帶方面發(fā)揮了巨大作用。但在實(shí)際應(yīng)用中，許多體系的密度函數(shù)變化率都比較大，這使得LDA的在計(jì)算物質(zhì)表面、晶體缺陷等方面誤差較大。為了更加精確的計(jì)算實(shí)際物質(zhì)體系，人們又創(chuàng)立了廣義梯度近似方法。這也是目前密度泛函計(jì)算中運(yùn)用最為廣泛的一類處理方法。

2.2 廣義梯度近似(GGA)

LDA方法中交換關(guān)聯(lián)項(xiàng)僅僅是電子密度的泛函，這使得在密度梯度變化率較大的體系中容易出現(xiàn)較大的偏差。GGA的近似處理正是在充分考慮到LDA近似算法的不足，將電子密度的梯度也作為交換關(guān)聯(lián)泛函的變量，以此對(duì)處在非均勻氣體條件下的體系進(jìn)行能量修正［12］。這就是GGA處理方法的基本思路，具體來(lái)看，就是使原項(xiàng)被改寫為包含電子密度和梯度函數(shù)的泛函形式，若在加上對(duì)自旋的描述，即可得到如下的交換關(guān)聯(lián)泛函:

同樣的，在GGA當(dāng)中，交換關(guān)聯(lián)項(xiàng)也可以拆解為交換能和關(guān)聯(lián)能兩項(xiàng)分別描述。那么，具體的問(wèn)題就是如何構(gòu)建這兩部分的合理表達(dá)式。關(guān)于這個(gè)問(wèn)題，目前主要有兩大觀點(diǎn)。其一是Beckc等人［13］主張的“一切皆合法”觀點(diǎn)。這個(gè)流派認(rèn)為具體的泛函形式原則上可以任意構(gòu)建，并不需要在泛函中考慮實(shí)際的物理意義，那么描述的合理性則完全由實(shí)際計(jì)算來(lái)評(píng)價(jià)，所以此流派的泛函項(xiàng)通常包含由計(jì)算和實(shí)驗(yàn)數(shù)據(jù)統(tǒng)計(jì)優(yōu)化而成的經(jīng)驗(yàn)參數(shù)。在化學(xué)相關(guān)領(lǐng)域應(yīng)用較為廣泛的GGA-PW91泛函就是根據(jù)這個(gè)思路構(gòu)建的。另一類觀點(diǎn)則是Perdew等人［14］主張的交換關(guān)聯(lián)泛函應(yīng)該從實(shí)際的物理意義出發(fā)，盡量回歸到純粹的量子力學(xué)計(jì)算理論，所有物理量的計(jì)算都僅從電子靜質(zhì)量、普朗克常數(shù)、光速等基本常量入手，所以計(jì)算泛函表達(dá)式也不應(yīng)過(guò)多包含經(jīng)驗(yàn)參數(shù)。在凝聚態(tài)物理學(xué)等相關(guān)領(lǐng)域的計(jì)算中經(jīng)常使用的GGA-PBE泛函即是以這種思想作為構(gòu)建依據(jù)的。

大量的算例證實(shí)，GGA通常更加符合實(shí)驗(yàn)結(jié)果，較LDA適用范圍也更加廣泛。除了對(duì)半導(dǎo)體和某些金屬晶格常數(shù)優(yōu)化等少許問(wèn)題還有待解決外，基本上已經(jīng)可以算作密度泛函計(jì)算中非常成熟的算法，所得結(jié)果能夠較好的與現(xiàn)有實(shí)驗(yàn)參數(shù)吻合。

2.3 雜化泛函近似

除了以上兩種基本計(jì)算形式，研究人員還針對(duì)不同體系下的計(jì)算需要做出了很多修正LDA和GGA算法。這使得泛函中包含了描述體系的信息越來(lái)越多，計(jì)算結(jié)果也與實(shí)際實(shí)驗(yàn)數(shù)據(jù)越來(lái)越接近，基本能夠以此對(duì)體系的實(shí)際狀態(tài)做出客觀的衡量。除了LDA和GGA這兩種基本泛函形式，還有考慮到Hatree-Fock形式交換作用影響的雜化密度泛函，各類雜化泛函尤其適用于有機(jī)化學(xué)領(lǐng)域，并在化學(xué)反應(yīng)機(jī)理計(jì)算方面獲得了很大的成功。以文獻(xiàn)中經(jīng)常出現(xiàn)的B3LYP雜化泛函為例:

式(2-6)中涉及到3個(gè)參數(shù)，并且在求算中用到了B88泛函和LYP泛函分別作為交換項(xiàng)和關(guān)聯(lián)項(xiàng)，即以雜化的方式將Hatree-Fock形式交換作用影響納入交換關(guān)聯(lián)項(xiàng)中。但需要指出的是，此類雜化泛函算法由于在類自由電子系統(tǒng)的處理上有所不足，因此對(duì)金屬和窄帶隙半導(dǎo)體材料的計(jì)算會(huì)出現(xiàn)較大的誤差，還有待改進(jìn)［15］。類似的三參數(shù)雜化泛函還有 B3PW91、B3P86等等［16］，這類泛函采取的思路大致相同，只是選取了不同的參數(shù)和泛函種類。

3 展望

近年來(lái)隨著各種表征技術(shù)的迅猛發(fā)展，各類電子顯微鏡和各類電子能譜都已經(jīng)較成熟地運(yùn)用于與材料學(xué)科相關(guān)的各個(gè)研究領(lǐng)域。但是傳統(tǒng)的研究方法仍然具有諸多不足。例如，使用掃描隧道顯微鏡(STM)對(duì)催化材料進(jìn)行表面分析時(shí)，隧道電流會(huì)破壞材料表面的原有結(jié)構(gòu)，同時(shí)也可能引起吸附原子的部分脫離，這些都將導(dǎo)致測(cè)量結(jié)果無(wú)法反應(yīng)材料表面的真實(shí)情況。此外，對(duì)于在研究反應(yīng)機(jī)理中非常重要的瞬態(tài)構(gòu)象信息，傳統(tǒng)的表征手段還不能做到及時(shí)捕捉，這使得很多反應(yīng)機(jī)理只能停留在各種理論推測(cè)的階段。另一方面，傳統(tǒng)的材料設(shè)計(jì)主要是通過(guò)實(shí)驗(yàn)手段來(lái)進(jìn)行，這往往會(huì)消耗較長(zhǎng)的時(shí)間和大量的財(cái)力，并且使得實(shí)驗(yàn)研究可能出現(xiàn)大量的盲目性。

為了解決這些問(wèn)題，越來(lái)越多的研究者把研究方向轉(zhuǎn)入對(duì)材料和反應(yīng)機(jī)理的理論研究，并取得了許多突破性進(jìn)展。DFT計(jì)算方法則是眾多理論研究方法中的一個(gè)重點(diǎn)，它的進(jìn)一步發(fā)展將會(huì)對(duì)突破目前材料學(xué)科的研究瓶頸起到重要作用。但DFT相比其他方法計(jì)算時(shí)間較長(zhǎng)，同時(shí)模型較難運(yùn)用于分子數(shù)較多的體系和含時(shí)系統(tǒng)，尤其是目前對(duì)于納米材料計(jì)算結(jié)果的可靠性也稍遜于模特卡羅方法、分子動(dòng)力學(xué)等方法，這些則是DFT領(lǐng)域發(fā)展中亟待解決的主要問(wèn)題和未來(lái)的主要發(fā)展方向。

［1］Hartree D R.The wave mechanics of an atom with a non-Coulomb central field ［J］.Rroc.Camb.Pil.Soc.，1928，24(1):89-110.

［2］Fock V，Petrashen M J.On the numerical solution of generalised equations of the self-consistent field ［J］.Phys.Z.Sowjetunion，1934，6:368-414.

［3］Hohenberg P，Kohn W.Inhomogeneous electron gas［J］.Phys.Rev.，1964，136(3B):B864-B870.

［4］Kohn W，Sham L J.Self-Consistent equations including exchange and correlation effects［J］.Phys.Rev.，1965，140(4A):A1133-A1138.

［5］Born M，Oppenheimer J R.Zur quantentheorie der molekeln［J］.Ann.Phys.，1927，389(20):457-484.

［6］Van Benthem K，Els?sser C，F(xiàn)rench R H.Bulk electronic structure of SrTiO3:experiment and theory［J］.J.Appl.Phys.，2001，90(12):6156-6164.

［7］Tomohito T，Katsuyuki M，Yuichi I.First-principles study on structures and energetics of intrinsic vacancies in SrTiO3［J］.Phys.Rev.B:Condens.Matter，2003，68(20):205213-1-205213-8.

［8］Arques M，Teles L K，Anjos V，et al.Full-relativistic calculations of the SrTiO3carrier effective mass and complex dielectric function ［J］.Appl.Phys.Lett.，2003，82(18):3074-3076.

［9］Kohn W，Beckc A D，Parr R G.Density functional theory of electronic structure［J］.Phys.Rev.，1965，140:B1133-B1138.

［10］Vignale G，Kohn W.Current-dependent exchange-correlation potential for dynamical linear response theory［J］.Phys.Rev.Lett.，1996，77(10):2037 ～2040.

［11］李震宇，賀偉，楊金龍.密度泛函理論及其數(shù)值方法新進(jìn)展［J］.化學(xué)進(jìn)展，2005 17(2):192-202.

［12］黃美純.密度泛函理論的若干進(jìn)展［J］.物理學(xué)進(jìn)展，2000，20(3):199-219.

［13］Webster F.Finite-size effects in DFT Calculation of The Phonon-Dispersion In crystalline Silicon［J］.Absteracts Of Papers Of The American Chemical Society，1993，206:201 ～inor.

［14］Perdew J P，Burke K，Wang Y.Generalized gradient approximation for the exchange-correlation hole of a many-electron system［J］.Phys.Rev.B，1996，54(23):16533-16539.

［15］Paier J，Marsman M，Kresse G.Why does the B3LYP hybrid functional fail for metals?［J］.J.Chem.Phy.，2007，127(2):024103-1-024103-10.

［16］Ventura O N.Transition states for hydrogen radical reactions:LiFH as a stringent test case for density functional methods［J］.Molecular Physics，1996，89(6):1851-1870.