CO2置換CH4水合物技術(shù)中主、客體分子間作用的DFT研究

陳 浩 顏克鳳* 李小森 陳朝陽 張 郁 徐純剛

(1.中國科學(xué)院 廣州能源研究所， 廣州 510640； 2.中國科學(xué)院 廣州能源研究所 中國科學(xué)院天然氣水合物重點實驗室，廣州 510640； 3.中國科學(xué)院 廣州天然氣水合物研究中心， 廣州 510640；4.廣東省新能源和可再生能源研究開發(fā)與應(yīng)用重點實驗室， 廣州 510640；5.中國科學(xué)院大學(xué)， 北京 100049)

引 言

天然氣水合物是以CH4、C2H6、CO2等氣體分子為客體，以H2O分子為主體，在低溫高壓下形成的一種籠型晶體化合物，其中主體分子以氫鍵相互連接形成不同的籠型孔穴，客體分子占據(jù)孔穴。天然氣水合物的晶體結(jié)構(gòu)主要包括I型(SI)、II型(SII)和H型(SH)。I型水合物晶胞是由2個512小孔穴和6個51262大孔穴總共46個H2O分子形成的體心立方結(jié)構(gòu)；II型水合物晶胞是面心立方結(jié)構(gòu)，包含136個H2O分子，由8個51264大孔穴和16個512小孔穴組成；H型水合物晶胞是簡單的立方結(jié)構(gòu)，包含34個H2O分子，晶胞中有3種不同的孔穴，即3個512孔穴、2個435663孔穴和1個51268孔穴。在自然界中天然氣水合物主要以SI晶體結(jié)構(gòu)存在，客體分子主要為CH4。由于具有儲量大、熱值高、無污染等優(yōu)點，天然氣水合物作為21世紀(jì)的一種新能源，其開采研究對于我國后續(xù)能源的可持續(xù)發(fā)展具有重要意義。目前，天然氣水合物的開采方法主要有熱激法、降壓法和化學(xué)試劑法，各個方法均有一定的局限，因此開采方法的優(yōu)化，以及新型開采方法的研發(fā)是天然氣水合物開采技術(shù)研究的重點。CO2置換開采水合物方法是一種新型的開采方法，在從水合物中置換出天然氣的同時注入CO2從而生成CO2水合物。該方法在開采出水合物資源的同時存儲CO2，可有效地減少碳排放，降低溫室效應(yīng)[1]。同時，CO2水合物的形成可有效防止水合物分解引發(fā)的地質(zhì)問題。因此，CO2置換開采天然氣水合物技術(shù)被認為是可用于水合物開采和減少碳排放的具有巨大應(yīng)用潛力的一種手段。

實驗和理論研究表明，注入CO2置換天然氣水合物中的CH4是可行的[2]。從熱力學(xué)上看，CO2水合物體系的平衡條件要比CH4水合物體系的平衡條件溫和[3]。從動力學(xué)上看，當(dāng)溫度為250 K時，CO2置換CH4水合物過程的吉布斯自由能ΔG是負值，這意味著置換反應(yīng)是自發(fā)進行的[4]。Yuan等[5]利用三維反應(yīng)器對CO2置換CH4水合物過程的有利條件進行研究，發(fā)現(xiàn)具有游離氣的水合物儲層適合置換開采。對CO2置換CH4水合物過程的微觀分析和組分測量表明，置換主要發(fā)生在水合物相中[6]，CO2水合物和CH4水合物同時存在于水合物中，置換過程不存在CH4水合物的解離[7]。有部分研究表明CH4水合物中只有51264籠中的CH4被置換，因此置換效率較低[8]。為提高置換效率，探究過程中的置換機理成為CO2置換CH4水合物技術(shù)研究的重點。目前，宏觀實驗和微觀測試的研究能夠獲得置換過程的熱力學(xué)和動力學(xué)特性，但還不足以解釋置換機理。計算機模擬是在微觀尺度上揭示置換機理的一種有效方法。近年來，研究者們將分子動力學(xué)(molecular dynamics，MD)模擬用于研究CO2置換CH4水合物過程的置換機理，模擬結(jié)果表明，置換過程可以分為3個步驟[9]:首先打破水合物籠子，然后CH4分子跑出籠子，同時CO2分子進入空籠并占據(jù)籠子進行置換。Liu等[10]利用MD模擬發(fā)現(xiàn)置換出的CH4分子在形成氣泡的過程中并沒有出現(xiàn)漆膜，而且CO2分子的存在能夠促使CH4分子形成氣泡，防止CH4重新生成水合物，從分子角度闡述了CO2的置換機理。MD模擬能從分子尺度解析置換的動態(tài)過程，量子力學(xué)(quantum mechanics，QM)模擬則能從量子尺度，通過從頭算(Ab initio calculation)以及密度泛函理論(density functional theory，DFT)的方法對水合物主客體分子間的電子效應(yīng)進行分析，從量子層面解釋置換機理[11]。Geng等[12]利用MD 模擬和QM模擬對CH4、CO2以及CH4- CO2水合物的穩(wěn)定性進行研究，通過分子尺度的徑向分布函數(shù)、均方位移、擴散系數(shù)以及量子尺度的結(jié)合能的對比分析，發(fā)現(xiàn)CH4- CO2混合水合物的穩(wěn)定性最高，此結(jié)果與Ota等[13]的實驗結(jié)果相一致。Liu等[4]利用從頭算方法對CO2取代CH4占據(jù)水合物籠子的結(jié)合能變化及熱力學(xué)特性進行研究，發(fā)現(xiàn)CO2占據(jù)51262籠在熱力學(xué)上的有利性。

由上述研究可知，水合物中主、客體分子對水合物的穩(wěn)定性具有重要的影響。為更深入地闡述CO2置換CH4水合物過程的本質(zhì)，本文采用DFT方法對水合物中主、客體雙分子聚體之間的相互作用進行模擬，通過對稱性匹配微擾理論(symmetry adapted perturbation theory，SAPT)進行能量分析，解析主、客體相互作用中靜電、色散、誘導(dǎo)和交換4種相互作用各自的貢獻。同時計算了波函數(shù)信息，分析了約化密度梯度函數(shù) (reduced density gradient，RDG)、獨立梯度模型(independent gradient model，IGM)、靜電勢等參數(shù)，探討主、客體對水合物籠子穩(wěn)定性影響的作用機理，為深層次地闡述置換機理提供理論依據(jù)。

1 模型及計算方法

本文利用DFT方法對水合物中主、客體分子聚體，即H2O- H2O、H2O- CH4、H2O- CO2雙分子聚體進行研究，針對它們的幾何結(jié)構(gòu)、相互作用能量、相互作用方式進行計算模擬。模擬通過幾何優(yōu)化確定雙分子聚體的最優(yōu)結(jié)構(gòu)，然后計算單點能，采用式(1)計算雙分子聚體的相互作用能。

ΔEa-b=Etotal-Ea-Eb

(1)

式中，ΔEa-b表示雙分子聚體中分子之間的相互作用能，a、b分別表示聚體中的分子；Etotal表示聚體的總能量；Ea和Eb分別表示雙分子聚體中各個分子的能量。模擬結(jié)構(gòu)優(yōu)化及能量計算都在Gaussian16程序C.01版本[14]上進行。

本文首次采用高階SAPT[21]對水合物中主、客體雙分子聚體進行能量分解，評估雙分子聚體相互作用中靜電、色散、誘導(dǎo)和交換4種相互作用對主、客體相互作用的貢獻，以及影響置換過程的最主要的相互作用方式。整個SAPT計算都在PSI4程序[22]上進行。同時在幾何優(yōu)化后的結(jié)構(gòu)的基礎(chǔ)上，使用Multiwfn程序[23]計算波函數(shù)信息。通過RDG分析(也稱為non-covalent interaction，NCI分析)[24]進行弱相互作用的可視化研究。采用IGM[25]對弱相互作用的強度進行分析。通過靜電勢分析獲得靜電作用在CO2置換CH4水合物過程中的影響方式。

2 結(jié)果與討論

2.1 主、客體雙分子聚體結(jié)構(gòu)優(yōu)化和相互作用能

水合物的晶體結(jié)構(gòu)主要由主體分子和客體分子組成的籠子構(gòu)建而成，其中主體分子與客體分子間的相互作用決定籠子的穩(wěn)定性。在CO2置換CH4水合物過程中，水合物籠子的主體H2O分子、客體CH4分子和置換的CO2分子之間的相互作用影響著整個置換過程。本文采用不同的泛函在6- 311+G(2d,2p)基組下進行結(jié)構(gòu)優(yōu)化，然后采用CCSD(T)/aug- CC- pVTZ水平計算相互作用能。模擬優(yōu)化得到的H2O- H2O、H2O- CH4、H2O- CO2雙分子聚體的幾何構(gòu)型如圖1所示，雙分子聚體相互作用能計算結(jié)果如表1所示。從表1中可以看出，經(jīng)B3LYP- D3泛函優(yōu)化后計算得到的雙分子聚體的能量相比于M062X- D3、WB97XD、B97D這3個泛函的計算結(jié)果，更接近于基準(zhǔn)MP2函數(shù)的計算結(jié)果，表明B3LYP- D3泛函更適合用于描述H2O- H2O、H2O- CH4、H2O- CO2之間的相互作用。由于基準(zhǔn)MP2函數(shù)計算時間較長[26]，因此在B3LYP- D3/6- 311+G(2d,2p)水平下進行CO2置換CH4水合物過程的幾何結(jié)構(gòu)優(yōu)化。

圖1 雙分子聚體的優(yōu)化結(jié)構(gòu)Fig.1 Optimized structures of dimers

表1 H2O- H2O、H2O- CH4、H2O- CO2雙分子聚體相互作用能Table 1 Interaction energies of H2O- H2O, H2O- CH4,H2O- CO2 dimers

在B3LYP- D3/6- 311+G(2d,2p)水平下進行的3種雙分子聚體結(jié)構(gòu)優(yōu)化的基礎(chǔ)上，分別利用B3LYP、B3LYP- D3、M062X、M062X- D3、WB97XD、B97D、MP2等方法在6- 311+G(2d,2p)基組下進行雙分子聚體相互作用能的計算，并與CCSD(T)/aug- CC- pVTZ水平下的計算結(jié)果進行對比，如表2所示。從表2可以看出，相比于其他計算水平，M062X/6- 311+G(2d,2p)水平下計算獲得的分子間相互作用能與CCSD(T)/aug- CC- pVTZ水平下的計算結(jié)果最接近，表明其計算精度較其他方法的要高。CCSD(T)/aug- CC- pVTZ水平主要用于計算小分子單體，對于大分子體系計算較為困難[27]，

表2 不同計算水平下的雙分子聚體相互作用能Table 2 Interaction energies of dimers for different calculation levels

因此，M062X/6- 311+G(2d,2p)水平是CO2置換CH4水合物過程的QM研究的較優(yōu)計算結(jié)合能參數(shù)。

從圖1、表1和表2可以看出，在H2O- H2O聚體中，不同H2O分子中的H原子與O原子由氫鍵連接形成相互作用，相互作用能為-21.05 kJ/mol。文獻[28]對氫鍵進行了系統(tǒng)分類，指出弱氫鍵相互作用能小于-10.46 kJ/mol時，相互作用由靜電和色散作用主導(dǎo)；“weak-to-medium”氫鍵相互作用能在-10.46～-58.60 kJ/mol時，相互作用由靜電作用主導(dǎo)；中氫鍵相互作用能在-46.04～-62.79 kJ/mol之間時，相互作用由靜電作用主導(dǎo)；強氫鍵相互作用能大于-62.79 kJ/mol時，相互作用由靜電和感應(yīng)作用主導(dǎo)。因此H2O分子中H原子與相鄰的H2O分子中O原子間的相互作用介于弱氫鍵相互作用與中氫鍵相互作用之間。在H2O- CH4聚體中，H2O分子中的H與CH4中的C形成微弱的吸引作用，相互作用能為-4.06 kJ/mol。在H2O- CO2聚體中，H2O中的O原子與CO2中的C原子相互作用，相互作用能為-12.18 kJ/mol。可以看出H2O分子之間的相互作用明顯高于H2O- CH4和H2O- CO2，因此在水合物中H2O分子與客體分子之間相互作用偏弱，更傾向于與其他H2O分子相互作用形成籠子。并且與H2O- CH4相比，H2O- CO2具有更高的相互作用能，當(dāng)形成籠型水合物時，則更利于籠子的聚集，得到的CO2水合物較CH4水合物更加穩(wěn)定[4]。

2.2 主、客體雙分子聚體相互作用中各項相互作用的貢獻

在化學(xué)體系中，弱相互作用為分子間或分子內(nèi)的非化學(xué)鍵作用，包括范德華作用、pi- pi堆積作用、氫鍵、二氫鍵、鹵鍵等。SAPT利用片段間的微擾將弱相互作用分解為靜電能(Eelst)、交換能(Eexch)、色散能(Edis)和誘導(dǎo)能(Eind)4部分，從而從能量角度更深入地理解弱相互作用的本質(zhì)。本文通過SAPT對CO2置換CH4水合物過程中主、客體雙分子聚體的能量進行分析，解析主、客體相互作用中靜電、色散、誘導(dǎo)和交換4種相互作用的貢獻。H2O- H2O、H2O- CH4、H2O- CO2雙分子聚體在不同SAPT計算水平下的相互作用能列于表3。將其與CCSD/aug- CC- pVTZ水平下的相互作用能比較發(fā)現(xiàn)，在SAPT下，計算級別的提升并不一定會伴隨計算精度的提升，并且加入δMP2項來考慮誘導(dǎo)和色散間的高階耦合作用，加入(CCD)項以優(yōu)化色散項并不會使得相互作用能的計算值更接近CCSD/aug- CC- pVTZ水平下的計算值。通過對比分析發(fā)現(xiàn)SAPT2+/aug- CC- pVDZ、SAPT2+/jun- CC- pVTZ和SAPT2/aug- CC- pVTZ 水平下計算的相互作用能更接近CCSD/aug- CC- pVTZ水平下的計算結(jié)果，具有較高的計算精度，因此它們在分析H2O- H2O、H2O- CH4、H2O- CO2雙分子聚體相互作用時更具優(yōu)勢。在Parker等[21]有關(guān)于效率計算的研究中指出SAPT2+/jun- CC- pVTZ水平耗時是SAPT2+/aug- CC- pVDZ和SAPT2/aug- CC- pVTZ水平的5倍以上。當(dāng)進行小分子聚體計算時，耗時差異較小，當(dāng)進行水合物籠體的相關(guān)研究時，耗時差異則較大。而在實際的水合物籠體計算中，SAPT2/aug- CC- pVTZ水平計算過程的臨時文件占用幾T的硬盤儲存空間，存在因硬盤空間不足而使得計算終止的問題。因此，考慮到計算精度與計算時長，水合物在SAPT分析中較適合采用SAPT2+/aug- CC- pVDZ水平計算，同時SAPT2+/aug- CC- pVDZ水平是文獻[21]提出的SAPT計算中的“銀”標(biāo)準(zhǔn)，能夠進行精準(zhǔn)的能量分解研究。

表3 不同SAPT計算水平下雙分子聚體相互作用能Table 3 Interaction energies of dimers for different SAPT calculation levels

在SAPT方法中，弱相互作用能可以被分解為具有不同物理意義的組分，即ΔE=Eelst+Eexch+Eind+Edis。其中，Eelst主要反映主客體之間的經(jīng)典靜電相互作用；Eexch表示在主客體中由電子的費米子行為引起的單體波函數(shù)以及反對稱的重疊交換排斥作用；Eind描繪了誘導(dǎo)作用，包括極化以及兩種單體之間的電荷轉(zhuǎn)移；Edis是由主客體的電子之間的庫侖相關(guān)性引起的色散貢獻。分子占據(jù)籠子時，主、客體相互作用能的差異明顯，為了闡明水合物中主、客體相互作用的電子性質(zhì)并探究影響主客體間穩(wěn)定性的最主要的相互作用形式，我們對H2O- H2O、H2O- CH4、H2O- CO2雙分子聚體在SAPT2+/aug- CC- pVDZ計算水平下進行SAPT分析，結(jié)果見表4。從表4中可以發(fā)現(xiàn)相互作用能的負值是由于靜電、誘導(dǎo)、色散3項的總和大于交換排斥項的貢獻產(chǎn)生的，因此總的相互作用主要為吸引作用。為了考察各吸引作用對總吸引作用的貢獻，更直觀深入地理解主客體相互作用的性質(zhì)，將表4中的數(shù)據(jù)轉(zhuǎn)化為每個吸引項的貢獻分數(shù)C，如式(2)所示。

表4 SAPT2+/aug- CC- pVDZ水平下計算的雙分子聚體相互作用能物理分量Table 4 The physical components of the interaction energy of the dimers calculated at the SAPT2+/aug- CC- pVDZ level

(2)

其中EX表示Eelst、Eind或Edis。SAPT衍生的各相互作用對總吸引作用的貢獻分數(shù)如圖2所示。從圖中可以看出，在H2O- H2O聚體中靜電作用占比為64.25%，即吸引作用主要是由靜電作用引起的，而誘導(dǎo)和色散作用的影響較小，占比分別為19.53%和16.21%。在H2O- CH4聚體中靜電作用占比為36.00%，低于色散作用的45.20%，即吸引作用主要是由色散作用引起的，而誘導(dǎo)作用占比僅為18.80%，對總的吸引作用影響較小。在H2O- CO2聚體中色散作用占比為29.62%，相比于H2O- CH4聚體中高的色散貢獻下降15.58%，靜電作用占比為59.38%，較H2O- CH4聚體增加23.38%。由此發(fā)現(xiàn)在雙分子聚體中，隨著相互作用能的增加，靜電作用有上升趨勢，色散作用則表現(xiàn)為下降，即靜電作用占比與相互作用能大小成正比，色散作用占比與相互作用能大小成反比。因此，可以推斷客體分子CO2和CH4與主體H2O分子作用時靜電和色散相互作用占比的差異導(dǎo)致了兩種客體分子形成水合物時結(jié)構(gòu)和能量的差別。

圖2 SAPT衍生的各相互作用對總吸引作用的貢獻分數(shù)Fig.2 The contribution percentage of each interaction derived from SAPT to the total attraction

綜上分析可知，在3種聚體的能量分解中H2O- H2O聚體的吸引作用主要由靜電作用貢獻，色散和誘導(dǎo)作用的影響較小。在H2O- CH4、H2O- CO2聚體的對比中發(fā)現(xiàn)，與H2O相互作用的分子由CH4轉(zhuǎn)變?yōu)镃O2時色散作用的影響減弱，靜電作用的影響加強，表明了靜電作用在CO2置換CH4中的重要地位。由此推斷，增加靜電作用是提升置換效率的有效方法，例如加入與H2O分子有較高靜電作用的試劑分子，可以提高置換效率。

2.3 主、客體分子的靜電勢分析

靜電勢V(r)(ESP)是由原子核和電子在周圍空間每一點r處產(chǎn)生的分子間或分子內(nèi)真實的物理特性，其值能夠通過衍射實驗以及計算的方法確定，計算公式為

(3)

式中，ZA表示位于RA處原子核A所帶的電荷數(shù)，ρ(r)表示分子的電子密度，r′表示距離r單位電子電荷增量距離。在任意點r上的V(r)反映了分子中每個核和每個電子在該點的貢獻[29]。在任何給定的區(qū)域，V(r)是正還是負取決于原子核的正電荷或電子的負電荷在特定位置的主導(dǎo)性。V(r)是正的區(qū)域為親核性，V(r)是負的區(qū)域為親電性。計算過程中當(dāng)電子密度為0.001 e/Bohr3(0.006 784 e/?3)時，在分子表面繪制靜電勢等值面，如圖3所示。V(r)的最正值和最負值分別用Vmax和Vmin表示。圖3(a)為主體H2O分子的靜電勢等值面圖，從圖中可以看出H2O分子中H原子靜電勢的表面電位為正，且Vmax值為189.12 kJ/mol和189.20 kJ/mol，O原子靜電勢的表面電位為負,Vmin值為-140.44 kJ/mol。由此可以得出H2O- H2O聚體中H原子與O原子之間的靜電勢差達到329.64 kJ/mol。圖3(b)為客體CH4分子的靜電勢等值面圖，可以看出H原子具有正的表面電荷，Vmax值為35.37 kJ/mol，Vmin點位于CH4的兩個H原子之間。由圖3(a)和圖3(b)得出，H2O分子的Vmax與CH4分子的Vmin之差以及CH4分子的Vmax與H2O分子的Vmin之差分別為197.45 kJ/mol和175.81 kJ/mol，因此H2O分子的Vmax位點更利于同CH4分子的Vmin位點相互作用，即H2O分子中的H與CH4分子中的C形成相互作用，結(jié)果與圖1一致。圖3(c)為客體CO2分子的靜電勢等值面圖，由圖發(fā)現(xiàn)碳原子表面點位為正，Vmax值為118.08 kJ/mol，氧原子表面靜電勢為負，Vmin值為-48.14 kJ/mol。由圖3(a)和圖3(c)得出，在H2O- CO2聚體中H2O分子的Vmax與CO2分子的Vmin之差以及CO2分子的Vmax與H2O分子的Vmin之差分別為237.26 kJ/mol和258.52 kJ/mol，因此H2O分子的O原子更利于與CO2分子的C原子形成相互作用，這與圖1中的幾何結(jié)構(gòu)優(yōu)化結(jié)果相同。在H2O- CO2聚體與H2O- CH4聚體的對比中，H2O- CO2聚體中兩種分子的靜電勢差值(258.52 kJ/mol)高于H2O- CH4聚體的相應(yīng)值(197.45 kJ/mol)，但低于H2O- H2O聚體的329.64 kJ/mol，因此H2O與CO2相互作用的強度高于H2O與CH4分子，而低于H2O分子與H2O分子間的相互作用，這與相互作用能計算結(jié)果(表1)的趨勢相匹配。圖4為 H2O- H2O、H2O- CH4、H2O- CO2雙分子聚體的靜電勢表面穿透圖，穿透位置與Vmax、Vmin的最大差值一致，對于H2O- H2O、H2O- CH4、H2O- CO2這3種聚體，其分子間靜電勢表面相互穿透距離分別為1.16 ?、0.45 ?、0.98 ?，相互作用強度與相互穿透距離的大小具有一致性，進一步驗證了靜電作用在H2O- H2O、H2O- CH4、H2O- CO2相互作用中的關(guān)鍵影響。

圖3 靜電勢等值面圖Fig.3 Isosurface diagrams of electrostatic potential

圖4 靜電勢表面穿透圖Fig.4 Surface penetration diagrams of electrostatic potential

2.4 RDG及IGM分析相互作用可視化

RDG分析基于電子密度及其降低的梯度，可定性鑒定和表征各種強度的弱相互作用。IGM理論則基于電子密度，能夠識別和量化由于相互作用而導(dǎo)致的電子密度凈梯度衰減，定義分子間相互作用的區(qū)域，從而直接表征分子間的相互作用。本文通過計算波函數(shù)信息，對主、客體分子進行RDG和IGM分析，定向研究主、客體分子對水合物籠子穩(wěn)定性的影響機理。RDG分析是通過生成雙分子聚體真實空間相應(yīng)區(qū)域的梯度等值面來定位它們，通過這些區(qū)域的sign(λ2)ρ值來區(qū)分相互作用的種類，其中sign(λ2)函數(shù)是電子密度Hessian矩陣的第二大本征值λ2的符號，用以反映相互作用類型；弱相互作用臨界點的ρ的數(shù)值和鍵的強度具有正相關(guān)性。sign(λ2)ρ值從0.05→0→-0.05逐漸減小，相互作用變化為強互斥作用(位阻作用)→范德華作用→強吸引作用(氫鍵、鹵鍵等)，在RDG圖中表現(xiàn)為顏色從紅色→綠色→藍色的變化，如圖5所示。當(dāng)sign(λ2)ρ>0時表現(xiàn)為互斥作用，sign(λ2)ρ<0時表現(xiàn)為吸引作用。從圖5中可以看出H2O- H2O聚體中H2O分子之間形成了較強的吸引作用，即氫鍵相互作用，而H2O- CH4、H2O- CO2聚體中分子間呈現(xiàn)出范德華相互作用。H2O- CO2聚體的范德華作用形成的sin(λ2)ρ峰值比H2O- CH4聚體的峰值更負，因此H2O分子與CO2分子之間的范德華作用高于H2O分子與CH4分子之間的范德華作用，這與客體分子和水合物籠子之間的相互作用趨勢一致[30]。

圖5 雙分子聚體的RDG圖Fig.5 RDG diagrams of the dimers

圖6為雙分子聚體間的IGM圖，圖中的sign(λ2)ρ=0附近有左右兩個峰，當(dāng)sign (λ2)ρ<0時，峰高表示吸引作用，當(dāng)sign (λ2)ρ>0時，峰高則表示排斥作用。從圖6可以看出，吸引作用的峰尖比排斥作用的峰尖高，表明雙分子聚體主要表現(xiàn)為吸引作用。由圖6(a)看出吸引作用的峰尖主要為藍色，即H2O- H2O聚體是以氫鍵相互作用為主。由圖6(b)和圖6(c)看出吸引作用的峰尖主要為綠色，即H2O- CH4、H2O- CO2聚體中分子間呈現(xiàn)范德華力相互作用，該結(jié)果與RDG分析一致。在分子間相互作用區(qū)域，IGM方法定義的δg函數(shù)值與相互作用強度具有正相關(guān)性。對比3種雙分子二聚體相互作用的δg值(圖6)發(fā)現(xiàn)δg, H2O- H2O>δg,H2O- CO2>δg, H2O- CH4。這與表2中相互作用能的大小具有一致性，進一步驗證了3種雙分子二聚體相互作用間的差異。

圖6 雙分子聚體的IGM圖Fig.6 IGM diagrams of the dimers

由上可知，RDG分析可表征主、客體分子之間的相互作用，IGM分析可獲得主、客體分子之間的相互作用區(qū)域，判斷相互作用大小，從而實現(xiàn)主、客體相互作用的定性及定量研究。本文研究從電子密度分析角度證實了H2O- H2O聚體是以氫鍵相互作用為主，H2O- CH4、H2O- CO2聚體中分子間是以范德華力相互作用為主。由此，主、客體雙分子聚體分子間的相互作用是影響水合物籠子穩(wěn)定性的關(guān)鍵因素。

3 結(jié)論

(1)對比不同DFT計算方法的精度，得到B3LYP- D3/6- 311+G(2d,2p)水平是CO2置換CH4水合物過程的QM研究的較優(yōu)結(jié)構(gòu)優(yōu)化參數(shù)，M062X/6- 311+G(2d,2p)水平是CO2置換CH4水合物過程的QM研究的較優(yōu)計算結(jié)合能參數(shù)。

(2)結(jié)構(gòu)優(yōu)化與靜電勢分析結(jié)果表明，在CO2置換CH4水合物過程中，H2O分子的H原子與CH4分子的C原子形成相互作用，H2O分子的O原子與CO2分子的C原子形成相互作用；主體H2O分子與客體分子相互作用的強度低于H2O- H2O相互作用，因此主體H2O分子傾向于相互作用形成籠子；H2O- CO2相互作用的強度高于H2O- CH4相互作用，因此，當(dāng)形成水合物籠子時，CO2水合物較CH4水合物更加穩(wěn)定。

(3)由靜電、色散、誘導(dǎo)和交換4種相互作用的貢獻分析得出，CO2置換CH4水合物過程中主、客體分子間的作用主要由靜電作用貢獻，色散和誘導(dǎo)作用的影響較小。置換過程中，客體分子由CH4轉(zhuǎn)變?yōu)镃O2時色散作用的影響減弱，靜電作用的影響加強。

(4)RDG分析和IGM分析從電子密度角度證實了CO2置換CH4水合物過程中主體H2O 分子之間以氫鍵相互作用為主，主、客體分子之間以范德華力相互作用為主。

由此，SPAT計算方法能較好地分析CO2置換CH4水合物過程中主、客體分子間的作用，研究結(jié)果表明靜電作用是置換的關(guān)鍵因素，提高與H2O的靜電作用是提升置換效率的有效方法。下一步我們將采用此QM方法用于天然氣水合物分解機理的微觀理論研究，與MD技術(shù)和微觀實驗結(jié)合分析CO2置換CH4水合物的動態(tài)過程，為天然氣水合物技術(shù)的發(fā)展提供理論指導(dǎo)。