趙慧霞
中北大學信息商務學院, 山西 晉中 030600
分子動力學模擬(即MD模擬),最早于1957年就得到了實現(xiàn),當時人們就實現(xiàn)了在微觀尺度下對納米材料的力學性能和物理性質做研究,并且成為當時研究的熱點,主要是在關于分子動力學模擬的方面,著重模擬關于受限于納米材料中的小尺寸材料,還有有關流體的動力學行為[1,2].隨后,Hummer等人在2001年研究了H2O分子受限于管徑較小的納米碳管中的微觀結構及分子動力學模擬的動力學特征,實驗模擬結果表明,H2O分子在碳納米管的中央形成了一條由H—H鍵形成的單分子水鏈.接著,在同一年,科學家Koga等人將H2O分子放到碳管中做研究,并且控制管徑大小的不同,同時將該裝置置于500 bar~5 000 bar的壓強下,實驗結果表明,H2O分子形成了環(huán)狀的結構.
在Koga的模擬實驗中,由于H2O分子受限于碳納米管中,受到碳納米管的熱力學屏蔽作用,會出現(xiàn)類似堆疊的環(huán)狀六聚體或者有序的螺旋狀結構等不同的相態(tài).Mashl R Jay等人在2003年運用MD模擬的方法進一步探究了H2O分子在納米碳管中的受限作用下的結構以及動力學行為[3],模擬結果表明,受限水分子與體相水不同,形成了六聚體冰相結構[4].
本實驗中用到的實驗原理為分子動力學的基本原理:模擬利用分子、離子、原子間的相互作用勢函數(shù),通過分析可以得到每個分子所受到的作用力,然后可以進一步得到體系中各個分子在微觀狀態(tài)下的運動規(guī)律,然后對系統(tǒng)中粒子微觀狀態(tài)下的動量和其位置對時間求平均,即可求得該系統(tǒng)粒子的擴散系數(shù)、粘度、空間分布、能量以及壓力等物理量的宏觀性質的微觀結構分析.我們還可以用此方法計算該模擬體系中各種動力學性質和該模擬體系的平衡性質.在體系統(tǒng)計分析的過程中,如果我們選擇的統(tǒng)計方法更為恰當,同時還可以得到整個系統(tǒng)的宏觀性質.通過對牛頓運動學方程積分可以得到分子動力學系統(tǒng)中一系列的粒子的力學性質,在今后的分子動力學模擬中可以通過該積分方程的曲線,很準確的描述模擬體系中分子的運動軌跡,同時可以推出系統(tǒng)內粒子的速度、位置隨時間的變化情況.根據(jù)求解牛頓第二定律的微分方程,即可得到模擬體系中粒子的運動軌跡[5,6].用公式(1)可以表示該粒子的運動軌跡
(1)
式(1)的物理意義為:質量為mi的分子、受到力Fi時的位置矢量ri隨時間變化的運動方程[7].式中,根據(jù)勢函數(shù)U求梯度后,即可得出第i個原子上的作用力可用公式(2)表示
Fi=-▽iU
(2)
該模擬體系的動能隨著溫度T的變化,其變化規(guī)律如公式(3)所示
(3)
在本實驗中,體系模擬所采用的是由Refson K等人研究發(fā)明的Moldy分子模擬程序包[8],在這之前,已經有大量分子動力學模擬檢驗過了它的可靠性[3].本文的模擬是在正則系綜的體系下進行實驗的,然后運用“Nose-Hoover”的算法控制模擬體系的溫度分別在75 K、150 K和300 K的溫度下進行分子動力學模擬[9,10].設定1×1×60個元胞作為模擬體系盒子的大小,碳納米管取長度為14.7 nm的碳管,在模擬的過程中將32個H2O分子置于管內進行實驗,勢能截斷半徑為半個模擬盒子的邊長,模擬中采用Ewald求和處理長程力,使用周期性邊界條件.
本文的研究對象H2O分子的轉動慣量較小,以至于非常容易發(fā)生高速旋轉,所以選取了比較小的時間步長,取0.5 fs 作為模擬時間步長,以保證后期的運動方程的解的穩(wěn)定性.本實驗總的模擬總步數(shù)為1 400 000步,其中前400 000步用于使體系達到平衡狀態(tài),后1 000 000步對已處于平衡狀態(tài)的體系進行抽樣模擬,然后每隔100步保存記錄一次研究體系的狀態(tài)數(shù)據(jù),接著利用程序對實驗模擬數(shù)據(jù)進行系綜平均分析,進而得到水的徑向分布函數(shù)(RDF)曲線圖.
除此之外,本文還模擬了在300 K的溫度下,32個SPC/E模型的H2O分子分別處于受限SWCNT(6,6)、(7,7)、(8,8)碳納米管中的運動狀態(tài).
表1 不同類型碳納米管的幾何尺寸及其原子數(shù)Tab.1 The geometry size and atomic number of different types of SWCNT
在本文中用到好多參數(shù),其中如上表1所示,為本模擬實驗中碳納米管模型的參數(shù).
徑向分布函數(shù)(radial distribution function),簡稱為RDF,一般用g(r)表示徑向分布函數(shù)的方程.徑向分布函數(shù)用來反映模擬體系分子微觀結構特征的物理量,還可以同時分析說明分子結構性質的相關參數(shù)以及直接反映物質分子的聚集程度.g(r)的物理意義:距離中心粒子r處微體積元的分子局部數(shù)密度相對于體系平均密度的比值.
在本實驗中,模擬了水分子分別在溫度為75 K、150 K、300 K的狀態(tài)下的運動狀態(tài),計算并研究了32個SPC/E模型水分子在受限碳納米管中的分子O—O鍵的徑向分布函數(shù).經分析徑向分布函數(shù)曲線,發(fā)現(xiàn)r在小于2.4 ?的范圍內,H—H鍵的徑向分布曲線值為零,這表明SWCNT(9,9)中H2O分子之間的距離至少為2.4 ?;在2.7 ?處徑向分布函數(shù)曲線圖出現(xiàn)了最高峰,這說明在距離被研究的中心H2O分子2.7 ?處出現(xiàn)其他H2O分子的可能性最大,說該處的區(qū)域H2O分子密度最大;而在r大于8.0 ?時,徑向分布函數(shù)曲線的值非常接近1,說明在此時的H2O流體與均勻流體的性質比較相似.根據(jù)以上三幅圖中的徑向分布曲線函數(shù)圖分析可得,徑向分布函數(shù)曲線的峰值位置在不同溫度下位置相同,但峰值會隨著溫度的變化有明顯的差異:在75 K時徑向分布函數(shù)曲線圖的峰值最高,在150 K時的峰值低于75 K時的峰值,在300 K時徑向分布函數(shù)曲線圖的峰值最低.由此可得到以下結論:在研究模型的碳納米管的管徑尺寸相同時,模擬體系的溫度越高,對應的徑向分布函數(shù)曲線的峰值越低.該現(xiàn)象說明,溫度的升高促進了分子的熱運動,從而降低了水的有序化程度.
圖1 (9,9)碳納米管中水分子在75 K溫度下的徑向分布曲線Fig.1 (9,9) radial distribution curve of water molecules in carbon nanotubes at 75 K temperature圖2 (9,9)碳納米管中水分子在150 K溫度下的徑向分布曲線Fig.2 (9,9) radial distribution curve of water molecules in carbon nanotubes at 150 K temperature
圖3 (9,9)碳納米管中水分子在300 K溫度下的徑向分布曲線Fig.3 (9,9) radial distribution curve of water molecules in carbon nanotubes at 300 K temperature圖4 H2O分子在300 K不同管徑碳納米管中的O—O徑向分布曲線Fig.4 O—O RDF of H2O molecules in 300 K carbon nanotubes with different diameters
本文還模擬了在300 K的溫度下,32個SPC/E模型的H2O分子分別處于受限SWCNT(6,6)、(7,7)、(8,8)碳納米管中的運動狀態(tài).如圖4所示,水分子在300 K不同管徑碳納米管中的O—O徑向分布曲線.實驗數(shù)據(jù)表明,O—O鍵的徑向分布函數(shù)曲線發(fā)現(xiàn)在(6,6)管中出現(xiàn)的峰值個數(shù)比較多,說明在其中水接近冰相結構,在(8,8)管中峰值最高,說明在(8,8)的碳納米管中水分子的有序化程度最高.