趙文麗 王永剛 張路路 岳大光 孟慶田
1)(山東農業大學信息科學與工程學院,泰安 271018)
2)(山東交通學院理學院,濟南 250357)
3)(山東師范大學物理與電子科學學院,濟南 250358)
基于一個最新的勢能面,運用切比雪夫波包方法對初始態為反應體系在1.0—2.0 eV 的碰撞能量范圍內進行了動力學研究.通過對角動量量子數J=60以下的所有分波進行計算,得到了反應幾率、積分散射截面和速率常數.計算中用到了耦合態近似方法和考慮科里奧利耦合效應的精確量子方法.通過對比發現,隨著角動量量子數以及能量的增加,科里奧利耦合效應的影響越發顯著,因而對于該反應體系,科里奧利耦合效應不可忽略.本文計算所得的積分散射截面和速率常數尚無實驗數據可以比較,對該反應的后續研究有一定的參考價值.
在燃燒化學[1]和天體化學[2,3]中,碳原子和氫分子間的反應是一種非常重要的化學反應,而CH2自由基作為中間產物,在這類碳氫反應中發揮著重要的作用[4].迄今為止,人們在實驗[5-7]和理論[8-19]上對亞甲基分子CH2進行了大量的研究,且主要集中于等幾個單態的勢能面性質[8-11]和相關的的動力學計算而對于三態的勢能面和相關的反應動力學研究相對較少,特別是正向反應少有涉及.本文的工作就是基于一個最新的CH2(X3A'')勢能面并對反應開展量子波包動力學研究.
在理論方面,Knowles 等[20]在1983年給出了一個態的勢能面,該勢能面能較好地反映體系在所有漸近區的性質.基于該勢能面,Murrell和Dunne[21]用準經典軌線(QCT)的方法研究了并計算了該反應的速率常數.1993年,Harding等[22]應用多參考組態相互作用(MRCI)[23,24]結合Dunning相關一致基組的方法得到了CH2基態的全維勢能面(HGS).隨后,該小組應用HGS計算了反應的速率常數,在高溫區域(1500—2000 K),理論結果和實驗符合得很好.基于 HGS,van Harrevelt等[25]應用量子波包和QCT的方法計算了 C H+H→C+H2的散射截面和速率常數.另外,Gamallo等[26]使用玻恩-奧本海默(BO)近似的量子力學方法對H(2S)+反應進行了動力學計算,并指出該反應只發生在解耦合的基態上.1996年,Guadagnini和Schatz[27]使用 QCT方法研究了 C(3P)+H2→H+CH 反應的產物分布,并指出產物 C H+H 的形成完全取決于C原子靠近H2的插入機制.
實驗上,Scholefield等[28]使用激光燒蝕石墨的方法獲得 C(3P)粒 子束,研究了 C(3P)與H2,HCl,HBr和CH3OH的反應,并指出產物CH的生成與反應物的轉動態呈現明顯依賴關系.1997年,Ehbrecht等[29]使用高能量(10 —180 eV)的C原子束與H2發生光化學反應,觀測到了生成CH產物不同 電 子態輻射 帶(2Δ -X2Π,2Σ--X2Π和并指出三種態可觀測的能量范圍不同.
在2019年,Zhang等[30]應用多參考組態相互作用方法(MRCI),采用包含Davidson修正的aug-cc-pVXZ(X=Q和5)基組計算得到的5285個能量點,將其外推到完備基組極限(CBS)并利用多體展開(MBE)的方法擬合構建了一個關于的全維解析勢能面,其中勢能最大值與最小值之間的跨度為17 eV,勢能面的總方均根偏差為0.0349 eV.基于該勢能面,該團隊運用QCT方法,計算了反應的積分散射截面和速率常數,與基于HGS勢能面的動力學計算結果和實驗結果符合得很好,從而驗證了該勢能面的精確性.
CH2自由基電子態性質復雜,相關的動力學研究涉及到BO近似和非絕熱過程,不過對該體系性質進行全面的研究并不是本文的目的.本工作主要是基于Zhang等[30]構建的勢能面,用切比雪夫量子波包方法對反應的反應幾率和散射截面等性質進行研究.本文所采用的方法,相對于傳統的含時方法有著更優越的標度律,只需要波包的一次傳播就能夠給出所有能量點處的反應幾率,并且波包在實空間的傳播沒有近似,因而在計算上更具有優勢[31-33].本文主要包含以下幾部分: 第二部分介紹所采用的理論方法以及相應的數值計算方法.第三部分給出理論計算結果和討論.第四部分總結得出相關的結論.除特別說明之外,文中的所有公式和計算都使用原子單位.
對于反應物滿足交換對稱的體系(C + H2),使用雅克比坐標(R,r,γ)能夠減少計算基組個數,在雅克比坐標下C + H2體系的哈密頓量的表達式為

式中 r 表示雙原子(H-H)之間的距離,R 表示碳原子與雙原子質心之間的距離(C-H2),μr=mH/2和μR=2mCmH/(mC+2mH)表 示約化質量,V(R,r,γ)為體系的勢能面,和是雙原子轉動角動量算符和軌道角動量算符.軌道角動量算符的平方可以表示為

在切比雪夫波包方法中,波包傳播的三項遞歸關系為[34,35]

其中,Hnorm=(H-H+)/H-是歸一化的哈密頓量,且通過光譜的上(下)限計算獲得[36].在 R(或 r)格點邊緣,用高斯型函數D來描述輸出波的邊界條件

其中 dx決定阻尼范圍,xd為阻尼起始點.初始波包|ψ0〉 表示為

N 是歸一化常數,k0為初動量,R0與δ 分別為初始波包的中心位置和寬度,| φi〉 為一個具體的振轉本征函數.
通常S-散射矩陣元計算耗時較長[37],為了避免這種情況,本文采用流計算法,應用這種方法獲得的初始態總反應幾率為[38]

式中 rf決定了產物通道的分界面,θ≡arccos[(E-H+)/H-]是切比雪夫角,ai(E)是 初始波包的振幅[39]

為了獲得給定初始態的積分散射截面(ICS),應該對所有J的反應幾率進行求和

式中Ec為碰撞能量,是宇稱為p、角動量量子數為 J 時,初始態的反應幾率.將相應的積分散射截面 συiji(Ec)對碰撞能量 Ec結合玻爾茲曼權重進行積分,得到初態的速率常數為

其中f是電子的簡并因子,kB和T分別為玻爾茲曼常數和熱力學溫度.

圖1 CH2等勢線,圖中等勢線間隔為0.1 eV(a)C沿著共線構型靠近H2分子;(b)C 沿著T構型插入H2Fig.1.Equipotential contour plot for CH2,the contour increments are 0.1 eV:(a)For bond stretching in C-H-H linear geometry;(b)for T-shaped insertion of C into H2 diatoms.

圖2 90°和180°的最小能量路徑Fig.2.The minimum energy paths as a function of RCHRHH at 90° and 180°.
本文采用切比雪夫量子波包方法研究C(3P)+H2→H+CH反應.由上述分析可知,該反應是一個有較大的勢壘和較小的勢阱的典型吸熱反應,為了計算得到精確的ICS和速率常數,需要較高的碰撞能量、大量的分波以及較大的基組數,尤其是考慮完全的CC效應時.對于 J=0 的情況,經過反復的收斂性測試,本文選取了最優的計算參數列于表1.應用這些計算參數,對于J=0,5,15,··,60的分波進行了CS 近似和完全的CC效應的計算,其他反應幾率通過上述反應幾率用三次樣條插值法插值得出[42],進而得出碰撞能量范圍在1.0—2.0 eV 的ICS,最后給出速率常數.

表1 波包計算中的數值參量(除特殊說明,均采用原子單位a.u.)Table 1.Parameters used in wave packet calculation(The atomic unit is used in the calculation unless otherwise stated).

圖3 不同分波的反應幾率隨著碰撞能量的變化Fig.3.The reaction probabilities vs.collision energy at different J.
圖3給出了J=0,J=10,J=30和J=50的精確量子計算(考慮完全CC效應)反應幾率隨碰撞能量的變化.從圖中可以看出,整體上反應幾率數值較小(< 0.22),這是因為C原子無論從哪個角度碰撞H2,過程中都有可能形成中間絡合物CH2,當能量較低、J較小時,絡合物還可能比較“長壽”,并且CH2也有很大的幾率衰變回反應物C+H2[27].圖3中不同分波的反應幾率顯示范圍是0—0.22,碰撞能量顯示范圍是1.0—2.0 eV.對于J=0,反應幾率逐漸從0開始增大,約為1.10 eV時,反應幾率大于0.1,這與前面MEP討論的C(3P)+H2→H+CH反應是一個吸熱反應是一致的,與該反應需要吸熱1.109 eV也是一致的.但是同時要注意到,在1.0—1.1 eV的范圍內,碰撞能量不足以越過該體系的勢壘,而反應幾率仍然不為0,這只能用量子隧穿效應解釋.隨著總角動量J的增加,離心勢壘逐漸增大,反應閾能逐漸增加,J=30的閾能約在1.2 eV,J=50的閾能約在1.65 eV.從整體上看,反應幾率隨著碰撞能量的增加而增加,這也符合一般的吸熱反應的反應幾率變化規律.在低能區域反應幾率呈現強烈的振蕩,這是由體系的勢阱(0.2 eV)所產生的共振引起的.當碰撞能量增加時,反應幾率的振蕩減弱,逐漸變得不明顯.
在分子反應動力學的量子方法中,CS近似被廣泛地應用于很多體系,而且這種方法已經被證明是相當成功的[43].但是,CS近似中忽略了科里奧利效應.實際上,有動力學研究表明在一些化學反應中,這種效應會產生非常重要的影響因而不能忽略[44,45].例如 Meijer和Goldfield[44]發現CC效應在H+O2反應中扮演著重要的角色,如果忽略CC效應將會使J > 0時的計算變得不可靠,這種偏差在高能區越發明顯.
圖4畫出了 J=5,10,20的 CS近似和CC兩種方法下的反應幾率.兩種計算結果都表現出共振特性,隨著J的增大,一些共振會消失,這是因為J越大,波函數計算中包含的分波更多.此外,CS的振蕩幅度比CC的大得多,這是反應中科里奧利耦合效應受到長程力支配的結果.整體上來講,在角動量J比較小的情況下(J=5),CC計算的反應幾率小于CS近似的數值;當J=10時,兩種方法的反應幾率基本符合;當J=20時CC計算的反應幾率大于CS近似的結果,J越大,二者差別越明顯.從J較大時的單一分波來看(J=20),在低能區域(小于1.3 eV),CS的計算結果比CC的大,在1.3—2.0 eV的范圍內,CC計算的結果大于CS近似的數值,并且碰撞能量越大偏差越大.理論認為,中間絡合物CH2的形成對于反應機制起著主導作用.當J較小時,CS近似下由于H2分子態的限制,C更容易與其發生碰撞,從而促進產物生成,所以CS的幾率較大.隨著J的增大,更多分波貢獻被考慮其中.一方面,CS近似下長程相互作用影響較大,長程勢會使得碰撞時間變長,原子的重新組合變慢從而使得反應幾率變小;另一方面,由于科里奧利耦合效應激發的振動模式更有利于破壞H-H鍵從而促進產物的形成.這種情況下,CC效應就會顯著影響反應過程,所以CC的幾率較大.

圖4 CC與CS反應幾率比較Fig.4.Comparisons between the CC and CS probability.
根據(8)式,ICS是把所有分波的反應幾率乘以各自權重,然后再對所有的分波求和得到的.在1.0—2.0 eV的碰撞能量范圍下,我們對所有J < 60分波進行加權求和得到了 C(3P)+H2→H+CH 反應的ICS.圖5給出了分別用CS近似和CC方法得到的ICS隨著碰撞能量的變化.從圖中可以看出,隨著碰撞能量的增加,兩種計算方法的ICS都會增加,與一般吸熱反應的ICS的變化趨勢一致.相比于CS近似的ICS,CC計算的結果振蕩得更加劇烈,但幅度比CS的小,并且CC的ICS隨碰撞能量的增加而增加得更快.當EC< 1.27 eV時,兩種方法的計算結果基本符合.當EC> 1.27 eV時,CC的ICS迅速增加,CS的ICS則緩慢增加,碰撞能量越高,二者的差距越大.例如,當EC=1.21 eV 時,σCC=0.135 ?2,σCS=0.130 ?2,二者數值幾乎相當;當EC=1.61 eV時,σCC=0.647 ?2,σCS=0.256 ?2,前者是后者的 2.53 倍;當EC=1.99 eV 時,σCC=1.022 ?2,σCS=0.305 ?2,前者是后者的3.35倍.考察勢能面的特征就會發現,勢能面的勢阱、鞍點和過渡態的位置對于中間絡合物的形成是非常有利的,即不利于產物的形成,這就導致了CS近似下ICS隨能量增加而緩慢攀升,然而,正如前文所述,科里奧利耦合效應就可以克服這一點,所以CC的ICS隨能量增加上升得非常快.

圖5 C+H2反應的積分散射截面Fig.5.The integral cross section of the C+H2 reaction.
基于積分散射截面,應用(9)式,我們計算出了 C(3P)+H2→H+CH 反應的速率常數,如圖6所示.溫度的取值范圍為1000—2000 K,對應ICS碰撞能量的取值范圍1.0—2.0 eV.圖中紅色實線表示CC計算的速率常數,黑色虛線表示CS近似的結果.從圖6中可以看出,在1000—1400 K,速率常數增加比較平緩,CC和CS的結果基本相一致;當T > 1400 K時,隨著溫度的增加,速率常數迅速增加,我們認為速率常數對于溫度如此敏感歸因于該反應為吸熱反應,有較高的反應閾能.另一方面,我們注意到文獻[46]報道的類似電子結構O(3P)與CH4反應的速率常數,該文獻認為計算過程中忽略勢壘交叉的影響會導致速率常數計算結果偏大.基于同樣的原因,我們預測本文速率常數計算結果可能比實際偏大.此外,類似于ICS的結果,隨著溫度的增加,CC速率常數比CS攀升得更快,例如在T=1400 K時,CC計算和CS近似的速率常數分別為 5.52 × 10—15和5.07 ×10—15cm3·molecule—1·s—1,二者的差值為4.50 ×10—16cm3·molecule—1·s—1,CS 比 CC 結果小 8.25%;在T=1700 K時,CC計算和CS近似的速率常數分別為 2.91×10—14和2.37×10—14cm3·molecule—1·s—1,CS比 CC結果小 18.58%;在 T=2000 K時,CC計算和CS近似的速率常數分別為9.43 × 10—14和6.98 × 10—14cm3·molecule—1·s—1,二者的差值為2.46 × 10—14cm3· molecule—1· s—1,CS 比 CC 結果小25.95%.可見,溫度越高,CC與CS的計算結果差別越大.

圖6 C+H2反應的速率常數Fig.6.The rate constant of the C+H2 reaction.