苗志文，甘曉松，許團委，晁 侃

(中國航天科技集團公司四院四十一所，西安 710025)

2016-10-26；

2016-11-15。

苗志文(1991—)，男，碩士，研究方向為固體火箭發(fā)動機內(nèi)流場分析。E-mailmu_qing44@163.com

固體火箭發(fā)動機噴管擴張段粒子沖刷流場分析

苗志文，甘曉松，許團委，晁 侃

(中國航天科技集團公司四院四十一所，西安 710025)

某翼柱形藥柱固體火箭發(fā)動機噴管擴張段出口部位在試驗后出現(xiàn)了與藥柱翼槽位置相對應(yīng)的沖刷痕跡，為了研究Al2O3粒子對噴管擴張段的沖刷規(guī)律，對噴管型面改進提供依據(jù)，對比了不同湍流模型、顆粒軌道模型對形成沖刷痕跡的影響，分析了發(fā)動機噴管擴張段兩相流場特征，確定了形成沖刷痕跡的粒徑范圍，判斷了沖刷痕跡的形成時間，提出了噴管型面改進方案。結(jié)果表明，噴管擴張段的沖刷痕跡形成于發(fā)動機工作的15 s時刻之前，主要由藥柱后翼燃燒產(chǎn)物中顆粒粒徑分布為10～16 μm區(qū)間的粒子造成，改進后的噴管型面可有效降低粒子對噴管擴張段的沖刷。

固體火箭發(fā)動機；翼柱形藥柱；粒子沖刷；噴管擴張段

0 引言

目前，針對Al2O3粒子在固體火箭發(fā)動機噴管中的運動規(guī)律已有較多研究，在不考慮粒子湍流擴散效應(yīng)的條件下，粒子在噴管擴張段存在極限流線，在極限流線和噴管壁面之間的區(qū)域為“無粒子區(qū)”[1-4]，“無粒子區(qū)”大小與粒子直徑呈正相關(guān)性[5-8]。李江[9]基于RTR技術(shù)研究了固體發(fā)動機噴管凝相粒子分布規(guī)律，根據(jù)試驗觀測結(jié)果，噴管擴張段整個截面上都有粒子分布，不存在真正的“無粒子區(qū)”，粒子會在湍流擴散的作用下進入“無粒子區(qū)”而形成“稀薄粒子區(qū)”。劉靜[10]使用顆粒隨機軌道模型對噴管內(nèi)的兩相流動進行了數(shù)值分析，統(tǒng)計結(jié)果表明，粒子的隨機擴散效應(yīng)使其在噴管中的分布更加均勻，且分布范圍更大，驗證了“稀薄粒子區(qū)”的說法。工程實踐中，“稀薄粒子區(qū)”的粒子對噴管擴張段造成沖刷的現(xiàn)象鮮有發(fā)生，公開文獻較少。近期，某發(fā)動機試驗后，噴管擴張段絕熱層出口部位出現(xiàn)了與燃燒室藥柱翼槽數(shù)量和周向位置相對應(yīng)的沖刷區(qū)域，沖刷痕跡清晰、規(guī)整。該現(xiàn)象表明，噴管擴張段尾部出口內(nèi)型面受凝相粒子流沖刷，這對發(fā)動機會產(chǎn)生兩方面影響。其一，擴張段局部絕熱層燒蝕加劇，有結(jié)構(gòu)失效致發(fā)動機工作異常的風(fēng)險；其二，噴管效率降低，致發(fā)動機性能不滿足要求。因此，有必要對噴管擴張段的粒子沖刷現(xiàn)象進行流場分析，提出噴管設(shè)計優(yōu)化的措施。

本文圍繞試驗現(xiàn)象，研究Al2O3粒子對噴管擴張段的沖刷規(guī)律，對比不同湍流模型、顆粒軌道模型對形成沖刷痕跡的影響，分析發(fā)動機噴管擴張段兩相流場特征，確定形成沖刷痕跡的粒徑范圍，判斷沖刷痕跡的形成時間，以為后續(xù)發(fā)動機噴管型面改進提供依據(jù)。

1 模型與仿真方法

1.1 物理模型

發(fā)動機整機長度約15 m，燃燒室為兩分段式，直徑2 m，采用前后翼柱形藥柱，翼槽數(shù)為12，采用潛入式噴管。試驗中噴管擴張段出現(xiàn)的沖刷痕跡與燃燒室藥柱翼槽的周向位置相對應(yīng)。選擇發(fā)動機工作0、15、28 s時刻為代表進行分析，這3個時刻的內(nèi)流場幾何模型均有明顯藥柱翼槽。根據(jù)對稱性選取內(nèi)流場幾何模型的1/4作為計算域，0 s時刻的發(fā)動機內(nèi)流場幾何模型如圖1所示。計算參數(shù)如表1所示，與實際計算的設(shè)置一致。表1中，氣相和凝相的入口質(zhì)量流量均為發(fā)動機實際工作總流量的1/4。假設(shè)所有的凝相粒子都是Al2O3粒子，粒子的初始溫度與氣相入口溫度相同，粒子直徑服從Rosin-Rammler分布。

工作時刻/s氣相質(zhì)量流量/(kg/s)凝相質(zhì)量流量/(kg/s)環(huán)境壓強/MPa083.0734.491592.2738.310.0952889.1937.03入口溫度/K環(huán)境溫度/K摩爾質(zhì)量/(g/mol)比定壓熱容/[J/(kg·K)]323529818.73169

1.2 網(wǎng)格劃分

發(fā)動機工作0 s時刻的計算域采用全結(jié)構(gòu)網(wǎng)格，網(wǎng)格單元數(shù)為135萬，對前后翼及噴管區(qū)域的網(wǎng)格進行加密。發(fā)動機工作15 s和28 s時刻的計算域采用混合網(wǎng)格劃分，在前翼和后翼劃分非結(jié)構(gòu)網(wǎng)格，其他區(qū)域劃分結(jié)構(gòu)網(wǎng)格，網(wǎng)格單元數(shù)分別為136萬和154萬。

1.3 計算方法

發(fā)動機穩(wěn)態(tài)氣相流場采用基于壓力的求解器，不考慮重力的影響，湍流模型選用Spalart-Allmaras模型和SSTk-ω模型，壓力-速度耦合方式采用Coupled算法。藥柱燃面設(shè)置為質(zhì)量流量入口邊界，噴管出口設(shè)置為壓力出口邊界，燃燒室和噴管壁面設(shè)置為無滑移壁面。

凝相粒子的運動采用顆粒確定軌道模型和顆粒隨機軌道模型計算，認為粒子從入口邊界均勻地釋放到計算域中，不考慮粒子的燃燒、蒸發(fā)、破碎過程，不考慮粒子間相互作用，粒子碰撞到壁面后發(fā)生彈性反彈。

1.4 工況設(shè)置

計算工況設(shè)置如表2所示，所有工況粒子的平均直徑均取為30 μm，最大直徑均取為120 μm。工況1～3為了探究選用不同的湍流模型、顆粒軌道模型對形成沖刷痕跡的影響，工況3～11為了確定形成沖刷痕跡的粒徑范圍，工況3、12、13為了判斷沖刷痕跡的形成時間。

表2 工況設(shè)置

2 結(jié)果分析

2.1 湍流因素對形成沖刷痕跡的影響分析

圖2是工況1～3噴管擴張段的粒子濃度云圖，圖3是工況1～3噴管出口的粒子濃度云圖。根據(jù)工況1的計算結(jié)果，噴管壁面附近存在無粒子區(qū)，完全沒有粒子碰撞到噴管擴張段壁面，噴管出口粒子濃度云圖的外緣基本為圓形。根據(jù)工況2的計算結(jié)果，噴管出口粒子濃度云圖的外緣呈現(xiàn)與藥柱翼槽位置相對應(yīng)的波浪形特征，外緣輪廓較模糊，噴管擴張段壁面附近粒子濃度很低。存在少量粒子碰撞到擴張段壁面，碰撞位置與藥柱翼槽的位置相對應(yīng)，碰撞痕跡不均勻，最大粒子濃度為7.2×10-4kg/m3。根據(jù)工況3的計算結(jié)果，噴管出口粒子濃度云圖的外緣呈現(xiàn)與藥柱翼槽位置相對應(yīng)的波浪形特征，外緣輪廓清晰。存在較多的粒子碰撞到擴張段壁面，碰撞位置與藥柱翼槽的位置相對應(yīng)，碰撞痕跡均勻、清晰，最大粒子濃度為1.03×10-3kg/m3，明顯高于工況2噴管擴張段的粒子濃度。

可見，燃氣的湍流運動對凝相粒子的運動軌跡有顯著影響。在工程計算中，Spalart-Allmaras模型[11]具有求解速度快、易收斂的優(yōu)點，主要適用于翼形、壁面邊界層等流動，不適合射流類自由剪切湍流問題。標(biāo)準(zhǔn)k-ω模型[12]中包含了低雷諾數(shù)影響、可壓縮性影響和剪切流擴散，適用于尾跡流動、混合層、射流以及受壁面限制的流動附著邊界層湍流和自由剪切流計算；SSTk-ω模型[13]綜合了k-ω模型在近壁區(qū)計算的優(yōu)點和k-ε模型在遠場計算的優(yōu)點，在廣泛的流動領(lǐng)域中有更高的精度和可信度。固體火箭發(fā)動機噴管擴張段的流動兼有壁面邊界層流動和自由剪切流的特點，相對而言，選用SSTk-ω模型更為合理。模擬凝相粒子的運動，顆粒確定軌道模型不考慮燃氣湍流的影響，削弱了粒子在流場中的擴散，粒子的軌跡相對集中；顆粒隨機軌道模型考慮燃氣湍流的影響，更加接近物理真實，粒子在燃氣脈動速度的影響下軌跡相對較發(fā)散。試驗中粒子遠離噴管軸線碰撞到擴張段壁面上，與粒子在燃氣湍流作用下的徑向運動有關(guān)。因此，選用顆粒隨機軌道模型更具有可信度。

表3統(tǒng)計了工況3條件下碰撞到噴管擴張段壁面的粒子來源，83.91%的粒子來源于藥柱后翼，僅有5.23%的粒子來源于藥柱前翼，10.86%的粒子來源于藥柱其他位置。所以，試驗中出現(xiàn)的沖刷痕跡主要由藥柱后翼產(chǎn)生的粒子造成。

表3 沖刷痕跡的粒子來源

2.2 噴管兩相流場特征分析

對比工況1～3的計算結(jié)果，工況3精度較高，更具有可信度。因此，根據(jù)工況3的計算結(jié)果，對噴管的兩相流場特征進行分析。

圖4是純氣相流場的燃氣速度云圖。

噴管擴張段除壁面邊界層處有速度突變，其他區(qū)域燃氣速度沿徑向分布較為均勻。

圖5是計算凝相后的燃氣速度云圖，受到粒子阻力的影響，噴管擴張段的燃氣速度整體小于純氣相燃氣的速度，與文獻[14]觀點一致。噴管出口位置燃氣速度分布不均勻，結(jié)合圖3(c)可發(fā)現(xiàn)與粒子濃度有關(guān)，粒子濃度低的位置燃氣速度較大，粒子濃度高的位置燃氣速度較小。噴管出口壁面附近粒子較稀薄，濃度小于10-2kg/m3，燃氣受到粒子的阻力較小，燃氣速度在2600～2980 m/s之間；噴管中心部位粒子聚集明顯，粒子濃度最大可達8.66×10-2kg/m3，燃氣受到粒子的阻力較大，燃氣速度在2100～2600 m/s之間。

圖6是噴管擴張段以粒子速度標(biāo)識的粒子軌跡圖。從圖6可觀察到，噴管中心部位粒子軌跡密集，粒子速度主要在1500～2400 m/s之間，噴管擴張段壁面附近粒子軌跡稀疏，粒子速度相對較高，在2400～2830 m/s之間。粒子軌跡外緣呈現(xiàn)出與藥柱翼槽位置相對應(yīng)的波浪形特征。圖7是以粒子直徑標(biāo)識的粒子軌跡圖，噴管擴張段壁面附近的粒子直徑小于噴管中心部位的粒子直徑。圖6和圖7的計算結(jié)果表明，在發(fā)動機噴管擴張段，小直徑粒子主要集中在靠近壁面的位置，粒子速度較大，大直徑粒子主要集中在噴管中部，粒子速度相對較小。

粒子在噴管擴張段的速度、直徑分布特征可從粒子的運動方程進行解釋。粒子在氣相中運動時，假設(shè)只受到氣體的曳力，則粒子的運動方程可寫為

(1)

方程左邊是粒子運動的加速度，方程右邊是單位質(zhì)量粒子受到氣相的曳力。

FD定義為

(2)

式中μ為氣體動力粘度；ρp為粒子密度；dp為粒子直徑；CD為曳力系數(shù)；Re為相對雷諾數(shù)。

Re定義為

(3)

假設(shè)粒子在運動時始終保持球形，則曳力系數(shù)CD[15]：

(4)

為得到粒子運動加速度和粒子直徑的關(guān)系，將上述公式整理為

(5)

由式(5)可看到，粒子的直徑越小，其運動加速度越大。

由于直徑較小的粒子具有更大的運動加速度，其速度的大小和方向更容易改變，所以小直徑粒子比大直徑粒子速度更大，如圖6所示；小直徑粒子更容易運動到遠離噴管軸線的位置，即粒子軌跡更發(fā)散，如圖7所示。

2.3 粒徑范圍對形成沖刷痕跡的影響分析

根據(jù)圖7所示的粒子直徑在噴管擴張段的分布規(guī)律，靠近噴管壁面的粒子直徑最小，可判斷試驗中的沖刷痕跡由小直徑粒子造成。為了進一步確定形成沖刷痕跡的粒徑范圍，在工況3的基礎(chǔ)上調(diào)整最小粒徑，對比噴管擴張段上粒子濃度的變化。圖8是工況3～11噴管擴張段的粒子濃度云圖，圖8(a)～(i)的最小粒徑分別為4、6、8、10、12、14、16、18、20 μm。

前文已指出，粒子直徑越小，粒子運動的加速度越大，粒子在噴管擴張段的軌跡更發(fā)散。所以，圖8中隨著最小粒徑的增大，粒子碰撞到噴管擴張段的范圍逐漸變小，粒子濃度逐漸減小。當(dāng)最小粒徑設(shè)置為4、6、8 μm時，粒子在擴張段的沖刷范圍相對試驗結(jié)果明顯偏大，故認為發(fā)動機實際工作時，在噴管擴張段直徑小于10 μm的粒子很少，可忽略其對壁面的沖刷。當(dāng)最小粒徑設(shè)置為10、12、14、16 μm時，粒子在噴管擴張段周向的沖刷位置與藥柱翼槽的位置相對應(yīng)，沖刷范圍與試驗結(jié)果基本一致，這4種工況下噴管擴張段壁面上最大的粒子濃度分別為1.03×10-3、5.5×10-4、1.1×10-4、3.2×10-5kg/m3，最小粒子直徑從10 μm增加至16 μm，噴管擴張段的最大粒子濃度減小了兩個數(shù)量級。當(dāng)最小粒徑設(shè)置為18、20 μm時，粒子的沖刷痕跡很模糊，沖刷位置無法與翼槽位置相對應(yīng)，擴張段邊緣處最大的粒子濃度在10-6kg/m3的量級，可認為粒徑為18、20 μm的粒子碰撞到噴管擴張段的幾率很小，可忽略不計。結(jié)合圖7的粒子直徑分布規(guī)律，當(dāng)粒子直徑大于20 μm時，其運動軌跡會更加遠離擴張段壁面，粒子碰撞到噴管擴張段壁面的幾率更小，可不予考慮。圖9對比了粒子沖刷核心區(qū)沿流動方向的某條特征線上工況3～9的粒子濃度分布，當(dāng)最小粒徑在4～10 μm之間粒子的沖刷范圍明顯偏大，當(dāng)最小粒徑在10～16 μm之間粒子的沖刷范圍與試驗結(jié)果吻合較好。綜上所述，試驗中出現(xiàn)的沖刷痕跡主要由直徑為10～16 μm的粒子形成。

2.4 沖刷痕跡形成時刻分析

發(fā)動機工作0 s時刻的翼槽最狹窄，粒子在翼槽處的聚集最為顯著，由前文分析0 s時刻可形成試驗中出現(xiàn)的粒子沖刷痕跡，但在該發(fā)動機實際工作過程中，藥柱翼槽存在的時間超過30 s，翼槽存在期間粒子是否都可碰撞到噴管擴張段？工況12、13即對發(fā)動機工作的15、28 s進行分析。

圖10是發(fā)動機工作15 s時刻噴管擴張段的粒子濃度云圖。從圖10可觀察到，存在少量粒子碰撞到噴管擴張段壁面，碰撞位置與翼槽位置相對應(yīng)，但粒子濃度僅在10-6kg/m3量級，可忽略不計。從噴管出口的粒子濃度云圖上已基本觀察不到與翼槽對應(yīng)的波浪形外緣，計算域?qū)ΨQ面上粒子濃度分布異常，出現(xiàn)非物理解，原因是計算時設(shè)置粒子碰撞到對稱面上發(fā)生彈性反彈。這一現(xiàn)象也說明少量粒子具有不可忽略的周向速度，周向速度的來源可能與部分翼槽燃面法向方向偏離發(fā)動機軸線有關(guān)。在發(fā)動機工作28 s，完全沒有粒子碰撞到噴管擴張段壁面。這樣的計算結(jié)果表明，出現(xiàn)試驗中的沖刷痕跡需要粒子在翼槽附近的聚集程度足夠大，該發(fā)動機的粒子沖刷痕跡出現(xiàn)于發(fā)動機工作的15 s時刻之前。

2.5 噴管擴張段型面改進

試驗中粒子碰撞到噴管擴張段形成沖刷痕跡，與藥柱翼槽有關(guān)，與粒徑分布有關(guān)，也與噴管擴張段型面有關(guān)。從喉部到出口，噴管擴張段與軸線的夾角逐漸減小。粒子由于其自身慣性，在噴管擴張段運動過程中速度方向變化幅度小于噴管擴張段與軸線夾角變化幅度，因此可能會碰撞到擴張段壁面。如果發(fā)動機工作在高過載環(huán)境下，粒子碰撞到噴管擴張段壁面的幾率會大大增加，此時粒子的速度在2000 m/s以上，粒子的破壞力不容忽視。因此，建議減小噴管擴張段型面的曲率變化，增大噴管擴張段出口半角。圖11是原噴管和改進后噴管的二維幾何型面。其中，噴管收斂段的幾何型面保持不變，噴管喉部直徑和出口直徑保持不變，僅減小噴管擴張段的曲率變化，增大噴管出口半角。

針對改進的噴管型面，建立了0 s時刻內(nèi)流場幾何模型，采用與工況3相同網(wǎng)格類型、相同的數(shù)值方法開展兩相流場計算。圖12是原噴管和改進噴管在粒子沖刷核心區(qū)沿流動方向的粒子濃度分布，可看到仍存在少量粒子會碰撞到噴管擴張段壁面，但碰撞區(qū)域的粒子濃度已顯著減小，最大的粒子濃度僅為1.53×10-4kg/m3。因此，減小噴管擴張段的曲率變化、增大噴管出口半角可有效降低粒子對擴張段壁面的沖刷。

3 結(jié)論

(1)對于計算凝相粒子的運動軌跡，SSTk-ω湍流模型比Spalart-Allmaras湍流模型精度更高，顆粒隨機軌道模型比顆粒確定軌道模型精度更高。

(2)在發(fā)動機噴管擴張段，壁面附近的燃氣速度大，粒子直徑小，粒子速度大；中心位置的燃氣速度小，粒子直徑大，粒子速度小。

(3)噴管擴張段的沖刷痕跡形成于發(fā)動機工作的15 s時刻之前，主要由藥柱后翼燃燒產(chǎn)物中顆粒粒徑分布為10～16 μm區(qū)間的粒子造成。

(4)減小噴管擴張段型面的曲率變化、增大噴管擴張段出口半角可有效降低粒子對擴張段壁面的沖刷程度。

[1] 劉君,黃琳,夏智勛.潛入式擺動噴管兩相內(nèi)流場數(shù)值模擬[J].固體火箭技術(shù),2003,26(4):18-21.

[2] 吳限德,張斌,陳衛(wèi)東,等.固體火箭發(fā)動機噴管內(nèi)氣粒兩相流動的CFD-DSMC模擬[J].固體火箭技術(shù),2011,34(6):707-710.

[3] 李耿,周為民,張鋼錘,等.翼柱型藥柱對發(fā)動機內(nèi)流場及噴管收斂段燒蝕影響分析[J].固體火箭技術(shù),2010,33(3):260-264.

[4] 陳福振,強洪夫,高巍然,等.固體火箭發(fā)動機內(nèi)氣粒兩相流動的SPH-FVM耦合方法數(shù)值模擬[J].推進技術(shù),2015,36(2):175-185.

[5] 于勇,劉淑艷,張世軍,等.固體火箭發(fā)動機噴管氣固兩相流動的數(shù)值模擬[J].航空動力學(xué)報,2009,24(4):931-937.

[6] 賀征,郜冶.星孔型裝藥發(fā)動機三維兩相流場的數(shù)值模擬[J].推進技術(shù),2004,25(2):118-121.

[7] 蘇鵬輝,艾邦成,潘宏祿.噴管兩相流耦合數(shù)值模擬與粒子熱增量預(yù)測[J].北京航空航天大學(xué)學(xué)報,2012,38(7):900-909.

[8] 李雅娣,陳林泉,蹇澤群.粒子尺寸分級的噴管兩相流場計算[J].固體火箭技術(shù),2003,26(3):32-34.

[9] 李江,蔡體敏,何國強,等.SRM噴管凝相粒子分布規(guī)律的實驗研究[J].推進技術(shù),1991,20(6):84-86.

[10] 劉靜,徐旭.隨機軌道模型在噴管兩相流計算中的應(yīng)用[J].固體火箭技術(shù),2006,29(5):333-353.

[11] Spalart P R, Allmaras S R. A one-equation turbulence model for aerodynamic flow[R]. AIAA 92-1439.

[12] Wilcox D C. Turbulence modeling for CFD[M]. La Canada: DCW Industries,1998.

[13] Menter F R. Two-equation eddy-viscosity turbulence models for engineering applications[J]. AIAA Journal,1994,32(8):1598-1605.

[14] 馮海云,李猛.火箭發(fā)動機噴管結(jié)構(gòu)優(yōu)化設(shè)計研究[J].計算機仿真,2013,30(10):85-102.

[15] Liu A B, Mather D, Reitz R D. Modeling the effects of drop drag and breakup on fuel sprays[R]. SAE Technical Paper 930072.

ParticleerosionflowfieldanalysisonnozzledivergentsectionofSRM

MIAO Zhi-wen, GAN Xiao-song, XU Tuan-wei, CHAO Kan

(The 41st Institute of the Fourth Academy of CASC, Xi'an 710025, China)

Particle erosion marks related to grain fin slots were observed on nozzle divergent section of a finocyl grain SRM after experiments. The main purpose of this research is to investigate the law of particle erosion at the nozzle divergent section and improve the nozzle profile. Effects of turbulence models and particle tracking models on particle erosion were compared. The features of Two-phase flow at the nozzle divergent section were analyzed. The diameter of particles that are relevant to erosion marks was estimated. The time interval of particle erosion was determined. Furthermore, the improved solution of nozzle profile was proposed. The results shows that the particle erosion marks on nozzle divergent section appear after 15 s' working process of the SRM. They were mainly caused by the particles among aft-fin combustion products ranging from 10 μm to 16 μm. The improved nozzle profile could decrease the particle erosion effectively.

solid rocket motor；finocyl grain；particle erosion；nozzle divergent section

V435

A

1006-2793(2017)05-0562-07

10.7673/j.issn.1006-2793.2017.05.005

(編輯：薛永利)