金納米管力學性能的分子動力學模擬*

蘇錦芳 宋海洋 安敏榮

(西安郵電大學理學院，西安 710121)

(2012年9月26日收到;2012年10月31日收到修改稿)

1 引言

納米材料由于其獨特的微結構[1,2]和特異力學性能如高強度[3]、超塑性[4]和高韌性[5,6]，納米技術正成為各國科技界所關注的焦點，被認為是21世紀最有前途的研究領域，其廣泛應用于物理、化學、生物、微電子等諸多學科.對于金屬納米材料，像 Mo[7,8]，Cu[9,10]，Ni[11]，Fe[12]，Al[13,14]，Mg[15,16]，Pt[17]，Ag[18]等力學性能已被廣泛研究，尺寸效應、溫度和應變率等會明顯地影響其拉伸力學性能，改變材料的變形機制[19-21].當材料尺寸減小到納米尺度時，它的硬度、強度和延展性都會發生很大的變化[22-24].納米結構的力學性能在未來的納米電子元件的發展中也占重要的地位.計算機模擬作為科學研究的重要手段，已被應用于多方面的學術研究，并取得豐碩的成果，計算機模擬方法具有重要的意義.目前，在納米材料力學性能的研究中，分子動力學模擬[25]是最為理想的模擬手段之一.

自從日本東京工業大學教授高柳邦夫領導的科研小組宣布成功研制出金的納米管以來，由于其獨特的性質而受到全世界科學家的關注，而對金納米管力學性能的研究還鮮有報道.金是良好的導體，但微觀的金納米結構卻不同于塊體金屬.類似于碳納米管[26-28]，金納米管也具有其獨特的性能，其在光電傳感器中具有較強的優勢[29].最近，Yang等[30]研究了手性結構的金納米管結合能及弦張力的影響，Zhang等[31]研究了表層與核心層的殘余應力和彈性模量對銅納米管的影響.因此對金屬類納米管力學性能的研究引起了世界各界學者的關注，而對金納米管力學性能的研究還鮮有報道.本文采用分子動力學模擬方法，研究了金納米管在拉伸和壓縮載荷下的力學性能，并研究了尺寸效應對金納米管拉伸變形過程中微觀結構的變化和楊氏模量的影響.

2 模擬方法

利用分子動力學方法研究了 〈100〉，〈110〉，〈111〉晶向金納米管在拉伸與壓縮載荷下的力學性能及半徑選取1.2—6.0 nm之間的金納米管的拉伸力學行為.在對三個晶向的研究中，金納米管的厚度均為0.6 nm，其對應的原子數分別為6580，6600，7660，不同半徑下的金納米管所選原子總數3480到21760不等.

圖1 金納米管的初始模型 (a)金納米管的截面圖;(b)金納米管的立體圖

初始模型如圖1所示，x，y，z方向均采用自由邊界條件，程序中采用“熱浴法”調節溫度，使系統控制在0.01 K，運動方程的積分采用velocity-verlet算法[25].分子動力學模擬結果是否準確取決于原子間相互作用勢函數的選取.為了能夠對金納米管進行準確的模擬，本文原子間相互作用勢函數選取鑲嵌原子勢(EAM)，系統的總勢能表示為

這里ξ代表有效配位數，rij代表兩個原子之間的距離，比截斷距離要小，r0是最近鄰原子距離，參數A，p，q和ξ分別由聚合能、晶格常數、體積彈性模量和剪切彈性常數的實驗數據決定.表1列出了當前模擬中Au的各項參數.

表1 Au的相關參數

在m面n方向的原子模擬的應力σmn由下式計算得出：

這里mi表示原子i的質量，Vi表示原子i指定的體積，NS表示在面積S范圍內所包含的粒子數，S定義為原子相互作用的范圍，rij為原子i與 j之間的距離.(2)式右邊部分第一項代表原子運動動力學的影響，并且是受溫度的影響;第二部分表示原子間的力和距離.

在模擬分析中采用開放性可視化工具OVITO(open visualization tool)[32]，用公共近鄰分析法CNA(common neighbor analysis)[33]分析了原子結構的轉變.

3 結果與討論

在模擬前先對系統進行弛豫，共弛豫50000步，使其能量最小化.模擬時使系統保持在0.01 K，時間步長選取為2 fs.

3.1 不同晶向金納米管力學性質的對比

應力-應變曲線關系能反映材料的基本力學性能.圖 2 為 fcc結構三個不同晶向 〈100〉，〈110〉和〈111〉金納米管的拉伸與壓縮載荷下的應力-應變曲線，研究的材料具有相同的長度和半徑.圖2中橫坐標代表應變，縱坐標為金納米管軸向所受的應力.在拉伸曲線達到第一個峰值時的應變為臨界應變，此峰值對應的應力為臨界應力(屈服強度)，其所對應的點為屈服點，屈服強度反映了材料對力作用的承受能力.

圖2 不同晶向金納米管拉伸與壓縮的應力-應變曲線

屈服點把應力-應變曲線分成了兩部分如圖2，一部分是屈服點之前的彈性階段，一部分是屈服點之后的塑性階段，這主要是由于不同階段的變形機制不同.在彈性階段，由于表面應力的存在使得拉伸載荷或壓縮載荷未超出金納米管的彈性限度;在屈服點處應力下降，這是因為此時部分晶格由于載荷的作用不再按規則的fcc結構排列，位錯開始發射，納米管在拉伸或壓縮載荷下積累的能量隨著位錯的移動而釋放，導致了應力的下降.應力應變曲線在達到屈服點后都有一個下降的過程，拉伸載荷下，〈110〉向的納米管在塑性階段其平均流動應力高于其他晶向，說明〈110〉向的金納米管具有較高的延展性，Coura等[34]在實驗上通過HRTEM觀察到的金納米線拉伸也有類似現象，由圖2的應力應變曲線也說明了不同晶向的金納米管其力學性能有顯著的不同.

為了進一步探索此刻的變形機制，圖3(a)和圖3(b)分別給出了在拉伸和壓縮載荷下，三個不同晶向 〈100〉，〈110〉，〈111〉的金納米管在屈服點處的原子結構.圖中可以明顯地觀測到在拉伸的兩端開始發射位錯，但是大部分原子仍處于原初始位置.模擬結果顯示，〈100〉，〈110〉，〈111〉晶向金納米管在拉伸載荷下的臨界應力分別為1.83，4.24，3.13 GPa，在壓縮載荷下的臨界應力分別為0.37，7.79，7.52 GPa.可以看出〈100〉晶向在拉伸和壓縮下的屈服強度都很小，這是因為金單晶各晶面之間的間距是不同的，各晶面的原子密度也不同.〈100〉面的金納米管其面原子間距大，面密度最小，位錯發射所需要克服的勢壘最大，易發生脆性斷裂.楊氏模量是固體材料中重要的性能參數，被視為衡量材料彈性變形的難易程度的指標，其值越大，使材料發生彈性變形的應力越大，即材料的剛度越大.楊氏模量相當于應力-應變曲線線性階段的斜率，可以從應變小于0.02的拉伸結果中線性擬合得出.圖3(a)中三個晶向的金納米管對應的楊氏模量分別為33.26，100.39，91.55 GPa.從模擬結果可看出，在加載負荷下〈100〉晶向的金納米管楊氏模量遠遠低于〈110〉向和 〈111〉向的金納米管，〈110〉和 〈111〉向的金納米管應力-應變曲線總體走勢一致，因此，〈110〉，〈111〉晶向的金納米管的力學性能優于〈100〉晶向.

圖 3 載荷下 〈100〉，〈110〉，〈111〉晶向金納米管剛剛發生位錯時的原子結構 (a)拉伸載荷下，從左到右分別為 〈100〉，〈110〉，〈111〉晶向的金納米管;(b)壓縮載荷下，從左到右分別為 〈100〉，〈110〉，〈111〉晶向的金納米管

圖4 拉伸載荷下〈111〉向金納米管不同應變的原子結構 (a)應變為0時的拉伸模型;(b)應變為0.020時的拉伸模型;(c)應變為0.055時的拉伸模型;(d)應變為0.288時的拉伸剖面圖

為了探索金納米管在載荷下的形變過程，圖4給出了拉伸載荷下〈111〉向不同應變下的原子結構，圖4(a)—(c)分別表示金納米管受到拉伸應變為0，0.020，0.055時的原子結構，分別對應于圖2中的〈111〉向的初始點、屈服點和拉伸應變0.055處，圖4(d)為拉伸應變為0.288時的原子剖面圖.圖4(a)為金納米管在未受到拉伸時的原子模型，此時原子全部處于原初始位置，圖4(b)為屈服點處的原子模型，可以看出金納米管此時發生了一定的形變，位錯開始發射.在彈性階段，應力-應變呈線性變化，在這一階段金納米管處于彈性形變.在塑性形變的過程中，金原子在發生位錯、滑移的同時不斷地重組，以使整個金納米管達到新的穩態.隨著拉伸的繼續，會出現頸縮現象，頸縮區域原子發生重構，形成晶格不完善的結構如圖4(c)，此時應變為0.055.隨著位錯的大量出現，原子不再按嚴格的規則排列，拉伸兩端的頸縮部分出現孔洞接著呈現單原子排列如圖4(d)，與實驗上觀察到的金納米線拉伸[34]的結果相符，直至局部斷裂.

3.2 截面尺寸對金納米管力學性能的影響

為了解截面尺寸對納米管拉伸力學性能的影響，圖5給出了7種不同半徑下的金納米管拉伸的應力-應變曲線，此處列出的金納米管的半徑分別為 1.2，2.2，3.0，4.0，4.4，5.0，6.0 nm.由 7 種不同尺寸金納米管的比較可看出，在拉伸載荷下，金納米管在半徑小于3.0 nm時，楊氏模量變化不大，在113 GPa到117 GPa之間上下浮動，屈服強度在4.46 GPa到4.75 GPa之間變化，尺寸效應對金納米管的楊氏模量和屈服強度影響均不明顯，在半徑大于3.0 nm時，金納米管的屈服強度明顯降低.這主要是因為直徑較小的金納米管曲率較大，其層間相互作用較大，使得其力學性能優越.金納米管表面原子占相當大的比例，這些表面原子配位數比較少，在空間方向失去原子價鍵形成斷鍵，有很高的表面剩余能，引起相應的表面張力.另外，表面原子的弛豫還會增加納米管的預壓應力，截面尺寸較小的金納米管會產生更多的預壓應力，需要較多的拉伸應力來克服預壓應力，因此截面尺寸小的金納米管屈服強度較大.隨著截面尺寸的增加，金納米管的預壓應力減少，其局部原子的活動自由度增大[35]，位錯發生的概率提高，從而使得金納米管強度降低，即金納米管承受載荷作用能力下降.

圖5 不同半徑下金納米管的應力-應變曲線

4 結論

本文采用分子動力學方法對金納米管的力學性能進行了研究和分析，在模擬中，采用EAM勢描述金原子之間的相互作用.根據計算結果，分析了金納米管沿不同晶向拉伸與壓縮的力學性能，及金納米管的半徑對其力學行為的影響.研究結果表明，不同晶向的金納米管在拉伸及壓縮過程中力學性能相差較大，〈100〉向的金納米管的力學性能相對其他晶向比較弱，其屈服強度和楊氏模量都遠遠小于其他晶向.當納米管的半徑小于3.0 nm時，金納米管的屈服強度沒有大的變化，而當半徑大于3.0 nm后，隨著半徑的增大，其屈服強度明顯降低.

[1]Liu Y L，GuiL J，Jin S 2012 Chin.Phys.B 21 096102

[2]Zhou G R，Teng X Y，Wang Y，Geng H R，Hur B Y 2012 Acta Phys.Sin.61 066101(in Chinese)[周國榮，滕新營，王艷，耿浩然，許甫寧2012物理學報61 066101]

[3]Lu L，Chen X，Huang X，Lu K 2009 Science 323 607

[4]McFadden S X，Mishra R S，Valiev R Z，Zhilyaev A P，Mukherjee A K 1999 Nature 398 684

[5]Zhang X，Wang H，Scattergood R O，Narayan J，Koch C C，Sergueeva A V，Mukherjee A K 2002 Appl.Phys.Lett.81 823

[6]Wang Y M，Chen M W，Zhou F H，Ma E 2002 Nature 419 912

[7]Jia M，LaiY Q，Tian Z L，Liu Y X 2009 Acta Phys.Sin.58 1139(in Chinese)[賈明，賴延清，田忠良，劉業翔2009物理學報58 1139]

[8]Wang J M，Hu W Y，LiX F，XiaoS F，Deng H Q 2010 Comput.Mater.Sci.50 373

[9]He A M，ShaoJ L，Wang P，Qin C S 2010 Acta Phys.Sin.59 8836(in Chinese)[何安民，邵建立，王裴，秦承森2010物理學報59 8836]

[10]Wang G H，Pan H，Ke F J，Xia M F，BaiY L 2008 Chin.Phys.B 17 259

[11]Huang D，Zhang Q，QiaoP Z 2011 Comput.Mater.Sci.50 903

[12]Wang Z G，Wu L，Zhang Y，Wen Y H 2011 Acta Phys.Sin.60 096105(in Chinese)[汪志剛，吳亮，張楊，文玉華2011物理學報60 096105]

[13]Pastor-Abia L，Caturla M J，SanFabian E，Chiappe G，Louis E 2011 Phys.Rev.B 83 165441

[14]An M R，Song H Y，Su J F 2012 Chin.Phys.B 21 106202

[15]Song H Y，LiY L 2012 Acta Phys.Sin.61 226201(in Chinese)[宋海洋，李玉龍2012物理學報61 226201]

[16]Song H Y，LiY L 2012 J.Appl.Phys.112 054322

[17]AoB Y，Xia J X，Chen P H，Hu W Y，Wang X L 2012 Chin.Phys.B 21 026103

[18]Liu Y H，GaoY J，Wang F Y，Zhu T M，ZhaoJ W 2011 Acta Phys.Chim.Sin.27 1341(in Chinese)[劉云紅，高亞軍，王奮英，朱鐵民，趙健偉2011物理化學學報27 1341]

[19]Gleiter H 1995 Nanostruct.Mater.6 3

[20]Qin Kun，Yang L M，Hu S S 2008 Chin.Phys.Lett.25 2581

[21]Deng C，Sansoz F 2010 Phys.Rev.B 81 155430

[22]SiegelR W 1996 Mater.Sci.Forum.235-238 851

[23]Sanders P G，YoungdahlC J，Weertman J R 1997 Mater.Sci.Eng.A 234-236 77

[24]Koch C C，Malow T R 1999 Mater.Sci.Forum.312-314 565

[25]Honeycutt R W 1970 Methods Comput.Phys.9 136

[26]LiR，Hu Y Z，Wang H，Zhang Y J 2008 Chin.Phys.B 17 4253

[27]Song H Y，Zha X W 2009 J.Phys.D：Appl.Phys.42 225402

[28]Song H Y，Zha X W 2009 Phys.Lett.A 373 682

[29]Yang G M，Xu G L，LiY T，Xu F，Chang Y B，Yang Y H 2010 Mater.Rev.24 81(in Chinese)[楊光明，徐國良，李月婷，徐鳳，常艷兵，楊云慧2010材料導報24 81]

[30]Yang J H，LiB，Zhang Q J，Chen L 2012 Phys.Lett.A 376 2707

[31]Zhang J，Wang C Y，Chowdhury R，AdhikariS 2012 Appl.Phys.Lett.101 093109

[32]StukowskiA 2010 Modelling Simul.Mater.Sci.Eng.18 015012

[33]Faken D，Jonsson H 1994 Comput.Mater.Sci.2 279

[34]Coura P Z，Legoas S B，Moreira A S，SatoF，Rodrigues V，Dantas S O，Ugarte D，GalvaoD S 2004 NanoLett.4 1187

[35]Horstermeyer M F，Baskes M I 1999 J.Eng.Mater.Technol.121 114