賀立超

(北京航空航天大學 北京 100191)

高超聲速飛行器的發(fā)展，對于未來國家安全、政治和經(jīng)濟利益有著巨大的戰(zhàn)略意義和潛在價值[1]。在高速飛行情況下，激波層和邊界層內(nèi)的氣體溫度急劇升高，離解后的O、N原子在表面熱防護材料的催化效應作用下發(fā)生復合反應并釋放熱量，因而不同的壁面催化特性會極大地影響飛行器表面熱流[2]。

采用數(shù)值計算方法預測高超聲速飛行器表面熱環(huán)境時，由于材料表面的催化過程較復雜，多數(shù)情況下假設材料表面處催化邊界條件為完全催化條件或者完全非催化條件[3]。數(shù)值計算表明：在飛行器再入過程中，完全催化作用的材料表面熱流約為完全非催化表面熱流的3 倍。準確預測催化反應速率，可以為高超聲速飛行器表面熱防護材料的設計提供依據(jù)。隨著高性能計算系統(tǒng)的飛速發(fā)展，反應分子動力學模擬方法已經(jīng)廣泛應用于不同材料的物理性質(zhì)的研究中，成為了研究表面催化反應過程的新方式[4-5]。

SiO2等硅基復合材料具有強度大、高溫穩(wěn)定性以及化學穩(wěn)定性好的特點，被廣泛應用在高超聲速飛行器的熱防護中[6-7]。雖然在先前的研究中，已經(jīng)有研究者通過RMD研究提高了對SiO2表面與O原子碰撞、吸附、催化等過程的認識，但是由于表面附近的環(huán)境條件相互耦合，O 原子與SiO2表面的反應過程的影響因素仍需要進一步詳細討論。因而，本文對SiO2表面溫度、入射O原子動能這兩種單一條件對α-SiO2表面與O原子催化重組反應機制的影響進行研究，以此揭示單一因素對壁面的影響。

1 氣固界面反應模型構建

本模擬使用了ReaxFFSiO反應力場，該力場計算結果與DFT 吻合得很好，能較好地再現(xiàn)多種SiO2表面結構，可以應用于大規(guī)模Si-O界面的分子動力學模擬。

為了研究O 原子與SiO2表面反應過程，首先建立了一個45 ?×45 ?×240 ?的模擬框，以距離模擬框x、y 方向的底面30 ? 高度的平面為中心，建立了一塊沿著（001）表面切割的α-SiO2平板，裸露在平板表面的原子為Si原子。α-SiO2沿z軸方向的尺寸為12 ?，沿x軸和y軸方向的尺寸為45 ?×45 ?，共使用了五層SiO2。為了避免模擬過程中SiO2平板受到O原子的沖擊向下移動或變形，距離平板底面0 ?～1 ? 范圍內(nèi)的原子被固定住，距離平板底面1 ?～10 ? 范圍內(nèi)的原子使用Langevin恒溫器控制溫度，Langevin恒溫器的阻尼常數(shù)為10 fs，距離平板底面10 ?～12 ? 范圍內(nèi)的原子自由移動。

為了保證O 原子持續(xù)與SiO2表面碰撞，使用了通量邊界條件。為避免在初始階段SiO2表面與O原子相互作用，在SiO2表面上方15 ?的平面上，每隔1 ps隨機位置產(chǎn)生一個O原子。設立距離表面15 ?以內(nèi)的區(qū)域為反應區(qū)域，高于表面15 ? 的區(qū)域設立為統(tǒng)計區(qū)域，如圖1所示。O原子z軸方向速度指向SiO2表面，以保證所有O 原子與表面相撞。每隔1 ps，對位于統(tǒng)計區(qū)域的原子統(tǒng)計種類和數(shù)目后刪除。一方面，可以避免與SiO2表面碰撞后的原子在模型內(nèi)大量堆積，在長時間的模擬中形成高密度氣態(tài)區(qū)域，改變表面附近環(huán)境條件；另一方面，統(tǒng)計區(qū)域內(nèi)原子密度的疏解也可以減少與SiO2表面反應后的產(chǎn)物在統(tǒng)計區(qū)域再次反應，進一步精確種類統(tǒng)計結果。

圖1 通量邊界模型的結構

使用兩種環(huán)境條件來分析不同條件對O 原子與SiO2表面反應過程：SiO2表面溫度的影響和入射O原子初始動能的影響。模擬過程對系統(tǒng)整體使用NVE 系綜，在模擬開始前，對SiO2平板進行能量最小化，并使用Nose/Hoover 恒溫器對SiO2平板進行50 ps 的充分弛豫，仿真模擬時長均為1 ns，時間步長為0.25 fs。在對SiO2表面溫度影響的研究中，O 原子入射方向垂直于SiO2表面，攜帶動能為0.05 eV。通過改變Langevin 恒溫器的溫度來獲取不同的表面溫度，溫度改變范圍500～2 100 K。在對O原子入射時攜帶的動能影響的研究中，O 原子入射方向垂直于SiO2表面，SiO2表面通過Langevin 恒溫器穩(wěn)定為1 000 K。O原子入射時的動能變化范圍從0.01～4 eV。

2 結果與分析

2.1 表面溫度的影響

通過對SiO2表面進行溫度控制，筆者對SiO2表面溫度的影響進行了研究。如圖2所示，經(jīng)過50 ps的弛豫后，系統(tǒng)溫度隨著模擬的運行保持穩(wěn)定；而系統(tǒng)勢能在平板弛豫過程中，不同溫度下的勢能均得到了釋放。隨著系統(tǒng)模擬的進一步運行，不同溫度下的勢能都有所升高，但是整體升高幅度較小。這是由于原子和表面碰撞后，碰撞、吸附和復合產(chǎn)生的熱能經(jīng)過熱傳導，通過Langevin 恒溫器重新進行了調(diào)整，而離開表面的碰撞后的產(chǎn)物，進入統(tǒng)計區(qū)域后被刪除，因而本模型很好地維持了系統(tǒng)溫度和勢能的穩(wěn)定平衡。

圖2 系統(tǒng)溫度和勢能隨時間的變化

如圖3所示，在模擬初始階段，O原子大量吸附在SiO2表面，隨著模擬過程繼續(xù)運行，較低溫度下O原子在SiO2表面吸附數(shù)目逐漸飽和，而較高溫度下O 原子仍在繼續(xù)附著。一方面，由于較高的表面溫度為吸附、復合等反應過程提供了足夠的活化能；另一方面，則是隨著溫度升高，α-SiO2表面逐漸呈現(xiàn)無定形SiO2的形態(tài)，提供了更多的反應位點供O原子吸附。

圖3 不同溫度下的表面O原子吸附特性

在O 原子附著數(shù)目逐漸穩(wěn)定后，更多與SiO2表面碰撞后的產(chǎn)物進入了統(tǒng)計區(qū)，進入統(tǒng)計區(qū)的O 原子數(shù)目初期呈現(xiàn)快速增長的趨勢。由于進入統(tǒng)計區(qū)域的O原子統(tǒng)計后被刪除，因而大幅度減少了氣相內(nèi)的反應，因而O 原子數(shù)目不會因為氣相內(nèi)的消耗而減少，O2分子增長速率也相對穩(wěn)定，由于O3分子鍵長更長，鍵能較弱，因而比O2分子不穩(wěn)定，生成數(shù)目較少。通量邊界模型刪除統(tǒng)計區(qū)域的原子，減少氣相內(nèi)反應產(chǎn)物的大量堆積，同樣可以減少O3、O4以及其他不穩(wěn)定的中間分子產(chǎn)物出現(xiàn)在統(tǒng)計區(qū)域內(nèi)的可能性。由圖4（a）可知，在相同時刻，隨著SiO2表面溫度的升高，產(chǎn)生的O2數(shù)目隨之增長。在獲取O2生成數(shù)目后，通過下式獲取了O2在SiO2表面的催化重組系數(shù)。

圖4 不同溫度下催化重組反應特性

不同時刻溫度下的催化重組系數(shù)如圖4（b）所示，隨著表面附著O原子數(shù)、O2生成數(shù)的穩(wěn)定，催化重組系數(shù)隨著模擬的進行逐步達到平衡。催化重組系數(shù)隨溫度升高而升高。

Eley Rideal（E-R）重組機制和Langmuir Hinshelwood（L-H）重組機制是O原子的重組過程占據(jù)主導地位的兩種機制。如圖5（a）所示，≡Si-O·位點和擴散而來的離解原子相撞，形成≡Si-O2，之后Si與O之間的化學鍵斷裂，最終結合成O2分子，這是E-R重組過程?！許i-O·位點得到了相鄰的≡Si-O2的末端O原子，重新形成了一個≡Si-O2，之后Si與O之間的化學鍵斷裂，并通過解吸附，O2分子離開表面，這是L-H 重組過程。通過自編程序的方式，通過對反應過程的軌跡文件進行后處理，對形成O2分子的兩個O原子的運動過程進行軌跡追溯，并對不同反應機制下產(chǎn)生的分子分別進行統(tǒng)計，最終獲取了不同溫度下通過E-R 重機制或L-H 重組機制催化復合形成的O2的比例。由于受到O原子連續(xù)不斷的轟擊，E-R重組機制始終占據(jù)著較大比例，在500 K時，通過E-R重組機制重組形成的O2所占比例高達83%，隨著溫度升高，O原子可以得到更多的活化能，L-H重組機制占據(jù)的比例隨之有所上升，通過E-R 重組機制重組形成的O2接近70%，具體詳見圖5（b）。

圖5 不同溫度下重組機制

2.2 O原子入射動能的影響

筆者通過改變O原子在入射時攜帶的動能（0.01～4 eV），以此來研究來流O 原子的動能對碰撞、催化重組過程的影響。隨著O 原子入射動能的提高，可以為O 原子在表面位點吸附或重組提供的活化能隨之增加，因而同一時刻下O 原子在SiO2表面吸附數(shù)目也隨之迅速升高，具體詳見圖6。O 原子入射動能0.05 eV的表面附著數(shù)目遠遠小于入射動能1 eV的表面附著數(shù)目。當入射動能超過0.5 eV 后，活化能逐漸不再是限制表面附著數(shù)目的關鍵因素，SiO2表面可以提供的位點逐漸趨于飽和，吸附數(shù)目不再繼續(xù)增長。

圖6 不同入射動能下的表面O原子吸附特性

如圖7（a）所示，在初始階段，入射原子大量附著在SiO2表面，隨著表面位點的減少，O原子吸附數(shù)目增長開始放緩。隨著表面分布的O 原子越來越密集，相鄰位點的O 原子會更容易復合，大量新入射的O 原子和附著在表面的O 原子碰撞的概率也會增加。如圖7（b）所示，隨著反應的進行，O2分子開始迅速生成。當入射動能較高時，有較多O3進入統(tǒng)計區(qū)域，這是由于入射原子攜帶的大量動能使相對不穩(wěn)定的O3更容易脫離SiO2表面的吸附進入氣相中。

圖7 O原子入射動能1 eV時SiO2表面附著數(shù)和產(chǎn)物變化

由圖8（a）所示，當入射動能小于1 eV 時，隨著O原子入射動能增加，重組生成的O2也隨之增加，這是由于隨著入射動能增加，SiO2表面O 原子附著更加密集，可以用于O 原子表面催化重組的活化能也進一步增加，O原子與位點碰撞概率與催化重組概率均提升。由圖8（b）所示，可見當入射O 原子動能小于1 eV 時，隨著O 原子入射動能的增加，催化重組系數(shù)也隨之升高。當入射O 原子動能大于1 eV 時，隨著O 原子入射動能增加，O2生成數(shù)、催化重組系數(shù)都隨之下降。這主要是由于較高的動能會使入射的O原子彈離表面或使吸附在SiO2表面的O 原子直接解吸附，而不是繼續(xù)發(fā)生催化重組反應。因此，O原子與表面位點碰撞時，過高的動能可以降低E-R催化重組反應效率。

圖8 不同O原子入射動能影響下的催化重組特性

3 結論

本文使用ReaxFF反應力場，通過建立α-SiO2表面與O 原子碰撞的通量邊界模型，從SiO2表面溫度和入射O 原子初始動能兩種條件出發(fā)，對O 原子碰撞SiO2表面的反應過程不同因素對催化重組的影響進行了研究，得到的主要結論如下。

（1）在模擬初期，入射的O原子會大量吸附在SiO2表面，這一過程隨著SiO2表面位點的飽和逐步放緩。

（2）由于更高的表面溫度可以為催化重組反應提供更充沛的活化能，導致L-H 催化重組機制也隨著表面溫度的升高而活躍，而SiO2表面的催化能力也隨著表面溫度的升高而增大。

（3）當入射O原子攜帶的初始動能低于1 eV時，催化重組系數(shù)隨著動能增加而升高，高于1 eV時，催化重組系數(shù)隨著動能增加而降低，這是由于入射O 原子過高的動能，使得O原子彈離表面或直接解吸附，減少了通過E-R催化機制形成的O2分子。