黨 宇，楊曉東，劉熠斌，馮 翔，楊朝合

(中國石油大學(xué) 重質(zhì)油國家重點實驗室，山東 青島 266580)

八面沸石(FAU)分子篩因具有較高的內(nèi)表面積、適宜的孔結(jié)構(gòu)、固體酸特性及高溫穩(wěn)定性而被用作催化裂化(FCC)反應(yīng)的催化劑[1-2]。目前，催化裂化過程正逐步適應(yīng)加工劣質(zhì)的重油，而劣質(zhì)重油中含有大量雜原子化合物。雜原子化合物通過誘導(dǎo)效應(yīng)來中和催化劑的酸中心，或以形成焦炭覆蓋酸中心的方式使催化劑中毒失活，嚴(yán)重阻礙其它烴類的反應(yīng)[3-4]。

原料分子中的雜原子主要包括S、N、O原子，不同的雜原子化合物對催化裂化反應(yīng)的影響程度不同。研究表明[5]，雜原子具有給電子共軛效應(yīng)且分散正電荷能力強，因此極易在酸作用下質(zhì)子化，引發(fā)開環(huán)或聚合反應(yīng)。與烴類分子相比，雜原子化合物在分子篩上吸附能力強，擴散性能差。考察典型雜原子分子在分子篩上的吸附和擴散行為，對研究劣質(zhì)重油非理想組分阻滯原料轉(zhuǎn)化的作用機理有重要意義。

分子模擬方法已被廣泛用于研究反應(yīng)物分子在多孔材料中的吸附和擴散性能[6-8]。該方法能夠準(zhǔn)確的計算客體分子和分子篩之間的相互作用，得到客體分子的宏觀、微觀性能數(shù)據(jù)，如吸附等溫線、吸附熱和分子在孔道中分布等。目前，眾多學(xué)者采用分子模擬方法計算了雜原子分子在FAU分子篩中的吸附行為，得到了與實驗一致的結(jié)果。沈喜周等[9]采用蒙特卡洛方法(GCMC)考察了噻吩和吡啶在FAU分子篩上的吸附行為，計算得到的吸附熱數(shù)據(jù)和實驗結(jié)果一致；同時，研究發(fā)現(xiàn)噻吩與分子篩孔道的相互作用強于吡啶。楊文平等[10]采用GCMC方法模擬了噻吩類硫化物在FAU分子篩中的吸附性能，結(jié)果表明，噻吩、2-甲基噻吩和2,5-二甲基噻吩在FAU上的單組分吸附等溫線均可用Langmuir吸附模型擬合；此外還發(fā)現(xiàn)，在573.15 K、0.01～1.0 MPa下，3種分子在FAU孔道中存在競爭吸附，對含硫化合物在催化反應(yīng)中的阻滯作用機理提供了微觀解釋。Fu等[11]通過GCMC方法模擬計算了苯和噻吩在FAU分子篩中的競爭吸附，結(jié)果表明，2種分子在FAU超籠中的吸附位隨吸附量的變化而改變；分子篩酸性質(zhì)(硅/鋁比)對2種分子的吸附能和分布有較大影響。

大量工業(yè)數(shù)據(jù)表明，氮化物較強的吸附性會嚴(yán)重影響催化劑的活性，給反應(yīng)過程帶來不利影響[12-16]。部分學(xué)者模擬了氮化物在FAU分子篩中的吸附行為。熊袖章等[17]采用GCMC方法考察了吡啶衍生物在FAU中的吸附行為，得到了吸附等溫線、吸附平衡常數(shù)和吸附熱數(shù)據(jù)，與實驗結(jié)果一致，并發(fā)現(xiàn)不同客體分子的吸附性能與其堿性和結(jié)構(gòu)有關(guān)。Injan等[18]通過DFT方法計算了吡啶在H-FAU分子篩團簇上的吸附，得到了和實驗一致的吸附熱結(jié)果。王迪等[19]考察了吲哚對減壓蠟油(VGO)在FAU催化劑上裂化反應(yīng)的影響，結(jié)果表明，隨著原料中氮含量的增加，反應(yīng)轉(zhuǎn)化率明顯降低、干氣選擇性變差、汽油和液化氣收率顯著下降；分析認(rèn)為是由于吲哚在反應(yīng)中優(yōu)先吸附導(dǎo)致催化劑活性中心減少所致。

有機含氧化合物會導(dǎo)致原油的化學(xué)性質(zhì)不穩(wěn)定，降低油品質(zhì)量，因此氧化物在分子篩中的吸附特性對催化劑性能的影響同樣重要，但關(guān)于油品中氧化物在FAU分子篩的吸附行為的研究還未見報道。同時，有關(guān)雜原子化合物在FAU分子篩孔道中的擴散行為的研究仍存在空白；對S、N、O這3種雜原子化合物在催化劑中吸附和擴散性能的定量對比研究還未見報道。

筆者采用GCMC方法和分子動力學(xué)(MD)方法，在催化裂化溫度(823 K)下，考察了劣質(zhì)重油中典型的五元環(huán)雜原子化合物噻吩、吡咯和呋喃在H-FAU分子篩中的吸附和擴散性能，得到了它們的吸附等溫線、吸附熱、吸附質(zhì)分布及其在孔道中的擴散系數(shù)。模擬結(jié)果為催化裂化過程中不同類型雜原子分子對反應(yīng)影響的研究提供理論指導(dǎo)。

1 計算模型與方法

1.1 分子篩模型

FAU分子篩屬于立方晶系，空間群為Fd-3m，晶胞參數(shù)a=b=c=2.503 nm，α=β=γ=90°，主要結(jié)構(gòu)是八面沸石籠，籠口孔徑為0.74 nm×0.74 nm。模擬過程中，按文獻[20-21]調(diào)整分子篩硅/鋁原子比為5.86，與實際使用分子篩的硅/鋁比一致，并引入H原子平衡電荷，構(gòu)建H-FAU分子篩模型。初始H-FAU分子篩模型的分子式為：H28O384Al28Si164，其中各原子電荷指定如下[20-21]：Si(+0.89)、Al(+0.73)、O(-0.33)、H(+0.083)。模擬中使用1(1×1×1)個單元晶胞，并在x，y，z3個方向添加周期性條件。對更大晶胞進行測試未發(fā)現(xiàn)明顯的尺寸效應(yīng)，表明以上晶胞大小合適。最終的 H-FAU 模型和孔道特征如圖1所示。

圖1 H-FAU分子篩模型及其孔道特征Fig.1 Model and pore structure of H-FAU zeolite(a)H-FAU zeolite model；(b)Pore structure of H-FAU zeoliteAl—Pink；Si—Yellow；O—Red；H—Gray

1.2 雜原子化合物模型

噻吩、吡咯和呋喃3種化合物的三維尺寸如表1 所示。由表1可知，3種分子的最小截面尺寸分別為：噻吩0.55 nm×0.21 nm、吡咯0.55 nm×0.22 nm、呋喃0.51 nm×0.20 nm，均小于FAU分子篩超籠籠口孔徑0.74 nm×0.74 nm。因此，理論上3種分子可以擴散進入FAU分子篩孔道中。采用分子模擬軟件Dmol3模塊對3種分子進行結(jié)構(gòu)優(yōu)化，并計算了各原子的ESP電荷，計算中使用了廣義梯度近似(GGA)的Perdew-Burke-Ernzerhof(PBE)泛函和雙數(shù)值軌道基組(DNP)[22-23]，結(jié)果如圖2所示。

表1 3種雜原子化合物分子的三維尺寸Table 1 The size of three heteroatom molecules

圖2 噻吩、吡咯、呋喃的結(jié)構(gòu)示意圖以及計算ESP電荷Fig.2 Molecular structures and calculated ESP charges of heteroatom molecules(a)Thiophene；(b)Pyrrole；(c)FuranH—White；C—Gray；S—Yellow；N—Blue；O—Red

1.3 參數(shù)設(shè)置

模擬軟件Sorption模塊中抽樣方法選用Metropolis方法，客體分子和分子篩的相互作用包括靜電作用和范德華作用，其中范德華作用采用Lennard-Jones勢函數(shù)，數(shù)學(xué)表達式為：

(1)

式(1)中，ULJ表示客體分子和分子篩之間相互作用；i和j表示不同原子；Rij表示原子間距，nm；Dij和(R0)ij為Lennard-Jones參數(shù)；qi和qj表示原子所帶電荷。模擬中采用Universal力場，靜電相互作用采取Ewald方法處理，非鍵相互作用的加和方法采用Atom based算法，范德華作用截斷距離設(shè)置為1.251 nm，正好小于晶胞邊長(2.503 nm)的一半。計算平衡步數(shù)為106步，生產(chǎn)步數(shù)為107步。

由于缺乏吡咯和呋喃的吸附熱數(shù)據(jù)，但噻吩和吡啶的吸附熱結(jié)果已有詳細(xì)數(shù)據(jù)。因此，為了驗證模擬參數(shù)的合理性，首先采用上述模擬參數(shù)計算了給定實驗條件[24-25]下噻吩和吡啶在FAU分子篩中的Henry常數(shù)，并采用范特霍夫(van’t Hoff)方程計算吸附熱：

Qst=RT2(?lnk/?T)0

(2)

式(2)中，Qst為吸附熱，kJ/mol；R為氣體常數(shù)，J/(mol·K)；k為Henry常數(shù)，T為反應(yīng)溫度，K。計算得到噻吩和吡啶在FAU上的吸附熱分別為93.45 kJ/mol和59.75 kJ/mol，與實驗值80.22～94.08 kJ/mol(噻吩)[25]和57.32～98.74 kJ/mol(吡啶)[24]一致。這表明上述模擬參數(shù)能夠準(zhǔn)確的描述雜原子化合物和分子篩的相互作用。

動力學(xué)性質(zhì)采用Forcite模塊計算，該模塊主要依據(jù)分子動力學(xué)(MD)方法，MD模擬中采用NVE系綜，溫度為FCC反應(yīng)溫度823 K。牛頓運動方程的積分采用蛙跳算法(Velocity-Verlet)，模擬步長1 fs，總模擬時長1 ns，其中前500 ps用于平衡體系。擴散系數(shù)Ds采用愛因斯坦方程計算：

(3)

式(3)中，Ds表示擴散系數(shù)，m2/s；N表示單位晶胞內(nèi)分子數(shù)；r(0)是客體分子質(zhì)心的初始位置坐標(biāo)，r(t)是客體分子在時間t時的質(zhì)心坐標(biāo)，[|r(t)-r(0)2|]是擴散分子均方根位移(MSD)的系綜平均。模擬中3種雜原子化合物在分子篩中的負(fù)載量均為每單元晶胞26個分子。

2 結(jié)果與討論

2.1 吸附等溫線分析

圖3為823 K、100～1000 kPa條件下，噻吩、吡咯和呋喃的單組分在H-FAU中的吸附等溫線。

圖3 823 K下噻吩、吡咯和呋喃在H-FAU分子篩中的吸附等溫線Fig.3 Adsorption isotherm of thiophene,pyrrole and furan in H-FAU zeolite at 823 K

從圖3可以看到，在整個壓力范圍內(nèi)3種分子的吸附量隨壓力的增大緩慢增加。低壓下3種分子吸附量的增速依次為：噻吩、吡咯、呋喃，表明 H-FAU 分子篩對噻吩的吸附能力最強，對呋喃的吸附能力最弱。采用朗格繆爾(Langumir)吸附等溫線模型(公式(4))對3種分子的吸附數(shù)據(jù)進行擬合，結(jié)果如表2所示。

(4)

表2 噻吩、吡咯和呋喃在FAU分子篩上吸附的Langmuir吸附模型擬合參數(shù)Table 2 Adsorption parameters of isotherms of thiophene，pyrrole and furan in FAU zeolite

式(4)中，N表示單位晶胞內(nèi)分子數(shù)；Nm為單元晶胞飽和吸附分子數(shù)；B為平衡常數(shù)，kPa-1。

由圖3和表2的結(jié)果可知，3種分子在FAU分子篩中的吸附數(shù)據(jù)能夠較好的采用Langmuir吸附方程擬合，相關(guān)系數(shù)(R2)均大于0.99。此外，3種分子的吸附等溫線屬于Ⅰ型等溫線，表明呈現(xiàn)典型的單分子層吸附。

表2中Langmuir方程擬合得到的3種分子的飽和吸附量大小順序依次為：噻吩、吡咯、呋喃，且噻吩的飽和吸附量是吡咯的1.17倍，是呋喃的1.21倍。這進一步證明，與吡咯和呋喃相比，噻吩與FAU分子篩有更強的相互作用。

為了考察實際反應(yīng)條件下3種雜原子化合物組分在催化劑中的實際吸附行為，進一步計算了 823 K 下3種組分的等摩爾比混合氣體在分子篩中的吸附等溫線，模擬結(jié)果如圖4所示。由圖4可知，混合氣體中3種組分的吸附量隨壓力的變化趨勢與單組分的吸附情況相似。值得注意的是，在低壓范圍內(nèi)，噻吩的吸附量隨壓力的增大而迅速增加，而吡咯和呋喃的吸附量變化較為平緩。同時，在各自分壓為700 kPa時3種分子的吸附均達到飽和，飽和吸附量的大小次序依次為噻吩、吡咯、呋喃，其中噻吩的飽和吸附量是吡咯的1.28倍，呋喃的1.18倍。這表明3種分子在H-FAU分子篩中存在比較明顯的競爭吸附，其中噻吩的吸附能力最好，呋喃相對較弱。由于3種雜環(huán)分子的結(jié)構(gòu)相似，均為五元環(huán)化合物，且計算的分子三維尺寸接近，因此吸附量和吸附變化趨勢的差異主要是受雜原子性質(zhì)的影響。這說明在823 K、100～1000 kPa條件下，雜環(huán)化合物對催化裂化反應(yīng)的影響由大到小依次為噻吩、吡咯、呋喃。

2.2 吸附熱分析

對于氣體在多孔材料中的吸附過程，吸附熱的其大小可用于定量描述客體分子和多孔材料之間相互作用大小。在GCMC模擬中，等量吸附熱ΔQ可由巨正則系綜中能量/粒子的漲落[26]計算得到：

圖4 823 K下，噻吩、吡咯和呋喃混合氣體在H-FAU分子篩中吸附等溫線Fig.4 Adsorption isotherm of thiophene,pyrrole and furan mixtures in H-FAU zeolite at 823 K

(5)

式(5)中，ΔQ為等量吸附熱，kJ/mol；T為體系溫度，K；R表示理想氣體常數(shù)，J/(mol·K)；UN為吸附相勢能；〈〉為系綜平均。823 K下，3種雜環(huán)分子在H-FAU分子篩中的吸附熱曲線如圖5所示。從圖5中可以看到，噻吩和吡咯的吸附熱均隨壓力的升高而增加。這說明：當(dāng)吸附量較低時，客體分子和分子篩孔道之間的相互作用占主導(dǎo)地位；隨著壓力的增大(即吸附量的增加)，客體分子之間的相互作用增強[27-28]，導(dǎo)致吸附熱增加。此外，從圖5中還可以看到，低壓下噻吩的吸附熱隨壓力變化明顯快于吡咯，說明噻吩的熱效應(yīng)更強。而呋喃的吸附熱曲線表現(xiàn)出與噻吩和吡咯完全不同的變化趨勢，其大小隨壓力的增大而降低。這與很多傳統(tǒng)小分子氣體的吸附呈相同的規(guī)律，表明在低壓下H-FAU分子篩對呋喃分子的捕捉效率較高，使得呋喃分子優(yōu)先占據(jù)活性位，放出較多熱量。但是隨著壓力增大，孔道表面活性位被占據(jù)，呋喃的吸附轉(zhuǎn)向其它普通的吸附位，導(dǎo)致吸附熱下降。這表明在實際反應(yīng)較高壓力下，呋喃對分子篩的影響要小于噻吩和吡咯，這一點也可以由飽和吸附時3種分子的吸附熱大小次序(噻吩、吡咯、呋喃)得到證明。

圖5 噻吩、吡咯和呋喃在H-FAU分子篩中的等量吸附熱Fig.5 Isoteric heat of heteroatom molecules in H-FAU zeolite

2.3 吸附質(zhì)分布

蒙特卡洛模擬可以得到一系列吸附質(zhì)分子的構(gòu)型數(shù)據(jù)，分析這些數(shù)據(jù)可以得到客體分子在分子篩孔道中各個吸附位點的出現(xiàn)概率，從而有助于我們進一步認(rèn)識不同雜原子化合物在H-FAU分子篩中的分布情況。

823 K時，3種雜環(huán)化合物單組分在H-FAU分子篩中概率密度分布如圖6所示。圖6中，不同顏色表示概率高低，紅色區(qū)域為高密度分布區(qū)，表示分子在該區(qū)域出現(xiàn)的概率較大；藍色區(qū)域為低密度分布，表示分子出現(xiàn)概率較低。由圖6可知，3種雜原子分子在分子篩中的分布類似，主要分布在超籠中，在方鈉石籠(SOD)中也有少許分布。與吡咯和呋喃比較，噻吩在超籠中的分布更加集中(紅色區(qū)域)，同時在SOD籠中的集中區(qū)域也更明顯。吡咯和呋喃的密度分布比較接近，表明它們在分子篩孔道中的分布類似。這與它們在分子篩中的飽和吸附量和吸附熱大小的規(guī)律一致。

圖6 噻吩、吡咯和呋喃在FAU分子篩中的分布概率密度圖Fig.6 Probability density distribution of heteroatom molecules in H-FAU zeolite(a)Thiophene；(b)Pyrrole；(c)Furan

基于上述吸附性質(zhì)的研究結(jié)果我們發(fā)現(xiàn)，3種雜原子化合物在分子篩中的飽和吸附量、吸附速率、吸附熱大小和孔道中分布密度由大到小的次序均依次為噻吩、吡咯、呋喃，表明在催化裂化條件下3種分子對催化劑分子篩的影響大小依次為噻吩、吡咯、呋喃。

2.4 吸附過程的動力學(xué)模擬

由于3種雜原子分子的最小截面尺寸遠(yuǎn)小于分子篩的孔徑，因此雜原子分子在孔道中的擴散性能主要受雜原子類型的影響。在雜原子分子擴散性能的研究中，擴散系數(shù)是描述擴散能力、研究傳質(zhì)過程和計算傳質(zhì)速率的重要參數(shù)。筆者采用分子動力學(xué)方法計算得到了823 K下噻吩、吡咯和呋喃在FAU分子篩中的擴散軌跡，并通過分析擴散軌跡得到了3種分子的均方根位移(MSD)隨時間的變化曲線，如圖7所示。

圖7 噻吩、吡咯和呋喃的MSD曲線Fig.7 MSD curves of heteroatom molecules in H-FAU zeolite

根據(jù)愛因斯坦方程(公式(3))計算得到了3種分子在H-FAU分子篩中的擴散系數(shù)，見表3。由表3可知，3種雜原子分子在H-FAU分子篩中的自擴散系數(shù)大小順序依次為：呋喃、吡咯、噻吩，即呋喃在分子篩孔道中的擴散阻力最小，噻吩的擴散阻力最大。這與3種分子的吸附模擬結(jié)果一致，噻吩的吸附能力強，導(dǎo)致噻吩和孔道的相互作用最大，擴散阻力也最大；而呋喃的吸附能力相對較弱，在孔道中受到束縛相對較小，擴散阻力小于噻吩和吡咯。

表3 噻吩、吡咯和呋喃的擴散系數(shù)Table 3 Diffusion coefficients of thiophene,pyrrole and furan

為了進一步考察3種雜原子化合物在H-FAU分子篩上的作用位點，首先計算了3種分子的電荷分布，結(jié)果如圖8所示。靜電勢能是指將單位正電荷從無窮遠(yuǎn)處移到分子周圍某處做的功，其中原子核貢獻正值，電子貢獻負(fù)值。從圖8中可以看到，3種分子的電荷分布具有一定差異，其中噻吩和呋喃的電荷主要分布在分子環(huán)區(qū)上、下兩側(cè)以及S、O原子附近，并且可以看到，噻吩的電荷分布更加離域，而呋喃的電荷則分布在O原子附近，分布相對局域。這主要是因為O原子的電負(fù)性(3.50)比S原子的電負(fù)性(2.44)強，導(dǎo)致O原子的吸電子誘導(dǎo)效應(yīng)更大[5]。吡咯分子電荷主要分布在分子環(huán)區(qū)兩側(cè)。這是因為分子中N—H鍵偶極矩的方向與p-π共軛方向相同所致。由于化合物分子和分子篩之間的相互作用本質(zhì)上是電荷的相互作用，因此根據(jù)電荷分布可以預(yù)測，3種雜原子化合物分子與分子篩相互作用主要通過分子環(huán)區(qū)和雜原子來實現(xiàn)。

圖8 噻吩、吡咯和呋喃分子中的靜電勢能分布圖Fig.8 Electrostatic potential density distribution of heteroatom molecules(a)Thiophene；(b)Pyrrole；(c)Furan；EP—Electronic potential

在擴散計算中，徑向分布函數(shù)(RDF)表示粒子在空間中的密度變化，是局域密度和平均密度的比，以g(r)表示，r為雜原子或雜原子分子質(zhì)心與H質(zhì)子間的距離。筆者采用RDF法，對3種雜原子化合物潛在作用位點(分子環(huán)區(qū)和雜原子)與H-FAU分子篩孔道結(jié)構(gòu)上的活性位(B酸位)的H質(zhì)子的間距分布進行統(tǒng)計分析，計算結(jié)果如圖9所示。

圖9中H-Center表示客體分子的質(zhì)心(Center)在H質(zhì)子周圍的分布；H-X(X=S、O、N)則表示客體分子上的雜原子在H質(zhì)子周圍的分布情況。由3種分子質(zhì)心的RDF圖(圖9(a)、(b)、(c))可知，主峰所對應(yīng)的距離均小于0.4 nm，表明存在分子環(huán)區(qū)的離域電子與孔道上的H質(zhì)子相互作用；主峰位置所對應(yīng)的距離在成鍵距離范圍內(nèi)說明，雜原子分子和H質(zhì)子之間可能發(fā)生電荷轉(zhuǎn)移。這需要進一步的量化計算。另一方面，由雜原子的RDF圖(圖9(d)、(e)、(f))可以看到，最強峰的位置同樣小于0.4 nm，表明雜原子和H質(zhì)子之間存在較強相互作用。值得注意的是，呋喃雜原子的RDF峰位置在0.3 nm附近，比噻吩和吡咯的RDF峰位置低，表明分子篩的活性位點(吸附位)對呋喃分子有較強的作用。這主要是因為呋喃分子中的電荷部分局域在電負(fù)性較強的O原子附近所導(dǎo)致的，同時這也解釋了為什么低壓下呋喃的吸附熱較噻吩和吡咯的大。綜上所述，3種雜原子分子和分子篩上H質(zhì)子之間存在2種作用形式：(1)分子環(huán)區(qū)的離域電子和H質(zhì)子作用；(2)分子上的雜原子和H質(zhì)子作用。

圖9 噻吩、吡咯和呋喃分子與H質(zhì)子作用的RDF圖Fig.9 RDF diagram of thiophene,pyrrole and furan molecules interacting with H protons(a)Thiophene,H-Center distance;(b)Pyrrole,H-Center distance;(c)Furan,H-Center distance;(d)Thiophene,H-S distance;(e)Pyrrole,H-N distance;(f)Furan,H-O distance

3 結(jié) 論

在催化裂化反應(yīng)條件823 K、100～1000 kPa下，3種雜原子化合物分子在分子篩上的飽和吸附量大小、與分子篩的相互作用強弱、對催化裂化反應(yīng)過程影響由強到弱的順序都依次為噻吩、吡咯、呋喃；而其在分子篩上的擴散能力強弱順序依次為呋喃、吡咯、噻吩。

噻吩和吡咯在FAU分子篩上的吸附熱隨壓力的升高而增加，噻吩的熱效應(yīng)比吡咯更強；呋喃表現(xiàn)出與噻吩和吡咯不同變化趨勢，其在FAU分子篩上的吸附熱隨壓力的升高而降低，直至吸附平衡時不再變化。

3種雜原子化合物分子在FAU分子篩的超籠和SOD籠結(jié)構(gòu)中均有分布，且噻吩在超籠中的分布較吡咯和呋喃更加集中。3種雜原子化合物分子與分子篩上的H質(zhì)子相互作用主要通過分子環(huán)區(qū)的離域電子和雜原子來實現(xiàn)。