亞大氣壓下AP/HTPB微尺度穩態燃燒的數值模擬

趙健鋒，余永剛

(南京理工大學 能源與動力工程學院，江蘇 南京 210094)

引 言

AP/HTPB復合推進劑廣泛應用于火箭發動機和底排裝置。相較于火箭發動機內的高壓環境，底排彈飛行時所處環境壓強一般低于大氣壓強。且在底排彈飛行過程中，隨著飛行高度的提升，環境壓強不斷降低。因此，為保證底排裝置最優的增程效果，對亞大氣壓環境下推進劑的燃燒特性進行分析十分必要。

由于AP/HTPB復合推進劑真實的燃燒過程十分復雜，為更清楚地了解其燃燒過程，研究人員進行了很多實驗與數值研究。其中黃蒙等[1]對HTPB/Al/AP/RDX各組元相互作用進行了實驗分析。更多研究人員采用數值計算，相繼提出了許多理論模型，目前BDP多火焰模型[2]使用較為廣泛。BDP多火焰模型將AP/HTPB的燃燒分為3個階段，分別受三步反應控制。除BDP模型外還有經簡化的BDP兩步反應模型和詳細的多步化學反應模型。

過去很多研究者對AP/HTPB的燃燒特性進行了實驗研究和數值分析。利用兩步反應模型，Cazan等[3]分析AP/HTPB的火焰結構；Knott等[4]研究了貝克萊數對AP/HTPB燃燒的影響；Massa 等[5]、 Gaduparthi等[6]對比了利用兩步和三步化學反應機理對AP/HTPB燃燒進行數值模擬的差異。

上述研究大都在2MPa壓強以上開展。目前關于亞大氣壓下AP/HTPB燃燒的研究多為實驗研究。Schoyer等[7]測量了在20～140kPa下AP/HTPB中加入草酸鹽和硝基胍后推進劑在不同初始溫度下的燃速。Krishnan等[8]在2～100kPa壓強范圍內研究了推進劑中加入氧化鐵和亞鉻酸銅的影響；Lee等[9]測量了含有不同添加成分的5種推進劑的燃速，得到20～100kPa下燃速壓強指數為0.6～0.88，高于高壓下的燃速壓強指數； Miller等[10]在20～100kPa內分別測量了加入催化劑和降速劑后推進劑的燃速，并得到了推進劑低壓熄火極限與實驗藥條直徑的關系；Tanaka等[11]在1～100kPa壓強范圍內測量了雙基和復合推進劑及加入催化劑后的燃速，催化劑并未影響AP/HTPB的低壓極限，但提高了燃速；Kuwahara[12]報道了在 20～100kPa四組 AP/ HTPB和AP/HTPE的燃速和火焰結構，得到壓強指數約為 0.6～0.98； HTPB含量不同，火焰顏色也不同。

上述關于亞大氣壓下 AP/HTPB 燃燒的文獻多關注于燃速的測量，以及火焰厚度及顏色的觀測，而關于亞大氣壓下AP/HTPB燃燒特性未給出詳細的數值分析結果。基于上述背景， 本研究利用BDP三步反應模型，針對亞大氣壓下AP/HTPB穩態燃燒進行數值模擬，同時與高壓環境下(4MPa)的燃燒特性對比，分析壓強對AP/HTPB燃燒特性的影響，為深入研究底排彈飛行過程中推進劑亞大氣壓燃燒狀況對底阻的影響奠定基礎。

1 物理模型

以HTPB黏合劑作為基體，將AP粒子填充進其中，經固化得到AP/HTPB復合推進劑。 Powling等[13]將這一結構簡化成三明治的夾層結構，該物理模型被廣泛用于AP/HTPB燃燒研究[3,5,6]。在微尺度條件下，利用三明治結構，建立了二維AP/HTPB模型，其結構如圖1所示。該模型為周期模型。

圖1 AP/HTPB二維三明治模型

文獻[3,5,6]表明，在2MPa壓強以上，AP氣相分解反應放熱量占AP/HTPB燃燒總放熱量的比例較大，但隨著壓強的下降，其所占總放熱量的比例不斷減少。

文獻[14]表明，純 AP 的 LPL(低壓熄火極限)約為2MPa，在該壓強以下，AP 的放熱已不能維持燃燒。

研究人員[8,11,12]得到 AP/HTPB復合推進劑的LPL低于20kPa，由此推斷在 20kPa～2MPa壓強范圍內，AP/HTPB 燃燒過程中氧化劑與黏結劑之間氣相燃燒放熱量已占總放熱量的很大一部分，且隨壓強的下降，氧化劑與黏結劑之間的燃燒放熱量占總放熱比例增加。

在數值模擬過程中，為簡化模型并根據引言中描述的實驗現象和推斷，本研究采取如下假設：

(1)亞大氣壓下AP和HTPB分解氣體之間的氣相反應優先于純AP氣相反應；

(2)固相分解均在表面進行，忽略里層可能發生的反應，固相物性參數不變；

(3)利用Arrhenius定律描述固相分解和氣相反應；

(4)氣體為理想氣體混合物，Le數、Pr數均為1；

(5)不考慮輻射影響。

2 數學模型

根據上述模型和假設，建立如下控制方程：

2.1 化學動力學方程

氣相的三步化學反應[2]：

其中X由AP顆粒分解產生，一般認為X是氨氣和過氯酸氣體，Y由HTPB基質分解產生的烴類氣體(CH2)，Z為O2+H2O+N2+HCl。β為化學計量比。第一步反應為AP與HTPB之間發生的氣相反應，第二步反應為AP自身的氣相分解反應。

化學反應速率R1、R2、R3遵循Arrhenius定律[5]分別表示如下：

R1=D1Pn1[X]3.3[Y]0.4exp(-E1/RuT)

(1)

R2=D2Pn2[X]exp(-E2/RuT)

(2)

R3=D3Pn3[Y][Z]exp(-E3/RuT)

(3)

式中：D為化學反應速率常數；P為氣體壓強；n為壓強指數；[X]、[Y]、[Z]分別代表組分X、Y、Z的質量分數；E為反應活化能；Ru為理想氣體常數。

2.2 固相方程

固相方程為：

SH+(λcT)=0

(4)

式中：λc為導熱系數；SH為能量源項，代表固相AP和HTPB分解吸熱,為分解熱、燃速、密度與燃面面積之積；Q為反應熱。

(5)

從宏觀上看，燃速是單位時間內推進劑燃燒的厚度，而從微觀上可看作是AP顆粒與HTPB基體分解速率的快慢，因此可用Arrhenius定律定義的燃速表示，A為指前因子。

(6)

2.3 氣相控制方程

連續性方程：

(7)

動量方程：

(8)

(9)

組分方程：

L(X,Y,Z)=(-R1-βR2,-R2-R3,R1-βR3)

(10)

能量方程：

(11)

上述算子L定義為：

(12)

理想氣體方程：

PV=nRuT

(13)

Sm為質量源項，由式(6)得到的燃速與面積和密度相乘得到。SM為動量源項，固體分解時生成氣體的相對初速度，由(20)得到：

(14)

SM=Sm·v

(15)

2.4 初始條件與邊界條件

模型左右兩側邊界為對稱邊界，梯度為零；底部模擬遠場邊界，設為恒溫；出口設為壓強出口，外部環境壓強與計算域內環境壓強一致。

(16)

T|y=-L=300K

(17)

Pout=Pin

(18)

在氣固兩相交界處保持溫度連續、質量守恒，能量平衡即：

T|0+=T|0-

(19)

(20)

(21)

式中：Pout為出口外壓強；Pin為出口處計算域內壓強；ρc為固相密度；ρg為氣相密度；r為燃速；v為生成氣體的相對初速度。

2.5 計算方法及參數設置

本研究使用CFD軟件FLUENT進行數值模擬 ，采用均勻結構化網格，網格尺寸為10-6m。采用二階迎風格式來離散對流項，用SIMPLE算法對壓強進行修正。模型已經過網格無關性驗證。動力學及物性參數設置如表1。認為固相反應動力學參數不發生變化，AP和HTPB的物性及動力學參數使用參考文獻的數據，見表1。

表1 AP/HTPB化學動力學及物性參數[5]

3 數值模擬結果

3.1 模型合理性驗證

Tanaka[11]通過將燃燒器連接真空泵和穩壓罐，使燃燒器內壓強穩定且低于大氣壓，用移動攝像機觀察燃燒情況，測得了AP質量分數為80%的AP/HTPB推進劑的燃速。將模擬得到的燃面溫度通過式(6)計算燃速，并與Tanaka實驗結果對比，結果如圖2所示，在壓強20～80kPa下基本與實驗結果相吻合，驗證了模型的合理性。

圖2 燃速模擬結果與Tanaka實驗結果對比

3.2 火焰位置與形狀

圖3～圖5為不同壓強下三步反應的放熱云圖。

對比圖3～圖5可知，第一步反應放熱區域較大，單位體積放熱率遠大于其他反應，因此亞大氣壓下第一步反應放出了AP/HTPB燃燒過程中絕大部分熱量，在很大程度上控制了AP/HTPB推進劑的燃燒特性。第一步反應放熱貼近燃面，核心放熱區近似橢圓狀。第三步反應放熱區遠離燃面，火焰呈現一維分布。第二步反應核心放熱區處于二者中間區域，近似圓形。

圖3 第一步反應放熱云圖

圖4 第二步反應放熱云圖

圖5 第三步反應放熱云圖

穩態燃燒三步反應放熱位置的相對高低，代表其反應的先后順序和反應速率。在亞大氣壓環境下，三步反應速率都較為緩慢，擴散過程進行得更充分，這是亞大氣壓下反應的特征，也是導致反應放熱區出現高低不同位置的重要原因。由于第一步反應速率相對較快，由AP分解得到的組分X直接參與了與組分Y的反應，這導致第一步反應放熱區域靠近燃面。圖4(b)表明，40kPa下，第二步反應核心放熱區域(紅色部分)中心的高度約為0.05mm，而圖7(c)組分X的質量分布表明，在40kPa下，組分X在0.05mm處的質量分數已不足30%。這表明確實大部分組分X直接參與了第一步反應，只有一小部分組分X參與了分解為組分Z的過程。

未參與第一步反應的組分X向上擴散，同時憑借第一步反應放熱，加速了反應速率，因此在第一步反應放熱區域上方，形成第二步反應放熱區。第二步反應是發生第三步反應前提，兩步反應存在著必然的先后關系，因此導致第三步反應放熱區處于最上方。

對比不同壓強下三步反應的放熱速率可以看出，隨著壓強的降低，放熱率逐步降低，第一、二步的放熱核心逐步上移。在80kPa下第一、二步反應放熱核心在y方向的中心坐標約為0.02mm和0.04mm，在20kPa時分別上升至0.03mm和0.06mm。這是由于壓強下降，反應速率變緩，擴散作用更加明顯使放熱區上移。由于擴散作用增強，將導致組分混合更加充分，使反應放熱區域面積增加。第三步反應放熱區面積隨著壓強的下降有減小的趨勢。前面分析到第一步反應提供了絕大部分熱量，隨著壓強的下降，反應速率變緩，整個計算域溫度降低，從而導致核心放熱面積減小。

圖6為高壓4MPa下的三步反應放熱云圖。由于前兩步反應的放熱面積較小，因此將最大縱坐標調整為0.1mm。相比于圖3～圖5中亞大氣壓下的放熱面積，高壓下放熱面積明顯縮小。第二步反應緊貼燃面，成為細帶狀，布滿AP顆粒表面；第一、三步反應分布在AP和HTPB交界上方，放熱中心化為兩個點。從圖6中可以看到，在4MPa的壓強下第二步反應的單位體積放熱率已經超過第一步反應，且在4MPa壓強下，第二步反應的放熱面積遠大于第一步反應，因此第二步反應釋放出了AP/HTPB燃燒過程中的絕大部分熱量。這表明在高壓下第二步反應對燃燒做出了主要貢獻，因此高低壓下控制AP/HTPB推進劑燃燒的放熱反應并不相同。

圖6 高壓4MPa下三步反應放熱云圖

圖7為4MPa和40kPa下組分X、Y的分布。造成AP/HTPB高低壓下燃燒現象不同的原因可以解釋如下：AP顆粒分解產生的組分X，通過前兩步反應參與燃燒。第二步反應對環境壓強更加敏感[2]，相較于第一步反應，壓強的變化使其反應速率變化更快。因此在亞大氣壓下，第一步反應放熱占據了主導。隨壓強升高，第二步反應以較大的“加速度”提高反應速率，以至超過第一步反應的速率。同時由于組分X同時參與兩種反應，因此存在著兩種反應對組分X的“爭奪”，在亞大氣壓下，當組分X與Y發生反應，勢必導致參與分解反應的組分X的濃度降低，加速減緩了組分X分解速率。同時由于化學反應速率的較慢，擴散效應增強也影響著放熱區域的形成，在擴散作用下，整個放熱區域較大。高壓下，化學反應速率較擴散速度快。如圖7(a)，只有十分貼近燃面處才有組分X的分布，這說明大部分組分X在與組分Y混合前已自行分解，導致放熱面積減小，最終導致在亞大氣壓下呈現擴散與混合過程共存的燃燒火焰，在高壓下為擴散火焰。

圖7 組分X、Y含量分布云圖

3.3 溫度分布

圖8為亞大氣壓下AP/HTPB上方氣相溫度分布情況。

圖8 亞大氣壓下AP/HTPB氣相溫度分布云圖

由圖8可見，氣相溫度分布特性基本一致，不隨壓強改變而變化，最高溫度隨著壓強的降低而下降。HTPB上方的溫度略低，AP顆粒上方的溫度較高，靠近燃面的等溫線呈現“凸”字形。這是由于第一、二步反應放熱均處于AP顆粒上方，在HTPB上方放熱較少，因此，HTPB上方氣相溫度較低，而AP顆粒上方氣相溫度較高。離燃面較遠處，在同一高度，溫度逐漸趨于一致。

圖9為高壓4MPa時氣相溫度分布云圖，圖中顯示AP上方的氣相溫度較低，而HTPB上方的氣相溫度較高。這是由于AP上方由純AP分解放熱，分解熱較低。而HTPB上方是黏結劑和氧化劑之間的氣相反應，反應熱較高。同時由于組分X不斷分解成組分Z，導致形成組分Z的“氣流柱”，帶走了從HTPB上方高溫區橫向導入AP上方低溫區的熱量。最終形成如圖所示的溫度分布。

圖9 高壓4MPa下AP/HTPB氣相溫度分布云圖

圖10為40kPa和4MPa環境壓強下燃面處固體的溫度分布情況。當AP質量分數為80%時，亞大氣壓下，HTPB與AP溫度較高，而交界面處溫度較低。從氣相圖中可看到，靠近燃面處的氣相溫度仍然是HTPB較低而AP較高，但由于單位質量AP從預熱到完全分解前需要吸收的熱量比HTPB要多，因此導致AP固相溫度稍低于HTPB。由式(6)計算可得交界面處燃速較慢，此處相對于整個燃面突出。4MPa時則完全相反，交界面處溫度較高，計算得到燃速較快，因而交界面將會凹陷，見圖11。

圖10 40kPa和4MPa下AP/HTPB固相溫度分布

圖11 40kPa和4MPa下AP/HTPB燃面形狀

4 結 論

(1)亞大氣壓下,AP/HTPB復合推進劑穩態燃燒,第一步反應靠近燃面，且該步反應放熱量在AP/HTPB推進劑燃燒過程中占主導地位，而其余反應的放熱量對燃燒影響較小。

(2)環境壓強影響了AP/HTPB復合推進劑火焰特性、反應位置。亞大氣壓下，反應放熱區較大，擴散過程與混合過程共存。4MPa下放熱區較小，以擴散火焰為主。

(3)放熱區域與放熱率影響AP/HTPB推進劑表面氣固兩相溫度分布，從而造成高低壓下不同的燃面形狀。亞大氣壓下AP與HTPB交界處相對于整個燃面突出，而高壓下交界處相對于整個燃面凹陷。