溫雄飛,趙 瑜,馬新建,趙 志,孟 浩
(上海航天動力技術研究所,上海 201109)
為改善固體火箭發動機工作性能,提高固體推進劑能量特性,通常在固體推進劑中加入金屬鋁粉[1],從而使推進劑燃燒產物中含有氣相與凝相顆粒兩種成分[2],在發動機燃燒室及噴管內形成典型的兩相流動。凝相顆粒的存在會對發動機噴管內流場特性產生較大影響[3],同時凝相顆粒流會沖擊噴管壁面,加速燒蝕而改變噴管型面和喉部面積[4],因此研究固體發動機噴管中的兩相流動情況,能夠為噴管設計與熱防護提供指導。國外Hwang等[5]應用拉格朗日方法與麥科馬克二步二階差分格式相結合的辦法來求解噴管中兩相流動問題;Thakre等[6]應用歐拉-拉格朗日方法結合多組分模型研究了HTPB復合推進劑的燃燒產物在石墨噴管中的流動與燒蝕問題;Greenfield[7]等針對超聲速微噴管兩相流問題,研究了兩相流對推力和噴管效率的影響。國內淡林鵬等[8]應用顆粒軌道模型研究了長尾噴管中顆粒與壁面在不同碰撞情況下的粒子運動軌跡。于勇等[9]應用顆粒軌道模型求解顆粒相,結合高精度空間離散格式,對噴管中不同尺寸顆粒運動規律及顆粒相對流場的影響進行了研究;楊育文等[10]對高含鋁固體推進劑低壓下固體火箭發動機噴管尾流場進行了數值模擬和實驗研究;童悅等[11]對擺動噴管進行了三維兩相數值模擬研究,重點分析了擺角對粒子濃度分布與噴管推力系數的影響。斜切噴管主要應用于推力矢量控制及導彈級間分離[12],非常規的幾何外形導致其流場特性要比一般直噴管更加復雜。國內學者近年來逐步開展了對斜切噴管的研究工作:劉君等[13]應用PHOENICS軟件,對有無斜切的雙噴管固體火箭發動機內流場進行了二維數值模擬研究;徐瑋等[14]應用有限體積法對某非軸對稱斜切噴管三維流場特性進行了數值模擬研究,并分析了斜切角度對發動機推力的影響;金賀龍等[15-16]采用非定常可壓縮N-S方程與Realizablek-ε湍流模型相結合的方法,對不同角度斜切噴管的流場特性與推力特性進行了數值模擬研究。
綜上,國內外諸多學者對于固體火箭發動機噴管內兩相流動的研究主要集中于常規直噴管,而針對斜切噴管兩相流場特性的研究工作相對較少[17]。文中基于固體火箭發動機斜切噴管,運用Realizablek-ε湍流模型和顆粒軌道模型,采用單元內顆粒源法進行兩相耦合計算方法,開展了斜切噴管內三維兩相流動數值模擬計算,分析了兩相流場特性及對推力的影響。
文中研究的斜切噴管結構如圖1所示,其主要由基本部分和斜切部分組成,基本結構參數有:基本部分擴張半角α,斜切部分擴張半角β,斜切角δ,噴管入口處直徑L1,喉部直徑L2,噴管擴張段軸線長度L3,噴管出口截面在中心軸線方向的投影距離L4。由于斜切噴管的非對稱結構,其推力由兩部分組成:一部分為沿噴管軸線方向的軸向力Fx,另一部分為垂直于噴管軸線方向的側向力Fy,其合力與噴管軸線的夾角為推力偏轉角φ。選取噴管喉部直徑d=0.06 m作為特征長度,該斜切噴管結構參數見表1。

表1 噴管結構參數表Table 1 Nozzle structure parameters
由于固體火箭發動機工作過程涉及復雜的物理和化學反應,為簡化計算,又能夠真實模擬斜切噴管內的兩相燃氣流動情況,引入如下假設條件:
1)忽略化學反應過程,推進劑燃燒生成的混合氣體滿足理想氣體狀態方程。
2)不考慮氣相、顆粒相以及固體壁面之間的輻射換熱現象。
3)不考慮兩相流對噴管型面的燒蝕和顆粒沉積。
1.3.1 控制方程
氣相控制方程[18]可表示為:
(1)
式中:A=[ρ,ρu,ρv,ρw,e]T,B,C,D為守恒型通量,BV,CV,DV為粘性通量,Sp為固體顆粒對氣相產生的源項。
顆粒軌道求解控制方程[18]為(單位質量下):
(2)
(3)
式中:vp表示顆粒速度;xp表示顆粒位置;u表示氣相的平均速度;u′表示氣相的瞬時湍流脈動速度;g為重力加速度。FD為顆粒的阻力函數:
(4)
式中:CD、ρp和dp分別表示曳力系數、顆粒密度和顆粒直徑。Re表示相對雷諾數,其計算式為:
(5)
式中:ρ表示氣相密度;μ表示動力粘度。
1.3.2 湍流模型
數值模擬選擇Realizablek-ε湍流模型,k方程和ε方程[19]分別為:
(6)
(7)
式中:Gk是由層流速度梯度引起的湍動能k的產生項;Gb是由浮力引起的湍動能k的產生項;YM是由于在可壓縮湍流中,過渡擴散而產生的波動;其他參數為相關特征常量和自定義源項。
1.3.3 發動機推力計算
固體火箭發動機推力為發動機內、外表面全部作用力的合力,其計算公式[3]為:
(8)
式中:Pin是作用于發動機內表面的壓強;Pout是作用于發動機外表面的壓強;dAs是表面微元面積。對于斜切噴管發動機,推力由式(8)對發動機內外表面的壓強積分得出。
1)氣相邊界條件:選取鋁粉含量18%的HTPB復合推進劑,采用質量流量入口條件,給定入口總溫總壓,流動方向垂直于入口邊界,出口采用壓力出口,壁面采用絕熱無滑移邊界條件。
2)顆粒相邊界條件:在入口處,取每個網格邊的中點作為顆粒的加入點。給定入口處顆粒的初始速度、溫度和顆粒質量流量。經大量實驗證實,固體火箭發動機中顆粒相尺寸分布一般服從對數正態分布或者Rosin-Rammler分布,文中取平均粒徑為70 μm的Rosin-Rammler分布,顆粒與壁面碰撞后反彈,切向和法向恢復系數均為0.8,顆粒出口邊界不加任何限制條件,達到出口即逃逸。
3)采用標準壁面函數法作為近壁處理方式;采用第一層網格y+值為35,增長率為1.2,共15層的邊界層處理方式。
4)顆粒相與氣相的耦合計算采用計算單元內顆粒源項算法(PSIC),即先計算氣相壓強、速度、溫度等參數,再計算顆粒軌道以及顆粒各項參數信息,最終將顆粒相代入氣相方程再次求解,將此過程循環迭代,直至收斂。計算過程采用二階迎風格式,收斂精度為10-6。
5)邊界條件數值:氣相質量流率為8.685 kg/s;氣相入口壓強為11.68 MPa;氣相入口溫度為3 556.3 K;氣相出口壓強為101 325 Pa;氣相出口溫度為293 K;顆粒相質量流率為3.722 kg/s;顆粒相入口速度為0.1 m/s;顆粒溫度為3 253 K。
采用上述Realizablek-ε湍流模型、顆粒軌道模型與單元內顆粒源法對固體火箭發動機內兩相流場特性開展研究的方法已被很多學者應用并取得了良好效果[20-21]。
首先選取顆粒相平均粒徑dp為70 μm的工況進行計算,得到如圖2~圖4所示噴管對稱截面的氣相及兩相流場的特征參數分布,其中X、Y為模型位置坐標,單位為m。分析圖2~圖4可知,顆粒相的存在會對斜切噴管的溫度場和速度場產生重要影響,但對壓力場影響較小,壓力場分布基本不變。

圖2 斜切噴管溫度場分布云圖Fig.2 Temperature field distribution contour of angle-cut nozzle

圖3 斜切噴管壓力場分布云圖Fig.3 Pressure field distribution contour of angle-cut nozzle

圖4 斜切噴管速度場分布云圖Fig.4 Velocity field distribution contour of angle-cut nozzle
對比兩相流場與純氣相流場情況下,斜切噴管中心軸線處與壁面附近的流場分布,通過圖5、圖6可以發現,顆粒相對噴管壁面附近的流場結構影響較小,兩相流場速度略小于純氣相流場速度,而溫度基本相同,但顆粒相對噴管軸線附近區域的流場結構影響較大。

圖5 斜切噴管距壁面0.02 m附近流場特征參數對比Fig.5 Comparison of flow field characteristic parameter curves along 0.02 m away from the angle-cut nozzle wall

圖6 斜切噴管軸線附近流場特征參數對比Fig.6 Comparison of flow field characteristic parameter curves along the angle-cut nozzle axis
通過圖7顆粒濃度分布可知顆粒相在軸線附近分布濃度較大,則噴管軸線附近區域受到顆粒阻力作用,燃氣流動速度明顯隨顆粒濃度增大而減小,存在一個低速區;同時由于顆粒相對氣相作用的熱增量,導致噴管軸線附近也存在一個高溫區,且此流場低速區和高溫區都起始于噴管喉部并沿噴管出口方向逐漸擴大。

圖7 dp=70 μm粒徑下顆粒濃度分布云圖Fig.7 Particle concentration contour for 70 μm particle diameter
結合流場顆粒濃度分布,分析該現象原因為:噴管軸線附近顆粒濃度較大,顆粒隨氣體一起流入噴管過程中,顆粒相速度滯后于氣相,阻礙燃氣流動,從而降低燃氣流速。此外,由于顆粒溫度高于氣體溫度,對氣相有部分熱量傳遞,且顆粒之間互相碰撞會將部分動能轉化為熱能,這兩部分熱增量導致該區域流場溫度升高。
圖8為顆粒相最大粒徑為120 μm,最小粒徑為10 μm,平均粒徑分別為70 μm,60 μm,50 μm,40 μm的Rosin-Rammler分布時,流場內顆粒運動軌跡及粒徑分布。

圖8 顆粒軌跡及粒徑分布對比Fig.8 Comparison of particle tracks and distribution for different particle diameters
由圖可知在噴管基本部分擴張段及斜切部分壁面附近存在一個無顆粒區,且隨著平均粒徑的增大,斜切部分無顆粒區域范圍變大。分析圖9可知,顆粒運動速度在經過喉部以后隨粒徑增大而減小。主要原因是顆粒直徑越小,其慣性越小,經過喉部后氣流對顆粒加速效果越顯著,體現出小顆粒較好的隨流性,在斜切部分的運動軌跡分布更加分散,無顆粒區范圍較小;當顆粒直徑變大時,慣性增大導致顆粒運動受到氣流擾動變小,運動軌跡更集中于軸線附近且比較平直,從而使無顆粒區域的范圍變大。

圖9 不同粒徑下噴管軸線顆粒速度分布對比Fig.9 Comparison of particle velocity along the angle-cut nozzle axis for different particle diameters
圖10給出純氣相以及不同粒徑時,沿出口方向斜切噴管軸線處氣相溫度及速度分布的對比曲線。

圖10 不同粒徑下噴管軸線附近流場特征參數對比Fig.10 Comparison of flow field characteristic parameter curves along the angle-cut nozzle axis for different particle diameters
分析曲線圖可知,顆粒相的存在導致噴管中氣相溫度明顯升高,噴管軸線處氣相溫度先減小后增大,即在斜切噴管基本部分擴張段起始處存在一個低溫區,粒徑越大溫度變化越小,而在斜切部分增大到與收斂段相同的溫度。顆粒相的存在使噴管軸線處兩相流場的氣相速度小于純氣相流場速度,且噴管軸線處氣相速度先增大后減小,存在一個高速區,粒徑越大速度變化越小。
分析原因認為,粒徑越大,慣性越大,對氣相阻力作用越大,使氣相速度增大較少,速度變化較小;并且粒徑越大,其速度越小,滯留時間越長,從而與氣相的作用時間越長,兩相之間的動量與熱量交換率增大,使氣相溫度降低較少,溫度變化較小;同時由于顆粒相存在熱慣性,將保持原來溫度,當顆粒運動到噴管斜切部分時,速度減小,與周圍燃氣熱量交換率增加,使氣相溫度升高。根據圖7可知由于通過噴管上下壁面的兩束濃度較大的顆粒流在斜切噴管基本部分擴張段起始區域附近交匯,導致在此區域顆粒相與氣相相互作用增強,從而動量與熱量相互轉化率提升,形成如圖11所示的低溫區與高速區。

圖11 dp=40 μm粒徑下兩相流場特征Fig.11 Characteristics of two-phase flow field for 40 μm particle diameter
定義Fx為噴管軸向推力;Fy為垂直于燃燒室軸線方向作用力,即噴管側向推力;φ為推力偏轉角;Pmax為流場最大壓強;Pmin為流場最小壓強。根據固體火箭發動機原理,發動機的推力等于發動機所有內、外表面作用力的合力,通過對發動機內外表面作用力進行積分計算可得推力F數值,計算結果如表2所示。

表2 不同粒徑下斜切噴管推力Table 2 Trusts of the angle-cut nozzle for different particle diameters
在顆粒相質量分數為30%工況下,通過表2中推力數據可知,顆粒相的存在使得噴管推力增大約1%~2%,且隨著顆粒直徑增大,推力減小,推力偏轉角減小。由圖12可知,顆粒相的加入使燃氣溫度滯后,出口處燃氣溫度升高,出口壓強增大,且隨著粒徑增大溫度和壓強提升幅度減小,計算所得的推力變化趨勢與壓強變化趨勢相同。

圖12 噴管出口中心軸線附近燃氣參數分布對比Fig.12 Comparison of gas parameter curves along the angle-cut nozzle outlet axis for different particle diameters
文中進行斜切噴管內三維兩相流動數值模擬,主要分析噴管內兩相流場特性以及粒徑大小變化對流場與發動機推力的影響,得到以下結論:
1)顆粒相的存在會對發動機斜切噴管的流場結構產生較大的影響,尤其是對噴管軸線附近區域流場結構影響較大,但對斜切噴管邊界處流場結構影響較小。噴管軸線附近存在一個燃氣流動速度較低,溫度較高的區域,且此區域起始于噴管喉部并沿噴管出口方向逐漸擴大范圍。
2)斜切噴管兩相流場中,斜切部分的壁面附近存在無顆粒區域。顆粒直徑越大,運動越集中于軸線附近且軌跡比較平直,無顆粒區域越大。
3)顆粒相直徑的變化會對發動機斜切噴管流場特征參數的分布產生影響,顆粒的加入使得沿噴管軸線方向氣相速度降低,溫度升高。隨著顆粒直徑增大,發動機噴管軸線處氣相速度先增大后減小,燃氣溫度先減小后增大;且顆粒直徑越大,溫度和速度變化越小。
4)兩相流狀態下的噴管推力大于純氣相推力,且隨著顆粒直徑的增加,斜切噴管推力減小,推力偏轉角減小。