張晨曦,林三春,李 易,唐 碩
(1.西北工業(yè)大學(xué) 航天學(xué)院,西安 710072; 2.陜西省空天飛行器設(shè)計重點實驗室(西北工業(yè)大學(xué)),西安 710072; 3.北京宇航系統(tǒng)工程研究所,北京 100076)
發(fā)展可重復(fù)使用的航天運輸系統(tǒng)對于滿足未來空間開發(fā)和降低發(fā)射成本等需求具有重要意義[1-2],而垂直回收則是目前運載火箭部分重復(fù)使用的主要方式。為了對再入過程中動壓及熱流峰值密度進行約束,一子級發(fā)動機在分離后于合適高度啟動,向飛行方向前方噴射氣體使火箭減速[3]。動力減速過程中,方向相反的發(fā)動機射流與來流形成反向射流流動,回流至飛行器附近后會顯著改變運載火箭的環(huán)境特征,因而對反向射流回流區(qū)流動特征以及動力減速過程中火箭氣動特性的研究就顯得尤為重要。
Klevanski等[4]以及Zilker[5]針對DLR擬開發(fā)的可回收火箭開展了力、熱環(huán)境數(shù)值模擬研究,研究表明箭體力、熱環(huán)境主要受反向射流的流動特征影響,而流動特征又由高度欠膨脹的反向射流決定。Rincón[6]對飛行速度較低情況下(包含亞音速)的動力減速過程進行了數(shù)值模擬研究,結(jié)果表明在低速飛行中來流弓形激波被推移至較遠的距離,故箭體表面熱流密度與阻力系數(shù)更小。Ecker等[7]開展了獵鷹9運載火箭反推階段CFD與結(jié)構(gòu)有限元耦合研究,研究表明反向射流產(chǎn)生了低密度氣體環(huán)境,大幅降低了箭體熱載荷。目前對可回收火箭動力減速段氣動特性研究相對較少,且均從特定工程問題展開,對其氣動特性規(guī)律及形成機理缺乏系統(tǒng)性的認知。
目前關(guān)于反向射流的研究主要集中于火星動力減速、超音速減阻以及熱防護[8-10]領(lǐng)域,這些領(lǐng)域內(nèi)射流與來流總壓比遠小于運載火箭動力減速的情況,相關(guān)研究[11-12]則表明總壓比是控制反向射流流動特性的關(guān)鍵參數(shù)。近年來有學(xué)者提出新的表征參數(shù),例如Rong[13]提出的參數(shù)RPA,RPA增加時來流弓形激波遠離鈍體表面,同時阻力減小。這類參數(shù)可以直接表征反向射流的流動特征,但無法較好地預(yù)示箭體氣動特性變化規(guī)律。
本文針對反向射流作用下火箭氣動特性變化規(guī)律等問題開展數(shù)值模擬研究,完成馬赫數(shù)4.0~6.0、高度30~50 km范圍內(nèi)單噴管構(gòu)型火箭倒飛狀態(tài)氣動特性分析,獲得飛行馬赫數(shù)、飛行高度以及迎角對升阻特性的影響規(guī)律以及影響機理,并提出反映氣動特性變化規(guī)律的表征參數(shù),以期對可回收火箭工程應(yīng)用提供指導(dǎo)。
受連續(xù)流超聲速風(fēng)洞試驗成本的限制,本文開展了火箭升阻力特性的數(shù)值模擬研究。研究對象以Langley SRP(supersonic retro-propulsion)實驗[14]單噴管構(gòu)型為基礎(chǔ)。為了在幾何尺度上接近運載火箭,將標模直徑放大至2.54 m。噴管幾何形狀按實驗?zāi)P偷缺壤糯?噴管出口直徑0.892 m。
圖1為本文研究所使用的網(wǎng)格,采用四面體非結(jié)構(gòu)網(wǎng)格,總網(wǎng)格數(shù)量約為500萬。所有物面附近劃分邊界層網(wǎng)格,根據(jù)流動特點對邊界層網(wǎng)格厚度進行調(diào)整,即保證邊界層網(wǎng)格可以覆蓋實際流動的附面層,也保證物面附近網(wǎng)格的貼體性。根據(jù)反向射流形狀設(shè)置空間錐形加密區(qū)域。考慮到無噴射狀態(tài)下脫體激波貼近火箭底部,在噴管出口附近設(shè)置額外的網(wǎng)格加密區(qū)域。

圖1 本文研究所采用的計算網(wǎng)格Fig.1 Computational mesh of study objects
數(shù)值模擬采用NNW-FlowStar軟件,該軟件平臺是中國新一代基于有限體積法CFD軟件,以面向?qū)ο蟮乃枷?設(shè)計開發(fā)了具有高可復(fù)用性、高并行可擴展性和高靈活多學(xué)科耦合性等特點的軟件框架[15]。控制方程為直角坐標系下三維可壓縮定常Navier-Stokes方程,包含質(zhì)量、動量與能量方程。目前國內(nèi)外關(guān)于反向射流的學(xué)術(shù)研究以單一組分氣體為主,主要考察在反向射流作用下箭體力、熱環(huán)境的變化,故本文的研究使用單一組分氣體(空氣)。湍流模型采用S-A(Spalart-Allmara)一方程湍流模型[16],S-A模型專門為涉及壁面有界流動的航空航天應(yīng)用而設(shè)計,對承受逆壓力梯度的邊界層流動的求解有著較好效果[17]。
對流項離散采用HLLE++格式,該格式兼具Roe格式的數(shù)值精度與HLLE+格式的穩(wěn)定性,具有可在2個網(wǎng)格內(nèi)捕捉到激波、無熵增、邊界層耗散小、精確捕捉接觸間斷、可大大降低紅玉現(xiàn)象的發(fā)生等優(yōu)點。時間項采用隱式LU-SGS方法。
數(shù)值模擬中所有物面均采用等壁溫?zé)o滑移壁面。流動的5個基本變量分別為密度ρ、壓力p以及速度(含3個分量)。發(fā)動機入口處的流動為亞聲速流動,流向速度使用Neumann條件,其余4個變量為Drichlet條件。來流為超聲速條件,5個流動基本變量為Drichlet條件。
以獵鷹9運載火箭為參考對象,通過SpaceX公司發(fā)布的數(shù)次任務(wù)中一子級遙測數(shù)據(jù)(由官方發(fā)射視頻中包含的高度、速度信息進行轉(zhuǎn)換)提取其高度-馬赫數(shù)關(guān)系,如圖2所示。

圖2 獵鷹9數(shù)次任務(wù)動力減速段高度-馬赫數(shù)關(guān)系曲線Fig.2 Height-mach number relationship curve of Falcon 9 first stage propulsive descent during multiple missions
動力減速段啟動高度一般在50~60 km,啟動時馬赫數(shù)分布于7.0~7.5之間。關(guān)機高度位于35~45 km之間,馬赫數(shù)在5.0附近。本文的研究選擇高度30~50 km的區(qū)間范圍,這一區(qū)間可以包括獵鷹9任務(wù)中主要的飛行任務(wù)段。高度60 km開始氣體性質(zhì)逐漸偏離連續(xù)介質(zhì)并導(dǎo)致求解誤差逐漸增加,且前期研究發(fā)現(xiàn)50 km高度時發(fā)動機射流已經(jīng)基本達到膨脹的極限,故認為反向射流的性質(zhì)在50~60 km區(qū)間內(nèi)不會發(fā)生顯著變化。
飛行馬赫數(shù)選擇為Ma=4.0,Ma=5.0,Ma=6.0,這一范圍可以包含絕大部分獵鷹9運載火箭動力減速段的速度范圍。馬赫數(shù)更高的情況下,反向射流流動特征沒有明顯的改變,且過高的馬赫數(shù)會對數(shù)值模擬的穩(wěn)定性產(chǎn)生一定的負面影響,故本文研究中最高飛行馬赫數(shù)選為6.0。
在數(shù)值模擬方法的驗證工作中,以NASA Langley的單噴管構(gòu)型實驗為驗證對象,選擇其推力系數(shù)CT=4、Ma=2.4工況以及Ma=4.6無射流工況進行數(shù)值模擬校驗。提取鈍體底面沿徑向方向的壓力系數(shù),實驗數(shù)據(jù)及數(shù)值模擬結(jié)果如圖3所示,圖3中“Exp”為實驗數(shù)據(jù),“CFD”為本文數(shù)值模擬數(shù)據(jù),后綴“Ma=4.6”為馬赫數(shù)。Ma=4.6工況的數(shù)據(jù)顯示在次坐標軸上。

圖3 RANS數(shù)值模擬與實驗結(jié)果對比Fig.3 Comparison between RANS numerical simulation and experimental results
由圖3可知,數(shù)值模擬方法得到的壓力系數(shù)分布與實驗測量值符合較好,表明本文所采用的數(shù)值模擬方法可以較好地計算反向射流作用下箭體氣動特性。
動力減速段中火箭底部(安裝發(fā)動機的部分)指向飛行方向,這與常規(guī)意義上火箭的飛行狀態(tài)恰好相反。本文的研究將火箭底部視為飛行器的頭部,對升阻力、迎角等定義則沿用飛行器的常規(guī)定義方式。升阻力通過來流動壓與箭體截面的乘積進行無量綱化得到對應(yīng)的升阻系數(shù)。
在本文的升阻特性研究中,飛行馬赫數(shù)設(shè)定為3個狀態(tài),即Ma=4.0,Ma=5.0,Ma=6.0,飛行高度設(shè)定為3個狀態(tài),即H=30 km、H=40 km、H=50 km,飛行馬赫數(shù)與飛行高度組合后共9個飛行條件。每個飛行條件選取4個迎角狀態(tài),迎角分別為0°、3°、6°、10°。
圖4為火箭各部件升力系數(shù)曲線,包括箭體側(cè)壁面和底部的升力系數(shù),分別以Body和Bottom表示。由圖4可知,有迎角時箭體底部產(chǎn)生的升力系數(shù)明顯小于側(cè)壁面,在工程應(yīng)用中可以忽略底部產(chǎn)生的升力。

圖4 火箭各部件升力系數(shù)曲線Fig.4 Lift coefficient curve of rocket components
9種飛行狀態(tài)中側(cè)壁面升力系數(shù)曲線均出現(xiàn)不同程度的非線性增長現(xiàn)象,且飛行高度越低,非線性增長現(xiàn)象越明顯。其中H=50 km時非線性現(xiàn)象尚不明顯;H=40 km時,迎角0°~6°區(qū)間內(nèi)線性度較好,而6°~10°時升力系數(shù)增量較大,使得升力系數(shù)曲線出現(xiàn)一定程度的非線性增長現(xiàn)象;H=30 km時則出現(xiàn)明顯的非線性增長現(xiàn)象。
圖5為全箭升力系數(shù)分布曲線圖,包含4種組合狀態(tài),即Ma=4.0,Ma=6.0、H=30 km,H=50 km,每種狀態(tài)包括無射流和有射流兩種情況。圖例中后綴N為無射流情況,均以虛線顯示,有射流情況均以實線顯示。

圖5 全箭升力系數(shù)曲線Fig.5 Lift coefficient curve for rocket
無射流、Ma=6.0工況在小迎角時出現(xiàn)了負升力系數(shù),無射流時倒飛火箭為典型的高超聲速鈍頭體情況,結(jié)合底部幾何特征即可得到負升力系數(shù)的形成原因,這里不再進行贅述。
有射流時Ma=4.0,H=30 km與Ma=6.0,H=30 km兩種情況的升力系數(shù)曲線增長情況相對接近,Ma=4.0,H=50 km與Ma=6.0,H=50 km更為接近,說明在有射流情況下,升力系數(shù)曲線受飛行高度影響較大,而受飛行馬赫數(shù)影響較小,這與無射流情況恰好相反。無射流時相同飛行馬赫數(shù)的曲線分布更為接近。
圖6為箭體底部與側(cè)壁面阻力系數(shù)曲線。H=40 km,H=50 km底部阻力系數(shù)分布情況相對接近且數(shù)值均較小,H=50 km時部分馬赫數(shù)情況下出現(xiàn)負阻力系數(shù)。相比于其他兩種高度,H=30 km時底部阻力系數(shù)明顯增加,特別是馬赫數(shù)增加時,大迎角情況下底部阻力有大幅提高。

圖6 火箭各部件阻力系數(shù)曲線Fig.6 Drag coefficient curve of rocket components
箭體側(cè)壁面的阻力系數(shù)曲線在分布上則更為集中且相似,迎角增加時箭體側(cè)壁面將產(chǎn)生升阻力。當(dāng)飛行高度降低時,側(cè)壁面阻力系數(shù)將增加。
圖7為全箭阻力系數(shù)曲線圖,包含4種組合狀態(tài),即Ma=4.0,Ma=6.0、H=30 km,H=50 km,每種狀態(tài)包括無射流和有射流兩種情況。圖例中后綴N為無射流,并以虛線顯示。

圖7 全箭阻力系數(shù)曲線Fig.7 Drag coefficient curve for rocket
無射流時阻力系數(shù)曲線均集中于3.1附近,這種情況為典型的超聲速鈍頭體流動。有射流時全箭阻力系數(shù)大幅減小。Peterson[18]在推力系數(shù)0~1.6的超聲速反向射流實驗研究中發(fā)現(xiàn)阻力系數(shù)隨推力系數(shù)的增加而減小,推力系數(shù)超過1.5時阻力系數(shù)逐漸趨于0.1。其他超音速反向射流研究亦觀察到阻力系數(shù)隨推力系數(shù)增加而減小的現(xiàn)象[19-20]。
上述研究指出,反向射流對箭體升阻特性產(chǎn)生了較大影響,這種影響必然通過改變箭體周圍氣動環(huán)境實現(xiàn)。因此需要對反向射流作用下回流區(qū)流場狀態(tài)進行研究,并為分析氣動特性改變機理提供依據(jù)。圖8為Ma=6.0對稱面馬赫數(shù)分布云圖。發(fā)動機中心處的流線為回流區(qū)中發(fā)動機射流與來流氣體的分界線。

圖8 Ma=6.0情況下馬赫數(shù)分布云圖及部分射流流線Fig.8 Mach number distribution contour map and partial jet streamlines for Ma=6.0 condition
H=50 km時發(fā)動機射流在極低壓力下以大角度擴張并形成較大直徑的反向射流,通過終端激波(terminal shock)后,發(fā)動機射流較高的徑向速度將導(dǎo)致射流回流區(qū)域半徑明顯大于終端激波半徑,并在箭體周圍形成低速流動區(qū)域。H=40 km的流動情況亦是如此,區(qū)別在于環(huán)境壓力有所提高,反向射流在徑向方向上的膨脹程度較低,導(dǎo)致回流區(qū)直徑也有一定程度的降低。
H=30 km時,發(fā)動機射流并沒有以大角度擴張,而是呈現(xiàn)出“棗核”形狀。盡管沒有明顯的正激波形式的射流終端激波,氣體在到達射流終端后仍然以相對緩慢的速度進行轉(zhuǎn)折并向后流動(即遠場來流方向)。
綜上所述,發(fā)動機射流回流后會將火箭箭體“包裹”,并顯著改變箭體所處的氣動環(huán)境。高空階段反向射流直徑較大,對箭體產(chǎn)生的包裹作用較強,低空階段反向射流直徑明顯減小,發(fā)動機射流回流后緊貼箭體表面。
上述關(guān)于升阻特性的研究發(fā)現(xiàn)反向射流會對氣動特性產(chǎn)生明顯影響,本文將對以下3點主要變化進行說明,并分析其改變的機理。
1)升阻特性對高度敏感,對飛行馬赫數(shù)相對不敏感。無射流時飛行高度變化直接影響環(huán)境靜壓,不會對飛行器外流場流動結(jié)構(gòu)產(chǎn)生明顯影響,故無量綱化的升阻系數(shù)變化較小。而反向射流對飛行器外流場產(chǎn)生顯著影響,環(huán)境靜壓的改變將直接影響反向射流的幾何尺度,進而對飛行器升阻特性產(chǎn)生明顯影響。
2)氣動阻力大幅降低。由上述研究可知,反向射流的存在使得箭體周圍形成低速流動區(qū)域,進而降低了箭體氣動阻力。部分工況下射流的引射作用還將進一步降低底部壓力,并產(chǎn)生負阻力。
3)升力系數(shù)非線性增長現(xiàn)象。對于通常的超聲速飛行器,其升力系數(shù)曲線一般隨迎角呈現(xiàn)線性增加的趨勢。而火箭在動力減速段中的升力系數(shù)曲線則出現(xiàn)了明顯的非線性增長趨勢。為了對這種現(xiàn)象的形成機理進行分析,這里以Ma=6.0,H=30 km工況作為示例進行研究,圖9為兩個迎角時壓力分布云圖。顯示時為了不被發(fā)動機內(nèi)過高的壓力影響,圖9中高于10 kPa的部分在顯示上達到飽和。

圖9 Ma=6.0壓力分布云圖Fig.9 Pressure distribution contour map for Ma=6.0
AoA=3°時反向射流出現(xiàn)一定程度的不對稱情況,箭體底部迎風(fēng)面出現(xiàn)較弱的脫體激波并產(chǎn)生一定范圍的高壓區(qū)。AoA=10°時來流弓形激波進一步向箭體迎風(fēng)面偏轉(zhuǎn),與箭體底部迎風(fēng)面的脫體激波相干涉,進一步增加底部附近的壓力并提高了此工況下的升力系數(shù)。由此可見,升力系數(shù)的非線性增長現(xiàn)象實質(zhì)是反向射流幾何尺度較小時,來流弓形激波靠近箭體而出現(xiàn)的升力系數(shù)額外增加。
若火箭動力減速段結(jié)束較晚或發(fā)動機總壓較小,則有可能形成直徑較小的反向射流,此時將會出現(xiàn)較為明顯的升力系數(shù)非線性增長現(xiàn)象,壓力分布的變化也將對縱向力矩系數(shù)產(chǎn)生影響。工程應(yīng)用中可從增加反向射流直徑與限制火箭飛行迎角兩個角度出發(fā)實現(xiàn)降低對氣動特性的影響。
上述研究表明反向射流對箭體氣動特性產(chǎn)生了明顯影響,升阻特性出現(xiàn)對飛行高度敏感的現(xiàn)象,故需要提出一種表征參數(shù)以實現(xiàn)對箭體氣動特性的快速預(yù)估。若單純以飛行高度作為箭體升阻特性的表征參數(shù),則無法體現(xiàn)反向射流對升阻特性影響的機理。反向射流的遮擋作用是箭體升阻特性改變的根本原因,從影響機理出發(fā)可以提出表征參數(shù)Dj/Db,即反向射流寬度(直徑)與箭體直徑之比,可直接反映反向射流對箭體的遮擋比例,如圖10所示。利用本文數(shù)值模擬計算得到的單噴管反向射流的主要幾何參數(shù),計算各工況情況下寬度比數(shù)值,見表1。

表1 反向射流表征參數(shù)Tab.1 Characteristic parameters of opposing jet

圖10 反向射流表征參數(shù)示意Fig.10 Schematic diagram of characteristic parameters for opposing jet
反向射流將導(dǎo)致火箭氣動特性對飛行高度較為敏感,結(jié)合火箭升阻特性曲線,發(fā)現(xiàn)Dj/Db可以較好的表征反向射流對箭體氣動特性的影響:
1)Dj/Db<1.0時,反向射流不能完全形成對箭體的遮擋,升力系數(shù)的增長將出現(xiàn)較為明顯的非線性現(xiàn)象。Dj/Db越小,全箭阻力系數(shù)越大,特別是小于0.5時,還存在較大的箭體底部阻力。
2)Dj/Db>1.0時,反向射流對箭體的遮擋作用明顯增強,全箭升阻系數(shù)將進一步減小,阻力系數(shù)逐漸接近于0,底部甚至出現(xiàn)負阻力系數(shù)。
表1中還加入了各工況下10°迎角情況下的升力系數(shù)以及零升阻力系數(shù)。從表1中可以觀察到,升力系數(shù)的整體變化規(guī)律是隨Dj/Db的減小而增加,但是在相同高度情況下,升力系數(shù)卻隨Dj/Db的減小而降低。在相同高度情況下Dj/Db減小時,受反向射流的影響(箭體周圍產(chǎn)生低速流動區(qū)域),箭體絕對升力有一定程度的增加,但參考速度增加更為明顯,反而導(dǎo)致升力系數(shù)出現(xiàn)一定程度的下降。而當(dāng)Dj/Db的數(shù)值接近0.5時,反向射流對箭體氣動特性的影響明顯變?nèi)?此時升力系數(shù)的變化程度開始明顯變小。
Dj/Db降低時,零升阻力系數(shù)逐漸增加,表明此參數(shù)可以較好反應(yīng)箭體零升阻力系數(shù)的變化規(guī)律。直至0.5附近時零升阻力系數(shù)開始下降,此時反向射流的直徑已明顯小于箭體直徑,反向射流對箭體氣動特性的影響明顯變?nèi)?且反向射流長度有一定程度的增加,這將減小來流弓形激波的強度,進而減小箭體零升阻力系數(shù)。
Dj/Db同時受到飛行高度以及飛行馬赫數(shù)的影響,但受飛行高度的影響最為明顯。在火箭動力減速過程中,Dj/Db將逐漸降低,反向射流對箭體氣動特性的影響也將逐漸減小。若發(fā)動機參數(shù)發(fā)生變化并導(dǎo)致反向射流直徑降低時,箭體氣動特性受反向射流的影響亦會降低。
1)反向射流會對箭體產(chǎn)生遮擋作用,回流區(qū)內(nèi)形成的低壓、低速流動會對箭體氣動特性產(chǎn)生明顯的影響。存在反向射流時,火箭的升阻力特性受飛行高度影響較大,而受飛行馬赫數(shù)的影響較低,與無射流情況存在明顯的差異。
2)飛行高度相對較低時,來流脫體激波距離箭體較近,在有迎角來流作用下,來流弓形激波在箭體底部附近的再附著將導(dǎo)致升力的額外增加,使升力系數(shù)出現(xiàn)明顯的非線性增長現(xiàn)象。
3)H=40 km及以上時,反向射流的遮擋作用以及引射作用將大幅減小箭體底部的阻力,阻力系數(shù)降至0.1以下,甚至在部分飛行工況下產(chǎn)生負阻力。
4)提出以寬度比Dj/Db作為箭體氣動特性變化規(guī)律的表征參數(shù)。Dj/Db<1.0時,反向射流無法對箭體形成完全的遮擋作用。其數(shù)值越小,反向射流對箭體升阻特性的影響越小。