离子液体吸收式制冷工质对气液相平衡研究及其分子动力学模拟

摘要：为了克服传统吸收式制冷工质对的种种缺陷，离子液体型制冷工质对被视作理想替代品，得到了广泛的开发与应用。采用静态实验和分子动力学模拟两种方法，对［EMIM］BF4/CH3OH离子液体二元体系的气液相平衡性质展开研究。研究结果显示，［EMIM］BF4/CH3OH离子液体二元溶液具有良好的气液相平衡性质，实验测得饱和蒸气压比其他醇类离子液体二元溶液低约21%；模拟结果与实验值有相同的数量级和变化趋势，相对误差普遍小于8%。研究结果为离子液体制冷工质对的筛选和进一步的理论循环系统研究提供了物性数据基础，另一方面也为离子液体基础物性研究提供了一种模拟预测的新方法。

关键词：离子液体；吸收式制冷；分子动力学模拟；气液相平衡

中图分类号：TB61+2"" 文献标志码：A"" 文章编号：1002-4026（2024）04-0084-09

开放科学（资源服务）标志码（OSID）：

Study on gas-liquid phase equilibrium properties of ionic liquid absorption

refrigeration working pairs and their molecular dynamics simulations

ZHANG Haoran，CHEN Wei*

（College of Electromechanical Engineering， Qingdao University of Science ＆Technology， Qingdao 266061， China）

Abstract∶To overcome the shortcomings of traditional absorption refrigeration working pairs， ionic liquid refrigeration working pairs have been widely developed and used as ideal substitutes. Herein， the gas-liquid phase equilibrium properties of the ［EMIM］BF4/CH3OH ionic liquid binary system were investigated using static experiments and molecular dynamics simulations. The results reveal that this binary solution shows favorable gas-liquid phase equilibrium properties， and its saturated vapor pressure is experimentally measured to be approximately 21% lower than that of other alcohol-based ionic liquid solutions. In addition， the simulation results exhibit the same order of magnitude and trend as the experimental results， and the relative errors are generally less than 8%. These findings provide physical property database for screening ionic liquid refrigeration working pairs and studying the theoretical cycle system as well as a new method for simulating and predicting the basic properties of ionic liquids.

Key words∶ionic liquid; absorption refrigeration; molecular dynamics simulation; gas-liquid phase equilibrium

近年来，随着我国经济快速发展，环境问题日益严峻，能源供需矛盾日益突出［1-2］。为解决我国能源供需矛盾，需要逐步转变我国能源发展结构。吸收式制冷是一种节能环保新技术，因对我国低品位能源回收利用有重要作用而受到广泛关注［3］。迄今为止，H2O/NH3和LiBr/H2O仍是吸收式制冷循环中最常用的工质对［4］。但这两种工质对仍存在诸多缺陷：LiBr/H2O系统高温易腐蚀、低温易结晶，增加了设备的维护运行成本；H2O/NH3系统同样存在着毒性大、难分离等缺陷。

离子液体（ionic liquid，IL）是一种在室温状态下呈液态的熔融盐，可以作为吸收剂与甲醇等醇类制冷剂组成制冷循环工质对，有着广泛的应用前景［5］。目前国内外已有大量针对离子液体型制冷循环工质对筛选的研究，肯定了其理论价值和发展潜力［6-9］。

工质对混合溶液的气液相平衡性质是筛选离子液体型吸收式制冷工质对的关键，其饱和蒸气压作为工质传质特性的本征参数，亟需展开研究。对离子液体二元溶液气液相平衡性质的研究主要采用实验研究方法。张晓冬等［10］采用静态法实验测量了［BMIM］［DBP］分别与水、甲醇和乙醇组成的3组二元溶液的饱和蒸气压。Krlikowska等通过沸点法测量了4种溴化盐离子液体与水组成的二元溶液的饱和蒸气压数据。影响二元溶液饱和蒸气压的主要因素有组份与温度，不同的离子液体与制冷剂组合得到的二元溶液在气液相平衡性质上的表现各有差异，因此，为了筛选性能优良的吸收式制冷工质对，其气液相平衡性质值得进一步研究。

本文以［EMIM］BF4/CH3OH离子液体二元溶液为研究对象，对其气液相平衡性质展开研究。通过静态法实验测量［EMIM］BF4/CH3OH二元体系在不同状态点参数下的饱和蒸气压，提出了一种新的分子动力学模拟方法，模拟计算其饱和蒸气压并与实验所得数值进行对比分析，为后续制冷系统理论循环研究提供了基础数值和理论依据。

1 气液相平衡实验研究

1.1 实验材料

无水甲醇CH3OH（纯度＞99.999%，CAS No.67-56-1），购自阿拉丁试剂有限公司；1-乙基-3-甲基咪唑四氟硼酸盐［EMIM］BF4（纯度＞99.999%，CAS No.43314-16-3），购自上海成捷化学有限公司。本文实验所有试剂配制和称重过程均使用高精度分析天平，其误差范围为±0.1 mg。

1.2 实验装置与原理

本文实验采用静态法来测定［EMIM］BF4/CH3OH二元体系溶液的饱和蒸气压。其中静态法的基本测量原理是将待测［EMIM］BF4/CH3OH二元溶液充入内部真空状态的反应釜内，保证反应釜密闭良好的状态下，将其放入恒温油槽中进行加热，通过内部搅拌器对反应釜内的二元溶液不断搅拌，使其充分混合加速两相接触，从而可以更快达到气液相平衡状态，最后记录下反应釜中平衡状态下的压力值。

图1为［EMIM］BF4/CH3OH气液相平衡实验装置原理图。本文实验装置主要分为两个部分：恒温控制部分和气液相平衡反应与测定部分。恒温控制部分包括恒温油槽和离心叶片搅拌机；气液相平衡反应部分主要由高压反应釜、真空泵、压力表、电动叶片搅拌器等组成。其中核心部分高压反应釜容积为50 mL，釜身和釜盖由密封圈和六角螺丝固定，保证了实验的密封性。

2 ［EMIM］BF4/CH3OH分子动力学模拟

2.1 力场参数与势函数

图2为模拟［EMIM］BF4/CH3OH二元体系气液相平衡所使用的［EMIM］+阳离子、BF4-阴离子和CH3OH分子构型，基础构型由LigParGen［12］平台计算得出。

Sambasivarao和Acevedo［13］于2009年开发并验证了用于分子模拟的68种离子液体组合的OPLS-AA（optimized potentials for liquid simulations）力场参数。Doherty等［14］对OPLS-2009IL的力场参数进行了重新修订，并针对其他离子开发了新的OPLS-AA参数，还研究了±0.8 e的电荷缩放作为模拟极化和电荷转移效应的手段，在0.8尺度下，蒸发热、表面张力和自扩散系数均有显著改善。本文模拟中使用的 ［EMIM］+阳离子与BF4-阴离子的力场参数来自于Doherty等开发的全原子力场，离子液体［EMIM］BF4中的各原子使用缩放电荷，总电荷数为±0.8 e。甲醇CH3OH的力场参数采用Jorgensen等［15］开发的经典OPLS-AA力场。表1为［EMIM］BF4和CH3OH的非键参数，其中q为原子电荷，s为Lennard-Jones半径，ε为势阱深度（1 kcal/mol=4.184 kJ/mol）。

选用力场所对应的势函数如下：

Ebonds=∑ikr，iri－r0，i2，（1）

Eangles=∑ikθ，iθi－θ0，i2，（2）

Etorsion=∑iV1，i21+cos φ+V2，i21－cos 2φ+V3，i21+cos 3φ+V4，i21－cos 4φ ，（3）

Enonbond=∑i∑jgt;iqiqje2rij+4εijσijrij12－σijrij6fij，（4）

E=Ebonds+Eangles+Etorsion+Enonbond。（5）

E为离子液体体系总能量，包括谐波键伸缩项能量之和Ebonds、角弯曲项能量之和Eangles、二面角扭转项能量的余弦级数Etorsion和非键的相互作用Enonbond，其中Enonbond又包含库伦静电相互作用项和12-6 Lennard-Jones势项。式中，kr与kθ为各相互作用的力常数，r0和θ0为平衡键长值和平衡键角度值，r、θ、φ分别代表键长、键角和二面角值，V代表傅里叶系数，q代表部分原子电荷，σij为Lennard-Jones半径即相互作用的势能正好为零时的两体距离，rij为原子间距，εij为势阱深度。式（4）中对Lennard-Jones系数采用几何合并法则，即σij =（σiiσjj）1/2，εij =（εiiεjj）1/2；为了对分子内和分子间相互作用应用相同的参数，1，4分子内相互作用减半，即i， j为1和4时，fij取0.5，其余情况fij =1。

2.2 二元体系气液相平衡模拟过程

所有体系模拟均使用分子动力学模拟软件LAMMPS进行。体系模型的建立采用的是Packmol分子动力学模拟软件，针对不同的参数要求建立不同原子数量和不同体积模拟盒子，并实现原子在模型中的随机分布。使用Moltemplate软件为模型中的原子添加特定的力场参数，最终形成LAMMPS中运行的data文件。

分子动力学模拟中压强的计算可以通过维里方程计算：

P=NVkBT+13V∑i，jrij·Fij，（6）

其中，P为系统的压强，N为模型盒子中的原子数目，V是模拟系统的体积，kB是玻尔兹曼常数，T为系统的温度，rij为原子间距，Fij为原子间相互作用力。在每个时间步中计算粒子之间的相互作用力带入公式即可得到压强的对应估计值。

要在有限的资源下计算出误差最小的结果，前期对不同大小的体系进行了预模拟，最终选择模型盒子大小在x、y方向均大于80 ，z方向均大于200 ，其中液相区域在z方向长度不小于60 ，在气液相平衡过程中盒子内原子数目约3～4万。

图3为离子液体二元溶液气液相平衡模拟体系。由图3可知，液相区域位于模拟体系中部，上下两端在模拟开始前为添加的真空区域。模拟开始后，中部二元液体区域中的甲醇分子会开始向上下两端的真空区域移动，随着模拟的进行，逐渐达到气液相平衡状态，此时便形成了气相-液相-气相的混合体系，上下两端气相区域的压力便是［EMIM］BF4/CH3OH二元体系的饱和蒸气压。

首先通过预期要模拟的状态参数预估液相体系盒子的大小和内部要放置的离子液体和甲醇分子的数目。通过Packmol和Moltemplate软件来建立模型并添加力场参数。然后对所建液相的体系盒子进行能量最小化处理，让模型内部达到平衡状态使其合理化；选择使用NPT系综来对体系进行弛豫，控制温度T为模拟状态点温度，压力P为0保持真空状态，压力控制方向设置为z方向，并保持x、y方向盒子的尺寸不发生变化。弛豫运行时间约500 ps，体系的能量和体积达到稳定状态后记录下z方向盒子的上下尺寸。随后在液相体系盒子z方向上下分别添加相同尺寸的真空层，将拼接好的模型在NVT系综下运行约6 ns。整个模拟过程中气液界面较为清晰，气相和液相区域未发生较大宏观位移。取上下气相区域离气液界面较远处，沿垂直于z方向划分，取长40～50 的长方体区域作为甲醇压力测定区域。每隔500步输出一次原子轨迹数据，每隔1 000步输出模拟体系的温度、能量等基础数据和测定区域内甲醇的数量和压力数据。当区域平均压力稳定即判定此时体系达到了气液相平衡状态，压力即为当前温度下二元体系的饱和蒸气压。

图4为在温度T=298.15 K时，3种不同力场下，分子数目相同的纯甲醇体系气相压力随时间变化的模拟结果。由图可知，全原子力场OPLS-AA模拟值约为18.9 kPa，最接近甲醇对应温度下的实验测量值16.8 kPa［16］。联合原子力场OPLS-UA和TraPEE-UA［17］的模拟误差均大于OPLS-AA。而后进行多组不同温度下纯甲醇溶液体系的气液相平衡模拟，得到的结果与实验测量值相对误差均在11%～15%区间，最终采用最小二乘法最终确定压力修正系数为0.87。同时模拟结果也显示，模拟温度对饱和蒸气压的影响规律与实际情况符合较好，间接表明全原子力场OPLS-AA可以准确反映模拟体系的变化状况，证明了所选力场的合理性与可行性。

3 结果与讨论

3.1 实验结果与分析

图5为［EMIM］BF4/CH3OH二元体系的压力-温度-摩尔分数（P-T-x）图。由图可知，［EMIM］BF4/CH3OH二元溶液的饱和蒸气压随着甲醇的摩尔分数x和温度T的升高而增大。整体来看，随着甲醇摩尔分数的增加，离子液体二元溶液的饱和蒸气压呈现出先慢后快的上升趋势；当甲醇摩尔分数较低时，饱和蒸气压增长趋势较为缓慢；随着甲醇摩尔分数的不断增大，二元体系的饱和蒸气压的增长速度逐渐变快。出现这一现象的原因是，离子液体［EMIM］BF4的摩尔质量要远大于甲醇的摩尔质量，约为甲醇的6倍，当甲醇的摩尔分数为0.4时，其质量分数约为0.1。由于质量占比较少，故当甲醇摩尔分数小于0.4时，其饱和蒸气压上升趋势缓慢；随着摩尔分数增大，质量分数也不断增大，二元体系的饱和蒸气压上升趋势逐渐增快。随着温度的增加，离子液体［EMIM］BF4/CH3OH二元体系的饱和蒸气压也呈上升趋势，在温度较低时上升缓慢，上升趋势不明显；随着温度的升高，饱和蒸气压上升速度逐渐增加。

图6是在温度T=333.15 K条件下，［EMIM］BF4/CH3OH离子液体二元体系饱和蒸气压实验数据与［HMEA］［HCOO］/CH3OH溶液、［HDEA］［HCOO］/CH3OH溶液［18］与［MMIM］DMP/CH3OH溶液［19-20］蒸气压对比。当甲醇摩尔分数x在0.55～0.9内，［EMIM］BF4/CH3OH二元体系的饱和蒸气压远小于［HMEA］［HCOO］/CH3OH、［HDEA］［HCOO］/CH3OH溶液的蒸气压，略低于［MMIM］DMP/CH3OH溶液，其蒸气压值整体较［MMIM］［PMP］/CH3OH二元溶液平均低约21%。由此可以看出离子液体［EMIM］BF4对甲醇来说是较好的吸收剂，相较于其他传统离子液体工质对，［EMIM］BF4/CH3OH二元溶液有较为优秀的气液相平衡性质，具有替代吸收式制冷传统工质对，成为新型离子液体工质对的潜质。

3.2 模拟结果与讨论

图7展示了在温度T=348.15 K，甲醇摩尔分数x=0.9时，离子液体［EMIM］BF4/CH3OH二元体系气相压力随时间的变化过程。由图可知，0～1 ns区间体系的蒸气压力快速上升，自2.5 ns后，饱和蒸气压的模拟数值波动较小，逐渐趋于稳定，此时可认为离子液体二元体系达到气液相平衡状态。本章在模拟计算饱和蒸气压数值时，统计的是最后4～6 ns压强的平均值。

图8展示了在不同温度条件下，［EMIM］BF4/CH3OH二元溶液达到气液相平衡状态时饱和蒸气压在不同甲醇摩尔分数下实验结果与模拟值的对比，图8（a）、（b）、（c）和（d）分别为228.15、318.15、348.15和378.15 K时的变化趋势对比。由图可知，通过实验数据和模拟结果的对比，在不同温度工况下，［EMIM］BF4/CH3OH二元体系气液相平衡饱和蒸气压的模拟结果与实验数值有相同的数量级，并且饱和蒸气压模拟结果随甲醇摩尔分数的变化趋势与实验结果基本一致，均随摩尔分数的增大而增大；且当温度升高，饱和蒸气压也随之增大。模拟结果与实验值的误差最大不超过16%，且随温度的增加，相对误差下降，模拟值与实验值的契合度较好。总体来说，鉴于模拟所得结果与实验值有相同的数量级和变化趋势，该分子动力学模拟方法可以较为精确地计算出不同工况下［EMIM］BF4/CH3OH二元溶液的饱和蒸气压力值，定性分析离子液体二元体系的气液相平衡状态。为离子液体基础物性实验研究提供了数据补充辅助验证，也为节省实验资源，减少人工时间成本，探究新的气液相平衡饱和蒸气压的计算方法提供了新的思路。

4 总结

本文讨论了［EMIM］BF4/CH3OH离子液体二元体系制冷工质对的气液相平衡性质。通过静态法实验测得了不同温度和甲醇摩尔分数下［EMIM］BF4/CH3OH饱和蒸气压；提出了一种新的预测气液相平衡模拟方法，基于OPLS-AA全原子力场建模，通过分子动力学模拟同样工况下的二元溶液的饱和蒸气压，最后与实验值进行对比分析。主要结论如下：

（1）实验下［EMIM］BF4/CH3OH二元体系的饱和蒸气压随着温度T的升高和甲醇摩尔分数x的增大而升高。与其他醇类离子液体制冷工质对相比，［EMIM］BF4/CH3OH饱和蒸气压比［MMIM］DMP/CH3OH平均低约21%。

（2）模拟得到的数据与实验结果有相同的数量级和变化趋势，相对误差普遍低于8%，最大不超过16%。证明了此种分子动力学模拟的方法可以很好地预测离子液体二元体系的饱和蒸气压。

参考文献：

［1］曹莉萍， 周冯琦. 能源革命背景下中国能源系统转型的挑战与对策研究［J］. 中国环境管理， 2017， 9（5）： 84-89. DOI： 10.16868/j.cnki.1674-6252.2017.05.084.

［2］JOUHARA H， KHORDEHGAH N， ALMAHMOUD S， et al. Waste heat recovery technologies and applications［J］. Thermal Science and Engineering Progress， 2018， 6： 268-289. DOI： 10.1016/j.tsep.2018.04.017.

［3］郭少欣， 张杰， 罗昊敏. 低品位工业余热采集方法优化研究［J］. 节能， 2020， 39（3）： 95-99.

［4］WANG J Z， ZHENG D X， FAN L H， et al. Vapor pressure measurement for the water + 1， 3-dimethylimidazolium chloride system and 2， 2， 2-trifluoroethanol + 1-ethyl-3-methylimidazolium tetrafluoroborate system［J］. Journal of Chemical amp; Engineering Data， 2010， 55（6）： 2128-2132. DOI： 10.1021/je900738e.

［5］KHAMOOSHI M， PARHAM K， ATIKOL U. Overview of ionic liquids used as working fluids in absorption cycles［J］. Advances in Mechanical Engineering， 2013， 5： 620592. DOI： 10.1155/2013/620592.

［6］张芳芳， 陈更， 郑飞飞， 等. 离子液体工质对吸收制冷循环研究进展［J］. 现代化工， 2019， 39（2）： 36-40.

［7］MEYER T， SCHULZE U. Experimental data on mutual mass diffusivities of binary mixtures of ethanol and［EMIM］［DEP］［J］. Chemical Engineering amp; Technology， 2019， 42（7）： 1414-1420. DOI： 10.1002/ceat.201900027.

［8］FERWATI M S， AHMAD A M， TAKALKAR G D， et al. Energy and exergy analysis of parallel flow double effect H2O-［MMIM］［DMP］absorption refrigeration system for solar powered district cooling［J］. Case Studies in Thermal Engineering， 2021， 28： 101382. DOI： 10.1016/j.csite.2021.101382.

［9］de ARAJO H V， MASSUCHETTO L H P， do NASCIMENTO R B C， et al. Thermodynamic performance analysis of a single-effect absorption refrigeration system operating with water and 1-ethyl-3-methylimidazolium-based ionic liquids mixtures［J］. Applied Thermal Engineering， 2022， 201： 117761. DOI： 10.1016/j.applthermaleng.2021.117761.

［10］ZHANG X D， HU D P， ZHAO Z C. Measurement and prediction of vapor pressure for H2O+CH3OH/C2H5OH+［BMIM］［DBP］ternary working fluids［J］. Chinese Journal of Chemical Engineering， 2013， 21（8）： 886-893. DOI： 10.1016/s1004-9541（13）60544-3.

［11］KRLIKOWSKA M， PADUSZYSKI K， ZAWADZKI M. （Vapor+ liquid） phase equilibria of an aqueous solution of bromide-based ionic liquids-measurements， correlations andapplication to absorption cycles［J］. Fluid Phase Equilibria， 2019， 494： 201-211. DOI： 10.1016/j.fluid.2019.05.003.

［12］DODDA L S， CABEZA DE VACA I， TIRADO-RIVES J， et al. LigParGen web server： an automatic OPLS-AA parameter generator for organic ligands［J］. Nucleic Acids Research， 2017， 45（W1）： W331-W336. DOI： 10.1093/nar/gkx312.

［13］SAMBASIVARAO S V， ACEVEDO O. Development of OPLS-AA force field parameters for 68 unique ionic liquids［J］. Journal of Chemical Theory and Computation， 2009， 5（4）： 1038-1050. DOI： 10.1021/ct900009a.

［14］DOHERTY B， ZHONG X， GATHIAKA S， et al. Revisiting OPLS force field parameters for ionic liquid simulations［J］. Journal of Chemical Theory and Computation， 2017， 13（12）： 6131-6145. DOI： 10.1021/acs.jctc.7b00520.

［15］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. DOI： 10.1021/ja9621760.

［16］VERCHER E， LLOPIS F J， GONZLEZ-ALFARO V， et al. Volumetric properties， viscosities and refractive indices of binary liquid mixtures of tetrafluoroborate-based ionic liquids with methanol at several temperatures［J］. The Journal of Chemical Thermodynamics， 2015， 90： 174-184.

［17］JARADAT A， AL-SALMAN R， OBEIDAT A. Molecular dynamics simulation of vapor-liquid equilibrium in 1-alkanol unary systems： a study of surface tension， density， and vapor pressure of TraPPE-UA force field［J］. Fluid Phase Equilibria， 2024， 577： 113967. DOI： 10.1016/j.fluid.2023.113967.

［18］LI X M， SHEN C， LI C X. Effect of mono- and di-ethanolammonium formate ionic liquids on the volatility of water， ethanol， and methanol［J］. Chinese Journal of Chemical Engineering， 2013， 21（10）： 1162-1171. DOI： 10.1016/s1004-9541（13）60582-0.

［19］ZHAO J， JIANG X C， LI C X， et al. Vapor pressure measurement for binary and ternary systems containing a phosphoric ionic liquid［J］. Fluid Phase Equilibria， 2006， 247（1/2）： 190-198. DOI： 10.1016/j.fluid.2006.07.007.

［20］赵杰， 梁世强， 陈坚， 等. 离子液体1， 3-二甲基咪唑磷酸二甲酯盐-甲醇浓溶液气液平衡［J］. 化学工程， 2010， 38（3）： 52-56.