返回艙高雷諾數再入過程底部流動穩定性1)

朱德華 沈清 楊武兵

(中國航天空氣動力技術研究院,北京 100074)

引言

飛船返回艙、再入彈頭、高速導彈等都存在飛行穩定性問題.早在20 世紀60 年代研制阿波羅返回艙時人們就認識到底部流動穩定性是一個影響飛行穩定性的關鍵因素.底部流動失穩會產生非定常振蕩的氣動力,會造成鈍體飛行器的非定常的振動,嚴重情形下會導致飛行失穩.高速鈍體繞流底部流動失穩機理的認知將會對飛行器的穩定性控制提供基礎理論支撐.

一直以來,由于傳統的數值模擬方法數值耗散過大,往往不能獲得底部流動失穩特征,而風洞試驗方面試驗模型需要支撐機構,模型支撐對底部流動影響很大,造成底部流動失真.因此人們在底部流動失穩現象觀察與機理認識等方面一直存在理論認識上的不足和偏差.

國外一直在嘗試通過風洞試驗的方法研究高速底部流動穩定性.在1968 年Moseley 等[1]就針對阿波羅返回艙的飛行穩定性問題進行了研究,認識到了底部流動不穩定是導致飛行不穩定的主要因素.他們通過在背風面添加肋板或穩定翼的方法來消除底部流動的不穩定行為,增強了返回艙的飛行穩定性.

近期國外在鈍體繞流底部附近的剪切層失穩方面取得了研究進展.美國NASA 的Danehy 等[2]在2006 年發現返回艙肩部剪切層會存在層流、轉捩以及湍流的形態特征,并且會對底部流動產生影響.

2014 年,Combs 等[3]再次針對返回艙外形進行了研究,其發現返回艙剪切層最不穩定的狀態是出現在較小的攻角下(0?和12?).在攻角0?和12?情形下剪切層為轉捩狀態,而在攻角24?情形下剪切層是層流狀態.Danehy 等[2]的研究由于試驗模型底部存在支撐機構,底部流動整體情況以及失穩現象尚未獲得.

針對鈍體外形的底部流動的試驗研究還有一些類似的代表性工作如文獻[4-5],他們發展的試驗測量技術使得認識底部流動失穩現象成為可能.

2009 年,美國CUBRC 的MacLean 等[6]針對返回艙外形進行了層流與湍流狀態的數值計算和試驗研究.試驗是在LENS-I,LENS-II 激波風洞進行的,其在返回艙表面布置了大量的測壓和測熱傳感器,希望用定量測量方式研究底部流動失穩.數值計算方法采用的是US3D 代碼中的DES 方法.數值計算也證實了模型支架(圓形尾支撐、方形尾支撐以及圓形側支撐)同樣會對底部流動結構產生很大影響,無法正確模擬底部流動失穩.

2011 年,Brock 等[7]利用RANS 及RANS/LES 方法對高速返回艙的底部流動進行了數值模擬,獲得了其底部流動的層流到湍流的數值結果.但是由于研究對象底部存在支撐,底部流動失穩特征被破壞.

針對鈍體外形底部流動的數值模擬研究還有一些代表性的工作如文獻[8-12],研究表明采用高精度的數值模擬方法可以獲得底部流場的精細結構,有助于對底部流動失穩機理的認知.

2012 年,朱德華等[13]針對高速阿波羅返回艙底部流動穩定性開展了研究,研究發現了在低雷諾數的情形阿波羅返回艙底部流動結構存在不穩定性.

在鈍體繞流的國內研究方面,2019 年涂佳黃等[14-15]在低速低雷諾數情形下的圓柱和方柱繞流的尾流干擾研究方面取得了一定的進展,但在高速高雷諾數情形下的鈍體繞流的尾流相關研究較少.

綜上所述,高速高雷諾數鈍體繞流底部流動穩定性的研究無論是在風洞試驗還是數值模擬方面均存在一定的困難,其底部流動失穩現象以及失穩機理仍然存在理論認識上的不足和偏差.隨著高精度數值模擬方法的進步以及大規模并行計算能力的提升,使得高速高雷諾數鈍體繞流底部流動現象認識和失穩機理的研究成為了可能.

本文將采用多區并行隱式大渦模擬方法針對返回艙類鈍體外形,在高雷諾數再入條件下,研究其底部流動現象以及失穩機理,為返回艙類鈍體外形外部擾動因素對底部流動的作用機理及其穩定性控制提供基礎理論支撐.

1 控制方程和數值方法

控制方程為曲線坐標系下三維可壓縮Navier-Stokes 方程[16],采用的是基于有限差分的高精度多區并行隱式大渦模擬方法求解該方程.對于空間離散,無黏項采用五階WENO 格式[17],黏性項采用六階中心差分格式.時間離散采用三階TVD Runge-Kutta 方法[18].

壁面處采用等溫無滑移邊界條件,由于坐標變換在體軸上出現奇性軸,軸上點的流動參數由奇性軸周圍點平均獲得,遠場為無反射邊界條件.

需要注意的是奇性軸附近點的計算需要利用奇性軸上的值,而奇性軸上的通量為零,無法進行特征分裂,本文采取的處理方法為將五階WENO 格式在軸方向降階為二階守恒型NND 格式[19],另外兩個方向保持五階WENO 格式.

2 數值模擬結果分析

計算條件選取參照的是Brock 等[7]的試驗條件.具體的計算條件為(Re以1 m 為特征長度)

計算模型選取的是縮比的類獵戶座返回艙外形,具體的尺寸如圖1 所示.

圖1 計算模型尺寸圖(mm)Fig.1 Sketch of Orion mode(mm)

高速鈍體繞流底部流動物理現象的認識可以為數值模擬研究提供正確的指導,因此本文利用已知的鈍體繞流研究經驗,初步給出了基于類獵戶座返回艙外形高速高雷諾數繞流的基本流動形態(圖2),為后續數值模擬研究提供一定的理論支撐.

圖2 返回艙外形基本繞流形態Fig.2 Base flow patterns of return capsule

基于上述高速高雷諾數鈍體繞流的基本流動形態,選取結構網格開展高精度數值模擬工作.類返回艙鈍體外形具有軸對稱特性,C 型網格劃分是比較容易實現的,但在高馬赫數,高雷諾數條件下其繞流流動結構復雜,需要在復雜流動現象區域進行局部網格加密,單一的C 型網格很難實現局部加密,因此本文選取C 型網格和H 型網格混合的方法來實現網格的局部加密.為了保證邊界層的模擬精度,計算網格在壁面處第一層網格尺度為10μm,=0.8,針對類獵戶座返回艙采用兩套網格進行數值模擬,一套為粗網格網格量約為1800 萬,一套為細網格網格量約為3300 萬(主要在流向和法向進行加密,展向網格點數均為121 個).針對典型流動結構進行局部加密的類獵戶座返回艙網格示意圖如圖3 所示.

圖3 類獵戶座網格示意圖Fig.3 Sketch of grid

2.1 鈍體繞流底部流動的網格效應

本節針對類獵戶座返回艙外形不同網格下的大渦模擬結果進行分析,揭示其底部流動穩定性的網格效應.

圖4 顯示的是粗細兩套網格同一時刻的0.1U∞速度云圖以及瞬時流線.

圖4 返回艙外形繞流流線分布圖Fig.4 Skin friction line around the return capsule

圖4 返回艙外形繞流流線分布圖(續)Fig.4 Skin friction line around the return capsule(continued)

從圖4 可以看出二者分離起始位置基本一致,肩部剪切層起始位置也基本一致,再附(鞍點S)位置隨著網格的加密向返回艙底部靠近,底部物面處的半鞍點S′與再附鞍點S的連線是不穩定的,在剪切失穩的作用下已經出現了振蕩特征.

圖5 顯示的是同一無量綱時間t=1.3 時刻數值模擬結果.其中返回艙中心對稱面是用密度梯度顯示的,尾跡發展區以及遠尾跡區的渦結構是用速度梯度第二不變量Q值來描述的,這里Q=3000.

圖5 返回艙外形繞流瞬時形態Fig.5 Instantaneous flow patterns around the capsule

從圖5 可以看出,細網格肩部剪切層失穩靠前,細網格發卡渦結構已經出現在了尾跡發展區,而粗網格發卡渦結構出現在遠尾跡區.發卡渦結構的出現和破碎行為預示著轉捩的發生,可以預測繼續加密網格會使發卡渦結構的出現略微提前但不會移動到分離區域,因此繼續加密網格來捕捉尾跡發展區的轉捩行為已非必要.當然遠尾跡區的湍流行為加密網格還是十分必要的,再壓縮波會受到遠尾跡區的湍流影響,從圖5(a)和圖5(b)的再壓縮波角度的對比中可以看出由于細網格轉捩更早,流動更早的向湍流發展因此其再壓縮波的角度更小.

最后給出類獵戶座返回艙外形不同網格分布下頭部對稱線上無量綱壓力(p=ˉp/ρu2∞)的分布對比曲線圖,如圖6 所示.從圖中可以看出細網格(grid2)相比粗網格(grid1)激波捕捉能力略有提升,駐點處的無量綱壓力值稍有增加,細網格下駐點壓力為0.935,無黏球頭繞流駐點壓力為0.93[20](其中M∞=6.41狀態為高精度差值結果).可以看出本文提出的奇性軸處理方法在駐點上的計算結果是精確的.

圖6 頭部對稱線上的無量綱壓力分布Fig.6 Pressure distribution on the symmetrical axis

本文的研究只限于底部以及尾跡發展區的流動穩定性,出于計算量考慮針對返回艙外形繞流研究采用的是細網格.后續將針對遠尾跡區繼續加密網格開展研究.

2.2 類獵戶座返回艙外形底部流動失穩特征

前文借鑒高速鈍體繞流底部流動研究經驗初步給出了類獵戶座返回艙外形基本流動形態,下面分析數值模擬結果與其進行相互驗證.

頭部區域主要出現的物理現象是激波和駐點.圖7 中返回艙表面顯示的是壓力云圖,可以看出,物面壓力分布呈現軸對稱形態,并且隨著計算時間推進頭部對稱性保持很好,證實了本文采用的奇性軸處理方式的合理性.

圖7 返回艙外形繞流瞬時流動結構及監測點分布Fig.7 Instantaneous flow patterns and monitor point around the return capsule

底部區域主要出現的物理現象是流動分離,在雷諾數超過臨界值后本質上是結構不穩定,其在擾動的存在下極易出現非對稱結構失穩,出現周期性以及非周期性的結構振蕩特征,并且流場結構失穩與流動失穩存在相互影響.

肩部區域主要出現的物理現象是剪切層.圖7中的中心對稱面是采用密度梯度來顯示的,可以看出,高馬赫數氣流經過弓形激波后變為亞聲速氣流,亞聲速氣流經過返回艙肩部由于其外形的收縮特征導致氣流出現膨脹加速,即出現膨脹波,氣流經過膨脹波后達到超聲速狀態,此時返回艙底部由于壓力梯度的存在出現流動分離,流動分離內部氣流速度很低為亞聲速狀態.這種流動形態的存在導致返回艙肩部出現明顯的剪切區域即剪切層,剪切層本質上是不穩定的,在擾動的作用下其會出現失穩特征.從圖中可以看出,剪切層經過短暫的平直段后開始出現基頻渦,此時剪切層對流馬赫數Mc=(u1?u2)/(a1+a2)近似為0.8 左右,存在明顯的壓縮性效應.

尾跡發展區以及遠尾跡區主要出現的物理現象是轉捩、湍流以及類卡門渦街.從圖7 中速度梯度第二不變量Q值的分布可以看出肩部剪切層失穩直接導致了在尾跡發展區出現了發卡渦結構并隨著流向距離的發展迅速破碎,尾跡流動開始出現轉捩并向湍流化發展.

綜上所述,數值模擬獲得的類獵戶座返回艙繞流符合高速高雷諾數鈍體繞流的基本特征.圖7 中黑色圓圈示意位置為壓力信號監測點,在下文中用于對其基本特征形成機制開展詳細分析.

2.2.1 底部流動結構失穩

類獵戶座返回艙外形的底部流動存在大面積的流動分離,其流動拓撲結構不同時刻具有不同特征.中心對稱線附近流動的拓撲結構主要是由流動分離形成的,底部物面處的半鞍點S′與再附鞍點S的連線本質上是不穩定的,在水平和垂直擾動的作用下其會出現拉伸和振蕩特征,導致局部的拓撲結構發生變化.

圖8 顯示的是時間間隔t=1.4 ～1.7 內其局部拓撲結構的變化規律.可以看出圖8(a)中拓撲結構主要是SNS連線結構,這是中心對稱線處SS連線結構受擾動作用形成的一種拓撲結構.圖8(b)中拓撲結構主要是SS連線結構,其在垂直擾動的作用下已經偏離中心對稱線.圖8(c)中拓撲結構又重新變為了SNS連線結構,至此形成了一個周期的振蕩過程,其無量綱振蕩頻率f≈1/0.3 ≈3.3.

圖8 底部流線的瞬時拓撲形態Fig.8 Topological structure patterns of base flow

為了確認上述流動結構失穩現象,任意選取中心對稱面上底部物面半鞍點S′與再附鞍點S的連線附近一點作為壓力信號監測點(參見圖7 監測點1),其無量綱壓力隨任選無量綱時間(t=ˉt/(L/u∞),其中L為特征長度為1 m)間隔內的分布如圖9 所示.

圖9 監測點1 任選時間間隔內的壓力分布Fig.9 Pressure distribution of monitor point 1

針對此壓力信號進行傅里葉分析獲得其頻譜特征,從圖10 中可以看出其存在無量綱頻率為3.3左右的振蕩頻率(1/t),這與拓撲結構振蕩頻率基本吻合.

圖10 基于監測點1 壓力信號的頻譜特征Fig.10 Spectrum on pressure at the monitor point 1

2.2.2 肩部流動失穩

肩部剪切是高速鈍體繞流普遍存在的物理現象,其在高雷諾數情形下極易出現失穩.任意選取肩部剪切層開始出現失穩以及失穩過程中位置作為壓力信號監測點(參見圖7 監測點2 和3),其壓力隨時間的分布如圖11 和圖12 所示.

圖11 監測點2 任選時間間隔內的壓力分布Fig.11 Pressure distribution of monitor point 2

圖12 監測點3 任選時間間隔內的壓力分布Fig.12 Pressure distribution of monitor point 3

針對上述兩點的壓力信號進行傅里葉分析獲得其頻譜特征,從圖13 和圖14 中可以看出二者均存在一個無量綱頻率為3.3 左右的拓撲結構振蕩頻率,這是底部分離造成的剪切層的整體振蕩頻率,而另一占優頻率約為10 左右,這是剪切層自身失穩產生的,對應著剪切層的基頻.隨著剪切層的空間發展,更多更高頻率的激發代表著基頻渦演化出了更多尺度的渦,剪切層失穩加劇.

從底部流動結構失穩和肩部剪切失穩的分析中可以看出,二者在各自主導區域內均具有各自的特征,相互作為擾動源推動各自的失穩歷程.

圖13 基于監測點2 壓力信號的頻譜特征Fig.13 Spectrum on pressure at the monitor point 2

圖14 基于監測點3 壓力信號的頻譜特征Fig.14 Spectrum on pressure at the monitor point 3

2.2.3 尾跡區耦合失穩

在尾跡發展區可以預見上述兩種失穩模式將會產生復雜的相互作用,出現耦合失穩,導致尾跡區出現轉捩、湍流以及類卡門渦街等復雜物理現象.任意選取尾跡發展區一點作為壓力信號監測點(參見圖7監測點4)分析其頻譜特征.

從圖15 可以看出,底部結構失穩與肩部剪切失穩相互作用后各自的占優頻率仍然存在,整個頻譜呈現出寬頻特征,更多高頻能量的激發意味著轉捩的發生.另外需要注意的是兩種失穩模式作用產生了約為0.95 的無量綱占優頻率,這可能對應著類卡門渦街失穩頻率.

為了驗證存在類卡門渦街物理現象,在遠尾跡區任意選取再附激波附近一點作為壓力信號監測點(參見圖7 監測點5)分析其頻譜特征(圖16).

圖15 基于監測點4 壓力信號的頻譜特征Fig.15 Spectrum on pressure at the monitor point 4

圖16 基于監測點5 壓力信號的頻譜特征Fig.16 Spectrum on pressure at the monitor point 5

與尾跡發展區的頻譜特征曲線對比可以看出,無量綱頻率約為0.95 的占優頻率仍然存在且能量更高.為了直觀的認識此類卡門渦街振蕩現象,將不同時刻遠尾跡區中心對稱面上的渦量分布顯示在圖17 中.從圖中可以看出,高頻的小渦被低頻流動結構帶動形成類卡門渦街形式的振蕩流動現象,振蕩的頻率約為1.0 左右與頻譜分析結果基本一致.

卡門渦街的振蕩頻率公式為f=S t×(U∞/D),無量綱參數S t在高雷諾數情形下約為0.27,U∞為來流速度,D為繞流物體直徑.將本文相關參數無量綱化代入可得f=1.06,這與分析獲得的頻率基本吻合,可以認為遠尾跡區的振蕩流動現象為類卡門渦街振蕩現象.

2.3 鈍體繞流底部流動失穩機制

圖17 遠尾跡湍流區渦結構振蕩特征Fig.17 Oscillation characteristics of vortex structure in the far wake turbulent region

綜上所述,高雷諾數狀態下的類獵戶座返回艙底部流動與低雷諾數狀態下的阿波羅返回艙底部流動一樣均存在流動不穩定問題,但二者在底部流動失穩模式和流動狀態方面存在明顯差異.高速高雷諾數狀態下的返回艙類鈍體繞流的底部流動穩定性主要存在兩類失穩模式,一類是肩部剪切區域具有明顯壓縮性特征的剪切失穩模式,另一類是底部分離區域具有振蕩特征的流動結構失穩模式.這兩種模式會在尾跡發展區以及遠尾跡區出現耦合失穩,促使尾跡發展區出現轉捩行為,遠尾跡區出現湍流行為以及類卡門渦街振蕩行為.而高速低雷諾數狀態下的返回艙類鈍體繞流底部流動穩定性主要存在的是流動結構失穩模式,遠尾跡區是層流振蕩行為.

困擾人們的鈍體繞流肩部高熱流問題,底部阻力問題,鈍體飛行時的振動問題等均嚴重依賴于以上兩類失穩模式及其耦合效應的正確模擬.

本文研究指出,在模擬高速高雷諾數狀態下的鈍體繞流問題時,底部分離的起始位置的正確捕獲十分重要,這決定了鈍體外形肩部剪切層的形成.另一方面是肩部剪切層失穩行為必須正確捕捉,其決定了分離渦再附點位置.尾跡發展區的轉捩行為以及遠尾跡區的湍流行為也需要準確模擬,其影響的是再壓縮波的產生以及強度.

需要注意的是,本文的研究中采用的是隱式大渦模擬方法,針對返回艙類鈍體外形的繞流數值模擬并未加入任何人為的擾動信息,此種模擬方法可以用來模擬靜音風洞或者真實飛行條件下的鈍體繞流問題.對于噪聲風洞的鈍體繞流問題,需要在此方法的基礎上考慮擾動信息的引入,后續將開展此項研究工作.

3 結論

本文采用大渦模擬方法,針對高雷諾數返回艙類鈍體外形再入過程的底部流動穩定性開展了研究.細致的刻畫了類獵戶座返回艙外形的底部流動形態以及穩定性特征.從肩部剪切失穩、底部流動結構失穩,尾跡發展區以及遠尾跡區的耦合失穩等多個角度分析了類獵戶座返回艙外形的底部流動失穩機制.研究發現了其底部流動穩定性主要存在兩類失穩模式即肩部剪切失穩模式以及底部流動結構失穩模式,二種模式存在耦合效應;同時在遠尾跡湍流區存在類卡門渦街形式的振蕩行為.這些結論為理解返回艙外部擾動因素對底部流動的作用機理及返回艙穩定性控制提供了基礎理論支撐.

具體結論如下:

(1)高速高雷諾數狀態下的鈍體繞流底部流動穩定性均主要存在兩類失穩模式,一類是肩部剪切區域具有明顯壓縮性特征的剪切失穩模式,另一類是底部區域具有振蕩特征的流動結構失穩模式.

(2)兩類失穩模式存在耦合效應,遠尾跡區的湍流流動在耦合效應的作用下會出現類卡門渦街振蕩行為.

(3)底部流動形態以及失穩特征的正確模擬有助于解決鈍體繞流肩部高熱流,底部阻力,鈍體飛行時的振動等問題.