超臨界二氧化碳微乳反膠束體系的分子動力學模擬進展*

銀建中 王 松 徐君臣 喻 文 徐琴琴

（大連理工大學化工機械學院)

超臨界二氧化碳微乳反膠束體系的分子動力學模擬進展*

銀建中**王 松 徐君臣 喻 文 徐琴琴

（大連理工大學化工機械學院)

超臨界CO2微乳液/反膠束體系,是在表面活性劑作用下形成的CO2包水 (W/C)或者水包CO2(C/W)型且熱力學性質穩定的透明液體。通過改善超臨界CO2的性質,可以極大地拓展超臨界流體技術的應用領域。其熱力學性質、傳遞性質以及相行為是應用的基礎。對此復雜體系由實驗量測所獲得的信息遠不能滿足實際需要,但分子動力學模擬則可以達此目的。總結了有關分子動力學模擬在超臨界CO2微乳液/反膠束體系研究中的進展,重點分析了模擬實施過程中應該注意的關鍵問題,對于分子力場構建、單元體設計、平衡控制等熱點問題給予分析和闡述,并就應該致力于研究的核心方向進行了初步探討。

超臨界二氧化碳 反膠束 微乳 分子動力學模擬 物化性質

0 前言

綠色化學技術是過程工業實現節能、減排、降耗并提升技術經濟性的有效途徑,也是國民經濟實現可持續發展,以高科技改造傳統產業的必然選擇。 綠色化學依賴于綠色溶劑,而超臨界流體(supercritical fluid,SCF)、室溫離子液體 (room temperature ionic liquids,RTILs)、二氧化碳膨脹液體 (CO2expanded liquids,CXLs)和水是目前公認的綠色溶劑,自然成為研究和關注的熱點課題。其中,超臨界CO2(supercritical carbon dioxide,scCO2)受重視的程度尤為突出[1]。 特別是基于SCF發展起來的scCO2微乳液或反膠束體系,更是近年來SCF技術領域備受矚目的優先探索課題[2-3]。 scCO2微乳液/反膠束是在表面活性劑作用下形成的CO2包水 (W/C)或者水包CO2(C/W)型且熱力學性質穩定的透明液體,借助于改善scCO2的性質,極大地拓展了SCF技術的實際應用領域。其熱力學性質、傳遞性質以及相行為是應用的基礎。對此復雜體系由實驗量測所獲得的信息遠不能滿足實際需要,運用分子動力學模擬 (molecular dynamics simulation,MD)則可以達此目的[4-8]。 因此,近年來很多學者均在致力于探索將MD用于scCO2微乳液/反膠束體系性質的研究[9-11]。

1 scCO2微乳/反膠束

scCO2微乳/反膠束的結構如圖1所示。當以scCO2為連續相時,水在表面活性劑作用下分散其中形成極性微水環境,因為聚團尺寸小于100 nm,故溶液是澄清透明的[12]。 微水環境的存在使得scCO2可以溶解極性化合物、生物大分子等,極大地拓展了scCO2的應用范圍。目前,學術界趨向于不對scCO2微乳或反膠束進行嚴格區分,故本文按照模擬類文獻習慣將其統稱為scCO2反膠束(reverse micelles,RMs)。

圖1 scCO2RMs結構 [12]

1987年,美國太平洋西北國家實驗室 (PNNL)的Smith教授等率先以琥珀酸二酯磺酸鈉 (AOT)為表面活性劑,在超臨界流體中觀察到RMs的形成[13]。 1996年,美國德州大學奧斯汀分校的Johnston教授等以氟化聚醚碳酸銨 (PFPE)為表面活性劑,首次成功制備了scCO2RMs，并萃取出分子量為6 700且具有生物活性的牛血清蛋白,開創了scCO2RMs研究的新局面[14]。 2000年,Nature報道美國匹茨堡大學的Beckman教授[15]等合成了非氟表面活性劑聚醚碳酸鹽,其突出優點是能在中等壓力 (8～16 MPa)下溶于scCO2，當壓力為17～35 MPa時,在水和scCO2的混合物中加入1% （體積百分比）該聚合物就可以形成穩定的微乳液。該成果極大地帶動了可溶于scCO2表面活性劑的合成以及對W/C型微乳液熱力學性質的深入探究工作。近年來,英國布里斯托爾大學Eastoe教授的研究最為活躍,通過設計合成了可溶于scCO2的表面活性劑,推動了scCO2RMs的應用[16-17]。

2 MD模擬scCO2RMs體系

2.1 MD模擬RMs體系

用MD方法模擬RMs的優勢在于獲得其微觀結構和自組裝動力學特性等信息,此外,還可以藉以研究RMs的形狀和局部熱力學性質等。 研究RMs的外形特征時主要用到回轉半徑、主慣性矩和偏心率等參數。其中，回轉半徑公式如式（1）所示：

式中Rg——回轉半徑;

mi——距離RMs質心 （r0） 為ri的i原子的

質量。

徑向分布函數 (radial distribution function,RDF)/對相關函數 (pair correlation function)則主要用于考察體系的微觀結構,并可作為研究自組裝動力學的重要依據,如式 (2)所示:

式中 ρa(r)——以粒子a為參照,與距離r相關的其它粒子密度值;

X(r)——與參照粒子距離為r、厚度為Δr殼層內其它粒子的物理量,常用質量或數量;

＜…＞——表示系統平均。

均方位移 (mean square displacement,MSD)則可用于研究分子的擴散特性,如式 （3）所示:

式中r(0)——初始時刻粒子的坐標;

r(t)——經過時間t后粒子的坐標;

N——需要計算的粒子的總數量;

Nt——0到t時刻之間的時間步。

從發展歷程看,用MD模擬RMs經歷了由半定量模擬到精確模擬、驗證性模擬到預測性模擬的發展過程[18]， 其適用的體系也越來越廣泛。 Brown等[19]最早嘗試用MD方法模擬RMs,用三個Lernnard-Jones作用位點模型 (頭基、尾基和鈉離子)描述表面活性劑分子AOT,用單位點 (singlesite)模型模擬主體相正辛烷,將水分子裝入表面活性劑分子形成的球形殼層中,總模擬時間僅為163.5 ps。 由于硬件條件所限,盡管模擬時間很短,體系可能遠未達到穩定的平衡狀態,但這種嘗試為MD用于RMs研究開了個好頭。 Linse等[20-21]在真空環境下創建一個內部包含水分子并加載電荷的剛性球形殼層(固定位置、固定結構),用以研究殼層內的電場梯度。Faeder等[22-23]也曾使用過此模型。與Brown[19]的模型相似,Linse等的力場模型無法用于研究RMs的動力學性質,由于模型過于簡化,對RMs結構特性的研究結論也有待進一步完善。

早期,用MD模擬RMs體系時,采用在模擬初始階段將極性分散相 (如水等)預組裝 (preassembled)到表面活性劑球形殼層內的處理方法,而且受到一定的關注。其基本思想是假設預組裝的初始結構和模擬體系達到平衡時的結構相近,顯然,如何預知其平衡結構顯得尤為重要,而這種預知只能依賴于實驗研究,如RMs的特征尺寸、水/表面活性劑分子摩爾比W0等[24-25]。 采用預組裝初始結構可以極大地縮短模擬時間,但不利于研究反膠束的自組裝行為。

2.2 MD模擬scCO2RMs的發展

相比通常的RMs體系而言,scCO2RMs體系因高壓和超臨界相態的復雜性,由實驗獲得結構信息和物化性質參數的難度更大。用MD模擬scCO2RMs系統,關鍵問題在于模擬體系和表面活性劑的選擇等，需要考慮單個表面活性劑分子作用位點數、表面活性劑在scCO2中的溶解度等。 例如,多數AOT類表面活性劑在scCO2中的溶解度很小[26]，需要考慮添加助表面活性劑 (如甲醇、乙醇等)。用簡化的模型描述原子數較多的表面活性劑分子,可能削弱MD模擬在分析微觀結構和自組裝動力學方面的優勢。因此,模擬體系的選擇要綜合考慮模擬精度和計算效率,實現合理匹配。

相比表面活性劑,CO2和水分子已有成熟的力場模型可以使用。 scCO2的分子結構并非嚴格直線[27]， 相比二位點模型[28]， 采用三位點EPM/EPM2模型[29]會得到更高的模擬精度,計算量也將相應增加。目前水分子模型主要包括TIP3P、TIP4P、TIP5P、SPC/E、SSD等[30-34]， 其中又以三位點TIP3P和SPC/E模型應用最廣。進行MD模擬時,因為鍵伸縮和鍵角振動頻率較高,計算時不利于增大積分步長,加之鍵伸縮和鍵角振動對體系總勢能的貢獻不大,所以對復雜體系的MD模擬,應考慮在滿足計算精度的前提下,通過約束(constraints)降低分子內部自由度來提高計算效率,目前常用的約束法則有SHAKE和RATTLE等[35-36]。

美國田納西大學的Cui和橡樹嶺國家實驗室的Cummings等采用預組裝初始結構模擬了以(C7H15)(C7F15)CHSO4-Na＋為表面活性劑形成的scCO2RMs體系[24,37-39]。 Cui和Cummings等[24]通過RDF曲線研究反膠束的微觀結構,探討模擬過程中膠束的自組裝動力學,并通過對比表面活性劑尾端烴鏈和氟化烴鏈上反式鍵 (trans-bonds)的數量,表明氟化烴鏈與CO2的親和力更強。為了提高計算效率,他們選擇了一個水和表面活性劑濃度明顯高于實驗值[40]的體系 (W0=～4.4,共2 614個分子)。 雖然在1 ns的平衡模擬過程中觀察到明顯的團聚傾向,但在模擬結束時仍未獲得穩定的單一聚團結構,其模擬結果的意義在于證明了從原子水平用MD模擬scCO2RMs的可行性。

在此基礎上,他們又以同樣方法模擬了水和表面活性劑濃度與實驗值[40]相似的體系 (W0=～35,共14 008個分子)[37]， 并對比了預組裝初始結構和隨機初始結構對模擬結果的影響,在大約1 ns平衡模擬結束時均未得到單一的聚團結構。通過對如圖2所示的微團數量-時間曲線分析,認為初始結構的形式對于長時間后體系的動力學作用 (long-time dynamics)不產生影響,并推論雜化表面活性劑(dichain surfactant)在scCO2中形成Eastoe實驗[40]中的穩定單一聚團結構需要10-7s級的模擬時間。

通過對比,Cui和Cummings等認為含水量W0是影響RMs外形和水核熱力學性質的重要因素之一[38]。如圖3所示,他們用與水分子偶極矩矢量關聯的角分布函數 (angle distribution function，ADF)曲線考察W0較大體系中RMs水核內水分子的空間指向,發現界面區內水分子指向受離子電場的強烈影響,而水核中心區的水分子空間指向更多地體現出水分子間氫鍵網絡的作用。

圖2 水/表面活性劑聚團數量隨時間的變化曲線(W0=～35)[37]

圖3 水核中水分子間偶極矩-偶極矩矢量角度分布函數曲線 [38]

圖4 水與表面活性劑分子結合率隨時間變化曲線（W0=～4.4）[38]

如圖4所示,通過觀察水分子與表面活性劑分子結合率隨時間變化曲線,Cui和Cummings等提出雜化離子表面活性劑在scCO2中形成RMs的兩個步驟:(1)陰陽離子對快速電離 (T0,10-10s級);(2)表面活性劑分子在氫鍵作用下緩慢聚集成RMs(T1),并且T1明顯大于T0[39]。 這種團聚機制與Eicke等[41]的結論一致。Cui和Cummings等是最早采用MD方法模擬雜化表面活性劑形成scCO2RMs體系的研究小組,他們的工作體現了scCO2RMs模擬研究的關鍵問題,并通過模擬驗證了文獻對scCO2RMs團聚過程的推測。

Senapati和Berkowitz等[42]用小角中子散射(SANS)和MD方法考察了雜化磷酸鈉表面活性劑[(C6F13(CH2)2O)(C6F13(CH2)2O)PO2-Na＋,DiF8)] 形成的scCO2RMs體系,并對實驗測試和模擬計算結果進行定量比較。 模擬過程中采用聯合原子模型描述表面活性劑分子，并使用了單位點 (single-site)CO2模型[43]， 體系中各物質的摩爾濃度與實驗值相同。對膠束形狀進行了分析,第一次定量證實scCO2RMs并非嚴格球形，模擬的水核半徑與SANS實驗值[44]吻合較好，表面活性劑頭基表面積明顯小于AOT類，和Eastoe關系式計算結果[45]一致。結構分析表明，水分子與Na＋之間存在很強的作用，且很多磷原子的水合層之間有重疊。此外，Berkowitz等用MD方法對PFPE和聚醚碳酸銨(CH3(OCH2CH(CH3))3OCH2COO-NH4＋,PE) 為表面活性劑的scCO2RMs體系進行研究[46-47]， 采用全原子力場參數和預組裝初始結構的模擬方法,分析了RMs尺寸、形狀、頭基表面積、尾基構型和疏水尾基在scCO2中溶解度等性質,結果與實驗吻合良好。 分析氟原子在PFPE尾端鏈堆積 (packing of tails)過程中的作用時,發現氟化表面活性劑尾端鏈的二面角分布范圍比氫化表面活性劑更廣,即氟化表面活性劑尾端鏈更易變形。但模擬結束時卻發現，采用氟化和氫化表面活性劑所形成的RMs分別有12%和32%的水核表面積暴露于scCO2中。Klein等[48]也得到了類似的模擬結果。

開發scCO2RMs目的之一在于用極性水核選擇性增溶生物大分子或極性化合物。Senapati等用MD方法模擬scCO2RMs增溶蛋白質的四元體系[11]。 模擬單元內分子均隨機放置,表面活性劑為PFPE,蛋白質分子為色氨酸籠 (Trp-Cage)人造蛋白,CO2分子仍以單位點模型描述。如圖5所示,體系總分子數近7 000個,在50 ns平衡模擬過程中觀察到蛋白質模型自發進入scCO2RMs水核中,并保持原生狀態 (native state),模擬結果與實驗結果吻合較好。

近年來,綠色化學研究進展很快。其中一個趨勢就是將scCO2和RTILs相結合,充分利用scCO2良好的輸運特性、RTILs優越的溶劑特性和物化性質的穩定性,在化學反應、材料制備、分離工程等諸多領域展現出了廣闊的應用前景和商業價值。2010年,Senapati等模擬了scCO2為連續相、RTILs為分散相的scCO2RMs體系,首次基于原子水平通過動力學模擬觀察到IL-in-scCO2(即二氧化碳包離子液體型微乳,I/C微乳)自組裝過程及RMs的微觀結構[10]。 采用優化的OPLS-AA/AMBER力場參數描述了表面活性劑正乙基全氟辛基磺酸胺(C2H5NHSO2C8F17,N-EtFOSA)和RTILs。 通過對結構和能量性質的分析,解釋了四甲基醋酸胍 ([TMG][Ac])scCO2RMs比1-丁基-3-甲基咪唑六氟化磷([Bmim] [PF6])形成的RMs更穩定的原因。徑向分布函數圖顯示，表面活性劑頭基與RTILs陰離子之間的相互作用強于其與陽離子間的相互作用。能量分析表明，RTILs陰離子與表面活性劑頭基之間存在較強的氫鍵作用。Senapati等[10]采用能量判據分析氫鍵作用,優點在于可以避開分子具體構型,而僅考慮分子間能量是否小于某一設定值 (判據)。Ladanyi等[49]在分析RMs水核內的氫鍵網絡時采用了構型判據,他們認為當水分子間O-O距離小于3.5×10-10m且O-O-H鍵角小于π/6時存在氫鍵作用,并以此為依據計算水核內氫鍵數量。

如圖6所示,楊曉寧等[9]模擬了水/AOK(即5,5-二甲基-4-酮-己氧基磺酸鈉琥珀酸酯)/scCO2形成的RMs體系,第一次提供了全原子水平上CO2包水型膠束的結構圖。創建單元時,將表面活性劑分子規則排布,水和CO2分子隨機排布，這種處理避免了預組裝初始結構對模擬結果的影響。還討論了在保證模擬有效性的同時如何縮短平衡模擬時間。這些結果對今后開展此類工作具有一定的借鑒作用。模擬過程中,他們以量子力學改進的CHARMM 27力場參數描述AOK的尾端雙鏈 (Berkowitz等[42]在確定表面活性劑尾端鏈上各偽原子電荷時也使用了量子力學方法)。隨著HyperChem等商業化量子化學計算軟件的出現和不斷完善,采用量子力學手段獲得滿足具體體系要求,且具有更加精確的分子間勢能參數的方法已被廣泛接受,這可能是未來改進特定體系力場參數的重要方向。

圖5 protein/water/PFPE/scCO2體系自組裝過程(50 ns時體系的快照)[11]

圖6 水/AOK/scCO2體系自組裝過程(動力學模擬36 ns后體系快照)[9]

如圖7所示,楊曉寧等分析徑向密度分布函數圖發現,AOK/水/scCO2體系存在四層結構,分別是水核區、水核與表面活性劑頭基間的界面區、表面活性劑尾基與CO2之間的界面區以及CO2主體相區域,這符合目前對于scCO2RMs結構的一般認識。他們還指出,少量AOK尾端能一定程度地插入水核區,Berkowitz等[46]也有類似發現。

圖7 水/AOK/scCO2體系密度徑向分布函數曲線(以RMs質心為參考點)[9]

目前,MD方法作為研究超臨界反膠束體系性質的有效手段正受到越來越多學者和研究機構的重視。盡管獲得有效信息需要耗費巨大的計算機資源,而且理論工作、合適的力場構建以及體系設計等方面的研究工作亟待加強,但隨著人們對超臨界微乳體系認識的深入、實驗研究和儀器測試技術的進步以及新體系的不斷出現,學術界對MD作用的期待和它所能夠解決的問題都將有更大的提升。

3 結語

在科學技術領域,計算機模擬技術已成為和理論、實驗共存的主要研究手段。分子模擬(molecular simulation,MS)技術在化學化工、物理等學科發揮著越來越重要的作用。作為MS中的主要方法,MD是目前研究復雜化學體系結構、熱力學相行為以及動力學性質的有效手段,自然也是scCO2RMs體系模擬中最常用的方法。盡管已有相關商業軟件問世,但受算法、并行運算速度等因素限制,用MD可以處理的體系最多數萬個分子,即使以周期鏡像法進行彌補,但對于scCO2RMs這樣的復雜體系,能夠計算的分子數與實際體系仍有很大差距。

目前用于scCO2RMs體系平衡模擬的時間最多達到10-7s級,遠小于實際體系熱力學穩定時所需的時間尺度 (實驗觀察的平衡時間)。 當包含目的成分在內時,體系將更加復雜 (四元以上體系),自然對MD技術提出更高要求。此外,所使用的分子模型在模擬過程中各局部電荷按常數處理,事實上這種電荷信息因粒子間相互作用始終是變化的,故提高模擬有效性還需對現有分子力場作進一步改進和完善。

即使如此,從上述文獻分析討論可知,MD方法用于研究scCO2RMs體系是值得關注的。從模擬計算所涉及的體系和取得的初步結果看,借助于MD深入認識和揭示scCO2RMs體系的物化性質是計算化學未來發展的一個重要方向。

[1]Joshi D K， Prausnitz J M.AIChE J， 1984,30:522.

[2]Ryoo W，Webber S E，Johnston K P.Ind Eng Chem Res,2003,42:6 348.

[3]銀建中,周丹,王愛琴.化學進展,2009,21:2 505.

[4]Rapaport D C.The art of molecular dynamics simulation[M].2nd ed.Cambridge:Cambridge University Press,1995.

[5]Marrink S J，Mark A E.JACS,2003,125:15 233.

[6]Jalili S， Akhavan M. Colloids and Surf, A:Physicochemical and Eng Aspects,2009,352:99.

[7]Tsuzuki S， Shinoda W， Saito H， et al.J Phys Chem，B,2009,113:10 641.

[8]Salles F，Jobic H，Devic T，et al.ACS Nano,2010,4:143.

[9]Wu B， Yang X N， et al.Colloids and Surf， A:Physicochemical and Eng Aspects,2010,367:148.

[10]Chandran A，Prakash K.JACS,2010,132:12 511.

[11]Chaitanya V S V，Senapati S.JACS,2008,130:1 866.

[12]Smith R D， Fulton J L， Jones H K.Separ Sci Technol,1988,23:2 015.

[13]Gale R W，Fulton J L，Smith R D.JACS,1987,109:920.

[14]Johnston K P，Harrison K L.Science,1996,271:624.

[15]Sarbu T， Styranec T， Beckman E J.Nature,2000,405:165.

[16]Eastoe J， Gold S， Rogers S， et al.Angew Chem Int Ed,2006,45:3 675.

[17]Eastoe J,Gold S,Steytler D C.Langmuir,2006,22:9832.

[18]張陽,楊基礎,于養信,等.化學進展,2005,17:955.

[19]Brown D，Clarke J H R.J Phys Chem，1988,92:2 881.

[20]Linse P， Halle B.Mol Phys， 1989,67:537.

[21]Linse P.J Chem Phys，1989,90:4 992.

[22]Faeder J，Ladanyi B M.J Phys Chem,B,2000,104:1 033.[23]Faeder J，Ladanyi B M.J Phys Chem，2001,105:11 148.[24]Salaniwal S，Cui S T.Langmuir,1999,15:5 188.

[25]Lu L，Berkowitz M L.JACS,2004,126:10 254.

[26]Fan X，Potluri V K.JACS,2005,127:11 754.

[27]Saharay M， Balasubramanian S.J Chem Phys， 2004,120:9 694.

[28]Vrabec J， Stoll J， Hasse H.J Phys Chem， B,2001,105:12 126.

[29]Harris J G，Yung K H.J Phys Chem，1995,99:12 021.[30]Jorgensen W L， Chandrasekhar J.J Phys Chem， 1983,79:926.

[31]Jorgensen W L，Madura J D.Mol Phys，1985,56:1 381.

[32]Mahoney M W， Jorgensen W L.J Chem Phys， 2000,112:8 910.

[33]Berendsen H C，Postma J P M.Intermolecular forces.Reidel:Dordrecht,1981:331-333.

[34]Liu Y，Ichiye T.J Phys Chem，1996,100:2 723.

[35]Ryckaert J P， Ciccotti G， Berendsen H J C.J Comput Phys，1977,23:327.

[36]Andersen H C.J Comput Phys，1983,52:24.

[37]Salaniwal S，Cui S T.Ind Eng Chem，2000,39:4 543.

[38]Salaniwal S， Cui S T， Cochran H D， et al.Langmuir,2001,17:1 773.

[39]Salaniwal S， Cui S T， Cochran H D,et al.Langmuir,2001,17:1 784.

[40]Eastoe J,Bayazit Z.Langmuir,1996,12:1 423.

[41]Eicke H F,Christen H.Chim Acta，1978,61:2 258.

[42]Senapati S，Keiper J S.Langmuir,2002,18:7 371.

[43]Higashi H，Iwai Y.J Supercrit Fluids,1998,13:93.

[44]Nave S，Eastoe J.Langmuir,2000,16:8 741.

[45]Eastoe J，Downer A.PCCP,2000,2:5 235.

[46]Senapati S，Berkowitz M L.J Phys Chem,B,2003,107:12 906.

[47]Lu L Y，Berkowitz M L.JACS,2004,126:10 254.

[48]Tobias D J， Klein M L.J Phys Chem， 1996,100:6 637.

[49]Chowdhary J，Ladanyi B M.J Phys Chem,B,2009,113:15 029.

Advances of Molecular Dynamics Simulation for Supercritical Carbon Dioxide Microemulsion/Reverse Micelles System

Yin Jianzhong Wang Song Xu Junchen Yu Wen Xu Qinqin

SupercriticalCO2microemulsion/reverse micellessystem,which istransparentliquid of thermodynamic stability and formed by water-in-CO2(W/C)or CO2-in-water (C/W)with the presence of surfactant.Improving the properties of supercritical CO2could expands the practical application fields of supercritical fluid greatly.The thermodynamic properties,transport properties and self-assembly dynamics mechanism of the system are the foundation of applications.However,the information obtained from experimental results is insatiable for practical needs that's the main reason for use of molecular dynamics simulation.The progress of supercritical CO2microemulsion/reverse micelles system investigation by molecular dynamics simulation is summarized here.Force field construction,cell creation,balance control and other issues should be noted along with simulation process are presented,and the core direction of the research is discussed preliminarily.

Supercritical CO2;Reverse micelles;Microemulsion;Molecular dynamics simulation;Physical and chemical properties

O 64

*國家自然科學基金 （20976026，20976028）資助項目。

**銀建中，男，1964年生，博士，教授。大連市，116024。

2011-09-12）