分子動力學模擬研究混合堿土金屬對無堿硼鋁硅酸鹽玻璃結構和性能的影響

王照猛，朱雪梅，陳淑勇，李從云，馬立云

(1.中建材蚌埠玻璃工業設計研究院有限公司，浮法玻璃新技術國家重點實驗室,蚌埠 233010 2.硅基材料安徽省實驗室，蚌埠 233010)

0 引 言

無堿硼鋁硅酸鹽玻璃被廣泛應用于電子玻璃、光致變色玻璃、玻璃纖維和放射性廢物的固化等領域[1]。TFT-LCD基板玻璃屬于無堿硼鋁硅酸鹽玻璃[2]。同時TFT-LCD液晶顯示對超薄基板玻璃提出了高化學穩定性、高熱穩定性、無堿無砷、低密度以及高彈性模量等性能要求[3]。因此，近年來人們對無堿硼鋁硅酸鹽玻璃的結構和性能進行了深入研究，并努力探索性能與結構之間的關系[4-5]。

然而，宏觀性質是由微觀結構決定的，為了更好地理解宏觀性質，研究相應玻璃體系的組成、結構和性能之間的關系至關重要。分子動力學(Molecular Dynamics，MD)模擬是研究玻璃材料微觀結構信息的一種有效方法，可以直接觀察到原子結構，但很難通過實驗方法獲得[6]。Du等[7]采用分子動力學模擬研究了硼硅酸鈉和鋁硼硅酸鹽玻璃的組成、表面結構和性能之間的關系，得到理論模擬數據與實驗數據相吻合，同時其在研究鈉-鈣-鍶-鋁硼硅酸鹽玻璃中也得到了類似的結果[8]。

堿土金屬氧化物是玻璃形成外體中游離氧的主要提供者，起斷鍵和電荷補償作用，同時并對玻璃網絡形成體的配位數有很大的影響。因此，堿土金屬對無堿硼鋁硅酸鹽玻璃的物理化學性能的影響具有多面性。此外，不同堿土離子因其離子場強度不同而改變了鋁硼硅酸鹽玻璃的結構。Wu等[9]研究發現，具有高場強的堿土金屬離子更容易促進非橋氧形成，同時，高電場強度的陽離子也會增加鋁硅酸鹽玻璃中五配位鋁的形成。

因此，本研究的目的是通過MD模擬研究混合Ca/Mg對無堿硼鋁硅酸鹽玻璃的結構和性能的影響。具體計算了陽離子的局部環境和中程結構信息，并討論了混合Ca/Mg對彈性模量和剪切模量的影響。

1 模擬過程

為了研究混合堿土金屬對無堿硼鋁硅酸鹽玻璃結構和力學性能的影響，B2O3-Al2O3-SiO2的摩爾分數恒定，改變Ca/Mg比設計6組無堿硼鋁硅酸鹽玻璃樣品。組分和密度如表1所示。

表1 玻璃的化學組成和密度Table 1 Compositions and densities of simulated glasses

采用LAMMPS軟件進行分子動力學模擬。在MD模擬中，經驗勢對于生成真實的結構和描述玻璃性能至關重要。本研究使用的是部分電荷對勢，結合了長程庫倫相互作用和短程相互作用的Buckingham勢。該方法在多組分玻璃的建模中具有良好的可移植性和計算效率。Buckingham勢Φij(rij)[10]的計算公式如下：

(1)

式中：rij為原子i和j之間的距離；qi、qj為原子i和j的電荷；Aij、ρij和Cij為勢參數[11](如表2)。

表2 Buckingham勢參數Table 2 Parameters of the Buckingham potential

所有的模擬都使用了由3 210個原子組成的立方體模擬盒。在立方體模擬盒中，通過實驗密度，如表1所示隨機生成初始原子位置坐標。初始結構在0 K下進行能量最小化，然后通過正則過程(原子數不變、體積不變、溫度不變，NVT系綜)在300 K下弛豫60 ps。之后，通過NVT系綜，系統在6 000 K熔化100 ps(目的是快速打破模型的初始結構，使原子充分混合均勻，使模型快速熔化。不考慮工業生產實際溫度是因為在工業溫度下很難打破模型的初始結構，模型熔化慢，增加模擬耗時)，在5 000 K下弛豫100 ps，然后以5 K/ps的冷卻速率降溫到300 K，冷卻時長940 ps。在300 K下，玻璃在等溫等壓過程(原子數不變、壓強不變、溫度不變，NPT系綜)下平衡20 ps，然后在NVT系綜下平衡10 ps。最后，在微正則系綜(原子數不變、體積不變、能量不變，NVE系綜)下平衡10 ps，消除內部應力，獲得玻璃的最終結構。

2 結果與討論

2.1 徑向分布函數

徑向分布函數是描述液體和非晶體材料微觀結構最常用的數學語言，是結構無序化程度的表征，反映的是從任意指定A原子為中心，半徑為r到r+dr的球殼內發現B原子的數目為n(B)，公式如下：

(2)

式中：gαβ為徑向分布函數，α、β代表兩個原子；r為距中心原子半徑;δr為球殼厚度;n(r)為球殼內粒子數;ρβ為理想晶體中B原子的密度；V為體積。

圖1為O-Si、O-Al和O-B的偏徑向分布函數，從圖中可以得出，O-Si的鍵長約為1.63 ?(1 ?=0.1 nm)，O-Al鍵長約為1.75 ?，O-B鍵長約為1.39 ?，與文獻值相符[12-14]，說明模擬選擇的勢函數是可靠的。混合堿土金屬Ca/Mg對O-Si的鍵長變化幾乎沒有影響，對O-Al的鍵長有較小影響。與O-Si、O-Al偏徑向分布函數不同的是，O-B偏徑向分布函數有雙峰，分別在1.39 ?和1.47 ?處，說明硼配位數存在三配位和四配位。1.39 ?處峰高于1.47 ?處峰，說明硼的配位數主要以三配位構型存在。從偏徑向分布函數的峰寬以及尖銳程度可以得出，Al-O結構單元比Si-O和O-B結構單元松散、穩定性差。

圖1 O-Si、O-Al和O-B的偏徑向分布函數Fig.1 Partial radial distribution functions of O-Si, O-Al and O-B

2.2 鍵角分布函數(BAD)

圖2為O-Si-O、O-Al-O和O-B-O的鍵角分布函數，O-Si-O鍵角分布在108.4°～108.6°范圍內，說明在硼鋁硅系玻璃中，[SiO4]四面體呈現為規則的四面體結構，且混合堿土金屬的變化對[SiO4]四面體構型幾乎沒有影響。與O-Si-O鍵角分布相比，O-Al-O鍵角分布在106.4°～107.3°范圍內，說明在硼鋁硅系玻璃中[AlO4]四面體呈現不規則的四面體結構，且混合堿土金屬的變化導致了[AlO4]四面體結構變形。與O-Si-O和O-Al-O鍵角分布函數相比，O-B-O鍵角分布函數分別在109.5°和120.0°附近呈現出單峰，說明硼配位數存在[BO3]平面結構和[BO4]四面體結構，與O-B偏徑向分布函數的結果相一致。隨著Ca/Mg比的降低，109.5°處峰高度降低同時120.0°處峰高度增加，說明[BO3]含量在增加，[BO4]含量在降低。O-Al-O的鍵角分布范圍比O-Si-O分布范圍寬，說明[AlO4]四面體結構的變形大，穩定性差，與偏徑向分布函數結果一致。在鋁硅酸鹽玻璃中[AlO4]四面體結構周圍存在堿金屬和堿土金屬離子來補償電荷，導致[AlO4]四面體結構易變形，與模擬的結果相一致。

圖2 O-Si-O、 O-Al-O和O-B-O的鍵角分布函數Fig.2 Bond angle distribution functions of O-Si-O, O-Al-O and O-B-O

2.3 配位數

圖3(a)為O-Si配位數的百分含量圖，從圖中可以看出，O-Si的配位數主要為四配位，存在少量的五配位，沒有發現三配位的存在。Stebbins[15]在1991年通過NMR發現在石英玻璃中存在Si的五配位結構。Ha和Garofalini[6]通過分子動力學方法研究表明，在硼鋁硅酸鹽玻璃中存在少量的Si的五配位結構。因此，在無堿硼鋁硅酸鹽玻璃中存在少量的Si的五配位結構。Horbach和Kob[16]在包含8 016個原子的SiO2玻璃系統中使用BKS勢，發現在相對較高的2 750 K溫度下只有99%的四配位硅。Hoang[17]在3 000個原子的SiO2玻璃系統中使用Morse勢，觀察到在低于700 K的溫度下，多達2.7%的Si原子處于五配位。Lee等[18-19]發現在模擬SiO2玻璃系統中，缺陷濃度高達24.7%，這表明在模擬結構中不完全四配位的系統是一個不可避免的現象。Vollmayr等[20]通過對粒子A和B的任意二元混合徑向分布分析表明，隨著冷卻速率的降低，系統的局部順序變得更加明顯，表明配位狀態與冷卻速率有關。

圖3(b)為O-Al配位數的百分含量圖，可以得出，O-Al的配位數以四配位主，存在少量的三配位和五配位。根據MD模擬結果表明，六倍配濃度很低，在統計誤差范圍內可以忽略不計，與以往MD模擬結果一致。隨著Ca/Mg比的降低，五配位含量先降低后增加，在Ca/Mg=6/9時達到最小值。說明增加Mg2+的濃度對五配位鋁的含量的影響不是呈線性關系，而是滿足一個最值。Stebbins[21]在2005年通過NMR發現在硼鋁硅酸鹽玻璃中存在高配位的鋁。Ha和Garofalini[6]通過分子動力學的方法研究表明，在硼鋁硅酸鹽玻璃中存在少量的高配位鋁，得出模擬與實驗的結構相一致。值得注意的是，由于三配位沒有通過實驗方法檢測到，因此將計算得到的三配位作為模擬方法生成了鋁結構缺陷。

圖3(c)為O-B配位數的百分含量圖，從圖中可以得出，O-B的配位數以三配位主，存在一定量的四配位。隨著Ca/Mg比的降低，三配位含量增加，四配位含量降低。這由于Mg2+比Ca2+電場強度大，給出游離氧的能力弱于Ca2+，導致玻璃組分中，隨著Mg2+的增加，三配位含量增加，四配位含量降低。Stebbins[22]在2010年通過控制溫度和退火速率研究硼鋁硅酸鹽玻璃中硼的配位，研究表明BO3/BO4比與玻璃中的網絡修飾離子和樣品中的B2O3含量有關。

圖3 O-Si、O-Al和O-B的不同配位數含量Fig.3 Content of different CNs of O-Si,O-Al and O-B

2.4 氧種類

在玻璃體系中存在四種類型的氧：橋氧BO(氧原子與兩個網絡中間體相連)、非橋氧NBO(氧原子與一個網絡中間體和外體相連)、自由氧FO(氧原子與兩個網絡外體相連)和TBO(氧原子與三個網絡中間體相連)。如圖4為硼鋁硅酸鹽玻璃體系中不同氧類型的百分含量圖。隨著Ca/Mg比值下降，BO呈下降的趨勢，NBO呈增加的趨勢。說明增加Mg2+的濃度對NBO的含量呈線性關系，由于Mg2+比Ca2+電場強度大，給出游離氧的能力弱于Ca2+，導致玻璃組分中，隨著Mg2+的增加，NBO呈增加的趨勢。與Wu等[9]研究發現具有高場強的堿土金屬離子更容易促進非橋氧形成的結果相一致。改變Ca/Mg比會引起BO、NBO和TBO轉變。模擬結果證實了FO原子的存在，其含量與BO、NBO和TBO相比非常低，可以忽略不計。在玻璃系統中，TBO已經在早期的核磁共振實驗中被發現。例如，在無定形的Al2O3-2SiO2中發現了大量的TBO，占總的氧原子30%[23]。在硅酸鹽玻璃中堿金屬和堿土金屬離子如Na+、Ca2+、Mg2+等會填充在網絡周圍，形成網絡補償體，破壞玻璃網絡結構，增加非橋接氧數目。非橋氧進一步增加玻璃體系的熱容和熵來改變離子的擴散和流動性，最終降低玻璃體系的粘度。因此，非氧數目是影響硅酸鹽玻璃物理化學性質的重要因素。

圖4 玻璃中氧類型及百分含量Fig.4 Oxygen speciation and content of glass

2.5 力學分析

在力學性質分析中，計算了硼鋁硅酸鹽玻璃的體積模量、剪切模量和楊氏模量。圖5為不同組分樣品的模量變化。可以得出，在低Mg成分中，隨著Mg2+濃度的增加，楊氏模量和剪切模量減小。在高Mg2+成分中，隨著Mg2+濃度的增加，楊氏模量和剪切模量增加。楊氏模量都大于65 GPa，滿足玻璃基板在運輸、操作和加工過程中的下陷要求。體模量無明顯規律變化。

圖5 玻璃的體模量、剪切模量和楊氏模量Fig.5 Bulk, shear and Young’s modulus of glass

3 結 論

(1)Ca/Mg比對O-Si和O-B的鍵長幾乎沒有影響；對O-Al的鍵長有較小影響。

(2)無堿硼鋁硅酸鹽玻璃體系中[SiO4]四面體呈現規則的四面體結構，[AlO4]四面體呈現不規則的四面體結構。

(3)O-Si和O-Al的配位數以四配位主，O-Al存在少量的高配位；O-B配位數以三配位主要，存在一定量的四配位。

(4)隨著Ca/Mg比的降低，[BO3]含量增加，[BO4]含量降低；橋氧(BO)含量降低，非橋氧(NBO)含量增加。

(5)楊氏模量和剪切模量出現先降低后增加的趨勢，在Ca/Mg=6/9時，楊氏模量達到最小值65 GPa。

綜合本文模擬分析結果表明，目前的分子動力學方法可以產生合理的硼鋁硅酸鹽結構，其微觀結構和力學性能與文獻報道的相一致。因此，本文研究結果對實驗設計無堿硼鋁硅酸鹽玻璃配方有一定的指導作用。