航空煤油冰點及黏溫關(guān)系的分子動力學(xué)模擬

陳雪嬌，侯 磊，李師瑤

中國石油大學(xué)(北京)油氣管道輸送安全國家工程實驗室/石油工程教育部重點實驗室/城市油氣輸配技術(shù)北京市重點實驗室，北京 102249

航空煤油冰點及黏溫關(guān)系的分子動力學(xué)模擬

陳雪嬌，侯 磊*，李師瑤

中國石油大學(xué)(北京)油氣管道輸送安全國家工程實驗室/石油工程教育部重點實驗室/城市油氣輸配技術(shù)北京市重點實驗室，北京 102249

航空煤油的冰點和黏溫特性是衡量煤油低溫流動性的重要質(zhì)量指標(biāo)。本文以某航空煤油的碳數(shù)分布、碳?xì)滟|(zhì)量比、芳香烴飽和烴質(zhì)量比等參數(shù)為依據(jù)，建立了5種烴混合的煤油體系模型，進行了分子動力學(xué)模擬。根據(jù)密度和擴散系數(shù)隨溫度的變化預(yù)測了煤油體系的冰點，模擬冰點值與實驗結(jié)果一致；通過計算徑向分布函數(shù)分析了冰點溫度附近煤油體系內(nèi)部的微觀結(jié)構(gòu)，發(fā)現(xiàn)同種分子間的聚集是影響煤油冰點溫度的主要原因，芳香烴對煤油體系冰點的影響大于鏈狀飽和烴的影響；采用Rouse Model模擬了煤油體系的黏溫關(guān)系，模擬黏溫指數(shù)與通過實驗擬合的黏溫指數(shù)相同，模擬黏度隨溫度的變化趨勢與實驗結(jié)果一致。

航空煤油；分子動力學(xué)；冰點；間接模擬；黏溫關(guān)系

0 引言

冰點是衡量航空燃料的重要質(zhì)量指標(biāo)，航空煤油經(jīng)過冷卻形成固態(tài)烴類結(jié)晶，然后在特定的條件下開始升溫，當(dāng)烴類結(jié)晶消失時的溫度即為航空煤油的冰點。航空煤油作為航空渦輪發(fā)動機的重要燃料，需要在低達(dá)-55 ～-60 ℃的氣溫下工作。若在高空氣溫下煤油黏度過大，噴嘴就無法把燃料充分霧化，導(dǎo)致燃燒不充分，影響發(fā)動機的正常運行；若高空氣溫接近煤油的冰點，煤油中就會出現(xiàn)結(jié)晶，堵塞油管，釀成重大事故。因此，航空煤油的冰點和黏溫特性對煤油的低溫流動性至關(guān)重要。

針對石油類混合物的研究趨勢是從宏觀定性向微觀定量轉(zhuǎn)變，而分子動力學(xué)模擬在微觀和定量方面有著強大的優(yōu)勢。分子動力學(xué)方法在石油領(lǐng)域的研究起步較晚，始于20世紀(jì)末期[1]。Mondello等[2-3]對線直鏈烴和支鏈烷烴的黏度進行了分子動力學(xué)模擬，分析了平衡分子動力學(xué)計算方法和非平衡分子動力學(xué)計算方法的模擬結(jié)果，驗證了Rouse Model在烴類物質(zhì)黏溫關(guān)系模擬方面的適用性。殷開梁等[4]采用分子動力學(xué)方法研究了以C16為主要組分的生物柴油的凝固過程，根據(jù)定容比熱容和擴散系數(shù)的突變點判斷生物柴油的凝點，得到的模擬值與實驗值基本一致。Tae等[5]模擬了2 700 K高溫時碳數(shù)為C12～C160正構(gòu)烷烴的黏度，分析了在極端高溫條件下正構(gòu)烷烴的黏溫關(guān)系。Makrodimitri等[6]用分子動力學(xué)方法模擬了C8～C96正構(gòu)烷烴混合物的黏度，并構(gòu)建了擴散系數(shù)與溫度間的關(guān)系式。分子動力學(xué)模擬技術(shù)以量子力學(xué)和經(jīng)典牛頓力學(xué)為基礎(chǔ)，能夠分析分子形態(tài)、分子間作用力等微觀性質(zhì)的變化。目前對單一體系烴類的相變點和黏溫關(guān)系的研究較為成熟，對多種烴類混合體系的相關(guān)研究較少。本文根據(jù)某航空煤油的碳數(shù)分布、碳?xì)滟|(zhì)量比、芳香烴飽和烴質(zhì)量比等參數(shù)建立了5種烴混合的煤油體系模型，在能量最小化后進行分子動力學(xué)模擬，根據(jù)密度和擴散系數(shù)隨溫度的變化預(yù)測煤油體系的冰點，通過計算徑向分布函數(shù)分析冰點溫度附近煤油體系內(nèi)部的微觀結(jié)構(gòu)，采用Rouse Model模擬煤油體系的黏溫關(guān)系。

1 煤油體系模型構(gòu)建與模擬方法

1.1 模型構(gòu)建

表1中列出了航空煤油的碳數(shù)分布及正構(gòu)、異構(gòu)烷烴組成，其中C11～ C13占97.499％，此外，芳香烴含量為3.9％，飽和烴含量為96.1％。根據(jù)C11～ C13烴分子中正構(gòu)、異構(gòu)烷烴比例及芳香烴、飽和烴比例，選出5種代表性分子，分別是2-甲基癸烷(C11H24)、正十二烷(C12H26)、2-甲基十一烷(C12H26)、2-甲基十二烷(C13H28)、1,4,5-三甲基萘(C13H14)[7-8]。每種分子在煤油分子動力學(xué)模型中的個數(shù)及比例見表2，可以看出，煤油分子動力學(xué)模型中的碳數(shù)分布、正構(gòu)/異構(gòu)烷烴比例、芳香烴/飽和烴比例與實際航空煤油的測試數(shù)據(jù)相近，煤油分子動力學(xué)模型中的元素組成可以反映實際航空煤油的組成。利用以上5種分子模型，使用PACKMOL和MOLTEMPLATE構(gòu)建煤油體系的全原子分子動力學(xué)模型，體系中共包含7 611個原子，模型的初始尺寸設(shè)置為42×42×42 ?，邊界條件設(shè)置為周期性邊界條件，力場選為OPLS-AA力場，該力場適用于模擬芳香烴、飽和烴等有機分子系統(tǒng)[9-10]。采用LAMMPS(Large-scale Atomic/Molecular Massively Parallel Simulator)大規(guī)模原子/分子并行模擬器進行整個模擬過程。

表1 煤油碳數(shù)分布測定結(jié)果Table 1 Test results of carbon number distribution in kerosene

表2 煤油分子動力學(xué)模型中分子的個數(shù)及比例Table 2 Molecule number and proportion of kerosene molecular dynamic model

1.2 模擬方法

用共軛梯度法對體系進行結(jié)構(gòu)優(yōu)化，分別在NPT和NVT系綜下進行173～273 K溫度范圍內(nèi)的2 ns的模擬退火，獲得合理的模型尺寸和分子間距，使模型的整體構(gòu)型穩(wěn)定合理。對每個模擬溫度點進行1 000 ps的NPT系綜下的系統(tǒng)平衡弛豫，再進行NPT和NVT系綜的模擬。在NPT系綜下，模擬時間為1 000 ps，時間步長為1 fs，每1 000步輸出一次分子的運動軌跡，對溫度控制采用Nosé Hoover方法，對壓力控制采用Berendsen法，對位能截斷采用球形截斷法，截斷半徑為18.5 ?，對庫侖作用和范德華作用均采用PPPM加和方法進行長程校正[11-12]。在NVT系綜下，模擬時間為2 000 ps，其余設(shè)置與NPT系綜下的設(shè)置相同。

在分子動力學(xué)模擬中，系統(tǒng)的平衡用溫度和體系總能量兩個參數(shù)來判定。當(dāng)體系溫度在模擬溫度的±5％范圍內(nèi)浮動且能量收斂參數(shù)ΔEconverge＜0.001[13]時，認(rèn)為系統(tǒng)達(dá)到平衡。ΔEconverge按下式計算：

式中，E(0)為體系在初始狀態(tài)時的總能量，kcal/mol；E(i)為第i步模擬時的體系總能量，kcal/mol；N為模擬總步數(shù)。

以273 K時的模擬為例，圖1為273 K時0～300 ps的溫度和能量隨時間變化的曲線，溫度在273 K的-5％～+5％(-259.35～286.65 K)范圍內(nèi)波動，能量滿足收斂條件，此時體系達(dá)到平衡。從300 ps開始每隔100 ps取一次密度的平均值，并計算每段平均值間的標(biāo)準(zhǔn)偏差，當(dāng)標(biāo)準(zhǔn)偏差小于0.001時[14]，認(rèn)為煤油體系模型達(dá)到平衡狀態(tài)，密度模擬結(jié)果可靠，取300～1 000 ps的密度平均值作為273 K時的密度值，具體數(shù)據(jù)見表3。273 K時的密度模擬值為0.776 995 g/cm3，與文獻(xiàn)[15]中的密度值0.784 6 g/cm3相差0.97％，說明煤油分子動力學(xué)模型的構(gòu)建和力場的選取都是合理的，能夠反映實際航空煤油的性質(zhì)。其余溫度點下的平衡判定方法和統(tǒng)計物理量的計算方法與273 K時的計算方法相同。

圖1 273 K時煤油體系模型的溫度和能量曲線Fig. 1 Temperature and energy curve of kerosene system at 273 K

2 模擬結(jié)果與分析

煤油體系模型經(jīng)過退火處理和1 000 ps的系統(tǒng)平衡弛豫后，進行1 000 ps的NPT系綜模擬用于計算體系的密度，再進行2 000 ps的NVT系綜模擬用于求取密度、均方位移、擴散系數(shù)、徑向分布函數(shù)、回轉(zhuǎn)半徑和松弛時間。

表3 273 K時的分段密度值和標(biāo)準(zhǔn)偏差Table 3 Segmented density value and standard deviation at 273 K

2.1 冰點預(yù)測

2.1.1 密度模擬與冰點預(yù)測

在NPT系綜下，體系的壓力和溫度固定，當(dāng)體系的溫度在模擬溫度的±5％范圍內(nèi)浮動且能量收斂參數(shù)ΔEconverge＜0.001時，通過模型中所有原子的質(zhì)量和模型體積的統(tǒng)計平均值求取密度。在NPT系綜下對173 ～223 K溫度區(qū)間內(nèi)的密度進行模擬，當(dāng)溫度達(dá)到冰點溫度時，煤油的密度會發(fā)生突變，根據(jù)密度-溫度關(guān)系曲線的轉(zhuǎn)折點判斷煤油體系的冰點[16-17]。

圖2表示173 ～223 K溫度區(qū)間內(nèi)煤油體系模型的密度模擬結(jié)果，每5 K計算一次，每個溫度下的密度取平衡后的運動軌跡的分段密度平均值，用這種分段取平均的方法可以消除不同的初始構(gòu)型對模擬結(jié)果的影響。隨著溫度的降低，在198 K附近，密度-溫度曲線出現(xiàn)明顯轉(zhuǎn)折點，表明在該溫度下煤油體系模型由液態(tài)逐漸轉(zhuǎn)為固態(tài)，預(yù)測煤油體系的冰點在198 K(-75 ℃)附近，這與航空煤油冰點低于-60 ℃的實驗結(jié)果是一致的。

圖2 173～223 K時煤油體系模型的密度-溫度曲線Fig. 2 Density curve of kerosene system versus temperature from 173 to 223 K

2.1.2 擴散系數(shù)模擬與冰點預(yù)測

在分子動力學(xué)模擬中，原子由起始位置不停地做無規(guī)則運動，原子位移平方的平均值為均方位移(Mean Square Displacement, MSD)。在煤油體系模型中，2-甲基十二烷和2-甲基十一烷等都是由多個原子組成的長鏈分子，若直接用每個原子的中心計算MSD會導(dǎo)致統(tǒng)計誤差增大，因此要用分子質(zhì)心來計算每種分子的均方位移[18]。均方位移的計算式如下：

Einstein擴散定律描述了均方位移和擴散系數(shù)間的關(guān)系，即

式中，D為擴散系數(shù)，cm2/s，反映的是體系內(nèi)分子遷移性和活動性的強弱，擴散系數(shù)發(fā)生突變的溫度點認(rèn)為是體系的相變點[19]，可以用來判斷煤油體系的冰點。當(dāng)計算時間t→∞時，均方位移曲線的斜率即為6倍的擴散系數(shù)[20]。

當(dāng)煤油體系接近冰點溫度時，體系內(nèi)分子的運動會受到阻礙，尤其是長鏈分子的遷移性和活動性會急速下降，導(dǎo)致擴散系數(shù)迅速變小[21]，但擴散系數(shù)-溫度曲線的轉(zhuǎn)折點要略高于煤油體系的冰點，因為擴散系數(shù)發(fā)生突變時，局部分子遷移性的降低還未影響到整個煤油體系，從宏觀來看，此時煤油中會出現(xiàn)局部的結(jié)晶，但是整體還未凝固[22]。

圖3表示NVT系綜下煤油體系模型中的2-甲基十二烷、2-甲基十一烷在173 ～223 K溫度區(qū)間內(nèi)的均方位移-時間曲線，二者的均方位移隨著溫度降低而不斷減小，說明溫度越低，分子的遷移和運動受到的阻礙越大。擬合每個溫度點下的均方位移-時間關(guān)系曲線的斜率，得到擴散系數(shù)-溫度關(guān)系曲線(圖4)，隨著溫度的降低，2-甲基十二烷和2-甲基十一烷的擴散系數(shù)D均在208 K處出現(xiàn)突變，說明此時分子排列形態(tài)發(fā)生了變化，煤油體系模型中分子的聚集阻礙了2-甲基十二烷和2-甲基十一烷長鏈烴的運動，模型中的分子運動趨于緩慢，體系中的分子趨于聚集。綜合密度模擬和擴散系數(shù)的模擬結(jié)果，可以判斷煤油中分子出現(xiàn)聚集的溫度在208 K附近，煤油體系的冰點在198 K附近。

圖3 均方位移-時間關(guān)系曲線Fig. 3 Mean squared displacements curves versus time

圖4 擴散系數(shù)-溫度關(guān)系曲線Fig. 4 Diffusion coeffcient versus temperature

2.2 徑向分布函數(shù)模擬與內(nèi)部結(jié)構(gòu)分析

徑向分布函數(shù)g(r)(radial distribution function, RDF)表示在距離參考粒子r處搜索到另一個粒子的概率，是系統(tǒng)在參考粒子處的區(qū)域密度與整個體系的平均密度的比，體現(xiàn)液體中分子聚集的性質(zhì)，反映液體內(nèi)部的結(jié)構(gòu)[23]。當(dāng)溫度降低時，煤油體系模型中的分子出現(xiàn)聚集，聚集程度的大小體現(xiàn)在徑向分布函數(shù)g(r)的峰值上，峰值越高，聚集在r處的分子越多。對于航空煤油，如果在小范圍內(nèi)聚集的分子突然大幅增加，會導(dǎo)致在煤油內(nèi)部出現(xiàn)結(jié)晶，使體系的黏度增大，影響航空煤油的正常使用。

圖5 不同種類分子間的徑向分布函數(shù)Fig. 5 Radial distribution function between different molecular types

圖6 同種分子間的徑向分布函數(shù)Fig. 6 Radial distribution function between same molecules

圖5和圖6表示NVT系綜下198 ～208 K溫度區(qū)間內(nèi)長鏈烴代表物質(zhì)2-甲基十二烷、2-甲基十一烷和芳香烴代表物質(zhì)1,4,5-三甲基萘的徑向分布函數(shù)曲線，橫坐標(biāo)代表其他粒子與中心粒子間的距離，即分子發(fā)生聚集的位置，縱坐標(biāo)代表分子聚集的程度。圖5表示不同分子間的徑向分布函數(shù)。以2-甲基十二烷為中心，2-甲基十二烷和2-甲基十一烷在不同溫度下的徑向分布函數(shù)數(shù)值均在1左右，表明這兩種不同種類的分子未出現(xiàn)明顯聚集。2-甲基十二烷和1,4,5-三甲基萘在不同溫度下的徑向分布函數(shù)數(shù)值在5 ～10 ?附近略大于1，表明長鏈烴2-甲基十二烷和芳香烴1,4,5-三甲基萘有輕微的小范圍的聚集，但并未對煤油體系造成影響。圖6表示同種分子間的徑向分布函數(shù)，2-甲基十二烷，2-甲基十一烷，1,4,5-三甲基萘在198～208 K溫度范圍內(nèi)都呈現(xiàn)了不同程度的聚集，隨著溫度的降低，徑向分布函數(shù)值增大，分子聚集程度增強。

表4 2-甲基十二烷，2-甲基十一烷，1,4,5-三甲基萘的聚集位置和徑向分布函數(shù)數(shù)值Table 4 Accumulation position and RDF value of 2-methyldodecane, 2-Methylundecane and 1,4,5-trimethylnaphthalene

3種分子聚集的位置和徑向分布函數(shù)值見表4，2-甲基十二烷和2-甲基十一烷均在1.11 ?處出現(xiàn)聚集，但1,4,5-三甲基萘的聚集發(fā)生在1.09 ?處，略小于兩個長鏈烴的聚集位置，這主要是由于芳香烴1,4,5-三甲基萘的尺寸更小，導(dǎo)致聚集位置要小于長鏈烴。徑向分布函數(shù)數(shù)值的大小關(guān)系為：1,4,5-三甲基萘＞2-甲基十一烷＞2-甲基十二烷，其中芳香烴1,4,5-三甲基萘的RDF約為長鏈烴的10倍，表明芳香烴的聚集更加集中，對分子聚集起到主要作用。因此，從RDF分析來看，同種分子間的聚集是影響煤油冰點溫度的主要原因，芳香烴對煤油冰點的影響要大于飽和鏈烴的影響。

2.3 黏溫關(guān)系模擬

在分子動力學(xué)模擬中，黏度模擬分為平衡動力學(xué)法(Equilibrium Molecular Dynamics, EMD)和非平衡動力學(xué)法(Non-Equilibrium Molecular Dynamics, NEMD)。平衡動力學(xué)方法又分為Green-Kubo、Einstein等直接模擬法，以及Rouse Model、Debye-Stokes-Einstein等間接模擬法[24]。對于長鏈分子或者低溫的系統(tǒng)，用Green-Kubo法和Einstein法等直接模擬法計算黏度十分耗費計算資源。Debye-Stokes-Einstein方程中有與水力邊界條件和分子形狀有關(guān)的經(jīng)驗系數(shù)，需要根據(jù)實驗推算。Rouse Model最初是用來估測符合Gaussian鏈性質(zhì)的聚合物的黏彈性，Gaussian鏈需滿足：

式中，R為分子鏈末端距，?；Rg為回轉(zhuǎn)半徑，?。

Rouse Model的計算式如下：

式中，η(Rg)為體系的黏度，Pa·s；ρ為體系的密度，kg/m3；R為摩爾氣體常數(shù)，8.314 J/(mol·K)；T為體系溫度，K；Rg為回轉(zhuǎn)半徑，m；M為摩爾質(zhì)量，kg/mol；D為擴散系數(shù)，m2/s。

Mondello[25]和Zhang[26]等用Rouse Model模擬了直鏈烴、帶有分支的鏈烴、多環(huán)芳香烴等碳?xì)浠衔锏酿ざ龋M發(fā)現(xiàn)碳數(shù)大于66的鏈烴才能滿足Gaussian的條件，導(dǎo)致使用Rouse Model模擬C6～C66間的相對較短的鏈狀烴的黏度時，模擬值比實驗值低，但碳數(shù)越大，偏差越小，當(dāng)碳原子數(shù)超過66時，模擬值與實驗值接近。Mondello等[25,27]用Rouse Model分別模擬了298 K時n-C10、n-C16的黏度，模擬值與實驗值的對比見表5，η(Rg)是用式5模擬計算出的黏度，η是實驗測得的黏度。結(jié)果表明，n-C16的黏度模擬值更接近實驗值。模擬值與實驗值偏差較大，n-C16模擬值的黏溫關(guān)系為η=32.67e-0.01T，實驗的黏溫關(guān)系為η=366.5e-0.01T，黏溫指數(shù)均為0.01，模擬黏度隨溫度的變化趨勢與實驗結(jié)果是一致的。因此，采用Rouse Model模擬計算煤油體系的黏度。

計算煤油模型的黏度時，Rouse Model中Rg、M、D三個參數(shù)需要取松弛時間最長的分子的統(tǒng)計平均值[28-29]。鏈烴的弛豫時間與主鏈的長短有關(guān)。在煤油體系模型中，正十二烷和2-甲基十二烷是兩個主鏈最長的分子，其主鏈上都有12個碳原子，兩種分子在293 K時的松弛時間分別為95 ps和102 ps，因此，用2-甲基十二烷的Rg、M、D值來計算煤油模型的黏度。

248～303 K溫度范圍內(nèi)的分子動力學(xué)參數(shù)及η(Rg)模擬值見表6，ρ為煤油體系模型整體的密度，Rg、M、D分別為2-甲基十二烷的回轉(zhuǎn)半徑、摩爾質(zhì)量、擴散系數(shù)。隨著溫度的降低，煤油模型的密度增大，擴散系數(shù)降低，這與密度和擴散系數(shù)的模擬趨勢是一致的，但是2-甲基十二烷的回轉(zhuǎn)半徑并未隨溫度發(fā)生明顯變化，這主要是由于在248～303 K的溫度區(qū)間內(nèi)，煤油模型整體還處于液態(tài)，2-甲基十二烷本身的形態(tài)和尺寸還未受到明顯影響。

圖7表示航空煤油的實驗黏溫曲線和模擬黏溫曲線。實驗黏溫關(guān)系式為η=3.3639e-0.0257T，模擬黏溫關(guān)系式為η=0.92879e-0.0257T，模擬值比實驗值低，這主要是由于煤油模型中的分子未滿足Gaussian鏈的要求。雖然模擬值與測量值間存在偏差，但是用Rouse Model計算的黏溫曲線與實驗擬合的黏溫曲線趨勢相同，黏溫指數(shù)均為-0.0257，說明Rouse Model可以用來模擬預(yù)測烴類及烴類混合物的黏溫關(guān)系。

圖7 航空煤油實驗黏溫曲線和模擬黏溫曲線Fig. 7 Experimental and simulated viscosity of aviation kerosene versus temperature

表5 298 K時模擬值與實驗值的的對比Table 5 Contrast between simulated and experimental value at 298 K

表6 248 ～303 K溫度區(qū)間內(nèi)分子動力學(xué)參數(shù)和η(Rg)模擬值Table 6 Molecular dynamics parameters and η(Rg) simulation value from 248 to 303 K

3 結(jié)論

(1)通過分析煤油分子動力學(xué)模型的密度-溫度曲線和擴散系數(shù)-溫度曲線的突變點，預(yù)測煤油體系的冰點在198 K(-75 ℃)附近，煤油體系中長鏈分子2-甲基十二烷和2-甲基十一烷集中聚集的溫度高于冰點溫度，在208 K(-65 ℃)附近。

(2)通過模擬冰點溫度附近煤油模型中2-甲基十二烷、2-甲基十一烷、1,4,5-三甲基萘的徑向分布函數(shù)，發(fā)現(xiàn)同種分子間聚集程度高，這是影響煤油冰點溫度的主要原因，芳香烴對煤油體系冰點的影響大于飽和鏈烴的影響。

(3)結(jié)合2-甲基十二烷的回轉(zhuǎn)半徑、摩爾質(zhì)量、擴散系數(shù)及煤油模型的密度，采用Rouse Model預(yù)測煤油體系的黏溫關(guān)系，模擬黏溫指數(shù)與通過實驗擬合的黏溫指數(shù)相同，模擬黏度隨溫度的變化趨勢與實驗結(jié)果一致。

[1]于共奇. 重質(zhì)油組分及其溶解機理的分子動力學(xué)研究[D]. 青島: 中國石油大學(xué)(華東), 2013. [YU G Q. Molecular dynamics study of the composition and the dissolution theory of heavy oil[D]. Qingdao: China University of Petroleum (East China), 2013.]

[2]MONDELLO M, GREST G S. Molecular dynamics of linear and branched alkanes[J]. Journal of Chemical Physics, 1995, 103(103): 7 156-7 165.

[3]STEVENS M J, MONDELLO M, GREST G S, et al. Comparison of shear fow of hexadecane in a confned geometry and in bulk[J]. Journal of Chemical Physics, 1997, 106 (17): 7 303-7 314.

[4]殷開梁. 分子動力學(xué)模擬的若干基礎(chǔ)應(yīng)用和理論[D]. 杭州: 浙江大學(xué), 2006. [YIN K L. The application and theory of molecular dynamics[D]. Hangzhou: Zhejiang University, 2006.]

[5]KWON T W, SONG H L. Molecular dynamics simulation studies of mid-size liquid n-alkanes, C12-C160[J]. Bulletin of the Korean Chemical Society, 2015, 36(4): 1 165-1 171.

[6]MAKRODIMITRI Z Α, HELLER A, KOLLER T M, et al. Viscosity of heavy n-alkanes and diffusion of gases therein based on molecular dynamics simulations and empirical correlations[J]. Journal of Chemical Thermodynamics, 2015, 91: 101-107.

[7]于維銘. 航空煤油替代燃料火焰?zhèn)鞑ニ俣扰c反應(yīng)動力學(xué)機理研究[D]. 北京: 清華大學(xué), 2014. [YU W M. Study on fame speed and chemical reaction mechanism for alternative fuels of aviation kerosene[D]. Beijing: Tsinghua University, 2014.]

[8]于維銘, 袁振, 鐘北京. 航空煤油火焰?zhèn)鞑ニ俣葘嶒炁c動力學(xué)研究[J]. 工程熱物理學(xué)報, 2014(10): 2 102-2 107. [YU W M, YUAN Z, ZHONG B J. Experimental and kinetic study on fame speed of aviation kerosene[J]. Journal of Physical Thermodynamics, 2014(10): 2 102-2 107.]

[9]SUN H. COMPASS: An ab initio force-field optimized for condensed-phase applicationsoverview with details on alkane and benzene compounds[J]. Journal of Physical Chemistry B, 1998, 102(38): 7 338-7 364.

[10]MIKAMI Y, LIANG Y, MATSUOKA T, et al. Molecular dynamics simulations of asphaltenes at the oil-water interface: From nanoaggregation to thin-film formation[J]. Energy & Fuels, 2013, 27(4): 1 838-1 845.

[11]張軍, 房體明, 王業(yè)飛, 等. 烷烴油滴在超臨界二氧化碳中溶解的分子動力學(xué)模擬[J]. 中國石油大學(xué)學(xué)報(自然科學(xué)版), 2015, 39(2): 124-129. [ZHANG J, FANG T M, WANG Y F, et al. Molecular dynamics simulation of dissolution of n-alkanes droplets in supercritical carbon dioxide[J]. Journal of China University of Petroleum(Editon of Natural Science), 2015, 39(2): 124-129.]

[12]任強, 代振宇, 周涵. 重油膠體結(jié)構(gòu)的介觀模擬[J]. 石油學(xué)報, 2013, 29(1): 86-94. [REN Q, DAI Z Y, ZHOU H. Mesoscale simulation on colloidal structures of heavy oil[J]. Acta Petrolei Sinica, 2013, 29(1): 86-94.

[13]王大放. 用分子動力學(xué)方法模擬高壓電場對碳酸鈣結(jié)晶過程的影響[D]. 青島: 中國石油大學(xué)(華東), 2010. [WANG D F. Molecular dynamics simulation on the impact of high voltage electrostatic feld on crystallization of calcium carbonate aqueous solution[D]. Qingdao: China University of Petroleum (East China), 2010.]

[14]RAPAPORT D C. The art of molecular dynamics simulation[M]. Cambridge: Cambridge University Press, 2004.

[15]鄭東, 于維銘, 鐘北京. RP-3航空煤油替代燃料及其化學(xué)反應(yīng)動力學(xué)模型[J]. 物理化學(xué)學(xué)報, 2015(4): 636-642. [ZHENG D, YU W M, ZHONG B J. RP-3 aviation kerosene surrogate fuel and the chemical reaction kinetic model[J]. Acta Physico-Chimica Sinica, , 2015 (4): 636-642.]

[16]YANG Q, YANG X P, LI X. D, et al. The curing and thermal transition behavior of epoxy resin: A molecular simulation and experimental study[J]. Rsc Advances, 2013 (3): 7 452-7 459.

[17]YANG Q, CHEN X, HE Z, et al. The glass transition temperature measurements of polyethylene: Determined by using molecular dynamic method[J]. Rsc Advances, 2016, 6 (15): 12 053-12 060.

[18]KWON T W, SONG H L. Molecular dynamics simulation studies of mid-size liquid n -alkanes, C12-C160[J]. Bulletin of the Korean Chemical Society, 2015, 36(4): 1 165-1 171.

[19]CHOI J, YU S, YANG S, et al. The glass transition and thermoelastic behavior of epoxy-based nanocomposites: A molecular dynamics study[J]. Polymer, 2011, 52, 5 197-5 203.

[20]陳正隆. 分子模擬的理論與實踐[M]. 北京: 化學(xué)工業(yè)出版社, 2007. [CHEN Z L. The theory and practice of molecular dynamics[M]. Beijing: Chemical Industry Press, 2007.]

[21]殷開梁, 閆朋克, 孫靜, 等. 生物柴油分子系統(tǒng)凝固過程的分子動力學(xué)模擬[J]. 計算機與應(yīng)用化學(xué), 2011, 28(5): 602-606. [YINK L, YAN K P, SUN J, et al. Molecular dynamics simulation on solidifcation process of bio-diesel molecular systems. Computers and Applied Chemistry, 2011, 28(5): 602-606.]

[22]LI M, LIU X Y, et al. Molecular dynamics simulation on glass transition temperature of isomeric polyimide[J]. Express Polymer Letters, 2009, 3(10): 665-675.

[23]陳敏伯. 計算化學(xué)--從理論化學(xué)到分子模擬[M]. 北京: 科學(xué)出版社, 2009. [CHEN M B. Computational chemistry—from theoretical chemistry to molecular dynamics[M]. Beijing: Science Press, 2009.

[24]MONDELLO M, GREST G S. Molecular dynamics of linear and branched alkanes[J]. Journal of Chemical Physics, 1995, 103(103): 7 156-7 165.

[25]MONDELLO M, GREST G S, GARCIA A R, et al. Molecular dynamics of linear and branched alkanes: Simulations and nuclear magnetic resonance results[J]. Journal of Chemical Physics, 1996, 105(12): 5208-5215.

[26]ZHANG L Q, GREEN M L. Analyzing properties of model asphalts using molecular simulation[J]. Energy & Fuels, 2007, 21(3): 1 712-1 716.

[27]MONDELLO M, GREST G S. Viscosity calculations of n-alkanes by equilibrium molecular dynamics[J]. Journal of Chemical Physics, 1997, 106(106): 9 327-9 336.

[28]ZHANG L Q. Physical and mechanical properties of model asphalt systems calculated using molecular simulation[D]. Rhode: University of Rhode Island, 2007.

[29]BALASUBRAMANIAN S. Shear viscosity of linear alkanes through molecular simulations: Quantitative tests for -decane and -hexadecane[J]. Molecular Simulation, 2012, 38(14): 1 234-1 241.

Molecular dynamics simulation of freezing point and viscosity-temperature relationship of aviation kerosene

CHEN Xuejiao, HOU Lei, LI Shiyao
National Engineering Laboratory for Pipeline Safety/ MOE Key Laboratory of Petroleum Engineering /Beijing Key Laboratory of Urban Oil and Gas Distribution Technology, China University of Petroleum-Beijing, Beijing 102249, China

Freezing point and viscosity-temperature relationship are important quality indexes for evaluating low temperature fuidity of aviation kerosene. A kerosene model composing of 5 alkane molecules was established based on a carbon number distribution, carbon/hydrogen mass ratio and aromatic hydrocarbons/saturated hydrocarbons mass ratio. The kerosene model was studied by molecular dynamics simulation. Freezing point was predicted according to the sudden change of density-temperature and diffusion coeffcient-temperature curves and the simulated freezing point is consistent with the experimental value. The radial distribution function was simulated to analyze the microstructure of kerosene modelled at temperatures close to freezing point. The results showed that aggregation between like molecules was the main reason that infuenced the freezing point. Aromatic hydrocarbons had greater infuence than saturated hydrocarbons on the freezing point. The Rouse Model was employed to study the viscosity-temperature relationship of the kerosene model. The simulated viscosity index was identical to test data, showing good agreement between simulated and experimental viscosity changes with temperature.

aviation kerosene; molecular dynamics; freezing point; indirect simulation; viscosity-temperature relationship

2016-11-04

國家自然科學(xué)基金面上項目(51674276)和國家自然科學(xué)基金重點項目(資助U1562212)聯(lián)合資助

10.3969/j.issn.2096-1693.2016.03.046

(編輯 馬桂霞)

陳雪嬌, 侯磊, 李師瑤. 航空煤油冰點及黏溫關(guān)系的分子動力學(xué)模擬. 石油科學(xué)通報, 2016, 03: 493-502

CHEN Xuejiao, HOU Lei, LI Shiyao. Molecular dynamics simulation of freezing point and viscosity-temperature relationship of aviation kerosene. Petroleum Science Bulletin, 2016, 03: 493-502. doi: 10.3969/j.issn.2096-1693.2016.03.046

*通信作者, houleicup@126.com