單層石墨烯力學性能尺寸效應的有限元分析

韋 琦, 黃 君, 黃立新,2*, 吳 宇

(1. 廣西大學土木建筑工程學院, 南寧 530004； 2. 廣西大學工程防災與結構安全教育部重點實驗室, 南寧 530004)

石墨烯是單層碳原子面材料，其碳原子間相互圍成蜂窩形的二維網狀結構，被認為是一維結構的碳納米管和三維結構的石墨片的最基本組成單元。石墨烯是已知的材料中最輕、最薄、最堅硬的材料，在力學、電學和熱力學等方面都具有極其優(yōu)異和獨特的性質[1-4]?；诘谝恍栽韀5]、分子動力學[6-8]和結構力學[9-11]等方法得到石墨烯的彈性模量約為1 TPa。石墨烯優(yōu)異的力學性能使其成為納米復合材料的理想納米增強材料。

研究[12-16]表明，石墨烯材料增強的納米復合材料的力學性能增強效果非常明顯。隨著石墨烯納米復合材料市場需求的迅速增長，石墨烯性能參數(shù)的研究是必不可少的。石墨烯在不同尺度下往往會出現(xiàn)不同的力學特性，表現(xiàn)出明顯的尺寸效應，這是石墨烯納米復合材料性能增強的影響因素之一。Xu等[17]研究了石墨烯與聚對苯二甲酸乙二醇酯基之間切向界面尺寸產生的力學性能的變化和邊界效應。通過原位拉曼光譜測定得出界面剪切應力具有尺寸效應，它的值隨石墨烯尺寸的增大而顯著減小，石墨烯達到宏觀尺寸時趨于定值。在力學性能方面也存在著尺寸效應。Chu等[18]采用了分子動力學模擬的方法，利用單軸拉伸來研究有限尺寸下石墨烯納米帶的力學性能尺寸效應。研究發(fā)現(xiàn)，鋸齒型石墨烯納米帶的尺寸效應顯著而扶手型石墨烯納米帶尺寸效應可以忽略不計。值得注意的是，用有限元法、分子靜力學方法、分子動力學方法對石墨烯尺寸效應的研究是無法比較的。最后，結果表明，方形鋸齒型石墨烯尺寸增大，彈性模量增大，極限破壞應力和應變減小。對于長方形鋸齒型石墨烯，寬度增加，彈性模量增大，極限破壞應力和應變減小。

采用桿件來模擬石墨烯之間碳-碳共價鍵鍵，桿件之間的連接為剛性連接。由此得到的桿件系統(tǒng)模擬單層石墨烯結構?，F(xiàn)基于桿件系統(tǒng)的有限元方法，研究單層石墨烯的尺寸效應。首先構建8個不同尺寸的有限元計算模型，然后分析單層石墨烯力學性能隨尺寸變化的情況。通過數(shù)值算例分析單層石墨烯力學性能的尺寸效應以及尺寸效應弱化的現(xiàn)象。

1 石墨烯有限元模型

1.1 石墨烯結構

石墨烯晶體屬于復式晶格，單層石墨烯呈片狀蜂窩型，如圖1所示。兩個碳原子組成一個元胞，每個碳原子會與3個相鄰的碳原子形成共價鍵σ鍵，這是自然界中最強的化學鍵。研究表明，配位數(shù)為3的石墨烯兩碳原子之間的距離為1.42110-10m。此外，每個鍵角為120°，因此幾何形狀上為正六邊形。

圖1 單層石墨烯結構示意圖Fig.1 Schematic diagram of monolayer graphene structure

1.2 桿件系統(tǒng)的有限元計算模型

將石墨烯的碳碳共價鍵理解為桿件，碳原子理解為桿件之間的連接點，桿件的變形可用來反映實際共價鍵的拉伸、壓縮以及面內鍵角的變化。在此基礎上，得到桿件系統(tǒng)，如圖2所示。用等效的桿件系統(tǒng)描述石墨烯結構，把石墨烯結構的受力行為用桿件系統(tǒng)進行模擬?；谧钚菽茉?，推導出桿件系統(tǒng)梁單元的剛度矩陣，由此得出有限元求解的基本步驟[19]。在大型有限元軟件ABAQUS平臺上，對桿件系統(tǒng)進行網格劃分，定義梁單元的材料屬性和尺寸參數(shù)，可以實現(xiàn)對石墨烯結構的力學性能分析。

圖2 單層石墨烯有限元模型示意圖Fig.2 Schematic diagram of finite element model of monolayer graphene

1.3 等效梁單元參數(shù)

在分子力學中，分子的總能量由分子動能及勢能疊加而成。分子勢能其表達式為

U=Unb+Ur+Uθ+Uφ+Uτ+Uel

(1)

式(1)中：Unb為范德華勢；Ur為共價鍵的拉伸勢；Uθ為共價鍵的鍵角彎曲勢；Uφ為二面角的扭轉勢；Uτ為平面震動勢；Uel為庫侖靜電勢。而在勢能計算中，因為非鍵勢的作用相對于其他勢能影響較小，因此有時可以忽略其作用，Cornell等[20]得出石墨烯的拉伸勢與鍵角變化勢等勢能函數(shù)。為了將梁單元等效碳-碳連接鍵，Li等[11]做出假定：用圓形截面梁單元模擬碳-碳鍵，從而得到梁單元的應變能函數(shù)。根據(jù)分子勢能函數(shù)以及結構應變能函數(shù)，可以推導出等效梁單元的抗拉剛度EA、抗彎剛度EI以及抗扭剛度GJ。圖3為圓形截面梁單元示意圖。

F為受力，下標x、y、z為受力方向；M為彎矩， 下標x、y、z為彎矩作用方向圖3 圓形截面梁單元示意圖Fig.3 Schematic diagram of beam element with circular section

在自定義圓形截面梁時，根據(jù)在分子模型下的能量變化的拉伸常數(shù)kγ，鍵角彎曲常數(shù)kθ，抗扭常數(shù)kτ和經典力學中能量變化的等價關系，由式(2)、式(3)和式(4)確定單元的相關參數(shù),包括彈性模量E，泊松比μ，直徑d。

(2)

(3)

(4)

Li等[11]得出參數(shù)kγ/2=469 kcal/(mol·?2),kθ/2=63 kcal/(mol·rad2),kτ/2=20 kcal/(mol·rad2)，其中，1 cal=4.186 J,從而得出彈性模量E，泊松比μ，直徑d。Li等指出式(4)經過計算無法得出使用的泊松比，但因為kτ對模型的彈性模量影響較小，因此歐拉梁的μ取0。

2 有限元建模

2.1 自定義梁單元建模

碳原子相對位置變化引起的勢能等于空間梁單元變形引起的勢能，由此可以導出分子力學的力常數(shù)與結構構件單元剛度系數(shù)的聯(lián)系。假設桿件的截面為圓截面，桿件之間為剛性連接，如圖4所示。

基于有限元法導出梁單元的抗拉剛度矩陣，抗彎剛度矩陣以及抗扭剛度矩陣。根據(jù)式(2)、式(3)和式(4)，從而得出梁單元的材料參數(shù)E、μ以及幾何參數(shù)d，有關參數(shù)值如表1所示。在得到一系列的桿件參數(shù)后，采用Fortran語言編寫有限元軟件ABAQUS的圓截面梁單元的UEL子程序,以備后續(xù)中導入自定義單元進行建模。

在編輯UEL子程序時，對圓截面梁單元的結點坐標、自由度以及單元剛度進行定義。因為最后是通過有限元軟件ABAQUS進行相應的力學模型分析，因此需要完成UEL子程序與ABAQUS語言的對接。

圖4 石墨烯有限元建模結點連接Fig.4 Node connection of finite element modeling of graphene

表1 梁單元材料參數(shù)和幾何參數(shù)

2.2 單層石墨烯模型建立

單層石墨烯呈正六邊形集合的蜂窩型結構，鍵角均為120°。通過AutoCAD輔助設計軟件建立單層石墨烯模型，設置正六邊形的邊長為0.142 1，這是單層石墨烯碳-碳鍵的長度，即0.142 1 nm。

在CAD建立完單層石墨烯的幾何模型后，導入ABAQUS有限元計算軟件進行相應的前處理。在前處理中，需要對單元進行截面屬性和材料屬性的定義，而這些工作在UEL子程序已經考慮。在劃分網格時，桿件的尺寸設置為0.142 1，這與碳-碳鍵的尺寸一致。對于載荷以及邊界條件，均在inp文件里進行相應的設置，并輸入單元編號和節(jié)點編號。在ABAQUS中建立作業(yè)后，通過編寫好的inp文件進行計算得出結果，并在ABAQUS中對計算結果進行后處理。

在建立模型時，假設單層石墨烯(正六邊形結構)無團聚，無褶皺。如圖5所示，單層石墨烯的長和寬分別為l和h。為了研究單層石墨烯力學性能的尺寸效應，建立有限元計算模型時，考慮模型尺寸的8個工況，表2是8個工況的模型尺寸參數(shù)。

圖5 單層石墨烯計算模型工況Fig.5 Monolayer graphene computational model

表2 8個工況的尺寸參數(shù)

如圖6所示，考慮扶手型與鋸齒型的情況，單層石墨烯的扶手型與鋸齒型的規(guī)定是根據(jù)石墨烯的邊界來確定的。在固定好坐標軸后，垂直于x軸的邊界形狀為鋸齒形狀，則研究該x軸方向性質的模型表述為鋸齒型石墨烯。垂直于x軸的邊界形狀為扶手型，則研究該x軸方向性質的模型表述為扶手型石墨烯。扶手型與鋸齒型的邊界情況并不相同，從而在inp文件中編輯的扶手型與鋸齒型的邊界條件、位移條件等也不同。

圖6 單層石墨烯片的手性定義Fig.6 Definition of chirality of monolayer graphene

在單層石墨烯模型中，碳-碳共價鍵理解為圓截面梁，碳原子分布在各個結點上。因此忽略碳原子的大小，僅以圓截面梁剛性連接。在對單層石墨烯進行尺寸效應的分析中，將會對不同類型的單層石墨烯從幾納米到幾十納米情況下進行研究，并且石墨烯的厚度取0.34 nm[21]。單層石墨烯模型理論上為正方形薄片，但因其扶手型和鋸齒型的尺寸并不相等，且不一定恰好等于整數(shù)。因此，未得到整數(shù)尺寸下的模型力學性能。在這種情況下，取能完整顯示六邊形的模型，并取最接近該尺寸的模型來研究。同時考慮石墨烯在鋸齒型方向與扶手型方向上力學性質的不同，從而全面地對單層石墨烯的尺寸效應進行分析。

2.3 力學性能分析

在小變形的基礎上，通過有限元計算，將研究單層石墨烯的彈性模量、剪切模量和泊松比。在單層石墨烯的邊界設置約束，使其位移在某些方向為零。同時在一條邊界施加適當?shù)男∽冃?。根?jù)有限元計算求得單層石墨烯的支座反力，從而求得石墨烯的材料參數(shù)。

如圖7所示，針對鋸齒型和扶手型石墨烯，在相應邊上施加拉伸小變形。由胡克定律得到應變，即

(5)

應力σ和應變ε分別為

(6)

(7)

式中：Fy為y方向支座反力的合力(鋸齒型)；Fx為x方向支座反力的合力(扶手型)；Lx和Ly分別為石墨烯x和y方向的初始長度；t為石墨烯的厚度；ΔLy為石墨烯y方向的變形(鋸齒型)；ΔLx為石墨烯x方向的變形(扶手型)。

將式(6)和式(7)代入式(5),求得彈性模量為

(8)

根據(jù)橫向線應變與縱向線應變的比值求得泊松比為

(9)

圖7 單層石墨烯受拉伸變形Fig.7 Monolayer graphene under tensile deformation

如圖8所示，針對鋸齒型和扶手型石墨烯，在相應邊上施加剪切小變形。由剪切胡克定律求得剪切模量為

(10)

剪應力和剪應變分別為

(11)

(12)

式中：Fx為x方向支座反力的合力(鋸齒型)；Fy為y方向支座反力的合力(扶手型)；Lx和Ly分別為石墨烯x和y方向的初始長度；t為石墨烯的厚度；ΔLx為石墨烯x方向的變形(鋸齒型)；ΔLy為石墨烯y方向的變形(扶手型)。式(11)和式(12)代入式(10),求得剪切模量為

(13)

3 結果與分析

在單層石墨烯的邊界設置約束，針對鋸齒型和扶手型石墨烯，在相應邊上施加小變形。借助于大型有限元計算軟件ABAQUS，計算單層石墨烯的變形和支座反力，根據(jù)上述公式計算出單層石墨烯的性能參數(shù)。

有限元計算采用以下數(shù)據(jù)和條件：單層石墨烯的厚度t=0.34 nm；單層石墨烯的長和寬尺寸考慮8個工況，如表2所示。圖7和圖8所示的計算模型中，8個工況的Lx和Ly分別取表2的l和h的尺寸。圖7(a)鋸齒型石墨烯的y方向的拉伸變形量ΔLy=0.001 nm;
圖7(b)扶手型石墨烯的x方向的拉伸變形量ΔLx=0.001 nm；圖8(a)鋸齒型石墨烯的x方向的剪切變形量ΔLx=0.001 nm;
圖8(b)扶手型石墨烯的y方向的剪切變形量ΔLy=0.001 nm。統(tǒng)一采用國際單位，根據(jù)公式計算得到的彈性模量和剪切模量的單位為Pa，然后換算成GPa。

3.1 彈性模量

為了考察彈性模量隨尺寸變化而變化，定義如下量：

(14)

基于圖7所示模型的有限元軟件ABAQUS計算，分別得到y(tǒng)方向支座反力的合力Fy(鋸齒型)和x方向支座反力的合力Fx(扶手型)，代入式(8)得出彈性模量，其結果如表3所示。

表3 單層石墨烯彈性模量

3.2 泊松比

定義

(15)

基于圖7所示模型的有限元軟件ABAQUS計算，分別得到鋸齒型石墨烯x方向的變形ΔLx和扶手型石墨烯y方向的變形ΔLy，代入式(9)得出泊松比，其結果如表4和圖10所示。

表4 單層石墨烯泊松比

圖10 單層石墨烯泊松比Fig.10 Poisson’s ratio of monolayer graphene

3.3 剪切模量

定義

(16)

基于圖8所示模型的有限元軟件ABAQUS計算，分別得到x方向支座反力的合力Fx(鋸齒型)和y方向支座反力的合力Fy(扶手型)，代入式(13)得出剪切模量，其結果如表5和圖11所示。

表5 單層石墨烯剪切模量

圖11 單層石墨烯剪切模量Fig.11 Shear modulus of monolayer graphene

無論是彈性模量、剪切模量，還是泊松比，在小尺度下都會存在顯著的尺寸效應，但尺寸大于一定值時，尺寸效應消失，數(shù)值將會收斂到一定的值。為了進一步研究單層石墨烯力學性能的尺寸效應，基于圖7所示的模型，尺寸從大約1 nm變化到80 nm，根據(jù)有限元計算得到彈性模量隨尺寸變化的情況，如圖12所示。從圖12中可以看出，小尺度下尺寸效應十分明顯，并且扶手型和鋸齒型在小尺寸下的彈性模量相差較大。隨著單層石墨烯尺寸逐漸增大，扶手型和鋸齒型單層石墨烯的彈性模量逐漸下降并且趨近，最后趨近于同一個值。

圖12 單層石墨烯彈性模量變化曲線Fig.12 Modulus of elasticity curve of monolayer graphene

4 結論

在片狀蜂窩型的單層石墨烯中，將碳-碳共價鍵視為桿件，碳原子視為桿件之間的連接點，由此得到桿件系統(tǒng)，其變形可用來反映實際共價鍵的拉伸、壓縮、以及面內鍵角的變化。基于桿件系統(tǒng)的有限元方法，計算出單層石墨烯的彈性模量、剪切模量和泊松比等力學性能參數(shù)。針對有限元計算模型，構建了8個從小到大的尺寸工況，求出單層石墨烯的力學性能參數(shù)，然后研究這些參數(shù)的尺寸效應。通過有限元計算分析，得到以下結論。

(1)小尺寸(0.923 nm×0.984 nm)的計算模型，單層石墨烯的彈性模量和剪切模量均達到最大值。當計算模型尺寸增大到1.969 0 nm×1.989 4 nm時，彈性模量和剪切模量均明顯下降?？梢姡谛〕叽绲那闆r下，單層石墨烯的力學參數(shù)變化幅度大，尺寸效應很顯著。

(2)隨著計算模型尺寸的增大，單層石墨烯的彈性模量和剪切模量逐漸下降，并且下降幅度逐漸減緩，最終彈性模量和剪切模量分別趨近于相應的值，尺寸效應出現(xiàn)弱化；而單層石墨烯的泊松比的數(shù)值迅速趨于某一區(qū)間范圍之內，并且出現(xiàn)上下波動的變化。

(3)相同尺寸的情況下，扶手型單層石墨烯的彈性模量和剪切模量總是大于鋸齒型單層石墨烯的相應量。除了小尺寸(0.923 nm×0.984 nm)的一個工況外，扶手型單層石墨烯泊松比的值大于鋸齒型的值。