單層MoS2力學性能的有限元分析

楊 璐, 齊 博, 楊宏旭, 王天韻

(沈陽工業(yè)大學 建筑與土木工程學院, 沈陽 110870)

單層MoS2是近些年發(fā)展起來的新型材料,屬于過渡金屬硫化物,是一種類石墨烯結構和性能的二維材料,主要采用微機械力剝離、離子插層和液相剝離法等為主的“自上而下”的剝離法,以及以化學氣相沉積法、水熱/溶劑熱合成法等為主的“自下而上”的合成制備方法.該材料具有優(yōu)異的電學、光學以及力學性能,具有廣闊的市場前景,可用作柔性光電器件、光電探測器、類腦邏輯器件和催化劑等[1-3].

單層MoS2的半導體性質允許其克服石墨烯的零帶隙,同時分享石墨烯在電子應用中的許多優(yōu)勢[4].為了確?；趩螌覯oS2的傳感器和納米材料在制作和使用期間保持其結構完整性,防止其出現(xiàn)機械失效,有必要對單層MoS2的力學性能進行更全面的了解.對于MoS2力學性能的研究大多數(shù)人采用了原子力顯微鏡(AFM)方法、密度泛函理論(DFT)、分子動力學方法(MD)得出其彈性模量、泊松比、破壞強度等力學參數(shù)[5-10].采用AFM試驗方法對所需儀器的要求較高,采用DFT和MD方法的模擬計算代價太大,經(jīng)常用來模擬少量分子和原子的系統(tǒng),以及相對短壽的現(xiàn)象.

部分學者還在致力于尋找其他高效準確的方法來研究納米材料的力學性能,Li等[11]提出在三維空間中用空間框架結構來等效石墨烯,采用有限元方法得出了石墨烯的彈性模量以及剪切模量,為其他共價鍵合材料的模擬提供了一個新思路;韋琦等[12]計算了石墨烯的彈性模量、泊松比、剪切模量等力學參數(shù);張續(xù)等[13]利用石墨烯的Morse勢能函數(shù)作為C—C鍵的本構方程,得出了含缺陷石墨烯的抗拉強度.

同為二維材料的MoS2與石墨烯有相似的結構,為了探究出更加高效的方法,采用有限元法模擬MoS2力學性能并探究其正確性.本文計算了MoS2的彈性模量,并得出單層MoS2彈性模量具有尺寸效應的結論,通過ABAQUS軟件模擬出了MoS2的破壞強度,為采用有限元法研究MoS2提供理論借鑒.然而采用有限元法的計算成本較低,只要采用足夠精確的化學鍵拉伸、彎曲、扭轉常數(shù)等參數(shù),在對MoS2的屈曲、振動、拉伸斷裂等線性或非線性范圍內的力學特征進行模擬時也能通過有限元法得到較為精確的結果.

1 有限元模型建立

1.1 MoS2模型建立

單層MoS2中間一層由Mo原子的一個平面組成,兩邊連著兩個S的原子平面,具體的原子結構如圖1所示.在有限元計算中采用桿系結構來模擬Mo原子和S原子之間的共價鍵,桿系結構之間采用剛性連接.本文采用AutoCAD完成三維建模,再導入到大型商用有限元模型ABAQUS中進行計算,模型如圖2所示,鍵長為0.241 nm[14].

圖1 MoS2 球棍模型

圖2 MoS2有限元模型

1.2 線性材料參數(shù)確定

在有限元模擬中采用圓形截面梁單元模擬Mo和S之間的共價鍵,共價鍵鍵長r0=0.241 nm.由于有限元方法需要材料屬性以及單元截面參數(shù),所以需要將微觀分子計算與宏觀結構力學建立聯(lián)系.在分子力學中分子總能量等于分子勢能與分子動能相加,總能量可表示為

U=∑Ur+∑Uθ+∑Uτ+∑Uω+∑Uvdw

(1)

式中:Ur為化學鍵拉伸能;Uθ為化學鍵彎曲能;Uτ為化學鍵扭轉能;Uω為平面外扭轉能;Uvdw為非鍵范德華相互作用.

一般來說,位能的作用主要來自前四項,為了方便研究,將平面內扭轉能與平面外扭轉能合并成一項,采用諧振子勢表示各勢能函數(shù),可得

(2)

(3)

(4)

式中:kr為化學鍵拉伸長度;kθ為化學鍵彎曲常數(shù);kτ為化學鍵扭轉常數(shù);Δr、Δθ和Δφ分別為化學鍵拉伸增量、鍵角轉動增量和鍵扭轉增量.

在結構力學中,框架結構中的梁也可以產生類似的拉伸變形和彎曲變形.在彈性小變形范圍內,拉伸應變能、彎曲應變能、扭轉應變能分別表示為

(5)

(6)

(7)

式中:A為截面面積,A=πd2/4,d為截面直徑;I為截面慣性矩;E為彈性模量;ΔL、Δα和Δβ分別為梁軸向拉伸變形量、梁端點的轉角、梁端點的扭轉角;L為梁長度;F、M、T分別為梁單元所受的拉伸作用力、彎矩、扭矩;G為剪切模量;J為截面極慣性矩.

將分子勢能與應變能等效,可得

(8)

(9)

(10)

進而得出

Li等[11]通過計算kτ得到石墨烯的泊松比不在合理范圍內,且kτ對彈性模量的影響很小.假設圓截面梁單元的泊松比為0,Nicolini等[14]將不同諧振子勢對比,認為諧振子勢參數(shù)在預測MoS2的結構和振動性能方面有較好的應用前景,其參數(shù)為

計算得到截面直徑d=0.317 nm,彈性模量E=452.859 GPa.

1.3 非線性材料參數(shù)確定

本文對理想單層MoS2的彈性性能進行了詳細描述,但在彈性范圍內研究MoS2納米帶在大變形下的拉伸斷裂等性能是不準確的.因此,采用結構力學的方法來研究MoS2的非線性特征,需要將MoS2的非線性本構等效成圓形截面梁的本構方程.

與研究單層MoS2納米帶彈性性能類似,建立同樣的空間框架結構.MoS2納米帶非線性模型的建立主要體現(xiàn)在Mo—S鍵的非線性上,采用Morse勢參數(shù)作為等效圓截面梁的本構方程[14].Morse勢的總方程可以表示為

E=Estretch+Eangle

(11)

Estretch=De[1-e-α(r-r0)]2

(12)

F=2Deα(1-e-αΔr)e-αΔr

(13)

F=2Deα(1-e-αr0)e-αr0

(14)

圖3 圓形截面梁應力應變曲線

張續(xù)[13]和王少培等[15]在采用Morse勢研究石墨烯時分別采用在應力最大值的應變以及應變?yōu)?%前的曲線斜率作為彈性模量值.經(jīng)過計算參數(shù)的取值均合理,表明材料屈服時應變的取值對結果影響不大.本文將應變?yōu)?%視為材料屈服.在應變達到16%時,此時Mo—S鍵間應力值達到最大,承受的力最大為2.273 nN.當原子間的距離增大時,原子之間的吸引力逐漸變大,2.273 nN即為原子間最大的吸引力.隨后原子間的作用力逐漸減小,在Mo—S鍵承受的力最大時,此時的Mo—S鍵可以認為發(fā)生斷裂,等效為Mo—S鍵的圓形截面梁發(fā)生斷裂.計算得到非線性階段的屈服應力和應變,通過Mo—S鍵的Morse勢能函數(shù)推導出結構力學中Mo—S鍵等效成的圓形截面梁的彈塑性本構.通過在ABAQUS中輸入圓形截面梁的彈塑性本構,計算其拉伸破壞時的破壞應力以及應變.

1.4 荷載施加

扶手型和鋸齒型的區(qū)別如圖4所示.沿著x軸方向的邊界形狀似鋸齒狀,所以研究MoS2在x軸方向上的性質表述為鋸齒型MoS2;沿著y軸方向的邊界形狀似扶手狀,所以研究MoS2在y軸方向上的性質表述為扶手型MoS2.

圖4 MoS2手性的定義

采用線性本構關系時,該方法是在彈性范圍內應變能與分子勢能等效,在施加力時只能施加小變形才能保證計算相對準確.采用非線性本構關系時,同樣的約束方式,在另一端可施加較大的位移荷載.

荷載施加示意圖如圖5所示.采用一端約束,另一端單向拉伸的方法.對于扶手型MoS2,在最右側原子上沿x方向施加0.001 nm的微小形變,在另一側施加x方向以及z方向上的約束.對于鋸齒型MoS2,在最上邊每個原子上沿y方向施加0.001 nm的微小形變,在另一側施加y方向以及z方向上的約束.

圖5 單層MoS2加載示意圖

2 線性本構下彈性模量的計算

2.1 模型尺寸確立

采用8種不同尺寸的單層MoS2,分別計算其在扶手型和鋸齒型方向下的彈性模量大小.由于x方向與y方向的長度不是嚴格意義上的相等,所以在保證結構完整的情況下,盡可能接近矩形大小,基于此8個工況的具體尺寸如表1所示.

表1 各工況MoS2尺寸以及單元信息

2.2 彈性模量計算

在單向拉伸試驗中,彈性模量是應力與應變的比值.對于扶手型MoS2,通過施加小變形,計算所對應方向支座反力Fx,并通過式(15)計算彈性模量大小.對于鋸齒型MoS2,通過施加小變形,計算所對應方向的支座反力Fy,并通過式(16)計算彈性模量大小.

(15)

(16)

式中:t為MoS2薄片厚度;Lx為x方向長度;Ly為y方向長度;ΔLx為施加在x方向上的小變形;ΔLy為施加在y方向上的小變形;A0和A1分別為橫向和縱向MoS2薄片橫截面面積大小,采用Lxt和Lyt計算.無論是對石墨烯還是MoS2的研究中,其厚度的選取存在著較大爭議,厚度的合理選取也是有限元模擬結果是否準確的關鍵,本文選取模型的厚度為0.61 nm[16].

2.3 數(shù)值模擬計算結果與分析

基于荷載施加方法,以及式(15)、(16)分別計算不同尺寸下扶手型和鋸齒型MoS2的彈性模量.為了更直觀地看出其變化,計算相鄰兩尺寸下MoS2納米帶彈性模量的相對誤差,其表達式為

(17)

表2 單層MoS2彈性模量

圖6 不同尺寸MoS2彈性模量變化曲線

由表2和圖6可以看出,扶手型MoS2和鋸齒型MoS2的彈性模量大小都隨尺寸的增加逐漸減小,并且相鄰兩工況之間的相對誤差也在逐漸減小.工況1兩種類型的MoS2彈性模量均達到了最大值,分別為202.78 GPa和173.07 GPa.工況8兩種MoS2彈性模量值均為8種工況下的最小值,分別為145.03 GPa和141.07 GPa.

扶手型和鋸齒型MoS2的相對誤差值都是逐漸減小的,最小值分別達到了1.08%和0.56%,可以看出,隨著尺寸的增大,尺寸效應逐漸弱化.在工況8中無論是扶手型MoS2還是鋸齒型MoS2的彈性模量均已趨于平穩(wěn).

從整體來看,在每種工況下,扶手型MoS2均比鋸齒型MoS2的彈性模量大,這表明扶手型MoS2抵抗變形的剛度較大,具有更好的抗拉性能.鋸齒型MoS2在每兩種工況之間的相對誤差比扶手型MoS2小.扶手型和鋸齒型MoS2在工況1與工況2之間的相對誤差達到最大,分別為14.42%和9.03%.扶手型和鋸齒型MoS2在工況7和工況8之間的相對誤差達到最小,最小值分別為1.08%和0.56%,表明鋸齒型MoS2的尺寸效應沒有扶手型MoS2明顯.

采用單向拉伸不同工況下MoS2納米帶的彈性模量可以清晰地看到,無論是扶手型還是鋸齒型MoS2的彈性模量都存在明顯的尺寸效應.

3 非線性本構下破壞強度的計算

3.1 破壞應力應變的計算

拉伸過程中一端固定,另一端施加較大的位移荷載,大到Mo—S鍵達到其破壞應力且發(fā)生斷裂,整個MoS2薄片無法再繼續(xù)承受拉力,通過模擬結果計算MoS2薄片的應力和應變.MoS2薄片的應力和應變計算公式為

(18)

(19)

3.2 數(shù)值模擬結果與分析

圖7 MoS2拉伸破壞過程應力應變曲線

扶手型MoS2的破壞應力達到18.9 GPa,此時對應的破壞應變?yōu)?2.5%.鋸齒型MoS2的破壞應力達到16.0 GPa,對應的破壞應變?yōu)?1.2%.

4 結 論

本文采用有限元方法,利用桿系結構模擬共價鍵,分別采用線性本構關系以及非線性本構關系進行單軸拉伸,得出了MoS2納米帶的彈性模量以及破壞應力應變.

1) 采用線性本構關系,單軸拉伸得到MoS2的彈性模量具有尺寸效應,隨著尺寸的增大,彈性模量不斷減小.扶手型MoS2的彈性模量為145.03～202.78 GPa,鋸齒型MoS2的彈性模量為141.07～173.07 GPa.扶手型MoS2的尺寸效應更為明顯.

2) 采用非線性本構關系,單軸拉伸得到扶手型MoS2的破壞應力為18.9 GPa,破壞應變?yōu)?2.5%.鋸齒型MoS2的破壞應力為16.0 GPa,破壞應變?yōu)?1.2%.扶手型MoS2的破壞應力和應變較鋸齒型MoS2更大,力學性能更好.