吳朋朋，楊月誠，高雙武，趙奇國

(第二炮兵工程大學(xué)，西安 710025)

0 引言

固體火箭發(fā)動機工作過程中，噴管承受著高溫、高壓、高速燃氣流的沖刷，環(huán)境非常惡劣。改進噴管是提高火箭發(fā)動機性能重要途徑之一。氣動性能設(shè)計、結(jié)構(gòu)強度設(shè)計都是噴管設(shè)計的重要內(nèi)容。

隨著推進技術(shù)的發(fā)展，運載火箭的助推級或第一級發(fā)動機正在采用越來越大面積比的噴管以提高高空性能。但大面積比噴管在地面試車以及發(fā)動機的啟動和關(guān)機過程中，經(jīng)常會產(chǎn)生噴管側(cè)向載荷，導(dǎo)致噴管受損。Wang T.S 等人[1－2]建立了三維非定常流動模型，對影響噴管側(cè)向載荷的因素進行了研究分析，得出了側(cè)向載荷與燃燒波、激波轉(zhuǎn)變和噴管出口部位的激波震動有關(guān)。

在很大的空氣彈性變形情況下，由于流動和結(jié)構(gòu)相互作用，將會引起側(cè)向載荷巨大的增長[3]。文中運用MpCCI(mesh-based parallel code coupling interface)耦合器作為計算流體動力學(xué)軟件FLUENT和有限元分析軟件ABAQUS的數(shù)據(jù)交換平臺，對流場與噴管結(jié)構(gòu)運動變形進行耦合計算。重點研究大面積比噴管在不同入口總壓時出現(xiàn)的側(cè)向載荷。分析了不同階段側(cè)向載荷的特點，展現(xiàn)了噴管流場與結(jié)構(gòu)之間的相互影響，并對噴管結(jié)構(gòu)的應(yīng)力、位移進行了分析。

1 流固耦合數(shù)理模型

文中以某大面積發(fā)動機比噴管為研究背景。該噴管總長1.6m，噴管面積比 ε=55.2。結(jié)構(gòu)計算只考慮燃氣流動對噴管結(jié)構(gòu)變形的影響。計算區(qū)域為噴管和外場，為捕捉非對稱現(xiàn)象，計算區(qū)域采用二維軸對稱模型旋轉(zhuǎn)360°得到，未做對稱性假設(shè)。

網(wǎng)格生成軟件GAMBIT用于構(gòu)造流場區(qū)域。對于流場區(qū)域，由于噴管結(jié)構(gòu)變形，需采用動網(wǎng)格技術(shù)來模擬，所以接近噴管內(nèi)壁的流場局部區(qū)域采用了四面體非結(jié)構(gòu)化網(wǎng)格，流場其它區(qū)域均采用六面體網(wǎng)格(見圖1(a))。因此文中對非對稱現(xiàn)象的研究將主要和流動性質(zhì)有關(guān)。有限元軟件ABAQUS用于構(gòu)造噴管結(jié)構(gòu)，單元類型為C3D8R(見圖1(b))。

圖1 流場與結(jié)構(gòu)計算網(wǎng)格模型

1)結(jié)構(gòu)參數(shù)

文中所選噴管結(jié)構(gòu)的主要參數(shù)為:密度 ρ=7850kg/m3，彈性模量 E=2.0 ×1011Pa，泊松比 γ =0.33，屈服強度 σs=1000MPa。

2)初始條件

文中以某固體火箭發(fā)動機噴管不同入口總壓(3～7MPa)，入口溫度均為1500K作為噴管流場非穩(wěn)態(tài)計算的初始值進行仿真計算。噴管結(jié)構(gòu)內(nèi)表面為耦合邊界。噴管流場取海平面的壓強與溫度作為流固耦合初始條件。噴管頭部壁面為絕熱固壁邊界。壓力遠場邊界條件取海平面的壓強和溫度。

3)網(wǎng)格無關(guān)性驗證

由于側(cè)向載荷較難算準(zhǔn)，側(cè)向載荷對網(wǎng)格密度有較強的依賴性。為選取合適的流場計算模型網(wǎng)格，文中進行了網(wǎng)格無關(guān)性驗證。選取了網(wǎng)格數(shù)為2303766的計算模型網(wǎng)格用于文中研究內(nèi)容(見圖1(a))。

4)數(shù)據(jù)監(jiān)測

數(shù)據(jù)監(jiān)測是獲得噴管內(nèi)流場壓強變化的有效方法。為了獲得噴管出口部位0°(1.59，0.49，0)和180°(1.59，－0.49，0)兩個對稱點的壓強變化并進行對比。文中選取這兩個對稱點作為監(jiān)測點，監(jiān)測噴管工作過程中噴管出口部位壓強變化。

5)計算方法

計算及程序編譯采用刀片服務(wù)器IBM BladeCenter高性能計算集群。流場計算中，湍流模型采用RNG k-ε模型。整個耦合過程中，采取時間步為10－5s。流場由FLUENT軟件計算，采用耦合隱式求解方法，整個計算具有二階精度。噴管結(jié)構(gòu)運動變形計算采用ABAQUS隱式求解模塊ABAQUS/Standard模塊。

2 結(jié)果分析

2.1 流場分析

耦合過程中，通過計算流體動力學(xué)軟件FLUENT計算結(jié)果，顯示了噴管在不同入口總壓(3～7MPa)，入口溫度均為1500K的條件下，噴管側(cè)向載荷隨時間變化的過程。

由圖2可以看出，隨著氣流向下游傳播，側(cè)向載荷數(shù)值先后經(jīng)歷了由分離模式轉(zhuǎn)換和噴管出口部位的激波振蕩，兩種不對稱狀態(tài)導(dǎo)致的較強側(cè)向載荷。

在0～0.5ms附近為噴管核心噴射流動狀態(tài)(core jet flow)[1]及噴管自由分離激波狀態(tài)(FSS)。此階段由兩個狀態(tài)組成，即在噴管核心噴射流動狀態(tài)之后，隨著流動的發(fā)展，噴管內(nèi)流場產(chǎn)生了不利于噴管流動的壓強梯度。噴管內(nèi)部喉部下游出現(xiàn)了激波。隨著流動的進一步發(fā)展，噴管流場呈現(xiàn)馬赫盤激波模態(tài)(Mach disk shock pattern)[4]。受噴管內(nèi)馬赫盤的強烈阻擋 ，氣流有繞開其流動的趨勢，從而使馬赫盤下游形成了低速區(qū)。由于康達效應(yīng)(Conda Effect)，靠近噴管壁面的附面層中流動是亞音速的，噴管內(nèi)部的流動往往會受外界擾動影響。當(dāng)附面層氣流的動能不能克服激波后的壓強梯度時，氣流與噴管壁面發(fā)生分離，激波進入噴管內(nèi)部，噴管內(nèi)部出現(xiàn)了激波和氣流分離。此時的激波分離后沒有再與噴管壁面附著，稱這種激波分離為自由分離激波。由圖2可以發(fā)現(xiàn)此階段出現(xiàn)的側(cè)向載荷較小。

圖2 不同入口總壓條件下噴管側(cè)向載荷隨時間變化曲線

0.5ms之后，由圖2可以發(fā)現(xiàn)在3～7MPa各個不同入口總壓條件下，均出現(xiàn)了第一次較大的側(cè)向載荷，此階段噴管內(nèi)部流場正在經(jīng)歷由自由激波分離模式向約束激波分離模式轉(zhuǎn)換。在此階段可以發(fā)現(xiàn)一個有趣的現(xiàn)象，隨入口總壓的增大，由激波轉(zhuǎn)變導(dǎo)致的側(cè)向載荷數(shù)值逐漸接近側(cè)向載荷峰值。當(dāng)入口總壓為6MPa、7MPa時，發(fā)現(xiàn)由激波轉(zhuǎn)變導(dǎo)致的側(cè)向載荷已經(jīng)成為此條件下的側(cè)向載荷峰值。

隨著流動的發(fā)展，側(cè)向載荷出現(xiàn)了回落，主要呈現(xiàn)持續(xù)較久的側(cè)向載荷波動。此階段噴管流場呈現(xiàn)馬赫盤激波模態(tài)中的約束激波分離模態(tài)(Mach disk flow:RSS)。在分離點與再附著點之間，流場形成了一個封閉的分離區(qū)。由于分離區(qū)的范圍有限 ，所以這種分離模式被稱為約束激波分離。

在約束激波分離模態(tài)之后，由圖2可以發(fā)現(xiàn)噴管側(cè)向載荷再一次出現(xiàn)了較大載荷波動。導(dǎo)致此階段出現(xiàn)側(cè)向載荷較大波動的原因是噴管出口部位的激波振蕩(lip RSS oscillation)。此階段約束激波已經(jīng)隨著流動的發(fā)展，達到噴管出口部位附近，由于流場區(qū)域需要不斷進行調(diào)整以適應(yīng)噴管擴張段流動和因為噴管結(jié)構(gòu)的振動，導(dǎo)致了約束激波在出口部位來回移動，但由于激波移動的不對稱性，導(dǎo)致了噴管出口部位流動再附著位置的不對稱，這樣就導(dǎo)致了嚴(yán)重的側(cè)向載荷。以入口總壓為5MPa，在6.05ms時的噴管內(nèi)流流場XY截面馬赫數(shù)云圖上(見圖3)，可以發(fā)現(xiàn)約束激波達到噴管出口部位附近，由圖3也可以清晰顯示此時的流場狀態(tài)。由圖2還可以發(fā)現(xiàn)入口總壓在3～5MPa條件下，側(cè)向載荷峰值均由噴管出口部位的激波振蕩導(dǎo)致。還發(fā)現(xiàn)入口總壓在4MPa條件下，由于噴管出口部位的激波振蕩導(dǎo)致的側(cè)向載荷達到了 8562N，遠遠高于在其它入口總壓條件下產(chǎn)生的側(cè)向載荷峰值。

由圖4可以看到在入口總壓為4MPa的條件下，在7～14ms期間噴管出口部位對稱點壓強呈現(xiàn)顯著的不一致性，這也進一步表明了此階段會出現(xiàn)較大的側(cè)向載荷。在圖4中入口總壓為4MPa與入口總壓為5MPa條件下的噴管出口部位對稱點壓強時間曲線進行對比，可以分析出入口總壓為4MPa的條件下將會導(dǎo)致更為嚴(yán)重側(cè)向載荷。結(jié)合圖2～圖4分析，可以得出此階段出現(xiàn)的較強側(cè)向載荷正是由于噴管出口部位流場不對稱性即噴管出口部位的激波不對稱性導(dǎo)致。

圖3 6.05ms時噴管內(nèi)流場馬赫數(shù)云圖

圖4 噴管出口部位對稱點壓強隨時間變化曲線

由上述分析，可以得出，對于此噴管設(shè)計使用中應(yīng)著力避免入口總壓為4MPa階段或者避免在噴管出口部位的激波振蕩。

圖5為噴管三維內(nèi)流場速度與溫度云圖，顯示了快速升溫階段中，6.05ms時噴管內(nèi)流場速度與溫度云圖。其中三維內(nèi)流場速度云圖為噴管XY平面截圖與三維視圖聯(lián)合顯示，三維內(nèi)流場溫度云圖為完全三維透視圖。圖5較好顯示了此時噴管內(nèi)流場的不對稱。

圖5 6.05 ms時噴管內(nèi)流場三維速度與溫度云圖

在噴管出口部位的激波振蕩之后，隨著流動的進一步發(fā)展，側(cè)向載荷回落至平穩(wěn)的狀態(tài)，此階段側(cè)向載荷數(shù)值很小，且噴管出口質(zhì)量流量平穩(wěn)。所以將此階段稱為噴管內(nèi)流場充分流動狀態(tài)。以入口總壓為5MPa為例，圖6顯示了此階段的XY截面內(nèi)流場速度云圖。

圖6 8.5ms時噴管的速度云圖分布

2.2 結(jié)構(gòu)分析

耦合過程中，通過有限元軟件 ABAQUS計算結(jié)果，顯示了入口總壓為4MPa的條件下噴管結(jié)構(gòu)應(yīng)力、位移分布圖。

圖7為11ms時噴管結(jié)構(gòu)Mises應(yīng)力圖?？梢园l(fā)現(xiàn)此時噴管結(jié)構(gòu)最大的Mises應(yīng)力約為 13MPa，遠小于噴管結(jié)構(gòu)屈服強度1000MPa。因此，此時噴管結(jié)構(gòu)不會發(fā)生破壞。

圖7 11ms時噴管結(jié)構(gòu)應(yīng)力圖

圖8 11ms時Y軸方向噴管結(jié)構(gòu)應(yīng)力圖

由圖 8、圖 9可以分析發(fā)現(xiàn)，噴管結(jié)構(gòu)應(yīng)力變化很大，將引起結(jié)構(gòu)振動。同時也可以發(fā)現(xiàn)此時噴管出口附近呈現(xiàn)了很明顯的應(yīng)力不對稱性，這也進一步說明了此時噴管內(nèi)流場出現(xiàn)明顯不對稱現(xiàn)象，產(chǎn)生了較大的側(cè)向載荷。

圖9 11ms時Z軸方向噴管結(jié)構(gòu)應(yīng)力圖

圖10 噴管外壁中點位移隨時間變化曲線

圖10為噴管出口外壁中點(1.6，0.52，0)的位移時間變化曲線。由圖可以發(fā)現(xiàn)，曲線較好反映了流場參數(shù)變化對結(jié)構(gòu)造成的影響。

3 結(jié)論

1)采用流固耦合方法，對固體火箭發(fā)動機噴管分離流動及側(cè)向載荷進行了研究。對噴管的受力、位移進行了分析，獲得了流場與噴管結(jié)構(gòu)耦合作用的機理。

2)通過對FLUENT軟件二次開發(fā)，使其具備了準(zhǔn)確反映側(cè)向載荷隨時間變化情況等功能。發(fā)現(xiàn)了較強側(cè)向載荷是由分離模式轉(zhuǎn)換和噴管出口部位的激波振蕩兩種不對稱狀態(tài)產(chǎn)生，提出了不同入口總壓條件下導(dǎo)致的噴管側(cè)向載荷大小不同，導(dǎo)致側(cè)向載荷峰值的原因也不盡相同，對于此噴管，要著力避免入口總壓為4MPa階段或者通過改進噴管避免在噴管出口部位的激波振蕩以避免較大側(cè)向載荷對噴管造成的損害。這些研究結(jié)論將對大面積噴管優(yōu)化設(shè)計具有重要的指導(dǎo)作用。

[1]Wang T-S. Transient three-dimensional startup side load analysis of a regenerative cooled nozzle，AIAA 2008－4300[R].2008.

[2]Wang T-S，Mike Guidos. Transient three-dimensional side load analysis of a film cooled nozzle，AIAA 2008－4297[R].2008.

[3]S J Zhang，T Fuchiwaki. Aeroelastic coupling and side loads in rocket nozzles，AIAA 2008 －4064[R].2008.

[4]Frey M，Hagemann G.Status of flow separation prediction in rocket nozzles[C]//34th AIAA/ASME/SAE/ASEE Joint Propulsion Conference Exhibit，1998.