李榮娟，孫民華，鄭大宇

(1.哈爾濱商業(yè)大學(xué)能源與建筑工程學(xué)院，哈爾濱150076;2.京都大學(xué)環(huán)境資源工程學(xué)院，京都日本6158540;3.哈爾濱師范大學(xué)物理與電子工程學(xué)院，哈爾濱150025)

玻璃轉(zhuǎn)變液體一直以來(lái)都是凝聚態(tài)物理研究的重點(diǎn)和難點(diǎn)問(wèn)題，多年來(lái)倍受國(guó)內(nèi)外很多學(xué)者的關(guān)注[1－6].近年來(lái)，人們對(duì)于玻璃轉(zhuǎn)變液體的結(jié)構(gòu)理解已經(jīng)有很大的進(jìn)步，但關(guān)于它的分子構(gòu)型還沒(méi)有形成統(tǒng)一定論[7－8].

最近幾年，分子動(dòng)力學(xué)(molecular dynamics簡(jiǎn)稱MD)模擬成為研究玻璃轉(zhuǎn)變的一個(gè)重要的理論工具.這是由于分子動(dòng)力學(xué)模擬不僅能像實(shí)驗(yàn)一樣進(jìn)行各種觀察，還可以描述微觀的分子運(yùn)動(dòng)細(xì)節(jié)，從而得到材料的結(jié)構(gòu)和性質(zhì).

勢(shì)能曲面(potential energy landscape簡(jiǎn)稱PEL)理論[9－10]是MD模擬中研究玻璃轉(zhuǎn)變行為的一個(gè)主要方法[11—15].這是因?yàn)?PEL方法能夠使玻璃轉(zhuǎn)變液體的分子位置在長(zhǎng)時(shí)間保持不變，分子就像被凍住一樣，從而減少了分析玻璃轉(zhuǎn)變液體結(jié)構(gòu)的難度.PEL是一個(gè)復(fù)雜的、多維的能量曲面，整個(gè)系統(tǒng)中勢(shì)能大小是根據(jù)系統(tǒng)中所有原子構(gòu)型的變化而定義的[16－17].參考文獻(xiàn)[3]給出了 PEL 的圖形解釋，從他們的圖形中我們能看到許多類似盆(basin)的形狀，每一個(gè)盆中都包含一個(gè)能量最小值(局部能量最小值)eIS－L，它對(duì)應(yīng)一個(gè)最小的分子構(gòu)型——內(nèi)在結(jié)構(gòu)(Inherent Structure簡(jiǎn)稱IS).研究IS的分子構(gòu)型變化對(duì)于理解玻璃轉(zhuǎn)變過(guò)程具有十分重要的意義.本文采用分子動(dòng)力學(xué)模擬的方法，利用PEL理論對(duì)玻璃形成液體進(jìn)行研究.同時(shí)，對(duì)一般液體(typical liquids－沒(méi)有經(jīng)過(guò)局域能量?jī)?yōu)化的液體)和IS液體中多面體的分子構(gòu)型隨溫度變化情況進(jìn)行分析.

1 勢(shì)函數(shù)和模擬方法

MD模擬玻璃轉(zhuǎn)變液體，所采用的勢(shì)函數(shù)是Lennard－Jones(L－J)勢(shì)，L－J勢(shì)函數(shù)的形式為:

方程(1)中的…α，β∈{A，B}，其中 表示第i(i=1，2，3，...，N)個(gè)粒子的位置坐標(biāo)指的是i原子與j原子之間的距離.參量ε和σ分別為L(zhǎng)－J勢(shì)壘的深度和L－J組元的尺寸.L－J勢(shì)的參數(shù):εAA=3.984 kJ/mol，σAA=0.341 nm，粒子的質(zhì)量為m，mA=39.948 amu，并且εAB/εAA=1.5，εBB/εAA=0.5，σAB/σAA=0.8，σAB/σAA=0.88，mB/mA=1，體積 V =(5.5σAA)3[18－20].

本文采用正則系綜(NVT)對(duì)204個(gè)A粒子和52個(gè)B粒子組成的二元液體系統(tǒng)進(jìn)行MD模擬.首先，將這些A和B粒子置于邊長(zhǎng)為L(zhǎng)的立方盒子中.其次，按照周期邊界條件，選取的交互截?cái)嗑嚯x rc=2.5σA[21]，時(shí)間步長(zhǎng)選取 0 .000 3 ps.研究的溫度范圍從600K降到4K，冷卻速率為1.3×1012K/s.為了獲得IS結(jié)構(gòu)，我們首先應(yīng)用標(biāo)準(zhǔn)的MD模擬的方法進(jìn)行高溫平衡;然后降溫平衡;最后經(jīng)過(guò)瞬時(shí)最陡下降(instantaneous steepest descent)的方法[22]使溫度迅速的降至2 K，模擬獲得IS能量.在每一個(gè)溫度T下，抽取了25個(gè)平衡構(gòu)型，在這25個(gè)平衡構(gòu)型中選取平均能量值eIS.

2 模擬計(jì)算結(jié)果及分析

2.1 結(jié)構(gòu)與動(dòng)力學(xué)分析方法

局部方向鍵體序參數(shù)(Local Bond－orientational Order)q6是用來(lái)描述原子團(tuán)簇的幾何形狀.其定義的形式如下:

其中:An中的n表示多面體的面，φf(shuō)n是多面體序的特征角，φn是鍵角分布函數(shù)N(θ)通過(guò)高斯擬合第一峰的位置，N(θ)是鍵角θ概率密度分布，θ是由一個(gè)中心原子和兩個(gè)近鄰原子組成的鍵角，hn是第一峰的高度，wn是高斯擬合第一峰的半高寬.

2.2 結(jié)果與討論

圖1顯示了玻璃形成液體IS中勢(shì)能最小eIS隨溫度的變化關(guān)系.從圖中我們很容易發(fā)現(xiàn)三個(gè)明顯的分區(qū)－平衡態(tài)(equilibrium liquid)，過(guò)冷態(tài)(supercooled liquid)和玻璃態(tài)(glass state).當(dāng)溫度高于390 K時(shí)，eIS幾乎不隨溫度而變化，此時(shí)體系處于平衡態(tài);當(dāng)溫度低于390 K時(shí)，eIS隨著溫度的降低迅速減小，液體系統(tǒng)的動(dòng)力學(xué)開(kāi)始減慢，此時(shí)體

圖1 玻璃形成液體IS中勢(shì)能最小值eIS與溫度T的關(guān)系

其中:Nb是團(tuán)簇原子之間的鍵數(shù)，Ylm(θi，i)是球諧函數(shù)，θi和φi分別是第i個(gè)鍵的極角和方向角.我們選取的團(tuán)簇一般是由13個(gè)最近鄰原子組成.選取方法:把一個(gè)原子作為中心原子，然后在它周圍找到與中心最近鄰的12個(gè)原子，這13個(gè)原子組成了一個(gè)團(tuán)簇.

M·Sun等人[23]提出了用來(lái)描述分子結(jié)構(gòu)堆砌的參數(shù)An，An的表達(dá)式如下:系跳出了平衡態(tài)進(jìn)入了過(guò)冷態(tài);當(dāng)溫度減小到36 K時(shí)，隨著溫度的繼續(xù)降低eIS值不再有明顯的變化，我們把能量開(kāi)始不再變化這一點(diǎn)定義為玻璃轉(zhuǎn)變點(diǎn)Tg，從Tg=36 K這點(diǎn)開(kāi)始體系進(jìn)入了玻璃態(tài).

我們采用局部方向鍵體序參數(shù)q6，對(duì)玻璃形成液體分子構(gòu)型結(jié)構(gòu)進(jìn)行統(tǒng)計(jì)分析.圖2描述二元玻璃形成液體系統(tǒng)IS中q6隨溫度的變化關(guān)系，此圖與參考文獻(xiàn)[23]中的表1對(duì)比，發(fā)現(xiàn)IS中存在二十面體，六面體和八面體等多面體的結(jié)構(gòu)轉(zhuǎn)變.但是在圖2中，很難發(fā)現(xiàn)玻璃形成液體系統(tǒng)的具體結(jié)構(gòu)趨向變化.為了更清晰的找到玻璃形成液體組成結(jié)構(gòu)趨向，采用An描述二元液體系統(tǒng)的IS結(jié)構(gòu)趨向.首先計(jì)算了IS中多面體序的鍵角分布函數(shù)N(θ)，對(duì)于一個(gè)13個(gè)原子的二十面體團(tuán)簇，理論上鍵角分布概率的峰值出現(xiàn)在63.4°，116.4°和180°.圖3描述了從24～380 K五個(gè)溫度下，玻璃形成體系IS中二十面體鍵角分布函數(shù)在各個(gè)角度(0°～180°)的百分比.從圖2中我們能夠觀測(cè)到N(θ)% 的最大峰值出現(xiàn)在 60°和 110°，而 180°的峰值幾乎觀測(cè)不到，這是因?yàn)?°和180°在計(jì)算中很難分清，所以我們沒(méi)有查找180°的鍵角.63.4°是二十面體中特殊存在的一個(gè)角度，因此我們的計(jì)算結(jié)果能夠說(shuō)明IS結(jié)構(gòu)中存在二十面體序.同時(shí)，我們得到的N(θ)%分布與D· Cicco[24]，和 M·Sun[23]所得到的結(jié)果相似.從圖2中我們還能夠觀測(cè)到，N(θ)%的第一峰的峰值隨溫度的減小而升高.這表明了當(dāng)液體系統(tǒng)被冷卻時(shí)，二十面體序?qū)⒃龆?

圖2 玻璃形成液體IS中的局部方向鍵序參數(shù)(q6)與溫度變化關(guān)系

圖3 IS中二十面體鍵角分布概率

圖4給出IS中六面體、八面體和二十面體結(jié)構(gòu)的鍵角分布鍵序參數(shù)An隨溫度的變化關(guān)系，溫度的變化范圍為380K～24K.從圖中可以觀察到六面體和八面體序隨著溫度的降低沒(méi)有明顯變化，而二十面體變化很明顯.二十面體隨著溫度的降低而增加，當(dāng)溫度接近玻璃轉(zhuǎn)變時(shí)，二十面體序增長(zhǎng)的趨勢(shì)幾乎消失.一般液體中與IS中二十面體相比較，前者比后者的少.根據(jù)鍵角分布鍵序參數(shù)An與溫度的變化關(guān)系我們能夠得出二十面體是形成玻璃轉(zhuǎn)變液體(非晶)的主要結(jié)構(gòu).玻璃轉(zhuǎn)變不僅是一般液體接近IS的過(guò)程，也是二十面體增加的過(guò)程.

圖4 玻璃形成液體中的六面體、八面體和二十面體結(jié)構(gòu)隨溫度的變化關(guān)系

3 結(jié)語(yǔ)

本文對(duì)IS中的二十面結(jié)構(gòu)體序進(jìn)行分析后得出以下結(jié)論:二十面體結(jié)構(gòu)隨溫度降低而穩(wěn)定增加，當(dāng)溫度下降到玻璃轉(zhuǎn)變溫度時(shí)，這種增加停止.玻璃轉(zhuǎn)變的IS比一般液體中的二十面體多.玻璃轉(zhuǎn)變過(guò)程是一般液體接近IS的過(guò)程，同時(shí)也是二十面體增加的過(guò)程.

[1] WEINTRAUB H，ASHBURNER M，GOODFELLOW P N，et al.Through the glass lightly[J].Science，1995，267(5204):1609－1618.

[2] EDIGER M D，ANGELL C A，NAGEL S R.Supercooled liquids and glasses[J].The Journal of Physical Chemistry，1996，100(31):13200－13212.

[3] DEBENEDETTI P G，STILLINGER F H.Supercooled liquids and the glass transition[J].Nature，2001，410:259 －267.

[4] MAYER P，BISSIG H，BERTHIER L，et al.Heterogeneous dynamics of coarsening systems [J].Physical Review Letters，2004，93(11):115701.

[5] BERTHIER L.Time and length scales in supercooled liquids[J].Physical Review E，2004，69:020201 －020204.

[6] TONINELLI C，WYART M，BERTHIER L，et al.Dynamical susceptibility of glass formers:contrasting the predictons of theoretical scenarios[J].Physical Review E，2005，71(4):041505－041525.

[7] FRECHERO M A，ALARCON L M，SCHULJZ E P，et al.Connections between structural jamming，local metabasin features，and relaxation dynamics in a supercooled glassy liquid [J].Physical Review E，2007，75(1):011502－011508.

[8] SPAEPEN F.Condensed－matter science:Five－fold symmetry in liquids[J].Nature(London)，2000，408:781－782.

[9] GOLDSTEIN M.Viscous Liquids and the Glass Transition:A Potential Energy Barrier Picuture[J].The Journal of Chemical Physics.1969，51(9):3728.

[10] STILLINGER F H，WEBER T A.Packing Structures and Transitions in Liquids and Solids[J].Science，1984，225:983 －989.

[11] BROOKS C L III，ONUCHIC J N，WALES D J.Statistical thermodynamics:Taking a walk on a landscape [J].Science，2001，293:612－613.

[12] SCIORTINO F，KOB W，TARTAGLIA P.Inherent structure entropy of supercooled liquids[J].Physical Review Letters，1999，83(16):3214－3217.

[13] BUCHNER S，HEUER A.Potential energy landscape of a model glass former:Thermodynamics，anharmonicities，and finite size effects[J].Physical Review E，1999，60(6):6507－6518.

[14] SASTRY S.The relationship between fragility，configurational entropy and the potential energy landscape of glass－forming liquids[J].Nature(London)，2001，409:164 －167.

[15] DOLIWA B，HEUER A.Hopping in a supercooled Lennard－Jones liquid:Metabasins，waiting time distribution，and diffusion[J].Physical Review E，2003，67(3):030501 －030504.

[16] STILLINGER F H，WEBER T A.Hidden structure in liquids[J].Physical Review A，1982，25:978 －989.

[17] RANDALL A，VIOLETTE L，STILLINGER F H.Multidimensional geometric aspects of the solid－liquid transition in simple substances[J].The Journal of Chemical Physics，1985，83(8):4079.

[18] KOB W.ANDERSEN H C.Testing mode－coupling theory for a supercooled binary Lennard－Jones mixture I:The van hove correlation function[J].Physical Review E，1995，51:4626 －4641.

[19] KOB W.ANDERSEN H C.Testing mode－coupling theory for a supercooled binary Lennard － Jones mixture.II.Intermediate scattering function and dynamic susceptibility[J].Physical Review E，1995，52:4134－4253.

[20] GLEIMN T，KOB W，BINDER K.How does the relaxation of a supercooled liquid depend on its microscopic dynamics? [J].Physical Review Letters，1998，81(20):4404 －4407.

[21] STODDARD S D，F(xiàn)ORD J.Numerical experiments on the stochastic behavior of a Lennard－Jones gas system[J].Physical Review A，1973，8(3):1504－1512.

[22] KOB W，DONATI C，PLIMPTON S J，et al.Dynamical heterogeneities in a supercooled Lennard－Jones liquid[J].Physical Review Letters，1997，79(15):2827 －2830.

[23] SUN M，LI C.A molecular dynamics simulation study of the polyhedral structure of liqud argon during glass transition[J].Modern physics Letters B，2009，23(8):1069－1075.

[24] DI C A，TRAPANANTI A，F(xiàn)AGGIONI S.Is there icosahedral ordering in liquid and undercooled metals[J].Physical Review Letters，2003，91(13):135505－135508.