雙級延伸噴管高空展開過程動力學耦合仿真研究

趙博文，田維平，董新剛，宋學宇，曹濤峰

(1.中國航天科技集團有限公司四院四十一所，西安 710025；2.中國航天科技集團有限公司第四研究院，西安 710025)

0 引言

延伸噴管技術是一種運用在火箭發動機上的新型技術，它改變了傳統大擴張比噴管的單級構型設計，而將噴管的擴張段分為了2～3節，最小的一節與噴管收斂段和喉襯固接在一起，稱為基礎噴管，其余的收攏在基礎噴管的周圍，稱為延伸錐。延伸噴管的類型眾多，作動筒式雙級延伸噴管作為目前延伸噴管的主要研究對象，其通常由基礎噴管、2節延伸錐、4個作動筒、能源系統、4個同步鎖緊機構等組成。當延伸噴管接到工作指令后，能源系統開始工作，通過展開機構(4個雙級作動筒)在發動機點火前或點火后的某個瞬間將延伸錐迅速展開到位并鎖緊，使得噴管的擴張比迅速增大，從而大大增加了上面級發動機的性能[1]。目前，國外已經有許多潛射及地地戰略導彈發動機，如美國的MX導彈Ⅱ、Ⅲ級發動機及慣性上面級發動機系統(Ius)第二級發動機奧巴斯6E，前蘇聯的SS-N-18、SS-N-20導彈Ⅱ、Ⅲ級發動機，俄羅斯白楊M導彈Ⅱ、Ⅲ級發動機，日本的入軌發動機KM-D，烏克蘭的SS-24Ⅱ、Ⅲ級發動機等都采用了延伸噴管技術[2-3]。大量的仿真和試驗都已表明，在發動機長度(或噴管外露長度)受限的條件下，采用延伸噴管技術可顯著提高噴管的擴張比，從而提高發動機的比沖、推力等性能，在高空發動機及上面級發動機上極具應用價值。鑒于延伸噴管技術在宇航和戰略導彈發動機上廣闊的應用前景，其必將成為下一代固體火箭發動機的關鍵技術之一。

延伸噴管按照工作過程，可分為發動機點火前展開式和點火后展開式。其中，點火前展開式適用于冷分離發動機，多用于運載領域。而在導彈武器等領域，由于導彈總體的控制需求，發動機往往采用熱分離方式進行級間分離，即延伸噴管采用先點火、后展開的程序，這意味著延伸錐在展開過程中噴管的流場已經建立，燃氣的流動勢必會影響延伸錐的展開運動，而延伸錐的展開反過來又會造成噴管流場的改變。因此，在延伸噴管展開過程中，噴管的流場與延伸錐的運動屬于耦合影響的過程，任何從單一角度來考慮和研究該問題都是不準確的，難以反映出延伸錐的真實展開情況。而且僅靠試驗手段，一般只能得到延伸錐的宏觀運動情況，以及流場的宏觀氣動力特性和局部觀測點處的流場參數，無法獲得詳細的流動演化情況，且試驗成本往往很高。因此，有必要從流場與動力學耦合的角度來對延伸噴管展開過程進行仿真研究[4]。

目前，與延伸噴管展開動力學或流場相關研究國內外均有報道。Venable等[5]說明了氣動計算和數值仿真對延伸噴管設計和研制的重要性，并介紹了通過一種地面試車裝置“CAEPE”來研究氣動力對延伸噴管影響的實驗原理和實驗方案； Gentil P[6]利用高空模擬試驗對折疊瓣式延伸錐上的氣動載荷變化及其力學性能進行了研究，驗證了該設計方案的可行性。王成軒[7]、閻德元[8]、尤軍峰[9]等對作動筒式延伸噴管進行了相關運動和力學分析，得到了簡化條件下延伸噴管的展開動力學模型，并通過數值計算得到了延伸錐的展開運動特性；董飛等[10]利用ADAMS仿真軟件對延伸噴管展開動力學進行了仿真計算，得到了在已知條件下延伸錐展開的位移、速度和加速度等運動參數變化情況；白宏偉[11]、劉華坪[4]等對帶有延伸噴管的固體火箭發動機動態熱分離過程級間區與尾流場進行數值計算，對發動機熱防護和展開系統的驅動力設計提供了參考。韓麗霞[12]、劉加明[13]等對通過數值仿真計算對延伸噴管展開前或展開狀態下的熱環境和溫度場進行了計算分析。

上述文獻通過計算和仿真得到的關于延伸噴管展開過程中延伸錐的運動規律或尾流場變化情況往往只適用于特定條件，與延伸噴管真實工作時的情況存在一定的差距。本文通過Fluent軟件中的用戶自定義函數(UDF)功能，將描述延伸噴管展開動力學模型的微分方程組通過C語言編譯到該軟件中，通過Fluent軟件流場模塊求解氣動參數，將每個時間步長離散點處計算得到的延伸錐壁面處的壓強值通過面積分求和得到延伸錐所受氣動力大小，將該氣動力在軸向進行投影，得到延伸錐所受軸向氣動力值，并傳遞給UDF函數中的動力學微分方程組，UDF函數通過數值求解延伸噴管動力學微分方程組，得到延伸錐的位移、速度和加速度等運動參數，并將位移值返回給Fluent流場模塊作為動網格劃分的依據，Fluent通過動網格方法在每一時間步上更新網格，判斷延伸錐是否運動到指定位置處，若延伸錐未達到預定位置，則再次計算延伸錐上的氣動力大小，并返回給UDF函數，如此往復迭代，從而實現延伸錐展開過程流場與動力學的耦合仿真。通過該方法，可仿真得到在燃氣流動影響下延伸噴管展開運動的真實情況，為延伸噴管展開時間的預估以及研究延伸錐的運動變化情況提供依據。

1 計算模型和計算方法

1.1 動力學計算模型

參考論文[7，14]中建立的較為典型的關于雙級套筒式延伸噴管展開動力學模型，本文以雙級作動筒式延伸噴管的1/4模型為研究對象，對其展開過程進行動力學建模，如圖1所示。

在該延伸噴管的固定支座中心建立直角坐標系，x軸與噴管軸線平行，y軸與其垂直。假設中間作動筒相對內作動筒的位移為f(t)，展開前固定支耳中心與可動支耳中心相互間的距離皆為l，令l+f(t)=F(t)，對Ⅰ、Ⅱ級延伸錐和內作動筒、中間作動筒和外作動筒分別建立動力學方程，并根據質點系運動定理，由于整個延伸噴管在y軸與z軸的合力為零，延伸錐只沿x軸運動，因此將整個延伸噴管的動力學微分方程在x軸進行投影，即得到式(1)：

Fg0-Ff01+0.25FgpxF(t)[F2(t)-s2]-1/2

={M1+2M2+(0.25M3+0.5M4)F2(t)[F2(t)-s2]-1+

s2{M0l0+M1[F(t)+l1]+M2[2F(t)+l2]}[F2(t)-s2]-1/2F2(t)×

[1-s2F-2(t)]-1/2}F″(t)+{2(M1+2M2)s2[F2(t)-s2]-1/2F-2(t)×

[1-s2F-2(t)]-1/2-s2{M1[F(t)+l1]+M2[2F(t)+l2]}×F-4(t)[1-s2×

F-2(t)]-1-s2{M0l0+M1[F(t)+l1]+M2[2F(t)+l2]}[F2(t)-s2]-1/2×

F-3(t)·[1-s2F-2(t)]-1/2·[2+s2F-2(t)][1-s2F-2(t)]-1+(0.25M3+

0.5M4){F(t)·[F2(t)-s2]-1-F3(t)[F2(t)-s2]-2}F′2(t)

(1)

式中Fg0為氣體作用于內作動筒的軸向力；Ff01為內作動筒與中間作動筒的摩擦力；Fgpx為兩級延伸錐在流場中受到總的軸向氣動力；f(t)、f'(t)、f''(t)分別為t時刻中間作動筒相對于內作動筒的位移、速度、加速度；l為展開前固定支耳中心與可動支耳中心相互間的距離；l0、l1、l2分別為內作動筒、中間作動筒、外作動筒各自質心與其支耳中心的距離；s為中間作動筒支耳距離固定支耳的豎直高度差；R0為固定支耳中心到噴管軸線的豎直高度差；M0、M1、M2、M3、M4分別為內作動筒、中間作動筒、外作動筒、Ⅰ級延伸錐、Ⅱ級延伸錐各自的質量。

圖1 雙級延伸噴管運動模型示意圖

1.2 能源系統計算模型

作動筒式延伸噴管一般將高壓氣瓶通過管路與作動筒連接，利用氣瓶中的高壓氣體作為展開能源，從而驅動作動筒帶動延伸錐一起進行展開運動。圖2展示了一種高壓氣瓶式展開能源系統的原理圖。

圖2 高壓氣瓶式能源系統原理圖

參考文獻[15]中關于貯氣瓶放氣過程和作動筒充氣過程的建模方法，將放氣和充氣過程都近似為絕熱過程，根據熱力學中開口系能量守恒定律可分別給出描述貯氣瓶壓強pb、減壓閥流量q以及作動筒壓強p三個參數變化的微分方程：

(2)

(3)

其中

(4)

式中k、R、pb0、Tb0分別為高壓氣體的比熱容比、氣體常數、初始壓強和溫度；Vb0為貯氣瓶容積；pi、Ti分別為貯氣瓶出口(閥門入口)壓強及溫度；p1、T1為閥門出口壓強和溫度，X=0.528為臨界壓比；Cv、A分別為閥門流量系數和流通面積；V0為初始時刻閥門出口管路與作動筒內容積之和；Ain、Amid、Aout分別為內作動筒、中間作動筒和外作動筒橫截面積。

令Z1(t)=F(t)，Z2(t)=F′(t)，通過移項及合并，使得式(1)～式(4)可用式(5)所示的統一形式來表示：

(5)

通過Runge-Kutta法求解(5)式中的常微分方程組可得到Z1(t)、Z2(t)、pb(t)和p(t)等參數值，將Z1(t)、Z2(t)代入相關表達式即可求得兩級延伸錐的位移、速度及加速度等運動參數。

通過數值計算驗證了該模型的可行性，因此在本文的耦合仿真計算中，采用上述延伸噴管展開動力學模型與能源系統模型的微分方程組。

1.3 流場計算模型

1.3.1 幾何模型

由于噴管內燃氣流動參數在軸向與徑向的變化遠大于在周向的變化，為了簡化計算對該問題采用二維軸對稱模型進行流場仿真計算，帶有雙級延伸噴管的發動機尾流場二維軸對稱幾何模型如圖3所示。尾流場計算區域大小設定為徑向高度為噴管出口直徑的3倍，軸向長度為噴管出口直徑的7倍。考慮到發動機內的燃燒過程，為使得仿真計算更加符合實際，對仿真的入口條件給定如圖4所示的基礎噴管入口總壓隨時間的變化規律。

圖3 雙級延伸噴管幾何模型圖

1.3.2 網格模型

在對本研究中雙級延伸噴管展開過程尾流場進行仿真計算時，網格劃分的質量和數量很大程度影響了仿真結果的精確度和計算效率，同時延伸錐在尾流場中的運動使得仿真中需要使用動網格技術。因此，如何構造初始網格，并保證在計算過程中的網格的數量和質量，是本研究的關鍵點之一。整個計算區域的網格劃分如圖5所示，噴管內流場和來流遠場采用結構化網格。

圖4 基礎噴管入口壓強變化曲線

(a)Grid distribution in stowed condition

(b)Grid distribution in deployed condition

為了適應動網格要求，尾流場采用分區網格劃分和非結構網格，對噴管出口處、噴管壁面以及延伸錐壁面附近進行網格加密，保證壁面附近的計算精度，計算域初始網格數量為58 767個。1.3.3 求解條件與邊界條件

在本研究中，流場計算基于雷諾平均N-S方程，采用密度基隱式求解器進行求解。求解格式為二階迎風格式并采用Roe-FDS通量格式。湍流模型為k-εRNG模型。氣體的密度設定為理想氣體，定壓比熱容采用基于溫度的多項式擬合方式計算，熱傳導率采用分子動力學理論進行計算，粘性系數使用Sutherland公式進行計算。

計算條件設定為飛行高空工況，用到的邊界條件如下所示：

(1)壓力入口?？倻卦O為3751 K，總壓按圖4所示升壓曲線變化，來流方向垂直于入口邊界。

(2)壓力出口。模擬海拔為25 km高度，出口壓強為2549.22 Pa，回流方向垂直于邊界，總溫為221.5 K。

(3)壓力遠場。在尾流區的左側與上側設立壓力遠場條件，總溫總壓同壓力出口，來流馬赫數為2。

(4)固體壁面。壁面設置為非滑移壁面，絕熱條件。

(5)軸線邊界。軸線處使用對稱邊界條件。

(6)連通邊界。相鄰計算域連接處的邊界，在Fluent中設為interface邊界進行數據插值計算。

2 仿真計算結果

在流場仿真計算中將0～0.2 s模擬為發動機點火進行熱分離階段，該時段內延伸錐呈收攏狀態。隨后，延伸噴管展開機構在0.2 s時刻接到展開指令進行展開運動。為了同時滿足耦合計算精度以及動網格更新的穩定性要求，通過計算將時間步長設置為0.000 1 s，總計算步數為5000步，以氣瓶初始壓強為1.2 MPa時計算結果為例分析流場變化。

在延伸噴管展開動力學與能源系統計算部分，通過對延伸噴管流場、展開動力學以及能源系統三者的耦合計算，得到了在不同氣瓶初始壓強下延伸錐在流場中的運動情況以及延伸錐受到的軸向氣動力變化情況。同時，對比分析了耦合計算結果與對應的忽略氣動力(即只考慮作動筒壓強對延伸錐展開運動的影響)條件下的計算結果，得到了氣動力對延伸錐展開運動的影響情況。計算中，將氮氣作為展開氣源，下面給出了與貯氣瓶和閥門等相關的主要參數：Tb0=253 K，Vb0=0.002 m3，Cv=0.06，A=0.000 51 m2，V0=0.000 49 m3，Ain=0.001 96 m2，Amid=0.002 375 m2，Aout=0.002 82 m2，作動筒內初始壓強p0=0.1 MPa。

2.1 流場計算結果

在模擬發動機在海拔25 km的工作高度下，通過耦合計算得到氣瓶初始壓強為1.2 MPa條件下發動機尾流場的速度、壓力及溫度隨時間變化的云圖，分析圖6～圖8可知：

(1)發動機在點火后0.005 s以內噴管尾流就已形成，在開始階段噴管入口壓強較低且噴管擴張比較大，噴管出口處產生斜激波并形成馬赫盤，影響燃氣的流動。隨著入口壓強的迅速增加，馬赫盤位置迅速向后移動消失。在0.2 s延伸錐開始運動時，噴管入口壓力上升到約為5.2 MPa，噴管出口產生膨脹波。

(a) t=0.005 s (b) t=0.2 s

(c) t=0.4 s (d) t=0.5 s

(a) t=0.005 s (b) t=0.2 s

(c) t=0.4 s (d) t=0.5 s

(a) t=0.005 s (b) t=0.2 s

(c) t=0.4 s (d) t=0.5 s

(2)在給定的高空工況以及氣瓶初始壓強為1.2 MPa的條件下，通過耦合仿真計算得到雙級延伸噴管的展開時間為0.197 0 s，即兩級延伸錐在0.397 0 s時刻展開到預定位置，此時由于噴管擴張比的增大，燃氣在噴管內的膨脹程度增加，噴管出口膨脹波的強度減弱。

(3)從溫度云圖可看出，在延伸噴管展開過程中，發動機后封頭與基礎噴管外壁組成的空腔內的溫度出現短時劇烈波動，瞬時溫度可達到1200 ℃左右，當延伸錐展開到位后該溫度值迅速下降。這是由于高空段外界壓強較小，高溫燃氣通過基礎噴管與延伸錐之間的狹縫出現短時回流，導致該空腔內溫度劇增，因此設計時，應對基礎噴管外壁、延伸錐壁面及其展開機構進行相應的熱防護措施。

2.2 動力學與能源系統計算結果

2.2.1 耦合仿真結果

在本研究中，利用耦合仿真分別計算了氣瓶初始壓強為0.8、1.0、1.2、1.5、1.7 MPa條件下延伸錐的展開情況。由于結構設計的原因，Ⅱ級延伸錐在展開過程中，其位移、速度及加速度值皆為Ⅰ級延伸錐對應值的2倍，因此下面僅給出Ⅰ級延伸錐在流場中其運動參數隨時間的變化情況。

表1給出了在不同氣瓶初始壓強(pb0)作用下延伸錐展開到預定位置所用的時間(t)。由表1可看出，在高空工況下，當氣瓶初始壓強低于0.8 MPa時，延伸錐將無法正常展開到預定位置。當氣瓶初始壓強大于1.0 MPa以后，隨著該壓強的增加，延伸錐的展開時間逐漸減小，但減小的幅度越來越小，表明延伸錐的展開時間受氣瓶初始壓強變化的影響程度隨該壓強的增大呈下降趨勢。

圖9分別對比了在不同氣瓶初始壓強作用下Ⅰ級延伸錐在流場中展開時的軸向位移、速度及加速度隨時間的變化情況。

表1 不同氣瓶初始壓強下延伸錐的展開時間

(a) Axial displacement (b) Axial velocity (c) Axial acceleration

由圖9可看出，在作動筒驅動力和流場氣動力的共同作用下，延伸錐在流場的展開過程中其軸向速度和軸向加速度先增大、后減小，表明當延伸錐越接近預定位置，其受到流場的阻礙作用越大。隨著氣瓶初始壓強的增加，延伸錐展開到預定位置處時速度也越大，延伸錐對接處的機構受到的沖擊也將越大。因此，在保證展開時間達到要求的前提下，氣瓶初始壓強的選取應盡可能小。

圖10對比了在不同氣瓶初始壓強作用下Ⅰ、Ⅱ級延伸錐所受的總軸向氣動力隨時間的變化曲線。可看出，在高空工況下，由于燃氣射流的引射作用，延伸錐在開始階段受到的流場氣動力為正向，但其值相對較小。隨著噴管入口壓強的增大以及延伸錐位置的移動，該軸向氣動力變為負向，且其值逐漸增大。當氣瓶初始壓強改變時，延伸錐受到流場的軸向氣動力變化情況大致相同。

圖10 Ⅰ、Ⅱ級延伸錐總軸向氣動力隨時間變化曲線

圖11、圖12對比了在不同氣瓶初始壓強作用下氣瓶內壓隨時間的變化情況、作動筒內壓隨時間的變化情況。由圖11和圖12可看出，在延伸噴管展開過程中，氣瓶內壓逐漸減小，作動筒內壓先增大、后減小。改變氣瓶初始壓強后，氣瓶與作動筒內壓的變化趨勢大致相同。

圖11 氣瓶內壓隨時間變化曲線

圖12 作動筒內壓隨時間變化曲線

表2給出了在0.5 s時刻時不同氣瓶初始壓強對應的各級延伸錐受軸向氣動力大小情況。可見，在0.5 s時刻時，當延伸錐展開到預定位置，且流場入口條件基本穩定后，Ⅰ、Ⅱ級延伸錐受到流場的氣動力大小基本保持不變，該計算結果可對發動機在高空段工作時延伸噴管鎖緊裝置的設計提供參考和依據。

表2 0.5 s時刻不同氣瓶初始壓強下對應的各級延伸錐受軸向氣動力大小

2.2.2 耦合計算結果與忽略氣動力計算結果對比

當忽略延伸錐在流場受到的氣動力時，則延伸錐展開運動的驅動力只有作動筒內氣體的壓力，即方程式(1)中的左端項變為Fg0-Ff01。通過將耦合仿真與忽略氣動力作用下二者的計算結果進行比較，可得到氣動力對延伸錐展開運動的影響規律和影響大小。

表3給出了在不同氣瓶初始壓強下，耦合計算與忽略氣動力計算得到的延伸錐展開時間之間的對比情況。在給定相同的氣瓶初始壓強下，耦合仿真計算得到的延伸錐展開時間比忽略氣動力計算得到的延伸錐展開時間長，二者間最大相差10.4%。這是因為在高空工況下，延伸錐在流場中主要受到負的軸向氣動力，該氣動力對延伸錐的展開過程起到了阻礙作用。隨著氣瓶初始壓強的增大，耦合計算與忽略氣動力計算得到的延伸錐展開時間二者間的差距逐漸減小，表明氣動力對延伸錐展開時間的影響隨氣瓶壓強的增大呈下降趨勢。

圖13分別給出了在不同氣瓶初始壓強下，耦合計算得到的Ⅰ級延伸錐的位移、速度、加速度以及氣瓶與作動筒內壓隨時間變化曲線與忽略氣動力計算得到的對應的曲線之間的對比情況?？煽闯?，忽略氣動力計算結果與耦合仿真計算結果二者間的差異主要體現在延伸錐的速度和加速度上，氣動力對延伸錐展開運動的影響較大，不可忽略。隨著氣瓶初始壓強的增大，延伸錐的運動及氣瓶與作動筒內壓的變化情況大致相同。

表3 不同氣瓶初始壓強下耦合計算與忽略氣動力計算關于延伸錐展開時間的對比

(a1) pb0=1.0 MPa,displacement (a2) pb0=1.2 MPa,displacement (a3) pb0=1.5 MPa,displacement

(a4) pb0=1.7 MPa,displacement (b1) pb0=1.0 MPa,velocity (b2) pb0=1.2 MPa,velocity

(b3) pb0=1.5 MPa,velocity (b4) pb0=1.7 MPa,velocity (c1) pb0=1.0 MPa,acceleration

(c2) pb0=1.2 MPa,acceleration (c3) pb0=1.5 MPa,acceleration (c4) pb0=1.7 MPa,acceleration

(d1) pb0=1.0 MPa, pressure of bottle (d2) pb0=1.2 MPa,pressure of bottle (d3) pb0=1.5 MPa,pressure of bottle

(d4) pb0=1.7 MPa,pressure of bottle (e1) pb0=1.0 MPa,pressure of cylinder (e2) pb0=1.2 MPa,pressure of cylinder

(e3) pb0=1.5 MPa,pressure of cylinder (e4)pb0=1.7 MPa,pressure of cylinder

3 結論

(1)雙級延伸噴管在高空展開時，發動機后封頭與基礎噴管外壁組成的空腔內由于高溫燃氣的短時迴流其內部溫度出現短時的劇烈波動，瞬時溫度可達1200 ℃左右。因此，須對基礎噴管外壁面和延伸錐及其展開機構進行相應的熱防護措施。

(2)在給定的高空工況下，耦合仿真計算表明，氣瓶初始壓強越大，延伸錐展開到位所用時間就越短，但延伸錐展開時間的減小幅度隨氣瓶初始壓強的增大呈下降趨勢。

(3)在高空工況下，當采用-20 ℃的氮氣作為展開能源，且氣瓶體積為0.002 m3時，則氣瓶初始壓強至少為1.0 MPa以上，才能保證延伸錐展開到位。若氣瓶初始壓強過小，則導致展開過程中作動筒壓強過小，延伸錐難以正常展開到位，若氣瓶初始壓強過大，延伸錐展開到預定位置處時的速度和加速度就越大，延伸錐減速時，對接處受到的沖擊就越大。因此，設計時需要合理地選擇展開能源輸入參數大小。

(4)高空工況下，關于延伸噴管展開時間，耦合仿真結果與忽略氣動力計算結果之間最大相差為10.4%，同時二者關于延伸錐展開速度和加速度的計算結果存在明顯差別，表明氣動力對延伸錐的展開運動影響較大，不可忽略。本文的耦合仿真方法對研究在流場和展開動力學相互影響下延伸噴管的展開情況提供了思路和依據。