一种模拟I型甲烷水合物包覆稳定II型氢气水合物的方法

一种模拟i型甲烷水合物包覆稳定ii型氢气水合物的方法
技术领域
1.本发明属于水合物储运氢气技术领域，具体涉及一种模拟i型甲烷水合物包覆稳定ii型氢气水合物的方法。


背景技术：

2.随着环境问题和能源危机的日益严重，氢能由于绿色环保、资源丰富等优点获得人们关注。作为氢气从生产到利用过程中的桥梁，储运环节是高效利用氢能的关键。水合物储氢技术具有成本低、易释氢、安全性高的优点。氢气水合物虽然储氢量可观，但是相平衡条件苛刻；此外，氢气分子尺寸小，易从水合物笼中逃逸。如何在温和条件下提高稳定性成为氢气水合物保存的难题。
3.在金属领域的实验和模拟工作中，在熔点低的金属外表面涂覆一层晶格匹配且熔点高的金属层，熔点低的金属就可以在过热的条件下以亚稳态的形式存在，因为两种金属形成连续无缺陷的边界，熔点低的金属缺乏能引发熔化过程的自由表面。类似的例子有银镀金、银镀镍等。基于以上原理，本发明提出利用i型甲烷水合物包覆稳定ii型氢气水合物，首先，甲烷水合物的稳定性高于氢气水合物；其次，i型水合物与ii型水合物在界面处通过5
12
63笼连接，虽然晶体结构存在差异，但仍能形成无缺陷的连续边界。此外，甲烷是燃料气，利用i型甲烷水合物包覆ii型氢气水合物能够提高水合物储气的整体热值。为了探究i型甲烷水合物包覆能否提高ii型氢气水合物的稳定性，本发明从分子层面上观测分析i型甲烷水合物外覆层存在下ii型氢气水合物的结构变化以及氢气分子在水合物内的稳定性，为实验提供指导。分子动力学模拟利用牛顿运动定律，从微观角度研究物质相互作用，其优点是可以实时输出分子的动力学信息，计算分子量较大的体系，广泛应用于各种水合物的模拟计算中。


技术实现要素：

4.针对现有技术中存在的问题，本发明提出一种模拟i型甲烷水合物包覆稳定ii型氢气水合物的方法，所述方法从微观角度上实时定量的记录体系分子的运动轨迹，探究i型甲烷水合物外覆层存在下不同储氢密度的ii型氢气水合物的结构变化以及氢气分子在水合物内的稳定性，为实际中氢气水合物在更温和条件下的储存和运输提供理论上的指导。
5.为达到上述目的，本发明采用以下技术方案：
6.本发明提出一种模拟i型甲烷水合物包覆稳定ii型氢气水合物的方法，所述方法包括以下步骤：
7.(1)采用materials studio建立i型水合物晶胞、ii型水合物晶胞及甲烷和氢气分子模型，并对i型水合物晶胞和ii型水合物晶胞进行2
×2×
2扩胞，将得到的i型水合物超晶胞、ii型水合物超晶胞、甲烷和氢气分子的坐标信息转换为分子模拟软件所需的坐标文件；
8.(2)修改步骤(1)得到的i型水合物超晶胞坐标文件，将步骤(1)得到的甲烷分子模型置于i型水合物超晶胞结构的5
12
62笼和5
12
笼中，同时切割超晶胞获得i型甲烷水合物超晶
胞坐标文件；修改步骤(1)得到的ii型水合物超晶胞坐标文件，将步骤(1)得到的氢气分子模型置于ii型水合物超晶胞结构的5
12
64笼和5
12
笼中，同时切割超晶胞获得具有不同储氢密度的ii型氢气水合物超晶胞坐标文件；
9.(3)采用gromacs软件，带入步骤(1)和步骤(2)得到的坐标文件，扩充气体分子数量，构建气相模型，按照i型甲烷水合物-ii型氢气水合物-i型甲烷水合物-气相的顺序进行叠加，调整ii型氢气水合物与i型甲烷水合物的界面连接方式，得到i型甲烷水合物包覆稳定ii型氢气水合物的初始构型；
10.(4)采用gromacs软件，带入步骤(3)得到的初始构型信息，设置模拟参数，依次通过能量最小化、位置限制性预平衡模拟后得到稳定构型；
11.(5)采用gromacs软件，带入步骤(4)得到的稳定构型信息，设置模拟参数，进行分子动力学计算，得到分子轨迹坐标信息。
12.(6)对步骤(5)得到的分子轨迹坐标信息进行构象分析和计算分析，得到i型甲烷水合物包覆下ii型氢气水合物结构变化的动态构象以及氢气分子在水合物内的稳定效果。
13.优选地，步骤(2)所述的储氢密度通过改变水合物5
12
笼和5
12
64笼中的氢气分子数量进行调整，如l2s1、l3s1、l4s1(l指5
12
64笼，s指5
12
笼)分别代表储氢密度为2.1wt％、2.7wt％、3.2wt％。
14.优选地，步骤(1)所述的2
×2×
2i型水合物晶胞的尺寸为2.406
×
2.406
×
2.406nm3；步骤(1)所述的2
×2×
2ii型水合物晶胞尺寸为3.460
×
3.460
×
3.460nm3；步骤(2)所述的切割和调整后的i型甲烷水合物超晶胞尺寸为2.505
×
3.650
×
2.406nm3；步骤(2)所述的切割和调整后的ii型氢气水合物尺寸为2.505
×
3.650
×
4.450nm3；步骤(3)所述的叠加后的初始模型由i型甲烷水合物-ii型氢气水合物-i型甲烷水合物-气相组成；步骤(3)所述的气相尺寸为2.505
×
3.650
×
2.400nm3；所述叠加后的初始构型尺寸为2.505
×
3.460
×
11.662nm3。
15.优选地，步骤(3)所述ii型氢气水合物与i型甲烷水合物的界面连接方式为通过5
12
63界面笼结构进行连接。
16.步骤(4、5)所述的模拟参数包括力场参数、分子间相互作用参数、温度参数、压强参数、长程作用力参数、键约束算法和模拟时长；其中，所述长程作用力参数采用particle-mesh ewald(pme)方法计算，采用lincs算法限制体系所有分子的键长键角；
17.优选地，所述力场参数包括水分子、氢气分子、甲烷分子的力场；其中水分子力场采用tip4p/ice模型，氢气力场采用三点模型，甲烷力场采用opls-ua模型。
18.优选地，所述分子的相互作用参数包括相同分子间相互作用和不同分子间相互作用，相同分子间相互作用参数为l-j参数，所述不同分子间相互作用参数采用lorentz-berthelot混合规则获得。
19.步骤(4)所述的能量最小化方法采用gromacs软件中的最陡下降算法，所述的位置限制性预平衡模拟包括nvt系综及npt系综，采用蛙跳(leap-frog)算法。
20.优选地，步骤(4)中位置限制性预平衡采用v-rescale耦合方法控制温度参数，温度设置为260k，时间常数为2ps，采用berendsen耦合方法控制压力参数，压力设置为10mpa，时间常数为2ps。
21.优选地，模拟盒子的三个方向都采用周期性边界条件。nvt限制性预平衡模拟时长
为1ns，npt限制性预平衡模拟时长为2ns。
22.上述的模拟参数都需预先设置进相应的nvt.mdp和npt.mdp文件中。
23.步骤(5)采用npt系综，所述的模拟参数中力场参数、分子间相互作用参数与步骤(4)中的模拟参数相同；采用nose-hoover耦合方法控制温度参数，温度设置为260k，时间常数为2ps，采用parrinello-rahman耦合控制压力参数，压力设置为10mpa，时间常数为4ps。模拟时长为800ns。
24.上述的模拟参数都需预先设置在相应的md.mdp文件中。
25.步骤(5)得到的分子轨迹坐标信息输出在trr文件中，为了便于体系中分子运动的可视化观察，利用gromacs自带gmx_trjconv命令对体系的周期性边界条件进行设置。另外，能量输出信息输出在edr文件中，可用自带gmx energy命令调出查看；最终构象信息输出在gro文件中。
26.步骤(6)所述构象分析采用vmd软件进行，打开md.gro文件，加载md.trr文件，对分子轨迹进行可视化处理，进而得到i型甲烷水合物包覆下ii型氢气水合物结构变化的实时动态构象。
27.步骤(6)所述的计算分析包括chill+识别算法以及氢气水合物内氢气分子的留存率，晶体结构的水分子能与周围水分子形成四个氢键，根据氢键两端水分子排列的差异将氢键分为重叠键和交错键，其中水合物内水分子的重叠键为4个，交错键为0个；冰ih内水分子的重叠键为1个，交错键为3个；冰ic内水分子的重叠键为0个，交错键为4个。通过统计固定区域内被识别为水合物的水分子数量对氢气水合物的分解情况进行分析。对于交错键与重叠键的区分采用下列公式计算：
[0028][0029]
其中，q
l
(i)是i原子的键取向序参数，q
l
(j)是j原子的键取向序参数，c(i,j)为相关系数，当0.25≥c(i,j)≥-0.35时，i与j原子之间的键被视为重叠键，而当c(i,j)≤-0.8时，i与j原子之间的键被视为交错键。
[0030]
氢气水合物内氢气的留存率采用下列计算公式：
[0031][0032]
其中n
h2 in hydrate
指停留在水合物内的氢气分子数目，n
total h2
指总的氢气分子数目。
[0033]
与现有的技术方案相比，本发明具有以下有益效果：
[0034]
本发明提出一种模拟i型甲烷水合物包覆稳定ii型氢气水合物的方法，所述方法从微观分子层面上直观分析i型甲烷水合物包覆下ii型氢气水合物的结构变化以及氢气分子在水合物内的稳定性，节省实验成本；并为实际中氢气水合物在更温和条件下储存和运输提供理论指导；应用i型甲烷水合物包覆ii型氢气水合物，一方面提高氢气水合物的稳定性，可以在更温和的条件下过热保存；另一方面，甲烷是燃料气，选择i型甲烷水合物作为包覆层能够提高水合物储气的整体热值。
附图说明
[0035]
图1为氢气水合物晶胞图。
[0036]
图2为所述方法中i型甲烷水合物包覆稳定ii型氢气水合物体系的初始构象图。
[0037]
图3为实施例1和2体系中氢气水合物区域被chill+算法识别为水合物水分子的数量随模拟时间的变化图。
[0038]
图4为实施例2、3和4体系中氢气分子在水合物内的留存率。
[0039]
具体实施方法
[0040]
为更好的说明本发明，下面结合附图及具体实施方式对本发明的内容作进一步详细说明，但本发明的实施方式不限于此。
[0041]
实施例1
[0042]
一种模拟无i型甲烷水合物包覆稳定储氢密度为2.1wt％的ii型氢气水合物的方法，所述方法包括以下步骤：
[0043]
(1)采用materials studio建立ii型水合物晶胞和氢气分子模型，将ii型水合物晶胞进行2
×2×
2扩胞，尺寸为3.460
×
3.460
×
3.460nm3，将得到的晶胞结构及氢气分子的坐标信息转换为分子模拟软件所需的pdb坐标文件。
[0044]
(2)修改步骤(1)得到的ii型水合物超晶胞坐标文件，将步骤(1)得到的氢气分子模型置于ii型水合物超晶胞结构孔穴中，5
12
64笼中放置2个氢气分子，5
12
笼中放置1个氢气分子。同时切割超晶胞获得储氢密度为2.1wt％、尺寸为2.505
×
3.650
×
4.450nm3的ii型氢气水合物超晶胞坐标文件；
[0045]
(3)采用gromacs软件，带入步骤(1)和步骤(2)得到的坐标文件，扩充气体分子数量，构建尺寸为2.505
×
3.650
×
2.400nm3的气相模型，按照氢气水合物相-气相的顺序进行叠加，得到尺寸为2.505
×
3.650
×
6.850nm3的体系初始构型；
[0046]
(4)采用gromacs软件，带入步骤(3)得到的初始构型信息，采用steepest descent minimization算法能量最小化后，依次进行蛙跳算法(leap-frog)的nvt和npt位置限制性预平衡模拟，得到稳定构型信息。模拟参数包括力场参数、分子间相互作用参数、温度参数、压强参数、长程作用力参数、键约束算法和模拟时长，都预先设置在相应的mdp文件中。其中，氢气分子采用三点氢气模型；水分子采用tip4p/ice模型；相同分子间相互作用参数采用l-j参数，不同分子间相互作用参数采用lorentz-berthelot混合规则得到；采用v-rescale耦合控制温度参数，设置为260k，时间常数为2ps；压强参数采用berendsen耦合，设置为10mpa，时间常数为2ps；长程作用力采用pme方法计算，选用lincs算法限制所有分子的键长键角，模拟盒子的三个方向都采用周期性边界条件。nvt预平衡模拟时长为1ns，npt预平衡模拟时长为2ns。
[0047]
(5)采用gromacs软件，带入步骤(4)得到的构型信息，并且在mdp文件中设置模拟参数，其中力场参数、分子间相互作用参数、长程作用力参数、键约束算法与步骤(4)中设置相同，采用nose-hoover耦合控制温度参数，设置为260k，时间常数为2ps；压强参数采用parrinello-rahman耦合，设置为10mpa，时间常数为4ps；模拟时长为400ns。计算得到的分子轨迹坐标信息输出在trr文件中，为了便于体系中分子运动的可视化观察，利用gromacs自带gmx_trjconv命令对体系的周期性边界条件进行设置。另外，能量输出信息输出在edr文件中，可用自带gmx energy命令调出查看；最终构象信息输出在gro文件中。
[0048]
(6)对步骤(5)得到的分子轨迹坐标信息和动力学信息进行构象分析和计算分析：构象分析采用vmd软件，打开gro文件并加载trr文件，对分子轨迹做可视化处理，进而得到
minimization算法能量最小化后，依次进行蛙跳算法(leap-frog)的nvt和npt位置限制性预平衡模拟，得到稳定构型信息。模拟参数包括力场参数、分子间相互作用参数、温度参数、压强参数、长程作用力参数、键约束算法和模拟时长，都预先设置在相应的mdp文件中。其中，氢气分子采用三点氢气模型，甲烷分子采用opls-ua联合原子模型；水分子采用tip4p/ice模型；相同分子间相互作用参数采用l-j参数，不同分子间相互作用参数采用lorentz-berthelot混合规则；采用v-rescale耦合控制温度参数，设置为260k，时间常数为2ps；压强参数采用berendsen耦合，设置为10mpa，时间常数为2ps；长程作用力采用pme方法计算，选用lincs方法限制所有分子的键长键角，模拟盒子的三个方向都采用周期性边界条件。nvt预平衡模拟时长为1ns，npt预平衡模拟时长为2ns。
[0061]
(5)采用gromacs软件，带入步骤(4)得到的构型信息，并且在mdp文件中设置模拟参数，其中力场参数、分子间相互作用参数、长程作用力参数、键约束算法与步骤(4)中设置相同，采用nose-hoover耦合控制温度参数，设置为260k，时间常数为2ps；压强参数采用parrinello-rahman耦合，设置为10mpa，时间常数为4ps；模拟时长为800ns。计算得到的分子轨迹坐标信息输出在trr文件中，为了便于体系中分子运动的可视化观察，利用gromacs自带gmx_trjconv命令对体系的周期性边界条件进行设置。另外，能量输出信息输出在edr文件中，可用自带gmx energy命令调出查看；最终构象信息输出在gro文件中。
[0062]
(6)对步骤(5)得到的分子轨迹坐标信息和动力学信息进行构象分析和计算分析：构象分析采用vmd软件，打开gro文件并加载trr文件，对分子轨迹做可视化处理，进而得到储氢密度为2.1wt％的ii型氢气水合物在i型甲烷水合物包覆的情况下结构变化的实时动态构象。计算分析包括chill+识别算法以及氢气水合物内氢气分子的留存率，晶体结构的水分子能与周围水分子形成四个氢键，根据氢键两端水分子排列的差异将氢键分为重叠键和交错键，其中水合物内水分子的重叠键为4个，交错键为0个；冰ih内水分子的重叠键为1个，交错键为3个；冰ic内水分子的重叠键为0个，交错键为4个。通过统计固定区域内被识别为水合物的水分子数量对氢气水合物的分解情况进行分析。对于交错键与重叠键的区分采用下列公式计算：
[0063][0064]
其中，q
l
(i)是i原子的键取向序参数，q
l
(j)是j原子的键取向序参数，c(i,j)为相关系数，当0.25≥c(i,j)≥-0.35时，i与j原子之间的键被视为重叠键，而当c(i,j)≤-0.8时，i与j原子之间的键被视为交错键。
[0065]
氢气水合物内氢气的留存率采用下列计算公式：
[0066][0067]
其中n
h2 in hydrate
指停留在水合物内的氢气分子数目，n
total h2
指总的氢气分子数目。
[0068]
图3显示了例2中氢气水合物区域内水合物水分子随模拟时间的变化。水合物水分子数量在800ns内维持稳定，表明氢气水合物未出现分解。对比例1，ii型氢气水合物在i型甲烷水合物的包覆下，结构稳定性明显提高。图4显示了模拟结束后例2中氢气水合物内氢气的留存率，数值高达100％，表明氢气分子全部停留在水合物内。
[0069]
实施例3
[0070]
一种模拟i型甲烷水合物包覆稳定储氢密度为2.7wt％的ii型氢气水合物的方法，所述方法除了步骤(2)中5
12
64笼中插入3个氢气分子，5
12
笼中插入1个氢气分子外均与实施例2相同。
[0071]
图4显示了模拟结束后例3中氢气水合物内氢气的留存率，数值为98.77％。
[0072]
实施例4
[0073]
一种模拟i型甲烷水合物包覆稳定储氢密度为3.2wt％的储氢水合物的方法，所述方法除了步骤(2)中5
12
64笼中插入4个氢气分子，5
12
笼中插入1个氢气分子外均与实施例2相同。
[0074]
图4显示了模拟结束后例4中氢气水合物内氢气的留存率，数值为94.86％。随着氢气水合物储氢密度从2.1wt％升高至3.2wt％(即案例2、3、4)，氢气在水合物内的留存率从100％降低至94.86％。通过实施案例2-4与案例1比较，可以清楚看出i型甲烷水合物的包覆能够将氢气分子稳定在水合物内，不过稳定效果受到氢气水合物储氢密度的影响。
[0075]
以上所述，仅用于说明本发明的技术方案而非限制，尽管参照较佳实施例对本实施例的技术方案进行了详细说明，但本发明的保护范围并不局限于此，本领域技术人员应当理解，基于本发明中的实施例，任何在技术范围内可轻易想到变化或替换，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围，因此本发明的保护范围应以所述权利要求的保护范围为准。