趙文麗 王永剛 張路路 岳大光 孟慶田

1)(山東農(nóng)業(yè)大學信息科學與工程學院,泰安 271018)

2)(山東交通學院理學院,濟南 250357)

3)(山東師范大學物理與電子科學學院,濟南 250358)

基于一個最新的勢能面,運用切比雪夫波包方法對初始態(tài)為反應(yīng)體系在1.0—2.0 eV 的碰撞能量范圍內(nèi)進行了動力學研究.通過對角動量量子數(shù)J=60以下的所有分波進行計算,得到了反應(yīng)幾率、積分散射截面和速率常數(shù).計算中用到了耦合態(tài)近似方法和考慮科里奧利耦合效應(yīng)的精確量子方法.通過對比發(fā)現(xiàn),隨著角動量量子數(shù)以及能量的增加,科里奧利耦合效應(yīng)的影響越發(fā)顯著,因而對于該反應(yīng)體系,科里奧利耦合效應(yīng)不可忽略.本文計算所得的積分散射截面和速率常數(shù)尚無實驗數(shù)據(jù)可以比較,對該反應(yīng)的后續(xù)研究有一定的參考價值.

1 引 言

在燃燒化學[1]和天體化學[2,3]中,碳原子和氫分子間的反應(yīng)是一種非常重要的化學反應(yīng),而CH2自由基作為中間產(chǎn)物,在這類碳氫反應(yīng)中發(fā)揮著重要的作用[4].迄今為止,人們在實驗[5-7]和理論[8-19]上對亞甲基分子CH2進行了大量的研究,且主要集中于等幾個單態(tài)的勢能面性質(zhì)[8-11]和相關(guān)的的動力學計算而對于三態(tài)的勢能面和相關(guān)的反應(yīng)動力學研究相對較少,特別是正向反應(yīng)少有涉及.本文的工作就是基于一個最新的CH2(X3A'')勢能面并對反應(yīng)開展量子波包動力學研究.

在理論方面,Knowles 等[20]在1983年給出了一個態(tài)的勢能面,該勢能面能較好地反映體系在所有漸近區(qū)的性質(zhì).基于該勢能面,Murrell和Dunne[21]用準經(jīng)典軌線(QCT)的方法研究了并計算了該反應(yīng)的速率常數(shù).1993年,Harding等[22]應(yīng)用多參考組態(tài)相互作用(MRCI)[23,24]結(jié)合Dunning相關(guān)一致基組的方法得到了CH2基態(tài)的全維勢能面(HGS).隨后,該小組應(yīng)用HGS計算了反應(yīng)的速率常數(shù),在高溫區(qū)域(1500—2000 K),理論結(jié)果和實驗符合得很好.基于 HGS,van Harrevelt等[25]應(yīng)用量子波包和QCT的方法計算了 C H+H→C+H2的散射截面和速率常數(shù).另外,Gamallo等[26]使用玻恩-奧本海默(BO)近似的量子力學方法對H(2S)+反應(yīng)進行了動力學計算,并指出該反應(yīng)只發(fā)生在解耦合的基態(tài)上.1996年,Guadagnini和Schatz[27]使用 QCT方法研究了 C(3P)+H2→H+CH 反應(yīng)的產(chǎn)物分布,并指出產(chǎn)物 C H+H 的形成完全取決于C原子靠近H2的插入機制.

實驗上,Scholefield等[28]使用激光燒蝕石墨的方法獲得 C(3P)粒 子束,研究了 C(3P)與H2,HCl,HBr和CH3OH的反應(yīng),并指出產(chǎn)物CH的生成與反應(yīng)物的轉(zhuǎn)動態(tài)呈現(xiàn)明顯依賴關(guān)系.1997年,Ehbrecht等[29]使用高能量(10 —180 eV)的C原子束與H2發(fā)生光化學反應(yīng),觀測到了生成CH產(chǎn)物不同 電 子態(tài)輻射 帶(2Δ -X2Π,2Σ--X2Π和并指出三種態(tài)可觀測的能量范圍不同.

在2019年,Zhang等[30]應(yīng)用多參考組態(tài)相互作用方法(MRCI),采用包含Davidson修正的aug-cc-pVXZ(X=Q和5)基組計算得到的5285個能量點,將其外推到完備基組極限(CBS)并利用多體展開(MBE)的方法擬合構(gòu)建了一個關(guān)于的全維解析勢能面,其中勢能最大值與最小值之間的跨度為17 eV,勢能面的總方均根偏差為0.0349 eV.基于該勢能面,該團隊運用QCT方法,計算了反應(yīng)的積分散射截面和速率常數(shù),與基于HGS勢能面的動力學計算結(jié)果和實驗結(jié)果符合得很好,從而驗證了該勢能面的精確性.

CH2自由基電子態(tài)性質(zhì)復雜,相關(guān)的動力學研究涉及到BO近似和非絕熱過程,不過對該體系性質(zhì)進行全面的研究并不是本文的目的.本工作主要是基于Zhang等[30]構(gòu)建的勢能面,用切比雪夫量子波包方法對反應(yīng)的反應(yīng)幾率和散射截面等性質(zhì)進行研究.本文所采用的方法,相對于傳統(tǒng)的含時方法有著更優(yōu)越的標度律,只需要波包的一次傳播就能夠給出所有能量點處的反應(yīng)幾率,并且波包在實空間的傳播沒有近似,因而在計算上更具有優(yōu)勢[31-33].本文主要包含以下幾部分: 第二部分介紹所采用的理論方法以及相應(yīng)的數(shù)值計算方法.第三部分給出理論計算結(jié)果和討論.第四部分總結(jié)得出相關(guān)的結(jié)論.除特別說明之外,文中的所有公式和計算都使用原子單位.

2 理 論

對于反應(yīng)物滿足交換對稱的體系(C + H2),使用雅克比坐標(R,r,γ)能夠減少計算基組個數(shù),在雅克比坐標下C + H2體系的哈密頓量的表達式為

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

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

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

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

N 是歸一化常數(shù),k0為初動量,R0與δ 分別為初始波包的中心位置和寬度,| φi〉 為一個具體的振轉(zhuǎn)本征函數(shù).

通常S-散射矩陣元計算耗時較長[37],為了避免這種情況,本文采用流計算法,應(yīng)用這種方法獲得的初始態(tài)總反應(yīng)幾率為[38]

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

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

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

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

3 結(jié)果和討論

3.1 勢能面和最小能量路徑

圖1 CH2等勢線,圖中等勢線間隔為0.1 eV(a)C沿著共線構(gòu)型靠近H2分子;(b)C 沿著T構(gòu)型插入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°.

3.2 計算細節(jié)

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

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

3.3 反應(yīng)幾率

圖3 不同分波的反應(yīng)幾率隨著碰撞能量的變化Fig.3.The reaction probabilities vs.collision energy at different J.

圖3給出了J=0,J=10,J=30和J=50的精確量子計算(考慮完全CC效應(yīng))反應(yīng)幾率隨碰撞能量的變化.從圖中可以看出,整體上反應(yīng)幾率數(shù)值較小(＜ 0.22),這是因為C原子無論從哪個角度碰撞H2,過程中都有可能形成中間絡(luò)合物CH2,當能量較低、J較小時,絡(luò)合物還可能比較“長壽”,并且CH2也有很大的幾率衰變回反應(yīng)物C+H2[27].圖3中不同分波的反應(yīng)幾率顯示范圍是0—0.22,碰撞能量顯示范圍是1.0—2.0 eV.對于J=0,反應(yīng)幾率逐漸從0開始增大,約為1.10 eV時,反應(yīng)幾率大于0.1,這與前面MEP討論的C(3P)+H2→H+CH反應(yīng)是一個吸熱反應(yīng)是一致的,與該反應(yīng)需要吸熱1.109 eV也是一致的.但是同時要注意到,在1.0—1.1 eV的范圍內(nèi),碰撞能量不足以越過該體系的勢壘,而反應(yīng)幾率仍然不為0,這只能用量子隧穿效應(yīng)解釋.隨著總角動量J的增加,離心勢壘逐漸增大,反應(yīng)閾能逐漸增加,J=30的閾能約在1.2 eV,J=50的閾能約在1.65 eV.從整體上看,反應(yīng)幾率隨著碰撞能量的增加而增加,這也符合一般的吸熱反應(yīng)的反應(yīng)幾率變化規(guī)律.在低能區(qū)域反應(yīng)幾率呈現(xiàn)強烈的振蕩,這是由體系的勢阱(0.2 eV)所產(chǎn)生的共振引起的.當碰撞能量增加時,反應(yīng)幾率的振蕩減弱,逐漸變得不明顯.

在分子反應(yīng)動力學的量子方法中,CS近似被廣泛地應(yīng)用于很多體系,而且這種方法已經(jīng)被證明是相當成功的[43].但是,CS近似中忽略了科里奧利效應(yīng).實際上,有動力學研究表明在一些化學反應(yīng)中,這種效應(yīng)會產(chǎn)生非常重要的影響因而不能忽略[44,45].例如 Meijer和Goldfield[44]發(fā)現(xiàn)CC效應(yīng)在H+O2反應(yīng)中扮演著重要的角色,如果忽略CC效應(yīng)將會使J ＞ 0時的計算變得不可靠,這種偏差在高能區(qū)越發(fā)明顯.

圖4畫出了 J=5,10,20的 CS近似和CC兩種方法下的反應(yīng)幾率.兩種計算結(jié)果都表現(xiàn)出共振特性,隨著J的增大,一些共振會消失,這是因為J越大,波函數(shù)計算中包含的分波更多.此外,CS的振蕩幅度比CC的大得多,這是反應(yīng)中科里奧利耦合效應(yīng)受到長程力支配的結(jié)果.整體上來講,在角動量J比較小的情況下(J=5),CC計算的反應(yīng)幾率小于CS近似的數(shù)值;當J=10時,兩種方法的反應(yīng)幾率基本符合;當J=20時CC計算的反應(yīng)幾率大于CS近似的結(jié)果,J越大,二者差別越明顯.從J較大時的單一分波來看(J=20),在低能區(qū)域(小于1.3 eV),CS的計算結(jié)果比CC的大,在1.3—2.0 eV的范圍內(nèi),CC計算的結(jié)果大于CS近似的數(shù)值,并且碰撞能量越大偏差越大.理論認為,中間絡(luò)合物CH2的形成對于反應(yīng)機制起著主導作用.當J較小時,CS近似下由于H2分子態(tài)的限制,C更容易與其發(fā)生碰撞,從而促進產(chǎn)物生成,所以CS的幾率較大.隨著J的增大,更多分波貢獻被考慮其中.一方面,CS近似下長程相互作用影響較大,長程勢會使得碰撞時間變長,原子的重新組合變慢從而使得反應(yīng)幾率變小;另一方面,由于科里奧利耦合效應(yīng)激發(fā)的振動模式更有利于破壞H-H鍵從而促進產(chǎn)物的形成.這種情況下,CC效應(yīng)就會顯著影響反應(yīng)過程,所以CC的幾率較大.

3.4 積分散射截面和速率常數(shù)

圖4 CC與CS反應(yīng)幾率比較Fig.4.Comparisons between the CC and CS probability.

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

圖5 C+H2反應(yīng)的積分散射截面Fig.5.The integral cross section of the C+H2 reaction.

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

圖6 C+H2反應(yīng)的速率常數(shù)Fig.6.The rate constant of the C+H2 reaction.

4 結(jié) 論