高超聲速滑翔飛行器駐點熱環境不確定性量化評估

吳 昊，楊 濤，豐志偉，葛健全，張青斌，黃 浩

（國防科技大學空天科學學院，湖南長沙 410073）

由于速度和突防方面的優勢，高超聲速滑翔飛行器受到了各軍事強國的廣泛關注。鄰近空間內長時間超高速飛行所帶來的嚴酷氣動加熱給飛行器防熱及結構設計帶來了諸多挑戰［1］。文獻［2］指出，飛行器長時間受氣動加熱作用后，熱變形可能會引起氣動特性的變化，從而導致飛行彈道和姿態的變化，影響飛行器的整體性能。目前的數值仿真以確定性仿真為主，即需要給出精確的來流條件、邊界條件和物理模型。而氣動加熱熱流受來流條件、幾何尺寸、邊界條件和物理模型的影響，這些不確定性因素是普遍存在的，當這些因素與實際情況存在偏差時，實際熱流就會偏離預示值。不確定性因素的來源主要可分為隨機不確定性（Aleatory Uncertainty）和認知不確定性（Epistemic Uncertainty）兩大類。前者是由物理系統及其環境固有的隨機性產生的；而后者是由于人的認識不足或者信息缺乏所造成的。氣動熱環境的不確定性量化評估對于高超聲速滑翔飛行器飛行試驗方案的確定具有重要參考意義。

Bettis和Hosder［3］對高超聲速再入流動的不確定度問題進行了分析，研究混合不確定性輸入變量對氣動熱的影響。Brian A.Lockwood等進一步研究了認知和隨機不確定性混合傳播的優化使用［4］，并發展出一種區間統計方法，結合Fay-Riddell駐點熱流公式證明其在流體仿真計算中的適用性；Lockwood在進一步研究中提出一種多變量輸入情況下省時高效的基于梯度信息的代理模型［5］，研究了理想氣體和不平衡氣體條件下高超聲速黏流中熱流的不確定性。相比于國外，國內對氣動熱不確定度量化評估的研究還處于起步階段。張偉等［6］運用基于非侵入式多項式混沌（NIPC）的方法和基于Sobol指標的方法對Apollo飛船返回艙開展了氣動熱不確定度量化分析和敏感性分析，得到了駐點和肩部熱流不確定度隨徑向距離變化的規律，為飛船返回艙氣動熱防護設計提供參考依據，證實了NIPC方法與Sobol敏感性分析方法在熱不確定性問題中的適用性，對進一步用于高超聲速飛行器的魯棒優化設計和可靠性分析提供了思路。

本文以某滑翔飛行器為研究對象開展氣動熱數值預示，選取來流速度、來流密度、來流溫度和端頭半徑4個不確定性變量，利用多項式混沌展開方法對氣動熱進行不確定度量化和敏感性分析，為飛行試驗方案的確定提供理論依據。

1 氣動熱分析方法

高超聲速滑翔飛行器氣動熱環境可分解為駐點、飛行器前緣和大面積等部位，下面分別給出具體的計算方法。

1.1 駐 點

球頭駐點熱環境預測采用壓縮性修正的Fay-Riddle公式，即

式（3）中，壓縮性系數Z與駐點壓力和溫度相關，采用壓縮性系數函數表插值計算。

1.2 飛行器前緣

飛行器前緣的熱環境預測采用后掠圓柱理論，即

式（4）中，qy為翼前緣駐點線熱流，qs為與翼前緣相同半徑的駐點熱流，λe為有效后掠角，λ為后掠角，α為攻角，指數n取1.2。

1.3 大面積

飛行器大面積的熱環境預測采用流線法，計算中考慮變熵效應。本文選用Zoby公式［7］計算熱流，該公式是基于Eckert參考焓的可壓縮平板摩擦系數關系得到的，計算結果與實驗數據和飛行數據均吻合較好。下面分別給出層流、湍流和轉捩區的熱流計算公式。

（1）層流熱流，即

式（5）中，qs為采用Fay-Riddle公式得到的駐點熱流。

（2）湍流熱流，即

式（6）中，c1，c2，c3，c4，m是速度剖面指數N的函數，分別為

在通常選用的1/7次速度分布中，N=7，但真實情況下N是變化的，由軸對稱噴管壁面實驗數據擬合出的N曲線［8］，即

式（8）中，Reθ是基于動量厚度的雷諾數，表達式為

本文根據式（8）和（9），采用迭代方法求解N與Reθ。

（3）轉捩區熱流，主要采用間歇因子方法計算轉捩區熱流，即

式（10）中，W為間歇因子，變化范圍是0＜W＜1，它與轉捩起始位置和轉捩結束位置有關，取

定義邊界層外緣當地雷諾數等于轉捩雷諾數Retrb的位置為轉捩起始點，其位置為Strb。文獻［9］總結碳-酚醛端頭的飛行試驗數據，給出了與邊界層外緣馬赫數相關聯的公式，即

轉捩結束點位置為

2 不確定性量化評估方法和敏感性分析

本文采用多項式混沌（Polynomial Chaos，PC）法作為均值隨機展開，對滑翔飛行器氣動熱環境進行不確定性量化分析。多項式混沌法基于不確定性的譜表示，相比傳統的蒙特卡洛仿真具有更高的效率。

2.1 多項式混沌不確定性量化方法

不確定性譜表示的基本思想是將響應值或隨機函數α*分解為由確定性和隨機性分量表示的形式，即

其中，αi是確定性分量，Ψi是與第i模態對應的隨機變量基函數，α*假設為獨立變量x和n維標準隨機變量ξ的函數。式（14）是一個無窮序列，實際使用時對其進行截斷處理并保留部分輸出模態作為輸出。對于p階多項式混沌展開（Polynomial Chaos Expansion，PCE）和n維隨機變量，可計算總階數Nt，即

PCE方法的核心是確定展開系數αi，通?？墒褂们秩牖蚍乔秩氲姆绞接嬎愕玫?。侵入式方法在理論上更為直接，但是對于復雜問題該過程具有耗時長、代價大等缺點。而非侵入式方法，則不需要修改原確定性模型，容易構建表示復雜模型的代理模型。

對于非侵入性多項式混沌（Nonintrusive Polynomial Chaos，NIPC），目前已發展了多種方法，其中在航空航天領域應用最廣泛的是配點NIPC法。該方法首先將隨機響應面或隨機函數近似表示為PCE；然后在隨機空間中選擇Nt個矢量，每個矢量即為一個采樣點，采用確定性模型評估這些采樣點，得到式（14）左邊；接下來，形成式（16）給出的Nt個方程的線性方程系統；最后求解隨機變量的譜模態。

需要說明的是，Nt是用于確定PCE系數的確定性采樣的最小數目。如果存在更多的采樣且線性無關，則系統為超定方程，可使用最小二乘方法求解。超過所需最小個數的采樣使用過采樣率（Oversampling Ratio，OSR）表示，定義為實際采樣與最小需要采樣的比值。一般而言，配點數量可以通過式（15）乘以OSR確定，PCE的精度依賴于配點數量。

2.2 Sobol靈敏度分析方法

Sobol靈敏度分析方法是基于方差的全局敏感性分析方法，通過計算一個或多個輸入變量對輸出方差的貢獻量，以此來評價該輸入變量在不確定性模型中的重要性。Sobol指標需要通過Sobol重積分分解（以下簡稱Sobol分解）計算［10］。

假設通過PCE方法得到多項式函數f(x)，x=(x1，x2，…，xn)，定義域為In，通過Sobol分解可將該函數分解為如下形式，即

其中，各項按照以下公式分解得到

得到式（17）后，可求得總方差，即

各階偏方差為

總方差和偏方差滿足

各階Sobol指標Si1…is的定義為

Sobol指標提供了敏感性的度量，即，每個輸入不確定變量Si的貢獻以及混合貢獻Si1，...，is。輸入參數的總組合定義為部分Sobol指標的和，即

即將所有與變量xi有關的Sobol指標相加，當n=4時，變量x1對總方差的總貢獻為

Sobol指標可用于給出每個輸入不確定性變量對總方差貢獻的相對排序，且能夠考慮輸入變量和輸出指標的非線性相關性。

3 仿真算例

以某高超聲速滑翔飛行器端頭為對象，開展氣動熱環境的不確定性量化評估和全局敏感性分析。其中，熱環境以輻射平衡熱流表示，邊界層外緣屬性按照激波后屬性計算［11］，邊界層黏性用Sutherland公式計算［12］，且壁面壓力可假設為邊界層外緣壓力。

飛行環境及端頭半徑如表1所示。分析中取4個隨機不確定變量：來流速度U∞、來流密度ρ∞、來流溫度T∞和端頭半徑Rn，由于飛行器加工過程中制造誤差和飛行環境的不確定性，來流條件的隨機波動是存在的。這些變量假設為關于某個均值的正態分布，變量分布和方差系數見表2。

表1 飛行環境及端頭半徑數值Tab.1 Aviation environment and the radius of curvature of the body

表2 隨機不確定性隨機變量分布Tab.2 Probability distribution of aleatory variables

由于該問題的不確定變量維數較小，可以使用最小二乘解來比較。對隨機不確定性在給定采樣規模下分析響應面提供累積分布概率（CDF）的上下邊界。當過采樣率取2時，能夠克服傳統配點法的不穩定性，最大程度減少各組配點對響應面的影響。圖1為過采樣率為2時熱流分布的P-box圖。相應的95%置信區間為表3所示的下邊界2.5%概率水平到上邊界97.5%概率水平。

圖1 過采樣率取2時的熱流分布P-box圖Fig.1 P-box plots for heat flow(OSR=2)

圖2 隨機變量的Sobol指標Fig.2 Sobol indices for aleatory variables

表3 熱流密度95%置信區間Tab.3 95%Confidence interval of heat flux for selected samples

4 總 結

本文以高超聲速滑翔飛行器為對象開展氣動熱環境數值預示，選取來流速度、來流密度、來流溫度和端頭半徑等4個不確定變量，采用NIPC方法對端頭半徑氣動熱進行不確定量化評估和敏感性分析，得到如下結論：

（1）利用非侵入式多項式混沌方法與大量蒙特卡洛采樣所得不確定度區間吻合，表明利用NIPC方法能夠在縮短不確定性量化評估時間的同時，得到可信的結果。

（2）給定變量不確定條件下，端頭熱流的均值大約為305.05W/cm2，不確定度不小于7.03%。

（3）采用基于Sobol指標的敏感性方法，對計算結果進行了敏感性分析，結果表明來流速度變化和來流密度變化對熱流不確定度貢獻最大，端頭半徑對熱流的不確定度有一定影響，而來流溫度變化幾乎不產生影響。在高超聲速氣動熱預測中應重點關注來流速度和來流密度的變化。