受限在纳米水环境空间中的甲烷分子

冯 梅，叶 超，陈艳燕，周晓艳，吴锋民，陆杭军

(浙江师范大学凝聚态物理研究所，浙江金华 321004)

0 引 言

能源是人类社会赖以生存与发展的物质基础.当今，煤、石油、天然气的储量正在日趋减少，为了防止能源枯竭，科学家们正通过多方面开发高效、清洁的新能源［1-5］.因此，早在20世纪60至90年代，世界主要发达国家就已相继开始了对甲烷水合物(又称可燃冰)的勘探与开发研究.作为非常规清洁能源的重要组成部分，甲烷水合物已经受到了世界各国的重点关注，欧美日韩等国纷纷展开对甲烷水合物的勘探开发工作，我国也已在甲烷水合物领域付出了巨大的努力，未来可燃冰有望成为天然气的“接班人”.虽然甲烷水合物的发现已经有200多年的历史，但是直到最近50年，科学家们才从微观层次上去深入研究其形成机理［6］.由于甲烷水合物的凝结核非常小，而现在的实验手段还无法动态跟踪原子的运动轨迹，因此很难通过实验在原子层次上研究甲烷水合物的形成机制与动力学行为［7］.当前，行之有效的方法是通过分子动力学模拟，模拟甲烷水合物的形成和演化，阐明其中的物理机制［8-9］.另外，利用热激化法开采可燃冰时，由于海底的多孔介质不是集中为“一片”，也不是一大块岩石，而是较为均匀地遍布着，很难收集甲烷气体，所以甲烷在小空间的水环境中的一些性质研究十分重要.

2008年，Hawtin等［10］进行了一系列的甲烷水合物的模拟，初步得到了其形成机理.2009年，Walsh等［11］模拟了甲烷水合物在μs时间尺度下的自然成核过程，并发现在一定的温度和压强条件下水分子在甲烷周围通过氢键形成的“包络”将其包围，形成了类冰状的结晶物质.但是，当甲烷分子与水分子受限在纳米空间时，它们又会表现出怎样的性质呢?2001年美国科学家Hummer等［12］通过计算机模拟发现，直径只有0.81 nm的碳纳米管能够充水，而且发现在碳纳米管中与水分子偶极矩方向形成一致，成为一维水分子链，水分子链偶极矩方向在热扰动作用下能够不断翻转，有趣的是有序水分子链能够步调一致地快速通过碳纳米管.Koga等［13］发现在直径为1 nm左右的一系列碳纳米管中，高压能使水发生相变.近年来，许多实验与理论研究结果都表明，水受限在纳米空间中，它的结构、性质与动力学行为都与宏观的水有很大的不同［14-20］.在国内，上海应用物理研究所的方海平研究组也在这方面做了一系列的深入研究［21-24］.

本文将利用分子动力学模拟，进一步研究甲烷与水受限在碳纳米管中的一些性质与动力学行为.

1 模型与方法

模拟系统如图1所示.图1(a)为模拟系统初始状态的径向截面图;图1(b)为模拟系统初始状态的轴向截面图.首先将一定数量的甲烷分子以随机的方式放入碳纳米管中，再在碳纳米管内部充满水.笔者构建了3个不同甲烷浓度的模拟系统:充入100个甲烷分子的系统称为系统1;充入200个甲烷分子的系统称为系统2;充入300个甲烷分子的系统称为系统3.

图1 初始模拟系统的径向剖面和轴向剖面示意图

本文采用GROMACS 4.0.7软件进行分子动力学模拟，碳纳米管半径R=2.67 nm，管长z=7.40 nm，在Z方向设置周期性的边界条件.水分子采用tip4p模型［25］.碳原子之间的莱纳德-琼斯势参数为εcc=0.439 kJ/mol.设定每一个步长为2.0×10-6ns，模拟时间为2.0 ns.静电相互作用的计算采用了Particle Mesh Ewald(PME)方法［26］.温度为300 K，并使用Nose-Hoover控温算法.其他相关参数都按照文献［27-28］中的设置.

图2是模拟系统的的剖面图，为了研究甲烷分子与水分子的轴向扩散系数，建立了以碳纳米管中心轴为z轴的柱坐标.与笔者以前研究工作中的方法［20］类似，将碳纳米管空间划分为4个子空间，半径依次为 r1=0.61 nm;r2=1.55 nm;r3=2.14 nm;r4=2.67 nm，并从内到外依次编号 1，2，3，4 层.根据分子的坐标位置，可以判断分子属于哪一层.假设该分子到碳纳米管轴中心的距离为r，当r＜r1时，分子属于第1层;当r1＜r＜r2时，分子属于第2层;当r2＜r＜r3时，分子属于第3层;当r3＜r＜r4时，分子属于第4层.

2 结果与分析

2.1 径向密度分布

径向分布函数是可观测的量，是表征流体局部结构的函数，可以用来研究物质的有序性，在液体态理论中有重要作用.如图3所示，笔者分别计算了水与甲烷径向的密度分布.本文中，密度分布ρ(r)是表示在距离碳纳米管中心轴距离为r处的密度分布，ρ0表示在中心轴(r=0)的密度分布，ρ(r)/ρ0为r处密度与中心轴处的密度比.

图2 径向截面分割图

图3 水分子的相对径向密度分布

图4 甲烷分子的相对径向密度分布(插图为相对径向密度分布的全图)

从图3可以明显地看出，管内的水分子相对密度分布在2.3 nm附近有一个极大值，随后迅速减小，在大于2.4 nm处水分子密度分布基本为0.这主要是由于碳纳米管壁的范德华势对水分子分布的影响.

从图4发现，管内的甲烷相对密度分布在2.3 nm附近也有一个峰值，同样在2.4 nm附近减小到0.分子动力学模拟结果说明甲烷分子在充满水的碳纳米管管壁附近会出现聚集现象.其中，模拟系统1中的甲烷分子相对密度分布在0～2.0 nm处保持稳定，而系统2和系统3在0～2.0 nm范围内的甲烷分子相对密度分布呈下降趋势.这是由于系统1中甲烷数量较少，在整个模拟过程中都富集到管壁附近.而在系统2和系统3中，由于管内甲烷浓度较大，一部分甲烷分子在管中央发生聚集.另一方面，由于管壁积聚的甲烷过多，出现不稳定，使得一部分距离碳纳米管壁较远的甲烷难以再向管壁富集，从而只能自发地在管中央处发生聚集.

通过比较图3与图4可知，甲烷分子和水分子在管内出现峰谷的位置基本相同.

2.2 扩散系数

分子扩散系数表示它的扩散能力，是物质的物理性质之一.笔者分区域分别计算了每一层水分子和甲烷分子在轴向(z轴)的扩散系数，结果如表1所示.

从分子动力学模拟结果分析可以看出，随着甲烷分子浓度的增加，水分子与甲烷分子的扩散系数不断增加.有趣的是系统1(甲烷分子100个，水分子个数为4 412)，第1层水的扩散系数为8.454×10-4，而当甲烷分子增加到200个时(系统2)，第1层水分子的扩散系数增加到1.414×10-3，是系统1的1.6倍.系统3中的第2层水分子的扩散系数是系统1的4倍.同样，甲烷分子的扩散系数也有类似的变化规律，随着甲烷浓度的增加，各部分的甲烷分子的扩散系数都在增大，特别是第1与第4层的甲烷扩散系数增加较快.在甲烷浓度较高时，甲烷分子更倾向于在管壁和管中央流动，这从甲烷分子的密度分布结果中可以得到证实.为什么水分子与甲烷分子的扩散系数有这样的变化规律呢?这与水分子之间的氢键变化有密切的关系.因此，笔者进一步分析了水分子之间的氢键变化.

表1 各部分的扩散系数

2.3 平均氢键

水分子与水合物的许多奇特性质与水分子之间的氢键有关.本文中，笔者分别计算了每一层水分子内部形成的氢键个数和相邻两层界面处形成的氢键个数.

表2 各部分形成的平均氢键个数

由于第4层水分子是最外层水分子层，所以第4层与外层形成的氢键数为零.从表2中的各个层面上的氢键数量可知，第4层内每个水分子形成的平均氢键个数最多.这是由于第4层水分子密度分布最大，受到碳纳米管壁的吸附作用影响，这个区域中的水分子容易形成有序的水分子层，这和之前由图2得出的结论相符.此外，第3层与外层形成的氢键数量较少，与图2比较，不难发现，在2.14 nm附近水分子密度分布出现一个极小值，导致在2.14 nm附近的水分子形成的氢键数量也随之减少.横向比较3个系统各个层面上形成的平均氢键和总体平均氢键，发现它们的数值随着甲烷个数的增长呈递减趋势，因为在一个饱和体系中甲烷分子的加入，意味着占据了原先水分子的位置，从而必定降低整个系统的水分子之间的平均氢键个数.

由此可知平均氢键个数在几个关键点的数量变化;同时还得知，该过程中表现出来的几个特点与水分子密度分布时相吻合.

3 结 论

利用分子动力学计算机模拟，系统研究了甲烷分子与水分子受限在碳纳米管中的性质.发现了甲烷分子与水分子的径向密度分布的规律，并计算了它们的扩散系数和氢键分布.分析发现，氢键的整体分布和数量符合水分子在管内的密度分布，特别是在2.10 nm到2.67 nm关键位置，水分子密度的变化和氢键数量的变化高度一致.甲烷分子在2.30 nm附近出现了高度富集现象，并且在甲烷含量较高时，一部分甲烷还可能在管中央聚集.随着甲烷浓度的增加，系统在管壁和中央的扩散系数会迅速提高，大大提高甲烷分子沿轴向的扩散能力.

该研究结果对甲烷、水与甲烷水合物受限在纳米孔道中奇特行为的理解有很大的帮助.为进一步理解小空间中的甲烷水合物的形成、分解、输运做了些有意义的前期研究.

［1］赵生才.德国气水合物研究计划简介［J］.天然气地球科学，2001，12(2):63-67.

［2］郑军卫.美国国家甲烷水合物多年研究计划简介［J］.天然气地球科学，2001，12(2):42-45.

［3］宋海斌，松林修.日本的天然气水合物地质调查工作［J］.天然气地球学，2001，12(2):46-53.

［4］Schüth F.Technology:hydrogen and hydrates［J］.Nature，2005，434(7034):712-713.

［5］Sugahara T，Haag J C，Prasad P S R，et al.Increasing hydrogen storage capacity using tetrahydrofuran［J］.Journal of the American Chemical Society，2009，131(41):14616-14617.

［6］Nelson W O，Rettig S J，Orvig C.Exoclathrate Al(C7H8NO2)3.cntdot.12H2O.A facial geometry imposed by extensive hydrogen bonding with the ice 1 structure［J］.Journal of the American Chemical Society，1987，109(13):4121-4123.

［7］Handa Y P.Effect of hydrostatic pressure and salinity on the stability of gas hydrates［J］.Journal of Physical Chemistry，1990，94(6):2652-2657.

［8］Frauenfelder H，Petsko G A，Tsernoglou D.Temperature-dependent X-ray diffraction as a probe of protein structural dynamics［J］.Nature，1979，280:585-563.

［9］Moon C，Taylor P C，Rodger P M.Molecular dynamics study of gas hydrate formation［J］.Journal of the American Chemical Society，2003，125(16):4706-4707.

［10］Hawtin R W，Quigley D，Rodger P M.Gas hydrate nucleation and cage formation at a water/methane interface［J］.Physical Chemistry Chemical Physics，2008，10(32):4853-4864.

［11］Walsh M R，Koh C A，Sloan E D，et al.Microsecond simulations of spontaneous methane hydrate nucleation and growth［J］.Science，2009，326(5956):1095-1098.

［12］Hummer G，Rasaiah J C，Noworyta J P.Water conduction through the hydrophobic channel of a carbon nanotube［J］.Nature，2001，414(6860):188-190.

［13］Koga K，Gao G T，Tanaka H，et al.Formation of ordered ice nanotubes inside carbon nanotubes［J］.Nature，2001，412(6849):802-805.

［14］Ball P.Water:Water—an enduring mystery［J］.Nature，2008，452(7185):291-292.

［15］Fang Haiping，Wan Rongzheng，Gong Xiaojing，et al.Dynamics of single-file water chains inside nanoscale channels:physics，biological significance and applications［J］.Journal of Physics D:Applied Physics，2008，41(10):103002.

［16］Hummer G，Rasaiah J C，Noworyta J P.Water conduction through the hydrophobic channel of a carbon nanotube［J］.Nature，2001，414(6860):188-190.

［17］Lu Hangjun，Li Jingyuan，Gong Xiaojing，et al.Water permeation and wavelike density distributions inside narrow nanochannels［J］.Physical Review B，2008，77(17):174115.

［18］Lu Hangjun，Nie Xuechuan，Wu Fengmin，et al.Controllable transport of water through nanochannel by rachet-like mechanism［J］.Journal of Chemical Physics，2012，136(17):174511.

［19］Wang Chunlei，Li Jingyuan，Fang Haiping.Ordered water monolayer at room temperature［J］.Rendiconti Lincei，2011，22(1):5-16.

［20］Zhou Xiaoyan，Wang Chunlei，Wu Fengmin，et al.The ice-like water monolayer near the wall makes inner water shells diffuse faster inside a charged nanotube［J］.The Journal of Chemical Physics，2013，138(20):204710.

［21］Fang Haiping，Wan Rongzheng，Gong Xiaojing，et al.Dynamics of single-ile water chains inside nanoscale channels:physics，biological significance and applications［J］.Journal of Physics D:Applied Physics，2008，41(10):103002.

［22］Lu Hangjun，Li Jingyuan，Gong Xiaojing，et al.Water permeation and wavelike density distributions inside narrow nanochannels［J］.Physical Review B，2008，77(17):174115.

［23］Lu Hangjun，Nie Xuechuan，Wu Fengmin，et al.Controllable transport of water through nanochannel by rachet-like mechanism［J］.Journal of Chemical Physics，2012，136(17):174511.

［24］Wang Chunlei，Li Jingyuan，Fang Haiping.Ordered water monolayer at room temperature［J］.Rendiconti Lincei，2011，22(1):5-16.

［25］Jorgensen W L，Chandrasekhar J，Madura J D，et al.Comparison of simple potential functions for simulating liquid water［J］.J Chem Phys，1983，79(2):926-935.

［26］Belhadj M，Alper H E，Levy R M.Molecular dynamics simulations of water with Ewald summation for the long range electrostatic interactions［J］.Chemical Physics Letters，1991，179(1):13-20.

［27］Kaminski G，Duffy E M，Matsui T，et al.Free energies of hydration and pure liquid properties of hydrocarbons from the OPLS all-atom model［J］.The Journal of Physical Chemistry，1994，98(49):13077-13082.

［28］Jorgensen W L，Maxwell D S，Tirado-Rives J.Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids［J］.Journal of the American Chemical Society，1996，118(45):11225-11236.