再生冷卻推力室熱機疲勞壽命預測研究

黃峻峰, 賀爾銘, 易金翔, 杜大華, 王紅建

(1.西北工業大學 航空學院, 陜西 西安 710072; 2.液體火箭發動機技術重點實驗室, 陜西 西安 710100)

可重復使用火箭發動機是發展綠色航天工程的重點研究工作之一[1-2]。美國私立太空企業SpaceX公司近年成功研制的“梅林”火箭發動機使得“獵鷹9號”和“重型獵鷹”運載火箭具有了一定的重復發射能力。英國開展能重復使用超過100次的“佩刀”發動機的研制工作，并計劃將其安裝在“云霄塔”太空飛機[3]。歐洲宇航院(ESA)的“未來發射器”計劃，基于先進3D打印技術研制的“普羅米修斯”低成本發動機及可重復使用飛行器樣機“西彌斯女神號”都在進行緊鑼密鼓的試驗[2]。俄羅斯針對10～500 kg小載荷商業發射市場，正研制可重復使用的“飛翼-SV”輕型運載火箭及相關發動機[4]。我國已成功研制出50噸級氫氧發動機YF-77和120噸級液氧煤油發動機YF-100，并作為新一代運載火箭的動力系統推動近些年我國航天事業的巨大發展[5-6]。與此同時，我國已展開了可重復使用液氧煤油發動機、液氧液氫發動機、液氧甲烷發動機關鍵技術攻關，對可重復使用運載火箭的研制制定了“三步走”的發展戰略[7]。

推力室作為液體火箭發動機的核心結構，承擔著推進劑混合燃燒、釋放熱量、產生推力及完成能量轉化的工作，自然成為了可重復使用關鍵技術研發和驗證的重點結構。在推力室工作過程中，其內部的燃氣壓力、溫度及流體物性變化十分劇烈，燃氣壓力最大可高達200個大氣壓，燃氣溫度最大可達到3 600 K，最大速度可超過6馬赫[8]。這些復雜的載荷循環會使推力室薄壁夾層結構產生低周疲勞和高溫蠕變效應，而且推力室在循環加載下受到的載荷往往是非對稱應力控制的循環加載，這會導致結構無法恢復的塑性應變逐漸累積產生棘輪效應，大大縮短推力室的使用壽命。

低周疲勞、高溫蠕變效應和棘輪效應這3種行為的共同作用會造成推力室薄壁夾層結構失效，大大降低了推力室的使用壽命。國內外專家學者先后對推力室結構的失效機理及壽命預測方法進行了大量有效工作。

羽中豪等[9]使用準二維傳熱計算方法、二維熱-固耦合方法和局部應變法研究了室壓、推力燃料混合比對推力室棘輪效應的影響。孫冰等[10]使用二維熱-固耦合分析方法對推力室內壁進行非線性平面應變分析，并預測了內壁結構的使用壽命。孫冰和宋佳文[11-13]使用整場耦合方法對液氧甲烷發動機推力室內部的燃燒和冷卻劑流動換熱特性進行三維仿真，獲得了推力室在不同試車階段的應變分布，并比較了瞬態加載和穩態加載下的最大殘余應變。Di等[14]研究了噴油器二維噴流導致推力室壁面周向溫度分布不均的現象，并發現棘輪效應是導致結構失效的主要原因。Pizzarelli[15]用代數模型對推力室使用壽命進行預測分析，預測準確率超過80%。Riccius等[16]發現由燃燒造成的推力室銅合金內壁粗糙度的增加能降低推力室28%的疲勞壽命。

本文針對某液體火箭發動機再生冷卻推力室的薄弱環節結構——流道喉部薄壁夾層結構，建立其局部二維平面應變有限元模型，通過傳熱仿真和非線性結構分析獲得循環載荷下的溫度分布和應力應變響應變化規律，研究了該結構產生的棘輪效應，并揭示了棘輪、蠕變及疲勞損傷對推力室結構使用壽命的影響。

1 物理模型與計算方法

1.1 物理模型

再生冷卻推力室采用了大量如圖1所示的帶銑槽的釬焊薄壁夾層結構，圖中從內至外依次為內壁、再生冷卻通道及外壁，相鄰冷卻通道之間由肋片隔開。在設計階段，為了使火箭發動機具有高推重比，在保證強度/剛度足夠的條件下盡可能使結構輕量化，導致發動機的薄壁夾層結構設計幾乎達到了材料的承載極限。

圖1 推力室薄壁夾層結構及“狗窩”失效示意圖

本文研究對象為某型液氧煤油火箭發動機推力室喉部結構，該結構的內壁和肋片采用NARloy-Z銅合金材料，外壁則采用1Cr21Ni5Ti不銹鋼材料，如圖2所示。

圖2 推力室結構及其喉部局部細節

推力室喉部截面是循環對稱結構，內壁和外壁之間均勻分布著270條再生冷卻通道。根據其軸對稱特性，選取如圖3所示的二分之一冷卻通道截面作為計算模型。

圖3 二分之一冷卻通道物理模型及有限元模型

使用八節點四邊形單元對模型進行網格劃分，共生成15 592個節點和4 901個單元。為了便于后續工作描述模型溫度變化和預測壽命，在模型上依次選取7個參考點，其中A,B點分別為外壁面上、下端，C,D點為冷卻通道肋片的上下頂點，E,F,G點位于靠近高溫燃氣通道的內壁面上。

1.2 控制方程

本文采用單向耦合方法，先通過傳熱分析計算整個工作循環中結構的溫度場分布，再在其基礎上進行瞬態結構分析得到結構的應力應變演化結果。

傳熱分析的控制方程為無內部熱源的導熱微分方程

(1)

式中:T為結構溫度；λ為推力室材料導熱率；ρ為材料密度；c為推力室材料比熱容；為哈密頓算子。

在進行結構分析時,結構的總應變εtotal由熱應變εth、彈性應變εel、塑性應變εpl和蠕變應變εcr組合而成

εtotal=εth+εel+εpl+εcr

(2)

在計算熱應變時使用的公式為

εth=aΔT

(3)

式中：a為材料的熱膨脹系數；ΔT為相對參考溫度的溫度變化量。

在計算彈性應變和塑性應變時,使用Von Mises屈服準則來判斷材料是否達到屈服,其表達式為

(4)

式中：f為屈服函數；s為偏應力張量；α為背應力張量,即屈服面中心；σy為屈服應力。達到屈服面后塑性應變產生,塑性應變沿著屈服面梯度方向增加

(5)

式中：dλ為塑性乘子,即等效塑性應變增量；σ為應力張量。

本文使用Chaboche隨動硬化模型來描述背應力張量的演化過程,在該本構關系中,總背應力由多級背應力分量疊加得到,如(6)式所示

式中:M為背應力分量的級數;Ci和γi為隨動硬化參數;dp為累積塑性應變,其表達式為

(8)

由于金屬材料力學性能受溫度影響較大,本文將Esposito等[17]得到的NARloy-Z材料在高低溫環境下的單軸拉伸試驗數據,使用背應力級數為4的Chaboche模型進行擬合,得到不同溫度下的隨動硬化參數，如表1所示。試驗曲線與擬合曲線如圖4所示,可以看到擬合曲線與試驗曲線[17]在彈性段和塑性段的吻合度都很高,說明隨動硬化模型參數合理準確,能夠準確描述該材料的彈塑性變形。

表1 NARloy-Z材料隨動硬化參數

圖4 NARloy-Z材料不同溫度下的單軸拉伸曲線

在計算高溫蠕變應變時,使用Norton蠕變模型來進行描述,其表達式為

(9)

表2 Norton蠕變模型參數

1.3 壽命預估模型

推力室喉部結構在每個循環下產生的總損傷Dtotal由低周疲勞損傷Df、棘輪損傷Dr和蠕變損傷Dc組成,可按照Miner線性累積損傷理論來計算每個循環下的總損傷[10]

Dtotal=4×(Df+Dr+Dc)

(10)

式中,4為安全因子。

在計算低周疲勞損傷時,以該次循環應力水平下的疲勞循環壽命Nf的倒數來定義低周疲勞損傷

(11)

本文使用基于平均應力的Morrow修正模型來計算每個循環下的低周疲勞損傷,其表達式為

(12)

每個循環下棘輪損傷可以用殘余應變增量與材料斷裂應變的比值來表示

(13)

式中：εend為當前循環結束時的殘余應變；εbegin為當前循環開始時的殘余應變；εu為材料單軸拉伸時的斷裂應變。

蠕變損傷定義為當前應力水平和溫度下的累計保載時間t與該應力下的斷裂時間tc之比

(14)

1.4 邊界條件設置

該型液體火箭發動機一個完整的循環加載分為啟動-熱試車-后冷-關機4個階段,各階段時長分別為5,300,20,600 s,環境溫度為295.15 K。在不同的試驗階段下,對模型邊界條件的設置是不同的。

1) 熱分析的邊界條件設置

將圖3中對稱邊界Γsym設置為絕熱邊界,即在該邊界法線方向上的溫度梯度為0

(15)

式中:n為邊界的法線方向;hsym為對稱邊界的對流傳熱系數。

而內壁與外部環境、冷卻通道和高溫環境直接接觸的邊界為耦合邊界,在該邊界上存在如(16)式所示的對流關系

hcon(Tbulk-Twall)=n·λTwall, ?x∈Γcon

(16)

式中:Γcon為耦合邊界;hcon為耦合邊界的對流傳熱系數;Tbulk為流體溫度;Twall為固體壁面溫度。在熱分析中各耦合邊界條件設置如表3所示,表中下標out、cool和hot分別表示外壁面、冷卻通道壁面和內壁面,CFD代表由流-熱耦合計算得到的壁面參數設置。

表3 熱分析邊界條件設置

2) 結構分析的邊界條件設置

在進行結構分析時,將傳熱分析得到的每一步結構溫度分布作為溫度體載荷。由于推力室薄壁夾層結構對稱特性,將模型兩邊設置為對稱邊界。將各階段下耦合邊界上的載荷分別施加在相應邊界上,各階段下邊界載荷如表4所示。

表4 結構分析邊界條件設置 MPa

2 計算結果與分析

2.1 熱分析結果

圖5為不同階段下的結構溫度分布,圖6為結構上7個參考點溫度隨時間變化曲線。可以看到,

圖5 各階段結構溫度分布

圖6 參考點溫度變化曲線

由于外壁材料的導熱率小于內壁材料,內壁上的參考點溫度會很快到達穩定,而外壁面上的參考點則會出現一定的延遲才能到達穩定。在整個工作過程中,下表面出現最高溫度且該部位溫度變化最大,溫度最大點出現在熱試車階段,此時內壁下表面G,F點處的最高溫度分別為862 K,851 K。

2.2 結構分析結果

圖7和圖8分別為各階段下結構的應力和應變分布。

圖7 各階段結構應力分布 圖8 各階段結構應變分布

在開機階段,結構總體溫度僅提升到309 K左右,雖然結構整體開始膨脹,但由于冷卻通道內煤油的壓力載荷較低,整體結構的壓應變和壓應力都處于較低水平。

在熱試車階段,外壁和肋片上半部的溫度普遍低于450 K,該區域膨脹時內部的壓應變較小。同時由于外壁的熱膨脹系數小于肋片的熱膨脹系數,肋片阻礙了外壁的膨脹,所以外壁C點處出現了最大壓應力434.23 MPa。高溫燃氣會在內壁下表面產生較大的壓力載荷進而阻止附近結構的膨脹,越靠近燃氣的區域受到的阻礙作用越大,因此內壁下表面結構內部的壓應力要小于內壁上表面。此時內壁下表面的G點處出現了最大壓應變2.307×10-2,這是因為該階段下冷卻通道內煤油作用在內壁上表面的壓力大于高溫燃氣對內壁下表面的壓力,使得內壁向燃氣側凸出并彎曲,由此造成G點附近區域受壓嚴重。

在后冷階段,內壁區域溫度在短時間內降低了很多,其中下表面溫度最高降幅達550 K,導致該區域出現了收縮趨勢,因此該區域的應力和應變普遍為正。冷卻通道內煤油對周圍壁面仍有壓力作用,而燃氣作用域內壁下表面的壓力在該階段因燃燒終止而消失,導致內壁向燃燒室凸出的趨勢更加明顯,F點出現了最大拉應變4.467×10-3,拉應力值為205.13 MPa。

盧一平沒有馬上回駁她。這時候,回駁就是爭論,就是口角,就是齟齬,就是擱淺。盧一平沉默了,沉默地盯視著郝桂芹。盯一會兒,有意無意地抬起手腕,看了眼手表。

在關機階段,由于結構在之前的試驗階段已經出現了塑性應變,在卸掉熱載荷和壓力載荷后這些不可恢復的塑性應變作為殘余應變仍存在結構中。關機階段結束時的結構溫度相較于后冷階段結束時的溫度只降低了約14 K,所以兩時刻下的應力和應變分布云圖相差不大。F點處的最大殘余應變為4.437×10-3。

2.3 循環載荷下結構棘輪效應研究

在獲得單次循環加載下結構的溫度變化和應力應變演化規律后,對10個連續循環加載下結構的應力應變演化過程進行了仿真。圖9為不同循環加載次數后的結構殘余應變分布,可以看到殘余應變集中出現在內壁區域,并且每個循環下殘余應變的分布云圖基本相同,只是殘余應變的絕對值隨著循環數的增加而逐漸增大。

圖9 多次循環加載后的殘余應變分布

圖10為內壁上D,E,F,G點在10次循環加載中的應力-應變曲線,可見每個點都是應力控制的非對稱循環加載,并且都產生了棘輪效應。而且在后續幾個循環中參考點的應力-應變曲線形狀趨勢基本相同,呈現了一種漸進穩定性,這就是材料在非對稱應力循環中的調整行為。

圖10 各點應力-應變曲線

為了定量評估各點的棘輪效應大小,采用棘輪應變作為評估指標。使用每個循環中初始殘余應變和循環結束時的殘余應變之差來確定棘輪(殘余)應變

εr=εend-εbegin

(17)

圖11為各點的棘輪應變發展曲線,可以看到在第二次循環后各點的棘輪應變變化與循環次數基本上呈線性變化,即體現了非對稱應力循環中的調整行為。隨著循環數的增加,位于內壁上表面的D,E點棘輪應變大小變化不大,而位于下表面的F,G點的棘輪應變逐漸增大。其中F點在每個循環下殘余拉應變率約為2.3×10-3/次,10次循環加載后殘余拉應變值為2.529×10-2;G點的殘余壓應變率約為2.5×10-3/次,10次循環加載后的殘余壓應變分別為3.193×10-2。

圖11 各點棘輪應變發展曲線

根據結構殘余應變分布,綜合各點的應力-應變曲線和棘輪應變發展曲線,發現內壁下表面的F點和G點會成為結構的潛在失效點。這是因為該兩點的殘余應變最大,而且該兩點在每個循環中的應力和應變變化范圍也較大導致其低周疲勞損傷也很大。

圖12為10次循環后結構的蠕變應變分布及位移云圖,可以看出高溫蠕變效應集中出現在溫度較高的內壁下半部,位于內壁下表面的F點和G點的蠕變效應明顯,進一步加速了這兩點所在區域的結構損傷。

圖12 10次循環后結構蠕變應變及徑向位移分布

圖13為10個循環下內壁上下表面的E點和F點的徑向位移發展曲線。隨著循環數的增加,內壁會逐漸向燃燒室凸出并減薄,這就是推力室結構典型的“狗窩”失效模式。

圖13 E點和F點的徑向位移發展曲線

2.4 推力室結構壽命預測

F點和G點在10次循環下的損傷計算結果如圖14所示?？梢钥吹紽點的損傷形式為棘輪損傷、蠕變損傷及低周疲勞損傷,其中棘輪損傷為最主要的損傷形式,而G點的損傷由蠕變損傷和低周疲勞損傷構成。并且由于調整行為,兩點每個循環下的損傷基本保持不變,其中F點的單循環損傷約為2.95×10-2,總損傷為31.12×10-2,G點的單循環損傷約為2.32×10-2,總損傷為24.124×10-2。根據(10)式計算壽命,F點、G點的失效循環數分別為33次及42次。

圖14 F,G兩點損傷計算結果

綜上，在經過33次循環加載后，推力室內壁下表面中心點會首先發生失效破壞。在導致該點的總損傷中，棘輪損傷、蠕變損傷和低周損傷分別占比52%,32%和16%，說明了棘輪效應是導致推力室結構失效的主要原因。

3 結 論

本文建立了某液氧煤油火箭發動機推力室喉部結構的周期對稱模型，通過三維傳熱分析獲得結構的循環載荷特征，基于各工況下結構應變及損傷分析初步揭示棘輪、蠕變及疲勞對推力室結構使用壽命的影響，得到以下研究結論：

1) 本文建立發動機推力室喉部結構的循環對稱模型，進行推力室流-熱-固載荷傳遞及結構損傷分析，文中壽命預測模型及方法計算效率高、具有重要工程參考價值；

2) 定量分析了推力室內壁結構“狗窩”失效的影響機理，在冷卻通道總損傷中棘輪損傷、蠕變損傷和低周疲勞損傷分別占比52%,32%和16%，因此棘輪效應是導致其失效的主要原因;

3) 結構的蠕變效應主要集中在內壁的下半部分，這是因為在熱試車階段該區域的溫度很高，進而引起應力松弛現象。

本文的研究方法及結果可為液體火箭發動機再生冷卻推力室結構優化設計及快速壽命預測提供重要的工程參考。后續工作可進一步考慮推力室工作時序效應、各損傷模式之間的耦合作用及相關的試驗驗證，將會使推力室結構壽命預測結果更加準確。