入射激波邊界層干擾分離流場結構研究

姚 瑤,高 波

(北京宇航系統工程研究所,北京 100073)

0 引 言

由于超聲速飛行的速度特點(Ma＞5,例如X-43,Ma＝9.8)以及其高度特點(H＞30 000 m),在飛行器表面會形成厚邊界層,另外由于飛行器的幾何外形特點,在超聲速氣流中會產生一系列的激波,激波打在飛行器物面上會發生激波反射,遭遇物面的邊界層會形成激波/邊界層干擾現象。在高超聲速飛行器流場的多個部位存在著形式多樣的激波/邊界層干擾現象。如在進氣道內流場中,存在著復雜的激波反射、激波相交以及激波/邊界層干擾現象,如圖1所示,可能造成結構破壞和快速疲勞,氣動彈性以及熱防護困難,進氣氣流畸變,甚至不啟動等等,對發動機性能及穩定性等有著重要的影響,在設計中需要充分考慮。

圖1 超燃沖壓發動機示意圖[1]Fig.1 Scramjet schematic[1]

1939年,Ferri[2]首次在試驗中觀察到機翼后緣附近的激波/邊界層干擾現象。隨后,Fage和Sargent[3]、Barry等[4]、Liepmann等[5]以及Johannesen[6]等均在試驗中觀察到了不同類型的激波/邊界層干擾現象。典型的二維激波邊界層干擾有入射斜激波反射、壓縮拐角、正激波、壓力突增以及迎風臺階[7]。Green[8]基于大量試驗數據認為壓縮拐角干擾與入射激波邊界層干擾流場具有相似性,并指出當拐角為入射激波對應的氣流偏轉角兩倍時,兩者的干擾近似一致。

隨著計算技術的發展,1971年MacCormack[9]首次開展了激波與層流邊界層干擾的數值模擬工作。Beam等[10]研究了復雜幾何條件的激波與層流邊界層干擾問題,采用有限差分法求解可壓縮的N-S方程。Bodonyi等[11]數值模擬跨聲速流的入射激波層流邊界層干擾,發現分離區和再附區的大小隨激波強度增加變化較大。Mailison[12]利用試驗方法研究壓縮拐角流動,分析表面壓力與熱流密度的分布規律以及確定初始分離出現的條件。結果表明分離流動的出現及分離區的大小與壓縮角的大小相關。Price等[13]研究了高超聲速湍流壓縮拐角激波邊界層干擾,隨著壓縮角的增加,壓力和熱流迅速上升,但壓力分布在壓縮面上緩慢地趨于平緩,而熱流在達到峰值后迅速下降。龔安龍等數值仿真雙錐分離流動,認為分離流動特性及流場物理參數對計算網格比較敏感,過稀的網格將導致分離區減小,壓力和熱流峰值增大[14]。趙云飛采用WCNS格式數值分析非定常激波/邊界層干擾對于邊界層轉捩的影響[15]。童福林等采用DNS數值模擬壓縮拐角激波/轉捩邊界層干擾,研究轉捩過程對干擾區內分離泡結構的影響[16]。董祥瑞等采用DES方法對激波邊界層干擾誘導的分離流動進行數值模擬,并采用單、雙微楔對齊進行控制,獲得其最佳安裝位置[17]。宋友富等發展可壓縮湍流模型CKDO,其對壁面壓力和摩擦阻力系數的捕捉能力優于其他模型[18]。陳宣亮等采用分離渦模擬方法,對X-37B機翼副翼偏轉拐角的激波邊界層干擾做非定常數值分析,捕捉該類干擾下的脈動壓力[19]。王嬌等探索Bump進氣道中鼓包誘導的錐形激波和機身湍流邊界層干擾問題,同樣采用數值仿真,認為鼓包最適合作為超聲速進氣道前緣壓縮面[20]。2016年,Kotov等考慮通過亞格子模型以及數字模擬兩方面的改進來實現湍流流動激波附近的精度提高[21]。

1958年,Chapman等[22]基于大量的試驗數據,發現分離流動的一些特點并不依賴于引起分離的形式。把這種不受下游幾何形狀影響也不受引起分離的形式的影響的干擾叫做“自由干擾”。并針對自由干擾給出簡單的理論分析——自由干擾理論(FIT)。

在分離區長度方面,Burggraf[23]基于Neiland[24]以及Stewartson和Williams[25]提出的漸進理論,針對激波層流邊界層干擾產生較強分離時的分離區長度給出數學模型。同時,對于上游影響長度,即無黏激波入射點至干擾起始點的長度,Miller等[26]以及Ramesh等[27]學者也給出了工程擬合公式。

2011年,Babinsky和Harvey基于前人的試驗結果描繪入射激波層干擾分離流場的流動圖畫和壁面壓力分布規律。典型流動結構包括分離泡、入射激波、膨脹波以及反射激波,壓增主要分為分離階段以及再附階段兩個部分[7]。

前人基于大量的試驗以及數值手段,已經給出了分離流動的基本流動結構以及典型流動參數的分布規律。但是尚沒有完整的理論求解模型,干擾的關鍵點位置不能明確,對于波系結構的定量描述也少有研究。工程實際中,迫切地需要快速獲得流動圖畫,定位波系位置,明確流動參數具體分布,捕獲高壓高熱流區,為整體流場分析提供幫助。

本文基于對入射激波致邊界層分離流場的分析建立完整的理論求解模型,該模型可以根據來流條件(來流馬赫數、入射激波強度以及邊界層參數)給出局部干擾流場結構以及流動參數分布。

1 干擾流場

本文用ρ、p、V、a、Ma、γ分別表征密度、壓力、速度、聲速、馬赫數以及比熱比,考慮的流體為量熱完全氣體。另外,激波角(亦即激波與上游流體的夾角)用β表示,氣流穿越激波或膨脹波的偏轉角用θ表示。

1.1 流場結構

入射激波與平板邊界層發生干擾,當入射激波足夠強時,引起邊界層分離,產生分離泡,如圖2所示。

分離泡起點S上游流體發生偏轉,產生壓縮波系,壓縮波系合并為分離激波。入射激波與分離激波發生正規相交,產生透射激波。透射激波穿越自由剪切層并反射為膨脹波,使得自由剪切層向壁面偏轉,到達再附點R,R下游流體轉為與壁面平行,產生壓縮波系,壓縮波系合并為再附激波[7]。因此,已知的主要結構包括入射激波、邊界層、分離激波、干擾點膨脹波、分離泡以及再附激波。工程實際中已知的參數通常包括來流馬赫數、干擾點雷諾數以及入射激波強度。由此需要獲得各區流動參數、各干擾波的方向和強度以及壓力沿壁面分布等。

圖2 入射激波致邊界層分離流場結構Fig.2 Sketchoftheflowinducedbyanincidentshock waveboundarylayerinteractionwithseparation

1.2 壁面壓力

入射激波致邊界層分離流場的壁面壓力分布如圖3所示,分離點S上游的i點(干擾起始點)壓力即開始增加,對應干擾起始點壓力pi,通常近似認為pi等于來流壓力。伴隨分離階段,壓力歷經第一個突增,分離階段下游接著為分離流動的典型壓力平臺段。激波穿越邊界層內的外部區域,但是內部亞聲速流不能承受較大的壓力梯度,因而局部形成一個近似等壓邊界。下游伴隨再附階段,壓力經歷第二次增加,達到最終壓力峰值pF。相比分離階段壓增,再附階段壓增比較平緩[7]。

圖3 入射激波致邊界層分離流場壁面壓力Fig.3 Wallpressuredistributioninducedbyanincidentshock waveboundarylayerinteractionwithseparation

Chapmen等[21]基于大量的試驗數據,發現分離流動的一些特點并不依賴于引起分離的形式。并把這種不受下游幾何形狀影響也不受引起分離的形式的影響的干擾叫做“自由干擾”,同時對于這種干擾給出簡單的理論分析,即自由干擾理論(FIT),得出平臺壓力pP:

式中:q0為來流動壓,Cfi為壁面摩擦力系數,無量綱函數F被認為是普適函數,和馬赫數以及雷諾數無關,由試驗確定,見表1。

表1 特定點處的F值Table1 ValuesofFatspecialpoints

1.3 分離區長度

Katzer[28]針對入射激波層流邊界層干擾,進行絕熱壁面情況下的數值模擬,其模擬工況范圍為1.4≤Ma∞≤3.4,1×105≤ReL≤6×105,1.2≤pF/p∞≤3.4。基于這些數值數據,給出了入射激波層流邊界層干擾的分離區長度工程擬合公式:

式中:P＝(pF－pinc)/p∞,為未受干擾時的平板無黏激波入射點xL處邊界層位移厚度:

pinc為引起初始分離的激波強度,通過自由干擾理論給出:

CfL為未受干擾時平板邊界層xL處的壁面摩擦系數:

常數Pinc＝1.85。

Krishnan等[29]將Katzer公式的范圍拓展到Ma∞≤6.85和pF/p∞≤5.89:

本文運用數值方法并結合文獻數據進行驗證,如圖4所示。

數值計算數據基本與擬合公式吻合。由于擬合公式適用范圍所限,若馬赫數較大,其誤差將超出工程允許的范圍。

圖4 分離區長度驗證[29]Fig.4 Validation of fitting formula by Krishnan et al[29]

2 理論模型

2.1 模型簡介

基于第1節對激波致邊界層分離流場的分析,將局部流場結構進行簡化,建立理論模型,如圖5所示。

首先將分離泡假定為一個三角形。壁面壓力分布大致分為如下幾個階段:干擾上游與無干擾平板邊界層分布基本一致,分離階段經歷第一個壓增,隨后為壓力平臺,再附階段經歷第二個壓增,到達壓力峰值。由此假設第一個壓增完全由分離激波(C2)引起,分離角θS為引起分離激波(C2)的氣流偏轉角,波后(2)區壓力為平臺壓力。第二個壓增完全由再附激波(C6)引起,再附角θR為引起再附激波(C6)的氣流偏轉角,波前(6)區壓力亦為平臺壓力。入射激波(C1)與分離激波(C2)發生正規相交,產生透射激波(C3)與(C4)以及滑移線(Σ),(R4)反射為膨脹波,(6)區氣流轉平產生再附激波(C6)。該模型忽略分離激波和再附激波根部的壓縮波系以及各激波在穿越邊界層時的彎曲,將分離激波及再附激波簡化為直接從壁面引出的兩道平直斜激波。

以來流(0)區氣流狀態作為參考,運用極曲線理論描述流場結構,如圖6所示。p0/p0＝1,θ0＝0°,因此(0)區狀態點位于原點。所有由(0)區穿越斜激波可能達到的狀態由極曲線Γ0表征。(0)區氣流穿越入射激波(C1)到達(1)區,因此根據區(1)相對于(0)區的壓力比以及氣流偏轉角(此處即為楔角θw,由于氣流向下偏轉,為負值)可以確定其在極曲線Γ0上的位置。所有由(1)區穿越斜激波可能達到的狀態由極曲線Γ1表征。同時,(0)區氣流穿越分離激波(C2)到達(2)區,由于假設(2)區壓力等于平臺壓力(即p2＝pP),并且(2)區氣流必定向上偏轉,由此可以確定其在極曲線Γ0上的位置,同時給出分離角θS的大小。所有由(2)區穿越斜激波可能到達的狀態由極曲線Γ2表征。入射激波(C1)與分離激波(C2)在H點發生正規相交,產生透射激波(C3),波后(3)區;以及透射激波(C4),波后(4)區。(3)區與(4)區壓力相等且氣流平行,因此兩者狀態點位于極曲線Γ1和Γ2的交點上。緊接著,(4)區氣流穿越膨脹波到達(6)區,由于假設(6)區壓力與(2)區一樣,同為平臺壓力(即p6＝pP),因此(6)區狀態點位于經過(2)區的等壓線(圖6中點劃線)與經過(4)區的等熵膨脹極曲線Γ4的(圖6中虛線)的交點處。所有由(6)區穿越斜激波可能達到的狀態由極曲線Γ6表征。(6)區氣流穿越再附激波(C6)到達(8)區,此時氣流轉為與壁面平行(即與來流一致,θS＝0)。因此(8)區狀態點位于極曲線Γ6與豎軸的交點處。

圖6 入射激波致邊界層分離極曲線示意圖Fig.6 Polar of the flow induced by an incident shock wave boundary layer interaction with separation

2.2 模型求解

首先,分離激波波后為平臺壓力,可由自由干擾理論獲得。分離角θS和分離激波(C2)波后參數可由激波關系式獲得。

入射激波(C1)與分離激波(C2)發生異側正規相交,產生透射激波以及滑移線,如圖7所示。

圖7 入射激波與分離激波異側相交Fig.7 Regular intersection of incident shock and separation shock

已知(0)區參數以及入射激波(C1),波后參數可由斜激波關系式獲得:

(1)區氣流穿越透射激波(C3):

(2)區氣流穿越透射激波(C4):

如圖8所示,(4)區氣流穿越分離泡頂點A處的膨脹波(RA)到達(6)區,運用膨脹波關系式:

(6)區氣流穿越再附激波(C6)到達(8)區,運用斜激波關系式:

簡化模型假設膨脹波波后(6)區的壓力與(2)區一樣,等于平臺壓力,即:

聯立式(7~18)就可以完整求解整個干擾流場波系結構、各區流動參數以及各干擾點位置。

圖8 再附階段示意圖Fig.8 Schematic of reattachment

最后求解分離區長度以及幾何關系可以進一步確定分離泡個關鍵點的位置坐標。

如圖9所示,在ΔASR中:

在ΔSHA中:

從總體結構來看,如圖10所示,干擾點H的坐標為:

由式(19~25)可確定各干擾點的具體位置。

圖10 整體結構關鍵點示意Fig.10 Schematic of key points in whole flow filed

3 結果與分析

3.1 流場結構

根據第2節建立的理論模型可以預測入射激波致邊界層分離流場,此處用數值計算結果進行驗證。

考慮工況Ma∞＝2.15、θw＝3.75°、ReL＝9.91×104,對比數值計算與理論結果的局部流場結構,如圖11所示。入射激波(C1)、(C1)與分離激波(C2)干擾后的透射激波(C4)以及反射膨脹波(R)等波系位置與形狀均與數值計算結果吻合較好,同時分隔流線也與數值結果基本重合。理論預測的分離點位置ST＝0.0633525 m,與數值計算結果的分離點位置SC＝0.0610716 m相差3.7%;理論預測再附點位置RT＝0.0964525 m,與數值計算結果的再附點位置RC＝0.097169 m相差僅0.74%。由于理論模型里忽略了邊界層的逐漸變形過程,而是將分離點S以及再附點R的氣流偏轉作為突變處理,因此由分離點S以及再附點R處分別直接引出分離激波(C2)以及再附激波(C6),與數值計算結果中兩處形成壓縮波系進而合并為激波不同。各區流動參數的預測也相當準確,如表2所示。

圖11 理論預測與數值計算局部流場結構對比Fig.11 Comparison of wave structure between theoretical prediction and numerical result

表2 理論預測與數值計算結果各區流場參數對比Table 2 Comparison of flow parameters

壁面壓力分布對比如圖12所示。由于理論模型忽略了激波根部的壓縮波系,因此壁面壓力分布呈臺階狀,不似實際情況的平滑過渡,但是關鍵位置如干擾起始點、壓力平臺以及壓力峰值均吻合較好,因此整體分布趨勢一致。

圖12 理論預測與數值計算壁面壓力分布對比Fig.12 Comparison of wall pressure distribution between theoretical prediction and numerical result

3.2 參數化分析

運用該理論模型可進行參數化分析,例如分離泡高度。

考慮工況Ma∞＝2.15,θw＝3.75°和ReL＝9.91×104,根據式(3)可知無干擾情況下入射點xL處的邊界層位移厚度為9.1873×10-4m。而理論模型預測的分離泡的高度為1.0524×10-3m,相比沒有考慮干擾情況增厚的14.5%。此處考慮來流馬赫數、飛行高度以及氣流偏轉角對分離泡高度的影響,工況列于表3。分離泡高度與分離角及分離區長度直接相關,分離角越大或者分離區長度越長則分離泡高度越大。而分離角與分離區長度主要由來流參數(包括來流馬赫數與壓力或者飛行高度)和壓縮角確定。

表3 分離角參數化分析的工況Table 3 Condition for parametric analysis

對于表3中的工況1,在相同的外壓縮角以及飛行高度下,來流馬赫數的影響如圖13所示。隨著馬赫數(2≤Ma∞≤6)的增大,分離泡的高度先略微降低隨后增大。這是由于隨著馬赫數的增大,平臺壓力緩慢減小(接近正比于Ma-1/2),分離角也隨之緩慢減小;而總壓增隨著馬赫數的增大迅速增大,導致分離區長迅速增長。兩者相互制約,小Ma時,分離區長度變化較小,分離角的變化起主導因素,因此分離區高度有輕微的降低;大Ma時,分離區長度變化劇烈,起主導因素,因此分離區高度隨馬赫數的增大逐漸增大。

圖13 分離泡高度隨來流馬赫數的變化Fig.13 Influence of incoming Mach number

對于表3中的工況2,在2°~15°范圍內,在相同的來流馬赫數以及飛行高度下,分離泡高度是入射激波氣流偏轉角的增函數,如圖14所示。根據自由干擾理論,平臺壓力基本不隨外壓縮角變化,因而分離角也基本不變。而總壓增隨著外壓縮角的增大而增大,因此分離泡高度也隨之單調增大。

圖14 分離泡高度隨外壓縮角的變化Fig.14 Influence of flow deflection for incident shock

對于表3中的工況3,保持來流馬赫數以及外壓縮角不變,在20 km到50 km范圍內,分離泡高度隨飛行高度的升高而增大,在30 km高度時略有減小,如圖15所示。由自由干擾理論結合激波關系式可知,來流壓力的降低并不會改變分離角的大小,但引起初始分離的pinc隨之顯著降低,造成分離區長度增大,因而分離泡高度增大。30 km高度處壓力并未平緩降低,而是出現小的增加,因此分離區高度也略有減小。

圖15 分離泡高度隨飛行高度的變化Fig.15 Influence of flight altitude

4 結 論

本文針對激波致邊界層分離流場建立了理論預測模型。將分離泡簡化為三角形,總壓增分為分離階段以及再附階段兩部分。假設分離階段的壓增完全由分離激波完成,分離角為引起分離激波的氣流偏轉角;再附階段的壓增完全由再附激波完成,再附角為引起再附激波的氣流偏轉角。分離激波與再附激波之間的壁面壓力一直為平臺壓力,并運用極曲線方法描述了該局部流場波系結構。該模型可以給出包括分離激波、反射膨脹波以及再附激波等波系結構的位置和形狀以及各區流動參數,并得到了很好的CFD驗證。運用該理論模型,發現分離泡高度隨來流馬赫數的增大先略有減小再增大,同時是外壓縮角以及飛行高度的增函數。該模型忽略分離激波和再附激波根部的壓縮波系以及各激波在穿越邊界層時的彎曲,將分離激波及再附激波簡化為直接從壁面引出的兩道平直斜激波。

本文所建立的理論模型是對真實流動結構進行了一定簡化,后續還有以下工作可考慮開展:

(1)本文模型忽略了分離激波和再附激波根部的壓縮波系以及各激波在穿越邊界層時發生的彎曲。可以將Henderson[30]的方法添加到模型中,將邊界層分層,建立更為詳細的理論模型。

(2)本文采用的是理想氣體模型,對于高超聲速流動,伴隨高溫真實氣體效應,可以采用考慮復雜化學反應的氣體模型,其對波系結構以及流動參數均會產生影響。