姚 瑤,高 波
(北京宇航系統(tǒng)工程研究所,北京 100073)
由于超聲速飛行的速度特點(Ma>5,例如X-43,Ma=9.8)以及其高度特點(H>30 000 m),在飛行器表面會形成厚邊界層,另外由于飛行器的幾何外形特點,在超聲速氣流中會產(chǎn)生一系列的激波,激波打在飛行器物面上會發(fā)生激波反射,遭遇物面的邊界層會形成激波/邊界層干擾現(xiàn)象。在高超聲速飛行器流場的多個部位存在著形式多樣的激波/邊界層干擾現(xiàn)象。如在進氣道內(nèi)流場中,存在著復(fù)雜的激波反射、激波相交以及激波/邊界層干擾現(xiàn)象,如圖1所示,可能造成結(jié)構(gòu)破壞和快速疲勞,氣動彈性以及熱防護困難,進氣氣流畸變,甚至不啟動等等,對發(fā)動機性能及穩(wěn)定性等有著重要的影響,在設(shè)計中需要充分考慮。
圖1 超燃沖壓發(fā)動機示意圖[1]Fig.1 Scramjet schematic[1]
1939年,Ferri[2]首次在試驗中觀察到機翼后緣附近的激波/邊界層干擾現(xiàn)象。隨后,Fage和Sargent[3]、Barry等[4]、Liepmann等[5]以及Johannesen[6]等均在試驗中觀察到了不同類型的激波/邊界層干擾現(xiàn)象。典型的二維激波邊界層干擾有入射斜激波反射、壓縮拐角、正激波、壓力突增以及迎風(fēng)臺階[7]。Green[8]基于大量試驗數(shù)據(jù)認(rèn)為壓縮拐角干擾與入射激波邊界層干擾流場具有相似性,并指出當(dāng)拐角為入射激波對應(yīng)的氣流偏轉(zhuǎn)角兩倍時,兩者的干擾近似一致。
隨著計算技術(shù)的發(fā)展,1971年MacCormack[9]首次開展了激波與層流邊界層干擾的數(shù)值模擬工作。Beam等[10]研究了復(fù)雜幾何條件的激波與層流邊界層干擾問題,采用有限差分法求解可壓縮的N-S方程。Bodonyi等[11]數(shù)值模擬跨聲速流的入射激波層流邊界層干擾,發(fā)現(xiàn)分離區(qū)和再附區(qū)的大小隨激波強度增加變化較大。Mailison[12]利用試驗方法研究壓縮拐角流動,分析表面壓力與熱流密度的分布規(guī)律以及確定初始分離出現(xiàn)的條件。結(jié)果表明分離流動的出現(xiàn)及分離區(qū)的大小與壓縮角的大小相關(guān)。Price等[13]研究了高超聲速湍流壓縮拐角激波邊界層干擾,隨著壓縮角的增加,壓力和熱流迅速上升,但壓力分布在壓縮面上緩慢地趨于平緩,而熱流在達到峰值后迅速下降。龔安龍等數(shù)值仿真雙錐分離流動,認(rèn)為分離流動特性及流場物理參數(shù)對計算網(wǎng)格比較敏感,過稀的網(wǎng)格將導(dǎo)致分離區(qū)減小,壓力和熱流峰值增大[14]。趙云飛采用WCNS格式數(shù)值分析非定常激波/邊界層干擾對于邊界層轉(zhuǎn)捩的影響[15]。童福林等采用DNS數(shù)值模擬壓縮拐角激波/轉(zhuǎn)捩邊界層干擾,研究轉(zhuǎn)捩過程對干擾區(qū)內(nèi)分離泡結(jié)構(gòu)的影響[16]。董祥瑞等采用DES方法對激波邊界層干擾誘導(dǎo)的分離流動進行數(shù)值模擬,并采用單、雙微楔對齊進行控制,獲得其最佳安裝位置[17]。宋友富等發(fā)展可壓縮湍流模型CKDO,其對壁面壓力和摩擦阻力系數(shù)的捕捉能力優(yōu)于其他模型[18]。陳宣亮等采用分離渦模擬方法,對X-37B機翼副翼偏轉(zhuǎn)拐角的激波邊界層干擾做非定常數(shù)值分析,捕捉該類干擾下的脈動壓力[19]。王嬌等探索Bump進氣道中鼓包誘導(dǎo)的錐形激波和機身湍流邊界層干擾問題,同樣采用數(shù)值仿真,認(rèn)為鼓包最適合作為超聲速進氣道前緣壓縮面[20]。2016年,Kotov等考慮通過亞格子模型以及數(shù)字模擬兩方面的改進來實現(xiàn)湍流流動激波附近的精度提高[21]。
1958年,Chapman等[22]基于大量的試驗數(shù)據(jù),發(fā)現(xiàn)分離流動的一些特點并不依賴于引起分離的形式。把這種不受下游幾何形狀影響也不受引起分離的形式的影響的干擾叫做“自由干擾”。并針對自由干擾給出簡單的理論分析——自由干擾理論(FIT)。
在分離區(qū)長度方面,Burggraf[23]基于Neiland[24]以及Stewartson和Williams[25]提出的漸進理論,針對激波層流邊界層干擾產(chǎn)生較強分離時的分離區(qū)長度給出數(shù)學(xué)模型。同時,對于上游影響長度,即無黏激波入射點至干擾起始點的長度,Miller等[26]以及Ramesh等[27]學(xué)者也給出了工程擬合公式。
2011年,Babinsky和Harvey基于前人的試驗結(jié)果描繪入射激波層干擾分離流場的流動圖畫和壁面壓力分布規(guī)律。典型流動結(jié)構(gòu)包括分離泡、入射激波、膨脹波以及反射激波,壓增主要分為分離階段以及再附階段兩個部分[7]。
前人基于大量的試驗以及數(shù)值手段,已經(jīng)給出了分離流動的基本流動結(jié)構(gòu)以及典型流動參數(shù)的分布規(guī)律。但是尚沒有完整的理論求解模型,干擾的關(guān)鍵點位置不能明確,對于波系結(jié)構(gòu)的定量描述也少有研究。工程實際中,迫切地需要快速獲得流動圖畫,定位波系位置,明確流動參數(shù)具體分布,捕獲高壓高熱流區(qū),為整體流場分析提供幫助。
本文基于對入射激波致邊界層分離流場的分析建立完整的理論求解模型,該模型可以根據(jù)來流條件(來流馬赫數(shù)、入射激波強度以及邊界層參數(shù))給出局部干擾流場結(jié)構(gòu)以及流動參數(shù)分布。
本文用ρ、p、V、a、Ma、γ分別表征密度、壓力、速度、聲速、馬赫數(shù)以及比熱比,考慮的流體為量熱完全氣體。另外,激波角(亦即激波與上游流體的夾角)用β表示,氣流穿越激波或膨脹波的偏轉(zhuǎn)角用θ表示。
入射激波與平板邊界層發(fā)生干擾,當(dāng)入射激波足夠強時,引起邊界層分離,產(chǎn)生分離泡,如圖2所示。
分離泡起點S上游流體發(fā)生偏轉(zhuǎn),產(chǎn)生壓縮波系,壓縮波系合并為分離激波。入射激波與分離激波發(fā)生正規(guī)相交,產(chǎn)生透射激波。透射激波穿越自由剪切層并反射為膨脹波,使得自由剪切層向壁面偏轉(zhuǎn),到達再附點R,R下游流體轉(zhuǎn)為與壁面平行,產(chǎn)生壓縮波系,壓縮波系合并為再附激波[7]。因此,已知的主要結(jié)構(gòu)包括入射激波、邊界層、分離激波、干擾點膨脹波、分離泡以及再附激波。工程實際中已知的參數(shù)通常包括來流馬赫數(shù)、干擾點雷諾數(shù)以及入射激波強度。由此需要獲得各區(qū)流動參數(shù)、各干擾波的方向和強度以及壓力沿壁面分布等。
圖2 入射激波致邊界層分離流場結(jié)構(gòu)Fig.2 Sketchoftheflowinducedbyanincidentshock waveboundarylayerinteractionwithseparation
入射激波致邊界層分離流場的壁面壓力分布如圖3所示,分離點S上游的i點(干擾起始點)壓力即開始增加,對應(yīng)干擾起始點壓力pi,通常近似認(rèn)為pi等于來流壓力。伴隨分離階段,壓力歷經(jīng)第一個突增,分離階段下游接著為分離流動的典型壓力平臺段。激波穿越邊界層內(nèi)的外部區(qū)域,但是內(nèi)部亞聲速流不能承受較大的壓力梯度,因而局部形成一個近似等壓邊界。下游伴隨再附階段,壓力經(jīng)歷第二次增加,達到最終壓力峰值pF。相比分離階段壓增,再附階段壓增比較平緩[7]。
圖3 入射激波致邊界層分離流場壁面壓力Fig.3 Wallpressuredistributioninducedbyanincidentshock waveboundarylayerinteractionwithseparation
Chapmen等[21]基于大量的試驗數(shù)據(jù),發(fā)現(xiàn)分離流動的一些特點并不依賴于引起分離的形式。并把這種不受下游幾何形狀影響也不受引起分離的形式的影響的干擾叫做“自由干擾”,同時對于這種干擾給出簡單的理論分析,即自由干擾理論(FIT),得出平臺壓力pP:
式中:q0為來流動壓,Cfi為壁面摩擦力系數(shù),無量綱函數(shù)F被認(rèn)為是普適函數(shù),和馬赫數(shù)以及雷諾數(shù)無關(guān),由試驗確定,見表1。
表1 特定點處的F值Table1 ValuesofFatspecialpoints
Katzer[28]針對入射激波層流邊界層干擾,進行絕熱壁面情況下的數(shù)值模擬,其模擬工況范圍為1.4≤Ma∞≤3.4,1×105≤ReL≤6×105,1.2≤pF/p∞≤3.4?;谶@些數(shù)值數(shù)據(jù),給出了入射激波層流邊界層干擾的分離區(qū)長度工程擬合公式:
式中:P=(pF-pinc)/p∞,為未受干擾時的平板無黏激波入射點xL處邊界層位移厚度:
pinc為引起初始分離的激波強度,通過自由干擾理論給出:
CfL為未受干擾時平板邊界層xL處的壁面摩擦系數(shù):
常數(shù)Pinc=1.85。
Krishnan等[29]將Katzer公式的范圍拓展到Ma∞≤6.85和pF/p∞≤5.89:
本文運用數(shù)值方法并結(jié)合文獻數(shù)據(jù)進行驗證,如圖4所示。
數(shù)值計算數(shù)據(jù)基本與擬合公式吻合。由于擬合公式適用范圍所限,若馬赫數(shù)較大,其誤差將超出工程允許的范圍。
圖4 分離區(qū)長度驗證[29]Fig.4 Validation of fitting formula by Krishnan et al[29]
基于第1節(jié)對激波致邊界層分離流場的分析,將局部流場結(jié)構(gòu)進行簡化,建立理論模型,如圖5所示。
首先將分離泡假定為一個三角形。壁面壓力分布大致分為如下幾個階段:干擾上游與無干擾平板邊界層分布基本一致,分離階段經(jīng)歷第一個壓增,隨后為壓力平臺,再附階段經(jīng)歷第二個壓增,到達壓力峰值。由此假設(shè)第一個壓增完全由分離激波(C2)引起,分離角θS為引起分離激波(C2)的氣流偏轉(zhuǎn)角,波后(2)區(qū)壓力為平臺壓力。第二個壓增完全由再附激波(C6)引起,再附角θR為引起再附激波(C6)的氣流偏轉(zhuǎn)角,波前(6)區(qū)壓力亦為平臺壓力。入射激波(C1)與分離激波(C2)發(fā)生正規(guī)相交,產(chǎn)生透射激波(C3)與(C4)以及滑移線(Σ),(R4)反射為膨脹波,(6)區(qū)氣流轉(zhuǎn)平產(chǎn)生再附激波(C6)。該模型忽略分離激波和再附激波根部的壓縮波系以及各激波在穿越邊界層時的彎曲,將分離激波及再附激波簡化為直接從壁面引出的兩道平直斜激波。
以來流(0)區(qū)氣流狀態(tài)作為參考,運用極曲線理論描述流場結(jié)構(gòu),如圖6所示。p0/p0=1,θ0=0°,因此(0)區(qū)狀態(tài)點位于原點。所有由(0)區(qū)穿越斜激波可能達到的狀態(tài)由極曲線Γ0表征。(0)區(qū)氣流穿越入射激波(C1)到達(1)區(qū),因此根據(jù)區(qū)(1)相對于(0)區(qū)的壓力比以及氣流偏轉(zhuǎn)角(此處即為楔角θw,由于氣流向下偏轉(zhuǎn),為負值)可以確定其在極曲線Γ0上的位置。所有由(1)區(qū)穿越斜激波可能達到的狀態(tài)由極曲線Γ1表征。同時,(0)區(qū)氣流穿越分離激波(C2)到達(2)區(qū),由于假設(shè)(2)區(qū)壓力等于平臺壓力(即p2=pP),并且(2)區(qū)氣流必定向上偏轉(zhuǎn),由此可以確定其在極曲線Γ0上的位置,同時給出分離角θS的大小。所有由(2)區(qū)穿越斜激波可能到達的狀態(tài)由極曲線Γ2表征。入射激波(C1)與分離激波(C2)在H點發(fā)生正規(guī)相交,產(chǎn)生透射激波(C3),波后(3)區(qū);以及透射激波(C4),波后(4)區(qū)。(3)區(qū)與(4)區(qū)壓力相等且氣流平行,因此兩者狀態(tài)點位于極曲線Γ1和Γ2的交點上。緊接著,(4)區(qū)氣流穿越膨脹波到達(6)區(qū),由于假設(shè)(6)區(qū)壓力與(2)區(qū)一樣,同為平臺壓力(即p6=pP),因此(6)區(qū)狀態(tài)點位于經(jīng)過(2)區(qū)的等壓線(圖6中點劃線)與經(jīng)過(4)區(qū)的等熵膨脹極曲線Γ4的(圖6中虛線)的交點處。所有由(6)區(qū)穿越斜激波可能達到的狀態(tài)由極曲線Γ6表征。(6)區(qū)氣流穿越再附激波(C6)到達(8)區(qū),此時氣流轉(zhuǎn)為與壁面平行(即與來流一致,θS=0)。因此(8)區(qū)狀態(tài)點位于極曲線Γ6與豎軸的交點處。
圖6 入射激波致邊界層分離極曲線示意圖Fig.6 Polar of the flow induced by an incident shock wave boundary layer interaction with separation
首先,分離激波波后為平臺壓力,可由自由干擾理論獲得。分離角θS和分離激波(C2)波后參數(shù)可由激波關(guān)系式獲得。
入射激波(C1)與分離激波(C2)發(fā)生異側(cè)正規(guī)相交,產(chǎn)生透射激波以及滑移線,如圖7所示。
圖7 入射激波與分離激波異側(cè)相交Fig.7 Regular intersection of incident shock and separation shock
已知(0)區(qū)參數(shù)以及入射激波(C1),波后參數(shù)可由斜激波關(guān)系式獲得:
(1)區(qū)氣流穿越透射激波(C3):
(2)區(qū)氣流穿越透射激波(C4):
如圖8所示,(4)區(qū)氣流穿越分離泡頂點A處的膨脹波(RA)到達(6)區(qū),運用膨脹波關(guān)系式:
(6)區(qū)氣流穿越再附激波(C6)到達(8)區(qū),運用斜激波關(guān)系式:
簡化模型假設(shè)膨脹波波后(6)區(qū)的壓力與(2)區(qū)一樣,等于平臺壓力,即:
聯(lián)立式(7~18)就可以完整求解整個干擾流場波系結(jié)構(gòu)、各區(qū)流動參數(shù)以及各干擾點位置。
圖8 再附階段示意圖Fig.8 Schematic of reattachment
最后求解分離區(qū)長度以及幾何關(guān)系可以進一步確定分離泡個關(guān)鍵點的位置坐標(biāo)。
如圖9所示,在ΔASR中:
在ΔSHA中:
從總體結(jié)構(gòu)來看,如圖10所示,干擾點H的坐標(biāo)為:
由式(19~25)可確定各干擾點的具體位置。
圖10 整體結(jié)構(gòu)關(guān)鍵點示意Fig.10 Schematic of key points in whole flow filed
根據(jù)第2節(jié)建立的理論模型可以預(yù)測入射激波致邊界層分離流場,此處用數(shù)值計算結(jié)果進行驗證。
考慮工況Ma∞=2.15、θw=3.75°、ReL=9.91×104,對比數(shù)值計算與理論結(jié)果的局部流場結(jié)構(gòu),如圖11所示。入射激波(C1)、(C1)與分離激波(C2)干擾后的透射激波(C4)以及反射膨脹波(R)等波系位置與形狀均與數(shù)值計算結(jié)果吻合較好,同時分隔流線也與數(shù)值結(jié)果基本重合。理論預(yù)測的分離點位置ST=0.0633525 m,與數(shù)值計算結(jié)果的分離點位置SC=0.0610716 m相差3.7%;理論預(yù)測再附點位置RT=0.0964525 m,與數(shù)值計算結(jié)果的再附點位置RC=0.097169 m相差僅0.74%。由于理論模型里忽略了邊界層的逐漸變形過程,而是將分離點S以及再附點R的氣流偏轉(zhuǎn)作為突變處理,因此由分離點S以及再附點R處分別直接引出分離激波(C2)以及再附激波(C6),與數(shù)值計算結(jié)果中兩處形成壓縮波系進而合并為激波不同。各區(qū)流動參數(shù)的預(yù)測也相當(dāng)準(zhǔn)確,如表2所示。
圖11 理論預(yù)測與數(shù)值計算局部流場結(jié)構(gòu)對比Fig.11 Comparison of wave structure between theoretical prediction and numerical result
表2 理論預(yù)測與數(shù)值計算結(jié)果各區(qū)流場參數(shù)對比Table 2 Comparison of flow parameters
壁面壓力分布對比如圖12所示。由于理論模型忽略了激波根部的壓縮波系,因此壁面壓力分布呈臺階狀,不似實際情況的平滑過渡,但是關(guān)鍵位置如干擾起始點、壓力平臺以及壓力峰值均吻合較好,因此整體分布趨勢一致。
圖12 理論預(yù)測與數(shù)值計算壁面壓力分布對比Fig.12 Comparison of wall pressure distribution between theoretical prediction and numerical result
運用該理論模型可進行參數(shù)化分析,例如分離泡高度。
考慮工況Ma∞=2.15,θw=3.75°和ReL=9.91×104,根據(jù)式(3)可知無干擾情況下入射點xL處的邊界層位移厚度為9.1873×10-4m。而理論模型預(yù)測的分離泡的高度為1.0524×10-3m,相比沒有考慮干擾情況增厚的14.5%。此處考慮來流馬赫數(shù)、飛行高度以及氣流偏轉(zhuǎn)角對分離泡高度的影響,工況列于表3。分離泡高度與分離角及分離區(qū)長度直接相關(guān),分離角越大或者分離區(qū)長度越長則分離泡高度越大。而分離角與分離區(qū)長度主要由來流參數(shù)(包括來流馬赫數(shù)與壓力或者飛行高度)和壓縮角確定。
表3 分離角參數(shù)化分析的工況Table 3 Condition for parametric analysis
對于表3中的工況1,在相同的外壓縮角以及飛行高度下,來流馬赫數(shù)的影響如圖13所示。隨著馬赫數(shù)(2≤Ma∞≤6)的增大,分離泡的高度先略微降低隨后增大。這是由于隨著馬赫數(shù)的增大,平臺壓力緩慢減小(接近正比于Ma-1/2),分離角也隨之緩慢減小;而總壓增隨著馬赫數(shù)的增大迅速增大,導(dǎo)致分離區(qū)長迅速增長。兩者相互制約,小Ma時,分離區(qū)長度變化較小,分離角的變化起主導(dǎo)因素,因此分離區(qū)高度有輕微的降低;大Ma時,分離區(qū)長度變化劇烈,起主導(dǎo)因素,因此分離區(qū)高度隨馬赫數(shù)的增大逐漸增大。
圖13 分離泡高度隨來流馬赫數(shù)的變化Fig.13 Influence of incoming Mach number
對于表3中的工況2,在2°~15°范圍內(nèi),在相同的來流馬赫數(shù)以及飛行高度下,分離泡高度是入射激波氣流偏轉(zhuǎn)角的增函數(shù),如圖14所示。根據(jù)自由干擾理論,平臺壓力基本不隨外壓縮角變化,因而分離角也基本不變。而總壓增隨著外壓縮角的增大而增大,因此分離泡高度也隨之單調(diào)增大。
圖14 分離泡高度隨外壓縮角的變化Fig.14 Influence of flow deflection for incident shock
對于表3中的工況3,保持來流馬赫數(shù)以及外壓縮角不變,在20 km到50 km范圍內(nèi),分離泡高度隨飛行高度的升高而增大,在30 km高度時略有減小,如圖15所示。由自由干擾理論結(jié)合激波關(guān)系式可知,來流壓力的降低并不會改變分離角的大小,但引起初始分離的pinc隨之顯著降低,造成分離區(qū)長度增大,因而分離泡高度增大。30 km高度處壓力并未平緩降低,而是出現(xiàn)小的增加,因此分離區(qū)高度也略有減小。
圖15 分離泡高度隨飛行高度的變化Fig.15 Influence of flight altitude
本文針對激波致邊界層分離流場建立了理論預(yù)測模型。將分離泡簡化為三角形,總壓增分為分離階段以及再附階段兩部分。假設(shè)分離階段的壓增完全由分離激波完成,分離角為引起分離激波的氣流偏轉(zhuǎn)角;再附階段的壓增完全由再附激波完成,再附角為引起再附激波的氣流偏轉(zhuǎn)角。分離激波與再附激波之間的壁面壓力一直為平臺壓力,并運用極曲線方法描述了該局部流場波系結(jié)構(gòu)。該模型可以給出包括分離激波、反射膨脹波以及再附激波等波系結(jié)構(gòu)的位置和形狀以及各區(qū)流動參數(shù),并得到了很好的CFD驗證。運用該理論模型,發(fā)現(xiàn)分離泡高度隨來流馬赫數(shù)的增大先略有減小再增大,同時是外壓縮角以及飛行高度的增函數(shù)。該模型忽略分離激波和再附激波根部的壓縮波系以及各激波在穿越邊界層時的彎曲,將分離激波及再附激波簡化為直接從壁面引出的兩道平直斜激波。
本文所建立的理論模型是對真實流動結(jié)構(gòu)進行了一定簡化,后續(xù)還有以下工作可考慮開展:
(1)本文模型忽略了分離激波和再附激波根部的壓縮波系以及各激波在穿越邊界層時發(fā)生的彎曲??梢詫enderson[30]的方法添加到模型中,將邊界層分層,建立更為詳細的理論模型。
(2)本文采用的是理想氣體模型,對于高超聲速流動,伴隨高溫真實氣體效應(yīng),可以采用考慮復(fù)雜化學(xué)反應(yīng)的氣體模型,其對波系結(jié)構(gòu)以及流動參數(shù)均會產(chǎn)生影響。