石英燈陣模擬非均勻氣動加熱的功率優化

朱言旦，魏東，劉深深，曾磊，＊，桂業偉

1.中國空氣動力研究與發展中心 空氣動力學國家重點實驗室，綿陽 621000 2.中國空氣動力研究與發展中心 計算空氣動力研究所，綿陽 621000

結構熱試驗是對飛行器相關材料和結構進行長時間熱響應考核的重要試驗方法［1］，而石英燈陣輻射加熱以其熱慣性小、便于控制、方便調整結構、對復雜結構適應性強等特點在高超聲速飛行器結構熱試驗中獲得廣泛應用［2］。

對于利用石英燈陣模擬均勻氣動加熱，目前國內外已進行了大量的研究，發展了包括理論解析方法［3－5］、蒙特卡羅方法［6］及其相關簡化方法［7］在內的數值模擬方法，并開展了相關試驗研究。對于利用石英燈陣模擬大面積的非均勻氣動加熱，目前一般做法是將試驗區域進行離散分區［8－9］，分區內仍舊利用石英燈陣模擬均勻熱流，分區之間利用隔板等隔斷措施來避免分區之間的干擾以達到各分區之間非均勻熱流分布的效果。雖然分區模擬的方法可以很大程度上簡化控制，并已 廣 泛 應 用 于 材 料［10］、結 構［11－12］甚 至 是 整機［13］的結構熱考核試驗中，但有限的分區導致石英燈陣得到的非均勻熱流分布呈區塊性非均勻分布，分區邊緣處熱流存在間斷，能夠達到的模擬加熱過程與實際飛行條件下加熱過程有較大差異，給試驗過程帶來較大誤差。

同時航空航天技術的發展對飛行器綜合一體化熱／能量管理系統提出越來越高的需求，這就要求在結構考核試驗中同時實現表面熱流加載和進入飛行器內部熱量的精確控制，而傳統的分區模擬通常是極限考核，越來越難以滿足精細化模擬的需求。為了利用石英燈陣得到整體連續性、光順性、局部分布特征更好的非均勻熱流分布，朱言旦等提出了一種通過異化調控各石英燈功率來實現非均勻氣動加熱的非分區模擬方法［14］，該方法基于遺傳算法通過反向設計優化石英燈陣中各燈的功率實現非均勻氣動加熱的石英燈陣非分區整場模擬，并對大面積二維非均勻氣動加熱進行了模擬，初步驗證了該非分區模擬非均勻氣動加熱方法的可行性，但對于梯度較大的位置熱流模擬結果偏差較大。

利用石英燈陣進行非均勻氣動加熱非分區模擬試驗時，通常需要根據試驗表面的溫度實時計算燈陣各燈功率，而基于遺傳算法的功率優化算法很難滿足實時性的需求，為進一步探索并改進非分區模擬非均勻氣動加熱的方法以滿足快速設計需求，并且對大梯度熱流模擬提出可能的改進方向，本文建立了對石英燈陣功率進行異化調控的線性分析方法，對石英燈陣直接模擬非均勻氣動加熱的優化問題進行了研究。

1 理論方法

1.1 石英燈陣熱流分布理論解析方法及驗證

求解石英燈陣熱流分布的方法有理論解析方法（THE）、蒙特卡羅方法（MCM）及其相關簡化方法等，文獻［15］對不同方法進行了對比分析，結果表明，理論解析方法和蒙特卡羅方法均可較好地進行簡單燈陣（平直石英燈＋平面反射屏）的優化設計，但理論解析方法具有實現簡單、計算效率高的優點。因此，為簡化計算，本文采用理論解析方法進行計算分析。理論解析方法［3］將單個石英燈簡化為一個表面均勻的黑體，通過理論推導獲得單個石英燈熱流分布。以單個石英燈中心在平板測試面上的投影位置為原點，石英燈軸向為x軸，垂向為y軸建立坐標系，如圖1所示，則單個石英燈在平板接收面上的輻照熱流可以表示為

式中：qs（x，y）為平板上坐標為（x，y）處的輻照熱流；Q0為單個石英燈發射功率；H 為單個石英燈中心距平板的距離；L為單個石英燈有效加熱長度，此處取為燈絲發熱長度。

石英燈陣熱流分布由單燈熱流分布疊加獲得。文獻［14］指出，在石英燈陣高度和間距較小時，線性疊加和直接模擬結果會有明顯差異。本文的研究對象為如圖1所示的由10支石英燈均勻并排排列組成的石英燈陣，石英燈有效加熱長度為200mm，外徑為10mm，間隔為20mm，石英燈陣高度為20mm。對于一維排列的石英燈陣，沿石英燈軸向的熱流分布不可控，故以石英燈垂直方向的熱流分布表征石英燈陣的熱流分布。每支石英燈加載功率為1 000W，蒙特卡羅方法中光束模擬量為1 000W－1。疊加計算結果與蒙特卡羅方法直接模擬結果對比如圖2所示。理論解析方法疊加得到的石英燈陣熱流分布結果與蒙特卡羅直接模擬結果符合良好，中間較為均勻的區域偏差在2.5%以內，而通常蒙特卡羅方法計算結果與試驗結果偏差在5%以內，因此可以利用理論解析方法進行石英燈陣功率優化配置分析。

Fig.1 Quartz lamp array

Fig.2 Comparison of THE and MCM computational results

1.2 石英燈陣功率優化的快速線性分析方法

石英燈陣功率優化問題本質上是一個已知試驗面上的熱流分布求解石英燈陣各燈功率的輻射反問題。由式（1）分析可知，基于理論解析方法的石英燈陣熱流與功率之間存在線性關系，文獻［16］的計算和試驗結果也獲得了同樣的規律。因此，對于由n支石英燈組成的石英燈陣，如果共有m個熱流測點，則可以列出如下方程組：

式中：xj為第j支石英燈的功率；aijxj為第j支石英燈在第i個測點上的熱流；bi為第i個測點上的熱流。

將式（1）寫成矩陣形式可得

式中：

為便于分析和求解，在式（6）兩邊分別左乘AT得

式（11）中，如果m≥n且系數矩陣B是滿秩的，則線性方程組有唯一解，可通過求解式（11）獲得石英燈陣各燈功率。

2 結果分析與討論

對于圖1所示的石英燈陣，給定每支燈的功率為1 000W，在y軸上±100mm范圍內均勻設置11個測點，基于理論解析方法得到了該功率下的熱流分布及11個測點的熱流，對于采用理論解析方法計算得到的測點熱流值分別疊加±0%、±2.5%、±5.0%、±7.5%、±10.0%以內的隨機偏差，由疊加偏差后的熱流分布通過快速線性分析方法計算得到了各燈的功率配置（P1～P10），如表1所示?？梢钥闯?，隨著熱流疊加偏差的增加，功率計算結果與預定功率偏差明顯增加。以不同狀態下的功率計算結果為加載功率，計算獲得石英燈陣的熱流分布如圖3所示。可以看出，隨著疊加偏差的增加，熱流分布結果與初始狀態的偏差逐漸增加，但偏差量在疊加的隨機偏差范圍之內。

表1 不同熱流偏差下的功率計算結果Table 1 Power computational results with different heat flux errors

在數值分析中，將矩陣的范數和其逆矩陣的范數之積稱為這個矩陣的條件數，條件數表示了矩陣計算對于誤差的敏感性，條件數越大矩陣計算對誤差越敏感，小的輸入誤差會帶來越大的結果偏差。計算得到不同測點范圍和測點數目時系數矩陣的條件數（2范數），如圖4所示。不同測點范圍和測點數目時，系數矩陣的條件數一直保持在較大的水平，說明石英燈陣功率優化問題的功率優化結果對熱流的輸入比較敏感，如果缺少有效的抑制手段，可能會對石英燈功率的裕度提出較高要求。對于相同的測點范圍，不同測點數目時條件數變化較大，但隨著測點數目的增加呈現收斂的趨勢。對于相同的測點數目，不同測點范圍時條件數變化較大。整體上看，當測點范圍基本等于燈陣尺度時，條件數相對較小且隨測點數目的增加變化幅度不大。綜上所述，測點范圍和測點數目對條件數影響較大，因此可以通過優化分析選擇合適的測點范圍和測點數量得到相對穩定的系數矩陣。

Fig.3 Optimized heat flux results with different heat flux errors

為驗證線性分析方法（LM）進行石英燈陣功率配置優化得到給定熱流分布的有效性，選定了3個典型熱流狀態進行分析，并同時采用遺傳算法（GA）進行結果對比。其中，Case 1以實現100kW／m2的均勻熱流為目標，Case 2以實現20～90kW／m2的線性變化熱流為目標，Case 3以實現30kW／m2階躍到80kW／m2的階躍熱流為目標。燈陣熱流分布優化設計的遺傳算法流程［14］如圖5所示。文中遺傳算法采用64位二進制編碼，交叉及變異概率分別取為0.8和0.005，目標函數如式（14）所示，各參考點權重相同。選取種群規模為500，優勢群體10%，進化到5 000代時停止。

Fig.4 Variation of condition number with test point number

Fig.5 Process of optimization design method for quartz lamp array heat flux distribution［14］

式中：x為長度為n的優化參數向量，此處為燈陣中各燈功率，約束條件為xi∈ ［0，Pmax］，Pmax為石英燈允許加載的最大功率；m為測點數目；wi為測點i的權重；q′i為測點i的目標熱流；qi為測點i的計算熱流，qi是x的函數。

在配有Intel i7處理器的計算機進行單機計算，采用線性分析方法求解一個狀態耗時約7.4×10－6s，采用遺傳算法求解一個狀態耗時約296.3s。不同狀態下各石英燈功率計算結果如表2所示?？梢钥闯?，采用線性分析方法與遺傳算法求解的結果基本趨勢具有一致性。Case 3中線性分析方法的結果個別石英燈的功率為0W，對應位置遺傳算法得到的功率值也比較小，這說明該位置的石英燈可能是不必要的，或者說該位置的石英燈間距可以適當增大?？紤]石英燈功率、數量和位置的綜合優化問題已超出線性分析方法的范圍，本文暫不予進行更多討論，但這也給后續更為復雜熱流的模擬優化提供了一條具有參考性的方向。

將表2中的優化功率分別作為石英燈陣的功率進行計算獲得不同狀態下的石英燈陣熱流分布，如圖6所示，圖中refer表示優化目標熱流。對于Case 1，線性分析方法與遺傳算法求解的熱流優化結果基本一致，中間區域基本模擬了100kW／m2的均勻熱流，兩端有較大的波動偏差。從圖2可以看出，石英燈陣兩端熱流衰減迅速，為了使兩端熱流盡量滿足100kW／m2的均勻熱流，兩端的石英燈優化功率會偏大，又由于石英燈陣的不連續性，熱流分布本身會有波動形起伏，這就導致優化結果兩端存在較大的波動偏差。而且，應該注意到的是，遺傳算法優化得到的熱流分布具有明顯的不對稱性，說明其并沒有收斂到全局最優解，這也可以在一定程度上反映出該優化問題對輸入的敏感性。對于Case 2，線性分析方法與遺傳算法求解的熱流優化結果基本一致，中間區域基本模擬了20～90kW／m2的線性變化熱流，兩端有較大的波動偏差，與Case 1的規律基本一致。對于Case 3，線性分析方法與遺傳算法的熱流優化結果和規律較為一致，基本得到了高低熱流平臺，但都有較大的波動，中間過渡區域比預定過渡趨勢緩和很多，說明利用石英燈陣模擬大梯度的熱流分布具有一定的難度。

綜上所述，對于選定的3個典型化狀態，線性分析方法與遺傳算法求解的熱流優化結果基本一致，一定程度上線性分析方法求解的熱流結果與比遺傳算法求解的熱流結果還要更符合目標熱流，兩種方法計算結果的中間區域均與目標結果符合較好，兩端則有較大偏差，石英燈陣在模擬大梯度熱流分布方面具有一定的困難。

表2 各燈功率計算結果對比（測點范圍：±100mm）Table 2 Comparison of calculation results for each lamp power（Test area：±100mm）

Fig.6 Comparison of optimized heat flux distributions（test area：±100mm）

由上面的分析可以知道，石英燈陣熱流分布的可控區域為中間部分區域，因此，下面將縮小測點分布范圍到－75～75mm，保持其他條件不變，進行計算分析。不同狀態下各石英燈功率計算結果如表3所示。縮小測點分布范圍后，線性分析方法與遺傳算法的結果基本趨勢是一樣的，但個別石英燈的功率依然有較大的差別。

從熱流結果來看（圖7），由于減小了測點分布范圍，即測點基本分布在熱流可控的區域，優化熱流結果較之前有了較大的改善，遺傳算法優化結果較線性分析方法優化結果波動較大，說明對于此類問題，常規的遺傳算法很難收斂到全局最優解，可以考慮采用多種群遺傳算法等全局尋優能力更強的策略。對于Case 3，優化結果在高低熱流平臺的優化結果有較大改進，但在大梯度區域依然存在較大的波動偏差，說明石英燈陣模擬大梯度熱流分布有其本身的缺陷，這跟熱輻射傳輸的物理規律是符合的。

表3 各燈功率計算結果對比（測點范圍：±75mm）Table 3 Comparison of calculation results for each lamp power（Test area：±75mm）

圖7 優化后熱流分布對比（測點范圍：±75mm）Fig.7 Comparison of optimized heat flux distributions（test area：±75mm）

3 結 論

1）建立的石英燈陣功率優化的快速線性分析方法基本可以對石英燈陣功率配置進行快速優化。對于比較平緩的熱流分布可以獲得較好的熱流模擬結果，對于大梯度熱流分布也能獲得與遺傳算法較為一致的熱流模擬結果和規律，具有較好的工程適用性。

2）模擬非均勻氣動加熱的石英燈陣功率優化問題的功率優化結果對熱流的輸入比較敏感。熱流測點的分布范圍和測點數量對問題系數矩陣的條件數影響較大，但可以通過優化選擇合適的測點范圍和測點數量得到相對穩定的系數矩陣。

3）對于大梯度熱流分布以及更為復雜的熱流分布，可能需要結合石英燈陣的石英燈數量和位置進行綜合優化，后續將開展相關研究以更好地適應工程需要。