側向噴流干擾流場建立與消退過程數值模擬

劉耀峰，薄靖龍

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

0 引言

側向噴流控制技術是利用發動機噴流提供的直接力，對飛行器的姿態與軌道進行控制，當舵面控制無法滿足需求時，可用于補充或替代舵面。側向噴流控制技術具有響應時間短，不受來流動壓影響等優點，正在應用到越來越多的飛行器上，如航天飛機、飛船、高機動導彈等［1］。

側向噴流干擾流場中包含復雜的波系、渦系結構、流動的分離與再附等復雜的流動現象。國內外對于這類復雜干擾現象的研究已持續了半個多世紀，研究工作包括理論分析、數值模擬與實驗等。但不論是在流動機理上還是在工程應用上仍然存在很多問題需要深入地探討。

非定常性是側向噴流干擾流動的一個重要特性，主要體現在噴流的非定常性、噴流與來流干擾流場建立及消退過程的非定常性以及分離區內激波/激波干擾、激波/剪切層干擾的非定常性等。非定常過程產生的動載荷和氣動力波動會對彈體結構以及控制精度產生影響。準確預測噴流干擾的非定常效應影響對總體及控制系統設計非常重要。

在以往的研究中，人們關注較多的是定常噴流干擾效應，其非定常干擾效應并未引起廣泛的重視，研究工作相對比較少。1997年在Arnold工程發展中心進行了瞬態噴流干擾效應的探索性實驗研究［2］，采用高速Kulite壓力傳感器測量得到瞬態表面壓力，并綜合應用加速度計和應變儀測量數據得到了瞬態力和力矩。1998年，日本科研人員［3－5］采用風洞試驗、飛行試驗和CFD計算等方法定性和定量評估導彈非定常側噴氣動干擾效應，結果表明側向噴流非定常干擾效應對彈體氣動力動態特性有較大的影響。Dash［6］在類AIT外形上采用數值方法研究了側向軌控噴流干擾流場的非定常特性，Ebrahimi［7］研究了攔截彈外形上噴流發動機開啟和關閉過程的瞬態效應，結果表明噴流開啟到穩定的過程中，力和力矩存在一個較高的峰值。劉耀峰［8］在帶有叉字形四片翼的尖拱頭部－圓柱組合體外形上采用顯式雙時間步方法計算了非定常側噴干擾流場，劉君等［9］利用多組分的N－S方程，研究了發動機非定常工作過程對導彈氣動力動態特性的影響。楊彥廣等［10］對發動機啟動和關車過程中典型時刻的流場進行了細致的刻畫，分析了流動參數和模型氣動力系數的動態變化過程，并給出了噴流壓比和攻角變化所帶來的影響。孫得川等［11］采用3階WENO方法計算了大氣層內帶有側向噴流的超聲速導彈非定常繞流流場。劉學強等［12］采用DES方法針對噴流在開啟關閉時的氣動特性進行了數值模擬。陳堅強等［13］研究了舵面運動與側向噴流間相互干擾的非定常動態相應問題。以上研究獲得了許多有價值的結論，但研究工作缺乏系統性與確定性，難以形成規律性的認識，對于噴流非定常干擾效應的機理以及如何應用于工程實踐尚需要開展進一步探索研究。

本文重點關注側向噴流干擾流場建立與消退過程中波系等復雜結構的形成和發展及氣動特性隨時間的變化情況。利用自編的適用于高超聲速側噴非定常干擾效應研究的CFD軟件［14］，并嘗試應用多重網格方法［15］加速內迭代收斂速度，數值研究了錐－柱－裙外形上高超聲速軌控側向噴流干擾流場建立及消退過程的非定常特性。獲得了詳細的瞬時噴流干擾流場結構特性，分析了法向力放大系數(Ky)、干擾力矩系數(Kmz)、法向力系數(Cy)、俯仰力矩系數(Cmz)隨時間的變化特性，并與定常結果進行了比較，結果表明高超聲速來流條件下噴流干擾流場的建立與消退過程存在強烈的非定常效應，工程設計中需要認真考慮。同時本文的研究工作可為進一步認識側向噴流非定常作用機制及控制系統設計提供基礎與參考。

1 數值模擬方法

1．1控制方程與離散方法

采用三維可壓縮層流Navier-Stokes方程為控制方程，通過有限體積方法進行離散，對流項離散采用二階精度的Roe格式［16］，黏性項離散采用中心差分格式，時間項離散使用雙時間步方法［17］，其中內迭代采用LU-SGS隱式方法［18］，并采用多重網格方法加快內迭代收斂速度［15］。根據以往研究經驗，S-A一方程湍流模型在復雜側噴干擾研究中的適用性較好，獲得的計算結果與試驗結果相符，本文湍流模型選用一方程S-A模型［19］。

計算過程中，需要對雙時間步方法的內迭代收斂過程進行判斷和控制，使得內迭代過程在保證時間精度的前提下盡快結束。合理的方法是給定一個內迭代收斂判據來確定內迭代次數。本文選用的收斂判據［20］為:

式中:x為網格單元總數，φ為流場特征變量，取壓力或密度。當Tol小于某一給定值時(如10－2)，認為內迭代過程收斂。對于有側向噴流干擾的復雜流場，內迭代過程常因達不到收斂判據而陷入死循環，本文在計算過程中在給定內迭代收斂判據的同時限定內迭代的最大迭代次數。

1．2邊界條件

1)入流邊界:來流為高超聲速，入流邊界處給定來流的靜壓、靜溫及馬赫數。

2)出流邊界:出口為超聲速時下游流場不影響上游流場，將所有參數數值外推。

3)物面邊界:采用無滑移絕熱壁條件。

4)對稱邊界:對稱面上法向速度為零，所有變量的法向梯度為零。

5)噴流邊界:直接使用噴管出口參數。

1．3初始條件

首先計算相同來流條件下無噴流定常流場，以此定常流場作為非定常流動計算的初場。

2 計算模型和計算條件

計算模型為圖1所示的錐－柱－裙外形［21］，模型總長為1461 mm，底部直徑370 mm，噴管位于模型0°子午線上，噴管中心線經過質心。圖1同時給出了計算用坐標系。計算網格采用分區對接方式，網格在壁面、噴口附近加密處理。網格單元總量為204萬，分為33塊，圖2給出了噴口附近局部網格。

圖1 幾何外形及氣動特性系數方向定義Fig．1 Model geometry and dynamic coefficient directions

圖2 噴口附近局部網格Fig．2 Local grids in the vicinity of nozzle exit

來流條件:來流馬赫數為8．0，單位雷諾數為1．54×107，其它參數詳見表1。噴流條件:本文選擇了一種簡化的發動機建壓曲線，如圖3所示，在零時刻打開噴流，不考慮噴流從啟動到穩定所需要的時間，假設噴流瞬時建立;在1．5 ms時關閉噴流，不考慮噴流完全關閉所使用的時間，認為噴流出口壓力從平臺壓力減小為零瞬時完成。這樣處理就把噴流本身的建壓過程和干擾流場的建立與消退過程分離開來。在計算過程中出口馬赫數和溫度保持不變，具體噴流參數見表1。

表1 計算參數Table 1 Computational parameters

圖3 噴管出口壓力隨時間變化曲線Fig．3 The jet pressure with time at nozzle exit

3 結果與討論

由于側向噴流干擾流場結構的復雜性，在使用雙時間步方法計算非定常側噴干擾流場時，內迭代收斂非常緩慢，往往需要數百步的計算，才能滿足對時間精度的要求，其計算量之大、計算時間之長在工程應用中是難以接受的。因此，需要采取合適的加速收斂措施，在不降低計算精度的前提下，最大程度地減小計算量和計算時間。

本文采用多重網格方法加速內迭代過程收斂速度，由于流場形成與發展過程中會出現激波系、低壓區等復雜結構，需要在傳統多重網格方法的基礎上采取改進措施以限制數值計算過程中負壓力、負密度等非物理現象的產生，提高多重網格方法在側向噴流干擾流場模擬中的魯棒性，詳見文獻［13］。

圖4 法向力系數隨時間的變化歷程Fig．4 Comparison of normal force coefficient process

本文多重網格技術應用全近似(FAS)格式，采用兩重網格V循環，圖4給出了本文算例分別采用單重網格和多重網格在不同最大內迭代次數條件下計算得到的法向力系數隨時間的變化歷程比較，從圖中可以看出單重網格不同最大內迭代次數條件下獲得的結果存在差異，多重網格最大內迭代次數取為50和75時的結果幾乎重合，多重網格最大內迭代次數為50時的結果與單重網格最大內迭代次數為600時的結果基本一致，可見采用多重網格方法可顯著地減少計算量及計算時間，表現出了良好的加速收斂效果。

3．1側噴干擾流場建立過程瞬時流場特性

圖5給出了噴流干擾流場建立過程中不同時刻壁面及對稱面壓力等值線圖，反映了干擾流場細節隨時間的變化特性。在0．05 ms時，噴流作用尚局限在噴口附近比較小的范圍內，可以觀察到在噴流前方有一道微弱的激波。0．21 ms時，噴流弓形激波的高度和強度有了明顯的增強，分離激波的影響范圍有明顯增加，裙尾已受到噴流干擾的影響。0．3 ms時，弓形激波的高度和強度繼續加大，同時噴流干擾影響到達下表面，形成了微弱的包裹效應。0．5 ms以后，噴流上游流場基本不再變化，分離區基本穩定;從流場結構上看，0．8 ms和1．5 ms時的流場基本一致。

圖5 建立過程彈體表面及對稱面等壓力云圖(M∞=8．0，α=0°)Fig．5 Pressure contours on body surface and symmetric plane at the jet startup process(M∞=8．0，α=0°)

圖6給出了不同時刻上表面對稱線上壓力分布，可以看出隨著時間的增加，噴流影響區域不斷擴展，緊靠噴流前的壓力峰值自0．05 ms以后逐漸降低，趨于穩定。緊靠噴口后的壓力降低，由于噴流的影響，裙尾前段壓力升高，后段壓力降低，噴口后的壓力分布在0．8 ms以后趨于穩定。

3．2側噴干擾流場消退過程瞬時流場特性

圖6 建立過程上表面對稱線壓力分布(M∞=8．0，α=0°)Fig．6 Pressure distribution of 0°centerline at the jet startup process(M∞=8．0，α=0°)

圖7給出了不同時刻壁面及對稱面壓力等值線圖，從圖中可以清晰地看到噴流關閉之后干擾流場結構隨時間的消退歷程。在噴流關閉之后0．05 ms內，弓形激波及相應流場結構與充分發展的有噴流流場相比無明顯的變化，在1．77 ms時，可以清晰地看到弓形激波的強度減弱、激波位置后移到噴口之后，1．98 ms時，激波越過裙尾，此時激波強度進一步減弱。2．4 ms以后，激波完全消失，流場結構與無噴定常狀態基本相同。

圖8給出了不同時刻上表面對稱線壓力分布。可以看出，隨時間推移，彈體上表面高壓區向后移動，其峰值量值不斷減小。而在1．77 ms～1．85 ms時峰值壓力有一定幅度回升，這是因為弓形激波達到裙尾，氣流在裙尾前再次壓縮。2．4 ms以后，上表面對稱線上壓力值已經接近無噴流值。

由于噴流干擾流場的建立和消退過程受到的主導因素不同，兩者花費的時間不同，噴流流場的建立過程花費了0．8 ms建立起了穩定的流場結構，噴流關閉后，流場結構的消退過程花費了約0．5 ms的時間。

3．3氣動特性隨時間的變化特性

通常采用力和力矩放大系數來表征噴流對氣動特性的影響，對于本文由于軌控發動機位于彈體質心，噴流本身產生的力矩理論上為0，傳統定義的力矩放大系數失效，本文選擇使用干擾力矩系數(Kmz)來表征軌控噴流對干擾力矩的影響，其物理含義是噴流干擾引起的噴流推力作用點位置的改變，負值代表前移，正值代表后移。法向力放大系數、干擾力矩系數定義見式(2)和式(3):

式中:下標jeton、jetoff分別表示有噴、無噴時總的氣動力/力矩系數，下標jet表示噴流本身產生的氣動力/力矩系數。

圖9和圖10分別給出了法向力系數(Cy)和俯仰力矩系數(Cmz)隨時間的變化曲線以及與有噴定常值和無噴定常值的比較，其中氣動力/力矩系數不包括發動機本身的貢獻。可以看出:在噴流啟動后，兩者均存在一定幅度波動，隨后趨于有噴定常值。在噴流關閉之后，法向力系數經歷小幅震蕩之后逐漸趨于無噴定常值，俯仰力矩系數經震蕩幅度大，由于弓形激波后退的影響，會產生較大的抬頭力矩。

圖7 消退過程彈體表面及對稱面等壓力云圖(M∞=8．0，α=0°)Fig．7 Pressure contours on body surface and symmetric plane at the jet shutdown process(M∞=8．0，α=0°)

圖11給出了噴流開啟過程中法向力放大系數隨時間的變化規律。從圖中可以看出:在噴流啟動后，法向力放大系數先增加后減小，而后再增大趨于定常解。這一變化趨勢可以從第3．1節干擾流場的建立過程中得到解釋，在0．21 ms時，由于噴流上游產生了很強的弓形激波和高壓區，而噴流下游的低壓區很小，從而出現了較大的附加法向力，而后下游的低壓區域繼續增大，附加干擾力逐漸減小，0．3 ms以后下游低壓區壓力逐漸升高，附加干擾力增大，因而法向力放大系數出現減小再增大的現象。圖中同時給出了定常收斂值(噴管出口壓力取Pj=32．21P∞，出口馬赫數與出口溫度與表1相同)和非定常時均值，可以看出非定常時均值與定常值存在差別。

圖8 消退過程上表面對稱線壓力分布(M∞=8．0，α=0°)Fig．8 Pressure distribution of 0°centerline at the jet shutdown process(M∞=8．0，α=0°)

圖9 法向力系數隨時間變化特性(M∞=8．0，α=0°)Fig．9 Variation of normal force coefficient with time(M∞=8．0，α=0°)

圖10 俯仰力矩系數隨時間變化特性(M∞=8．0，α=0°)Fig．10 Variation of pitching moment coefficient with time(M∞=8．0，α=0°)

圖12給出了干擾力矩系數隨時間的變化規律，并與定常值及非定常時均值進行了比較。可以看出，在噴流開啟后干擾力矩系數正負值交替出現，說明在噴流開啟后噴流推力作用點不穩定，時而前移時而后移，在0．3 ms時距離質心位置最遠。比較定常值和非定常時均值，發現兩者存在差別。

圖11 法向力放大系數隨時間變化特性(M∞=8．0，α=0°)Fig．11 Variation of normal force amplification factors with time(M∞=8．0，α=0°)

圖12 干擾力矩系數隨時間變化特性(M∞=8．0，α=0°)Fig．12 Variation of interference moment coefficient with time(M∞=8．0，α=0°)

4 結論

利用基于三維雷諾平均N－S方程、有限體積離散及雙時間步方法的側噴干擾非定常計算方法，數值研究了錐－柱－裙外形上高超聲速軌控側向噴流干擾流場建立及消退過程的非定常效應，給出了詳細的噴流瞬時干擾流場結構特性。分析了法向力系數、俯仰力矩系數、法向力放大系數及干擾力矩系數隨時間的變化規律，并與定常值及非定常時均值進行了比較，通過本文研究，可得到如下結論:

1)高超聲速來流條件下噴流干擾流場的建立與消退過程具有強烈的非定常效應，瞬態值和時均值及定常值差別明顯。

2)噴流干擾流場建立過程中，噴流開啟到穩定狀態之間法向力存在一個較高的峰值。

3)噴流完全關閉之后，在較長一段時間內流場中還存在噴流干擾力的殘余影響，俯仰力矩系數變化劇烈。

本文給出了一定條件下軌控噴流干擾流場建立與消退過程的非定常特性及引起的劇烈的氣動力變化特性，并得到了一些有意義的結論，但側噴干擾效應非常復雜，其建立與消退過程的非定常效應與研究外形、來流條件、噴流條件、發動機布局等密切相關，不同條件下的干擾效應差異很大，為深入了解這些條件的影響，還需要開展更為詳實的研究。

［1］ 李素循．近空間飛行器的氣動復合控制原理及研究進展［J］．力學進展，2009，39(6):740－755．［Li Su-xun．Progress in aerodynamics of combination control for vehicles at high speed［J］．Advances in Mechanics，2009，39(6):740－755．］

［2］Chamberlain R，McClure D，Dang A．CFD analysis of lateral jet interaction phenomena for the THAAD interceptor［C］．AIAA 2000－0963，The 38th AIAA Aerospace Sciences Meeting and Exhibit，Reno，USA，January 10－13，2000．

［3］Shinichiro T，Etsuroh S，Akihiko Y，et al．Sidejet aerodynamics interaction effect of the missile，Part 1-Estimation of missile sidejet interaction force by modeling in pressure field［C］．AIAA 1998－4273，The AIAA Atmospheric Flight Mechanics Conference and Exhibit，Boston，August 10－12，1998．

［4］Shinichiro T，Etsuroh S，Akihiko Y，et al．Sidejet aerodynamics interaction effect of the missile，Part 2:Prediction of interaction effect by the force measurement［C］．AIAA 1998－4346，The AIAA Atmospheric Flight Mechanics Conference and Exhibit，Boston，August 10－12，1998．

［5］Shinichiro T，Etsuroh S，Akihiko Y，et al．Sidejet aerodynamics interaction effect of the missile，Part 3－Flight test results［C］．AIAA 1998－4347，The AIAA Atmospheric Flight Mechanics Conference and Exhibit，Boston，August 10－12，1998．

［6］Dash S M．Lateral jet aerodynamic interaction simulations for dynamic pressure loads［C］．AIAA 2000－2036，The 6th AIAA/CEAS Aeroacoustics Conference，Lahaina，Hawaii，June 12－14，2000．

［7］Ebrahimi H B．Numerical investigation of jet interaction in a supersonic freestream［C］．AIAA 2005－4866，The 17th AIAA Computational Fluid Dynamics Conference，Toronto，Canada，June 6－9，2005．

［8］ 劉耀峰．CFD軟件研發和戰術導彈噴流/繞流復雜流動數值模擬［D］．北京:北京理工大學，2005．［Liu Yao-feng．Development of CFD software and complex flows numerical simulation for interaction flow field between external flow and jet flow of tactical missiles［D］ ．Beijing:Beijing Institute of Technology，2005．］

［9］ 劉君，楊彥廣．帶有橫噴控制的導彈非定常流場數值模擬［J］．空氣動力學學報，2005，23(1):25－28．［Liu Jun，Yang Yan-guang．Numerical simulation of lateral jet control induced by impulse rocket motor for a supersonic missile［J］．Acta Aerodynamica Sinica，2005，23(1):25－28．］

［10］ 楊彥廣，劉君．高超聲速主流中側向噴流非定常干擾效應研究［J］．空氣動力學學報，2004，22(3):295－301．［Yang Yan-guang，Liu Jun．Unsteady characteristic research of lateral jet in hypersonic external flow［J］．Acta Aerodynamica Sinica，2004，22(3):295－301．］

［11］ 孫得川，賈曉洪．帶側向噴流的導彈非定常流場模擬［J］．固體火箭技術，2006，29(1):25－27．［Sun De-chuan，Jia Xiao-hong．Numerical simulation on nonsteady flow-field around missile with lateral jet［J］．Journal of Solid Rocket Technology，2006，29(1):25－27．］

［12］ 劉學強，李青，柴建忠，等．用DES方法數值模擬噴流開啟/關閉時的流場［J］．宇航學報，2007，28(6):1520－1523．［Liu Xue-qiang，Li Qing，Chai Jian-zhong，et al．The numerical simulaiton for lateral jet flow during the jet opening/closing using DESmethods［J］．Journal of Astronautics，2007，28(6):1520－1523．］

［13］ 陳堅強，陳琦，謝昱飛，等．側向噴流與舵面運動相互干擾的數值模擬研究［J］．宇航學報，2014，35(5):515－520．［Chen Jian-qiang，Chen Qi，Xie Yu-fei，et al．Numerical study on the interaction of lateral jet and rudder movement［J］．Journal of Astronautics，2014，35(5):515－520．］

［14］ 薄靖龍．側向控制噴流非定常干擾效應研究［D］．北京:中國航天空氣動力技術研究院，2013．［Bo Jing-long．Investigation of unsteady lateral control jet interaction effects［D］．Beijing:China Academy of Aerospace Aeronautics，2013．］

［15］Jameson A．Solution of Euler equations by a multigrid method［J］．Applied Mathematics and Computions，1983，13:327－356．

［16］Roe P L．Approximate Riemann solvers，parameter vectors and different schemes［J］．Journal of Computational Physics，1981，43:357－372．

［17］Jameson A．Time dependent calculations using multigrid with applications to unsteady flow past airfoils and wings［C］．AIAA 1991－1596，The 17th AIAA Computational Flutid Dynamics Conference，Honolulu，Hawaii，June 24－26，1991．

［18］Yoon S，Jameson A．Lower-upper Symmetric-Gauss-Sediel method for the Euler and Navier-Stokers equations［J］．AIAA Journal，1988，26(9):1025－1026．

［19］Spalart P R，Allmaras SR．A one-equation turbulence model for aerodynamic flows［C］．AIAA 1992－0439，The 30th AIAA Aerospace Sciences Meeting and Exhibit，Reno，NV，January 6－9，1992．

［20］ 袁先旭，張涵信，謝昱飛，等．非定常數值模擬方法的發展及其在動態繞流中的應用［J］．空氣動力學學報，2004，22(4):432－437．［Yuan Xian-xu，Zhang Han-xin，Xie Yu-fei，et al．The development of an unsteady numerical methods and its application on dynamic vehicle flows［J］．Acta Aerodynamica Sinica，2004，22(4):432－437．］

［21］John M S，Bernard J J．Acoustic associated with fluid/structure coupling for interceptor missile fow fields［C］．AIAA 2000－2035，The 6th AIAA/CEASAeroacoustics Conference，Lahaina，Hawaii，June 12－14，2000．