基于響應面法的喉栓式噴管型面優化設計①

成 沉，鮑福廷，劉 旸，許 昊

(西北工業大學航天學院，西安 710072)

0 引言

為提高導彈的機動性與靈活性，使得系統滿足多任務需求，作為動力裝置的固體火箭發動機要具備更強的推力調節能力。在噴管中嵌入一根喉栓，就能在不改變傳統固體火箭發動機結構形式情況下，實現推力大小的調節［1－3］。這種推力調節技術因其推力調節比大，理論上能實現無級調節而受到越來越多的關注。

噴管優化設計是喉栓式固體火箭發動機設計中的一個重要方面。由于在噴管中嵌入的喉栓會擾亂傳統噴管中的氣流路徑，使噴管擴張段產生復雜的激波系，甚至引發流動分離等現象，造成比沖損失。所以，喉栓式噴管的優化設計需綜合考慮喉栓的型面與噴管型面的相互作用，減小因喉栓插入而引起的比沖損失。

本文通過CFD模擬對喉栓式噴管的噴管型面及喉栓型面進行綜合氣動優化設計，探尋兩者的設計參數對比沖損失的影響，為喉栓式噴管的設計提供參考。由于CFD本身計算量較大，若直接進行優化必須進行大范圍的搜索，需消耗巨大的計算量。本文采用響應面法，在設計空間內對有限的設計試驗點進行CFD計算，并在此基礎上建立響應面模型，再結合序列二次規劃法(NLPQL)，對噴管型面及喉栓型面進行綜合氣動優化設計。結果表明，影響比沖效率的因素主要為噴管的設計參數，在進行喉栓式噴管的優化設計時，噴管部分完全可采用傳統噴管的優化設計方法，再針對喉栓進行匹配性優化。應用響應面法，能以相對較小的計算量進行快速優化設計，滿足工程優化需求，提高工程優化效率。

1 CFD模型

1．1 數學模型

采用可壓縮流動的二維軸對稱守恒型N-S方程。對流項采用二階迎風格式進行離散，粘性項按照中心差分格式進行離散。湍流模型采用RNG k-ε模型，增強壁面函數。

1．2 邊界條件

邊界處理是流場數值模擬的重要環節之一，當喉栓從完全打開位置向著發動機喉部運動時，喉部面積減小，燃燒室壓力升高，發動機入口質量流率隨燃燒室壓力升高而增加。因此，入口邊界條件定義為自適應的質量入口，在計算過程中讀入發動機燃燒室壓力，自動調整入口質量流率，使其滿足藥柱的燃速變化規律［3］。單位面積質量流率為

式中 Ab為藥柱截面積;Ain為噴管入口面積;ρp為推進劑密度;pc為燃燒室壓強;n為燃速指數;a為燃速系數。

發動機噴管出口面上為超音速出口，不需給定任何邊界條件，全部氣流參數二階外推即可。

2 響應面法

響應面法是試驗設計與數理統計相結合的優化方法，當實驗結果與已知參數間的函數關系為隱式時，在試驗測量、經驗公式或數值分析的基礎上，對指定的設計點集合進行連續的試驗，可設計實驗結果與參變量間的函數關系，并求得參變量的系數，最終建立響應與參變量間的函數關系，在設計空間構造測定量的全局逼近。響應面方法計算簡單，因而可僅采用簡單的代數形式，依靠目標函數本身的性質確定最優解。通過回歸模型的選擇，可擬合復雜的響應關系，具有良好的魯棒性［5］。本文用設計試驗點的CFD數值計算結果替代實驗結果生成響應面，可大大減少尋求最優解的優化計算次數。

2．1 響應面初始樣本點選取方法

幾乎所有的基于響應面的優化算法，第一步就是構造初始響應面，而初始響應面構造的首要問題就是初始樣本點的選取，科學合理的樣本點在樣本空間的分布應是均勻的。本文采用中心復合實驗設計(CCD)，生成初始試驗樣本點。CCD的主要特點是可評估因素的非線性影響。

2．2 響應面回歸模型

選擇適合的響應面回歸模型是至關重要的。這關系到響應面模型是否能真實反映目標函數與設計變量之間的關系。

常用的響應面模型有多項式(polynomial regression surrogate，PRS)、徑向基函數 (radial basis functions，RBF)、支持向量回歸(support vector regression，SVR)、多元自適應回歸樣條(multivariate adaptive regression spline，MARS)以及 Kriging等。本文選擇多項式模型及Kriging模型分別生成響應面，并進行了對比分析。

多項式模型:通常選擇一階多項式為響應面模型，但一階模型很難反映真實的響應情況;而選擇大于二階的多項式雖然有較高的擬合精度，但它由于包含較多的項，需付出較大的計算代價，尤其是在多變量情況下，擬合響應面要花費的計算時間將是無法承受的。對于氣動優化問題，相對來說二階模型形式比較靈活，對真實響應近似程度較好。本文基于CFD的計算結果，采用完全二次多項式構造響應面模型。其一般公式為

式中 k為變量個數;ε為統計誤差。

Kriging模型:Kriging方法是一種通過已知點來預測未知觀察點的一種插值方法。Kriging方法利用方差的變化來表達空間的變化，且可保證由空間分布得到的預測值的誤差最小。

由于Kriging方法生成的響應面經過所有的試驗樣本點，無法直接通過試驗樣本點與響應面計算值的擬合程度來評估響應面模型有效性。文中采取額外隨機抽取測試樣本點的方法，用測試樣本點與響應面計算值的擬合程度來更好地評估響應面模型有效性。

2．4 優化算法

本文采用序列二次規劃法在響應面上進行尋優。這種算法假設目標函數是連續可微的。基本思想是將目標函數以二階拉氏方程展開，并把約束條件線性化，使其轉化為一個二次規劃問題。二階方程通過Quasi-Newton公式得到了改進，而且加入了直線搜索，提高了算法的穩定性。

3 優化算例及結果分析

3．1 優化模型

本文優化的噴管型面采用雙圓弧法設計，而喉栓型面為頭部帶小圓弧的圓錐形，噴管及喉栓構型如圖1所示。一般設計喉栓式噴管時，要滿足噴管長度限制的要求。為了滿足發動機的內彈道要求，還需保證噴管的最大、最小等效喉部截面積(喉栓插入噴管中所形成的最小截面面積)一定。所以，本文在噴管長度和喉栓直徑一定、且保證等效喉部面積恒定的條件下，為了達到比沖損失最小目標，選擇了圖1所示的設計參數，對噴管型面及喉栓頭部型面進行了優化。

圖1 噴管及喉栓型面構型圖Fig．1 The contour of the pintle nozzle

其中，噴管設計參數為噴管收斂段與水平方向的夾角α，喉部小圓弧半徑R1，擴張段大圓弧半徑R2，噴管出口半徑Re;喉栓設計參數為頭部圓錐半角θ，頭部小圓弧半徑r。

本文以與喉栓式噴管等效喉部截面積相等的最優標準噴管(無喉栓的傳統噴管)的比沖為標準，評價喉栓式噴管的比沖損失。比沖損失定義如下:

式中 Isp為喉栓式噴管的比沖;Iss為最優標準噴管的比沖。

最優標準噴管的比沖是基于響應面法對相同的噴管設計參數，在同樣條件下，對標準噴管進行優化后計算得到的。

3．2 結果及分析

基于CCD方法，生成了45個初始試驗樣本點(DOE points)，并隨機抽取了38個測試樣本點(Verification points)進行CFD計算。分別用多項式模型及Kriging模型生成響應面。圖2和圖3分別顯示了2種模型的擬合效果。圖中橫坐標為樣本點CFD計算值，縱坐標為響應面預測值;方塊表示初始試驗樣本點，圓圈表示測試點，樣本點越接近中線，說明擬合程度越好。可看出，對于測試樣本點，2種模型的預測效果都不錯，多項式模型的均方根誤差(Root Mean Square Error)RMSE=0．009%，Kriging 模型的 RMSE=0．005%。均方根誤差RMSE(%)的定義如下:

式中 yi為CFD計算的準確值;y^i為響應面模型預測值;N為評估測試點數。

圖2 多項式模型擬合圖Fig．2 The response surface created by the PRS

圖3 Kriging模型擬合圖Fig．3 The response surface created by the Kriging

可見，2種模型的精度都足以滿足要求，而Kriging模型初始試驗樣本點全部通過響應面，多項式模型中存在偏差較大的初始試驗樣本點。相比而言，Kriging模型表現較好，較適用于喉栓式噴管的優化設計問題。本文優化時，選擇Kriging模型。用序列二次規劃法在Kriging響應面上優化的結果如表1所示。在響應面上尋優的結果，使比沖損失從9．01%降到了4．52%。對此優化點進行CFD計算，得到的比沖損失為4．52%。響應面誤差為1．6%，精度滿足要求。對優化點整化之后，再進行 CFD計算，得到最終的比沖損失為4．50%。圖4和圖5是優化前后流場的對比。可看出，優化后喉栓頭部的亞音速區明顯減小，從喉栓頭部延伸出的斜激波強度減弱。從圖6可看出，噴管出口處的徑向速度明顯減小。由于只有軸向速度對比沖有貢獻，而徑向速度越大，比沖損失越大。優化后，喉栓對流場的擾動作用減弱，流動損失自然減小了。

表1 優化結果(Kriging模型)Table 1 Optimization results(Kriging model)

圖4 優化前流場馬赫數云圖Fig．4 Mach number contour before the optimization

圖5 優化后流場馬赫數云圖Fig．5 Mach number contour after the optimization

圖6 噴管出口徑向速度Fig．6 The radial velocity of the outlet

3．3 關鍵參數分析

各設計參數的敏感性如圖7所示。其中，噴管設計參數擴張段大圓弧半徑R2、噴管出口半徑Re及喉栓設計參數頭部圓錐半角θ對比沖損失的影響較大。

圖7 參數敏感性Fig．7 The sensitivity of parameters

以優化點為基礎，再進一步針對這幾個參數進行單獨的優化，得到的響應面如圖8所示?？煽闯?，R2小于140 mm時，比沖損失較小;超過140 mm時，比沖損失會急劇增加。比沖損失隨Re的增大幾乎呈線性減小，在允許范圍內，Re越大越好。θ對比沖損失的影響相對前兩個參數較小，小于70°時，隨θ的增大，比沖損失增大的較明顯;超過70°時，比沖損失幾乎不變。

再對這幾個參數進行單獨尋優，組成新的優化點，比沖損失下降為3．63%。優化點為 α=44°，R1=3 mm，R2=120 mm，Re=26 mm ，θ=10°，r=1．3 mm。

3．4 優化方式改進

總體而言，噴管設計參數比喉栓設計參數的影響要明顯一些，而喉栓的設計參數對比沖損失的影響不是很大，影響比沖效率的主要因素還是噴管的設計參數。

圖8 關鍵參數響應面Fig．8 The response surface of the key parameters

基于上述分析，改變優化方式，先對噴管進行優化，再對喉栓進行優化。這樣試驗點從45個，減小到了34個，減小了計算量。優化結果與綜合優化結果的對比如表2所示，分開優化得到的優化結果與綜合優化的差別不大，這種優化方式能在達到同樣優化效果的前提下，減小計算量。

在進行喉栓式噴管的優化設計時，噴管部分完全可采用傳統噴管的優化設計方法，再針對喉栓進行匹配性優化。

表2 改進的優化結果Table 2 Improving optimization results

4 結論

(1)喉栓式噴管的型面優化問題需綜合考慮喉栓型面和噴管型面的相互影響，以得到最優氣動型面。

(2)對于文中的問題，Kriging模型比完全二階多項式模型擬合程度高，預測能力好，適用于喉栓式噴管的優化設計。

(3)基于響應面法進行喉栓式噴管的型面優化設計，計算量相對較小，計算結果的精度和可靠度相對較高，可快速進行氣動優化設計，具有一定的工程實用價值。

(4)影響比沖效率的因素主要為噴管的設計參數，在進行喉栓式噴管的優化設計時，噴管部分完全可采用傳統噴管的優化設計方法，再針對喉栓進行匹配性優化。

［1］ 張淑慧，胡波，孟雅桃．推力可控固體火箭發動機應用及發展［J］．固體火箭技術，2002，25(4)．

［2］ Randall Smith-kent，Hai-tien Loh，Pawel Chwalowski．Analytical contouring of pintle nozzle exit cone using computational fluid dynamics［R］．AIAA 1995-2877．

［3］ 李娟，李江．喉栓式固體火箭發動機噴管性能影響研究［J］．彈箭與制導學報，2007，27(3):154-160．

［4］ 鄒林君．基于Kriging模型的全局優化方法研究［D］．武漢:華中科技大學，2011．

［5］ 熊俊濤．基于響應面方法的氣動優化設計［D］．西安:西北工業大學，2005．

［6］ 王筱蓉，周長省，鞠玉濤，等．固體火箭發動機特型噴管的型面設計［J］．彈道學報，2008，20(4):77-80．