王智聰, 孫永麗, 幸曉鳳

(太原理工大學(xué) 物理與光電工程學(xué)院, 太原030024)

1 引 言

當(dāng)冷卻速度足夠快, 液體可避免結(jié)晶而在過冷至低于某一溫度(玻璃轉(zhuǎn)變溫度)時, 內(nèi)部原子大規(guī)模的運(yùn)動在時間尺度上被凍結(jié), 從而形成非晶結(jié)構(gòu)[1,2].玻璃是一種典型的非晶態(tài)材料, 人類發(fā)現(xiàn)玻璃到應(yīng)用玻璃已經(jīng)有近千年, 但其轉(zhuǎn)變機(jī)制仍然是凝聚態(tài)物理的難題之一[3-5].許多關(guān)于玻璃轉(zhuǎn)變的理論模型, 如自由體積模型[6]、 熱力學(xué)統(tǒng)計模型[7]、 能量勢壘理論[8]、 固體模型理論和模態(tài)耦合理論[9,10]等等, 并不能完整而明確地解釋玻璃轉(zhuǎn)變現(xiàn)象.玻璃態(tài)的原子結(jié)構(gòu)和其動力學(xué)性質(zhì)是密不可分的, 玻璃態(tài)形成過程中, 其原子結(jié)構(gòu)的變化必然會引起其動力學(xué)性質(zhì)的顯著變化.目前的相關(guān)研究[11-14]多是關(guān)于對擴(kuò)散系數(shù)D、 非高斯參數(shù)α2(t)、 粘度η等動力學(xué)性質(zhì)隨時間變化趨勢的研究, 也有壓強(qiáng)對非晶體結(jié)構(gòu)和動力學(xué)性能的研究[15,16], 原子尺寸與動力學(xué)性質(zhì)的關(guān)系并未過多涉及.本文研究了原子尺寸對LJ液體玻璃轉(zhuǎn)變過程中動力學(xué)性質(zhì)的影響.

2 模擬研究的條件和方法

本文采用Lammps經(jīng)典分子動力學(xué)軟件[17]進(jìn)行模擬.初始狀態(tài)為8781個A原子和2195個B原子隨機(jī)排列而成的立方原胞, 約化單位下A, B的質(zhì)量均為1.采用周期性邊界條件和等溫等壓(NPT)系綜, 使用LJ勢[18], 勢函數(shù)為

(1)

其中α,β∈{A, B}, 相互作用勢參數(shù)σ和ε分別為長度和能量單位.對于氬, 其對應(yīng)的長度為3.41?, 能量為120KkB.設(shè)置勢函數(shù)參數(shù)如下:εAA= 1.0,εBB= 0.5,εAB= 1.5,σBB= 1.0,σAB=(σAA+σBB) / 2, 這與文獻(xiàn)[19, 20]中所研究的參數(shù)ε相同.設(shè)置σAA分別為0.9、 0.86、 0.82和0.78進(jìn)行模擬.初始溫度為T*= 0.9, 對應(yīng)于T= 109 K, 時間步長為0.01 ps, 具體模擬過程為: 初始結(jié)構(gòu)在溫度為0.9時弛豫足夠的時間確保體系達(dá)到平衡, 之后采用梯度降溫法以1.88×1012K/s的冷卻速率每隔0.05降溫至0.3, 每降溫一次平衡一次并輸出相關(guān)數(shù)據(jù)以供后續(xù)分析使用.

3 模擬計算結(jié)果與分析

3.1 玻璃轉(zhuǎn)變溫度Tg分析

冷卻過程中, 原子尺寸不同的體系平均原子體積隨溫度的變化曲線的斜率是連續(xù)的, 并沒有出現(xiàn)突變, 如圖1, 這是由于降溫速率較快, 原子重組速度較慢,無法在冷卻條件限定的時間內(nèi)完成充分重組.此時, 液態(tài)的結(jié)構(gòu)好像被“凍住” ,從而發(fā)生玻璃轉(zhuǎn)變形成非晶結(jié)構(gòu).玻璃轉(zhuǎn)變溫度Tg是判斷液體非晶形成能力的重要依據(jù)[21].通過體積判別法[22,23]能夠得到原子尺寸為0.9、 0.86、 0.82、 0.78時所對應(yīng)的Tg的值分別為0.6037、 0.5904、 0.5808、 0.5772, 在所研究尺寸范圍內(nèi)Tg隨原子尺寸的減小而減小.

圖1 不同原子尺寸下LJ體系平均原子體積隨溫度的變化關(guān)系Fig.1 Relationships between average atomic volume and temperature of LJ system under different atomic sizes

圖2 (a) σAA= 0.9時不同溫度下MSD隨時間的變化; (b)T*= 0.45時不同原子尺寸下體系的MSD隨時間的變化Fig.2 (a) The MSD of LJ system at different temperatures when σAA= 0.9; (b) The MSD of LJ system with atomic sizes of 0.78 - 0.9 when T*= 0.45

3.2 均方位移MSD分析

均方位移MSD是用來描述體系中的原子在時間t內(nèi)平均位移的平方, 通過計算MSD來研究原子尺寸對體系動力學(xué)性質(zhì)的影響.圖2(a)是尺寸為0.9時不同溫度下MSD隨時間的變化.在高溫階段, MSD對時間呈線性關(guān)系, 此時體系為液態(tài).隨著溫度的不斷下降, MSD逐漸出現(xiàn)平臺區(qū)且越來越明顯, 說明原子運(yùn)動受到的束縛越來越大.t< 0.2 ps, MSD隨時間的增加很快, 原子受到的束縛較小做自由振動.0.2 10 ps, MSD隨時間增加較快, 原子脫離籠子的束縛做長程擴(kuò)散運(yùn)動, MSD斜率隨溫度的下降而不斷減小, 粒子運(yùn)動越來越慢, 動力學(xué)不均勻性越來越強(qiáng).圖2(b)為溫度0.45時原子尺寸不同的體系的MSD隨時間的變化, 在自由振動階段, MSD幾乎重合, 說明原子尺寸對做自由振動的原子影響微乎其微; 在長程擴(kuò)散階段,尺寸越小, 原子擴(kuò)散受到的限制越大; 在β弛豫階段, 原子尺寸大的體系MSD高于尺寸小的體系, 說明尺寸越小, 原子之間的相互作用力越強(qiáng), 動力學(xué)不均勻性越強(qiáng).

3.3 擴(kuò)散系數(shù)D分析

通過對MSD曲線求斜率可以得到不同原子尺寸下體系的擴(kuò)散系數(shù)D, 圖3為lnD隨1/T*的變化關(guān)系, 可以得到原子尺寸越小擴(kuò)散系數(shù)D越小.在高溫階段,不同尺寸下體系lnD與1/T*的關(guān)系均可以用Arrhenius方程[25]表示, 即

圖3 不同原子尺寸下體系的擴(kuò)散系數(shù)D隨溫度的變化及擬合曲線Fig.3 The variation of diffusion coefficient D with temperature under different atomic sizes and fitting curves

(2)

其中,D0為指前因子,E為擴(kuò)散激活能,KB為玻爾茲曼常數(shù),T*為溫度.通過擬合得到不同尺寸下的D0分別為6.07、 6.37、 6.91、 6.93;E分別3.19 KB、3.29 KB、 3.41 KB、 3.46 KB, 即原子尺寸越小, 躍遷過程中克服周圍原子對其阻礙所需要的擴(kuò)散激活能越大.當(dāng)溫度低于0.6附近時, 擴(kuò)散規(guī)律逐漸偏離高溫時的規(guī)律, 這表明在降溫過程中, 整個體系的微觀結(jié)構(gòu)隨溫度發(fā)生了改變, 使得其動力學(xué)性質(zhì)發(fā)生不連續(xù)變化.

3.4 非高斯參數(shù)α2(t)

非高斯參數(shù)α2(t)是用來表征動力學(xué)不均勻性的重要方法, 其峰值αmax表示動力學(xué)不均勻程度.圖4為四種原子尺寸不同的體系分別對應(yīng)的不同溫度下非高斯參數(shù)α2(t)隨時間的變化關(guān)系, 由此得出四種體系中, 當(dāng)t< 0.2 ps, α2(t) = 0, 原子運(yùn)動在振動階段符合高斯分布; 當(dāng)0.2 τβ, 體系進(jìn)入α弛豫階段, 溫度較高時的非高斯參數(shù)α2(t)不斷減小直到趨于零, 表明原子擴(kuò)散動態(tài)特性滿足高斯分布, 溫度較低時的非高斯參數(shù)α2(t)未減小到零, 這是由于低溫時體系的弛豫時間不夠充分.同時, 經(jīng)過擬合四種體系β弛豫中不同溫度下非高斯參數(shù)α2(t)隨時間的變化趨勢, 得到四種體系的非高斯參數(shù)α2(t)與時間t的關(guān)系均滿足冪律函數(shù)[26]:

α2(t)=α0×(t-t0)γ

(3)

其中指前因子α0分別為0.155、 0.177、 0.204、 0.212,β弛豫開始時刻t0= 0.2,指數(shù)γ分別為0.556、 0.577、 0.559、 0.609.圖5是四種體系不同溫度下的非高斯參數(shù)峰值αmax統(tǒng)計圖, 隨著溫度下降, 粒子運(yùn)動變慢, 動力學(xué)不均勻性增強(qiáng), 低溫時原子尺寸越小, 原子之間的引力越強(qiáng), 體系的動力學(xué)不均勻性越強(qiáng).

3.5 粘度η分析

液體粘度是決定其玻璃形成能力的關(guān)鍵動力學(xué)參數(shù), 因此計算體系粘度對研究玻璃轉(zhuǎn)變過程中動力學(xué)性質(zhì)的變化至關(guān)重要.

根據(jù)Green-Kubo公式[27], 即

(4)

圖4 不同原子尺寸下不同溫度的非高斯參數(shù)α2(t)隨時間的變化及擬合曲線Fig.4 The variation of non-Gaussian parameter α2(t) with time at different temperatures under different atomic sizes and fitting curves

圖5 不同原子尺寸下非高斯參數(shù)峰值αmax隨溫度的變化Fig.5 The peak values of non-Gaussian parameter αmax at different temperatures under different atomic sizes

其中,V是系統(tǒng)體積,KB是玻爾茲曼常數(shù),Pxy是壓力張量的偏對角線分量, 〈…〉是求系綜平均.壓力張量表達(dá)式為

(5)

η=η0exp[D*T0/(T-T0)]

(6)

其中,η0是溫度趨于無窮大時的粘度值,T0是粘度趨于無窮大時的溫度值, D*代表動力學(xué)強(qiáng)度, 是衡量液體脆性的指標(biāo).四種不同體系的粘度值均隨溫度的下降而增大并在各自的玻璃轉(zhuǎn)變點(Tg)附近發(fā)生突變[29], 尺寸越小的體系的粘度值越大且突變越明顯，其動力學(xué)不均勻性越強(qiáng).

圖6 不同原子尺寸下粘度η隨溫度的變化及擬合曲線Fig.6 The variation of viscosity η with temperature under different atomic sizes and fitting curves

4 結(jié) 論

本文基于Lammps分子動力學(xué)軟件, 研究了原子尺寸對LJ液體玻璃轉(zhuǎn)變過程中動力學(xué)性質(zhì)的影響, 得到如下結(jié)論:

(1)在研究范圍內(nèi), 四種原子尺寸不同的體系快速冷卻凝固均形成非晶體, 且玻璃轉(zhuǎn)變溫度Tg和體系的擴(kuò)散系數(shù)D均隨著原子尺寸的減小而減小.研究均方位移MSD得到在β弛豫階段原子尺寸越小, 體系動力學(xué)不均勻性越強(qiáng), 體系溫度越低, MSD斜率越小, 即整個降溫過程中, 體系動力學(xué)性質(zhì)發(fā)生了不連續(xù)變化.

(2)β弛豫階段非高斯參數(shù)α2(t)隨時間的變化關(guān)系符合冪律函數(shù), 與原子尺寸大小無關(guān).非高斯參數(shù)峰值隨溫度的降低以及原子尺寸的減小而不斷增加, 并隨著溫度的降低向長時間區(qū)域移動, 證實了玻璃轉(zhuǎn)變過程中動力學(xué)不均勻性的存在以及不均勻性程度與溫度和原子尺寸大小有關(guān).

(3)對體系粘度的計算發(fā)現(xiàn)LJ液體粘度與溫度和體系內(nèi)原子尺寸大小緊密相關(guān), 粘度隨著溫度的降低及粒子尺寸的減小而增大, 與非高斯參數(shù)峰值隨溫度及原子尺寸大小的變化趨勢相同, 即體系的動力學(xué)不均勻性隨著原子尺寸的減小而增大.