原《天然气化工—C1化学与化工》，经国家新闻出版署批复更名为《低碳化学与化工》。

  王芳，女，西南石油大学副教授，硕士生导师。致力于CO2捕集绿色新技术及转化利用研究，以第一作者或通讯作者身份发表SCI论文27篇，研究成果获四川省自然科学奖二等奖1项、中国化工学会科学技术“基础研究成果奖”二等奖2项，授权发明专利6项。

  西南石油大学周莹教授课题组2016年入选四川省青年科技创新研究团队，2021年入选“天府峨眉计划”团队。研究方向最重要的包含太阳能燃料技术、绿电转化综合利用和油气田开发利用功能材料。承担国家及省部级科研项目40余项，承担省科协首批科技智库调研课题等软科学项目5项，获省部级科研奖励11项。在Nature Commun、Angew Chem Int Ed等期刊发表SCI论文220余篇，被Nature/Nature Chemistry等正面引用9400余次，H指数56， ESI高被引论文19篇；授权中国发明专利30件，授权美国发明专利1件。培养博士/硕士研究生近100人次，留学深造10余人次。

  （1. 西南石油大学 新能源与材料学院，碳中和研究院，四川成都610500 ；2. 西南石油大学 油气藏地质及开发工程国家重点实验室，四川成都610500）

  基金项目：四川省重点研发项目（2022YFSY0050）；四川省自然科学基金（2023NSFSC0108）。

  摘要：CO2 置换开采天然气水合物法可在开采CH4 的同时将CO2 以固态水合物形式封存海底，对于解决能源和环境问题意义重大，但目前该技术尚不成熟，置换率较低。以构建CO2 及混合气体置换水合物团簇及晶体超胞结构为基础，通过计算置换能、相互作用能及态密度建立了CO2 占据孔穴状态、相互作用能和置换能三者之间的关联，研究了置换反应的可行性和CH4开采率受限的微观机理，同时进一步通过分子动力学预测了晶体结构的稳定性。根据结果得出，CO2 优先置换水合物51262 大笼中的CH4，相比于CH4 与51262 大笼间的相互作用能，其相互作用能增大2.26 eV，可稳定水合物结构，置换能为2.16 eV。当CO2开始占据512 小笼时，其相互作用能比51262 大笼小0.17 eV。以n(CO2):n(N2) = 1:3 的混合气体全部置换CH4 具有热力学可行性，且晶体结构在3 MPa、263 K下可稳定存在，水分子的均方位移小于CH4 水合物，表明混合气体置换CH4 后水合物晶体具有动力学稳定性。该结果为CO2 及混合气体置换天然气水合物技术的基础研究和应用提供了科学依据和理论指导。

  全球工业现代化进程的加快速度进行发展，引发了一系列环境污染和能源短缺危机。中国作为目前全球CO2 排放量最高的国家[1-2]，碳捕集、利用与封存（CCUS）技术起步较晚，实现2030 年前“碳达峰”及2060 年前“碳中和”的目标艰巨，发展CO2 封存与转化利用技术迫在眉睫。近年来CO2 置换天然气水合物法作为新兴CO2 储存技术，自OHGAKI等[3]提出后，在环境治理和清洁能源利用领域受到国内外高度关注[4-6]。该方法既能将CO2 转化为固态水合物稳定封存于海底，又可从水合物中开采天然气，相比于注热法、降压法和注化学剂法[7-9]等开采手段，在减轻温室效应的同时可稳固海床，减小地质灾害和避免环境污染，实现CO2 减排和CH4 清洁能源获取，可谓“一石三鸟”，对实现“碳中和”具备极其重大意义。

  已有研究[10-11]发现，单位物质的量的CO2 水合物生成释放的热量（-57.98 kJ/mol）能够供单位物质的量的CH4 水合物分解吸收（54.49 kJ/mol），可推断此置换反应易发生。王林军等[12]和HUANG等[13]证明了CO2 置换CH4 水合物反应的可行性，并通过拉曼光谱观察到CO2 选择占据水合物孔穴的51262 大笼。由于sⅠ型CH4 水合物中51262 大笼与512 小笼孔穴的数量比为3:1，因此理论上CH4 置换率最高达75%。目前，CO2 置换的回收率为64%[14]，虽然加热法和置换法耦合可提高其回收率，但加热法中的注热升温等辅助手段会带来经济性和环保等额外问题[15]，因此，怎么样提高置换率是置换法需解决的关键问题。

  CO2 和N2 混合气体（烟道气）作为置换气具有废气利用和提高置换率的特点，近年来非常关注[16-18]。范德华半径较小的N2 分子易占据512 小笼，与占据51262 大笼的CO2 分子组合可有效提升孔穴占据率，来提升CH4 置换率。PARK等[18]在274.15 K、3.5 MPa实验条件下通过V(CO2):V(N2) = 1:4 的混合气将置换率提升至85%；操原[19]通过混合气辅以注热法将置换率提升至91.6%；目前已有理论计算表明，N2可促进CO2 向水合物更深层次扩散以进一步促进置换反应进行[19]。虽然客体分子的范德华半径与笼径的比值大小可判断其占据态，但针对置换反应过程中客体分子与水笼的相互作用、晶体结构的稳定性和置换反应难易程度间的关联仍缺乏认知，对置换反应微观机理尚不明确，限制了置换率的进一步提高。

  基于以上研究，本文以sⅠ型CH4 水合物这一自然界存在最为广泛的水合物为理论计算模型，采用密度泛函理论（DFT）的第一性原理及分子动力学（MD）计算方式，深入研究团簇与晶体中混合气体的置换反应，通过对比置换前后晶体结构和电子结构、客体分子与水笼相互作用以及置换能，明确客体分子占据态、相互作用能和置换趋势三者的关联，并在此基础上预测CO2 和N2 混合气体完全置换的热力学可行性。

  sⅠ型CH4 水合物结构参数如表1 所示[20-22]。sⅠ型CH4 水合物为体心立方晶系结构，空间群Pm3n，晶格常数1.188 nm，其中，有两种不一样的笼型孔穴，单胞含有2 个512 小笼和6 个51262 大笼，512 小笼和51262 大笼的平均半径分别为0.395 nm和0.433 nm，512 小笼整体的结构由20 个H2O分子构成，51262 大笼整体的结构由24 个H2O分子构成。氧原子位置由X射线]，氢原子分布遵循“冰规则”[23]。

  两种类型的水合物笼状结构模型如图1 所示，sⅠ型CH4 水合物晶体结构如图2 所示。几何中心为客体分子，512 小笼可容纳N2、CH4 等小分子，51262 可容纳分子半径略大的CO2 分子，客体分子填充于孔穴中支撑晶体结构并稳定存在。水合物结构与客体分子的填充率有关，填充率越大结构越稳定，当客体分子和主体分子孔穴的直径比达到77%以上时，其状态最为稳定。

  本文所有第一性原理计算均使用VASP软件[24-25]，对笼型团簇及晶体超胞结构可以进行结构优化（晶胞类型和原子位置）和能量收敛计算。计算设置如下：采用GGA/PBE交换泛函[26]和PAW平面波赝势，Gaussian路径修正，展宽0.2 eV。平面波函数截断能400 eV，自洽迭代收敛精度1.0 × 10-5 eV/atom，所有计算均在倒易空间中进行。采用Monkorst-Pack对布里渊区k点积分，k点网格设置为1 × 2 × 2，原子间相互作用力小于0.5 eV/nm，原子位移小于0.0002 nm，最大应力为0.1 GPa。分子动力学使用Materials Studio软件的Forcite软件包进行计算，采用dynamic模块，选取263 K和3 MPa作为温压条件，在等温等压（NPT）下进行模拟。对结构最优化处理后，使用Nose-Hoover控制温度，在正则系综下进行50 ps的等温等容（NVT）计算，随后解除对氧原子限制，使用Andersen控制压力进行100 ps的等温等压计算。

  晶体的稳定性取决于客体分子与水笼相互的作用[27]，相互作用能越负，表明客体分子与水笼的结合能力越强，水合物越稳定。置换能（Eexchange，eV）可判断CO2 置换CH4 反应的难易，置换能越小，置换反应越易发生[29-30]。客体分子与H2O分子间主要是通过范德华力作用，相互作用能（Eguest-cage，eV）计算公式[27-28]如式(1)所示。

  式中，Etotal代表总能量，eV；Eguest代表客体分子能量，eV；Ecage代表为笼型结构的所有H2O分子总能量，eV。

  在CO2 占据51262 大笼后会出现两个CO2 同时占据1 个51262 大笼的双占据情况[34]，第3 个CO2 进入其中1 个被CO2 占据的空穴并与之形成双占据的Eguest-cage计算公式如式(2)所示。

  式中，m和n分别代表客体分子CO2 与N2 的数量；Eproducts代表生成物总能量（CH4 分子与CO2 水合物的能量加和），eV；Ereactant代表CO2 分子与CH4 水合物的能量加和， eV。

  CO2 置换CH4 水合物团簇的反应如图3 所示，H2O分子通过氢键连接，512 小笼有30 个氢键，51262大笼有36 个氢键，具有高度Td对称性的CH4 的占据位置受氢键影响较小。由共价键键合的非极性线 中氧原子与水中氢原子形成氢键，该氢键作用对CO2 在团簇中占据位置的影响不可忽视。

  由图3 可知，512 团簇由12 个五边型构成，CO2分子核间轴指向位于对立面五边形中心时，团簇结构最为稳定，其占据态与FERNANDO等[32]计算结果一致。51262 团簇中包含氢键数量较多，相对六边形对CO2 分子起到拉伸作用，产生“悬浮”于两平行面的效果。侧面的氢键作用导致CO2 与水平面呈10°左右倾斜角度，收敛结果与UDACHIN等[31]通过X射线衍射得出的占据态一致（水平倾斜角度为6.5°~14.4°）。

  上述置换反应中，CO2 置换512 团簇中的CH4 所需0.66 eV的能量，而置换51262 团簇中的CH4 释放0.44 eV的能量，因此CO2 置换51262 团簇中的CH4 反应更容易自发进行，具备热力学置换反应的可行性。

  为进一步研究水笼和客体分子的电子轨道以及水合物稳定作用的机理，通过态密度（DOS）和分波态密度（PDOS）研究不同分子间轨道叠加和键合情况[33]，客体分子及占据客体分子的团簇的PDOS如图4 所示。

  由图4(a)~(c)可知，512 小笼和51262 大笼中CH4分子H 1s轨道相较于孤立CH4 向能量更负方向挪动，表明CH4 分子使水笼更稳定[34-35]，CH4 分子与H2O分子间的相互作用是CH4 水合物结构稳定存在的关键。对于512 小笼CH4 团簇，在费米能级附近，在能量-2 eV左右出现了特征峰，主要由H2O中O 2p轨道电子贡献，碳原子和氢原子外层的电子贡献较小；在能量-4 eV左右出现了微小的峰，因CH4 中C 2p和H 1s轨道叠加形成了σ键。51262 大笼CH4 团簇轨道贡献情况与512 小笼CH4 团簇类似，即价带顶和导带底主要由H2O分子的O 2p轨道贡献，但是峰值呈弥散状态。此外，512 小笼在能量为-4 eV处出现非常明显的O 2p和H 1s峰叠加，51262 大笼在能量为-5 eV处出现了峰偏移，表明CH4 与512 小笼中H2O相互作用更强，CH4 更易稳定存在于512 小笼中。由图4(d)~(f)可知，51262 大笼中的CO2 分子的O 2p轨道相较于孤立CO2 向能量更负的方向挪动，并且51262大笼相较于512 小笼移动范围更大，CO2 分子在51262 大笼中更稳定。51262 大笼CO2 团簇中，在费米能级附近，在能量为-8 eV附近出现特征峰叠加，主要由CO2 中O 2p轨道电子贡献，在能量为-7 eV处出现特征峰且发生偏移，表明CO2 与51262 大笼中H2O相互作用更强，CO2 更易稳定存在于51262 大笼中。

  为深入探究置换反应前后水合物的热力学稳定性，采用不一样浓度（物质的量分数，下同）的CO2及N2 占据不一样的水笼孔穴以实现置换率的递增，以在提升置换率的同时通过客体分子的置换能和相互作用能判断客体分子的占据态及结构稳定性，进而通过置换后水合物的晶体结构稳定性及变化趋势判断置换反应的可行性。

  不同置换率下的占据态如图5 所示，除去边界效应后晶体中占比权重为1 的51262 大笼和512 小笼各两个。构建上述方案研究了混合气在晶体内完整孔穴结构的置换趋势，通过改变客体分子的占据态和占据率以模拟置换率，并构建置换率为100%的混合气体水合物模型，进一步开展了完全置换下的水合物稳定性研究。

  不同置换率的晶体稳定结构如图6 所示，对各结构均进行了结构优化和能量收敛计算。在相应置换率下的客体分子置换后，晶体结构仍然能够保持稳定，H2O分子相对于起始位置变化甚微，这是由于CO2 和N2 与H2O之间的强化相互作用导致的。

  晶胞参数的平均变化率指置换反应中晶胞结构参数中相较于起始结构的变化趋势，体现了置换反应前后的结构稳定性，如表3 所示。不同置换率下的水合物晶体参数平均变化率均保持在2.00%以内，与H2O分子产生更强相互作用的CO2 和N2 分子在置换后使得结构有所收缩。当置换率为100%时，CO2 和N2 的水合物晶体结构出现轻微膨胀，b0 和c0的参数增大，0.51%的平均变化率足以证明完全置换后结构仍保持高稳定性。由此可知，CO2 置换后，sI型CH4 水合物结构保持稳定性，并且CO2 与N2 混合气体置换反应对水合物结构的变化影响较小，从微观结构的变化解释了CO2 和N2 置换CH4 水合物的可行性，这也是CO2 置换法开采CH4 水合物使得海床稳固的原因之一。

  置换率提升至12.50%，首先形成的占据态将影响后续CO2 置换CH4 水合物的反应，反应式如式(6)所示。

  置换率逐步提升至18.75%，此时除边界效应外，晶体中完整51262 大笼的CH4 已全部置换，后续CO2 进入512 小笼置换CH4 的反应式如式(8)所示。

  相互作用能和置换能计算结果如表4 所示。CH4 与512 小笼的相互作用能为-0.93 eV，远小于51262 大笼，表明CH4 与512 小笼间的作用力更强，CH4更趋向于占据sⅠ型水合物的512小笼孔穴，这与张庆东[5]基于NPT系综下MD模拟水合物，计算团簇的结论（CH4 与512 小笼的相互作用能为17.96 kJ/mol，与51262 大笼的相互作用能为-16.08 kJ/mol）一致，并且符合HUANG等[13]通过拉曼成像得出的CH4 占据小笼孔穴结论。当CO2 置换率为6.25%时，CO2 与51262 大笼的相互作用能（-2.24 eV）更负，CO2 优先占据51262 大笼孔穴发生置换反应，相比CH4，其相互作用能增大2.26 eV，可稳定水合物结构，CO2 与sⅠ型CH4 水合物51262 大笼的结合能力比512 小笼更强，当CO2 开始占据512 小笼时，其相互作用能比51262大笼小0.17 eV，且相互作用能的变化趋势与CO2 占据不同水笼的结论[33]一致。当CO2 置换率提升至12.5%，水合物的稳定性逐渐增强，后进入51262 大笼的CO2 置换能为0.12 eV，相较于首个CO2 置换能降低2.04 eV，首个置换的CO2 极大提高了水笼稳定性，并促进了置换反应的进行。逐步提升置换率为18.75%时，第3 个CO2 的相互作用能为-2.02 eV，相较于双占据态，其相互作用更强，当51262 大笼完全被CO2 占据后，CO2 趋向于占据512 小笼孔穴，这与双占据态间的斥力有关。但第3 个CO2 的置换能为3.96 eV，所需的能量较高，难以持续发生置换反应，使置换率提高受限，是CO2 置换CH4 水合物过程中置换率较低的原因之一[32]。

  N2 是目前报道可稳定存在于水合物中最小的客体分子[7]，因此，为进一步提升CH4 置换率，引入客体分子N2 参与置换反应。N2 与512 小笼相互作用能为-2.23 eV，稳定性增强，且置换能小于第3 个CO2 置换512 小笼中CH4 的置换能，因此N2 置换有利于进一步提升置换率。在置换率达到100%后，虽-0.20 eV的平均相互作用能低于低浓度置换率时的作用能，但是可稳定水合物结构，这可能是CO2 和N2 的高浓度置换率使得相邻孔穴之间产生了影响，0.65 eV的平均置换能使得置换反应更易进行，使用n(CO2):n(N2) = 1:3 实现100%置换率具备热力学可能性。

  综上所述，混合气体CO2 和N2 置换CH4 的最佳路径如图7 所示，相较于置换51262 大笼中的CH4，置换512 小笼在热力学上反应难度更大。在CH4 水合物置换反应中CO2 优先置换51262 大笼中CH4，且占据51262 大笼后，后续置换反应的置换能降低，将有利于后续置换反应的持续进行，N2 的加入有利于置换512 小笼中CH4 以提升整体的CH4 水合物置换率，但是3.68 eV的置换能使得置换反应更难进行。基于以上分析，利用CO2 和N2 混合气体置换sI型CH4 水合物可以在一定程度上完成100%的置换率。

  不同客体分子占据态（占据率为6.25%）的PDOS如图8 所示。由图8(a)可知，在费米能级附近，能量在1.5 eV附近出现了O 2p特征峰，-4 eV左右出现了CH4 H 1s与H2O O 2p的轨道叠加的峰，说明氢键利于水合物的稳定[36]。除特征峰位置能量稍有变化，由图8(b)可知，CH4 在51262 大笼中的轨道贡献与512 小笼类似。CO2 进入水合物置换CH4后，客体分子与水笼的氢键类型由C—H…O变为C—O…H，氢键作用主要由CO2 O 2p和H2O 1s轨道贡献。由图8(c)和图8(d)可知，51262 大笼中CO2 O 2p轨道在价带顶处发生峰的劈裂，且轨道分布更加弥散，表明CO2 与H2O相互作用更强，有利于CO2 进入51262 大笼发生置换。

  为进一步研究客体分子与H2O分子间的电子转移情况，通过差分电荷密度图计算各个最优占据态的差分电荷密度，结果如图9 所示。由图9 可知，相比于CH4，N2 与水笼的电子聚集程度更高，这说明N2 与H2O分子之间的电子转移相较于CH4 更多，CO2 与水笼的电荷转移量（0.041 e）大于CH4 与水笼的电荷转移量（0.012 e），电子转移量越多，电子越聚集，证明CO2 与水笼中的H2O分子相互作用更强。

  不同客体分子占据态（占据率 6.25%）的PDOS如图10 所示。由图10(a)和图10(b)可知，随CO2 的增加，CO2 C 2p 和O 2p轨道分布更加弥散，与H2O H 1s在-9~-6 eV能量区间的PDOS叠加程度更强，CO2 的增加有利于增强晶体结构的稳定性。由图10(c)可知，晶体结构中的N2 相较于孤立的N2 2p向能量更低的方向挪动，归因于N 2p与H2O的相互作用使N2 稳定存在于512 小笼孔穴中，这为N2 易占据512 小笼孔穴提供理论依照。N2 与512 小笼之间的电荷转移量（0.024 e）大于CH4 与512 小笼之间的电荷转移量（0.012 e），所以N2 与512 小笼相互作用较强。由图10(a)和图10(c)可知，N2 对相邻占据51262 大笼的CO2 的PDOS分布影响甚微，相邻水笼之间的相互影响较小，水合物的稳定性主要来自于客体分子与水笼之间的相互作用。

  基于DTF计算得出最佳孔穴占据态的水合物结构，选取温度为263 K、压力为3 MPa，对3 种类型的水合物进行了动力学模拟以研究水合物晶体结构的动力学稳定性。

  在263 K、3 MPa条件下的水合物模型的最终构象如图11 所示。与初始结构相比，氧原子的位置移动较小，并且重叠度较高，氢键排列整齐且保持笼型结构，晶体结构稳定存在与实验结果贴合[37-40]。

  径向分布函数可反映水合物中的微观结构特征，go-o(r)表示一个氧原子与另外一个氧原子在距离为r出现的概率，为系统的区域密度与平均密度的比值。由于水合物晶体中H2O分子的排列顺序具有一定的周期性，因此稳定的水合物晶体的go-o(r)具有多个峰[5]。263 K、3 MPa下水合物最终构象的go-o(r)如图12 所示。CH4 水合物的go-o(r)曲线峰值较低，峰谷较高，而当CO2 替换51262 大笼的CH4 之后，峰值增加，峰谷降低，表明H2O分子排列的无序性增强，笼形结构的有序性得到了提高，说明CO2 能够稳定水合物[41]，这与前文所提到的CO2 与水合物的相互作用大于CH4 的结果吻合；而N2 替换CH4后在0.25 nm~0.30 nm的峰值出现显著地增强，同样表明N2 的加入增强了客体分子与水笼的相互作用，使H2O分子有序性得到提高，而在0.30~0.60 nm的峰值和峰谷与CH4 水合物基本等同，表明混合气能够稳定水合物晶体结构，水合物晶体结构在此压力条件下稳定存在。

  均方位移（MSD）能够反映水合物体系内氧原子的空间位置相对于其初始位置的偏离程度，判断水合物能否稳定存在。MSD较小时，氧原子的偏离程度较低，说明水合物晶体结构变化较小，能够稳定存在[5]。

  在263 K、3 MPa下，CH4 水合物的氧原子MSD曲线 及混合气体水合物的MSD跟着时间呈现先增大后减小的趋势，在0.6 × 10-20~1.2 × 10-20 m2 的范围波动，表明在模拟过程中H2O分子在其位置产生振动，结合上述结论，CH4 水合物仍以晶体结构存在。对于CO2 置换51262大笼CH4 后的水合物来说，其MSD在0.4 × 10-20 m2的小范围内波动，表明CO2 置换后的水合物结构更稳定，使氧原子的振动值大幅度的降低，表现出良好的晶体结构稳定性，进而随着N2 置换512 小笼后，其氧原子的振动值增大，在0.5 × 10-20~0.7 × 10-20 m2 的范围内波动，小于CH4 水合物中氧原子的振动范围，其分解速率大于CO2 和CH4 混合气体水合物，但是小于CH4 水合物，证明CH4 水合物的置换率达100%具有理论可行性。

  本文基于第一性原理与分子动力学，研究了CO2占据孔穴状态、相互作用能和置换能三者之间的关系，探究了置换反应的热力学可行性和CH4 开采率受限的微观机理，并对CO2 和N2 混合气体提高水合物中CH4 置换率达100%的热力学可行性和动力学稳定性进行了预测，得到如下结论。

  （1）随着部分CO2 优先占据CH4 水合物51262大笼，水合物晶体稳定性增加，同时反应的置换能降低，促进了CO2 置换反应持续进行。

  （2）N2 与CH4 水合物512 小笼更强的相互作用使其更容易置换512 小笼孔穴中的CH4，可通过N2稳定512 小笼孔穴和CO2 占据51262 大笼的方式进一步提升CH4 置换率。

  （3）以n(CO2):n(N2) = 1:3 的混合气分别置换51262 大笼和512 小笼中的CH4 具有热力学可行性，客体分子与水笼的相互作用能使水合物的结构保持稳定，较小的置换能能轻松实现CH4 置换率达100%，且置换后水合物中H2O分子的均方位移减小，晶体结构在263 K、3 MPa下能够保持动力学稳定。

  引用本文：王芳，陈佳豪，黄泽皑，等. CO2 及混合气体提高CH4 水合物置换率机理研究[J]. 低碳化学与化工, 2023, 48(3): 154-164.

  本文为原创作品，文字及图片版权均属《低碳化学与化工》编辑部所有，如需转载、摘编本文，务请联系本编辑部。欢迎按本微信原文转发分享！

  特别声明：以上内容(如有图片或视频亦包括在内)为自媒体平台“网易号”用户上传并发布，本平台仅提供信息存储服务。

  伊周周日热推：综艺《跳进地理书的旅行第二季》；电视剧《19层》......

  学龄期的孩子，每年都要做两次的眼睛检查哦！假期到了，快带你家孩子安排上#儿童的护眼小常识

  六个月宝宝体检时一直骂骂咧咧。还想攻击人家医生奈何没长牙 没成功。妈妈体检能有八百个心眼刚打过针以为又要再打一针紧皱眉头。