N-丁基吡啶四氟硼酸鹽/水二元體系的分子動力學模擬研究

王 玉,朱光來,王晨晨,徐建強,馬趙鵬

(1.安徽師范大學 物理系 原子與分子物理研究所,蕪湖 241002; 2.安徽達爾智能控制系統股份有限公司 研發中心,蕪湖 241003)

1 引 言

室溫離子液體是指在室溫或接近室溫下呈液態、由有機陽離子和無機陰離子所組成的離子化合物，簡稱為離子液體[1].離子液體具有良好的熱穩定性和化學穩定性、極低的蒸汽壓、可循環利用等優點[1-3].由不對稱的陰、陽離子組成的結構特點，決定了離子液體具有這些獨特的性質，研究離子液體的微觀結構有助于了解離子液體的物理化學性質、溶劑特性等[3].

實驗與理論研究告訴我們在離子液體及其溶液內部呈現靜態或動態的非均質結構，存在著不同形態聚集體[4-7].離子液體與分子溶劑的混合體系中，組分比例不同會引起溶劑的熱力學性質和局部微觀結構的變化，對溶質分子的擴散、遷移等輸運性質產生重要影響.而水作為生活中的常見分子溶劑，研究它與離子液體混合物的性質和結構具有重要的科學意義和應用價值[8].分子模擬方面，目前報道較多的是對烷基咪唑型離子液體與乙腈、水等混合物的徑向分布函數、擴散系數、黏度的研究[9，10].對于吡啶離子液體，其熱力學性質方面的實驗研究已有較多報道[11-13]，而其與分子溶劑的混合物的性質與結構的分子模擬研究還有待進一步擴展[14，15].

本文以N-丁基吡啶四氟硼酸鹽[BPy]BF4與水作為研究對象，用分子動力學模擬方法研究不同比例下二元體系的密度和不同組分的擴散系數，通過分析各組分間的徑向分布函數、空間分布函數以及氫鍵作用，重點研究二元體系的微觀結構.

2 力場與模擬方法

2.1 力場參數

[BPy]BF4和水分子結構如圖1示意，其力場采用Lopes等人在OPLS_AA/AMBER力場基礎上發展起來的全原子力場模型[16]，函數形式如下：

(1)

圖1 [BPy]BF4和水分子結構示意圖

其中Vtotal表示系統的總能量，包含了鍵長、鍵角、二面角、van der Waals(VDW)作用以及靜電作用五個部分.其中非鍵作用力由Lennard-Jones (L-J)項和庫侖項組成.鍵拉伸和鍵彎曲以簡諧勢表示.水分子的力場參數采用TIP4P四點結構模型[17].

2.2 模擬方法

分子動力學模擬采用GROMACS軟件包[18]，結果分析也基于GROMACS提供的分析手段.模擬過程采用NPT系綜，溫度和壓力維持在298 K和0.1 MPa.模擬步長為2 fs，每0.2 ps取樣一次，范德華相互作用的截斷半徑以及庫侖相互作用的截斷半徑均取為1.5 nm.模擬過程采用標準的周期性邊界條件，用V-rescale方法維持體系溫度，用Parrinello-Rahman方法維持體系壓力，遠程靜電作用采用PME方法，鍵長約束采用LINCS算法.每個體系均先計算2 ns以使得體系構象達到平衡，然后在此基礎上再運行4 ns用來進行統計分析.

3 結果與討論

3.1 密度與自擴散系數

如表1 所示，二元體系的密度隨著離子液體比例的增加呈遞增的規律，模擬的數值與實驗值(純水0.997 g/cm3，純離子液體1.2141 g/cm3)相比較[11]，誤差均在允許的范圍內.

表1 模擬得到體系的密度以及各組分的擴散系數

體系的自擴散系數D由Einstein關系(式(2))擬合得到：

(2)

式中t為時間，ri為組分i質心的矢量位移，本工作中各組分的自擴散系數D是通過線性擬合系統平衡后1ns時間內的均方位移(MSD)得到.從模擬結果中不難看出，隨著離子液體比例的增大系統各組分的擴散系數均有明顯的遞減變化，這說明[BPy]BF4對該混合體系的擴散性質有很大的影響.另一方面來看，混合體系各組分的擴散系數隨著水分子比例的減少而減小，這一現象預示著水分子的存在可能會有效降低[BPy]BF4離子液體的黏度.

3.2 微觀結構

3.2.1徑向分布函數

溶液的局部組成和微觀分布情況可以通過分析它們的徑向分布函數(RDFs)來描述，該函數的物理意義是指以某個粒子為中心，在距離為r處另一粒子的分布幾率.圖2～圖6分別給出了二元體系中不同組分間的徑向分布函數隨離子液體摩爾分數xi的變化情況.

圖2 陽離子上丁基鏈末端碳原子CT4之間的徑向分布函數

圖2表示陽離子丁基鏈末端碳原子(CT4)間的徑向分布函數，其峰值隨著離子液體比例的增加并沒有明顯的變化；圖3表示陽離子吡啶環的HA1與陰離子的F原子間的徑向分布函數，其峰值呈現出隨離子液體比例的增加遞增的規律.上述結果表明離子液體的濃度對陽離子丁基鏈組成的非極性區域的聚集程度沒有多大影響，但對于吡啶環與陰離子組成的極性區域的影響卻很大.這主要由于極性區域的主要作用是陰、陽離子間較強的庫侖作用，而非極性區域則是烷基鏈之間的VDWs相互作用.

圖3 陽離子吡啶環上HA1原子與陰離子F原子間的徑向分布函數

圖4和圖5分別表示水分子中HW與HA1、F原子間的徑向分布函數，雖然二者均隨著離子液體比例的增加而遞增，但圖4中HW-F的峰值更高，且變化趨勢明顯，這表明混合體系中水分子與陰離子間的相互作用更強些，水分子更易靠近陰離子，原因可能是由于水分子與陰離子間的由氫鍵作用形成了網狀聚集體.

圖4 水分子HW原子與陰離子F原子間的徑向分布函數

圖5 水分子HW原子與陽離子吡啶環HA1原子間的徑向分布函數

圖6所示的是水分子中OW原子間的徑向分布函數，很高的峰值說明水分子因較強的氫鍵作用聚集.而且其峰值隨離子液體比例的增加也呈遞增的趨勢，可能是由于離子液體增加導致體系的黏度增大，擴散變慢，促使由氫鍵作用而形成的水分子小團簇更加穩定.

圖6 水分子OW原子間的徑向分布函數

3.2.2空間分布函數

空間分布函數(SDFs)能夠更加直觀、清晰地反映混合體系中各組分的空間分布情況.圖7和圖8分別表示離子液體的摩爾分數為0.1和0.5時的二元體系中不同組分間的空間分布函數.

圖7 離子液體摩爾分數為0.1的混合體系中不同組分間的空間分布函數

圖8 離子液體摩爾分數為0.5的混合體系中不同組分間的空間分布函數

由圖7(a)和圖8(a)可以看到，離子液體摩爾分數xi為0.1和0.5時，陰離子周圍的水分子分布幾率圖形非常相似，在四個F原子周圍的都近似呈均等分布，這是由于F原子與水分子之間因OW-HW…F氫鍵作用形成了聚集體，但xi為0.5時的isovalue值要高許多，說明這時水分子在陰離子周圍的幾率更大，這一結果與HW-F原子間的徑向分布函數隨離子液體比例增大而增加的規律相一致.由圖7(b)和圖8(b)可以看到，在xi為0.1時，陽離子周圍的水是環繞整個陽離子頭部分布的，這是因為水分子可以與陽離子吡啶環上的氫原子形成不同強度的氫鍵，而xi為0.5時，水主要分布在一側，此時水在陰離子周圍的分布更有序；如圖7(c)和圖8(c)所示，陰離子在陽離子周圍主要分布在頭部，這主要是由于陰、陽離子之間存在較強的庫侖作用所致.

3.2.3氫鍵作用

通過觀察圖3和圖4所示的F-HA1、F-HW之間的徑向分布函數可以發現，兩者第一個波峰的位置在0.170-0.253 nm范圍內，均短于H、F原子間的Van der Waals半徑之和(0.267 nm)[19]，這表明陽離子與陰離子、水分子與陰離子之間分別有類似C-HA1…F和OW-HW…F氫鍵形成.進一步比較它們的峰位可以看出，OW-HW…F氫鍵(最大峰對應位置約為0.172 nm)要強于C-HA1…F氫鍵(約為0.253 nm)，并且這兩種氫鍵作用都隨著水分子含量的減少而增強.

為了探究水分子自身之間形成的氫鍵性質隨離子液體變化的趨勢，我們計算了不同離子液體比例下二元體系中每個水分子周圍形成的平均氫鍵數目以及氫鍵的存活壽命，如圖9所示.由圖可知，隨離子液體的增加，水分子數目明顯減少，體系中每個水分子周圍形成的平均氫鍵數目減少，而氫鍵壽命則增大，這與文獻報道的水分子氫鍵的生存周期隨溶質濃度增大而增大的規律類似[20].此時離子液體各組分的擴散系數呈遞減的趨勢，說明離子液體比例越高各組分擴散的速度越慢，使得水分子間的運動也變得緩慢，增加了氫鍵斷裂的難度，其氫鍵壽命得以延長.

圖9 水分子的平均氫鍵數目與氫鍵壽命隨離子液體比例的變化

4 結 論

利用分子動力學模擬方法研究分析了吡啶離子液體[BPy]BF4與水組成的二元體系的微觀結構和各組分之間的相互作用，主要結論如下：水的存在會有效地降低離子液體的黏度，各組分的自擴散系數隨著水分子比例的減少而變??；隨著離子液體比例的增大，除了丁基鏈末端碳原間的徑向分布函數無明顯變化外，其它組分原子間的徑向分布函數的峰值均隨之增加 說明水與陰離子、水與陽離子吡啶環之間的相互作用逐漸增強，并且水與陰離子間的相互作用更強一些，水分子在微觀分布上更容易趨近陰離子，這與三者之間形成的氫鍵有密切關系；空間分布函數則直觀地反映出陰離子周圍的水的分布幾率近似呈均等分布，陰離子主要分布在陽離子的頭部周圍，而且水分子在陰離子周圍的分布幾率隨離子液體比例增大而增加；通過分析水分子間的氫鍵還發現，每個水分子周圍的平均氫鍵數目隨離子液體比例的增加而減少，而氫鍵壽命則隨之增大.