04實(shí)驗(yàn)四-基礎(chǔ)電腦模擬奈米結(jié)構(gòu)及能態(tài).doc

實(shí)驗(yàn)四 基礎(chǔ)電腦模擬奈米結(jié)構(gòu)及能態(tài)一、 目的：本實(shí)驗(yàn)透過理論計(jì)算化學(xué)套裝程式，讓學(xué)生了解電腦模擬的步驟以及學(xué)習(xí)利用電腦模擬的能力。這個(gè)實(shí)驗(yàn)的設(shè)計(jì)是為了讓學(xué)生學(xué)習(xí)如何應(yīng)用電子的結(jié)構(gòu)模型來研究材料的物理與化學(xué)性質(zhì)。透過這堂實(shí)驗(yàn)，大家將會(huì)學(xué)習(xí)到如何建構(gòu)分子的幾何結(jié)構(gòu)、如何計(jì)算電子性質(zhì)、以及如何模擬光學(xué)性質(zhì)，同時(shí)我們也將學(xué)習(xí)到如何在個(gè)人電腦上以視覺化圖形來顯示計(jì)算的結(jié)果。二、 原理：量子化學(xué)中的兩個(gè)主要分子軌域理論近似方法分別為：半經(jīng)驗(yàn)分子軌域理論以及全始算之Hartree-Fock分子軌域理論。這兩個(gè)近似法都是以電子之間無交互作用或是電子彼此獨(dú)立的設(shè)計(jì)來架構(gòu)計(jì)算的出發(fā)點(diǎn)，為了改進(jìn)這樣設(shè)計(jì)的誤差，因此便引進(jìn)電子相干的理論。而另一方面，另外一種研究材料的電子與幾何結(jié)構(gòu)的近似法則是使用密度泛函數(shù)理論的概念。在這個(gè)理論裡，電子的交換與相干作用被視為是電子密度的泛函數(shù)，由於密度泛函數(shù)理論的計(jì)算規(guī)模較低，因此更適用於計(jì)算大型的材料系統(tǒng)。根據(jù)量子力學(xué)觀點(diǎn)，一個(gè)分子在某一穩(wěn)定狀態(tài)（stationary state）的能量、電子狀態(tài)及一切物理性質(zhì)皆可藉由薛丁格方程式求得。描述一個(gè)分子系統(tǒng)的特定狀態(tài)的薛丁格方程式可表示為 = E（2-1）其中為系統(tǒng)的漢米頓算子（Hamiltonian operator），對(duì)其右邊的波函數(shù)運(yùn)算可以得到一特徵值（eigenvalue）E，即系統(tǒng)的能量。為假設(shè)存在的一個(gè)描述系統(tǒng)狀態(tài)的函數(shù)，稱之為波函數(shù)（wave function，或狀態(tài)函數(shù)state function）。一個(gè)完整的漢米頓算子必須包含系統(tǒng)的動(dòng)能與位能，可以將其表示為（2-2）系統(tǒng)中粒子的動(dòng)能為（2-3）系統(tǒng)中帶電粒子的位能為（2-4）對(duì)於一個(gè)分子系統(tǒng)，動(dòng)能項(xiàng)部分包含電子與原子核的動(dòng)能，位能項(xiàng)則包含以核與核、電子與電子，及電子與核之間的作用力，除此之外，還必須包含電子自旋和各個(gè)軌域電子運(yùn)動(dòng)之間的相互作用力，而且當(dāng)電子運(yùn)動(dòng)速度比較高時(shí)還得考慮其相對(duì)論效應(yīng)。由這種完整的分子漢米頓算子所構(gòu)成的薛丁格方程式太過複雜，即使是現(xiàn)在的大型工作站也難以對(duì)較大的分子進(jìn)行計(jì)算。為簡化計(jì)算過程，其中至少忽略了二個(gè)效應(yīng)，一個(gè)是相對(duì)論效應(yīng)，雖然氫原子中電子的速度小於光速的1%，但對(duì)於重原子之內(nèi)層電子的相對(duì)質(zhì)量校正是相當(dāng)可觀的，另一個(gè)是忽略核與電子間的電子自旋和軌域運(yùn)動(dòng)相互作用（spin-orbit）效應(yīng)。如此一來，一個(gè)分子系統(tǒng)的漢米頓算子可以表示為（2-5）其中MA為原子核的重量，me為電子的重量，Z為原子核所帶的電荷，h為蒲朗克常數(shù)（Plancks constant），e為電子所帶的電荷，RAB為核與核之間距離（RAB =| rA - rB |），rAi為核與電子間距離（rAi =| rA ri |），rij為電子與電子間的距離（rij =| ri rj |）。漢米頓算子的動(dòng)能項(xiàng)部分為式2-5中核與電子動(dòng)能的和（2-6）位能項(xiàng)則為式2-5中電子-核、核-核及電子-電子間靜電位能的和（2-7）量子化學(xué)中的一些方程式，為了要簡化其呈現(xiàn)的方式，時(shí)常將一些常數(shù)項(xiàng)目設(shè)定成基本單位。在原子單位的表達(dá)方式中，長度是以波耳半徑（Bohr radius）為單位：= 0.529177249 （2-8）一般座標(biāo)可藉由除上a0以轉(zhuǎn)換成波耳的表示法，能量則以hartree為單位，其大小定義為兩個(gè)電子相距1波耳距離的庫倫排斥能1 hartree = 27.2114 eV 627.510 kcal/mol（2-9）質(zhì)量則是以一個(gè)電子的質(zhì)量為單位。如此一來先前的漢米頓算子則可改寫為（2-10）由於原子核的質(zhì)量是電子質(zhì)量的數(shù)千倍，因此相對(duì)於電子來說，原子核的移動(dòng)是非常緩慢的，電子瞬間的交互作用幾乎感受不到原子核位置的改變，所以可將原子核與電子的運(yùn)動(dòng)視為相互獨(dú)立，故將分子內(nèi)電子的分佈情形視為僅與原子核間的瞬間距離有關(guān)，而與原子核的運(yùn)動(dòng)無關(guān)，這就是所謂的波恩-歐本海莫近似法（Born-Oppenheimer approximation），主要的假設(shè)便是視原子核為固定不動(dòng)，電子就像是在一群固定的原子核間運(yùn)動(dòng)一樣，因此接下來的薛丁格方程式便僅需描述電子運(yùn)動(dòng)與分佈情況，也就是成為一個(gè)只與電子狀態(tài)有關(guān)的薛丁格方程式。ee = Eee（2-11）而簡化後的漢米頓算子則變?yōu)椋?-12）其中r表示所有電子的座標(biāo)，而R表示所有原子核的座標(biāo)。描述原子核動(dòng)能的項(xiàng)被忽略是因?yàn)檫h(yuǎn)小於，另外在計(jì)算過程中亦被省略，因?yàn)樵谠雍凸潭ú粍?dòng)的情況下，核與核間的靜電位能可被視為是一個(gè)常數(shù)，它只對(duì)能量特徵值（eigenvalue）變動(dòng)一固定的量。式2-12中，電子的動(dòng)能為（2-13）電子-電子之間的靜電位能為（2-14）電子-核之間的靜電位能為（2-15）因此，漢米頓算子可表示為式2-16，由此可求出系統(tǒng)中所有與電子有關(guān)的能量Ee，接下來在系統(tǒng)總能量（total energy）（指的是固定核位置下的總能量，即忽略核的動(dòng)能）的部份中，必須再把視為常數(shù)的核-核間靜電斥力當(dāng)作一附加項(xiàng)，加在只與電子狀態(tài)有關(guān)的薛丁格方程式所求出的電子能量Ee之後，即可計(jì)算出系統(tǒng)的總能（式2-17）。（2-16）（2-17）在量子力學(xué)上，波函數(shù)的平方、代表粒子分布的機(jī)率密度，具有統(tǒng)計(jì)上的意義，但是真實(shí)在計(jì)算上除了類氫的單電子系統(tǒng)（如H2+、He+、Li+等）可得到正確的波函數(shù)解之外，根據(jù)測(cè)不準(zhǔn)原理，多電子分子系統(tǒng)中電子為瞬間相對(duì)運(yùn)動(dòng)，存在於漢米頓算子中的各電子間靜電排斥位能項(xiàng)，無法分離變數(shù)而分別求得，我們僅能求得合理的近似解，使其結(jié)果僅可能地接近薛丁格方程式的真實(shí)解，因此大部份的計(jì)算都是採用分子軌域法來求其近似解。現(xiàn)在考慮電子波函數(shù)，因?yàn)殡娮硬ê瘮?shù)中R為定數(shù)，所以我們不再標(biāo)出R這個(gè)註記，僅以r表示所有電子的座標(biāo)，例如，。在此將電子波函數(shù)表示為。由於電子是自旋量子數(shù)為的費(fèi)米子（Fermion），因此電子波函數(shù)必須滿足一個(gè)重要的條件，就是所謂的反對(duì)稱原理（antisymmetry principle）。此原理說明在一個(gè)系統(tǒng)中的費(fèi)米子中，任兩個(gè)費(fèi)米子座標(biāo)（空間和自旋）互相交換時(shí)其波函數(shù)必需互為反對(duì)稱的。也就是（2-18）鮑利不相容原理（Pauli exclusion principle）就是反對(duì)稱原理的結(jié)論。電子的自旋軌域（spin orbital，(x)）就是單電子之空間和自旋座標(biāo)的函數(shù)，而空間軌域（spatial orbital，(r)）則只是單電子之空間座標(biāo)的函數(shù)。我們可將自旋軌域?qū)懗煽臻g軌域和自旋函數(shù)的乘積（2-19）或（2-20）對(duì)於一個(gè)空間軌域(r)，可以形成二個(gè)自旋軌域，一個(gè)對(duì)應(yīng)自旋，一個(gè)對(duì)應(yīng)自旋。所以一個(gè)空間軌域可填入兩個(gè)電子。分子軌域理論是利用自旋軌域所組成的函數(shù)式來描述一組電子的分布。此自旋軌域包含空間函數(shù)p(x,y,z)與自旋函數(shù)g(, )，因此在描述2n個(gè)電子的分子系統(tǒng)時(shí)，就有n個(gè)自旋軌域函數(shù)來組成合適的多電子函數(shù)以解薛丁格方程式，求出最穩(wěn)定的能量描述2n個(gè)電子之分子系統(tǒng)，處於封閉殼層（closed-shell）基態(tài)（ground state）的多電子波函數(shù)可由Slater行列式來表示：（2-21）上式中每一個(gè)單電子波函數(shù)可以分成由兩個(gè)部份的乘積所組成，一部份與電子的空間座標(biāo)有關(guān)，另一部份則是與電子的自旋座標(biāo)有關(guān)，電子的自旋狀態(tài)分別以、表示，整個(gè)分子系統(tǒng)只需要n個(gè)不同的空間座標(biāo)函數(shù)來描述。在實(shí)際計(jì)算中，分子軌域乃採用以n個(gè)已知的單電子軌域函數(shù)線性組合（linear combination）而成，即s = 1,2,3,b（2-22）s是單電子的基底函數(shù)（basis function），所謂基底函數(shù)指的是原子軌域波函數(shù)或一組描述原子軌域的波函數(shù)組，這個(gè)方法就是一般常用的MO-LCAO（molecular orbital-linear combination of atomic orbitals），有了由原子軌域線性組合的分子軌域之後，只要再根據(jù)變分原理（variational principle），利用矩陣形式來解出時(shí)間無關(guān)薛丁格方程式（time-independent Schrdinger equation），即可求出系統(tǒng)的能量。Hartree-Fock方程式對(duì)於一個(gè)封閉殼層（closed-shell）組態(tài)的分子，我們選用Slater行列式的自旋軌域反對(duì)稱積（antisymmetrized product of spin-orbital）作為描述其電子的波函數(shù)。根據(jù)變分原理，電子在基態(tài)時(shí)的波函數(shù)，其能量的期望值應(yīng)趨於最低，對(duì)此Slater行列式的多電子波函數(shù)，其近似解的能量期望值應(yīng)符合此原理。為了簡化問題我們將分子系統(tǒng)的電子漢米頓算子分離成單電子跟雙電子算子兩部份 = 1 + 2（2-23）（2-24）（2-25）1：表示一個(gè)電子在原子核整體的靜電位能場(chǎng)下的能量算子，一般稱為核漢米頓算子（core Hamiltonian operator，core）。2：表示各個(gè)電子間的靜電斥力位能。此時(shí)電子總能量的期望值E為：（2-26）將Slater行列式的波函數(shù)代入式2-26中可得到：（2-27）（2-28）（2-29）（2-30）式2-28代表一個(gè)電子在之分子軌域中core漢米頓算子的能量期望值，式2-29與式2-30代表著電子1跟電子2在分子軌域與的位能積分項(xiàng)，式2-29稱之為庫倫積分（coulomb integral），式2-30稱之為交換積分（exchange integral）。式2-27在經(jīng)過變分法處理後，使E產(chǎn)生最小值的計(jì)算過程中，Hartree-Fock方程式必須滿足下式：（2-31）此即為所謂的Hartree-Fock方程式，式2-31中i為軌域能量，有效能量算子即稱為Fock算子，其中（2-32）庫倫算子（coulomb operator）與交換算子（exchange operator）定義如下：（2-33）（2-34）其中f為任意函數(shù)且積分是定義為全域積分。式2-32的Fock算子可以視為一個(gè)電子在其分子軌域中的有效算子，將電子所受到其它電子靜電斥力的效應(yīng)，以分子軌域近似表示為電子密度分佈所產(chǎn)生的平均場(chǎng)位能，藉以分別解出分子軌域最適當(dāng)?shù)牟ê瘮?shù)形式。在求解的過程中，首先必須假設(shè)一組嘗試解的分子軌域，依此假設(shè)軌域去計(jì)算初步近似的Fock算子，再以此Fock算子解得一組更適當(dāng)?shù)姆肿榆売?，如此反覆?jì)算，即可求出幾乎與真實(shí)系統(tǒng)完全一致的Fock算子與一組最佳解的分子軌域，進(jìn)而求出系統(tǒng)的總能量，這就是所謂的自洽場(chǎng)理論（self-consistent field theory, SCF）。若將分子軌域依式2-22以N個(gè)原子軌域函數(shù)作線性組合，則原Hartree-Fock方程式可以轉(zhuǎn)化為求各分子軌域線性展開係數(shù)的聯(lián)立方程式，稱為Roothaan方程式。（2-35）兩邊同乘上r*並積分後，如此一來Hartree-Fock方程式可改寫為：,r = 1,2,b（2-36）,（2-37）為了解Roothaan方程式，我們首先必須以基底函數(shù)來表示矩陣元素Frs，因此將式2-32的Fock算子代入式2-36中，整理後得到（2-38）其中Hrs為core漢米頓算子所對(duì)應(yīng)的矩陣元素，而Ptu與雙電子排斥積分項(xiàng)(rs|tu)（two-electron repulsion integral）定義如下 ,t = 1,2,b ,u = 1,2,b（2-40）（2-41）Ptu被稱為電子密度矩陣（electron-density matrix elements）或電荷密度（charge-density）、鍵級(jí)矩陣（bond-order matrix elements），由於電子的機(jī)率密度為軌域波函數(shù)的平方和，若代入式2-22的基底函數(shù)則可得到封閉殼層分子的電子機(jī)率密度為（2-42）在式2-37中，由於Frs本身定義包含有Ptu，而Ptu為ctj與cuj乘積的和，因此並無法直接求得其解，故一般仍須採用反覆疊代的方式來求解。因此首先必須假定出一組近似解csi，以此計(jì)算初步近似的Frs矩陣以及總能E，如式2-43所示，再對(duì)Roothaan方程式進(jìn)行求解，即可得到更接近真實(shí)解的一組近似解，並得到另一組分子軌域係數(shù)csi，以自洽場(chǎng)過程計(jì)算求得一組可以被接受的近似解。（2-43）在整個(gè)計(jì)算過程中，對(duì)於所有需要的積分運(yùn)算均不做任何的近似與省略，故被稱為全始算（ab initio）。密度泛函數(shù)理論（Density Functional Theory）在1964年時(shí)Hohenberg和Kohn的理論證明了基態(tài)時(shí)分子的能量、波函數(shù)以及其他分子的電子性質(zhì)都可由電子的機(jī)率密度(x,y,z)唯一決定。換句話說，分子在基態(tài)時(shí)的能量E0為的泛函數(shù)，即E0 = E0()（2-44）在傳統(tǒng)的全始算理論中，只要給予一個(gè)系統(tǒng)所含的原子和電子數(shù)，就可以藉由解薛丁格方程式而知道這個(gè)系統(tǒng)的性質(zhì)，在薛丁格方程式中有3N個(gè)變數(shù)（N =總粒子數(shù)），如果系統(tǒng)所含的粒子數(shù)很多，這個(gè)方程式將變的極為複雜。跟傳統(tǒng)全始算理論不同的是，不論粒子數(shù)多寡，粒子密度函數(shù)都只有三個(gè)空間變數(shù)，如此問題就變的簡單多了。因?yàn)槿绻阑鶓B(tài)的電子密度(x,y,z)，則由Hohenberg-Kohn理論，理論上可以求得基態(tài)分子的所有性質(zhì)。這與以往必須先求出波函數(shù)，然後才能計(jì)算分子性質(zhì)的方法不同。其後在1965年，Kohn及Sham提出了目前學(xué)界DFT計(jì)算應(yīng)用上最常採用的版本(Kohn-Sham DFT; W. Kohn的貢獻(xiàn)也使其於1998年獲得諾貝爾化學(xué)獎(jiǎng))，即對(duì)一個(gè)給定的原子、分子或固態(tài)真實(shí)系統(tǒng)，可透過找尋一個(gè)虛擬的無電子交互作用的系統(tǒng)(non-interacting system)來產(chǎn)生真實(shí)系統(tǒng)的電荷密度分布，進(jìn)而再得到真實(shí)系統(tǒng)的基態(tài)能量：（2-45）其中i稱為Kohn-Sham orbitals、i = 1,2,n，Exc即為所謂的交換-相干能（exchange-correlation energy）。正確的基態(tài)電子密度函數(shù)可以由i組成：（2-46）這是因?yàn)樵跓o電子間交互作用的系統(tǒng)中(解真實(shí)系統(tǒng)本徵態(tài)薛丁格方程式ee = Eee的困難處即在於Hamiltonian中有電子間的交互作用項(xiàng))，本徵態(tài)多體波函數(shù)將變得非常簡單，即是解單體薛丁格方程式得到單體波函數(shù)(orbitals)後再取其中任意Ne(電子數(shù)目)個(gè)建構(gòu)出的Slater determinant；其電荷密度分布也可直接表達(dá)如(2-46)的形式。此虛擬系統(tǒng)(Kohn-Sham system)的單體波函數(shù)稱作Kohn-Sham orbital，可以由下列單電子本徵態(tài)薛丁格方程式求得（2-47）其中，即是Kohn-Sham虛擬系統(tǒng)的Hamiltonian，（2-48）為庫倫算子定義於式2-33中Vxc稱作交換-相干位能（exchange-correlation potential）以Exc泛涵數(shù)的變分表示（2-49）如果Exc已知，則由式2-49可求得Vxc，經(jīng)由這一項(xiàng)來處理所有有關(guān)電子交換與電子相干的效應(yīng)。將式2-48代入式2-47並使用式2-33與式2-46可以得到（2-50）首先假設(shè)一個(gè)起始的，由式2-49可以求得Vxc，再將之代入式2-50，解出一組i，在由式2-46得到一個(gè)新的。如此重複這相同的步驟一直反覆計(jì)算到函數(shù)收斂為止，最後再利用式2-45求得基態(tài)能量E0。所以在解Kohn-Sham方程式的過程和解Hartree-Fock一樣都是屬於SCF的方法。在Hartree-Fock方法中，是以近似的波函數(shù)（approximation wave function），即假設(shè)波函數(shù)是單一的Slater Determinant，去解正確的漢米頓算子。透過基態(tài)能量對(duì)每個(gè)組成的orbital變分找極值的原理，找出一組真實(shí)系統(tǒng)波函數(shù)的近似解。而在Kohn-Sham-DFT中，是透過找尋一個(gè)虛擬系統(tǒng)，其電荷密度分佈與真實(shí)系統(tǒng)相同，再由此密度可單一決定並解出真實(shí)系統(tǒng)的其餘物理性質(zhì)(根據(jù)Hohenberg-Kohn Theorem)，因此省去直接解真實(shí)系統(tǒng)波函數(shù)的困難。然而雖然原理上(in principle)若擁有正確的交換-相干位能，所解得的虛擬系統(tǒng)的基態(tài)(ground state)電荷密度將與真實(shí)系統(tǒng)基態(tài)電荷密度完全相同，但理論上交換-相干項(xiàng)位能（effective exchange-correlation potential）的正確形式依然是未知，需要做近似，因此衍生出非常多樣化的近似方法(LDA,GGA,B3LYP,等)。在DFT的理論方法中，通常在交換-相干能的計(jì)算上都把Exc分開成交換泛函數(shù)Ex及相干泛函數(shù)Ec兩項(xiàng)，即Exc() = Ex() + Ec()（2-51）在目前各種DFT的計(jì)算程式中均會(huì)提供數(shù)種的交換泛函數(shù)與相干泛函數(shù)，將任一組交換泛函數(shù)與任一組相干泛函數(shù)組合為完整的交換-相干泛函數(shù)以適用於各種不同的分子系統(tǒng)。例如，在Gaussian 03的計(jì)算程式中分別提供了Slater、Xa、Becke 88等數(shù)種交換泛函數(shù)以及VWN、VWNV、LYP、Perdew Local、Perdew 86、Perdew/Wang 91等多種相干泛函數(shù)。三、 儀器與材料硬體電腦叢集伺服器(Xeon 2.4GHz CPU*14)，P4個(gè)人電腦(P4 2.4GHz CPU)軟體HyperChem Pro 7.5，Gaussian 03，GaussView 03四、 實(shí)驗(yàn)步驟（一）分子幾何結(jié)構(gòu)的建置-在HyperChem的視窗環(huán)境中，選擇要加入的原子型態(tài)(圖一)，接著在工作視窗中放置原子(圖二)，並利用滑鼠拖曳的方式在兩個(gè)原子間拉出連結(jié)的直線(鍵結(jié))。 圖一圖二圖四圖三圖六（二）幾何結(jié)構(gòu)預(yù)最佳化-選擇以下的半經(jīng)驗(yàn)分子軌域理論方法進(jìn)行幾何結(jié)構(gòu)最佳化：AM1、PM3。(圖五)-在能量到達(dá)收斂之後，將此結(jié)構(gòu)儲(chǔ)存為PDB (ENT)的檔案格式。(圖六) 圖五（三）利用DFT方法的完全幾何結(jié)構(gòu)最佳化-在GaussView的視窗環(huán)境中，開啟之前所儲(chǔ)存的預(yù)最佳化幾何結(jié)構(gòu)檔案。-選擇計(jì)算的工作類型為Opt-Freq。(圖七)圖七圖八-選擇理論方法與基底函數(shù)組(basis set)分別為B3LYP與6-31G(d) (也可嘗試其他的基底函數(shù)組)。(圖八)-開始最佳化。-在能量到達(dá)收斂之後，檢查其振動(dòng)頻率是否都為正值。（四）電子相干計(jì)算-在GaussView的視窗環(huán)境中，開啟步驟三所完成的完全最佳化結(jié)構(gòu)。-選擇計(jì)算的工作類型為Opt-Freq。-選擇理論方法與基底函數(shù)組分別為MP2與6-31G(d)。-開始最佳化。-在能量到達(dá)收斂之後，檢查其振動(dòng)頻率是否都為正值。-將所得之能量與DFT的結(jié)果作比較