多網格下解析法計算石英燈輻射熱流對比研究

關鍵詞：多網格； 解析法； 輻射熱流； 有限元法； 石英燈

中圖分類號：V216.4 文獻標識碼：A DOI：10.19452/j.issn1007-5453.2024.10.008

石英燈加熱裝置是高速飛行器熱結構地面熱強度試驗[1]中廣泛采用的加熱裝置。加熱裝置設計的合理性直接關系到結構表面熱載荷施加的準確性，可以通過計算石英燈陣列在結構表面產生的輻射熱流分布情況來評估熱載荷的均勻性、峰值水平以及其他分布特征等。

常用的石英燈陣列輻射熱流計算方法包括解析法、蒙特卡羅法，以及采用基于上述方法的有限元軟件工具等。采用有限元軟件進行輻射熱流計算的缺點主要有：（1）需要對大量的石英燈元件進行實體建模，需要設定受熱結構表面與石英燈的復雜輻射換熱關系，需要明確受熱結構和石英燈的相關材料參數；（2）對于大規模石英燈陣列的輻射換熱通常計算時間較長；（3）部分有限元軟件中限制輻射換熱單元總數，在進行大規模加熱陣列計算時，對石英燈燈絲長度方向和環向的網格劃分數量往往只能采取較為粗略的網格劃分，在一定程度上影響計算結果的精度。孔凡金等[2]對比了區域法、離散坐標法、蒙特卡羅法等幾種經典的輻射傳熱計算方法以及常用的輻射傳熱計算軟件，討論了幾種輻射熱流場計算模型的適用性，認為可根據計算目的和精度要求選擇合適的計算方法，以實現對石英燈陣列熱流場的快速準確模擬。

蒙特卡羅法通過光子軌跡追蹤能夠有效模擬光子在石英玻璃管內、反射屏作用下的復雜折、反射過程，因此計算精度較高，在輻射熱流計算中應用十分廣泛，國外很早便建立了采用蒙特卡羅法預測石英燈輻射熱流的計算方法[3-5]，被國內學者廣泛引用。宋偉浩等[6]利用蒙特卡羅法對紅外燈陣出射光線進行了隨機抽樣和追跡，建立了受熱平面上光線落點坐標數據庫，計算了受熱平面輻射照度分布情況。夏吝時等[7]采用蒙特卡羅法對平板形石英燈加熱器中水冷反光板面積、水冷反光板與燈陣間距離、熱源疏密程度、熱源陣列與材料受熱面間距等因素對輻射熱場中典型隔熱材料受熱面溫度分布均勻性和熱流密度進行了模擬計算。朱言旦[8-9]基于蒙特卡羅法發展了石英燈陣輻射熱流分布預測方法及計算程序，以燈陣中各燈功率為優化參數的石英燈陣熱流模擬優化設計方法，對某飛行器結構部件迎風面氣動加熱進行了模擬。這些采用蒙特卡羅法的研究工作需要的輸入參數較多、計算結果精度直接依賴于樣本數量大小、計算時間相對較長等不足。

解析法基于斯忒藩-玻耳茲曼定律、蘭貝特定律等傳熱學基本定律，求解時往往需要對研究對象進行簡化，使之滿足基本定律的應用條件，優點是顯式計算過程收斂性好，通過模型簡化可以實現較快的計算速度，同時可以根據需求進行編程擴展，缺點則是不能考慮光線在傳播過程中復雜的折、反射現象。宋偉浩等[10]建立了紅外燈陣輻照度分布數學模型和受熱平面溫度分布數學模型，針對燈間距離參數對受熱平面溫度分布進行了仿真與試驗。李昕昕等[11]針對解析法不適用于帶反射屏的石英燈輻射熱流計算、蒙特卡羅法計算時間過長的問題，將蒙特卡羅法與解析法相結合，提出了基于解析法的“等效面光源”求解構想，首先采用蒙特卡羅法獲得反射屏作用下的典型熱流場分布結果，再利用輻射傳熱定律建立關于面光源各離散點溫度的線性方程組，求解得到等效面光源的溫度分布，之后采用解析法對等效面光源作用下的熱流場進行計算，實現了帶反射屏的石英燈加熱裝置的熱流快速計算。張肖肖等[12]針對平板形、圓柱形、圓錐形石英燈加熱陣列，提出了基于解析法的熱流快速計算方法，構建了解析傳熱模型，建立了圓柱和圓錐加熱陣列中單根石英燈管照射范圍的判斷依據，得到了典型燈管排布下熱流分布規律。這種計算方法將燈絲半圓柱面作為平面處理，模型得到簡化，計算速度顯著提高，但在精度方面有所下降，需要進一步探討半圓柱面多網格劃分時的解析計算方法。

本文在參考文獻[12]研究工作的基礎上，建立適用于半圓柱面多網格劃分時的熱流密度解析計算方法，分析不同網格數量劃分下熱流分布差異，并與有限元軟件熱流計算結果進行對比，以評估該解析計算方法的正確性，同時探討解析法在熱流計算工程應用中的相關注意事項。

1 多網格時解析法計算理論

考慮燈絲面向試驗件表面的半圓柱面（見圖1），將燈絲沿長度方向N等分，對燈絲半圓柱面按照圖2 采用數量為n的網格進行等分，將每一個弧面視為平面微元處理，把半圓柱面發出的總熱流平均到每段平面微元上，則每段微元對外發出熱流

式中， r 為燈絲半徑；L 為燈絲長度；N為燈絲沿長度方向均分數量；σ 為斯忒藩-玻耳茲曼常數；T為燈絲溫度。

從圖3可以看到，各微元中心位置與接收面上不同位置的接收點之間的相對關系呈現出更為復雜的變化，n=1 時，可認為所有接收點均可受到平面微元的照射，而當n>1 時，對于非水平位置的平面微元，在接收面上均存在一個其微元對應半球空間能夠產生照射的極限位置，即微元所在平面與接收面的交點位置（如圖3中的C和C'、D和D'）。

將不同θ位置的各段平面微元對G點產生的熱流密度累加，便得到了在單根石英燈作用下G點的熱流密度。對于多根石英燈組成的陣列，可以按照上述過程逐根計算每根石英燈對各接收點產生的熱流密度再累加，便可以得到不同網格數量劃分下的熱流密度場分布情況。

2 計算結果

計算長度300mm、直徑2.6mm、距接收面高度50mm的單根石英燈燈絲在300mm×400mm區域的熱流分布情況，接收面長度和寬度方向每10mm選取一個節點，燈絲沿長度方向等分為100 段，燈絲溫度取2000K，計算得到不同網格數量下的熱流分布云圖如圖5 所示，圖中白色虛線為燈絲在接收面上的投影。可以看到當n 從1 增加到4 時，熱流分布以及熱流峰值水平呈現出相對顯著的變化，而當n 從4繼續增加后，熱流分布以及熱流峰值水平相對變化不明顯。

提取上述云圖中y=200 以及x=150 兩條直線上各節點以及坐標為（150，200）的中心點的熱流密度計算結果，如圖6 所示。可以看到，隨著n 的增大，熱流密度峰值逐漸減小，當n 從1 增加到2 時，熱流密度峰值變化尤為明顯，當n 從4繼續增加后，峰值水平逐漸趨于穩定；與中心點處相反，坐標為（0，200）和（300，200）的左右邊緣處熱流密度，當n=1 時取得最小值；而坐標為（150，0）和（150，400）的上下邊緣處熱流密度隨著n 的增大，呈現出與中心點處相似的變化趨勢。這是由于當n=1 時燈絲所有輻射能量集中于一個矩形窄帶內，根據蘭貝特（Lambert）定律，此時燈絲對正下方區域的定向輻射強度達到最高值，而對左右邊緣處的輻射作用受較小的可見面積影響相對最弱；而隨著n 的增加，燈絲的輻射能量均勻分散于沿環向的各矩形窄帶內，燈絲對正下方區域的輻射強度逐漸降低，但由于兩者間距離最短，所以正下方區域的熱流密度仍然達到當前分布的最高值，而對左右邊緣處的輻射作用受到環向窄帶對其可見面積增大的影響相對n=1 有了一定增強，但受到可對其照射的窄帶數量降低的影響，增強幅度較小；隨著n 的增加，環向窄帶對上下邊緣處的可見面積相對n=1 進一步降低，因此上下邊緣處受到的輻射作用當n=1時達到最強。

3 有限元計算結果對比

3.1 熱流計算結果對比

目前，商用有限元軟件在石英燈輻射熱流計算中被廣泛采用，將解析法計算結果與相同參數條件下的有限元計算結果進行對比，通過結果的一致性可以驗證解析法的正確性。采用與第2 節中相同尺寸、高度和溫度的燈絲以及相同大小的接收面，接收面網格尺寸為10mm×10mm，燈絲沿長度方向等分為100個網格，按照面面輻射法在有限元軟件計算不同環向網格數量下的熱流密度分布，計算模型如圖7 所示。同樣提取橫豎中軸線上的熱流分布結果，與解析法計算結果進行對比，如圖8 所示。從圖8 中可以看到，當n=1 時（此時燈絲有限元模型與解析法模型均為矩形窄帶），解析法與有限元軟件計算結果呈現出顯著差異，峰值熱流水平相對誤差達到61.7%；隨著n 的增加（此時燈絲有限元模型與解析法模型均為半圓柱面，并且環向網格數相同），兩種方法的熱流計算結果的差異逐漸縮小，當n>4時，計算結果一致性相對保持穩定；與解析法計算結果不同，有限元軟件計算得到的中心點處熱流密度隨著n 的增大呈現出小幅的增大。

因此，燈絲對外輻射總能量差異應是解析法與有限元熱流計算結果差異的最主要原因。α 隨n 的變化如圖9 所示，從圖9 可以看到，隨著n 的增大，輻射總能量之比逐漸趨近于1；當n=1 時，解析法中輻射總能量是有限元軟件中的1.57 倍，在相同燈絲溫度情況下，有限元軟件中采用燈絲直徑寬度的矩形窄帶來計算石英燈輻射傳熱會顯著降低結構溫度響應的精度；當n=3 時，輻射總能量之比達到1.047，n=4 時，輻射總能量之比達到1.026，因此建議在有限元建模時對燈絲半圓柱至少采取3～4個的網格進行劃分。

為驗證解析法的有效性，將有限元計算中的燈絲溫度結合式（9）進行調整，以保證二者在計算過程中燈絲對外輻射總能量一致，進而通過對比驗證解析法計算結果的正確性。將有限元計算中的燈絲溫度T調整為式（10）中的T

式中，T為調整后的有限元計算中燈絲溫度；T 為給定的石英燈燈絲溫度。

3.3 調整燈絲溫度后的熱流計算結果對比

對有限元模型中的燈絲溫度按照式（10）調整后，重新計算并對比解析法與有限元計算熱流密度結果，如圖10所示。從圖10可以看到，經過燈絲溫度調整、保證輻射總能量一致后，有限元計算結果與解析法計算結果的一致性相對調整前得到顯著改善，n=1 時中心點熱流密度相對誤差僅有3.23%。其他產生誤差的原因可能是有限元軟件中熱流密度為單元計算結果，此處的節點熱流密度結果是將其相鄰單元的結果進行平均得到，而解析法計算結果則是該節點位置的理論計算結果，兩者存在一定偏差。

4 結論

針對石英燈輻射熱流計算，本文提出了適用于半圓柱面多網格劃分時的解析計算方法，對比了與有限元軟件熱流計算結果的差異，確定了差異產生原因并建立了調整方法，驗證了解析計算方法的正確性，得到了以下結論：

（1）采用解析法計算石英燈熱流密度分布時，采用不同的網格數量劃分燈絲半圓柱面會對計算結果產生影響，當n≥4后，熱流密度峰值水平逐漸趨于穩定。

（2）n=1 時采用解析法計算石英燈熱流密度雖然計算精度相對較差，但由于燈絲網格數量較少且不存在照射范圍判斷問題，在同等計算條件下能夠獲得最快的計算速度，因此在大規模石英燈陣列的輻射熱流預估中仍然有其工程應用價值；當n 達到3 以上時，解析法計算結果可以滿足工程中對石英燈熱流仿真計算精度的要求。

（3）在相同燈絲溫度情況下，有限元軟件中采用燈絲直徑寬度的矩形窄帶來計算石英燈輻射傳熱會顯著降低結構溫度響應的精度，因此在有限元建模時對燈絲半圓柱面至少采取3～4個的網格進行劃分。

（4）燈絲對外輻射總能量差異是解析法與有限元熱流計算結果差異的最主要原因，將燈絲溫度按照總能量比進行調整后，解析法與有限元計算結果的一致性得到大幅改善，證明了本文提出的多網格時熱流密度解析計算方法的正確性。

（5）采用本文提出的解析法計算的熱流密度值可作為熱流載荷邊界施加于有限元分析模型中，從而省去在有限元軟件中對石英燈建模的過程，同時可以不受部分有限元軟件對輻射換熱單元總數的限制，能夠滿足大尺度結構與大規模石英燈陣列輻射傳熱下結構溫度響應快速計算的需求。