V 字形鈍前緣激波反射遲滯現象

王 軍，張志雨，劉 愿，錢戰森，李祝飛

（1.中國科學技術大學近代力學系，合肥 230027；2.中國航空工業空氣動力研究院，沈陽 110034）

吸氣式寬速域飛行器從地面起飛逐漸加速至高超聲速以及返航過程中，經歷了多個氣動環境迥異的飛行階段，給氣動外形設計提出了巨大挑戰［1］。特別是，寬速域飛行條件下的激波干擾更加復雜多變，甚至與飛行階段的歷史效應緊密相關，認識其變化規律對工程設計尤為重要。

三維內轉式進氣道［2-3］具有較好的寬速域適應能力，受到廣泛關注，這類進氣道在V 字形唇口根部位置產生的復雜激波干擾，容易導致嚴酷的氣動力/熱問題［4-6］。為了深入認識V 字形唇口部位的復雜流動，肖豐收等［7-8］最早提煉出由半徑比R/r（根部倒圓半徑R和前緣鈍化半徑r）和半擴張角β表征的V 字形鈍前緣模型，在來流馬赫數Ma∞=6條件下，發現不同幾何參數下V 字形根部主要產生異側激波規則反射（Regular reflection，RR）、馬赫反射（Mach reflection，MR）以及同側激波規則反 射（Regular reflection from the same family，sRR）3 種類型。進一步的研究表明，V 字形鈍前緣的激波反射類型對流場非定常振蕩特性［9-13］、壁面氣動力/熱載荷［13-17］、下游流場演化［18-21］以及熱防護優化設計［22-23］等都起著關鍵作用。肖豐收等［7］和蒙澤威等［14］對V 字形鈍前緣的研究表明，反射類型在不同來流馬赫數下的變化明顯，激波干擾結構在寬速域條件下將變得更加復雜。值得注意的是，當飛行器經歷不同速域時，可能產生激波反射類型雙解以及轉變遲滯現象，一旦出現雙解，極容易引發流場結構突變，并給壁面氣動載荷預測帶來不確定性。張志雨［13］通過數值模擬初步證實了存在V 字形鈍前緣激波反射遲滯現象，但對遲滯過程中反射類型轉變機制的認知，仍不明晰。

本文采用數值模擬方法，并輔以風洞實驗，通過改變Ma∞，研究V 字形鈍前緣激波反射遲滯現象，以期闡明遲滯轉變機制以及遲滯現象對壁面氣動載荷影響，為內轉式進氣道V 字形唇口設計提供有價值的參考。

1 模型與方法

1.1 V 字形鈍前緣模型

V 字形鈍前緣模型如圖1 所示，由兩側斜掠直前緣和根部倒圓區域組成，模型半擴張角β=18°，倒圓半徑R=15 mm，前緣鈍化半徑r=15 mm，以及斜掠直前緣長L=21r。其中，x方向為流向，y方向為橫向，z方向為展向，φ為V 字形根部倒圓區域的周向角。

圖1 V 字形鈍前緣模型示意圖Fig.1 Schematic of the V-shaped blunt leading edge model

1.2 風洞實驗方法

實驗在中國空氣動力研究與發展中心的FL-31 下吹式常規高超聲速風洞［13］中進行，來流馬赫數Ma∞=6，靜溫T∞=55.5 K，靜壓p∞=633.4 Pa。此外，風洞配備了高速紋影系統來捕捉激波干擾結構，紋影的拍攝速率為16 kHz，曝光時間為62.5 μs。

1.3 數值模擬方法

數值求解基于雷諾平均的三維Navier-Stokes控制方程，無黏通量采用Roe 格式進行差分分裂，對流項采用二階迎風格式離散，黏性項采用二階中心差分格式離散，使用k-ωSST 湍流模型，空氣采用量熱完全氣體假設，分子黏性系數由Sutherland公式計算。計算域、邊界條件以及駐點附近網格分布如圖2 所示，鑒于來流和幾何的對稱性，為提升計算效率，采用1/4 計算域進行模擬，來流采用壓力遠場條件、壓力出口、對稱邊界（y=0 和z=0 平面）和無滑移絕熱固壁邊界條件。為了便于同本文風洞實驗結果對照，來流條件參照FL-31 風洞參數進行設置。

圖2 計算域、邊界條件和駐點附近網格Fig.2 Computational domain, boundary conditions and surface mesh near the stagnation point

計算域采用結構化六面體網格離散，并在近壁面區域進行了網格加密。以殘差下降4 個數量級，或者殘差不再變化同時駐點位置流場參數保持穩定作為收斂判據。為了進行網格無關性驗證，如圖3 所示，采用了4 套不同疏密程度的網格，提取壁面中心線上的壓力進行對比，其中壓力p采用相同來流條件下正激波波后的總壓p0進行無量綱化?？梢钥闯?，采用第3 套和第4 套網格獲得的壓力結果幾乎重合，最大的p/p0相差不超過0.4%，表明本文采用的數值模擬方法能夠滿足需求。因此，采用第3 套網格（ζ × ξ × η=450×300×90）進行后續研究，總網格量約為1 200 萬個，壁面網格第一層高度為1 μm，保持壁面y+＜1。上述模擬方法，在筆者前期的研究中已得到廣泛地驗證［8，16，18-21］，在本文2.1 節將進一步結合風洞實驗考核該方法的可靠性。

圖3 壁面中心線壓力對比Fig.3 Comparison of surface pressure along the centerline

在計算過程中，首先，分別在Ma∞=5.7 和6.7條件下，對應求解得到RR 和MR 的穩定流場。接著，在收斂解基礎上，改變0.1 個馬赫數繼續計算，直至激波干擾類型發生轉變，進而得到反射類型轉變臨界Ma∞。本文中，p∞和T∞始終保持和FL-31風洞來流參數一致，Ma∞逐步增大過程為5.7→6.5，Ma∞逐步減小過程為6.7→5.9。研究發現，遲滯轉變區間約為Ma∞=5.9～6.5，第2 節將重點對此區間中V 字形鈍前緣激波反射雙解現象以及轉變遲滯過程進行介紹。

2 結果與討論

針對R/r=1、β=18°的V 字形鈍前緣構型，首先介紹增大和減小Ma∞過程中激波反射類型演變過程；然后，分析導致激波反射遲滯現象的根源，并基于流動結構特征建立轉變邊界；最后，評估遲滯現象對壁面壓力特性影響。

2.1 改變Ma∞引起的激波反射遲滯現象

圖4 和圖5 分別給出了遲滯區間中Ma∞=5.9和6.5 時，V 字形鈍前緣根部反射類型對應為RR和MR 的三維流場。圖6 給出了Ma∞在增大和減小過程中x-z對稱面流場的馬赫數云圖，并疊加了壓力等值線，為了揭示遲滯現象對壁面壓力分布影響，對應補充了壁面無量綱p/p0云圖。

初始Ma∞=5.9 時，如圖4（a）和圖4（b）所示，直前緣脫體激波（Detached shock，DS）在V 字形根部倒圓區域相交，發生了同側激波規則反射RR，在DS 相交點（Intersection point，IP）產生兩道透射激波（Transmitted shock，TS）。如圖4（b）和圖4（c）中流線所示，由于V 字形根部倒圓區對氣流的幾何約束，來流通過激波偏折并匯聚到倒圓區根部中心位置，氣流滯止導致壓力顯著提升，沿直前緣方向產生較大的逆壓梯度，誘導駐點附近氣流向兩側偏折，形成大范圍的流動分離，分離區內產生大尺度渦和分離激波（Separation shock，SS）。兩側的SS 在駐點上游與TS 發生規則反射后繼續相交于點IP′。以Ma∞=5.9 的流場為初場并逐漸增大Ma∞，如圖6（a～d）所示，當Ma∞=5.9→6.4 時，激波反射類型一直為RR，DS 的脫體距離隨Ma∞增大而減小，TS 長度逐漸縮短，交點IP 向下游流場移動并不斷靠近IP′。

圖4 Ma∞= 5.9 時RR 反射類型三維流場Fig.4 Three-dimensional flow field of RR at Ma∞= 5.9

當Ma∞=6.5 時，原有的RR 平衡狀態遭到破壞，激波干擾結構通過調整后重新達到穩定格局。如 圖5（a）和 圖5（b）所 示，DS 與 馬 赫 桿（Mach stem，MS）相交，發生了馬赫反射MR，由三波點（Triple point，TP）發出的TS 以及剪切層Σ，入射到直前緣與倒圓區域相連接的位置。由于TS 前后存在較大的逆壓梯度，導致出現小范圍的流動分離，并產生SS。SS 與TS 直接相交，形成規則反射。從圖5（b）和圖5（c）中流線可以看出，包裹在剪切層中的氣流由兩側向中心匯聚，并在駐點附近對撞后向上游偏折，使得MS 下游形成大尺度的反轉渦對（Counter-rotating vortex pair，CVP），MS 后的流動分為兩部分，一部分為對撞導致的逆流，另一部分是直接穿過MS 的氣流，兩部分氣流最后都沿橫向溢流。以Ma∞=6.5 的流場為初場并逐漸減小Ma∞，如圖6（e～h）所示，當Ma∞=6.5→6.0 時，激波反射結構一直為MR，Ma∞減小直接導致DS脫體距離增大，MS 長度逐漸減小，TP 向下游移動，兩側的分離區范圍逐漸增大。當Ma∞=5.9 時，如圖6（a）所示，MS 突然消失，激波反射類型由MR轉變為RR。

圖5 Ma∞= 6.5 時MR 反射類型三維流場Fig.5 Three-dimensional flow field of MR at Ma∞= 6.5

圖6 改變Ma∞引起的激波反射遲滯現象Fig.6 Shock reflection hysteresis phenomena by variation in Ma∞

圖7 給出了來流Ma∞=6 時FL-31 風洞的實驗紋影和數值紋影。由于風洞實驗中存在擾動等因素，在實驗中也觀察到V 字形鈍前緣激波反射存在雙解現象，相應的兩種波系干擾結構均與數值模擬吻合良好。

圖7 風洞實驗和數值模擬中的激波反射雙解現象Fig.7 Dual-solution phenomena in the wind tunnel experiments and numerical simulations

前文的流場分析已經表明，當來流Ma∞=5.9～6.5 時，Ma∞增大或減小過程中，來流條件和幾何構型相同，激波反射類型卻分別對應為RR 和MR，其雙解特性與經典二維理論認識的差異，以及對應的轉變機制如何，將在2.2 節進行分析。

2.2 RR-MR 遲滯轉變機制

由于x-z對稱面流場具有準二維特征，并且能夠體現激波結構的主要特征，利用二維理論可以近似分析該平面內的流場參數變化。圖8 給出了來流Ma∞=6 條件V 字形鈍前緣x-z對稱面主激波結構對應的激波極曲線，考慮到流動對稱特性，圖中只展示了其中一側分支。從圖8 可以看出，來流首先經過激波DS 壓縮之后，對應由狀態（1）點轉變為狀態（2）點，當RR 時，氣流通過TS 后對應于激波極曲線上的狀態（3，4）點，當MR 時，氣流通過TS 后對應于激波極曲線上的狀態（3′，4′）點。

結合圖8 中的反射類型示意圖，理論上可以將MR 分為3 類［33］：當狀態（3′，4′）的凈偏折角為正時，對應為直接馬赫反射（Direct-Mach reflection，DiMR）；當狀態（3′，4′）的凈偏折角為零時，對應為固定馬赫反射（Stationary-Mach reflection，StMR）；當狀態（3′，4′）的凈偏折角為負時，對應為逆馬赫反射（Inverse-Mach reflection，InMR），在前人的研究中InMR 被認為是定常激波反射中的反常結構［34］。經典二維激波反射類型轉變理論［30］認為雙解存在于von Neumann 準則θvN和脫體準則θD之間，雙解中 的MR 通 常 為DiMR 或StMR。然 而，V 字 形 鈍前緣的雙解位于von Neumann 轉變邊界之下，特別地，如圖6（e～h）和圖7 所示，V 字形鈍前緣的MR 結構中，由三波點TP 發出的滑移線Σ 呈擴張狀，狀態（3′，4′）點的氣流的凈偏折角為負，對應出現了InMR 結構。實際上，本文構型中InMR 中的MS 之所以能穩定存在，主要是依靠其下游的反轉渦對CVP 維持。

圖8 Ma∞=6 時雙解區激波極曲線Fig.8 Shock polar diagram for the double solution at Ma∞=6

V 字形鈍前緣雙解區內的InMR 反射類型以及遲滯轉變特性，都與經典二維理論認識存在顯著差異，有必要探究其激波反射遲滯現象的轉變機制。圖9 給出了激波反射類型分別為RR 和MR 時的x-z對稱面的流場結構示意圖。2.1 節研究表明，流場中主激波結構上的特征位置（如IP、IP′、TP 和CVP）與激波干擾類型的轉變密切相關，對這些特征位置的定量分析，有助于理解遲滯轉變機制。為了便于下文對這些特征位置的分析，將特征點到駐點S 的水平距離定義為參數d。圖10 給出了特征點位置參數d隨Ma∞變化，其中，散點符號表示數值結果數據，曲線為理論分析或由散點通過二次多項式擬合得到的結果。

圖9 不同激波反射類型示意圖Fig.9 Sketches of different shock reflection structures

圖10 遲滯過程中不同反射類型特征點位置變化Fig.10 Positions of characteristic points of different shock reflections during the hysteresis process

當激波反射類型為RR 時，如圖9（a）所示，兩側DS 的交于點IP，兩側SS 與TS 分別干擾后進一步相交于點IP′。如圖10 所示，當增大Ma∞時，DS脫體距離δ逐漸減小，導致點IP 向下游移動，即點IP 到 駐 點S 的 距 離dIP，S逐 漸 減 ??；點IP′的 位 置 主要受大分離區范圍以及SS 影響，數值模擬結果顯示大分離區范圍以及SS 隨Ma∞變化并不明顯，如式（1）和圖10 所示，可以將點IP′到駐點S 的距離dIP′，S視 為 定 值。當 點IP 和 點IP′重 合 時，即 滿 足dIP，S=dIP′，S，波系結構將不再穩定，將會引起激波干擾類型從RR 向MR 轉變。因此，確定點IP 的位置，是準確預測RR 向MR 轉變的關鍵。

由幾何關系可知，點IP 到駐點S 的距離dIP，S主要與DS 的脫體距離δ有關，dIP，S和δ滿足式（2）關系。為了求解δ，可以將斜掠直前緣看作一段后掠圓柱，DS 將沿其流向逐步發展，當直前緣長度L足夠長時，DS 最終達到充分發展狀態，即激波的脫體高度δ到達定值δf。Zhang 等［32］發現δ主要與垂直于后掠圓柱方向的氣流分量Ma∞sinβ有關，如式（3）所示，根據無黏圓柱脫體激波高度近似理論［35］，可以求解得到充分發展狀態的脫體距離δf。由于本文采用的構型，直前緣長度限制，δ未達到充分發展狀態，但基于δf變化規律，根據數值模擬結果進行修正，可以近似獲得δ，如式（4）和圖10 所示。聯立式（1～4），進而可獲得式（5）表示的轉變邊界fRR→MR（Ma∞）。在圖10 中，dIP，S和dIP′，S變化曲線在Ma∞= 6.48 時相交，從理論上確定了臨界轉變邊界fRR→MR位置，而數值計算得到的fRR→MR位置約為Ma∞= 6.5，兩者偏差僅約0.3%。

當激波反射類型為MR 時，如圖9（b）所示，DS、TS 和MS 交于TP，此時，匯聚到V 字形根部氣流對撞后向上游偏折，進而產生CVP。當減小Ma∞時，δ逐漸增大，TP 逐漸向下游移動，同時MS長度減小。由圖6 可知，MS 長度減小至一定值后突然消失，反射類型由InMR 轉變至RR。實際上，InMR 中MS 產生或消失都與CVP 存在很強的關聯。本文采用對稱面上逆流最高位置（鞍點，見圖5（c）至駐點S 的距離dCVP，S表征CVP 尺度變化，采用TP 至駐點S 的距離dTP，S表征MS 變化，如圖10所示，dTP，S和dCVP，S變化趨勢基本一致，這進一步印證了MS 主要受CVP 尺度變化影響。

結合圖10 中dCVP，S、dTP，S和dIP，S變化過程，不難推斷出，InMR 轉變至RR 的充分條件為，CVP 尺度減 小 至 不 再 影 響DS 直 接 相 交，即dCVP，S=dIP，S，此時MS 將 消 失，反 射 類 型 轉 變。InMR 中dCVP，S和dTP，S均難以理論獲得，從數值模擬結果中，可以擬合獲得dCVP，S和dTP，S隨Ma∞的變化規律，如式（6）和 圖10 所 示，其 中，dTP，S可 以 近 似 看 作dCVP，S與CVP 導致的MS 脫體距離之和。聯立式（2～4）和式（6），進而可獲得式（7）表示的轉變邊界fMR→RR（Ma∞）。圖10 中，dCVP，S和dIP，S變化曲線在Ma∞=5.85 時相交，從理論上確定了轉變邊界fMR→RR位置；而數值計算得到的fMR→RR位置約為Ma∞=5.9，兩者偏差僅約0.8%。

以上分析表明，遲滯現象將V 字形鈍前緣激波干擾問題變得更加復雜，反射類型雙解的出現不僅使流場結構、流場參數明顯改變，而且不同反射類型對應產生的透射激波、剪切層和超聲速射流等與壁面相互作用，將會導致其氣動載荷特性出現極大的差異。

2.3 遲滯現象對壁面壓力特性影響

圖11 給出了不同Ma∞條件下RR 時V 字形根部倒圓段壁面中心線上p/p0分布和局部壓力峰值附近流場，結合圖6（a～d）可以看出，復雜激波干擾導致其壁面壓力分布極不均勻。當φ=0°～12°時，流場上、下兩側的脫體激波DS 之間，以及分離激波SS 之間的相互干擾，使得經過多道激波壓縮后的氣流從兩側向流場中心附近匯集，超聲速氣流直接沖擊駐點，形成滯止激波BS，導致駐點位置處出現極其嚴酷的中心壓力峰值（Central peak）。當Ma∞= 6.4 時，BS 形態受到上游波系影響發生變化，導致中心壓力峰值分布一定程度改變。當φ=12°～72°時，由于處于大分離區位置，壁面p/p0相對較低。

圖11 RR 反射類型壁面壓力分布Fig.11 Surface pressure along the centerline of RR

圖12 給出了不同Ma∞條件下MR 時壁面中心線上p/p0分布和局部壓力峰值附近流場。結合圖6（e～h）可以看出，由TP 發出的TS 和Σ 入射倒圓段的壁面，TS 與近壁面附近的SS 發生規則反射，產生透射激波，在下游氣流再附，產生再附激波。因此，從倒圓區域的起始位置處（φ=72°）開始，p/p0開始迅速升高，在φ=36°～42°位置的p/p0出現最外側峰值（Outermost peak）。進一步地，包裹在剪切層中的氣流由兩側向中心匯聚，在駐點附近發生對撞，在φ=0°位置的p/p0出現中心峰值（Central peak）。

圖12 MR 反射類型壁面壓力分布Fig.12 Surface pressure along the centerline of MR

圖13 給出了不同Ma∞條件下RR 和MR 時的壁面不同位置壓力峰值變化。當Ma∞增大或減小時，中線峰值在Ma∞=5.9～6.5 范圍出現遲滯環。當RR 時，中心峰值隨著Ma∞升高而逐漸增大，在Ma∞=6.4 時有所減小；當MR 時，中心峰值和最外側峰值都隨著Ma∞降低而逐漸減小。由于不同激波反射類型對應的壓力峰值產生機制存在差異，RR 工況對應中心壓力峰值約為MR 工況對應中心壓力峰值的2～3 倍。此外，根據前期將壁面條件設置為等溫壁的研究結果［13，16-17］，RR 工況對應的最大熱流值也遠高于MR 工況對應的最大熱流值，更容易造成燒蝕破壞。

圖13 遲滯過程中壓力峰值變化Fig.13 The maximum surface pressure during the hysteresis process

V 字形鈍前緣激波反射的雙解現象，將會導致流場結構、壁面載荷的幅值及其位置發生突變。在實際工程設計中，可以通過調整V 字形唇口幾何參數，規避設計馬赫數附近的遲滯現象，避免對飛行器造成的不利影響。

3 結 論

本文采用數值模擬并輔以風洞實驗，在Ma∞=5.7～6.7 范圍，研究了R/r=1，β=18°的V 字形鈍前緣激波反射類型的演變，主要得到了以下結論：

（1）隨著Ma∞的逐步增大或減小，在Ma∞=5.9～6.5 區間，V 字形后掠前緣上的脫體激波DS能夠產生規則反射RR 和馬赫反射MR 兩種類型，并出現轉變遲滯現象。進一步地，Ma∞=6 的風洞實驗證實激波反射類型存在RR 和MR 雙解。

（2）基于對流場特征結構的認識，建立了V 字形鈍前緣隨Ma∞變化時RR?MR 的轉變邊界。以RR 為初場，逐漸增大Ma∞，當DS 的交點向下游移動至與分離激波SS 的交點重合時，RR 轉變為MR；以MR 為初場，逐漸減小Ma∞，當反轉渦對的尺度減小至允許DS 直接相交時，MR 轉變為RR。

（3）激波反射類型雙解的出現使得流場結構、氣動載荷峰值及其分布特性存在突變。在雙解區內，RR 工況的壓力最大值出現在駐點（φ=0°）附近，MR 工況的壓力最大值出現在倒圓段兩側（φ=36°～42°），RR 工況的壁面壓力最大值為MR 工況的2～3 倍。

鑒于V 字形鈍前緣三維流動的復雜性，未來需要將激波干擾遲滯轉變準則拓展至寬馬赫數和寬幾何參數范圍，并借助細致的風洞實驗進行深入驗證。