Cu66Ti34非晶合金凝固過程的分子動力學模擬

邱克強, 厲 虹, 孫 晶, 尤俊華, 任英磊, 李慶豐

(沈陽工業大學 材料科學與工程學院，沈陽 110870)

塊體Cu基非晶合金自從被報道以來[1]，因其具有優異的力學性能和相對較低的材料成本而受到關注。目前通過成分設計可以獲得多種 Cu基塊體非晶合金[2?6]，其中Cu50Zr50二元合金可以形成直徑為2 mm的非晶合金[7?9]?？紤]Ti具有較低的成本和與Zr類似的性質，因此，Cu-Ti二元非晶合金已成為研究熱點。但是，目前的實驗結果表明，Cu-Ti二元合金還不具有制備成塊體非晶合金的能力[10?11]。

目前，基于分子動力學(MD)方法對Cu-Zr二元系非晶合金原子團簇結構演化的模擬比較多[12?13]，而從凝固過程出發，對Cu-Ti二元非晶合金性能的模擬還很少，本文作者基于 WADLEY等[14]和ZHOU等[15]提出的普適嵌入原子模型(GEAM 勢)，從動力學角度模擬過冷液態合金在不同溫度下原子的均方位移(MSD)、偶關聯函數和焓隨溫度的變化關系，進一步給出定壓比熱容與溫度的函數關系式，從而判定玻璃轉變的發生。

1 計算方法

1.1 勢函數模型

GEAM勢的表達形式如下：

式中：Fni、Fi、Fe和n為模型參數，可由結合能和體模量計算得到；ρe為平衡時的電子密度。

兩體作用勢可寫成如下形式：

式中：r為原子間的距離；re為平衡時最近鄰原子間的距離；A、B、α和β為模型參數；κ、λ為有關截斷半徑的附加參數。當計算合金時，合金勢采用以下形式[16]：

式中：φab(r)為合金的兩體勢；φaa(r)和φbb(r)分別為a和b組分的兩體勢；fa(r)和 fb(r)分別為a和b組分的電子密度函數。Cu66Ti34二元合金的勢參數如表1所列。

1.2 模擬方法

模擬過程中采用二元合金Cu66Ti34，將2 000個原子，按照比例置于立方體盒子中，并讓系統在三維周期邊界條件下以4×1013K/s冷卻速率進行冷卻，時間步長選定為2.5×1015s。模擬過程中采用NTP系綜，首先讓系綜在1 600 K等溫運行80 000步以使系綜達到平衡態，然后以給定冷卻速率逐步冷卻到300 K，最后在300 K弛豫40 000步，其中每100 K讓系統等溫運行2 000步，用于系統的數據統計，從而算出體系的各種性質參數。

表1 Cu66Ti34 二元合金的勢參數Table 1 Potential parameters of binary Cu66Ti34 alloy

2 計算結果及討論

2.1 偶關聯函數

偶關聯函數被廣泛用來描述液態和非晶態的結構特征，是一個重要的形態譜參數，通常定義如下[17]：

式中：gi,j(r)表示以a原子為中心、在距離為r到r+Δr的球殼范圍內發現原子的概率；L為模擬元胞的邊長；Ni和Nj分別是原子i和j的個數；nαj是指在 i原子周圍r到r+Δr范圍內發現j原子的個數；Δr是計算的步長。

圖1所示為Cu66Ti34在4×1013K/s冷卻速率下的偏偶關聯函數和總偶關聯函數的圖像。從圖1可以看出，高溫時偶關聯函數都顯示出典型的液態結構，而在低溫300 K時，gCu-Cu(r)、gTi-Ti(r)、gCu-Ti(r)和gtot(r)第二峰均發生分裂，標志著非晶結構的形成。由于總的偏偶關聯函數 gtot是所有原子偏偶關聯函數三維信息的平均[18]，因此，總偶關聯函數的第二峰分裂不如偏偶關聯函數明顯。與高溫液態結構的 g(r)相比較，低溫時第一峰的寬度變窄、高度增加、變銳，這表明隨著溫度的降低，原子配位數增加，短程有序性增強；另一方面，非晶結構的偶關聯函數在第二峰過后和液態結構的無太大區別，偶關聯函數都是比較平緩的曲線，說明仍然是長程無序的。從3種偏偶關聯函數來看，gCu-Ti(r)的第一峰都高于 gCu-Cu(r)和 gTi-Ti(r)的，這說明合金中異類原子之間的相互作用較強，容易形成非晶[19]。在600 K時，總偶關聯函數就開始出現比較微弱的第二峰分裂，隨著溫度的降低，此現象表現更加明顯。

2.2 均方位移(MSD)

圖1 Cu66Ti34二元合金不同溫度下偶關聯函數曲線Fig.1 Pair correction function curves of Cu66Ti34 binary alloy at different temperatures: (a)gCu-Cu; (b)gTi-Ti; (c)gCu-Ti; (d)gtot

粒子位移平方的平均值稱為均方位移，均方位移隨時間的變化表征了液態金屬粒子的擴散行為。為了便于觀察Cu66Ti34合金中Cu和Ti原子在急冷到各個溫度下MSD的圖像變化，以900 K為分界，將圖像分開展示，圖2(a)和(c)所示分別為Cu66Ti34合金中Cu和 Ti原子在急冷到 1 500、1 400、1 300、1 200、1 100、1 000和900 K時MSD(圖中以D表示)的變化；圖2(b)和(d)所示分別為Cu和Ti原子急冷到800、700、600、500、400和300 K時MSD隨時間的變化情況。在800 K之前，它們都表現出液態的特點，尤其是在1 200 K之前，MSD曲線的斜率較大，說明原子具有較強的擴散能力。從圖2(b)和(d)可以看出，800 K時，MSD曲線都變得比較復雜，從 MSD的數值看到，Cu和 Ti在800 K時MSD的最大值均已小于1×10?2nm2，合金熔體已經很粘稠。600 K時，在0.25 ps后，曲線的斜率降低，此時可認為擴散受阻，同時，也是動力學明顯呈現凍結狀態的起始溫度。到500 K時，原子擴散能力很弱，但局域結構仍可以表現出隨時間的波動性。到300 K時，原子幾乎沒有擴散，說明原子已經找到它們的平衡位置。通過圖2(a)與(c)，圖2(b)與(d)中 MSD的數值比較，還可以看出 Cu的擴散能力比Ti的強。

為了探討原子的擴散性與玻璃轉變溫度之間的關系，圖3給出了偶分布函數第一谷的最小值與第一峰最大值之比R[20?21]的分布規律，即R=g(r)min/g(r)max。在冷卻速率為4×1013K/s條件下，測得玻璃化轉變溫度Tg約為600 K，略低于Cu50Ti50的玻璃化轉變溫度641 K[22]，考慮到Cu50Ti50具有較高的固相溫度，因此預測 Cu66Ti34非晶合金玻璃化轉變溫度在600～641 K是符合實際的。

2.3 熱力學性質分析

在冷卻速率為4×1013K/s條件下，Cu66Ti34合金熔體凝固過程中焓隨溫度的變化如圖4所示，總體上，合金的生成焓隨溫度的降低而降低。這說明液態金屬過冷溫度越低，合金中所含熱量越少，越有利于非晶的形成。模擬得到的數據在400～900 K之間有波動，說明在這個溫度范圍內，此體系應處于過冷態，結構變化比較復雜，玻璃化轉變在這個溫度范圍內發生。能量變化較為連續，只是在一定溫度以下曲線斜率減小，這是非晶形成的標志[23]。當凝固過程中沒有發生結晶時，這一過程表現在能量、體積曲線上， 也是一條連續變化的曲線，但在Tg處應有一個微小的斜率變化[24?25]。由于曲線具有連續性，因此，可以對數據進行整體擬合，得到如下關系式：

圖2 在不同溫度下Cu66Ti34合金中Cu原子與Ti原子的均方位移曲線Fig.2 MSD curves of Cu and Ti atoms at different temperatures: (a)MSD for Cu at 900?1 500 K; (b)MSD for Cu at 300?800 K;(c)MSD for Ti at 900?1 500 K; (d)MSD for Ti at 300?800 K

圖3 Cu66Ti34合金偶關聯函數第一谷的最小值與第一峰最大值之比Fig.3 Ratios of minimum to maximum of the first peak of PCF for Cu66Ti34 alloy

圖4 Cu66Ti34合金凝固過程中焓隨溫度變化的關系Fig.4 Enthalpy of Cu66Ti34 vs temperature during solidification process

金屬液體處于Fragile和Strong兩種極限之間，比較接近Fragile極限，這個系統的特征是在Tg點時比熱容有大的突變。計算定壓比熱容有兩種方法[26]：一是利用能量的漲落來計算，二是用熱力學求導的方法計算，即cp=dh/dT。由于第一種方法需要的原子數目較多，故采用第二種方法計算。即

當T的取值范圍為300～1 600 K時，可得到圖5所示曲線。從圖5可以看出，函數在800～900 K之間存在最大值。這條曲線雖不能看出 cp的階躍性變化[25]，但可以估計玻璃化轉變的發生。BAILEY等[27]認為函數的最大值點即是 Tg點。對 cp求極值，得到Tg=892 K，這樣求出的Tg值比較具體。這與用偶分布函數第一谷的最小值與第一峰最大值之比求得的 Tg點相比，處于溫度更高的位置。一方面，BAILEY[27]所考慮的是極值點，而通過熱分析獲得的cp具有階躍變化特征，這個點是偏向溫度更低的范圍。另一方面，是由于 cp的增量是由液態金屬平移自由度所做的貢獻[28]，這說明玻璃化轉變與原子的運動有關，而原子的運動本身具有滯后性，從而導致了玻璃化轉變相對于時間和溫度具有滯后性，進一步說明Tg只是動力學參量。綜合研究從動力學和熱力學性質分別得到的玻璃化轉變溫度可知，用動力學性質計算的Tg與實驗結果更接近。

圖5 Cu66Ti34合金凝固過程比熱容cp隨溫度的變化Fig.5 Specific heat capacity cp vs temperature during solidification process of Cu66Ti34 alloy

3 結論

1) GEAM勢函數能夠較好地反映Cu66Ti34合金玻璃化轉變的動力學性質，通過偶分布函數第一谷的最小值與第一峰最大值之比獲得玻璃化轉變溫度為 600 K，與相近成分Cu50Ti50的實驗值接近。

2) 在模擬條件下根據均方位移數值，可以知道Cu的擴散能力要比Ti的強，且都在600 K開始呈現明顯凍結狀態。

3) 合金的生成焓隨溫度的降低而降低，比熱容與溫度具有二次函數關系?？梢岳胏P函數圖像來估計玻璃化轉變的發生。由于Tg的動力學性質，采用動力學方法計算的Tg與實驗結果更接近。

[1]LIN X H, JOHNSON W L. Formation of Ti-Zr-Cu-Ni bulk metallic glasses[J]. Journal of Applied Physics, 1995, 78(11):6514?6519.

[2]XU D H, DUAN G, JOHNSON W L. Unusual glass-forming ability of bulk amorphous alloys based on ordinary metal copper[J]. Physical Review Letters, 2004, 92(24):245503?245506.

[3]LEE S W, HUH M Y, FLEURY E, LEE J C.Crystallization-induced plasticity of Cu-Zr containing bulk amorphous alloys[J]. Acta Mater, 2006, 54(2): 349?355.

[4]DAI C L, GUO H, SHEN Y, LI Y, MA E, XU J. A new centimeter-diameter Cu-based bulk metallic glass[J]. Scripta Materialia, 2006, 54: 1403?1408.

[5]SHEN Y, MA E, XU J. A group of Cu(Zr)-based BMGs with critical diameter in the range of 12 to 18 mm[J]. Journal of Materials Science and Technology, 2008, 24(2): 149?152.

[6]JIA P, GUO H, Li Y, XU J, MA E. A new Cu-Hf-Al ternary bulk metallic glass with high glass forming ability and ductility[J].Scripta Materialia, 2006, 54: 2165?2168.

[7]TANG M B, ZHAO D Q, PAN M X, WANG W H. Binary Cu-Zr bulk metallic glasses[J]. Chinese Physics Letters, 2004, 21(5):901?903.

[8]XU D H, LOHWONGWATANA B, DUAN G, JOHNSON W L,GARLAND C. Bulk metallic glass formation in binary Cu-rich alloy series-Cu100?x(x = 34, 36, 38.2, 40 at.% )and mechanical properties of bulk Cu64Zr36glass [J]. Acta Materialia, 2004, 52(9):2621?2624.

[9]WANG D, LI Y, SUN B B, SUI M L, LU K, MA E. Bulk metallic glass formation in the binary Cu-Zr system[J]. Applied Physics Letters, 2004, 84 (20): 4029?4031.

[10]SATO S, TAKAHASHI Y, SANADA T, WAGATSUMA K,SUZUKI S. Small-angle X-ray scattering characterization of precipitates in Cu-Ti alloys[J]. Journal of Alloys and Compounds,2009, 447(1/2): 846?850.

[11]GHOSH G. First-principles calculations of structural energetics of Cu-TM (TM =Ti, Zr, Hf)intermetallics[J]. Acta Materialia,2007, 55: 3347?3374.

[12]OGATA S, SHIMIZU F, LI J, WAKEDA M, SHIBUTANI Y.Atomistic simulation of shear localization in Cu-Zr bulk metallic glass[J]. Intermetallics, 2006, 14(8/9): 1033?1037.

[13]LO Y C, CHOU H S, CHENG Y T, HUANG J C, MORRIS J R,LIAW P K. Structural relaxation and self-repair behavior in nano-scaled Zr-Cu metallic glass under cyclic loading:Molecular dynamics simulation[J]. Intermetallics, 2010, 18(5):954?960.

[14]WADLEY H N G, ZHOU X W, JOHNSON R A, NEUROCK M.Mechanisms, model and methods of vapor deposition[J].Progress in Material Science, 2001, 46(5): 329?377.

[15]ZHOU X W, WADLEY H N G, JOHNSON R A, LARSON D J,TABAT N, CERZO A, PETFORD-LONG A K, SMITH G D W,CLIFTON P H, MARTENS R L, KELLY T F. Atomic scale structure of sputtered metal multilayers[J]. Acta Materialia, 2001,49(19): 4005?4015.

[16]JOHNSON R A. Analytic nearest-neighbor model for fcc metals[J]. Physical Review B, 1988, 37(8): 3924?3931.

[17]WASEDA Y. The structure of non-crystalline materials[M]. New York: Mcgraw-Hill, 1980: 292?300.

[18]梁永超, 劉讓蘇, 朱軒民, 周麗麗, 田澤安, 劉全慧. 液態Mg7Zn3合金快速凝固過程中微觀結構演變機理的研究[J]. 物理學報, 2010, 59(11): 7930?7940.LIANG Yong-chao, LIU Rang-su, ZHU Xuan-min, ZHOU Li-li,TIAN Ze-an, LIU Quan-hui. Simulation study of evolution mechanisms of microstructures during rapid solidification of liquid Mg7Zn3alloy[J]. Acta Physica Sinica, 2010, 59(11):7930?7940.

[19]朱志雄, 張 鴻, 劉超峰, 齊衛宏, 易丹青, 李志成. Ni-Al合金凝固過程的分子動力學模擬[J]. 中國有色金屬學報, 2009,19(8): 1409?1416.ZHU Zhi-xiong, ZHANG Hong, LIU Chao-feng, QI Wei-hong,YI Dang-qing, LI Zhi-cheng. Molecular dynamics simulation for solidification process of Ni-Al alloys[J]. The Chinese Journal of Nonferrous Metals, 2009, 19(8): 1409?1416.

[20]WENDT H R, ABRAHAM F F. Empirical criterion for the glass transition region based on Monte Carlo simulations[J]. Physical Review Letters, 1978, 41(18): 1244?1246.

[21]ABRAHAM F F. An isothermal-isobaric computer simulation of the supercooled-liquid/glass transition region: Is the short-range order in the amorphous solid FCC[J]. Journal of Chemical Physics, 1980, 72(1): 359?365.

[22]NEUDECKER M, MAYR S G. Dynamics of shear localization and stress relaxation in amorphous Cu50Ti50[J]. Acta Materialia,2009, 57: 1437?1441.

[23]劉 俊, 趙九洲, 胡壯麒. Cu-50%Ni合金快速凝固過程中原子團簇演變的分子動力學模擬[J]. 金屬學報, 2005, 41(2):219?224.LIU Jun, ZHAO Jiu-zhou, HU Zhuang-qi. Molecular dynamic simulation of atomic clusters evolution Cu-50%Ni in alloy during rapid solidification[J]. Acta Metallurgica Sinica, 2005,41(2): 219?224.

[24]CAHN R W, HAASEN P. Physical metallurgy[M]. North Holland: Elsevier Science B V, 1996.

[25]?ZDEMIRKART S, TOMAK M, ULUDO?AN M, ?A?IN T.Structural, thermodynamical, and transport properties of undercooled binary Pd-Ni alloys[J]. Materials Science and Engineering A, 2006, 435/436: 736?744.

[26]楊 弘, 陳 民. 深過冷液態Ni2TiAl合金熱物理性質的分子動力學模擬[J]. 物理學報, 2006, 55(5): 2418?2421.YANG Hong, CHEN Min. A molecular dynamics simulation of thermodynamic properties of undercooled liquid Ni2TiAl alloy[J].Acta Physica Sinica, 2006, 55(5): 2418?2421.

[27]BAILEY N P, SCHI?TZ J, JACOBSEN K W. Simulation of Cu-Mg metallic glass: Thermodynamics and structure [J].Physical Review B, 2004, 69 (14): 144205?144215.

[28]陳 光, 傅恒志. 非平衡凝固新型金屬材料[M]. 北京: 科學出版社, 2004: 99.CHEN Guang, FU Heng-zhi. New metal materials of nonequilibrium solidification [M]. Beijing: Science Press, 2004: 99.