基于分子动力学模拟的功能化石墨烯杨氏模量研究

黄立新，邓雅丹 ，周 昆，黄 君

(1. 广西大学土木建筑工程学院, 广西 南宁,530004；2. 广西大学工程防灾与结构安全教育部重点实验室, 广西 南宁，530004)

在采用化学法制备石墨烯材料的过程中，不可避免地导致石墨烯功能化以及造成原子缺失，官能团与缺陷对石墨烯性能的影响较为复杂，近年来，研究者针对功能化石墨烯的力学性能开展了大量研究工作，如 Zheng等[4]研究了功能化对石墨烯杨氏模量、剪切模量及起皱性能的影响，发现石墨烯杨氏模量强烈依赖于其官能团引入程度及分子结构，功能化造成石墨烯剪切模量和临界起皱应变降低；Zhang等[5]研究了官能团对GO/聚合物界面结构、动力学和力学性能的影响，结果表明，氢键密度的增加有助于提高界面强度、强化界面结构并降低破坏应变；Yang等[6]基于分子动力学模拟及有限元分析，研究了氧官能团对单层石墨烯/聚乙烯复合材料的影响；Zaminpayma等[7]借助分子动力学模拟，研究了多种功能化石墨烯纳米带的力学和电学性质，结果表明，引入官能团的石墨烯纳米带杨氏模量随温度的升高而降低，同时官能团的存在还造成石墨烯纳米带电阻减小、电流强度增大；Melro等[8]研究了石墨烯与环氧树脂的相互作用，计算了功能化石墨烯层和环氧复合材料界面的剪切强度、拔出力并测定了杨氏模量，结果表明，功能化对复合材料机械性能不利，但能增加石墨烯和环氧树脂之间的相互作用；Haghighi等[9]研究了缺陷和官能团对碳纳米管和石墨烯增强环氧纳米复合材料力学性能的影响,结果表明，缺陷和官能团降低了纳米填料和纳米复合材料的弹性模量，二者对环氧树脂纳米复合材料力学性能的影响程度还与填料尺寸有关。官能团、缺陷对功能化石墨烯力学性能的影响机制复杂，尚有许多问题需要解决，因此，本文基于分子动力学模拟，采用分子模拟软件建立石墨烯模型，重点探究了羧基、羟基含氧官能团以及缺陷对功能化石墨烯杨氏模量的影响规律，以期为功能化石墨烯材料弹性性能的研究提供参考。

1 模拟方法及模型的建立

1.1 分子动力学模拟

基于分子动力学原理，利用Materials Studio(MS)分子模拟软件对材料开展弹性性能研究，该软件可用于创建分子结构模型并进行相关计算分析，采用COMPASS力场[10]模拟原子力场，COMPASS力场的vander Waals相互作用势为Lennard-Jone-9-6势函数[11]，该力场的总能量Utotal可以表示为：

Utotal=Udiagonal+Ucross-coupling+Uij+Uelec

(1)

式中，Udiagonal为对角的交叉偶合项，Ucross-coupling为非对角的交叉偶合项，Uij为vander Waals能，Uelec为库伦能，Uij和Uelec的计算公式分别为：

(2)

(3)

(4)

(5)

在进行分子动力学模拟时，设应力和应变分别为σ和ε，材料杨氏模量为E，则有:

σ=Eε

(6)

利用最小二乘法对小变形应力和应变样本进行直线拟合[12]，得到杨氏模量E的计算公式[13]为：

(7)

1.2 模型的建立

1.2.1 石墨烯模型

2018年10月22日，省委常委会审议通过《关于支持检察机关公益诉讼工作的意见》。省委书记娄勤俭强调，检察机关公益诉讼工作要坚持司法为民，要紧盯生态环境、资源保护等社会关切、群众关注的重点领域，加大公益诉讼办案力度，推动有关部门及时履职、解决问题，切实增强人民群众的获得感、幸福感、安全感，让公益诉讼真正体现公益、惠及民生。这为检察机关开展公益诉讼工作进一步提供了制度保障，并指明了工作方向。

图1所示石墨烯晶胞结构，其平衡键长为0.1421 nm[14]。利用MS分子动力学模拟软件导入石墨烯晶胞建立相应石墨烯模型，然后进行分子动力学模拟计算，具体操作步骤如下：

(1)在导入石墨烯晶胞的基础上，创建水平方向尺度Lx和竖直方向尺度Ly分别为6.04730、6.28454 nm的石墨烯超晶胞，采用周期性边界条件消除边界效应对模型性能带来的影响，再对石墨烯超晶胞进行切片，把所得石墨烯片放入空胞中，从而构建出石墨烯模型(图2)；

(2)借助MS软件中Forcite模块对石墨烯晶胞结构进行几何优化，使结构总能量最小，即对应其结构的势能面能量最低；

(3)在COMPASS力场中应用NPT系综，采用压力模式，借助恒压器Souza-Martins对模型进行调节，并设置时间常数为10 ps，恒温器采用NHL，时间步长设置为1 fs；

(4)将零应力设置为 5个平衡周期，并且每个周期的时间步长设置为100。晶胞的大小和形状均在每个循环开始时重置，以达到上个周期的平均规模；

(5)分别在X、Y方向对石墨烯晶胞施加应力，其它方向为零应力。所施加的初始应力为σ0，然后以增量Δσ施加一系列应力，并将每个应力值的平衡步数设置为500，每500 步生成一帧，每一帧统计一次应力、应变值，在此过程中借助恒压器Souza-Martins计算应力、应变值并绘制相应关系曲线。

图1 石墨烯晶胞模型

图2 石墨烯模型

1.2.2 功能化石墨烯模型

在石墨烯晶胞模型基础上，利用MS软件脚本Script进行石墨烯功能化的模拟计算。官能团随机引入石墨烯表面，设置1个官能团与1个碳原子相连，且2个官能团不能与1对相邻的碳原子相连。图3所示分别为羧基和羟基官能团接枝在石墨烯上，其中红球为氧原子、白球为氢原子。本研究中以Script脚本中羟基、羧基官能团的个数n与石墨烯碳原子总数N之比来衡量石墨烯的功能化程度w，即：

(8)

(a)羧基 (b)羟基

1.2.3 含缺陷石墨烯模型

空位缺陷是石墨烯中最常见的结构缺陷[15]，在此设计具有不同形状空位缺陷的石墨烯模型如图4所示，其中图4(a)～4(c)分别为含横向长条形、纵向长条形以及孔洞形空位缺陷的石墨烯模型，图4(d)是含有12个随机空位缺陷的石墨烯模型。本研究中利用碳原子缺失率Q来衡量缺陷对石墨烯杨氏模量的影响，碳原子缺失率是缺失碳原子数m与理想石墨烯碳原子总数M之比，即：

(9)

(a)横向长条形 (b)纵向长条形

(c)孔洞形 (d)随机空位缺陷

2 模拟与分析

2.1 理想石墨烯的杨氏模量

对图2所示石墨烯模型施加X方向的应力，设置初始应力值σ0为10 GPa，依次递增增量Δσ为20 GPa，根据MS软件模拟计算所得应力、应变值作出应力-应变关系曲线如图5所示。此外，将模拟所得应力、应变值带入式(7)，即可求得石墨烯杨氏模量E为1.158 TPa，与已有文献所报道的相应值(1.04～1.24 TPa)[16-19]基本吻合，表明本研究石墨烯模型建构正确。

图5 应力与应变的关系

2.2 官能团对功能化石墨烯键长及杨氏模量的影响

通过分子动力学模拟，含羧基或羟基石墨烯的键长测量结果如图(6)所示。由图6可知，官能团的引入改变了石墨烯的键长，在含羧基石墨烯中，羧基接枝处附近原子键O1A1、O1B1和O1C1的键长依次为0.1523、0.1545、0.1712 nm(图6(a))，相比理想石墨烯平衡键长(0.1421 nm)分别增加了7.2%、8.7%和20.5%，键长变化十分明显，这是因为与被接枝碳原子相邻的sp2杂化轨道转换成sp3杂化，同时该碳原子被拉离平面，从而导致周围键长发生改变；羟基接枝在石墨烯表面时的状况与羧基相似，但状态平衡后被接枝碳原子的外移程度不大，因此对键长的影响也相对较小，含羟基石墨烯中的最大键长相比理想石墨烯平衡键长的增幅也仅有5.8%(图6(b))。

(a)含羧基石墨烯 (b)含羟基石墨烯

功能化程度分别为0.67%、3.33%、6.67%的含羧基或羟基官能团石墨烯分子动力学模拟结果如图7所示。由图7可以看出，随着功能化程度的提高，石墨烯结构破坏程度也不断加剧，其中功能化程度为0.67%的石墨烯结构仅在接枝官能团处发生形变，其它部分则较为完整(图7(a)、图7(d))，而功能化程度为6.67%的石墨烯表面十分崎岖，结构明显变形(图7(c)、图7(f))，并且羧基官能团的引入对石墨烯结构的影响较羟基官能团更大。

(a)0.67%羧基 (b)3.33%羧基 (c)6.67%羧基

(d)0.67%羟基 (e)3.33%羟基 (f)6.67%羟基

通过分子动力学模拟，获得不同手性石墨烯杨氏模量与石墨烯功能化程度的关系如图8所示。由图8可见，随着石墨烯功能化程度不断提高，不同类型石墨烯的杨氏模量均逐渐降低，并且相比含羟基石墨烯，含羧基石墨烯的杨氏模量降速更快，即羧基量变化对石墨烯弹性性能产生的影响更大，此外，对于官能团相同的石墨烯，其手性为锯齿型时的杨氏模量降速较手性为扶手椅型时更慢，即扶手椅型石墨烯的杨氏模量受石墨烯功能化程度变化的影响更明显。

(a) 含羧基石墨烯 (b) 含羟基石墨烯

在采用化学法制备石墨烯时，所制石墨烯中往往同时存在多种官能团，为此，本研究利用分子动力学模拟，考察了功能化石墨烯中羧基和羟基的含量变化对石墨烯杨氏模量的影响，结果如图9所示，图中石墨烯均为扶手椅型。由图9可知，石墨烯杨氏模量随其功能化程度的增大而不断减小，当石墨烯功能化程度相同时，只含羟基的石墨烯杨氏模量最大，只含羧基的石墨烯杨氏模量最小，同时含有羟基和羧基的石墨烯杨氏模量则介于上述二者之间，并且随所含羟基比例的增加而增大，随所含羧基比例的增大而减小，该结果再次表明羧基量变化对石墨烯弹性性能的影响更大。需要指出的是，对锯齿型石墨烯进行分子动力学模拟也获得了相同的结果，此处不再赘述。

(a) 羧基与羟基含量之比为1∶1 (b) 羧基与羟基含量之比为1∶2或2∶1

(c) 羧基与羟基含量之比为1∶3或3∶1 (d) 羧基与羟基含量之比为2∶3或3∶2

2.3 空位缺陷对功能化石墨烯杨氏模量的影响

借助MS软件对羧基或羟基功能化程度为3%且含有空位缺陷的扶手椅型石墨烯模型进行分子动力学模拟，施加应力沿X方向，所得结果分别如图10和图11所示。从图10和图11中可以看出，所有石墨烯表面官能团和缺陷位置附近都出现一定程度的变形，并且在相同模拟条件下，由纵向长条形空位导致的石墨烯结构形变最明显，这应归因于在X方向的应力作用下，纵向缺陷更易使模型发生变形，而其它3种空位缺陷在该方向应力作用下稳定性相对较好，造成的结构变形也较小。功能化扶手椅型及锯齿型石墨烯的杨氏模量与石墨烯中碳原子缺失率的关系分别如图12和图13所示。从图12和图13中可以看出，随着碳原子缺失率的增加，功能化石墨烯的杨氏模量均不断减小；在碳原子缺失率不变的条件下，含相同官能团的石墨烯杨氏模量随其功能化程度的增加而减小，并且含羧基石墨烯的杨氏模量均小于含等量羟基的石墨烯相应值。由此可见，空位缺陷及官能团的存在均导致石墨烯杨氏模量较理论值有所降低，同时，二者降低石墨烯杨氏模量的作用具有叠加效应。图12和图13中的受力方向分别为X和Y方向，分别对比分析图12(a)和图12(b)、图13(a)和图13(b)发现，图12(a)和图13(b)中石墨烯杨氏模量随原子缺失率变化的幅度相对较小，这表明，当受力方向与缺陷分布方向相同时，石墨烯结构较稳定；当受力方向与缺陷分布方向垂直时，石墨烯结构更易被破坏。

(a)横向长条形空位 (b)纵向长条形空位 (c)孔洞形空位 (d)随机空位

(a)横向长条形空位 (b)纵向长条形空位 (c)孔洞形空位 (d)随机空位

(a) 横向长条形空位 (b) 纵向长条形空位

(c) 孔洞形空位 (d) 随机空位

(a) 横向长条形空位 (b) 纵向长条形空位

(c) 孔洞形空位 (d) 随机空位

3 结论

(1)随着石墨烯功能化程度的增加，石墨烯的杨氏模量逐渐减小且羧基对石墨烯杨氏模量的影响比羟基大；当石墨烯功能化程度相同时，同时含有羟基和羧基的石墨烯的杨氏模量随所含羟基比例的增加而增大，随所含羧基比例的增大而减小。

(2)空位缺陷及官能团的存在均导致石墨烯杨氏模量较理论值有所降低，二者降低石墨烯杨氏模量的作用具有叠加效应；当受力方向与缺陷分布方向相同时，石墨烯结构较稳定；当受力方向与缺陷分布方向垂直时，石墨烯结构则易被破坏。