多環芳烴在超臨界環己烷中的溶解:分子動力學模擬研究

趙向勃, 付 峰, 陳 權, 高曉明, 閆 挺

(1. 延安大學 化學與化工學院, 延安 716000; 2. 廈門海關技術中心, 廈門 361026)

1 引 言

近年來，盡管太陽能、氫能、風能等新能源發展迅速，但傳統化工原料仍然在人類的社會發展中發揮著重要的作用. 原油經過長期的過度開采及使用，輕質油含量在逐漸減小. 而重油和超重油儲量已占據我國現有石油儲量的50%以上[1]. 重油具有粘度高、雜原子含量高的特點，傳統的加氫脫碳工藝在重油加工中存在諸多局限性[2-4]. 脫碳過程會產生大量殘焦，造成資源浪費[5]. 而重油中的S、N等雜原子在加氫過程中會使催化劑失活. 超臨界流體(SCFs)技術是一種新型的加氫脫碳工藝，可以在不使用催化劑的情況下使得重油改質，抑制焦炭的生成[6]. 超臨界流體已經廣泛用于重質油輕質化、減黏等方面的研究[7-11]. 超臨界水(SCW)具有低的介電常數和良好的擴散性能，使得重油在SCW中有著良好的溶解性能，再結合潛在的供氫能力，重油在超臨界流體中的輕質化與減黏主要集中在SCW領域中[12-14]. 然而，近年來研究發現SCW在重油熱解的過程中的供氫效果并不明顯，SCW供氫能力有待深入了解[15]. 同時，Zhu等[16，17]基于量子力學的理論計算否定了SCW在涉及碳氫自由基的反應中給予氫自由基的可能性，所以越來越多的人認為SCW在重油熱解過程中主要作為溶劑發揮了物理作用.

然而，水和烴類的混合體系表現出復雜的相行為. 對其復雜相行為的檢測和預測一直是研究的重點[18]. 目前，水和眾多低分子量碳氫化合物的相行為已經有了較大的研究[19]. 然而，重烴混合物的研究較少，因為由于混合物的不透明、乳化現象均可能干擾相邊界測量. 此外，烴類高溫高壓下的熱解過程采樣困難，對這類混合物的相行為測量具有挑戰性. 在之前的研究中，根據Van Konynenburg和Scott分類方案[20]，Athabasca瀝青+水的偽二元混合物表現出Ⅲb型相行為[21，22]，并測量了瀝青在水中的溶解度. 結果表明在含水的瀝青混合物中加入溶劑時，相對于水和瀝青，溶劑通常具有較低的密度和較低的臨界溫度，相對于瀝青具有較高的水溶性. 預測富烴液相中的混合相行為和水含量會發生顯著變化.

隨著對重油-超臨界流體體系相結構、反應機理和反應動力學的研究，研究者認識到不同的相結構及其形成的擴散環境對超臨界流體中重油的非催化熱解起著至關重要的作用. Morimoto等[23]基于漢森溶解度參數(HSP)提出了SCW的靜態介電常數大于2.2，HSP氫鍵組分小于10 MPa0.5是SCW與重油良好混相的必要條件. He等[24，25]利用基團貢獻法，分析了273 - 644 K、26 MPa壓力下原油和水混合物的相平衡，在壓力為25 MPa、溫度為640 - 660 K的條件下，研究了軟瀝青、瀝青質和水的混合物的相平衡，揭示了單相和兩相區域的存在. 而Guo等[26]通過混相指數解釋了SCW中重油的改質，發現混相指數的降低導致軟瀝青產量增加的同時瀝青質和殘焦的收率減少. 在恒定的混相指數下，溫度增加可以導致軟瀝青產率的增加但是殘焦的產率也會有所提升，即在不結焦的高溫條件下最大限度的降低混相指數可以提高軟瀝青的產率. 所以重油-超臨界流體體系的混相可以揭示重油-超臨界流體體系中的相結構和溶解行為. 分子動力學模擬(MD)可以研究超臨界狀態下的復雜系統[27，28]. Jian等[29]通過MD模擬研究了SCW的結構和擴散性能. 從室溫到超臨界溫度的升高對水中氫鍵的影響比沿著超臨界等溫線密度的變化要明顯得多. 計算證實，隨著溫度的升高和密度的降低，每個分子的平均氫鍵數和團簇大小都在減少. Xin等[30]通過分子模擬研究了瀝青質在超臨界流體中的溶劑化反應，結果表明排斥Van der Waals相互作用對瀝青質與水分子相互作用的性質和程度有主要影響.

SCW是重油改質中常用的流體，但由于多環芳烴與SCW的相互作用表現為斥力，不利于多環芳烴的溶解. 而超臨界環己烷與多環芳烴之間的相互作用為引力. 因此，本文利用分子動力學方法研究了不同環數的多環芳烴以及多環芳烴混合物在超臨界環己烷中的溶解過程，揭示了多環芳烴在超臨界環己烷中的溶解機理. 首先采用MD模擬方法研究了純多環芳烴或多環芳烴混合物在超臨界環己烷中的溶解行為，考察了多環芳烴的尺度、多環芳烴混合物的組成對溶解過程的影響. 通過對多環芳烴在超臨界環己烷中的短程溶劑結構和溶劑化自由能的研究，表征了多環芳烴與超臨界環己烷環境的相互作用. 然后，通過計算多環芳烴間的內聚能密度來揭示了多環芳烴本身之間的相互作用. 本文將重油在超臨界流體中改質中超臨界流體的篩選，操作范圍的確定提供理論指導.

2 模擬與理論方法

首先MD仿真的可靠性取決于所用的力場，COMPASS力場在烴類/SCW體系中的有效性已被證實[31，32]，因此本文所用力場采用COMPASS力場應用于多環芳烴與超臨界環己烷體系，所采用的計算均基于Materials Studio 2019軟件的Forcite模塊. 由于重油當中所含多環芳烴的結構變化較多，因此本文只考慮沒有雜環摻入和烷基取代基存在的多環芳烴. 在本次研究中，選擇2環的萘(NAP)為輕多環芳烴，選擇6環的苯并芘(Bghip)為重多環芳烴如圖1所示. 這些多環芳烴的詳細信息如表1所示.

表1 多環芳烴模型性質

圖1 多環芳烴分子結構Fig. 1 Molecular structures of different PAHs.

在模擬中，首先構建充滿超臨界環己烷分子的立方體溶劑盒來模擬超臨界環己烷環境，如圖2所示，環境溫度為634 - 674 K，密度為0.1 g/cm3- 0.3 g/cm3，當模擬多環芳烴的溶解時，將含有純多環芳烴或多環芳烴混合物的油滴放在溶劑盒的中心. 在研究多環芳烴的近程溶劑結構和溶劑化自由能時，將一個多環芳烴分子放入隨機溶劑盒中. 溶劑盒的剩余空間裝滿了特定狀態下的環己烷分子.

圖2 超臨界環己烷溶劑盒Fig. 2 Schematic diagram ofSC-cyclohexane cube box.

多環芳烴在超臨界環己烷中溶解的典型動態模擬過程如下所述：首先，構建純多環芳烴的本體結構，對溶劑盒進行結構優化并通過NPT系綜的退火操作獲得了純多環芳烴的穩定本體結構. 從多環芳烴本體結構中取出一個半徑為15 ?的油滴分子，將其放置在一個立方體溶劑盒的中心，立方體溶劑盒邊長為80 ?. 此外，溶劑盒中充滿了環己烷分子，其密度由環己烷分子的數量來調節. 環己烷密度從0.1 g/cm3到0.3 g/cm3不等. 然后，基于NVT系綜進行持續時間為500 ps的分子動力學模擬，以保證動力學平衡. 在模擬過程中，溫度由NHL恒溫器調節. 靜電法和范德華求和法分別設置為van der Waals based和atom based.

3 結果與討論

3.1 多環芳烴在超臨界環己烷中的溶解

3.1.1溶劑密度對多環芳烴溶解度的影響

在674 K的溫度下，采用MD模擬了NAP和Bghip的多環芳烴液滴在0.10 -0.30 g/cm3三種超臨界環己烷密度下的溶解情況. 模擬時間在500 ps時的結果如表2所示.

表2 多環芳烴在不同密度超臨界環己烷中的溶解情況

溫度674 K下不同密度的壓力在0.005 GPa - 0.008 GPa之間，根據觀察，在不同密度的環己烷中，NAP油滴分子能夠迅速地均勻分散在溶劑盒中. 隨著環己烷密度的增大，多環芳烴與環己烷的體系的混相指數逐漸減小. 同時，環己烷密度的變化對NAP和Bghip的溶解有重要影響. 當環己烷密度為0.10 g/cm3時，NAP分子快速溶解于環己烷相中，球形油滴形態嚴重偏離最初的形態. 增加環己烷密度至0.30 g/cm3時，油滴擴散的效果較好，但是相比于0.10 g/cm3的超臨界環己烷體系仍然存在NAP分子的聚集. 與此同時，對于Bghip，在0.10 g/cm3時，Bghip分子迅速溶解于超臨界環己烷中. 但是相比于NAP分子溶解性稍差，在0.30 g/cm3模擬結束時，溶劑盒中只發現少量解離的Bghip分子.

從圖中可以看出NAP油滴可在超臨界環己烷中完全分離，但仔細觀察表2的圖像可以發現，在溶劑箱中仍存在NAP分子的區域富集，表明NAP發生了團聚. 由于富集區域的位置在整個模擬過程中隨機變化，因此在超臨界環己烷中形成的NAP團聚體是不穩定的. 而當重多環芳烴溶解到環己烷相時，Bghip在溶劑盒中出現穩定的團聚，其主要以二聚體、三聚體等類似于瀝青質π-π堆積的形式存在.

3.1.2溫度對多環芳烴溶解度的影響

在溫度為634 - 674 K、超臨界環己烷密度為0.10-0.30 g/cm3時，采用MD模擬方法研究了NAP和Bghip液滴在超臨界環己烷中的溶解情況，模擬時間500 ps情況下的結果如表3和表4所示.

表4 苯并芘在不同溫度下的溶解行為

當環己烷密度為0.30 g/cm3時，溫度升高20/40 K對NAP液滴溶解的影響有輕微的提高. NAP分子從油滴中溶解在環己烷相中的數量基本相同. 在環己烷密度為0.10 g/cm3時，在634和654 K溫度下也能觀察到類似的現象，但隨著溫度的進一步升高，情況發生了較大的變化. 在674 K溫度下，大量的NAP分子會以單體的形式溶解到環己烷相中.

當環己烷密度為0.30 g/cm3時，溫度升高20/40 K對Bghip液滴的溶解與NAP分子同樣具有明顯的影響. 大量的Bghip分子從油滴中解離并溶解在環己烷相中. 在環己烷密度為0.10 g/cm3時，在634 K和654 K溫度下也能觀察到類似的現象. 但隨著溫度的進一步升高，發生了的變化. 在674 K溫度下，大量的Bghip分子以單體、二聚體或更大規模的聚集體的形式溶解到環己烷相中，此時Bghip液滴的剩余部分可以看作是一個大尺度的團簇，類似于重油熱解過程中的結焦前驅體.

從表3和表4的圖像可以看出，較低環己烷密度與適宜的溫度有利于NAP、Bghip在超臨界環己烷中的溶解.

3.1.3多環芳烴混合物在超臨界環己烷中的溶解

上文對純NAP分子和Bghip分子的溶解行為進行了研究，而多環芳烴具有不同尺度，組成差異較大的特點，采用MD模擬方法研究了多環芳烴混合物即Bghip80NAP20(下標代表油滴中苯并芘與萘的比例)混合油滴在超臨界環己烷中的溶解行為，如表5所示.

表5 Bghip80NAP20油滴在不同密度超臨界環己烷中的溶解行為

相比較于純Bghip與NAP體系，NAP與Bghip混合物體系在超臨界環己烷環境下的溶解有重要影響. 模擬500 ps后，Bghip80NAP20的油滴在溶劑箱中有明顯的聚集成團簇狀. 在不同溫度下，溫度越高Bghip80NAP20的含量相較于0 ps時的含量有著明顯的降低，說明溫度升高會促進多環芳烴混合物在環己烷中的溶解，Bghip80NAP20油滴大部分完全分離，只有小部分的聚集. Bghip初始含量有了明顯的降低，促進了多環芳烴混合物在環己烷相中的溶解，油滴發生完全解離. 環己烷相中以NAP為主要成分的多環芳烴分子數量增加，從而減小了油滴的尺寸. 導致NAP和Bghip在溶劑箱中均以單體或聚集體的形式擴散. 通過對比表456的圖像，輕PAHs比重PAHs在超臨界環己烷中的溶解性更好.

3.1.4多環芳烴在超臨界環己烷中的溶解過程

從表6中可以看出，環己烷密度的增加對多環芳烴混合物在超臨界環己烷中的溶解有不利的影響. 當環己烷密度為0.10 g/cm3時，溶劑箱中僅存多環芳烴單體和團聚體，而當環己烷密度為0.20 g/cm3或0.30 g/cm3時，溶劑箱中雖然也存在多環芳烴單體和團聚體，但是相比于0.10 g/cm3來說，溶解性似乎更差. 并且在低環己烷密度下還是在高環己烷密度下進行模擬，都能觀察到類似的現象，即輕多環芳烴優先溶解到環己烷相中. 以環己烷密度為0.10 g/cm3時的模擬為例. 僅經過125 ps的模擬，就有相當數量的Bghip80NAP20分子優先溶解到環己烷相中. 而0.20 g/cm3和0.30 g/cm3在125 ps明顯可以看出Bghip80NAP20分子存在明顯的聚集現象.

表6 674 K下Bghip80NAP20在超臨界環己烷中的溶解過程

3.2 多環芳烴與超臨界環己烷相互作用

3.2.1徑向分配函數

多環芳烴與超臨界環己烷的混溶是由溶質-溶劑與溶質-溶質相互作用決定的. 超臨界環己烷中多環芳烴的近程溶劑結構和溶劑化自由能反映了體系中溶質-溶劑間的相互作用. 為了解釋多環芳烴的近程溶劑結構，計算了674 K溫度下、0.10 - 0.30 g/cm3環己烷密度中gBghip- H(r)和gNAP- H(r)的徑向分布函數(RDF)如圖3和圖4所示.

圖3 674 K下0.10 - 0.30 g/cm3密度范圍中超臨界環己烷中Bghip徑向分配函數Fig. 3 gBghip-H(r)in SC-cyclohexane at temperature of 674 K and cyclohexane densities from 0.10 to 0.30 g/cm3

圖4 674 K下0.10 - 0.30 g/cm3密度范圍中超臨界環己烷中NAP徑向分配函數Fig. 4 gNAP-H(r)in SC-cyclohexane at temperature of 674 K and cyclohexane densities from 0.10 to 0.30 g/cm3

在三種環己烷密度下，0.10 g/cm3的gBghip-H(r)與其他兩者分布不同. 在Bghip分子的油滴半徑內，gBghip-H(r)的值保持為0. 在2.1 ?以上半徑處，gBghip-H(r)的值逐漸增大，而在0.20 g/cm3和0.30 g/cm3時在10.0 ?附近開始接近1，而后持續保持. 而Fujisawa和Kajimoto等[33，34]根據光譜表征，得出如果溶質和溶劑之間存在強吸引的非鍵相互作用，則容易在溶質周圍形成溶劑化殼. 在0.10 g/cm3時候gBghip-H(r)發現了對應溶劑殼存在的峰，則可以認為多環芳烴與超臨界環己烷間具有較強的相互作用，0.10 g/cm3的環己烷可作為催化劑促進Bghip的溶解，其詳細的機理值得進一步研究. 而0.20 g/cm3和0.30 g/cm3則沒有產生溶劑峰，說明高密度下重多環芳烴與超臨界環己烷相互作用較弱. 而在三種密度下的gNAP- H(r)的RDF分布一致，在2.3?處開始增加，而后在10.0 ?附近開始趨向于1，繼而保持. 說明輕多環芳烴與超臨界環己烷的相互作用較弱，其結果與SCW有著相似的結果.

3.2.2溶劑化自由能

計算了NAP和Bghip的溶劑化自由能，如圖5所示. 634 K時，NAP溶解于超臨界環己烷相中，其溶劑化自由能絕對值為4.7989 kJ/mol. Bghip分子的溶劑化自由能絕對值為21.0463 kJ/mol，其絕對值明顯大于NAP分子. 當環己烷密度為0.30 g/cm3時，NAP的溶劑化自由能絕對值為7.1729 kJ/mol，高濃度的環己烷密度對多環芳烴分子在超臨界環己烷環境中的溶解度有消極影響. 將Bghip和NAP在三種不同環己烷密度的超臨界環己烷中溶解，結果表明，環己烷密度的變化影響多環芳烴的溶解度. 隨著環己烷密度的增加，溶劑化自由能的絕對值顯著增加，在0.30 g/cm3的環己烷密度下，溶劑化自由能絕對值約為7.729 kJ/mol和16.8161 kJ/mol. 另外，在相同溫度不同環己烷密度下，Bghip中的溶劑化自由能絕對值相對于超臨界環己烷環境中的輕多環芳烴的溶劑化自由能絕對值較大. 而溶劑化自由能絕對值越小說明溶解性越好，所以在低密度下和高溫下超臨界環己烷對NAP和Bghip有著良好的溶解性.

(a)不同溫度下的溶劑化自由能(a)Solvation free energy at different temperatures

盡管多環芳烴的尺度和環己烷的熱力學狀態發生了變化，但溶劑化自由能的計算值始終為負值. Bghip分子的溶劑化自由能的絕對值總是比NAP分子更大，證實了在超臨界狀態下重質多環芳烴比輕多環芳烴更加穩定，同樣也意味著NAP在超臨界環己烷中的溶解度應該比Bghip效果好.

3.2.3內聚能密度

由前文的MD模擬結果可知，多環芳烴與超臨界環己烷環境下的混相指數同時由另一個因素決定，即多環芳烴分子間的相互作用. 為了了解多環芳烴與超臨界環己烷環境之間的相互作用，本文計算了NAP和Bghip在超臨界環己烷環境下對應溫度和壓力下的內聚能密度(Cohesive energy density，CED)，計算得到的CED如圖6所示.

(a)不同密度下多環芳烴內聚能密度(a)CED of PAHs at different densities

超臨界環己烷環境下中多環芳烴的CED越高，代表著多環芳烴的溶解越困難. 在相同環己烷密度下6環的Bghip的CED和2環NAP的CED相差不大，并且伴隨溫度的升高兩者的CED在逐漸減小，說明在超臨界環己烷環境下對NAP和Bghip都有著不錯的溶解能力，并且溫度的升高對于NAP和Bghip的溶解性有所提高，這與其前文的MD結果相對應. 而對于不同密度來說，0.30 g/cm3的Bghip分子和NAP分子是0.10 g/cm3的Bghip分子和NAP分子的8倍左右. 這也證實Bghip分子比NAP分子更加穩定且不容易被溶解.

隨著環己烷溫度的增加，Bghip在超臨界環己烷環境之間的吸引作用增強. 相比之下，Bghip的內聚能密度對環己烷密度變化更為敏感. 在674 K溫度下，環己烷密度從0.10 g/cm3增加到0.30 g/cm3，導致Bghip的CED大大增加了2.64E7 j/cm3. 當環己烷密度為0.30 g/cm3時，隨著溫度從634 K升高到674 K，Bghip的CED只降低了2.2E5 j/cm3，表明多環芳烴之間的相互作用增加. 并且隨著環己烷密度的增大超臨界環己烷對NAP和Bghip的溶解能力減弱，而溫度的增加，會對超臨界環己烷環境下NAP和Bghip的溶解有著積極的影響.

4 結 論

本文采用分子動力學模擬方法研究了多環芳烴在超臨界環己烷下的溶解度行為，有助于了解重油在超臨界環己烷中的溶解行為. 結果表明在不同溫度與密度下超臨界環己烷對于具有2個芳香環的NAP(輕多環芳烴)均具有良好的溶解性. 而溫度對于含有6個芳香環的BghiP的油滴(重多環芳烴)在超臨界環己烷中的溶解影響較小. 在由多環芳烴混合物組成的油滴溶解過程中，輕環芳烴優先溶解到環己烷相，重環芳烴則集中在粒徑減小的油滴中. 同時，溫度升高以及密度的將降低會導致多環芳烴的內聚能密度降低，增加了多環芳烴與超臨界環己烷之間的相互作用. 因此，環己烷密度越小以及溫度的升高可以促進多環芳烴或多環芳烴混合物在超臨界環己烷中的溶解.