火箭發(fā)動機水擊壓力數(shù)值模擬分析

徐 峰，劉英元，陳海峰

（西安航天動力試驗技術(shù)研究所，陜西 西安710100）

0 引言

液體火箭發(fā)動機啟動、關(guān)機過程是十分重要而復(fù)雜的瞬變過程。此過程中所有部件幾乎都是從靜止?fàn)顟B(tài)（或工作狀態(tài)）迅速變?yōu)楣ぷ鳡顟B(tài)（或靜止?fàn)顟B(tài)），在很短的時間內(nèi)各種參數(shù)都要發(fā)生劇烈變化。實踐證明發(fā)動機的大部分故障是發(fā)生在啟動和關(guān)機這樣的瞬變過程中。開、關(guān)機過程中的水擊是造成液體推進(jìn)劑泄漏或發(fā)動機出現(xiàn)故障、降低發(fā)動機可靠性的原因之一。

1 推進(jìn)劑供應(yīng)管路流動分析

如何建立數(shù)學(xué)模型，理論模擬發(fā)動機熱試車時的工作狀態(tài)成為現(xiàn)今試驗技術(shù)的研究課題。對于定型的發(fā)動機試驗，發(fā)動機構(gòu)造基本固定，相同型號的發(fā)動機熱試車水擊壓力的數(shù)值變化不是很大，但對于研制階段的液體火箭發(fā)動機，發(fā)動機型號變化較大、管路復(fù)雜、接口較多，水擊壓力在開關(guān)機時的沖擊很容易造成推進(jìn)劑管路接口處的破裂、推進(jìn)劑的泄露，從而影響發(fā)動機正常工作。目前獲得水擊壓力的常用方法有分析法(能量法)、圖解法和特征線法3種。目前發(fā)動機試驗常采用圖解法，用于判斷開關(guān)機時泵前水擊壓力的峰值與谷值，這種方法軟件處理響應(yīng)時間短、圖形清晰明了，但無法綜合分析影響水擊壓力數(shù)值的全面因素。

結(jié)合現(xiàn)有計算機技術(shù)，經(jīng)分析論證采用有限差分特征線法作為數(shù)值求解的主要方法。采用龍格庫塔法動力學(xué)方程，能夠獲得發(fā)動機關(guān)機時水擊產(chǎn)生的最高壓力，分析發(fā)動機系統(tǒng)關(guān)機水擊產(chǎn)生的壓力和流量瞬變特性。

現(xiàn)以某液體火箭發(fā)動機為例，采用特征線法建立推進(jìn)劑供應(yīng)系統(tǒng)一維流動數(shù)學(xué)模型。液體火箭發(fā)動機推進(jìn)劑輸送管路總長度較短，但構(gòu)成復(fù)雜，管路系統(tǒng)中有直管、彎管、異徑管和3通等分支管。在同一系統(tǒng)中的管道材料、壁厚和直徑可有不同，故水擊波的傳播速度沿管道而變，這造成了復(fù)雜的水擊波反射和疊加。實際推進(jìn)劑供應(yīng)管路如圖1所示，推進(jìn)劑從恒定壓力的儲囊壓入管道，流量采用閥門控制。因試驗需要，引出一路測壓導(dǎo)管，其中位置P為測壓點。

圖1 推進(jìn)劑供應(yīng)管路Fig.1 Propellant feeding pipe

對于圖1所示管路，水擊一維不穩(wěn)定流動的運動方程和連續(xù)方程分別為：

式中：h為水頭，表示壓力；x為距離；t為時間；v為流速；g為重力加速度；λ為摩擦系數(shù)；D為管道直徑；a為管道壓力波的波速；E為管道材料彈性模量；K為流體體積彈性模量；δ為管道壁厚；ρ為流體密度。

水擊運動方程和連續(xù)方程是一階線性偏微分方程，用特征線法轉(zhuǎn)換為常微分方程，可得沿兩條特征線的常微分方程組：

1)沿c+特征線

2)沿c-特征線

式中：c為特征值，d為合定段長度。對公式 (4)、(5)，分別令dh=hpi-hi-1,dv=vpi-vi-1和dh=hi+1-hpi,dv=vi+1-vpi，差分離散化可得沿c+、c-特征線的單一管路瞬時壓力計算式分別為：

式中：Δx為距離步長，且Δt=Δx/c；下標(biāo)i-1、i+1分別表示第i-1和i+1個節(jié)點；下標(biāo)pi表示節(jié)點i在下一個時間步長時的變量。此處Δt為時間步長。

當(dāng)管道包括端點N+1個分點的v，h的初值已知時，用公式 (5)、 (6)可求得網(wǎng)格內(nèi)部結(jié)點的hpi，vpi，管道兩端的hp，vp，則需由公式(5)、 (6)及管道兩端的邊界條件確定。管道兩端的基本邊界條件如下：

1)上游水擊壓力hp1恒定為h0；

2)管道封閉端流量為零，有v恒等于0；

此處的τt為閥門相對開度隨時間變化、vmax為全開時恒定流動時閥門處最大流速、h1為此時閥門兩端壓力差。

分析管道內(nèi)劇烈的壓力瞬變造成的液柱分離，假設(shè)液體中無空氣釋放，當(dāng)管道中某個計算截面的壓力低于液體汽化壓力時，液體汽化形成蒸汽穴；截面保持恒定的蒸汽壓力，截面積為A。當(dāng)空穴存在時，令Qin，Qout分別為Δt內(nèi)的平均流入和流出截面的流量，則空穴體積由Σ(Qout-Qin)給定 (求和自汽穴開始出現(xiàn)后進(jìn)行)?？昭ㄏr，由兩液柱接觸形成的壓頭增加量用下式表示：

2 水擊壓力數(shù)值模擬

采用Flowmaste仿真軟件對實際推進(jìn)劑供應(yīng)系統(tǒng)的管道、閥門、3通及異徑管等各部件建立數(shù)學(xué)模型。在相應(yīng)的模塊計算中，管道、彎管及異徑管等常規(guī)部件的流動阻力從數(shù)據(jù)庫提取，下游處的電磁閥和5通連接管的過流特性數(shù)據(jù)則用流體模擬通用軟件CFD獲得。取物性參數(shù)20℃時推進(jìn)劑的粘度 μ=0.419×10-1Pa·s，ρ=1445 kg/m3，飽和壓力ps=0.102 MPa，K=1.33 GPa，不同材料管道楊氏彈性模量分別為E1=70 GPa，E2=115 GPa。不同管段的a根據(jù)管道直徑和壁厚由公式 (3)給出，計算條件為：

1） 儲囊出口處恒定壓力p=6 MPa；

2） 穩(wěn)態(tài)工作推進(jìn)劑體積流量Q=1.8×10-1m3/s；

3） 閥門全開至全關(guān)時間T=0.3 ms。

因管路分段多，取瞬態(tài)計算的Δt=0.001 ms，使時間步長滿足瞬態(tài)Courant在基于CFD軟件群的思想設(shè)計的FLUENT流動計算軟件下，Courent是指時間步長和空間步長的相對關(guān)系。一般來說，隨著Courent的越小，收斂速度越快，但穩(wěn)定性越差。在實際設(shè)置的時候，先設(shè)置較小Courent，如果穩(wěn)定性越好，則適當(dāng)增加Courent）的穩(wěn)定性條件Δt/Δx≤1/a。對不同時間測壓點和閥門處水擊壓力進(jìn)行試驗計算，計算的壓力波傳遞周期和最高壓力值及其出現(xiàn)時刻與試驗結(jié)果均吻合；試驗和計算結(jié)果都出現(xiàn)了明顯的汽穴現(xiàn)象,水擊最小壓力為0.102 MPa。受實際條件和現(xiàn)有設(shè)備的影響，實際結(jié)果還需以后進(jìn)行驗證。

3 測壓導(dǎo)管對水擊壓力測量的影響

目前常溫壓力傳感器用于低溫測量時用測壓導(dǎo)管引出，測壓導(dǎo)管的長度在0.6～2 m、內(nèi)徑6～8 mm范圍內(nèi)。通過水擊計算，可以發(fā)現(xiàn)測壓導(dǎo)管端點的最高壓力大于閥門處的最高壓力。由圖2～3可知：測壓點處水擊壓力最大為18.5 MPa，關(guān)閉閥處水擊壓力最大為16.1 MPa，且壓力波周期也不一致，在12.5 ms的計算時間內(nèi)兩者相差1個周期。試驗測量壓力最大為約19 MPa，壓力波的周期與支管端點的計算周期一致。實驗結(jié)果說明實測壓力數(shù)據(jù)反映的是測壓導(dǎo)管端部的壓力瞬變特性，與發(fā)動機實際工作時管路內(nèi)壓力瞬變特性有差異。因該系統(tǒng)尺寸相對較小 (管路直徑為毫米級，總長度小于300 mm)以及發(fā)動機空間結(jié)構(gòu)的限制，無法直接測量閥門處的壓力，因此若根據(jù)測壓導(dǎo)管端部壓力瞬變數(shù)據(jù)設(shè)計推進(jìn)劑供應(yīng)管路，則存在較大誤差。為此，用實例計算分析測壓導(dǎo)管端點壓力與閥前壓力的關(guān)系。

為簡化計算，在推進(jìn)劑管路實驗基礎(chǔ)上構(gòu)建簡單分支管計算模型，如圖2所示。

圖2 簡單分支管路Fig.2 Simple branch pipeline

圖中節(jié)點3和4分別為閥門和支管封閉端。計算時流體為推進(jìn)劑N2O4，Q=1.28×10-1m3/s，上游恒定工作壓力p=6 MPa，主管總長L=L2+L3=600 mm，取波速a=920 m/s，使用電磁閥，關(guān)閉時間T=0.3 ms。

經(jīng)過試驗驗證：與無支管相比，有支管時因封閉端壓力波反射會使閥門處水擊壓力的升高；總體來說，支管越短，對閥門處水擊壓力的影響就越小，且支管端部壓力越接近閥門處的真實水擊壓力；支管內(nèi)的流體汽化對支管封閉端壓力的影響較大，使之明顯異于閥門處壓力，且支管直徑越小，發(fā)生汽化的可能性和強度也大。

由圖2可知：當(dāng)L4/L3=0.5，d4/d3=1時分支管對閥門處水擊壓力影響最大，產(chǎn)生的最高水擊壓力為12 MPa，較無支管時最高水擊壓力 (9.4 MPa)增加了28%。分析表明：測壓導(dǎo)管對原管路水擊壓力的影響較大，且附加測壓導(dǎo)管后最高壓力測量值均大于實際最高壓力；當(dāng)主管路無汽化時，測壓導(dǎo)管端部也可能出現(xiàn)汽化，使壓力瞬變復(fù)雜化，增大測量壓力與實際壓力偏差；對小尺寸管路系統(tǒng)，應(yīng)盡可能避免使用小直徑長測壓導(dǎo)管，如確有需要，應(yīng)合理選擇長度和直徑，盡可能減小測量值與真實值的偏差。

4 結(jié)束語

本文利用特征線法對推進(jìn)劑供應(yīng)系統(tǒng)一維流動的理論數(shù)學(xué)模型，用Flowmaster模擬軟件進(jìn)行關(guān)機水擊數(shù)值模擬，并通過理論分析測壓導(dǎo)管不同點的水擊壓力，確定出測壓導(dǎo)管對原管路水擊壓力的影響。由于流體計算的復(fù)雜性和發(fā)動機管路結(jié)構(gòu)的多樣性，本文僅對一種實例進(jìn)行了簡單的理論分析，是否對各種發(fā)動機管路均有好的模擬結(jié)果，還需要對此法進(jìn)一步考核和研究。

[1]胡平信,劉國球.液體火箭發(fā)動機的技術(shù)發(fā)展和展望[J].導(dǎo)彈與航天運載技術(shù),1998,27(2):12－17.

[2]趙萬明.液氧/煤油發(fā)動機試車主要參數(shù)測量方法研究[J].火箭推進(jìn),2006,32(10):23-26.

[3]黎勤武,張為華,王振國,等.空間發(fā)動機系統(tǒng)關(guān)機過程水擊現(xiàn)象理論分析[J].推進(jìn)技術(shù),1998,19(3):27－30.

[4]聶萬勝,戴德海,夏鵬.姿控推進(jìn)系統(tǒng)發(fā)動機關(guān)機的管路瞬變特性[J].推進(jìn)技術(shù),2003,24(1):628－630.

[5]WYLIE E B,STREETER V L.Fluidtransients[M].New York:McGraw2Hill Inc.,1978.

[6]JENNINGS P C,PEDDIESON J,MUNUKUTLA S.New exact solution for fluid transients[J].Mechanics Research Communications,2005,32(3):59-72.

[7]侯詠梅.水頭損失對水擊計算的影響分析[J].河南科學(xué),2003,21(3):304-307.

[8]李治勤.流速及管道特性對水擊的影響[J].太原理工大學(xué)學(xué)報,2000,3(2):156-158.