碳納米錐力學(xué)特性的分子動力學(xué)研究*

李明林 林凡 陳越

(福州大學(xué)機(jī)械工程及自動化學(xué)院，福州 350108)

(2012年6月6日收到;2012年8月6日收到修改稿)

1 引言

自從1991年Iijima[1]發(fā)現(xiàn)碳納米管以來，各種碳納米結(jié)構(gòu)吸引了眾多科技工作者的廣泛關(guān)注和深入研究，包括碳納米管[2]、石墨烯[3]、碳納米錐[4]等.這主要是因?yàn)檫@些碳納米結(jié)構(gòu)存在著諸多奇特且異常顯著的物理化學(xué)性能，如高強(qiáng)度、低密度、高導(dǎo)電率、高導(dǎo)熱率和大的比表面積等特性.基于這些卓著的物理化學(xué)性能，這些碳納米結(jié)構(gòu)在各個工程領(lǐng)域存在著廣泛而巨大的應(yīng)用潛力，例如納米機(jī)電系統(tǒng)[5].碳納米管和石墨烯由于其在制備能力方面所取得的突破性進(jìn)展而較早受到廣泛而深入的研究.緊隨其后的是碳納米錐由于其獨(dú)特的錐形結(jié)構(gòu)正受到越來越多的科技工作者的關(guān)注.

碳納米錐的研究最早見于1994年Ge和Sattler[6]對其結(jié)構(gòu)的預(yù)測.隨后，Krishnan等[7]通過實(shí)驗(yàn)證實(shí)碳納米錐存在5種錐角.而Iijima等[8]在1999年將單壁碳納米錐的每小時產(chǎn)量提高到了10 g.隨著對碳納米錐制備方法的成功掌握，其熱學(xué)[9]、電學(xué)[10]和力學(xué)[11]特性的研究也逐漸鋪展開來.然而，至今針對碳納米錐力學(xué)特性的研究文獻(xiàn)尚不多見.2004年Wei和Srivastava[12]通過連續(xù)介質(zhì)彈性理論和分子動力學(xué)模擬推測出單壁碳納米錐的楊氏模量為其等價碳納米管的cos4θ.同年，Jordan和Crespi[13]利用分子動力學(xué)模擬發(fā)現(xiàn)碳納米錐的封閉頂部受球體緩慢擠壓后不發(fā)生化學(xué)鍵的斷裂，而是發(fā)生完美的手性翻轉(zhuǎn).這類似于“翻口袋”現(xiàn)象，即錐體的外表面變?yōu)閮?nèi)表面，而內(nèi)表面翻轉(zhuǎn)為外表面.2007年，Tsai和Fang[14]利用同樣的方法研究顯示碳納米錐的錐角越大，其結(jié)構(gòu)越穩(wěn)定.同年，Wei等[15]通過分子動力學(xué)模擬研究等長度的單壁碳納米錐的拉伸彈性和塑性特性.研究結(jié)果顯示，錐角越大，碳納米錐的彈性應(yīng)變極限越小，而其斷裂強(qiáng)度則越大.由此，他們得出錐角為零的等長度碳納米管具有最小的斷裂強(qiáng)度和最大的彈性應(yīng)變極限.同年，Liew等[16]則通過分子動力學(xué)模擬研究單壁碳納米錐臺的受壓曲屈行為.研究結(jié)果顯示，錐臺的高度與底半徑之比固定的情況下，錐頂半徑越小其剛度越大.2012年，F(xiàn)irouz-Abadi等[17]通過分子動力學(xué)模擬研究懸壁的單壁碳納米錐的橫向共振特性.其研究結(jié)果表明，一般情況下旋轉(zhuǎn)位移角度(disclination angle)和錐高越小其共振頻率越大，但240°角除外.此外，研究也顯示懸臂的單壁碳納米錐的共振頻率幾乎不受環(huán)境溫度的影響.

目前國內(nèi)針對碳納米錐的研究雖然也取得了一些有價值的成果，但關(guān)于碳納米錐的力學(xué)特性分析卻少有報(bào)道.2007年沈海軍等[18]基于Brenner勢分子動力學(xué)模擬方法研究平板與碳納米錐作用過程的能量、力和構(gòu)形的演化.本文針對由相當(dāng)數(shù)量或質(zhì)量的碳原子所組成的各種錐角和長度不等的碳納米錐和碳納米管，通過經(jīng)典分子動力學(xué)模擬方法，研究其在恒溫條件下的拉壓力學(xué)行為，以期獲得較為完整的拉壓載荷-應(yīng)變關(guān)系曲線及其承載能力.

圖1 碳納米結(jié)構(gòu)圖 (a)碳納米錐;(b)碳納米管

2 模型和模擬方法

碳納米錐和碳納米管的結(jié)構(gòu)如圖1所示.本文主要研究由約580個碳原子組成的4組碳納米錐：(112.88°，1.33 nm)C1，(83.62°，2.0 nm)C2，(60.0°，2.68 nm)C3，(38.94°，3.57 nm)C4，括號內(nèi)分別表示錐角α和錐高h(yuǎn).與碳納米錐形成對照組的碳納米管為4.0 nm長的(15，0)單壁碳納米管(carbon nanotube，CNT).這些納米結(jié)構(gòu)的碳原子數(shù)分別為580，580，576，576，576.隨著錐角的減小，碳納米錐的錐頂原子的個數(shù)分別為5，4，3，2.

本文模擬采用大規(guī)模原子/分子并行模擬器(large-scale atomic/molecular massively parallel simulator，LAMMPS[19])的經(jīng)典分子動力學(xué)模擬方法.經(jīng)典分子動力學(xué)本質(zhì)上是一種粒子方法，因其求解對象是基于牛頓第二定律的粒子動力學(xué)控制方程，即

其中，mi和ri分別表示第i個粒子的質(zhì)量和空間位置，V是系統(tǒng)的經(jīng)驗(yàn)作用勢，而?表示空間梯度.對粒子位置和速度的積分運(yùn)算采用Verlet方法.碳納米結(jié)構(gòu)的原子間相互作用勢通常采用Brenner半經(jīng)驗(yàn)勢，也叫Tersoff-Brenner勢[20].Brenner勢是對Tersoff勢的進(jìn)一步修正，使之能更恰當(dāng)?shù)胤从矯—C化學(xué)鍵的形成和斷裂過程，已成功用于預(yù)測和驗(yàn)證碳納米結(jié)構(gòu)的各種力學(xué)特性.由于第二代Brenner勢 (reactive empirical bond order potential，REBO勢)[21]主要考慮原子間的短程相互作用，本文增加了考慮原子間長程相互作用的Lennard-Jones勢，即

其中，rs和rl分別表示短程作用勢和長程作用勢的截?cái)喟霃剑謩e取值0.2 nm和0.85 nm;VR和VA分別表示短程作用勢的斥力項(xiàng)和引力項(xiàng)，bij是原子間的反應(yīng)經(jīng)驗(yàn)鍵序;ε和σ分別是一對原子間的平衡勢能和平衡距離，對C—C原子而言，ε=4.2038×10-3eV，σ=0.34 nm[22].

考慮到納米結(jié)構(gòu)受拉/壓過程的準(zhǔn)靜載恒溫條件，以及連續(xù)的加載和力學(xué)響應(yīng)過程，本文對碳納米結(jié)構(gòu)力學(xué)特性的分子動力學(xué)模擬實(shí)驗(yàn)步驟如下：1)固定碳納米結(jié)構(gòu)底部2層原子，并通過共軛梯度法對此結(jié)構(gòu)進(jìn)行能量最小化，能量優(yōu)化指標(biāo)為10-10eV，力優(yōu)化指標(biāo)為10-10eV/?A，隨后通過10 ps的趨衡模擬，使其達(dá)到恒溫條件下的平衡;2)針對碳納米錐和碳納米管的軸向拉伸，均采用逐步定向移動其自由端2層原子，并每隔1000步記錄此2層原子所受合力的軸向分量，由于橫截面的應(yīng)力沿軸向位置變化，本文主要統(tǒng)計(jì)軸向載荷-軸向應(yīng)變關(guān)系曲線;3)針對碳納米錐和納米管的軸向壓縮，均采用逐步移動剛性石墨烯平板進(jìn)行壓縮.軸向載荷以平板所受合力的軸向分量表示.系統(tǒng)環(huán)境溫度采用Berendsen溫控方法進(jìn)行恒溫控制，其指數(shù)衰減系數(shù)為0.1 ps.

為確定合適的模擬時間步長、拉壓速率和環(huán)境溫度，首先在較低的溫度(0.1 K)和加載速率 (0.001 nm/ps)條件下，分別采用 0.1，0.2，0.5，0.8和1.0 fs時間步長觀察碳納米錐C2的拉伸載荷-應(yīng)變曲線.模擬結(jié)果表明，時間步長越大，C2的受拉載荷極限和應(yīng)變極限越小.在0.1，0.2和0.5 fs下載荷極限均值約為279.66 nN(相對偏差低于0.6%)，載荷單位由實(shí)驗(yàn)的eV/?A變換為nN(1 eV/?A≈1.602 nN).在0.8和1.0 fs下載荷極限約為208.26 nN(與0.5fs相差達(dá)約25%).因此，后續(xù)模擬的時間步長均采用0.5 fs.緊接著，在0.1 K溫度和0.5 fs步長條件下，分別取用0.0005，0.001，0.003，0.005，0.007，0.01，0.02，0.05，0.07 和 0.1 nm/ps速率拉伸C2.結(jié)果表明，拉伸載荷極限隨著拉伸速率的增加而減小.在[0.0005，0.01］nm/ps范圍，載荷極限的均值約為278.88 nN(相對偏差低于0.5%);當(dāng)拉伸速率大于0.02 nm/ps時，載荷極限偏離上述均值均大于2.5%.因此，綜合考慮模擬精度和時間消耗，后續(xù)模擬均采用0.01 nm/ps的加載速率.最后，在0.5 fs時間步長和0.01 nm/ps加載速率條件下，分別觀察 0.1，1，10，30，100，300 和 500 K 溫度下 C2 的拉伸載荷極限.結(jié)果表明，隨著溫度升高，拉伸載荷極限均出現(xiàn)不同程度的降低.在[0.1，10］K低溫區(qū)間，載荷極限的均值約為276.88 nN(相對偏差低于0.7%).在[30，500］K區(qū)間，載荷極限加速降低至最低158.3 nN.從30 K起，載荷極限已偏離低溫均值約13.2%至240.3 nN.雖然載荷極限在0.1至10 K的低溫環(huán)境變化不大，但溫度升高時熱擾動對載荷-應(yīng)變曲線的影響漸趨明顯.因此，為減少熱擾動對結(jié)構(gòu)本征力學(xué)性能的影響，后續(xù)模擬實(shí)驗(yàn)均在0.1 K下進(jìn)行.綜上所述，針對碳納米結(jié)構(gòu)力學(xué)特性的分子動力學(xué)模擬，合適的模擬參數(shù)確定為：0.1 K的系統(tǒng)溫度，0.5 fs的時間步長和0.01 nm/ps的加載速率.

3 結(jié)果與討論

3.1 拉伸力學(xué)特性分析

受軸向拉伸時，等量碳原子組成的碳納米錐和納米管的載荷-應(yīng)變關(guān)系曲線如圖2所示.軸向應(yīng)變定義為納米結(jié)構(gòu)軸向長度的變化量與其原長的比值的百分?jǐn)?shù).與文獻(xiàn)[15］的分步拉伸-松弛的實(shí)驗(yàn)結(jié)果不同，圖2表明碳納米錐的軸向載荷隨著應(yīng)變的增加趨近于線性遞增，未見明顯非線性變化.碳納米管的軸向載荷與應(yīng)變則以非線性變化為主.對于碳納米錐，其拉伸應(yīng)變極限隨著錐角的增大而增加;但其拉伸載荷極限由C4增加到C2達(dá)到最大值后，C1的載荷極限反而陡然降低.對于碳納米管，其拉伸載荷極限與C1相當(dāng)，而應(yīng)變極限則介于C2和C3之間.這些納米結(jié)構(gòu)的拉伸載荷極限和應(yīng)變極限數(shù)值分別列于表1.拉伸實(shí)驗(yàn)的系統(tǒng)環(huán)境溫度通過Berendsen方法控制在0.1 K，以盡量減少熱擾動對結(jié)構(gòu)本征力學(xué)性能的影響.碳納米錐C2受拉時的系統(tǒng)溫度隨時間的變化曲線如圖3所示，嵌入圖為曲線的部分放大圖形.可見，在C2受拉過程中，系統(tǒng)溫度除了在破壞階段附近出現(xiàn)激烈波動以外，在其他過程基本保持恒定.因?yàn)橄到y(tǒng)的熱力學(xué)溫度T主要取決于系統(tǒng)動能Ke，即T=2/3(Ke/N/kB)，N為系統(tǒng)原子數(shù)，kB是Boltzmann常數(shù).通過觀察C2拉伸模擬的構(gòu)形演化，系統(tǒng)溫度的突變主要是由C—C原子鍵斷裂引起的碳原子飛逸脫離速度的激增產(chǎn)生的.

圖2 拉伸軸向載荷-應(yīng)變率曲線圖

碳納米錐和碳納米管在拉伸時顯現(xiàn)出截然不同的構(gòu)形演變特點(diǎn).圖4(a)給出C1碳納米錐拉伸初始時的構(gòu)形和拉伸斷裂之前的構(gòu)形.可見，碳納米錐出現(xiàn)受拉橫向膨脹.然而并未如文獻(xiàn)[15］所述，碳納米管的拉伸應(yīng)變最大而其拉伸極限載荷最小.這可能是由于碳納米管的拉伸強(qiáng)度與其手性和直徑有關(guān)，而本文只考慮手性為(15，0)的碳納米管.由圖2碳納米管的載荷-應(yīng)變關(guān)系曲線可以看出碳納米管在拉伸過程中經(jīng)歷了彈性、屈服和強(qiáng)化等明顯的塑性材料拉伸變形階段.在達(dá)到其拉伸強(qiáng)度極限之后，碳納米管出現(xiàn)了明顯的階段性頸縮變形(圖2CNT曲線的階梯狀部分).最后，碳納米管形成了藕斷絲連式的單列碳原子鏈[23](如圖4(b)所示).對于碳納米錐，類似單列碳原子鏈卻是只出現(xiàn)在C3的拉伸實(shí)驗(yàn).之前的模擬結(jié)果表明，碳納米錐C2在高溫(500 K)條件下也能形成單列原子鏈.然而，關(guān)于通過拉伸碳納米錐形成單列碳原子鏈的機(jī)理及其影響因素仍有待進(jìn)一步的研究.

表1 拉伸載荷極限和拉伸應(yīng)變極限

圖3 C2拉伸過程的溫度-時間變化曲線，嵌入圖為局部放大圖

圖4 拉伸構(gòu)形演變特點(diǎn) (a)膨脹的C1構(gòu)形;(b)頸縮的CNT構(gòu)形

3.2 壓縮力學(xué)特性分析

等量碳原子組成的碳納米錐和納米管在受石墨烯剛板軸向壓縮時的載荷-應(yīng)變關(guān)系曲線如圖5所示.壓縮應(yīng)變定義為石墨烯剛板和納米結(jié)構(gòu)固定端之間距離的變化量與納米結(jié)構(gòu)原長的比值的百分?jǐn)?shù).總體而言，碳納米錐的受壓載荷極限和應(yīng)變極限均隨著錐角的增加而增大.然而，錐角最大的C1受壓載荷極限卻最小.對于C1，C2，C3和碳納米管，受壓載荷極限定義為圖5曲線中載荷的最大值;對于C4，受壓載荷極限則取其載荷-應(yīng)變曲線的第一個峰值(點(diǎn)1).表2列出了各納米結(jié)構(gòu)的受壓載荷極限和應(yīng)變極限數(shù)值.顯然，碳納米管在這些納米結(jié)構(gòu)中的受壓力學(xué)特性既不突出也不遜色.

表2 壓縮載荷極限和應(yīng)變極限

圖5 壓縮軸向載荷-應(yīng)變率曲線圖

觀察其構(gòu)形演化過程，碳納米錐的受壓載荷極限均出現(xiàn)在其屈曲的臨界狀態(tài).然而，碳納米管所受載荷還沒達(dá)到最大值之前，其結(jié)構(gòu)已先后出現(xiàn)末端原子外擴(kuò)和徑向屈曲，如圖6(a)和圖6(b)所示.在經(jīng)受最大壓載之后，碳納米管承載能力隨之下降.此時碳納米管正在另一區(qū)域逐步形成新的徑向屈曲結(jié)構(gòu)，如圖6(c)所示.待到新結(jié)構(gòu)穩(wěn)定之后，碳納米管的承載能力有所恢復(fù)(見圖5曲線上點(diǎn)2至3段).隨著壓縮的繼續(xù)，兩個徑向屈曲結(jié)構(gòu)通過碳納米管的扭轉(zhuǎn)變形逐漸融合成穩(wěn)定的側(cè)向屈曲變形結(jié)構(gòu)，如圖6(d)所示.在這過程中碳納米管的承載能力迅速減弱至圖5曲線上點(diǎn)4位置.注意到圖5的C4載荷-應(yīng)變曲線存在兩個峰值，通過觀看C4受拉時構(gòu)形的演化過程發(fā)現(xiàn)，第一個峰值對應(yīng)C4頂端2個碳原子受壓側(cè)向滑移的臨界狀態(tài)(圖7(a))，而第二個峰值則是C4發(fā)生明顯側(cè)向失穩(wěn)(圖7(b))的臨界載荷.雖然C3受壓也是發(fā)生側(cè)向屈曲，但由于其頂端原子側(cè)移并不明顯且?guī)缀跖c側(cè)向屈曲同時發(fā)生，因此其載荷-應(yīng)變曲線只出現(xiàn)一個峰值.此外，C2和C1受壓屈曲變形均為頂端原子沿軸向向內(nèi)塌陷并逐漸擴(kuò)大，進(jìn)而最終發(fā)生類似“翻口袋”的完美手性翻轉(zhuǎn)[13]，如圖8(a)和8(b)所示.

圖6 CNT階段受壓的構(gòu)形圖 (a)頂端原子外擴(kuò);(b)首個徑向屈曲;(c)雙向徑向屈曲;(d)側(cè)向屈曲

圖7 C4的受壓構(gòu)形圖 (a)頂端原子側(cè)向滑移;(b)側(cè)向屈曲

圖8 碳納米錐受壓頂端內(nèi)陷構(gòu)形圖 (a)C1;(b)C2

4 結(jié)論

本文通過分子動力學(xué)方法模擬對等量碳原子組成的碳納米錐的拉伸和壓縮實(shí)驗(yàn)，獲得其受拉/壓載荷-應(yīng)變關(guān)系曲線、載荷極限和應(yīng)變極限、構(gòu)形演化等特性，并與碳納米管進(jìn)行對比研究.研究結(jié)果表明：1)碳納米錐的受拉/壓載荷極限隨著錐角的增大先增大后減小，而其受拉/壓應(yīng)變極限均隨著錐角的增大而增大，其中，83.62°錐角的碳納米錐的受拉/壓載荷極限最大，而112.88°錐角的C1受拉/壓應(yīng)變極限最大，這些碳納米錐的力學(xué)特性毫不遜色于等量原子組成的碳納米管;2)在受拉過程中，與碳納米管不同的是碳納米錐并沒有呈現(xiàn)出明顯的屈服和強(qiáng)化階段特點(diǎn)，且其受拉載荷-應(yīng)變曲線在斷裂之前近似于線性變化;3)在受壓過程中，與碳納米管呈現(xiàn)出豐富的徑向屈曲/扭轉(zhuǎn)/側(cè)向屈曲組合變形不同，碳納米錐呈現(xiàn)出較為單一的軸向手性反轉(zhuǎn)(83.62°和112.88°錐角的C2和C1)或側(cè)向屈曲(38.94°和60°錐角的C4和C3).本文的研究結(jié)果表明，碳納米錐的力學(xué)特性可媲美碳納米管而應(yīng)用于廣泛的納米科技工程領(lǐng)域，例如納米傳感器或納米復(fù)合材料.

[1］Iijima S 1991Nature354 56

[2］Mahar B，Laslau C，Yip R，Sun Y 2007Sens.J.IEEE.7 266

[3］Hill E W，Vijayaragahvan V，Novoselov K 2011Sens.J.IEEE.11 3161[4］Adisa O O，Cox B J，Hill J M 2011J.Phys.Chem.C 115 24528

[5］Allen M J，Tung V C，Kaner R B 2010Chem.Rev.110 132

[6］Ge M，Sattler K 1994Chem.Phys.Lett.220 192

[7］Krishnan A，Dujardin E，Treacy M M J，Hugdahl J，Lynum S，Ebbesen T W 1997Nature388 451

[8］Iijima S，Yudasaka M，Yamada R，Bandow S，Suenaga K，Kokai F，Takahashi K 1999Chem.Phys.Lett.309 165

[9］Yang N，Zhang G，Li B 2008Appl.Phys.Lett.93 243111

[10］Lu X，Yang Q，Xiao C，Hirose A 2006Appl.Phys.A：Mate.Sci.Proc.82 293

[11］Zhang S，Yao Z，Zhao S，Zhang E 2006Appl.Phys.Lett.89 131923

[12］Wei C Y，Srivastava D 2004Appl.Phys.Lett.85 2208

[13］Jordan S P，Crespi V H 2004Phys.Rev.Lett.93 255504

[14］Tsai P C，F(xiàn)ang T H 2007Nanotechnology18 105702

[15］Wei J X，Liew K M，He X Q 2007Appl.Phys.Lett.91 261906

[16］Liew K M，Wei J X，He X Q 2007Phys.Rev.B 75 195435

[17］Firouz-Abadi R D，Amini H，Hosseinian A R 2012Appl.Phys.Lett.100 173108

[18］Shen H J，Shi Y J 2007J.Atom.Mole.Phys.24 883(in Chinese)[沈海軍，史海進(jìn)2007原子與分子物理學(xué)報(bào)24 883］

[19］Plimpton S 1995J.Comp.Phys.117 1

[20］Brenner D W 1990Phys.Rev.B 42 9458

[21］Brenner D W，Shenderova O A，Harrison J A，Stuart S J，Ni B，Sinnott S B 2002J.Phys.：Cond.Matt.14 783

[22］Liew K M，Wong C H，He X Q，Tan M J 2005Phys.Rev.B 71 075424[23］Fu C X，Chen Y F，Jiao J W 2008Sci.ChinaE：Tech.Sci.38 411(in Chinese)[付稱心，陳云飛，焦繼偉2008中國科學(xué)E輯：技術(shù)科學(xué)38 411］