基于傳熱增廣模型的軌跡優化與防熱結構分析

張騰飛, 龔春林, 粟 華,*, 薛鵬飛

(1. 西北工業大學陜西省空天飛行器設計重點實驗室, 陜西 西安 710072;2. 空間物理重點實驗室, 北京 100076)

0 引 言

高超聲速飛行器在未來國際局勢中有非常重要的戰略作用[1],這類飛行器具有突防能力強、運輸效率和打擊精度高、殺傷威力大等優勢,但高速飛行造成的熱防護問題是其設計過程中必須著重考慮的因素之一[2-4]。

軌跡優化在飛行器初步方案設計階段,對方案可行性論證、設計參數分析、指導制導控制系統設計等具有重要意義。飛行器軌跡優化問題是一類包含非線性約束的最優控制問題,數值方法已成為解決復雜軌跡優化問題的主要方法[5]。根據是否直接對性能指標進行尋優將其分為間接法和直接法[6-8],隨著偽譜法的興起[9-12]以及數值優化技術的不斷改進,使得在軌跡優化問題中使用更加精細的模型成為可能[13-17]。

在軌跡優化中考慮熱約束是為了保證飛行過程中結構不遭受破壞,而結構破壞與否一般是依據是否超出對應材料的耐溫極限進行判斷的。飛行過程中的熱流密度僅是計算結構溫度的邊界條件,而現有的軌跡優化大多只通過約束熱流密度來體現熱防護要求[18-21]。需要明確的是,過高的熱流密度不是結構破壞的原因,是否超出結構材料的耐溫極限才是關鍵,熱流密度本質上屬于不必要的“軟約束”。隨著高超聲速飛行器設計中對熱防護的重視程度越來越高,在初步軌跡優化設計中就要考慮結構溫度約束。如此,約束熱流密度的方式便存在明顯的不足。

首先,最大熱流密度無法與結構最高溫度建立等價的簡單映射關系,事實上大多數軌跡優化是通過經驗或表面熱輻射平衡條件給出最大熱流密度約束,這顯然是不夠準確的。在一些著重考慮熱約束的研究中[22-23],將熱流密度積分值最小作為優化目標,但這樣得到的軌跡是十分保守的。因此,為了解決考慮結構溫度約束的軌跡設計問題,采用約束熱流密度峰值的軌跡優化模型需基于熱流密度上限越大,結構溫度峰值越高的設想,并且經過多次軌跡優化與傳熱分析的迭代,才能給出保證結構安全又不十分保守的熱流密度上限與飛行軌跡。其次,將傳熱分析與軌跡優化分步進行,則不能考慮到飛行軌跡與結構傳熱之間的相互影響,得到的軌跡優化結果必然不是最優的。綜上,基于一般模型求解會導致軌跡與熱防護結構設計之間反復迭代、互相妥協的問題。此外,為了對后續詳細設計提供參考,方案設計階段需要對不同防熱結構尺寸方案進行可行性分析,而由于熱流密度約束調節范圍有限,導致結構尺寸的可行域較小,不利于詳細設計階段總體層面的協調。

總之,一般約束熱流密度的軌跡優化模型不便于解決考慮結構溫度約束的軌跡優化設計問題,需建立耦合傳熱過程的軌跡優化模型,以便對結構溫度直接進行約束。文獻[24-25]盡管考慮了運動和熱傳導的耦合關系,但研究主要針對最優軌跡周期性巡航的規律研究,并沒有對防熱結構作詳細分析,且所建立的傳熱模型未考慮材料屬性隨溫度的變化導致誤差較大。

本文通過空間離散的方法,將一維傳熱方程轉化為一組一階微分方程,考慮到了材料屬性隨溫度的變化,使得計算精度更高,有利于得出更為準確的分析結論。使用該模型對運動方程進行增廣,增廣之后的系統狀態方程引入了結構溫度這一狀態量,在軌跡優化中可加入更為直接的結構溫度約束,通過一次軌跡優化便可得到滿足結構溫度約束的最優軌跡,而不用反復驗證。

以高超聲速飛行器CAV-H的再入軌跡優化問題為例,首先,對于基準熱防護結構方案,計算結果表明:基于傳熱增廣軌跡優化模型的求解過程更為便捷和高效,且得到的最大航程較一般模型的結果有所提升,最優軌跡在防熱結構性能的利用以及姿態變化等方面也更為合理。其次,基于不同軌跡優化模型確定了防熱結構尺寸的可行域,結果表明,增廣模型對應的可行域遠大于一般模型,這為設計提供了更大的選擇空間。最后,在結構尺寸可行域內計算最優軌跡,分析總結各層結構尺寸對最大航程的影響規律以及以不同形式再入對應的結構傳熱機制的影響,對防熱結構設計提供一定規律性的參考意見。

1 運動模型

為簡化問題和突出重點,僅考慮地球為非旋轉均質圓球的二維情況。根據文獻[26],給出簡化運動方程組:

(1)

式中:r和θ分別為飛行器質心到地心的距離(簡稱地心距)和航程角;v和γ為速度和當地速度傾角;L和D為升力和阻力;m為飛行器質量;μ為引力常數。

氣動力計算如下:

(2)

式中:qd為動壓;ρa為大氣密度(使用USSA1976大氣模型);Sref為參考面積;CL和CD為升阻力系數,是攻角α和速度v的函數。

2 熱傳導模型

2.1 一維傳熱模型

2.1.1 一維傳熱問題

本文不考慮主動熱防護。典型的被動熱防護結構方案一般為高溫防熱層+隔熱層+承力層的多層方案[27-28],其中防熱層一般采用陶瓷基復合材料或C/C復合材料;隔熱層一般采用氣凝膠、柔性隔熱氈等材料;內層承力結構為金屬材料,在計算傳熱時可以忽略,即將三層結構簡化為兩層,并將內表面視為絕熱壁面。

若非詳細研究結構的熱傳導問題,一維模型已經足夠精確[29-30],即只考慮溫度沿厚度方向的分布。一維傳熱問題如圖1所示,圖中qs和qr分別表示表面熱流密度和熱輻射密度,Tω為外表面溫度,δf和δs表示兩層結構的厚度。

一維熱傳導的控制方程和邊界條件分別為

(3)

(4)

式中:T為溫度;ρ、c、k分別為材料的密度、比熱和熱導率,后二者為溫度T的函數;x為結構厚度方向的位置。qs和qr的計算如下:

qs=Csρ0.5v3

(5)

(6)

式中:Cs為與飛行器外形有關的常數;σ和ε分別為表面輻射率和斯特藩-玻爾茲曼常數;T∞為環境溫度。

2.1.2 單元劃分

式(3)是一拋物線型二階非線性偏微分方程,需要對其進行離散化處理。一般熱傳導問題的求解一般是將熱傳導方程在時間和空間上都進行離散,而這里則需要保留溫度對時間的偏導數項,從而與運動方程保持一致。

首先,沿x方向布置節點:要求在內外表面以及交界面上必須布置節點,而在每層結構中,則沿厚度方向均勻布置若干節點。假設防熱層和隔熱層分別布置n和m個節點(交界面處的節點是共用的)。接著,作每兩個相鄰節點連線的中垂線,連同內外表面共同將整個結構劃分為N=n+m-1個單元,每個單元對應一個節點,稱為該單元的控制點。單元劃分結果如圖2所示,Δf和Δs分別表示兩層結構內部的單元厚度。

2.1.3 有限差分模型

對于劃分好單元的防熱結構,通過在單元內部對控制方程進行積分,并利用差分代替導數以消除其中關于x的偏導數項。

首先,對于而非交界面處的單元i,對式(3)在單元內部沿x方向進行積分,得

(7)

式中:xi-和xi+表示單元i左右兩側邊界的位置。上面處理得到的內外邊界處的單元厚度是內部單元厚度的一半,在積分時需要注意積分區間長度的差異。對于式(7)左端項,認為單元內部溫度相同,則式(7)可進一步簡化為

(8)

其次,對于交界面處的單元n,由于材料不同,需要分段積分,即:

(9)

將式(9)進行簡化如下:

(10)

認為交界面處沒有熱阻,則等號右端中間兩項可以消去。進而將式(8)和式(10)表示如下:

(11)

式中:

(12)

(13)

(14)

2.2 傳熱模型驗證

對于式(11),給定熱防護結構參數和熱流邊界條件后,便可以通過數值積分的方法進行求解。下面通過算例對本文所建一維傳熱模型進行分析驗證,具體參數如表1所示。結構整體初始溫度為280 K。

表1 傳熱問題算例

表1中,兩種材料的熱導率和比熱與材料溫度的關系如下所示:

(15)

(16)

熱流密度和環境溫度變化規律(總時間為2 500 s)如下:

(17)

(18)

由于內外表面溫度具有代表性,后面分析主要針對這兩處的溫度。計算發現,隨著隔熱層節點個數的增加,各節點溫度迅速收斂,且收斂速度逐漸變慢;由于防熱層厚度較小,熱導率較大,導致溫度梯度很小,因此較少的節點(但至少有一個內部節點,即防熱層節點數目取3即可)就可得到十分準確的外表面溫度變化規律。

(19)

圖3給出了內外表面節點溫度相對誤差隨隔熱層節點數的變化情況,如前面所述,外表面溫度誤差幾乎為0,而隨著熱量向內表面傳遞,溫度誤差也逐漸積累增大(內表面溫度相對誤差最大,其余各處的相對誤差均在圖中兩條曲線之間)。可以看出,當隔熱層節點個數增大到8之前,外表面溫度相對誤差持續降低,但在此之后基本不再提升,而且節點數為8對應的相對誤差已經足夠小(不超過2%)。因此,最終選取離散模型的節點個數為:n=3,m=8。

3 基于增廣模型構建軌跡優化問題

3.1 系統狀態方程

由式(1)和式(11)即可構成傳熱增廣方程,簡化表示如下:

(20)

式中:[XT]=[rθvγT1T2…T10],共同構成狀態量;F,G分別為4維和10維向量值函數,由式(1)和式(11)的右端項組成;U為控制量,在這里指攻角α。兩式通過qs,T∞關聯在一起。增廣后的系統狀態方程組維數由4增加到14,而且考慮到各變量之間數值上差異較大,都使得數值求解變得困難,因此對方程進行無量綱化處理是必要的。

(21)

需要說明的是,由于關于溫度的一組等式中,右端包含很多量綱復雜的材料屬性,不對每個變量單獨進行無量綱化處理,而是對其整體除以Tref/tref。

3.2 路徑約束與優化目標

一般軌跡優化會對過載、動壓以及熱流密度添加路徑約束以體現結構強度和防熱的要求。除了過載和動壓,本文直接對結構溫度進行路徑約束以代替熱流密度約束,即:

(22)

(23)

(24)

式中:qdmax和nmax分別為動壓和過載上限;最大節點溫度Timax根據對應材料的耐溫極限確定,如對于外表面處的T1max即材料PM1000的耐溫極限:1 570 K,而對于內表面處的T10max根據內部設備要求,取為620 K。

優化目標為最大航程角(在末端高度一定的情況下與最大航程是等價的),即:

minJ=-θ(tf)

(25)

4 算例分析

4.1 問題描述

本文以CAV-H高超聲速再入飛行器的最大航程再入問題為例,構建基于傳熱增廣模型的軌跡優化問題進行求解與分析。首先對飛行器相關參數以及再入軌跡要求作以說明。CAV-H飛行器的相關參數如表2所示。

表2 CAV-H相關參數

氣動模型根據文獻[31]給出的數據擬合如下所示:

(26)

CAV-H飛行器被運載器運送到指定高度和速度狀態之后與分離,開始再入[28]。再入過程一般分為返回段、能量管理段和攻擊段。其中,返回段占據大部分時間,最大航程主要取決于這一階段,本文也主要針對這一階段的最大航程軌跡優化問題進行分析。動壓約束和過載約束分別為70 000 Pa和2.5。具體軌跡要求如表3所示(表中高度H即地心距減地球半徑Re=6 378.145 km)。

表3 再入軌跡要求

4.2 基于不同模型的計算結果對比

針對上述軌跡優化問題,使用偽譜法進行求解,作為對比,同樣求解了相同條件下的約束熱流密度的軌跡優化問題。如前面所述,采用一般模型進行軌跡優化之前無法給定合適的熱流密度約束值,只能使用大致確定的約束值進行初步軌跡優化,然后對結構溫度進行安全性驗證,再根據具體情況對熱流密度約束值進行更新,經過迭代計算,最終確定保證結構溫度滿足要求且熱防護性能無過多剩余的熱流密度約束值。

由于數值計算噪聲的影響,通過梯度算法進行迭代較為困難,本文通過以下爬山法迭代確定熱流密度約束:根據表面熱輻射平衡條件大致確定初始熱流密度約束為10×105W·m-2、經過多次嘗試,確定收斂速度最高的初始迭代步長約為-2×105W·m-2。迭代過程中,當節點溫度與溫度上限之差的符號發生改變時改變搜索方向并將步長減半,由于無法同時使得內外表面溫度均恰好達到溫度上限,迭代以內表面溫度作為收斂依據(低于溫度極限不超過0.1 K)。經過17次迭代計算,得到收斂結果。上述計算與基于增廣模型的計算對比如表4所示。盡管增廣模型變量維度大大增加,但由于飛行狀態與結構溫度之間的耦合關系并不強,一次軌跡優化耗時并沒有過久,且不需要迭代計算,整體計算效率并沒有劣勢。另外,采用一般模型計算時對迭代算法有一定要求,一般很難順利進行迭代和快速收斂,整體來看,采用傳熱增廣模型求解考慮結構溫度約束的軌跡優化問題是十分高效的。

表4 不同模型計算結果對比

一般模型熱流約束為4.047×105W·m-2,內表面最高溫度達到溫度上限,而外表面的最高溫度則低于對應的溫度極限,說明外層結構的防熱性能存在剩余;基于增廣模型的結果,內外表面溫度均達到耐溫極限,末端航程角為106.15°,相比于前者的103.56°提升了2.5%。

考慮到再入過程中速度隨時間是遞減的,而且為便于表示動壓、熱流邊界以及最優升阻比攻角,圖5和圖6中分別給出最優軌跡的高度和攻角隨飛行速度的變化規律。可以看出,基于一般模型的最優軌跡在高度第一次降低到極小值點時達到熱流上限,此后高度基本持續下降;而由增廣模型求解的最優軌跡,取消了對熱流的限制,高度第一次可降低到動壓約束邊界處,由于“再入走廊”的放寬,最終航程有所增加。盡管基于增廣模型的軌跡存在更大程度的跳躍,但對于無動力再入來講嚴格的平衡滑翔條件并非必要的,跳躍式再入是合理的[2],事實上,盡管高度波動較大,但傾角始終保持在±2.5°以內。

從姿態(攻角)變化的角度來看,基于增廣模型結果不存在攻角的大幅變化,是更為穩定的。另外,在無約束的情況下,CAV-H為實現最大航程再入,在整個中期滑翔階段攻角是在最大升阻比攻角附近作幅度逐漸減小的波動變化,這與基于增廣模型的求解結果是類似的,不同的是后者波動略有滯后,且振幅較小,基本保持在最大升阻比攻角曲線之上,主要是為了控制高度下降,降低熱流密度峰值從而滿足外表面溫度約束。反觀限制熱流密度的軌跡,由于再入過程飛行器的勢能沒有充分轉化為動能,向上跳躍的能力有限,后續則維持最大升阻比攻角飛行。

熱流密度、熱輻射密度以及結構溫度隨時間變化曲線如圖7所示。圖7(a)中,增廣模型對應的熱流密度最高達到7.27×105W·m-2,遠超一般模型對熱流密度所限制的4.051×105W·m-2。盡管峰值較高,但與高度類似,熱流密度變化也有較大波動。對熱流密度曲線積分得到增廣模型軌跡總的氣動加熱密度為5.888×108J·m-2,而一般模型則為5.981×108J·m-2。由此可見,“熱流密度上限越大,結構溫度峰值越高”的設想并不嚴謹。另外,一般模型的結果中,熱流密度與熱輻射密度之間存在唯一的分界點,在該點之前外界向結構內部傳熱,之后結構才向外散發熱量;而增廣模型的結果中,熱流密度在再入過程前期就存在小于熱輻射密度的情況,這使得內表面溫度的上升趨勢略有減緩,從圖7(b)可以看出,最終達到溫度限制的時間也有所推遲,前者飛行時間為2 808.2 s,后者為2 855.8 s。

由前面的分析可知,約束熱流密度的模型不僅計算過程繁瑣,而且會造成飛行性能的損失;而采用增廣模型較為便捷且結果更加合理。

4.3 防熱結構分析

在初步方案設計階段,對于一些關鍵設計參數,一方面需要大致確定其可行域,為后續的詳細設計階段提供參考;另一方面則需要分析與總結其對總體性能的影響規律。下面就利用傳熱增廣軌跡優化模型對防熱結構尺寸(這里指厚度)這一關鍵設計參數進行分析。

4.3.1 防熱結構尺寸可行域分析

由于防熱層尺寸一般較小且要考慮諸如結構變形、維持外形等問題,這里給定防熱層厚度的設計范圍為0.012～0.030 m。另外,考慮內部空間的制約,要求整體結構厚度不超過0.13 m,防熱結構尺寸基本設計空間如圖8灰色區域所示。首先,針對不同防熱層厚度,基于增廣模型對不同隔熱層厚度方案進行軌跡優化,根據是否可行大致確定隔熱層最小厚度構成圖8中的尺寸下界,黃色區域表示基于增廣模型的防熱結構尺寸可行域。其次,對一般模型,計算確定了熱流密度約束值最低為2.5×105W·m-2,基于此約束下最優軌跡的熱流密度規律,確定不同防熱層厚度下的最小隔熱層厚度,圖8中紫色區域則表示基于一般模型的可行域。由于熱流密度約束可調節的范圍十分有限,即便對于熱流密度峰值最小的軌跡,仍需較厚的隔熱層,導致基于一般模型的尺寸可行域較增廣模型有很大的收縮。

以黃色區域內(紫色區域外)一點:[0.025 m,0.075 m]為例,其軌跡高度與結構溫度變化規律如圖9所示。相比于基準方案,為滿足結構溫度要求,這里高度“跳躍”的次數和幅度都有所增加,外表面降溫的次數也更多。盡管航程稍短,但通過一般模型判定為不可行的方案,使用增廣模型計算卻是可行的。

4.3.2 防熱結構尺寸對飛行性能的影響

4.3.3 不同再入形式對結構傳熱的影響

可以看出,能夠滿足一定航程的防熱結構尺寸方案并不唯一,而在相同的航程下,不同的結構尺寸對應的再入軌跡形式是有所差別的。下面選取尺寸可行域內航程相同的兩點:a=[0.012 m,0.099 m],b=[0.030 m,0.089 m],分析不同形式的再入軌跡對防熱結構傳熱的影響。對應的軌跡以及結構內外表面溫度變化規律如圖12所示(a和b對應的曲線分別用1和2表示)。

可以看出,軌跡1只在前期進行了2次跳躍,后面一直保持滑翔飛行,而軌跡2則持續跳躍,但最終的航程角是相同的。兩條不同形式的再入軌跡,均受到結構溫度約束的影響:對于軌跡1來講,由于防熱層較薄,為了保證防熱層溫度不超出限制,需要在飛行前期進行跳躍,減緩外表面溫度持續升高,而較厚的隔熱層使得內表面溫度升高較慢,允許外表面溫度長時間保持較高溫度。對于軌跡2來講,防熱層較厚,相比之下飛行初期外表面溫度升高的趨勢十分緩慢,但由于隔熱層太薄,即使外表面溫度較低,再入飛行中后期需要始終保持跳躍以周期性降低外表面溫度,延緩內表面溫度的升高。圖13給出了兩組結果對應的結構溫度分布隨飛行過程的變化情況對比(左側為a點對應的結果,右側為b點對應的結果),δx表示距離內表面的厚度。

可以看出,對于防熱層更薄的a,從30%航程開始一直到70%航程的飛行歷程中,防熱層結構溫度都保持在較高值,而對于b,由于防熱層很厚,結構溫度遠低于極限溫度,造成飛行前中期的防熱結構性能浪費,但盡管外層結構溫度較低,而由于長時間的飛行,熱量持續地向內傳遞,也使得內表面溫度最終達到極限。另外,在再入結束時刻,a對應的溫度分布比較均勻,而b對應的結構溫度分布則存在較大梯度,說明b在再入飛行結束時,結構有更多的熱冗余,存在一定的安全隱患。因此,對于長時間的再入飛行,在一定的結構尺寸約束下,應盡可能增大隔熱層厚度,防熱層只需滿足在再入初期飛行軌跡的協助下,外層結構溫度不超出極限即可。

4.3.4 軌跡與防熱結構相互影響規律總結

上面首先對防熱結構尺寸可行域進行了計算和對比,基于增廣模型的可行域較一般模型有很大擴展,有利于避免由于結構尺寸可行域過小導致設計過程不能繼續進行的隱患,體現出飛行軌跡對結構傳熱過程的良好調節作用。

其次,利用增廣模型分析了防熱結構尺寸改變對飛行性能的影響規律,防熱結構足夠厚時,熱約束變為消極約束,而在熱約束為積極約束的結構尺寸區域內,隨著厚度的增大航程角有所提升,其中隔熱層厚度影響較大而隔熱層厚度影響很小;結合結構總尺寸以及總質量密度隨各層尺寸的變化趨勢,發現在有限的結構尺寸和質量內,盡可能厚的隔熱層是使得航程最大的防熱結構尺寸方案。

最后,從更為細致傳熱過程的角度分析了相同航程下,較厚的隔熱層方案相比于較厚的防熱層方案的優勢,即后者對結構防熱性能利用不充分,而且最終時刻的熱冗余較大。

5 結 論

本文采用空間差分的方法處理一維傳熱方程,考慮了材料屬性隨溫度的變化,計算精度更高。將所建傳熱模型與運動模型組成傳熱增廣的軌跡優化系統狀態方程,用于構建高超聲速飛行器的軌跡優化問題。通過對算例的數值仿真,得到以下結論:

(1) 與約束熱流密度的軌跡優化相比,基于傳熱增廣模型的軌跡優化不需要進行軌跡優化與傳熱分析反復迭代驗證,且能夠體現軌跡對結構溫度的調節能力,航程也更遠;

(2) 傳熱增廣軌跡優化模型變量之間耦合關系并不強,因此維數的增大并沒有十分影響計算效率,能夠滿足方案初步設計階段對防熱結構方案可行性與設計參數快速分析的要求;

(3) 使用增廣模型可以得到更大的尺寸區域內的最優軌跡,有利于詳細設計階段總體層面的協調;

(4) 內表面溫度限制是制約較長時間再入過程的主要因素,在有限的結構尺寸和質量下,應使用盡可能厚的隔熱層,防熱層只需保證在再入前中期外層結構溫度不超出約束即可。

本文提出的傳熱增廣模型仍有改進空間,下一步工作中可建立飛行器表面不同防熱區域處的一維傳熱模型,為飛行器全表面的被動熱防護結構提供設計參考;另外,對于熱防護系統為主被動防熱相結合的飛行器,還可在傳熱模型中考慮主動冷卻,從而為主動冷卻控制規律設計以及主被動防熱規模的權衡提供參考。