密度泛函理論(五):實際計算中的問題

（絕緣體，金屬和半導體中價帶，導帶和費米能級的關係）

在前面的文章中，我們介紹了利用Bloch理論、傅立葉變換和平面波展開的方法，在倒易晶格中求解各個k點處波函數與電子能量的求解方法。在實際應用中，由於計算量和波函數函數積分的問題，我們還需要考慮能量截斷和模糊處理(smearing)。然後，本文將簡要介紹密度泛函理論在實際計算中的算法和第一性原理-分子動力學方法。

能量截斷

在每個k點處的波函數可以用平面波展開，原則上，這些平面波可以有無窮多個，但是在實際計算中，平面波的個數必須為有限個，因此，我們需要對平面波的能量進行截斷。高能量平面波的分量對波函數的貢獻很小，因為高能量平面波的展開係數很小，因此可以對平面波的能量進行截斷。

在前面的文章中，我們計算過平面波的動能為：

截斷能Ecut定義為：

在對每個k點的平面波展開時，我們只考慮能量低於Ecut的分量。截斷能的具體取值與計算的系統有關，一般來說，我們需要重複計算越來越大的截斷能，知道體系的能量趨於穩定，進而確定一個體系的截斷能。

截斷能的取值與密度泛函理論的計算量有直接關係，根據截斷能的定義，我們可以估算展開波函數所需要的平面波與體系大小和截斷能的關係，如下式所示：

通常地，每個原子的波函數約需要上百個平面波展開，對於一般的體系，一般需要1000~100000個平面波對波函數進行展開，這就意味著，需要計算與平面波個數對應的矩陣的對角化，尤其是當體系比較大時，計算量將大大增加。

模糊處理

對於一個孤立的原子，該原子的電子會佔據特定的離散的軌道，比如: 1s, 2s, 2p, 2d...等。當很多原子形成固體時，非常靠近原子核的電子不參與成鍵，仍然保持在原有的能量軌道上，而參與成鍵的價電子都試圖佔據同一能級，n個原子將會有n個相同能級的軌道，因此，多原子體系中同一能級有更多的狀態，進而形成能帶。如下圖所示:

費米能級定義為，體系處於基態時，電子所佔據的能量最高的能級。費米能級與導帶和價帶的位置關係決定了固體材料的導電性。例如，金屬材料的價帶很寬，導致費米能級位於導帶中，因此，金屬大多具有很好的導電性。而絕緣體和半導體在費米能級附近，電子的狀態密度很低，存在較大的能隙，因此他們導電性較差。

對於絕緣體或者半導體而言，因為態密度在價帶頂部逐漸衰減，電子密度也逐漸衰減，直至能隙下端，因此在每個k點對波函數進行積分時可以直接計算。但是對於金屬來說，由於部分電子會在導帶上分布，態密度在接近費米能級時不會逐漸衰減，而且，電子在費米能級附近的狀態佔據數會直接由1變為0，因此，電子密度在費米能級處會急劇衰減，由於電子密度是波函數在第一布裡淵區的某種積分，電子密度急劇變化時，波函數在利用平面波展開時的精度就會大大下降。因此，需要對金屬費米能級附近的能級的狀態的佔據數進行模糊處理，如下圖所示：

一般而言，狀態佔據數的模糊方法有：高斯模糊 (Gaussian smearing)， 費米模糊 (Fermi smearing)和Methfessel-Paxton模糊。

高斯模糊是Fu和Mo在1983年提出的，其主要原理是利用高斯函數虛構了一個電子溫度，等效對體系稍微進行加熱，電子溫度為0.1 eV，此時在費米能級附近的電子的能量範圍就會稍微拓寬，不會出現驟降。

費米模糊利用的原理類似，不同的是費米模糊採用Fermi-dirac函數進行有限溫度的虛構，如下式所示：

費米模糊中的電子溫度為KBT, 也就是25 meV約為300 K。

Methfessel-Paxton模糊由Methfessel和Paxton在1989年提出，是在實際中計算中使用最多的模糊方法。其主要原理是對step函數進行多項式和高斯函數展開。

實際算法

密度泛函理論的計算需要重複處理非線性特徵值問題，對於給定計算能力的計算機。計算效率就取決於實際使用的算法。對於電子密度和能量的計算其實耗費很少的計算資源，而K-S矩陣的對角化基本上佔據了70%的計算資源，合理的選擇計算算法對提高計算效率很重要。

電子能量最小化指保持離子的位置不變，使電子達到基態。對於整個體系的能量最小化，需要在電子能量最小化後，進行離子能量的最小化。電子能量最小化的方法包括：直接對角化方法(Direct diagonalization)，迭代戴維森算法(Iterative Davidson method)和RMM-DIIS方法。直接對角化方法不適合很大的矩陣的對角化，也就是有很多平面波(10000~100000)需要展開的時，使用直接對角化方法效率會很低。

迭代戴維森算法的計算效率比較慢，但是有很好的可靠性，也就是不同的初始迭代條件，最終都會收斂到基態能量。RMM-DIIS算法的計算效率很快，也可以進行並行運算。但是，如果RMM-DIIS算法在某些初始迭代條件下，會出現無法收斂到基態能量的情況，也會丟失特徵值。因此，在實際DFT計算中，往往首先採用迭代戴文森算法，然後使用RMM-DIIS算法。

離子能量的最小化指，在進行電子能量最小化後，計算每個離子的受力，然後將離子能量最小化，在連續進行電子能量最小化和離子能量最小化後，整個體系的能量將達到最小化。

Hellmann-Feynman力: 在進行離子能量最小化的時候，離子受到的力來自其他離子和電子，Hellmann-Feynman指出，離子的受力等於哈密頓量關於原子位置的偏導數的平均值，因此離子的受力如下：

在實際的密度泛函理論計算中，離子的受力計算效率很高，相比電子能量的最小化而言其計算時間幾乎是可以忽略的。得到離子的受力後，可以採用常規的共軛梯度法，近似牛頓法或者阻尼分子動力學方法進行最小化。

第一性原理-分子動力學方法

在經典分子動力學方法中，我們採用經驗化的勢函數，這種勢函數有兩個缺點：其一、不同的情形下我們需要採用不同的勢函數和參數來計算相同的體系的相互作用。其二、體系的電子和磁性特性無法被計算。如果我們能夠把密度泛函理論和分子動力學耦合起來，這些缺點將會被克服。

兩種常見的第一性原理-分子動力學方法包括：Born-Oppenheimer分子動力學方法(BOMD)和Car-Parrinello分子動力學方法(CPMD)。

BOMD方法可以看做是前面講到的離子能量最小化方法的推廣，其主要步驟如下：

1、固定原子核的位置，求解電子能量最小化；

2、計算Hellmann-Feynman力；

3、按照經典牛頓運動力學移動原子，例如速度verlet算法；

4、重複以上步驟，直至系統達到穩定。

這個方法在實際計算中的效率很慢，因為每一步分子動力學都需要對電子能量最小化，而電子能量最小化非常耗費計算資源。

CPMD方法將變分最小化應用到電子基態的尋找中，而變分最小化方法可以由模擬退火方法實現。這種方法將電子和離子的自由度都耦合到經典坐標體系中。CPMD方法在計算金屬時，由於價電子的移動會出現電荷波動。儘管如此，CPMD方法仍然是很成功的一種第一性原理分子動力學方法。

分子模擬交流學習，請掃碼關注