CuB2O4晶體在第一原理的研究

全文

(1)國立臺灣師範大學物理學研究所碩士學位論文. 指導教授：陳穎叡博士. CuB2O4 晶體在第一原理的研究 Ab initio study of the anti-ferromagnetic, non-collinear CuB2O4 crystal. 研究生 : 王俊凱撰. 中華民國一 ○ 八年一月.

(2) 摘要研究 CuB2O4 反鐵磁性材料在第一原理中的能帶架構以及材料特性。CuB2O4 在第一原理中，使用六個單元的方式呈現週期性的結構分布，在這個晶格中 a=11.357、c=5.542 Å , k-mesh 使用 4x4x4, nbands 為 400 並且使用 VASP GGA+U 的方法，其中 Cu(A)&Cu(B)的 U 分別為 U1 = 6.4eV, U2 = 7.0eV。在計算中，我們發現 CuB2O4 的 density of state(DOS)在導帶的地方會出現兩個非常明顯的尖點。一個是由 Cu(A)&O(1) 組成，另外一個則是由 Cu(B)&O(2)&O(3)&O(4)所貢獻，因此我們希望可以分析出兩者的不同。關鍵詞:第一原理、晶體結構、能帶分析 Ab initio study of anti-ferromagnetic non-collinear CuB2O4 crystal , CuB2O4 crystallizes in first-principles calculations, with six formula units were performed using periodic density function theory. The cell of dimensions a=11.357 and c=5.542 Å , kmesh=4x4x4, nbands=400 and use GGA+U on VASP where U1=6.4eV, U2=7.0eV for Cu(A)&Cu(B).In the case, we found that CuB2O4’s density of state(DOS) has two peaks in conduction band. One is contributed by Cu(A) & O(1), the other one is contributed by Cu(B) & O(2) & O(3) & O(4). As the result of, we can sort atoms to two types. A type is Cu(A) and O(1) which donate states to peak(1), the other type is Cu(B) and O(2)&O(3)&O(4) which donate states to peak(2). And we could do analyze their difference. Key Words: Ab initio、VASP、first-pronciple、CuB2O4、non-collinear calculation.

(3) 目次 Chapter 1 緒論………………………………………………………………………1 Chapter 2 密度泛函理論(DFT)和計算方法………………………………………..2 2.1. 密度泛函理論(Density function theory)……………………………………..2. 2.1.1 The Hohenberg-Kohn theorems………………………………………….2 2.1.2 The Kohn-Sham equations……………………………………………….4 2.1.3. Exchange-correlation energy……………………………………………..6. 2.2 GGA+U method………………………………………………………………7 Chapter 3 CuB2O4 的文獻探討…………………………………………………….8 3.1 CuB2O4 在基態時的鐵磁性研究(Ferromagnetism in CuB2O4)…………...8 3.2 CuB2O4 在 VASP 中的鐵磁性計算研究…………………………………..11 Chapter 4 CuB2O4 在 VASP 中的材料計算研究及討論………………………...15 4.1 CuB2O4 在 VASP 中的自洽(self-consistent field method, SCF)材料計算..15 4.2 CuB2O4 在 VASP 中使用 GGA+U 方法的計算…………………………..16 4.2.1. 在 VASP 測試 U 值的大小…………………………………………….16. 4.2.2 CuB2O4 在 VASP 中的收斂(relaxation)計算…………………………19 4.3 CuB2O4 材料在導帶上的特徵分析………………………………………..24 Chapter 5. 結果與討論……………………………………………………………..29. 參考文獻……………………………………………………………………………..31.

(4) Chapter 1. 緒論. 第一原理計算(First principle calculation)是以量子力學為基礎的計算方法。根據量子力學理論來計算物質原子核與電子間的交互作用，進而得到材料的各種性質。隨著 Kohn 和 Sham 提出了一種有效的近似方法，解決了多體量子系統長久以來計算結果不理想的問題後，第一原理計算方法以及密度泛函理論(DFT)才逐漸受到重視。隨著科技的進步，高速電腦(HPC)的計算能力有跳躍式成長，大幅縮短了理論計算所需的時間，也被允許進行更加精確的材料計算。近幾年來，銅-氧(Cu-O)化合物因為低溫時的磁性特徵而開始被科學家所關注。根據 M. Boehm 論文中的實驗結果指出，CuB2O4 材料會在 T <10K 時，磁性結構由 commensurate 轉變為 incommensurate 的結構。在這篇文章中，我們希望可以提供詳細的 CuB2O4 材料分析，因此我們使用 VASP 5.4.1 版本來進行第一原理計算。我們將 CuB2O4 計算後的結果配合 J.Y Lin 教授的光譜實驗結果來與 M. Boehm 以及 R.V. Pisarev 的實驗數據來做比對。我們將焦點集中在分析 CuB2O4 在導帶兩個尖銳且清晰的 peaks，我們發現在第一個 peak 上主要是 Cu(A) & O(1)的貢獻，而在第二個 peak 上主要是 Cu(B) & O(2、3、4)所貢獻的 states，並藉由 VESTA 軟體畫出兩個 peaks 的局域電子雲分布圖，來表現出我們的計算成果。. 1.

(5) Chapter 2. 密度泛函理論(DFT)和計算方法. 2.1 密度泛函理論(Density function theory) 密度泛函理論是一種研究多體量子系統(many-body system)的方法，此方法透過量子力學的計算，解出物質系統的電子總波函數及基態總能量，在物理及化學上常用來研究原子、分子、凝態物理材料的性質。基於 Hohenberg-Kohn 的理論，電子密度函數與這些系統的組成結構，具一對一的對應關係，而系統的基態能量，以至於基態的所有性質，皆僅為電子密度的泛函。密度泛函理論由此開始，但卻直到 Kohn 和 Sham 提出了一種有效的近似方法來描述，才開始被普遍應用。在 Kohn-Sham 方法中，電子密度仍以單電子波函數為元素所構成，因此，解出對密度泛函變分所產生的 Kohn-Sham equation，即得到單電子波函數，而各單電子波函數又產生總電子密度函數，此函數又進到 Kohn-Sham equation 裡面的電子交互作用中，如此周而復始，以迭代的方式進行直至收斂。. 2.1.1 The Hohenberg-Kohn theorems 密度泛函理論起源於 Thomas 和 Fermi 在 1927 年提出的方法，將能量分成動能及位能兩項來求得𝐸𝐺.𝑆. [𝜌(𝑟)]的近似：藉由庫倫定律求得電荷在空間中的靜電位能，利用自由電子氣的密度泛函來近似動能項，如此一來我們便具備初步的泛函可以求出了𝐸𝐺.𝑆. [𝜌(𝑟)]。但是 Thomas-Fermi 的近似過程中過分簡化了電子跟電子間的交互作用，導致計算結果不理想，因此 1928 年 Dirac 在 Thomas-Fermi 模型泛函中增加了交換關聯項。密度泛函理論直到 Hohenberg-Kohn 定理提出後才有了堅實的理論基礎。 Hohenberg-Kohn 定理證明了多體系統中與基態相關的兩件事。理論一：薛丁格方程式中的基態能量是電荷密度的泛函。 𝐸𝐺.𝑆. [ψ(𝑟1 , 𝑟2 , … … , 𝑟𝑛 )] = 𝐸𝐺.𝑆. [𝜌(𝑟)] 2.

(6) 理論二：得到最低總能量的電荷密度分布，該能量就會是基態的總能量，且該電荷密度分布就是基態的電荷密度分布。意即，若將非基態的電荷密度𝜌(𝑟) 帶入，會得到E[𝜌(𝑟)] > 𝐸𝐺.𝑆. [𝜌𝐺.𝑆. (𝑟)]。而獲得基態電荷密度𝜌𝐺.𝑆. (𝑟)之後，我們便可以預測許多系統性質。以下是 Hohenberg-Kohn 定理的證明：如果有兩個相異 Hamiltonians 分別為 H1 和 H2，他們之間只有位能項不同，且形式如下: {. H1 = T + V𝑒𝑒 + 𝑉1𝑒𝑥𝑡 H2 = T + V𝑒𝑒 + 𝑉2𝑒𝑥𝑡. 其中V𝑒𝑒 為電子間的庫倫靜電力，𝑉1𝑒𝑥𝑡 、V2𝑒𝑥𝑡 為系統其他的位能，例如電子與原子核的庫倫作用力，或外加電場。若H1 與H2 的基態狀態函數相同，兩位能之間的差異必為 𝑉1𝑒𝑥𝑡 (𝑟) − 𝑉2𝑒𝑥𝑡 (𝑟) = 常數函數欲證此定理，我們假設有兩個相異的位能𝑉1 (r)和𝑉2 (𝑟)，兩位能差𝑉1 (r)-𝑉2 (𝑟) 不等於定值，且他們的基態電子密度分布相同。擁有這兩個相異位能的 Hamiltonians 分別為𝐻1 和𝐻2 ，其基態波函數分別為𝜓1 和𝜓2 ，即: 𝐻1 𝜓1 = 𝐸1 𝜓1 , 𝐻2 𝜓2 = 𝐸2 𝜓2 根據變分定理: E1 <⟨ψ1 |H1 |ψ1 ⟩=⟨ψ2 |H2 |ψ2 ⟩+⟨ψ2 |H1 -H2 |ψ2 ⟩=E2 + ∫ ρ(r){V1 (r)-V2 (r)}dr E2 <⟨ψ1 |H2 |ψ1 ⟩=⟨ψ1 |H1 |ψ1 ⟩+⟨ψ1 |H2 -H1 |ψ1 ⟩=E1 - ∫ ρ(r){V1 (r)-V2 (r)}dr 將兩方程式相加，會得到一個矛盾的結果 𝐸2 + 𝐸1 < 𝐸1 + 𝐸2 由此我們可以推論出𝑉1 (r)-𝑉2 (𝑟)必須為常數函數，也就是說，當我們擁有相同的基態電子密度時，以此為基態的 Hamiltonian 已經被決定了。而 Hamiltonian 的形式決定其基態電子密度，是當然之事，其反向論述亦成立，意即 Hamiltonian. 3.

(7) 與基態電子密度分布之間有一對一的對應，而系統的基態能量為電荷密度的泛函。根據 Hohenberg-Kohn 的理論一，當我們用一個臆測的電子密度𝜌′ 帶入已知 Hamiltonian 的泛函 : H = T + V𝑒𝑒 + 𝑉𝑒𝑥𝑡 ，可以得到一個基態波函數ψ′ 。利用變分原理來試驗這個波函數ψ′，會得到 Hohenberg-Kohn 的理論二: ⟨𝜓′|𝐻′|𝜓′⟩ = ∫ 𝜌′(𝑟)𝑉′(𝑟)𝑑𝑟 + 𝐹𝐻𝐾 [𝜌′] = 𝐸(𝜌′) ≥ 𝐸(𝜌0 ) 其中 FHK 泛函項的定義為: 𝐹𝐻𝐾 [𝜌] = 𝐸[𝜌] − ∫ 𝜌(𝑟)𝑉(𝑟) 其中，𝐻′和𝑉′與原來的𝐻和V之不同，皆因𝜌′與𝜌之不同而引起。. 2.1.2 The Kohn-Sham equations 一個多電子系統的基態總能泛函，由下面三項泛函所構成： E[ρ] = 𝑇𝑚 [𝜌] + 𝐸𝑒𝑒 [𝜌] + 𝐸𝑒𝑥𝑡 [𝜌] 其中𝑇𝑚 [𝜌]是動能、𝐸𝑒𝑒 [𝜌]是電子對電子交互作用能，這兩者密度泛函的形式是未知；而𝐸𝑒𝑥𝑡 [𝜌]則是外界對此粒子密度分佈的相互作用所產生之能量，其大小為 ∫ 𝑉𝑒𝑥𝑡 (𝑟)𝜌(𝑟)𝑑3 𝑟，其中𝑉𝑒𝑥𝑡 (𝑟)依個別問題而給定，此項泛函之形式為己知。原本在 Hohenberg-Kohn 的論述中，意欲放棄使用總波函數，轉而取用粒子密度作為自變數來使總能量的泛函簡化，但 Kohn-Sham 方法為了寫出較精確的動能，又重新引入單粒子波函數的寫法，並藉由說明粒子密度、總波函數及單粒子波函數的確切關聯，使得最後得到的是單粒子波函數的方程式。也因為有了單粒子波函數的寫法，針對𝑇𝑚 [𝜌]與𝐸𝑒𝑒 [𝜌]，Kohn 與 Sham 分離出了各粒子在獨自運動時的總動能𝑇𝑠 [𝜌]，以及電磁學的庫倫位能𝐸𝐻 [𝜌]，至於 𝑇𝑚 [𝜌]與𝐸𝑒𝑒 [𝜌]剩下的部分，則合併成為交換關聯能( Exchange correlation energy ) 𝐸𝑥𝑐 [𝜌] ： E[ρ] = 𝑇𝑠 [𝜌] + 𝐸𝑥𝑐 [𝜌] + 𝐸𝐻 [𝜌] + 𝐸𝑒𝑥𝑡 [𝜌] 除了 𝐸𝑒𝑒 [𝑟] 之外，各泛函項可以用以下公式表示，分別為: 獨自運動時的總動能𝑇𝑠 [𝑟]: 4.

(8) 𝑇𝑠 [𝜌] =. 1 ∑ 𝜓 ∗ (𝑟)∇2 𝜓𝑖 (𝑟)𝑑 3 𝑟 2 𝑖 𝑖. 其中 ρ(r) = ∑ |𝜓𝑖 (𝑟)|2 𝑖. 並滿足歸一化條件 N = ∫ 𝜌(𝑟)𝑑𝑟 = ∑ ∫ 𝜓𝑖∗ (𝑟)𝜓𝑖 (𝑟)𝑑 3 𝑟 𝑖. 𝜓𝑖 (𝑟)是單獨粒子的波函數，不是系統的總波函數 ψ(𝑟1 , 𝑟2 , … … , 𝑟𝑛 )，這些 𝜓𝑖 (𝑟) 被稱作是 Kohn-Sham 軌域。靜電位能 𝐸𝐻 [𝜌] 的部分則是： 𝐸𝐻 [𝜌] = ∫ 其中 𝑉𝐻 (𝑟) =. 𝜌(𝑟 ′ )𝜌(𝑟) 3 ′ 3 𝑑 𝑟 𝑑 𝑟 = ∫ 𝑉𝐻 (𝑟)𝜌(𝑟)𝑑3 𝑟 |𝑟 ′ − 𝑟|. 𝜌(𝑟 ′ )𝜌(𝑟) ，稱作 |𝑟 ′ −𝑟|. Hartree energy，即庫倫靜電作用。最簡單的則是外. 加場影響下的位能，其電荷密度的泛函可寫為 𝐸𝑒𝑥𝑡 [𝜌] = ∫ 𝑉𝑒𝑥𝑡 (𝑟)𝜌(𝑟)𝑑 3 𝑟 最後，即便不知道公式的情況下，我們一樣可以定義單一粒子的交換關聯能 𝜀𝑥𝑐 [𝜌(𝑟)]（exchange-correlation energy per particle），如下： 𝐸𝑥𝑐 [𝜌] ≡ ∫ 𝜀𝑥𝑐 [𝜌(𝑟)]𝜌(𝑟)𝑑 3 𝑟 在此𝜀𝑥𝑐 [𝜌(𝑟)]仍是未知的，Kohn-Sham 方法裏的第二大重點，就是針對交換相關效應從泛函近似成普通的函數，如下: 𝐿𝐷𝐴 𝜀𝑥𝑐 [𝜌(𝑟)]〈未知〉～𝜀𝑥𝑐 (𝜌(𝑟)) 〈近似〉. 我們原本需要知道整個𝜌(𝑟) 函數分佈才知道空間中各點的𝜀𝑥𝑐 大小，在此近 𝐿𝐷𝐴 似法中，只要給定位置𝑟0 ，代入𝜌(𝑟)得到該處的𝜌0 ，就能得到該處的𝜀𝑥𝑐 值。因. 此，𝜀𝑥𝑐 大小只跟局域位置的電荷密度大小有關，故被叫作局域密度近似（Local Density Approximation），簡稱 LDA。在實用上可分為兩種： local density approximation (LDA) 與 generalized 5.

(9) gradient approximation (GGA)，其分別在於 GGA 除了電荷密度大小還考慮電荷密度在空間中的梯度分佈，因此對於材料電性和磁性的預測與實驗結果更為吻合。. 2.1.3 Exchange-correlation energy 在 Kohn-Sham 方程的架構下，多電子系統基態問題可以在形式上轉化成有效的單電子問題。這種計算與 Hartree Folk 相似，然而我們只有在找出準確且便於表達的交換相關能(Exchange-correlation energy)形式才有實際的意義。交換相關能(𝐸𝑥𝑐 [𝜌])一般來說與密度函數 ρ(r)有關，因此交換相關能在意義上是非局域的，因此很難精確地描述。形式上交換相關能可以表示成: 𝐸𝑥𝑐 [𝜌] = ∫. 𝜌(𝑟 ′ )𝜌𝑥𝑐 (𝑟. 𝑟′) 3 ′ 3 𝑑 𝑟𝑑 𝑟 |𝑟 − 𝑟′|. 1. 在這裡， 𝜌𝑥𝑐 (𝑟, 𝑟 ′ ) = 𝜌(𝑟 ′ ) ∫0 𝑑𝜆[𝑔𝜆 (𝑟, 𝑟 ′ ) − 1]，而𝑔𝜆 (𝑟, 𝑟 ′ )是非均勻電子氣的函數，且𝜌𝑥𝑐 (𝑟, 𝑟 ′ )滿足 ∫ 𝑑𝑟′ 𝜌𝑥𝑐 (𝑟, 𝑟 ′ ) = −1 在計算時利用 Kohn-Sham 提出的局域密度近似(LDA)，利用密度為𝜌(𝑟)的均勻電子氣的交換相關能來得到非均勻電子氣的交換相關能𝐸𝑥𝑐 [𝜌]。對於一個變化曲線平坦的密度函數而言，用一均勻電子氣的交換相關能密度函數𝜀𝑥𝑐 [𝜌]代替非均勻電子氣可以表示成 𝐿𝐷𝐴 [𝜌] E𝑥𝑐 = ∫ 𝜀𝑥𝑐 [𝜌(𝑟)]𝜌(𝑟)𝑑 3 𝑟. 或利用 generalized gradient approximation (GGA)，均勻電子氣可以表示成 𝐺𝐺𝐴 [𝜌] E𝑥𝑐 = ∫ 𝜀𝑥𝑐 [𝜌(𝑟), ∇𝜌]𝜌(𝑟)𝑑 3 𝑟. 則 Kohn-Sham 方程中的交換相關位能可以近似成 𝑉𝑥𝑐 (𝑟, 𝜌) = 𝜀𝑥𝑐 (𝑟, 𝜌) + 𝜌. 𝑑𝜀𝑥𝑐 (𝜌) 𝑑𝜌. 將𝑉𝑥𝑐 表示成 r 的函數，並將交換相關能分成交換和關聯兩部分(𝜀𝑥𝑐 = 𝜀𝑥 + 𝜀𝑐 )可 6.

(10) 4. 以得出交換位能𝑉𝑥 (𝑟) = 3 𝜀𝑥 (𝑟)。. 2.2. GGA+U method 當我們計算過渡金屬氧化物的電子結構時，我們必須考慮軌域上電子間的庫. 倫作用力，但是 GGA 方法並未能正確的描述此現象，因此我們需要增加額外的函數來修正。 D 軌域及 f 軌域上的電子會對其他電子產生強大的庫倫作用力，這時系統總能量我們會使用 U(庫倫靜電位能)和 J(Hund’s rule exchange energy)來表示: 𝐸𝐺𝐺𝐴+𝑈 = 𝐸𝐺𝐺𝐴 + ∑. 𝑈−𝐽 𝑇𝑟(𝜌𝑎 − 𝜌𝑎 𝜌𝑎 ) 2. 此時的 ρa 表示有電子佔據的軌道矩陣(atomic orbital occupation matrix)。在過渡金屬氧化物中，若是我們沒有加入 U 的計算容易低估了物質的 energy gap 和 magnetic moment，因此適時的使用+U 的計算方法，有可能策略性地彌補 LDA 物理性質的問題，而目前我們利用+U 所採取的手段，屬於擬合參數的性質。. 7.

(11) Chapter 3. CuB2O4 的文獻探討. 3.1 CuB2O4 在基態時的鐵磁性研究 (Ferromagnetism in CuB2O4) CuB2O4 因為其特別的磁性結構而受到大家的關注，根據實驗發現 CuB2O4 在常溫下為 Commensurate 的磁性結構，然而將溫度降低到 T=10K 附近時 CuB2O4 會開始出現 incommensurate(每個單位晶格的分布具有 spin wave)的磁性結構，這樣有趣的材料性質，我們希望可以藉由第一原理計算模擬方法來更進一步分析、探討 CuB2O4 這種材料特性背後的原因。 CuB2O4 為一 tetragonal 結構(space group:I4̅2𝑑)，其中 a=11.565 Å、c=5.620 Å，每一個單位晶格(unit cell)中包含了 6 formula units，在材料中發現有兩種不一樣位置的 Cu 原子，一種記為 Cu(A) Wyckoff position 為 4b，位置在一個幾乎平行於 ab 平面之正方形的中心，而正方形的四個頂角上都有一個 O 原子存在；另外一種記為 Cu(B) Wyckoff position 為 8d，同樣位於以 O 原子為四個頂角的正方形中心，但 Cu(B)的正方形卻幾乎垂直於 ab 平面，如 FIG.1.。. FIG. 1. CuB2O4 分子結構，黑色為 Cu 原子，灰色為 O 原子，白色為 B 原子 (Boehm et al., 2003) 8.

(12) 根據由 R.V. Pisarevn 在穩定磁場下的諧振波光譜: FIG.2 及可見光吸收光譜 FIG.3 實驗，我們可以從中得知 CuB2O4 的材料特性。從 FIG.2 中，我們可以看到 CuB2O4 在能量大約為 1.41 eV 時會開始 transition，另外，在光譜 FIG.2(a)和 FIG.2(c)中我們看到 X 方向值和 Z 方向值的差異，有可能是因為兩種 Cu 具有不一樣的電子結構特性所導致的結果，而且在 FIG.2(c)裡面，我們可以將大圖(T=6K) 和小圖(T=30K)比對，可以發現 CuB2O4 會因為溫度的不同而呈現不一樣的諧振波光譜曲線。這是當實驗溫度低到一個程度時，材料出現了 incommensurate 的磁性結構導致光譜結果改變的現象。再藉由可見光吸收光譜的實驗結果: FIG.3，我們可以發現 CuB2O4 在 X 方向 absorption edge 出現在 1.403eV 的位置；Z 方向 absorption edge 則出現在 1.577eV 的位置，若將兩個方向的光譜曲線與 FIG.2 對照可以發現兩者非常的相近，因此這個實驗的結果對我們了解 CuB2O4 材料性質具有很大的幫助。藉由 R.V. Pisarevn 的實驗，我們可以知道 CuB2O4 在 X 和 Z 方向具有兩種非常不同的電子結構分佈，而且此分佈跟溫度有確切的關聯。而這樣的性質，我們希望能藉由第一原理計算來呈現出更詳細的電子分佈狀態，並跟實驗數據來做結果的對照。. 9.

(13) FIG. 2. (a)當一穩定磁場中沿著 x 軸作用時，(010)CuB2O4 的二次諧振波光譜 (MFISH: magnetic field induced second harmonic)。(b)的插圖中可以看到二次諧振波與磁場大小的關係，(c)是從(110)作用的光譜。 ( R.V. Pisarev et al., PhysRevLett.93.037204). FIG. 3. CuB2O4 在 H=0，k and E denote wave vector and polarization of the incoming light wave at frequency ω 時的吸收光譜(a)作用在(100)面上，(b)和(c)為作用在(010) 面上。( R.V. Pisarev et al., PhysRevLett.93.037204) 10.

(14) 3.2 CuB2O4 在 VASP 中的鐵磁性計算研究然而在第一原理的理論計算中，我們無法有效計算 incommensurate phase 的磁性結構，但是我們相信，藉由了解 commensurate phase 磁性結構的計算結果，我們還是可以在一定程度上了解這個材料的內涵，因此我們希望以 VASP GGA+U 的計算方式來跟 CuB2O4 的實驗結果來做比較，企圖了解 CuB2O4 的材料特性。在 Boehm 論文中，就利用了 CuB2O4 在 commensurate phase 的反鐵磁性結構材料計算來跟自己的實驗結果:TABLE I 做比對。Boehm 利用同樣的原子位置，如: TABLE II，來嘗試三種不一樣的磁性結構: FIG.4、FIG.5、FIG.6，如: TABLE ⅠII，利用這三種不同的磁性分佈來進行 CuB2O4 的材料計算，最後選擇𝜒 2 較小的結構 1 磁性結構: FIG.4，來進行後續的探討。而這組磁性結構也是我們在之後做第一原理計算時的參考設定值以及與 Boehm 的結果比較根據，然而請注意，這組磁性結構卻是破壞了原先原子結構的對稱性。我們的計算也是從此一結構出發進行結構(原子幾何位置以及磁矩向量)的優化，所以更進一步的對稱性破壞是可以預期的，我們也會注意這種破壞的程度。. TABLE Ⅰ.. CuB2O4 在 T=12 K 時的自旋常數(Boehm et al., 2003)，μＢ單位為 Bohr. magneton 11.

(15) TABLE ⅠI. CuB2O4 材料中，各原子的分佈 x、y、z 座標。B 為 isotropic Debye-Waller factor。. TABLE ⅠII. CuB2O4 在 commensurate phase 下，結構中較有可能的磁性常數. 12.

(16) FIG. 4. CuB2O4 在 commensurate phase 時的反鐵磁性結構，TABLE ⅠI.中的結構 1，其中 Wyckoff 位置為 4b 的 Cu(A)的磁性平行 ab 平面；Wyckoff 位置為 8d 的 Cu(B) 磁性平行 c 軸方向。. FIG. 5. CuB2O4 在 commensurate phase 時的反鐵磁性結構，TABLE ⅠI.中的結構 2，其中 Cu(A)的磁性平行 a 軸，Cu(B)磁性平行 c 軸方向。 13.

(17) FIG. 6. CuB2O4 在 commensurate phase 時的反鐵磁性結構，TABLE ⅠI.中的結構 3，其中 Cu(A)和 Cu(B)的磁性皆平行 ab 平面，且兩種 Cu 的夾角約為一直角。. 14.

(18) Chapter 4. CuB2O4 在 VASP 中的材料計算研究及討. 論 4.1 CuB2O4 在 VASP 中的自洽(self-consistent field method, SCF)材料計算為了更深入的了解 CuB2O4 的材料性質，我們利用 VASP 軟體來進行第一原理的材料計算。CuB2O4 為一 tetragonal 結構(space group:I4̅2𝑑)，在計算中的每一個單位晶格(unit cell)包含了 6 formula units，晶格的 a=11.357 Å ， K-points 的設定值為 4x4x4，使用的 DFT 為 PE(Perdew-Burke-Ernzerhof)，進行 non-collinear 的 GGA 方法計算，在 INCAR 中的 MAGMOM 的設定值依據 Boehm 文章中的 structure 1 來設定，數值如 TABLE ⅠI.。首先我們不使用+U 的方式進行 SCF 自洽計算，得到結構如 FIG.7，計算出來的 density of state(DOS)結果如 FIG.8。. FIG. 7. CuB2O4{6 formula, K-points (4,4,4), space group:I4̅2𝑑}的晶格結構以及原子的分布情形。. 15.

(19) FIG. 8. CuB2O4{6 formula, K-points (4,4,4), PE, nbands = 400}收斂後的 DOS 分布圖。. 4.2 CuB2O4 在 VASP 中使用 GGA+U 方法的計算 4.2.1 在 VASP 測試 U 值的大小從 FIG.8 中我們可以發現，在一般的 GGA 計算中，CuB2O4 材料的 band gap 不見了!這個結果與實驗比較可以發現是一個非常明顯的錯誤，而這個錯誤非常容易在計算鐵磁性材料時出現。在使用了 GGA+U 方法後，我們的計算結果出現 band gap，並且在 DOS 的導帶出現了兩個非常狹窄清晰的 peaks，我們相信這兩個 peaks 對於 CuB2O4 的材料特性有關，也是我們這次研究的重點目標。然而我們必須先決定我們要加入多大的 U，才有辦法在接下來的計算中繼續進行下去，因此我們進行了 U 的測試計算。在 CuB2O4 的材料中，我們知道它會 16.

(20) 有兩種不一樣特性的 Cu 原子，所以我們在加 U 的過程中，對兩種不一樣的 Cu 原子給了不同大小的 U，其中 Cu(A)的 U 我們表示為 U1，而 Cu(B)的 U 我們表示為 U2。為了可以找出較為理想的 U 值，我們利用 R.V. Pisarev 光譜實驗數據，如 FIG.3，來得到 CuB2O4 的 band gap 大小；以及根據 J.Y. Lin 的 O-K-edge 光譜實驗，如 FIG.9，來得知 CuB2O4 中兩個 peaks 的能量差(△E)，並將所有的測試結果與實驗值做比較，結果如 TABLE ⅠII.，我們可以發現在 VASP 的計算過程中，隨著我們慢慢的將 U 加大，band gap 與 peak 的△E 也會隨之增加，因此我們決定先以 peak 的△E 為主要依據，然後再從中挑出 band gap 比較小的那組 U 作為我們接下來計算所要使用的 U 值。因此最後我們選擇了(U1,U2) = (6.4,7.0)這組 U 值作為我們接下來計算的 U。另外在測試 U 值的同時，我們必需要對每一組 U 值的結果進行 DOS 分析，如 FIG.10，來確保加上 U 以後，每個軌域的 DOS 分佈是維持一致性的分布，沒有因為 U 值增加導致某些 DOS 的特徵位置(如:兩個 peaks)，出現了不一樣的軌域的貢獻。. 17.

(21) Peak2 Peak1. FIG.9. CuB2O4 O-K-edge 光譜 (Exp. graph by J.Y. Lin). U1_U2 exp Gap 1.403 Peak's delta E 0.600. 3_3.5 1.390 0.418. 4.5_5 1.642 0.488. 5_5.5 1.727 0.549. 5.5_5 1.798 0.366. 6.4_7 1.991 0.610. 7_6.5 2.104 0.366. 7.5_8 3.132 0.610. 8_7.5 2.346 0.366. TABLE ⅠII. CuB2O4 在加入(U1,U2)後的 gap 以及兩個 Peaks 的分布變化，我們可以發現隨著 U 的增加，band gap 會一起增加，但是當 U1 > U2 時，兩個 Peaks 的間距並沒有上升符合實際的數值；反觀當 U1 < U2 時，兩個 Peaks 的間距也會隨之慢慢增加。. 18.

(22) FIG.10. 在使用 GGA+U 計算之後，我們可以看到 CuB2O4 的 band gap 以及導帶上兩個 peak 的差值皆隨著 U(eV)的增加而慢慢變大，差值大小如 TABLE ⅠII.。. 4.2.2 CuB2O4 在 VASP 中的收斂(relaxation)計算在確定了 U1 及 U2 後，我們重新進行 CuB2O4 材料的收斂(relaxation)計算，得到 CuB2O4 的 DOS(FIG.11)以及 band structure 的能帶結構(FIG.12)，根據計算的結果，我們可以詳細分析在 DOS 中導帶上的兩個 peak 到底是由哪些原子的哪個軌域上的電子貢獻這些 state，藉由得知哪個軌域的貢獻，讓我們可以知道在導帶上這兩個 peaks 的電子雲分布方向；而 band structure 則可以讓我們清楚 CuB2O4 材料的能帶以及 gap 的分布狀態。. 19.

(23) FIG. 11. CuB2O4 的 density of state，可以看到在 fermi energy 右邊的導帶上出現了兩個尖銳的 peak，其中 246 表示為若 spin=1，到 fermi energy 總共有 246 個電子 (實際為 spin=1/2，故基態時總共有 492 個電子)。. FIG.12. 20.

(24) CuB2O4{6 formula, K-points (4,4,4), PE, nbands = 400, (U1,U2)=(6.4,7.0)}的 Band Structure 能帶結構分佈，右邊為 space group:I4̅2𝑑的 Brillouin zone 以及該結構上的特殊 K 點。在收斂的過程中，我們利用分析 LDOS 在價帶上的分布情形，來幫助我們觀察 CuB2O4 晶格是否在收斂時保持晶格結構的對稱性。由於 CuB2O4 為一 tetragonal 結構(space group:I4̅2𝑑)，所以我們可以發現 Cu 的 LDOS 分布在價帶中 dyz 以及 dxz 的曲線會完全重合一致，如 FIG.13；而 O 原子的 px 與 py 軌域在價帶的的 LDOS 分布也會完全相同，如 FIG.14，所以我們可以確定 CuB2O4 的結構對稱性沒有在收斂的過程中被破壞。在確認 CuB2O4 的結構收斂完畢後，我們得到 CuB2O4 晶體的邊長 a，以及高 c 的收斂值，並將其與 G. Petrakovskii 的實驗數據做比較得到 TABLE ⅠV；另外我們也將 CuB2O4 晶體中每一種原子單獨出來，如 TABLE V，並將每個原子的 x、y、z 座標與 Bohem 使用的值(TABLE II)做比較，得到 TABLE VI。除了 CuB2O4 的 x、y、z 軸的座標以外，我們也將每個 Cu 原子在收斂計算後的 spin component 從 OUTCRAR 中拿出來，與 Bohem 的結構一(TABLE III)做數值上的比較，如 TABLE VII。. 21.

(25) Cu(A). Cu(B). FIG. 13. 兩種 Cu 原子在價帶的 LDOS 分佈圖，其中因為結構對稱性，Cu(A)與 Cu(B) 的 dxz 和 dyz 兩條軌域 LDOS 分佈都是重疊的狀況。. FIG. 14. 22.

(26) 四種 O 原子在價帶上的 LDOS 分佈圖，可以看到每一種 O 原子的 px 和 py 兩條軌域 LDOS 分佈都是重疊的狀態。. A (Å). C (Å). GGA+U. 11.56513754. 5.66133755. EXP(JMMM1999). 11.484. 5.620. △/EXP. 0.707%. 0.736%. TABLE ⅠV. CuB2O4 的結構(space group:I4̅2𝑑)收斂後的邊長 a、c 與 G. Petrakovskii 的實驗數據比較，以及誤差百分比。. X. Y. Z. Cu(A). 0.000000. 0.000000. 0.500000. Cu(B). 0.081791. 0.250000. 0.125000. O1. 0.158580. 0.069219. 0.494946. O2. 0.253324. 0.250000. 0.625000. O3. 0.250000. 0.083231. 0.875000. O4. 0.074225. 0.187930. 0.798658. B1. 0.185315. 0.147711. 0.697342. B2. -0.000035. 0.250000. 0.625000. TABLE V CuB2O4 的結構(space group:I4̅2𝑑)收斂後每一種原子的 x、y、z 座標位置。. Cu(A). △x/x. △y/y. △z/z. 0.0. 0.0. 0.0. 23.

(27) Cu(B). 1.457%. 0.0. 0.0. O1. 1.654%. 4.877%. 0.613%. O2. 0.128%. 0.0. 0.0. O3. 0.0. 0.915%. 0.0. O4. 4.542%. 0.566%. 1.224%. B1. 1.265%. 1.171%. 0.664%. B2. 98.833%. 0.0. 0.0. TABLE VI CuB2O4 的結構(space group:I4̅2𝑑)收斂後每一種原子的 x、y、z 座標位置與 TABLE II 的誤差百分比，其中我們可以看到在 B2 的誤差特別顯著是因為 B2 本身的值非常小導致。. μx. μy. μz. μdft. μexp. Cu(A1). 0.495. 0.496. 0.129. 0.712518. 0.863. Cu(A2). -0.496. -0.495. -0.129. 0.712518. 0.863. Cu(B1). -0.074. -0.080. 0.735. 0.743035. 0.2. Cu(B2). -0.074. -0.081. -0.735. 0.743143. 0.2. Cu(B3). 0.081. 0.074. 0.735. 0.743143. 0.2. Cu(B4). 0.080. 0.074. -0.735. 0.743035. 0.2. TABLE VII CuB2O4 的結構(space group:I4̅2𝑑)收斂後每一個 Cu 原子的 spin component 與 Bohem 的實驗數值。(單位: Bohr magneton). 4.3 CuB2O4 材料在導帶上的特徵分析根據 FIG.11 我們可以看到在 DOS 導帶上有兩個很尖的 peaks，其中第一個 (較靠近左邊的)peak 是由 4 個電子所貢獻的 state，而第二個 peak 則是由 8 個電 24.

(28) 子所貢獻的 state。為了可以更進一步了解這兩個 peaks 的特性，我們分析每一種原子的 local density of states(LDOS)分布情形，結果我們發現第一個 peak 的主要組成貢獻是由 Cu(A)與 O(1)所提供，而第二個 peak 則是由 Cu(B)以及 O(2)、O(3)、 O(4)所提供的 state。我們現在將焦點放在這兩個 peaks 上，並將原子每個軌域的 DOS 貢獻都畫出來後發現，Cu(A)的 LDOS 在兩個 peak 中只有貢獻給第一個 peak，而且是 dxy 以及 dx2 這兩個軌域所提供的 state；而 Cu(B)的 LDOS 也只有貢獻給第二個 peak，但是這些 state 則是由 Cu(B)的 dyz、dxz、dx2 以及 dz2 這四個軌域一起貢獻的結果，如 FIG.15。更有趣的是，CuB2O4 材料中的 O 原子也呈現了兩種不一樣的 LDOS 分布關係，其中 O(1)在兩個 peaks 中幾乎將 states 都貢獻給了第一個 peak，而貢獻的軌域為 px 以及 py，並且在 LDOS 的分布圖中可以發現兩者的圖形呈現一個完全重疊的情況；O(2)、O(3)、O(4)則幾乎將 states 都貢獻在第二個 peak，不同於 O(1) 這三個 O 原子都提供了 px、py、pz 三條軌域的 states，但是 O(2)、O(3)的 px、py 都只提供在第二個 peak，而 O(4)則是將所有的 px、py、pz 軌域都貢獻在了第二個 peak，如 FIG.16。藉由 FIG.15 以及 FIG.16 我們可以分析第一個 peak 是 Cu(A)與 O(1)所貢獻的 states，而且 Cu(A)提供的是 dxy 以及 dx2 這兩個軌域；O(1)提供的是 px 以及 py 兩個軌域，因此我們可以知道第一個 peak 的電子雲分布應該是在 xy 平面上，而這個結果也與 J.Y. Lin 教授的 O-K-edge 光譜實驗(FIG.9)相符合。在第二個 peak 中，則是由 Cu(B)以及 O(2)、O(3)、O(4)所提供的 state，而且 Cu(B)提供的四個軌域可以再分成 dyz、dx2 及 dxz、dz2 兩組；而 O(2)、O(3)、O(4) 則提供 px、py、pz 三個軌域。因此第二個 peak 的電子雲分布在 dyz、dx2 這組會斜躺(接近平行)在 yz 平面上，在 dxz、dz2 這組會斜躺(接近平行)在 xz 平面上，，而這個結果也跟 FIG.9 的實驗結果相符。. 25.

(29) FIG. 15. 我們將 Cu 原子 LDOS 的 d 軌域投影出來，並跟 total DOS 做比較發現可以分出兩種類的 Cu 原子，其中左邊的 peak 主要都是 Cu(A)的 dxy 以及 dx2 這兩個軌域所貢獻得 states；而右邊的 peak 主要是 Cu(B)的 dyz、dxz、dx2 以及 dz2 四個軌域貢獻。. 26.

(30) FIG.16. 我們將 O 原子 LDOS 的 p 軌域投影出來，並跟 total DOS 做比較發現可以分出四種類的 O 原子，其中左邊的 peak 主要都是 O(1)的 px、py 兩個軌域所貢獻得 states；而右邊的 peak 主要是 O(2,3,4)氧原子的 px、py、pz 三個軌域貢獻。. 最後我們在 VASP 軟體中計算兩個 peaks 的 PARCHG(partial charge)，並將其利用 VESTA 軟體呈現出來。為了更便利的觀察，我們在這邊將原本一個 unit cell 有 6 個 formula 增加到 12 個 formula，其結構如 FIG.17，我們可以觀察到電子雲的分布情形，如 FIG.18，除了呈現出我們已知可以將原子歸類成 Cu(A) & O(1) 以及 Cu(B) & O(2、3、4)兩類的情形之外，在原子間的電子分布呈現出一群一群的分布情形，而 Cu(A) & O(1)每一群的電子雲分布方向平行於 xy 平面，如 FIG.19(a)；Cu(B) & O(2、3、4)每一群的電子雲分布方向接近平行於 z 軸方向，如 FIG.19(b)，這個結果跟我們從 J.Y. Lin 教授所作的實驗數據是非常符合的結果。 27.

(31) 這兩個 DOS 峰分別來自兩種 CuO plaquettes，而從這些 plaquettes 分散不連結的幾何分布，或直接觀察他們這兩個峰的 charge density 分布，皆可見這兩群能帶如此平坦(dispersionless)的原因。. FIG.17. CuB2O4{12 formula, K-points (4,4,4), PE, nbands = 400, (U1,U2)=(6.4,7.0)} 經過 VASP 收斂後的晶格結構以及原子的分布情形。. FIG. 18. CuB2O4 晶格中的局域電子密度分佈圖，左邊為 type A [Cu(A)&O(1)]在晶格中的分佈情形，可以對應到小圖中的藍線；而右邊為 type B [Cu(B)&O(2,3,4)]電 28.

(32) 子在晶格中的分佈情形，可以對應到小圖中的黑線，中間的小圖是 JY Lin 實驗室做出來的 O-K-edge 光譜。. FIG.19 CuB2O4 晶格中的局域電子密度分佈圖，左邊(a)為 type A 分佈主要是 ab 平面方向；而右邊(b)為 type B 可以看出 type B 的分佈幾乎都接近垂直 ab 平面。. Chapter 5. 結果與討論. 藉由 CuB2O4 材料的鐵磁性計算，我們了解到若要使用第一原理的模擬方法來驗證實驗結果時，我們必須要使用加入 Hubbard U 的方式來打開鐵磁性材料的 band gap，這時我們才有辦法得到與實際結果相符合的材料性質。在這一次的 CuB2O4 材料模擬計算中，我們使用了 VASP 5.4.1 版本的 GGA+U 方法來進行計算，我們藉由 R.V. Pisarev 的光譜實驗以及 J.Y Lin 教授的 O-K edge 實驗結果來對照我們計算數據，最後我們得到 U1=6.4 eV、U2=7.0eV 作為我們計算所使用的 U 值。我們使用 tetragonal 結構(space group:I4̅2𝑑)，每一個單位晶格(unit cell)包含 6 formula units；a=11.357 Å ；K-points = 4x4x4；U1=6.4 eV、U2=7.0eV 進行 noncollinear 的 GGA+U 方法進行收斂(relaxation)計算，得到收斂後的 CuB2O4 晶格 a=11.565 與 JMMM1999 論文中的數據 11.484 只有大約 0.7%的誤差值存在；若 29.

(33) 是將每種原子的位置來與 Bohem 的實驗數據來做比對時，除了 B 原子因為本身的 x 軸位置座標非常小(x= -0.000035)，導致誤差百分比極大之外，其他原子的誤差百分比皆在合理的範圍之內，因此我們可以初步的判定計算的結果非常符合實際的材料狀態。為了可以更進一步的判斷計算結果是否有破壞掉材料的對稱性，我們觀察計算後的 LDOS 分布情形，發現 Cu 原子的 dyz 與 dxz 兩個軌域以及 O 原子的 px 與 py 軌域的 LDOS 分布完全相同，因此我們可以確定收斂後的 CuB2O4 晶格對稱性並沒有被破壞掉。經過收斂後的 CuB2O4 的 DOS 分布在導帶上出現了兩個 peaks，這與 J.Y. Lin 教授所做的 O-K edge 實驗上黑線與藍線出現兩個 peaks 的現象是符合的。為了知道這兩個 peak 的特性我們對 Cu 原子和 O 原子進行了 LDOS 的分析，發現在第一個 peak(藍線)上主要是 Cu(A) & O(1)的貢獻，而在第二個 peak(黑線)上主要是 Cu(B) & O(2、3、4)所貢獻的 states。最後，我們利用 VASP 來計算出兩個 peaks 的局域電子密度分布(partial charge) 並使用 VESTA 軟體繪圖。我們可以很明顯的發現電子雲分布都是以 Cu 原子為中心與鄰近 O 原子所形成的區域分布，而且第一個 peak 的電子雲分布上主要為平行 ab 平面(即 xy 平面)方向；而第二個 peak 的電子雲分布主要為垂直 ab 平面。這個現象我們可以從 LDOS 分布中，第一個 peak 是由 Cu(A)的 dxy 以及 dx2 以及 O(1)的 px、py 軌域貢獻，因此電子雲的分布為平行 ab 平面(即 xy 平面)；而第二個 peak 主要是 Cu(B)的 dyz、dx2 及 dxz、dz2 與 O(2,3,4)的 px、py、pz 軌域貢獻，因此電子雲的分布為垂直 ab 平面。藉由這些計算結果與實驗比對之後，我們對於 CuB2O4 這樣的鐵磁性材料特性有了更多的認識可以在以後發展更多相關的應用。. 30.

(34) 參考資料: 1. Andreas Görling. Density-functional theory beyond the Hohenberg-Kohn theorem. Phys. Rev. A 59,3359 (1999) 2. A.I.Liechtenstein, V.I.Anisimov, and J.Zaanen. Density-functional theory and strong interactions: Orbital ordering in Mott-Hubbard insulators. Phys. Rev. B 52, R5467(R)(1995) 3. C.F. von WEIZSÄ CKER. Z. Phys. 96, 431 (1935) 4. John P. Perdew and Mel Levy. Physical Content of the Exact Kohn-Sham Orbital Energies: Band Gaps and Derivative Discontinuities. Phys. Rev. Lett. 51, 1884(1983) 5. Larry Spruch. Pedagogic notes on Thomas-Fermi theory (and on some improvements): atoms, stars, and the stability of bulk matter. Rev. Mod. Phys. 63, 151 (1991) 6. M. Boehm, B. Roessli, J. Schefer, A. S. Wills, B. Ouladdiaf, E. Lelièvre-Berna, U. Staub, and G. A. Petrakovskii. Complex magnetic ground state of CuB2O4. Phys. Rev. B 68, 024405(2003) 7. P.A.M. Dirac. Note on Exchange Phenomena in Thomas Atom. Proceedings of the Cambridge Philosophical Society, vol. 26, issue 03, p. 376 (1930) 8. R. V. Pisarev, I. Sänger, G. A. Petrakovskii, and M. Fiebig. Magnetic-Field Induced Second Harmonic Generation in CuB2O4. Phys. Rev. Lett. 93, 037204(2004) 9. Vladimir I. Anisimov, Jan Zaanen, and Ole K. Andersen. Band theory and Mott insulators: Hubbard U instead of Stoner. Phys. Rev. B 44, 943(1991) 10. V.I.Anisimov ,F.Aryasetiawan, and A.I.Liechtenstein. First-principles calculations of the electronic structure and spectra of strongly correlated systems: the LDA+ U method. J.Phys:Condens.Matter, Volume 9, Issue 4, pp. 767-808 (1997) 31.

(35)