Mg熔體凝固過程中的分子動(dòng)力學(xué)模擬*

孫 慧，堅(jiān)增運(yùn)，姜冰清，劉翠霞

(1.西安工業(yè)大學(xué) 材料與化工學(xué)院，西安 710021；2.西安工業(yè)大學(xué) 組織部，西安 710021)

diffusion coefficient

近年來，研究金屬熔體中原子團(tuán)簇與熔體結(jié)晶相在結(jié)構(gòu)上的相關(guān)性對(duì)認(rèn)識(shí)金屬的凝固規(guī)律、控制金屬凝固過程有著重要的意義[1-3].金屬熔體的結(jié)構(gòu)和性質(zhì)直接影響到凝固過程的形核和晶體生長(zhǎng)的特性[4-6]，但是由于金屬熔體中原子團(tuán)簇的結(jié)構(gòu)無法進(jìn)行實(shí)驗(yàn)測(cè)定，所以人們一直沒有找到金屬熔體中原子團(tuán)簇與熔體中結(jié)晶相在結(jié)構(gòu)上面的關(guān)系，目前的形核理論關(guān)于金屬熔體中臨界晶核、亞臨界晶核與熔體中結(jié)晶相在結(jié)構(gòu)上的關(guān)系僅僅建立在推測(cè)基礎(chǔ)上，文獻(xiàn)[4]曾指出研究非晶中亞臨界晶核的結(jié)構(gòu)及形成過程對(duì)完善和發(fā)展結(jié)晶動(dòng)力學(xué)具有重要意義，但是對(duì)一般金屬來說，目前實(shí)驗(yàn)所能獲得的冷卻速率遠(yuǎn)遠(yuǎn)達(dá)不到其發(fā)生玻璃化轉(zhuǎn)變所需的冷速，計(jì)算機(jī)模擬技術(shù)的飛速發(fā)展為上述問題的解決提供了可能，分子動(dòng)力學(xué)模擬是研究金屬結(jié)構(gòu)的一個(gè)有效方法，它能獲得比目前實(shí)驗(yàn)大得多的冷卻速率，可以從原子尺度直接模擬出熔體、非晶體及晶體的微觀結(jié)構(gòu)[7-10].

國(guó)內(nèi)外涌現(xiàn)出了大量的關(guān)于分子動(dòng)力學(xué)模擬金屬熔體凝固過程中微觀結(jié)構(gòu)演變規(guī)律的報(bào)道[11-13]，由文獻(xiàn)[14]采用分子動(dòng)力學(xué)對(duì)液態(tài)金屬的凝固過程進(jìn)行了研究，拉開了以分子動(dòng)力學(xué)模擬研究物質(zhì)微觀變化的序幕.文獻(xiàn)[15]采用分子動(dòng)力學(xué)方法對(duì)NaCl溶液進(jìn)行了模擬研究，發(fā)現(xiàn)自擴(kuò)散系數(shù)與溶液的溫度和濃度有關(guān)，溫度越高，自擴(kuò)散系數(shù)越大；文獻(xiàn)[16]采用EAM勢(shì)對(duì)凝固過程進(jìn)行了模擬研究，發(fā)現(xiàn)凝固以后的結(jié)構(gòu)和冷速有關(guān)，當(dāng)金屬熔體在較慢的冷速下冷卻，形成晶體，在較快的冷速下形成非晶，同時(shí)發(fā)現(xiàn)分析玻璃化轉(zhuǎn)變溫度對(duì)于分析晶體團(tuán)簇結(jié)構(gòu)有著重要的作用.文獻(xiàn)[17]采用分子動(dòng)力學(xué)方法對(duì)Pb-K合金的凝固過程進(jìn)行了模擬研究，通過與實(shí)驗(yàn)數(shù)據(jù)的對(duì)比發(fā)現(xiàn)模擬所得結(jié)果與理論結(jié)果相符合，證明了分子動(dòng)力學(xué)模擬結(jié)果的可靠性.其他一些研究者采用EAM勢(shì)以及相應(yīng)的微觀結(jié)構(gòu)分析方法對(duì)液態(tài)金屬Al凝固以后的組織進(jìn)行了研究.以鎂為基體的材料是目前工程應(yīng)用中最輕的金屬材料，具有較高的比強(qiáng)度和比剛度，故鎂金屬是輕質(zhì)金屬材料的理想基體.隨著各個(gè)學(xué)科的飛速發(fā)展，以鎂為基體的合金材料的優(yōu)勢(shì)逐漸體現(xiàn)出來，然而由于鎂金屬的化學(xué)活性較高，并且容易燃燒同時(shí)成孔性差，所以對(duì)于金屬鎂熔體的性能的實(shí)驗(yàn)研究有一定局限性.分子動(dòng)力學(xué)模擬方法不僅可以對(duì)實(shí)驗(yàn)難以達(dá)到和無法觀察的細(xì)節(jié)進(jìn)行研究，同時(shí)也能為傳統(tǒng)實(shí)驗(yàn)提供科學(xué)的預(yù)判.目前采用分子動(dòng)力學(xué)模擬的方法對(duì)鎂熔體的團(tuán)簇種類和團(tuán)簇?cái)?shù)量與冷速的關(guān)系以及亞臨界晶核、臨界晶核晶體在結(jié)構(gòu)上的關(guān)系也并沒有進(jìn)行深入研究，而搞清熔體結(jié)構(gòu)與結(jié)晶相結(jié)構(gòu)之間的關(guān)系是正確認(rèn)識(shí)形核生長(zhǎng)特性和規(guī)律的前提.

本文采用分子動(dòng)力學(xué)模擬方法研究了Mg熔體以不同冷卻速率凝固后的微觀結(jié)構(gòu)的演變規(guī)律，并通過理論計(jì)算確定出Mg熔體凝固后所得的晶態(tài)結(jié)構(gòu)為臨界晶核的臨界冷速，最后對(duì)不同初始溫度和不同保溫溫度下Mg熔體的擴(kuò)散系數(shù)進(jìn)行了研究.

1 模擬過程及分析方法

模擬軟件采用LAMMPS軟件，模擬對(duì)象為13 500個(gè)鎂原子，模擬采用三維周期性邊界條件，選擇這個(gè)邊界條件主要是能夠避免邊界效應(yīng)或表面效應(yīng)對(duì)研究模型力學(xué)作用的影響，適合本模擬的研究模型，為了使模擬結(jié)果準(zhǔn)確，故選擇時(shí)間步長(zhǎng)為2 fs，步長(zhǎng)過大會(huì)使部分原子溢出，影響模擬結(jié)果的準(zhǔn)確性.所選擇的勢(shì)函數(shù)為1990年ACKLAND G J等的勢(shì)函數(shù)，由于模擬系統(tǒng)可以和外部的時(shí)刻進(jìn)行熱交換，因此模擬系統(tǒng)的溫度能夠保持恒定，采用Nose-Hoover控溫控壓算法來控制系統(tǒng)的溫度壓力.在常壓下，將系統(tǒng)在980 K下運(yùn)行100 000步使其達(dá)到液態(tài)平衡；然后將Mg熔體在1 400 K下運(yùn)行100 000步使其轉(zhuǎn)化為液態(tài)并弛豫保溫500 000步，以不同冷速降溫至50 K并保溫弛豫，采用徑向分函數(shù)、H-A鍵型分析和晶體團(tuán)簇法分析Mg熔體的微觀結(jié)構(gòu).擴(kuò)散系數(shù)的模擬在1 400 K下運(yùn)行50 000步使Mg熔體熔化并達(dá)到平衡，然后以1.0×1011.0K·s-1的冷速冷卻至800 K，溫度每降低100 K，弛豫保溫5 000步測(cè)量均方位移.

2 模擬結(jié)果與討論

2.1 熔體熱歷史對(duì)Mg熔體均質(zhì)形核過冷度影響的分子動(dòng)力學(xué)模擬

2.1.1 熔點(diǎn)的模擬

采用1990年ACKLAND G J等的勢(shì)函數(shù)，以13 500個(gè)Mg原子作為模擬體系進(jìn)行熔點(diǎn)的測(cè)量，繪制溫度與時(shí)間的關(guān)系如圖1所示.從圖中可以看出，隨著時(shí)間的推移，整個(gè)體系中的溫度并不是一個(gè)定值，而是不斷變化的.但是這種變化又有一定的趨勢(shì)，總是在890 K附近波動(dòng)，金屬M(fèi)g理論熔點(diǎn)為921 K，模擬結(jié)果比實(shí)際值僅相差31 K，這在分子動(dòng)力學(xué)模擬中基本可以忽略不計(jì)，說明此勢(shì)函數(shù)下模擬結(jié)果可以達(dá)到穩(wěn)態(tài)，因此認(rèn)為金屬M(fèi)g的模擬熔點(diǎn)Tm=890 K.

圖1 采用ACKLAND G J等的勢(shì)函數(shù)模擬Mg的熔點(diǎn)

2.1.2 初始液態(tài)平衡性驗(yàn)證

模擬體系的初始狀態(tài)是否為液態(tài)平衡體系對(duì)于模擬結(jié)果的可靠性和準(zhǔn)確性有著至關(guān)重要的影響，因此在開始進(jìn)行一切模擬運(yùn)算前都要先進(jìn)行平衡性驗(yàn)證.圖2是不同原子數(shù)的液態(tài)鎂在980 K時(shí)的徑向分函數(shù).

圖2 不同原子數(shù)體系中液態(tài)鎂在980 K時(shí)的徑向分布函數(shù)

采用原子數(shù)分別為2 048、13 500和50 000個(gè)Mg原子為研究對(duì)象，將其在980 K(處于研究對(duì)象熔沸點(diǎn)之間的任意溫度值，為了方便后續(xù)實(shí)驗(yàn)的說明，這里取一個(gè)較低溫度值)的初始溫度下運(yùn)行50 000步使其轉(zhuǎn)化為液態(tài)并弛豫保溫50 000步，所得體系中不同原子數(shù)的徑向分布函數(shù)如圖2所示.從圖中可以看出，雖然體系中原子數(shù)不同，但是其徑向分?jǐn)?shù)的液態(tài)結(jié)構(gòu)是一致的.因此可以認(rèn)為Mg熔體經(jīng)過100 000步的保溫弛豫過程能夠使液態(tài)Mg達(dá)到平衡.

2.1.3 初始溫度對(duì)形核過冷度影響的模擬

將Mg熔體分別在七個(gè)不同的初始溫度(980 K，1 080 K，1 200 K，1 280 K，1 350 K，1 380 K，1 400 K)下運(yùn)行50 000步使其轉(zhuǎn)化為液態(tài)并弛豫保溫50 000步得到弛豫組織.再將弛豫組織以1.0×1011.0K·s-1的冷速快速冷卻到200 K，繪制能量隨溫度變化曲線如圖3所示.由圖可知，Mg熔體以1.0×1011.0K·s-1的冷速冷卻到200 K.一方面，初始溫度越高，升溫完成后，體系能量越高；另一方面，體系降溫過程中能量均隨著溫度的降低逐漸減小，當(dāng)溫度分別為619 K、614 K、606 K、590 K、584 K、584 K、584 K時(shí)，能量出現(xiàn)突變，急劇減小，使原本平滑的曲線出現(xiàn)明顯的拐點(diǎn)，說明當(dāng)溫度足夠高時(shí)，能量曲線出現(xiàn)轉(zhuǎn)折，故在此處系統(tǒng)開始凝固.根據(jù)熱力學(xué)理論，在恒溫恒壓的條件下，系統(tǒng)總是自發(fā)地趨向低能量的穩(wěn)定態(tài)，當(dāng)過熱溫度低于固液兩相的理論凝固溫度時(shí)，體系中固相的自由能總是低于液相的自由能，此時(shí)體系就會(huì)發(fā)生凝固結(jié)晶這說明金屬熔體在降溫的過程中發(fā)生了相變，形成了晶態(tài)凝固組織，此時(shí)拐點(diǎn)所對(duì)應(yīng)的溫度為Mg熔體凝固過程中的實(shí)際相變溫度.

圖3 不同初始溫度下Mg熔體凝固過程能量 隨溫度的變化曲線

為了研究過冷度隨初始溫度的變化關(guān)系，結(jié)合上述能量變化圖獲得初始溫度對(duì)應(yīng)的過冷度見表1.

表1 不同初始溫度的Mg熔體冷卻 到200 K時(shí)的過冷度

從表1可以看出，隨著過熱溫度的升高，實(shí)際相變溫度Tg大致上是逐漸減小的，其對(duì)應(yīng)的過冷度ΔT逐漸增大，當(dāng)溫度達(dá)到1 350 K以后，Tg和ΔT都趨于穩(wěn)定，不再變化，此時(shí)金屬的結(jié)晶溫度維持不變，相應(yīng)的過冷度也將確定不變.從圖4可以看出在過熱溫度較低時(shí)，體系內(nèi)存在一些晶體質(zhì)點(diǎn)，原子鍵對(duì)不易破壞，此時(shí)體系較為平衡，所以過冷度較小，但是隨著過熱溫度的升高，金屬原子的擴(kuò)散能力增強(qiáng)，體系的平衡狀態(tài)被逐漸打破，原子形核率增大，過冷度增大，曲線出現(xiàn)一次轉(zhuǎn)折；當(dāng)過熱溫度持續(xù)增大時(shí)，體系內(nèi)原子完全掙脫束縛，在體系內(nèi)任意運(yùn)動(dòng)，為了達(dá)到新的平衡狀態(tài)，故整個(gè)系統(tǒng)的自由能應(yīng)該在最小值附近，故而過冷度穩(wěn)定在306 K.

圖4 Mg熔體凝固過程中過冷度隨 初始溫度的變化曲線

2.1.4 冷速對(duì)均質(zhì)形核過冷度的影響的分子動(dòng)力學(xué)模擬及檢驗(yàn)

將Mg熔體在初始溫度為1 400 K下運(yùn)行50 000步使其轉(zhuǎn)化為液態(tài)并弛豫保溫50 000步得到弛豫組織.再將弛豫組織分別以十個(gè)不同冷速：1.0×1011.0K·s-1,，1.0×1011.5K·s-1，1.0×1011.8K·s-1， 1.0×1012.0K·s-1，1.0×1012.3K·s-1，1.0×1012.5K·s-1，1.0×1013.0K·s-1，1.0×1013.5K·s-1，1.0×1014.0K·s-1，1.0×1015.0K·s-1(選擇這些冷速主要是因?yàn)槟壳皩?shí)驗(yàn)條件下所能獲得的最大冷速為109K·s-1，在這個(gè)冷速條件下許多晶化能力強(qiáng)的金屬和合金不能實(shí)現(xiàn)玻璃化轉(zhuǎn)變，這限制了對(duì)非晶形成機(jī)理及結(jié)晶動(dòng)力學(xué)的研究，故將冷速區(qū)間細(xì)化可以更好的找出臨界冷速點(diǎn))冷卻到200 K(此溫度點(diǎn)可以在室溫下的任一溫度，為了計(jì)算效率更高更快，所以選擇溫度為200 K)并弛豫保溫后體系原子的勢(shì)能的變化曲線如圖5所示，在不同冷速下，隨著溫度的降低，體系勢(shì)能總體呈現(xiàn)降低的趨勢(shì).當(dāng)冷速降低到1.0×1012.5K·s-1后體系原子的勢(shì)能變化曲線上會(huì)出現(xiàn)拐點(diǎn)，這標(biāo)志著體系內(nèi)部發(fā)生相變，凝固組織為晶態(tài)；當(dāng)冷速大于1.0×1012.5K·s-1時(shí)，體系原子的勢(shì)能隨著溫度的降低平緩下降不出現(xiàn)突變點(diǎn)，這表明體系內(nèi)部沒有發(fā)生相變，凝固組織為非晶態(tài).因此，可以初步確定Mg熔體形成非晶的冷速應(yīng)大于1.0×1012.5K·s-1.

圖5 Mg熔體在不同冷速下原子勢(shì)能隨溫度的變化曲線

為了研究過冷度隨冷卻速度的變化關(guān)系，結(jié)合上述能量變化圖獲得不同冷速對(duì)應(yīng)的過冷度見表2.從表2可以看出，初始溫度為1 400 K時(shí)Mg熔體以不同冷速冷卻后，相變溫度Tg大致上是逐漸減小的，其對(duì)應(yīng)的過冷度ΔT逐漸增大.隨著冷速的增大，熔體凝固后的結(jié)晶溫度減小，說明過冷度受冷卻速率Rc的影響較大，當(dāng)冷速增加時(shí)，整個(gè)體系受到破環(huán)的溫度就越低，由于受冷卻速率的影響，產(chǎn)生了時(shí)間的滯后性，使得金屬的過冷度越大.

表2 1 400 K時(shí)Mg熔體在不同冷速下的過冷度

2.1.5 模擬所得過冷度結(jié)果的可靠性驗(yàn)證

采用文獻(xiàn)[18]中的理論模型計(jì)算Mg熔體凝固時(shí)的過冷度，當(dāng)熔體在降溫過程中發(fā)生均質(zhì)形核時(shí)，液相結(jié)構(gòu)中將會(huì)產(chǎn)生晶核，這一現(xiàn)象可用下面的公式來表示：

(1)

式中：NA為阿佛加德羅常數(shù)；r*為臨界晶核半徑；Vm為摩爾體積.

可通過下式計(jì)算出臨界晶核的原子半徑r*：

(2)

式中：σ為金屬的固/液界面能；ΔSf為金屬的熔化熵；ΔTv為凝固后的均質(zhì)形核過冷度.

查出固/液界面能σ，而均質(zhì)形核過冷度ΔTv可由下式確定：

(3)

其中Iv為

(4)

式中:Av為均質(zhì)形核率固定值(Av=1 043±1 m-3·s-1)；ΔGA為界面擴(kuò)散激活能；σ(T)為界面能；ΔHv為金屬結(jié)晶潛熱.如果固/液界面能知道，過冷度ΔTv就能計(jì)算出來.

對(duì)于公式中的參數(shù)值，可以通過查閱有關(guān)資料獲得，并將結(jié)果繪制成表3.模擬結(jié)果是否可靠，需要同其他理論或者方法作為對(duì)比.凝固時(shí)的最大均質(zhì)形核過冷度，需要將理論計(jì)算結(jié)果與模擬結(jié)果進(jìn)行對(duì)比.選1.0×1011.0K·s-1和1.0×1011.5K·s-1兩個(gè)冷速條件下的過冷度進(jìn)行比較，其中，各項(xiàng)數(shù)據(jù)繪制見表3.

表3 鎂在不同冷速下的均質(zhì)形核過冷度

圖6是模擬所得過冷度與理論對(duì)比的情況，通過固/液界面能理論模型獲得冷速變化時(shí)對(duì)應(yīng)的過冷度，在圖6中，實(shí)線表示由界面能理論所得過冷度與冷速之間的變化關(guān)系，藍(lán)色星狀符號(hào)表示模擬所得結(jié)果.從圖中對(duì)比可發(fā)現(xiàn)，模擬結(jié)果與理論吻合的很好，表明本次模擬的結(jié)果是準(zhǔn)確可靠的.

2.2 熔體熱歷史對(duì)Mg熔體凝固組織微觀結(jié)構(gòu)影響的分子動(dòng)力學(xué)模擬

2.2.1 徑向分布函數(shù)分析結(jié)果

為了進(jìn)一步確定Mg熔體凝固后形成非晶的冷速，采用徑向分布函數(shù)進(jìn)行分析.徑向分布函數(shù)曲線如圖7所示.

圖6 模擬結(jié)果與界面能理論結(jié)果的對(duì)比曲線

圖7 Mg熔體在1 400 K下以不同冷速降溫 至200 K后的徑向分布函數(shù)

從圖7可知，曲線10峰型圓鈍，線型平緩，具有液態(tài)的典型特征，曲線1峰型明顯而尖銳，線型波動(dòng)性大，說明熔體凝固后形成晶態(tài)結(jié)構(gòu).曲線2到曲線9，峰型逐漸變鈍，劈裂成兩個(gè)峰的第二峰逐漸消失，線型逐漸變得平緩，說明此時(shí)體系中既有非晶態(tài)的特征又具有晶態(tài)的特征.分析結(jié)果表明，冷速小于等于1.0×1011.0K·s-1時(shí)，熔體冷卻降溫后的結(jié)構(gòu)為晶態(tài)組織；冷速大于等于1.0×1012.5K·s-1時(shí)，降溫后的結(jié)構(gòu)為非晶態(tài)組織，冷速在1.0×1011.0K·s-1到1.0×1012.5K·s-1之間時(shí)，熔體降溫后的結(jié)構(gòu)為晶態(tài)和非晶態(tài)混合組織.

2.2.2 鍵型指數(shù)分析

Mg熔體以不同冷速降溫到200 K弛豫保溫后的鍵對(duì)分析結(jié)果如圖8所示，可以看出，當(dāng)冷速為1.0×1011.0K·s-1時(shí)，表征 fcc 結(jié)構(gòu)1421鍵型的分?jǐn)?shù)為0.5，隨著冷速的增加，1421鍵型的分?jǐn)?shù)先增加后降低，當(dāng)冷速超過1.0×1012.0K·s-1時(shí)降低的趨勢(shì)逐漸平緩；表征 hcp結(jié)構(gòu)1422 鍵型的分?jǐn)?shù)隨著冷速的增加總趨勢(shì)也是降低的，但沒有1 421鍵型降低得快，當(dāng)冷速超過1.0×1012.0K·s-1時(shí)降低的趨勢(shì)逐漸平緩.表征非晶結(jié)構(gòu)的1551，1541，1431 鍵型的分?jǐn)?shù)一開始為0，隨著冷速的增加才開始出現(xiàn)并升高，表明 Mg 熔體以 1.0×1011.0K·s-1冷速降溫凝固后完全形成晶體結(jié)構(gòu)；冷速大于1.0×1011.0K·s-1時(shí)，表征非晶結(jié)構(gòu) 1551，1541，1431 鍵對(duì)的分?jǐn)?shù)隨冷速的增加而增大，這說明隨著冷速的增大形成非晶的分?jǐn)?shù)增多.這與徑向分布函數(shù)分析所得到的結(jié)論基本相符合.

圖8 最終平衡構(gòu)型中各鍵對(duì)分?jǐn)?shù)隨冷速的變化

2.2.3 最大晶體團(tuán)簇中的原子數(shù)及其結(jié)構(gòu)

表4為Mg 熔體以不同冷速冷卻凝固后最大晶體團(tuán)簇中的原子數(shù).

表4 不同冷速下最大原子團(tuán)簇中的原子數(shù)Nc及 fcc原子數(shù)Nfcc，bcc原子數(shù)Nbcc，hcp原子數(shù)Nhcp

從表4中可以看出，當(dāng)冷速為1.0×1011.0K·s-1時(shí)，晶體團(tuán)簇由fcc，hcp結(jié)構(gòu)原子構(gòu)成，其中hcp結(jié)構(gòu)原子占絕大多數(shù)，fcc結(jié)構(gòu)原子只有很少一部分.隨著冷速的增大，最大原子團(tuán)簇的原子總數(shù)，以及fcc，hcp結(jié)構(gòu)原子數(shù)量都呈減少趨勢(shì).

為了能直觀真實(shí)的顯示最大晶體團(tuán)簇的微觀結(jié)構(gòu),采用周期性邊界條件對(duì)研究原子進(jìn)行了鏡像處理,結(jié)果如圖9所示.可以看出，Mg熔體以大于1.0×1011.0K·s-1冷速凝固后形成的最大晶體團(tuán)簇均是由 hcp 和 fcc 晶態(tài)結(jié)構(gòu)組成的混合偏聚體，其中絕大多數(shù)為hcp結(jié)構(gòu)，混有少量的fcc 結(jié)構(gòu)；當(dāng)冷速處于1.0×1011.0K·s-1到1.0×1012.3K·s-1之間時(shí)，凝固后形成的最大晶體團(tuán)簇中原子大小均勻排列緊密；以大于1.0×1012.3K·s-1的冷速冷卻凝固后的最大晶體團(tuán)簇，原子開始雜亂排布；隨著冷速的增加,最大團(tuán)簇的原子個(gè)數(shù)逐漸減少.

圖9 不同冷速下最終平衡構(gòu)型中最大晶體團(tuán)簇結(jié)構(gòu)(紫色表示hcp 晶態(tài)原子，黃色表示fcc 晶態(tài)原子)

2.3 Mg熔體擴(kuò)散系數(shù)的分子動(dòng)力學(xué)模擬

2.3.1 凝固過程中不同保溫溫度下的擴(kuò)散系數(shù)

使用LAMMPS軟件進(jìn)行分子動(dòng)力學(xué)模擬研究并不能直接得到我們想要的擴(kuò)散系數(shù)，直接得到的是一個(gè)與擴(kuò)散系數(shù)相關(guān)的量：均方位移MSD.根據(jù)菲克第一定律以及能斯特-愛因斯坦( Nernst-Einstein)方程，可以通過做MSD與時(shí)間的關(guān)系曲線的斜率來得出擴(kuò)散系數(shù).在1 400 K下對(duì)Mg熔體進(jìn)行加熱，運(yùn)行50 000步，使之熔化并達(dá)到平衡；然后以1.0×1011.0K·s-1的冷速冷卻至800 K，溫度每降低100 K，弛豫保溫500步，則可得到其均方位移.圖10是模擬所得不同保溫溫度下的均方位移與時(shí)間的關(guān)系曲線.

由圖10可以看出，一方面，在同一保溫溫度下，隨著保溫時(shí)間的延長(zhǎng)均方位移的數(shù)值會(huì)增大，但是增大的速度會(huì)減緩，表現(xiàn)在圖中就是均方位移曲線的斜率逐漸減小，即在同一保溫溫度下擴(kuò)散系數(shù)會(huì)隨著保溫時(shí)間的延長(zhǎng)而減??；另一方面，隨著保溫溫度的降低，熔體中原子的均方位移逐漸減小，曲線越來越平緩，所以原子的擴(kuò)散系數(shù)也就越來越小.

圖10 不同保溫時(shí)間下的均方位移

為了能直觀的反映擴(kuò)散系數(shù)與保溫溫度的變化規(guī)律，利用Einstein關(guān)系式，對(duì)不同溫度下的擴(kuò)散系數(shù)D0進(jìn)行統(tǒng)計(jì)見表5，從表5中可以看出隨著保溫溫度的降低，熔體原子的擴(kuò)散系數(shù)逐漸降低.

2.3.2 不同初始溫度下凝固過程中的擴(kuò)散系數(shù)

在上面已經(jīng)提到Mg熔體在980 K下就能完全熔化達(dá)到初始液態(tài)平衡結(jié)構(gòu).這里為了研究不同初始溫度對(duì)Mg熔體擴(kuò)散系數(shù)的影響，只需將熔體溫度分別從1 000 K、1 200 K和1 400 K、以1.0×1011.0K·s-1的冷速降溫到200 K，分別測(cè)定其凝固過程的MSD，從而進(jìn)一步分析凝固過程中擴(kuò)散系數(shù)的變化規(guī)律即可.圖11為不同初始溫度下的均方位移.

表5 不同保溫溫度下的均方位移及對(duì)應(yīng)的擴(kuò)散系數(shù)

圖11 不同初始溫度下凝固過程中的均方位移

從圖11中可以看出，初始溫度越高，開始降溫時(shí)的均方位移增加越快，隨著降溫過程的進(jìn)行，溫度越低均方位移的變化幅度減小，擴(kuò)散系數(shù)也隨之減小.當(dāng)時(shí)間足夠長(zhǎng)時(shí)，均方位移不再隨著變化，出現(xiàn)一個(gè)穩(wěn)定的值，此時(shí)，擴(kuò)散系數(shù)也不在隨著初始溫度的改變?cè)俑淖?

3 結(jié) 論

1) 當(dāng)冷速小于等于1.0×1011.0K·s-1時(shí)，Mg熔體凝固后形成晶態(tài)組織，具有fcc結(jié)構(gòu)；當(dāng)冷速大于等于1.0×1012.5K·s-1時(shí)，Mg熔體凝固后形成非晶態(tài)組織；當(dāng)冷速在1.0×1011.0K·s-1到1.0×1012.5K·s-1之間時(shí)，Mg熔體凝固后形成晶態(tài)結(jié)構(gòu)和非晶態(tài)結(jié)構(gòu)組成的混合組織，由fcc足夠長(zhǎng)時(shí)，均方位移不再變化，出現(xiàn)一個(gè)穩(wěn)定值，此時(shí)，擴(kuò)散系數(shù)不再隨著初始溫度的改變而改變.

2) 在一定的冷速范圍內(nèi)，Mg熔體凝固后得到的組織是晶體團(tuán)簇和非晶組成的混合體，最大晶體團(tuán)簇尺寸隨著冷速的增加而減小.

3) Mg熔體中晶態(tài)結(jié)構(gòu)原子團(tuán)簇尺寸小于臨界晶核尺寸時(shí)，用H-A鍵型指數(shù)法雖能確定出一定數(shù)量晶態(tài)結(jié)構(gòu)原子鍵對(duì)的存在，但徑向分布函數(shù)反映不出晶態(tài)結(jié)構(gòu)的特征.

4) 在同一保溫溫度下擴(kuò)散系數(shù)隨著保溫時(shí)間的延長(zhǎng)而減小，隨著保溫溫度的降低，熔體中原子的均方位移逐漸減小，原子的擴(kuò)散系數(shù)越來越??；不同初始溫度下隨著初始溫度的升高，開始降溫時(shí)的均方位移增加較快，隨著降溫過程的進(jìn)行，溫度越低均方位移的變化幅度減小，擴(kuò)散系數(shù)也隨之減小，當(dāng)時(shí)間位移不再隨著變化，出現(xiàn)一個(gè)穩(wěn)定的值，此時(shí)，擴(kuò)散系數(shù)也不在隨著初始溫度的改變?cè)俑淖?

參考文獻(xiàn)：

[1] LEE B S,BURR G W,SHELBY R M,et al.Observation of the Role of Subcritical Nuclei in Crystallization of a Glassy Solid[J].Journal of Science,2009,326(3):980.

[2] 堅(jiān)增運(yùn),相敏,朱滿,等.Al凝固特性隨冷卻速率變化規(guī)律的分子動(dòng)力學(xué)模擬[J].西安工業(yè)大學(xué)學(xué)報(bào),2015,4(35):311.

JIAN Zengyun,XIANG Min,ZHU Man,et al.Molecular Dynamic Simulation of the Dependence of Solidification Behaviour of Aluminum on the Cooling Rate[J].Journal of Xi’an Technological University,2015,4(35):311.(in Chinese)

[3] 張愛龍,劉讓蘇,梁佳,等.冷卻速率對(duì)液態(tài)凝固過程中微觀結(jié)構(gòu)演變影響的模擬研究[J].物理化學(xué)學(xué)報(bào),2005,21(4):347.

ZHANG Ailong,LIU Rangsu,LIANG Jia,et al.A Simulation Study for Effects of Cooling Rate on Evalution of Microstructures During Solidification Process of Liquid Metal Ni[J].Acta Physico-Chimica Sinica,2005,21(4):347.(in Chinese)

[4] GIBSON J M.Viewing the Seeds of Crystallization[J].Journal of Science,2009,326(5955):942.

[5] 堅(jiān)增運(yùn),鐘亞男,許軍鋒,等.冷速對(duì)Au熔體凝固組織影響的分子動(dòng)力學(xué)模擬[J].西安工業(yè)大學(xué)學(xué)報(bào),2016,36(6):468.

JIAN Zengyun,ZHONG Yanan,XU Junfeng,et al.Molecular Dynamic Simulation of the Influence of Cooling Rate on Solidification Structure of Au Melt[J].Journal of Xi’an Technological University,2016,36(6):468.(in Chinese).

[6] QI D W,WANG S.Icosahedral Order and Defects in Metallic Liquids and Glasses[J].Physical Review B,1991;44(2):884.

[7] LI Y D,HAO Q H,CAO Q L,et al.Two-order-parameter Description of Liquid Al Under Five Different Pressures[J].Physical Review B,2008,78(17):174202.

[8] LIU H R,LIU R S,ZHANG A L,et al.A Simulation Study of Microstructure Evolution During Solidification Process of Liquid Metal Ni[J].Chinese Physics B,2007,16(12):3747.

[9] LI J,LIU R,ZHOU Z,et al.Microstructural Transitions During the Rapid Cooling Proess of Liquid Metal Al[J].Journal of Materials Science & Technology,1998(14):461.

[10] WOLDE P R T,DAAN F.Numerical Calculation of the Rate of Crystal Nucleation in a Lennard-Jones System at Moderate Undercooling[J].Journal of Chemical Physics,1996,104(24):9932.

[11] MONDAL K,MUTRY B S.Prediction of Maximum Homogeneous Nucleation Temperature for Crystallization of Metallic Glasses[J].Journal of Non-Crystalline Solids,2006,352(50):5257.

[12] SCHENK T,HOLLANDMORITZ D,SIMONET V,et al.Icosahedral Short-Range Order in Deeply Undercooled Metallic Melts[J].Physical Review Letters,2002,89(7):46.

[13] QIN J,GU T,YANG L.Structural and Dynamical Properties of Fe78Si9B13 Alloy During Rapid Quenching by First Principles Molecular Dynamic Simulation[J].Journal of Non-Crystalline Solids,2009,355(48):2333.

[14] ALDER B J,WAINWRIGHT T E.Phase Transition for a Hard Sphere System[J].Journal of Chemical Physics,1957,27:1208.

[15] RAHMAN A.Correlations in the Motion of Atoms in Liquid Argon[J].Physics Review,1964,136(2):405.

[16] PARK S,SCHULTEN K.Calculating Potentials of Mean Force from Steered Molecular Dynamics Simulations[J].Journal of Chemical Physics,2004,120(13):5946.

[17] NOSE S.A Unified Formulation of the Constant Temperature Molecular Dynamics Methods[J].Journal of Chemical Physics,1984,81(1):5.

[18] JIAN Z,KURIBAYASHI K,JIE W Q.Solid-liquid Interface Energy of Metals at Melting Point and Undercooled State[J].Materials Transactions,2002,43(4):721.