粒子插入法計算一氧化碳在乙醇中的溶解度

宣愛國，侯俊燁，吳元欣，閆志國，嚴亮

1．武漢工程大學化工與制藥學院，湖北 武漢 430074；2．綠色化工過程省部共建教育部重點實驗室（武漢工程大學），湖北 武漢 430074

粒子插入法計算一氧化碳在乙醇中的溶解度

宣愛國1，2，侯俊燁1，2，吳元欣1，2，閆志國1，2，嚴亮1，2

1．武漢工程大學化工與制藥學院，湖北 武漢 430074；2．綠色化工過程省部共建教育部重點實驗室（武漢工程大學），湖北 武漢 430074

采用等溫等壓（NpT）系綜結合Widom粒子插入法計算了298、323 K下1～6 MPa壓力范圍內一氧化碳在乙醇中的溶解度，探討了壓力、溫度對一氧化碳溶解度的影響.模擬過程中，一氧化碳采用DREIDING力場，乙醇采用OPLS-ua力場，模擬結果與文獻實驗值吻合良好.結果表明溫度恒定時一氧化碳在乙醇中的溶解度隨著壓力的增加而增大，且溶液中一氧化碳的摩爾分數與壓力的關系滿足亨利定律；壓力相同時，一氧化碳的溶解度在298、323 K下的差別很小，說明在一定范圍內，一氧化碳在乙醇中的溶解度受溫度影響不大.

粒子插入法；蒙特卡羅模擬；一氧化碳；乙醇；溶解度

0 引言

氣液相平衡數據在化工生產與理論研究中有著重要作用，但其實驗測定較為困難.近年來，隨著計算機技術的高速發展，分子模擬已經成為除實驗測定、模型計算外的一種全新而獨特的研究手段［1］.分子模擬建立在統計熱力學與分子力學理論基礎之上，在計算熱力學數據的過程中只依賴于描述原子間相互作用的分子力學力場，直接從微觀狀態分布來求解宏觀量，在計算相平衡方面具有巨大優勢.

計算體系氣液相平衡常用且計算效果較好的方法包括Widom粒子插入法［2］與Gibbs系綜蒙特卡羅（Gibbs ensemble Monte Carlo，GEMC）［3-4］方法.前者通過計算插入的虛擬粒子與體系粒子之間作用的超額化學勢來計算氣體的溶解度或亨利系數，而后者則分別模擬氣、液兩相，通過加入粒子在兩相間的交換移動達到氣液相平衡.與GEMC方法相比，粒子插入法在實施過程中更為便捷，且在計算亨利系數方面更為直接.在2004年的第二屆工業流體模擬競賽（the Second Industrial Fluid Property Simulation Challenge，2nd IFPSC）［5］中，Widom方法被用于計算N2、O2、CO2和CH4四種氣體在乙醇中的亨利系數.此后，隨著力場的不斷完善，粒子插入法逐漸被用于多種體系［6-9］氣液相平衡的計算，且都取得了令人滿意的結果.

一氧化碳、乙醇是羰基化法制備碳酸酯［10-11］的重要原料和溶劑.研究碳酸酯相關合成體系的相平衡可為產物分離方法的選擇、分離設備平衡級數的計算及尺寸設計等提供依據.本文采用NpT系綜結合粒子插入法，通過Monte Carlo算法計算了不同條件下一氧化碳在乙醇中的溶解度及亨利系數，將計算結果與文獻值進行比較，發現粒子插入法可以取得較好的結果.

1 模擬方法

1．1 粒子插入法計算原理

一氧化碳在乙醇中的溶解度較低，溶解了一氧化碳的乙醇可視為無限稀釋溶液.溶液中的相互作用可直接視為乙醇分子與一氧化碳分子之間的相互作用，而不考慮液相中氣體分子之間的相互作用.故可采用Henry定律來描述一氧化碳-乙醇體系的相平衡關系.

根據溶液熱力學，亨利系數H2，1可由溶質在無限稀釋溶液中的超額化學勢μ2ex，∞來進行計算［12］：

式中，kB是Boltzmann因子，ρ1為溶劑密度，T為體系溫度.由此可知，計算亨利系數的關鍵在于準確地計算溶質的超額化學勢μ2ex，∞.采用Widom粒子插入法計算化學勢的基本步驟為：在模擬進行過程中凍結一模擬構型，在模擬盒中隨機插入一個粒子，計算此粒子與體系中其他粒子相互作用的勢能之和，然后將插入粒子撤出，恢復原來的模擬構型，繼續進行模擬.采用NpT系綜模擬氣液相平衡時，插入的溶質粒子的超額化學勢μ2ex可表示為［13］：

表示對體系能量進行等溫等壓系綜平均.當體系足夠大時，即得到溶質在無限稀釋溶液下的超額化學勢μ2ex，∞.

1．2 分子力場

模擬過程中，溶劑分子乙醇采用OPLS-ua力場，乙醇分子由甲基、亞甲基及羥基氫、氧原子四個作用位點構成，電荷分布在亞甲基及羥基氫、氧原子上，詳細力場參數見文獻［14］.溶質分子一氧化碳采用DREIDING力場，碳、氧原子上分別負載有0．020 3 e的正、負電荷［15］，其余力場參數見文獻［16］.

體系中分子i和j之間的相互作用勢能Eij可表示為：

式中，ELJ為Lennard-Jones勢能，Echarge為電荷作用部分勢能，分別表示為：

式中，na，nb分別為分子a和b的作用位點數；qi，qj分別為作用位點i及j的電荷量；rij為作用位點i與j之間的距離；ε0為真空介電常數.

不同作用位點之間的相互作用采用Geometric混合規則：

式中，εi，σi分別代表作用位點i的LJ勢能能量參數與尺寸參數.

2 模擬結果與討論

采用NpT系綜結合Widom粒子插入法計算一氧化碳-乙醇體系的氣液相平衡.模擬盒采用周期性邊界條件，初始邊長設為30 nm，模擬溫度分別為298、323 K，壓力范圍為1～6 MPa.采用球形割去法進行勢能截斷，截斷半徑為1 nm，截斷距離外的范德華作用通過長程矯正法進行矯正.靜電作用通過Ewald加和法進行計算.模擬過程共進行300 000次Monte Carlo（MC）循環，其中，前250 000次為平衡階段，后50 000次為統計計算階段.模擬過程的MC移動及其發生的概率分別為：體積波動，1%；偏倚構型增長，32%；質心平移，34%；繞質心的旋轉，33%.各移動的接受概率均設為0．5.本文的模擬工作均使用Towhee-7．0．2［17］軟件包完成.

2．1 壓力對氣體溶解度的影響

圖1（a）、圖1（b）分別表示溫度298、323 K時一氧化碳在乙醇中的溶解度的模擬結果與文獻結果［18］之間的比較.

圖1 一氧化碳在乙醇中的溶解度隨壓力的變化Fig．1Effect of pressure on solubility of carbon monoxide in ethanol

從圖中可知，溫度一定時，一氧化碳在乙醇中的溶解度隨壓力的增大而增加，且溶質在溶液中的摩爾分數與壓力呈很好的線性關系，滿足亨利定律.在模擬計算采用的壓力范圍內，溶解度的模擬值與文獻值非常接近.在1～4 MPa壓力范圍內，溶解度模擬值總是略大于文獻實驗值，但總體均在可接受范圍內：298 K下，模擬值與文獻值的相對偏差在4．2%～12．9%之間；323 K下，二者的相對偏差在5．5%～18．3%之間.壓力大于4 MPa后，由圖中擬合模型可知，溶解度的模擬值與模型值較1～4 MPa范圍內的值相比更為接近，298 K下模擬值與模型值的相對偏差范圍降至2．1%～3．3%；323 K下二者相對偏差范圍則降至0．47%～3．9%.

2．2 溫度對氣體溶解度的影響

圖2（a）、圖2（b）分別為298～323 K范圍內一氧化碳在乙醇中溶解度的文獻結果與模擬結果.實驗結果表明，壓力一定時，一氧化碳在乙醇中的溶解度幾乎不隨溫度的變化而發生改變.壓力一定時，模擬得到的溶解度值在298～323 K的溫度范圍非常相近，最大偏差為0．28%，且兩個溫度下回歸所得的溶解度模型與文獻模型有著很相似的趨勢，模擬結果證實了實驗結果的準確性，同時也證明使用粒子插入法可以預測一定范圍內溫度對一氧化碳-乙醇體系氣液相平衡的影響.

圖3為298 K、323 K下一氧化碳-乙醇體系亨利系數的模擬擬合值與文獻值之間的對比.

圖3 一氧化碳在乙醇中的亨利系數模擬值與文獻值之間的對比Fig．3Comparison of Henry’s law constants between simulated and experimental data

模擬與實驗結果均表明一氧化碳-乙醇體系的亨利系數在298 K與323 K下幾乎可視為相等；且模擬所得的亨利系數與文獻結果非常接近，在298、323 K下，二者的相對偏差分別為2．9%、2．1%，說明采用與NpT結合的粒子插入法模擬一氧化碳-乙醇體系亨利系數是可取的，其結果可用于預測本研究范圍內其他壓力下一氧化碳在乙醇中的溶解度.

3 結語

本研究采用NpT系綜結合Widom粒子插入法計算了298、323 K下1～6 MPa范圍內一氧化碳在乙醇中的溶解度，計算結果及擬合模型均與文獻結果吻合良好，證實了分子模擬可有效地檢驗與預測實驗數據.同時，模擬結果與文獻結果均表明：溫度恒定時，一氧化碳在乙醇中的溶解度隨壓力的增加而增大，且二者之間的關系符合亨利定律；一定范圍內，溫度對一氧化碳溶解度的影響并不明顯，一氧化碳在乙醇中的亨利系數幾乎不隨著溫度的變化而改變.

致謝

感謝國家自然科學基金委和乙烯工程下游產品開發及過程強化湖北省協同創新項目組提供的資助.

符號說明

E-勢能，J；Hi-溶液中組分i的亨利系數，MPa；kBBoltzmann常數，J·K-1；N-模擬分子總數；n-分子作用位點數；p-壓力，MPa；q-作用位點的電荷量，e；R-氣體常數，J· mol-1·K-1；rij-作用位點i與j之間的距離，nm；T-絕對溫度，K；V-體積，m3；xi-溶液中組分i的摩爾分數．

希臘字母：ε-分子能量參數，J；ε0-真空介電常數，F·m-1；μi－組分i的化學勢；ρ－數密度，m－3；σ－分子尺寸參數，nm．

上角標：ex－超額熱力學性質；∞－無限稀釋值．

［1］朱宇，陸小華，丁皓，等．分子模擬與化學工程［J］．化工學報，2004，55（8）：1213－1223．

ZHU Yu，LU Xiao－hua，DING Hao，et al．Molecular simulation in chemical engineering［J］．CIESC Journal，2004，55（8）：1213－1223．（in Chinese）

［2］WIDOM B．Some topics in theory of fluids［J］．J．Chem Phys．，1963，39（11）：2808－2812．

［3］PANAGIOTOPOULOS A Z．Direct determination of phase coexistence properties of fluids by Monte Carlo simulationinanewensemble［J］．Molecular Physics，1987，61（4）：813－826．

［4］PANAGIOTOPOULOS A Z，QUIRKE N，STAPLETON M，et al．Phase equilibria by simulation in the Gibbs ensemble alternative derivation，generalization and application to mixture and membrane equilibria［J］．Molecular Physics，1988，63（4）：527－545．

［5］CASE F，CHAKA A，FRIEND D G，et al．The second industrial fluid properties simulation challenge［J］．Fluid Phase Equilibria，2005，236（1）：1－14．

［6］吳曉萍，劉志平，汪文川．分子模擬研究氣體在室溫離子液體中的溶解度［J］．物理化學學報，2005，21（10）：1138－1142．

WU Xiao－ping，LIU Zhi－ping，WANG Wen－chuan．Molecular dynamics simulation of gas solubility in room temperature ionic liquids［J］．Acta Phys．－Chim．Sin．，2005，21（10）：1138－1142．（in Chinese）

［7］WU Chuan－jie，LI Xiao－feng，DAI Jian－xing，et al．Prediction of Henry’s law constants of small gas molecules in liquid ethylene oxide and ethanol using forcefieldmethods［J］．FluidPhaseEquilibria，2005，236（1）：66－77．

［8］HUANG Yow－lin，MIROSHNICHENKO S，HASSE H，et al．Henry’s law constant from molecular simulation：a systematic study of 95 systems［J］．International Journal of Thermophysics，2009，30（6）：1791－1810．

［9］潘楊．TDI相關合成體系相平衡的測定與分子模擬［D］．武漢：武漢工程大學，2013．

PAN Yang．Determination and molecular simulation of TDI related synthesis system［D］．Wuhan：Wuhan Institute of Technology，2013．（in Chinese）

［10］田琦峰，張光緒，吳廣文，等．碳酸二苯酯的合成研究進展［J］．精細石油化工進展，2001，2（9）：37－41．

TIAN Qi－feng，ZHANG Guang－xu，WU Guangwen，etal．Researchprogressinsynthesisof diphenyl carbonate［J］．Advances in Fine Fetrochemicals，2001，2（9）：37－41．（in Chinese）

［11］劉定華，徐明，陳振松，等．銅鹽催化乙醇氧化羰基化反應性能的研究［J］．高校化學工程學報，2013，27（3）：425－430．

LIU Ding－hua，XU Ming，CHEN Zhen－song，et al．Study on reaction performance of copper salt in oxidative carbonylation of ethanol［J］．Journal of ChemicalEngineeringofChineseUniversities，2013，27（3）：425－430．（in Chinese）

［12］SHING K S，GUBBINS K E，LUCAS K．Henry constants in non－ideal fluid mixtures：Computer simulation and theory［J］．Molecular Physics，1988，65（5）：1235－1252．

［13］ALLEN M P，TILDESLEY D J．Computer simulationsofliquids［M］．Oxford：ClarendonPress，1987：2－240．

［14］JORGENSEN W L．Optimized intermolecular potentialfunctionsforliquidalcohols［J］．J．Phys．Chem．，1986，90（7）：1276－1284．

［15］SWEATMAN M B，QUIRKE N．Modelling gas adsorption in slit－pores using Monte Carlo simulation［J］．Molecular Simulation，2001，27（5－6）：295－321．

［16］MAYO S L，OLAFSON B D，GODDARD III W A．DREIDING：A genericforcefieldformolecular simulations［J］．Phys Chem．，1990，94（26）：8897－8909．

［17］MARTIN M G.MCCCS T owhee:a tool for Monte Carlo molecular simulation[J].Mloecular Simulation,2013,39(14/15):1212-1222.

［18］TONNER S P，WAINWRIGHT M S，TRIMM D L．Solubility of carbon monoxide in alcohols［J］．Journal of Chemical and Engineering Data，1983，28（1）：59－61．

Solubility simulation of carbon monoxide in ethanol using Widom insertion method

XUAN Ai－guo1，2，HOU Jun－ye1，2，WU Yuan－xin1，2，YAN Zhi－guo1，2，YAN Liang1，2
1．School of Chemical Engineering＆Pharmacy，Wuhan Institute of Technology，Wuhan 430074，China；2．Key Laboratory of Green Chemical Process（Wuhan Institute of Technology），Ministry of Education，Wuhan 430074，China

Widom insertion method was used to calculate the solubility of carbon monoxide in ethanol．The simulation was performed under isothermal－isobaric conditions of 298，323 K and 1－6 MPa．Effects of pressure and temperature on solubility of carbon monoxide in ethanol were discussed．During the simulation，the universal force field DREIDING was chosen for solute carbon monoxide，while the unitedatom force field OPLS－ua was chosen for solvent ethanol．The simulated data show good agreement with experimental results．The results indicate that solubility of carbon monoxide in ethanol increases with increasing pressure when the temperature is specified；the relationship between mole fraction of carbon monoxide and pressure can be described by Henry’s Law；however，the solubility of carbon monoxide obtained at 298 K and 323 K shows little difference when the pressure is the same，which means the dissolution of carbon monoxide in ethanol is not sensitive to temperature．

Widom insertion method；Monte Carlo simulation；carbon monoxide；ethanol；solubility

TQ02

A

10.3969/j.issn.1674-2869.2014.012.003

1674－2869（2014）012－0011－05

本文編輯：張瑞

2014-11-08

宣愛國（1957-），女，湖北松滋人，教授，博士，碩士研究生導師．研究方向：化工過程模擬與優化．