乙腈在石墨烯极板表面的分子动力学模拟*

姚眉捷，王永朋，梁斐，刘书乐

中山大学材料科学与工程学院，广东 广州 510275

超级电容器，又称为电化学电容器，是一种主要依靠双电层结构或氧化还原赝电容电荷来储存电能的新型储能装置，它的性能介于物理电容器和二次电池之间，兼有物理电容器的功率密度大和二次电池的能量密度高的特点［1-10］。按照反应机理来划分主要分为双电层电容器和法拉第赝电容电容器，其结构主要由电极材料、隔膜和电解质组成［9，11］。对于双电层电容器，它一般通过活性电极表面对电解液中粒子的可逆性物理吸附来储存电荷，进而实现电能的储存。其中，离子液体是常用的电解质［4，12-15］。

离子液体在服役工况温度下具有独特的物理化学性质，如热稳定性好、挥发性较低和工作温度宽等优点，但离子液体高黏度的特点是双电层超级电容器功率性能较差的主要原因［16-17］。通常可以在离子液体中加入有机溶剂改善离子传输，进而提高超级电容器的功率密度［18］。乙腈ACN（Acetonitrile）作为最简单的腈类分子，其具有高介电常数而且可以和水以任意比例互溶，是常用电容器的溶剂，同时也是液相色谱中的流动相溶剂［19-20］。Chaban 等［21］在先前的模拟工作中报道了乙腈的加入能够显著提高咪唑类离子液体的电导率，其原因在于乙腈溶剂稳定了离子对，而更大尺寸的离子簇被破坏。Frolov 等［22］针对碳纳米管与几种不同室温离子液体在乙腈混合物中的相互作用做了大量分子动力学模拟工作，研究中发现离子液体的所有成分都在碳纳米管表面富集，甚至在未充电的碳纳米管表面也形成了几个不同的层。Wang 等［23-24］在乙腈分子和二氧化硅表面所构成的界面体系的模拟中，通过修整二氧化硅表面极性，发现乙腈分子的结构和动力学会随着表面的极性的改变而发生显著变化。

综上，基于乙腈和离子液体混合溶剂作为电解质，使用碳纳米管等作为电极的双电层电容器的模拟研究已经有一部分展开。但是到目前为止，尚未有和乙腈在双电层电容器碳极板表面的结构动力学相关的详尽研究。在服役工况下的双电层电容器中，极板通常带有一定量的电荷，且通常极板上会有微小的杂质［22，25］或多孔结构形成的团簇［6，26］，这些因素对极板表面乙腈结构和动力学的影响值得深入探究。从乙腈在平板石墨烯电极表面的体系出发，本文采用分子动力学模拟手段，依次对极板带电和极板修饰两种条件下乙腈在极板表面的结构和动力学性质进行了模拟。

1 分子动力学模拟

模拟体系由两块石墨烯极板和乙腈分子构成，如图1（a）所示。每块石墨烯极板由4 层碳原子构成，石墨烯极板的尺寸为5.15 nm × 4.25 nm。在最初的体系设定中，两块石墨烯极板均为电中性，它们之间的距离为12 nm。总共有1 787 个乙腈分子限制在两组石墨烯极板之间，这样体系中乙腈的平均密度等于343 K 和一个大气压下乙腈的体相密度。模拟中，在x、y方向施加周期性边界条件；而由于界面的存在，非周期性边界条件施加在z方向。所有的分子动力学模拟都是在温度为343 K 的正则系综（NVT） 下使用Nose-Hoover 热浴［27-28］进行。

图1 平板石墨烯电极-乙腈和团簇修饰石墨烯电极-乙腈体系示意图Fig.1 Schematic diagram of planar graphene-acetonitrile system and cluster modified graphene-acetonitrile system (white,cyan and blue spheres represent hydrogen,carbon and nitrogen atoms)

在模拟的力场参数选择方面，乙腈选用的是Nikitin 和Lyubartsev 开发的全原子模型［29］，该模型能够很好的给出乙腈的密度、扩散系数等物性。石墨烯极板上的碳原子力场采用Amber96 力场中sp2杂化碳原子的参数［30］。计算模拟过程中，体系中非键相互作用的截断半径是1.4 nm。长程的静电相互作用采用PPPM（Particle-particle Particlemesh）［31］算法处理。同时考虑到系统在z方向是非周期性的，在PPPM 计算中我们使用了Berkowitz等［32-33］提出的修正项。全部的分子动力学模拟均使用LAMMPS［34］软件运行。

在电中性平板石墨烯电极的基础上，我们对极板做了两种不同设置以模拟不同条件下乙腈在极板附近的结构以及动力学性质。第一种设置是在每块极板上修饰两个宽度为1 nm 的纳米团簇，团簇由碳原子构成，两个团簇之间的距离为1.22 nm，如图1（b）所示。我们将团簇对应的碳原子层数从一层逐渐增加到五层，这对应于团簇的高度从0.334 nm 到1.67 nm。第二种设置是让石墨烯极板带一定量电荷，这样可以模拟在不同充放电状态下乙腈的结构和动力学性质。在带电体系中，左边的石墨烯极板带正电，右边的石墨烯极板带等量电荷的负电，这样体系整体保持电中性。极板的电荷密度从0 开始，以2.4µC/cm2为间隔逐渐增加到12.0µC/cm2，总共模拟了5个带电体系。

对于所有的模拟体系，我们首先预平衡5 ns（时间步长为1 fs），之后再运行10 ns，并收集轨迹做数据分析，其时间步长为2 fs。

2 分析与讨论

2.1 平板与团簇修饰体系

2.1.1 密度分析我们计算了平板体系和团簇修饰体系的乙腈氮原子的密度（数密度）分布。不同体系氮原子的密度分布，如图2 所示。其中，z= 0 nm对应于石墨烯极板最靠近乙腈的顶层碳原子的位置。对于平板体系和修饰体系，密度分布的第一个峰出现在z= 0.295 nm 左右（图2）。对于平板体系的最靠近极板的第一个峰，其峰值约为体相密度的5倍，这个密度值比较接近之前亲水二氧化硅表面乙腈的密度［23，25］。这说明尽管本研究中使用的石墨烯电极表面并未被任何极性基团修饰，但由于极板与乙腈之间存在范德华相互作用的吸引力，因此乙腈在平板电极表面仍然能够有比较强的吸附。随着乙腈分子逐渐远离极板，极板对于乙腈的吸引逐渐减弱，当z＞3 nm 后，乙腈密度值逐渐趋向于体相密度。当石墨烯极板表面被团簇修饰后，第一个峰的位置虽然没有发生变化，但是峰值较平板体系已经减小了近2/3。很明显，这种现象出现的原因是修饰团簇对乙腈存在空间位阻效应，并且从图2中第一个峰也可以看出，无论修饰团簇高度如何，其修饰体系的第一个峰的数密度都几乎不变，这表明修饰团簇对石墨烯表面附近的空间位阻效应不受其团簇高度的影响。随着乙腈向体相延伸，团簇高度与乙腈的数密度呈现出强烈的同步效应。当修饰高度为一层时，这个体系在z= 0.685 nm 左右（第二个峰）出现最高峰，当修饰高度为二层时，这个体系在z= 1.015 nm 左右（第三个峰）出现最高峰，依次类推。简而言之，修饰团簇层数为n时，它在n+1 个峰会出现最大值（与其他修饰体系和平板体系对比）。

图2 平板电极（黑色曲线）与不同层数团簇修饰电极体系（彩色曲线）的乙腈氮原子数密度Fig.2 ACN nitrogen atom number densities of planar electrode(black curve)and modified electrode with different number of layers of clusters(colored curve)

2.1.2 取向分析为了表征乙腈分子在界面上的分子取向，我们计算了乙腈分子从氮原子指向氰基碳原子的N—C 向量在距离石墨烯极板表面不同距离内的区域的取向分布。因为乙腈分子的N—C向量与其偶极向量相平行，所以可用其来表征乙腈偶极在界面上的排布。图3给出了不同修饰层数下乙腈N—C 向量的取向在不同区域的分布，其中θ是N—C 向量和z轴，也就是垂直于碳极板的法向的夹角。为了更好地分析界面不同区域的取向分布，我们将原来界面区域分为三个层，其z轴坐标范围依次为：第一层，z= 0.296～1.015 nm；第二层，z=1.015～1.775nm；第三层，z= 1.775～2.62 nm。如图3（a）所示，在最靠近石墨烯极板的第一层，几乎所有的取向分布都在三个角度值（cosθ=-1.0，0.1，1.0）附近有明显的分布，这三个角度分别对应于N—C 向量垂直于表面指向碳极板、平行于极板表面（略微倾斜一定角度）、垂直于表面指向乙腈体相。cosθ= ±1.0 处的峰值说明存在排列取向相反的乙腈分子，乙腈可能在碳极板表面形成了类似磷脂双分子层的双层结构［23，25］。

对于平板体系，cosθ= 1.0 处的峰值并不算太高，大多数第一层的乙腈分子的N—C 向量仍然倾向于平行于表面。随着修饰团簇的层数的增加，cosθ= ±1.0 处的峰值出现了明显的增加，说明由于团簇的空间位阻效应，有更多的乙腈分子倾向于垂直于极板排列。

如图3（b）和3（c）所示，随着乙腈分子逐渐远离极板，在第二层和第三层的乙腈分子的取向分布相对变得更均匀。对于一层修饰和二层修饰的体系，乙腈分子在第二层和第三层的取向分布逐渐趋近于一条水平线，也就是各个角度的分布都存在，这非常类似于体相中乙腈分子的取向分布。从图3（c）可以看出，对于五层原子团簇修饰的体系，即便在距离表面较远的第三层，乙腈分子的取向分布仍然显示出强烈的不均匀性，说明修饰团簇扰动了其周围的乙腈分子的排布。

图3 平板（黑色）与团簇修饰体系（彩色）的乙腈N—C向量取向分布Fig.3 ACN N—C vector orientation distributions for planar(black)electrode system and cluster modified electrode system(colored)

2.1.3 扩散系数在体相液体中，分子的扩散一般使用一个扩散系数来表征，即使用经典的Einstein 公式计算均方位移MSD（Mean Squared Displacement）对时间的斜率或者通过Green-Kubo 公式对速度自关联函数进行积分。由于在界面上分子的扩散是各向异性的，因此无法如体相那样用一个单一的扩散系数进行表征，而需要考虑扩散张量的对角元［35-37］。在界面体系中扩散张量仅包含三个对角元Dxx，Dyy，Dzz，其中Dxx=Dyy通常表示为D||，用于表征平行于表面的扩散，而Dzz用于表征垂直于表面的扩散。在本项研究中，我们主要计算了在固液界面中乙腈分子平行于表面的二维扩散系数，可通过计算乙腈分子的二维均方位移得到。二维均方位移的计算公式为

在获得固液界面中乙腈分子的二维均方位移后，我们通过式（2）计算固液界面中乙腈分子的二维扩散系数D||［37-38］，其公式为

图4给出了不同修饰表面体系中平行于表面的扩散系数随碳极板距离变化的趋势。为了更好地得到极板附近的扩散性质，我们分别计算了六个亚层中的扩散系数。亚层的定义如下：第一亚层z= 1.28～1.56 nm；第二亚层z= 1.56～1.95 nm，第三亚层z= 1.95～2.35 nm，第四亚层z= 2.35～2.72 nm，第五亚层z=2.72～3.21nm，第六亚层z=3.21～3.60 nm。这些亚层的定义大体上遵循图2 中各密度峰的范围。如图所示，扩散系数随着修饰层数以及距离碳极板表面的距离都呈单调变化的趋势。首先，对于所有的表面层，随着修饰团簇的尺寸逐步增加，扩散系数呈逐渐下降的趋势，这说明团簇的修饰限制了表面上乙腈分子的活动空间；其次，随着乙腈分子逐渐远离碳极板，扩散系数逐渐增加，说明修饰团簇对于乙腈分子的扩散有阻碍作用。即便对于无修饰的平板表面，其最靠近极板的乙腈分子受到极板的吸引力，其扩散系数也低于远离表面。

图4 平板体系（黑色）与修饰体系（彩色）扩散系数Fig.4 Diffusion coefficient of planar electrode system(black)and cluster modified electrode system(colored)

2.2 平板带电体系

2.2.1 密度分析我们分别对左右极板计算了带电石墨烯极板体系的腈基氮原子密度分布，如图5所示。乙腈的第一个峰出现在z= 0.28 nm，这说明由于极板的范德华力和静电相互作用的吸引，使得乙腈能够在极板上有一定的吸附。值得注意的是，第一个峰在z=0.47 nm 处出现了一个肩峰，这与之前Hu 等［18］研究乙腈在二氧化硅表面的密度的特征类似，说明在石墨烯表面很可能存在两层取向相反的乙腈分子。进一步地，我们观察到第二个峰出现在z=0.73 nm处，随着极板表面正电荷密度的增加，峰值呈现不断扩大的趋势，如图5（a）所示。与此对应的，在带负电荷极板，我们发现了一些更让人出乎意料的结果。如图5（b）所示，对于靠近右侧极板的第一个峰，均出现了双重峰的现象，且当电荷密度0、2.4、4.8µC/cm2时，肩峰均在主峰的左侧；而随着极板电荷量的增加，主峰和肩峰的位置发生了对换。这说明随着极板负电荷的增加，对带负电的乙腈氮原子的排斥增强，使更多乙腈氮原子远离极板表面。一方面，肩峰的出现如前面所述，石墨烯的表面很可能存在两层取向相反的乙腈分子；另一方面肩峰的出现和极板表面的电荷量存在极强的耦合性，且当电荷密度为-9.6 µC/cm2和-12.0 µC/cm2时，发现其左边的肩峰几乎可以忽略不计，我们猜想由于电荷密度的增大所带来的排斥作用，导致乙腈分子无法在z= 0.28 nm 左右处驻留。密度曲线继续向体相延伸，在第二重峰时肩峰特征已经消失，此时随着电荷的增加，峰的位置向体相中间偏移。整体上来看，峰偏移的同时也有峰值增大的趋势。

图5 不同电荷密度下氮原子数密度图Fig.5 ACN nitrogen atom number densities under different charge densities

通过对正负电极板不同电荷密度下表面乙腈密度的分析，发现密度峰的位置和峰高随电荷密度都有一定程度的改变，这说明乙腈的取向随电场强度增加发生了变化。在接下来的分析中，我们将通过取向分析揭示这一变化。

2.2.2 取向分析为了更进一步地揭示带电体系界面处的结构，我们计算了带电体系在极板界面第一层乙腈分子的N—C 向量与垂直于极板表面的z轴夹角的分布，如图6 所示。当极板所带正电荷密度逐渐增加，在cosθ= 1.0 处的概率密度越来越大，这说明乙腈分子倾向于垂直于极板表面排列，同时乙腈的氮原子更加靠近表面。值得我们注意的是，当电荷密度为9.6 µC/cm2和12.0 µC/cm2时，取向分布在cosθ=1.0 近乎是单调上升的趋势。这正好和我们在2.2.1节的推测呼应：电荷密度的增大不仅让乙腈分子取向变化并且也使分子远离极板驻留。其中最直接的原因是乙腈的分子取向发生改变，例如在正极板附近，远离腈基的带正电的一端被极板强烈排斥，导致分子发生旋转，进而取向发生改变。相应地，当极板所带的负电荷密度逐渐增加时，在cosθ= -1.0 处的概率密度越来越大，这说明乙腈分子倾向于垂直于极板表面排列，同时乙腈的带正电的甲基更加靠近表面，能够被带正电的表面吸引。

图6 不同电荷密度下第一层乙腈分子的取向分布Fig.6 Orientation distribution of the first layer of ACN under different charge densities

2.2.3 扩散系数图7 给出了不同带电极板体系中平行于表面的扩散系数在距离碳极板表面不同距离的亚层中的分布图，类似2.1.3节扩散系数的分析，我们在这里计算了最靠近极板表面的前三个亚层的扩散系数。对于所有的体系，随着乙腈逐渐远离表面，扩散系数都出现了一个逐渐增加的趋势。这说明在带电极板上，由于受到极板静电力的吸引，乙腈分子的平移受到了一定程度的限制。另外一个很明显的趋势就是，随着极板电荷密度的增强，乙腈的扩散系数也出现了下降，这同样可以归因为极板吸引力的增强。比较有趣的是，对于带正电荷的极板，当极板从不带电荷增加到电荷密度为2.4µC/cm2时，扩散系数反而出现了一个增加。这与带负电的极板非常不同，乙腈在所有带负电极板表面的扩散系数均低于不带电荷极板的扩散系数。这种现象的可能原因是在该电荷密度的情况下，乙腈的分子取向从原来的腈基氮靠近极板转变成远离极板，此时静电作用强烈排斥乙腈分子，间接导致其扩散速度变快，而对于其他带电体系，乙腈分子取向虽然也发生改变，但是由于距离极板较远，排斥力减小，故扩散速度也相应变慢。

图7 不同电荷密度下乙腈的扩散系数Fig.7 Diffusion coefficient of ACN under different charge densities

3 结 论

在本文中，我们采用分子动力学模拟研究了超级电容器中团簇修饰以及极板电荷密度对极板表面乙腈溶剂的结构和动力学的影响。对于团簇修饰，模拟发现：团簇会对乙腈的吸附起一定的阻碍作用。且随着修饰团簇高度的增加，这种阻碍作用就越强。这种阻碍效应会降低乙腈分子在极板表面的吸附，并使分子取向从平行于极板表面转变为垂直于极板表面。同时，团簇的存在会使界面上分子的扩散动力学减慢。当极板处于带电状态下时，极板电荷密度的增强不仅增强了乙腈分子的吸附，同时使得乙腈分子的取向发生了明显改变，进而影响了乙腈在极板表面的扩散。模拟工作揭示了乙腈在复杂极板表面上的结构和动力学性质与平板表面非常不同，这将有助于未来深入理解乙腈离子液体混合溶液在复杂极板表面上的电化学性能。