周晨初，李舒欣，陳宏玉，王 丹，任孝文

(西安航天動(dòng)力研究所 液體火箭發(fā)動(dòng)機(jī)技術(shù)重點(diǎn)實(shí)驗(yàn)室， 陜西 西安 710100)

0 引言

液體火箭發(fā)動(dòng)機(jī)在工作過程中，由于液體局部加速度過大引起管路中壓強(qiáng)急劇變化形成水擊[1]。水擊所產(chǎn)生的壓強(qiáng)通?？梢赃_(dá)到管道穩(wěn)態(tài)壓強(qiáng)的幾倍到幾十倍不等，經(jīng)常造成發(fā)動(dòng)機(jī)及其試驗(yàn)臺(tái)系統(tǒng)結(jié)構(gòu)破壞。火箭發(fā)動(dòng)機(jī)推進(jìn)劑供應(yīng)管路是一個(gè)復(fù)雜的管網(wǎng)系統(tǒng)，尤其近年來隨著發(fā)動(dòng)機(jī)高度集成和快速響應(yīng)的發(fā)展，其關(guān)機(jī)水擊問題更加突出[2]。

Menabrea最早開展水擊理論研究，指出管道具有彈性和水流可以壓縮，計(jì)算水擊時(shí)應(yīng)考慮波速的影響[3]。李文勛建立了瞬變流的基本微分方程，奠定了水擊分析的理論基礎(chǔ)[4]。WYLIE E B等人探討了水擊連續(xù)性方程能夠在恒定流條件下恒成立的問題，猜測可能由于沒有考慮摩擦阻力或管道傾斜對(duì)水流狀態(tài)的影響[5]。侯詠梅對(duì)水擊連續(xù)性方程不滿足恒流條件的問題做了初步探討，考慮管道的摩擦阻力與傾斜度的影響，并對(duì)基本方程進(jìn)行了重新推導(dǎo)[6]。水擊問題可由一維瞬變流動(dòng)基本方程描述，它們是一組非線性雙曲型偏微分方程，一般情況下不存在解析解，只能采用數(shù)值計(jì)算方法求解[7-8]。隨著計(jì)算機(jī)技術(shù)的發(fā)展，國內(nèi)外學(xué)者采用數(shù)值仿真的方法，對(duì)水擊問題進(jìn)行了大量研究，其中包括：特征線法[5，9-13]、有限元法[9-12]、無網(wǎng)格法[13-15]和譜方法等[8,16-18]。由于發(fā)動(dòng)機(jī)管網(wǎng)復(fù)雜，采用特征線法在時(shí)間步長協(xié)調(diào)上比較困難，邊界條件的處理也比較復(fù)雜，而且非線性迭代解算收斂較慢，限制了特征線法在復(fù)雜系統(tǒng)(如補(bǔ)燃循環(huán)發(fā)動(dòng)機(jī)) 動(dòng)態(tài)特性分析方面的應(yīng)用[12]。譜方法的最大優(yōu)點(diǎn)是“無窮階收斂性”，對(duì)于方程的解充分光滑情況下，所求得的近似解將以基函數(shù)個(gè)數(shù)的任意次冪速度收斂于精確解[18]，但水擊仿真過程中由于存在參數(shù)間斷或大梯度變化時(shí)，采用譜方法容易造成Gibbs振蕩現(xiàn)象。一維有限體積法將計(jì)算域沿流體流動(dòng)方向劃分為若干不重復(fù)但可重疊的控制體(狀態(tài)單元和速度單元)，認(rèn)為控制體內(nèi)物理量不變，將描述連續(xù)介質(zhì)流動(dòng)力學(xué)守恒方程應(yīng)用于每個(gè)控制體，可得到流體的一維方程[13]。該方法物理意義明確，物理量的守恒性好，處理復(fù)雜邊界連接的問題方便，在流體力學(xué)數(shù)值模擬中得到了廣泛應(yīng)用。本文主要基于一維有限體積法開展發(fā)動(dòng)機(jī)水擊特性研究。

1 仿真模型

發(fā)動(dòng)機(jī)閥門關(guān)閉時(shí)，由于流體慣性，貯箱到閥門之間的推進(jìn)劑會(huì)繼續(xù)向閥前流動(dòng)，并逐漸滯止，當(dāng)所有推進(jìn)劑完全滯止時(shí)，流體處于最大壓縮狀態(tài)，管壁處于最大拉伸狀態(tài)，閥前壓強(qiáng)達(dá)到最大；下一時(shí)刻，由于閥前壓強(qiáng)遠(yuǎn)大于貯箱壓強(qiáng)，推進(jìn)劑壓力波回傳，假設(shè)液體和管壁完全是彈性的，不考慮能量損耗，當(dāng)管壁恢復(fù)到正常狀態(tài)，閥前壓強(qiáng)恢復(fù)到平衡態(tài)；當(dāng)管壁進(jìn)一步壓縮，將會(huì)在閥前形成低壓區(qū)，壓力波反射回貯箱入口時(shí)，管壁處于最大壓縮狀態(tài)，閥前壓強(qiáng)達(dá)到最??；此時(shí)由于管壁處于最大壓縮狀態(tài)，管內(nèi)液體處于舒張狀態(tài)，當(dāng)管壁狀態(tài)恢復(fù)后，液體壓力恢復(fù)，形成一個(gè)完整水擊周期。發(fā)動(dòng)機(jī)水擊壓強(qiáng)、頻率與管網(wǎng)組成和閥門關(guān)閉時(shí)間有關(guān)。采用一維有限體積法建立管路、閥門的一維有限體積模型，對(duì)發(fā)動(dòng)機(jī)關(guān)機(jī)水擊進(jìn)行仿真分析。

1.1 管道模型

由于水擊頻率一般較高，為精確捕捉水擊頻率和峰值特征，需要采用分布參數(shù)描述管道水擊過程。從描述黏性不可壓縮流體的N—S方程出發(fā)，采用一維有限體積法，沿管路軸線方向進(jìn)行離散，離散方法如圖1所示。忽略熱傳導(dǎo)、重力場和管壁彈性變形的影響，得到等直徑管可壓縮流一維瞬變流動(dòng)的連續(xù)性方程、動(dòng)量守恒方程和能量守恒方程。

圖1 一維有限體積法交錯(cuò)網(wǎng)格劃分Fig.1 Grid division of one—dimensional finite element

(1)

(2)

式中：ρ為介質(zhì)密度；u為介質(zhì)流速；p為介質(zhì)壓力；E為單位體積介質(zhì)總能；l為分段管長；i=0,1,…，n-1為管道軸線方向狀態(tài)控制單元編號(hào)，認(rèn)為狀態(tài)控制單元內(nèi)流體壓強(qiáng)、溫度等狀態(tài)量保持不變；需要指出的是式中未考慮局部流阻影響。

(3)

式中：A為管路截面積；W為流體動(dòng)量；f為單位質(zhì)量介質(zhì)受到管壁摩擦阻力；j=0,1,…，n為管道軸線方向速度控制單元編號(hào)，且速度控制單元與狀態(tài)控制單元交替存在，認(rèn)為速度控制單元內(nèi)流體速度保持不變。

為使方程封閉，需要補(bǔ)充控制體比能量計(jì)算方程

流體狀態(tài)方程

P=P(ρ,e)

式中e為介質(zhì)比內(nèi)能。

采用真實(shí)流體的立方型狀態(tài)方程[19-20]，能夠較為準(zhǔn)確地模擬純物質(zhì)液相及氣相的壓力、內(nèi)能及密度之間的關(guān)系。流體動(dòng)量計(jì)算方程

式中V為控制單元體積。

作用在管壁上的摩擦阻力系數(shù)計(jì)算公式

式中fλ為無量綱摩擦損失系數(shù)；D為管路直徑。

1.2 閥門模型

由于閥門進(jìn)出口容腔長徑比較小，按流容單元處理，認(rèn)為單元內(nèi)流體狀態(tài)混合均一，閥芯節(jié)流按噴孔模型處理，得到閥門一維有限元模型。

進(jìn)口容腔模型

(4)

出口容腔模型

(5)

節(jié)流孔模型

(6)

考慮到閥門動(dòng)作時(shí)間對(duì)關(guān)機(jī)水擊壓強(qiáng)沖擊影響較大，補(bǔ)充閥芯開度方程A=A(t)。

2 水擊特性仿真分析

2.1 模型驗(yàn)證

利用現(xiàn)有常規(guī)推進(jìn)劑液體火箭發(fā)動(dòng)機(jī)試驗(yàn)臺(tái)系統(tǒng)對(duì)上述管道、閥門模型水擊仿真準(zhǔn)確性進(jìn)行驗(yàn)證。試驗(yàn)臺(tái)系統(tǒng)由貯箱、供應(yīng)路截止閥、管路、推力室和供應(yīng)系統(tǒng)主閥組成，壓力測點(diǎn)位于氧化劑和燃料主閥前，如圖2所示。其中氧化劑管路長13 m；燃料供應(yīng)管路長20.8 m；管路直徑均為32 mm；閥門動(dòng)作時(shí)間12 ms。

圖2 試驗(yàn)臺(tái)系統(tǒng)組成Fig.2 Diagram of test system

試驗(yàn)臺(tái)系統(tǒng)關(guān)機(jī)水擊仿真與試驗(yàn)對(duì)比如圖3和圖4所示。圖3和圖4中，閥門關(guān)閉時(shí)產(chǎn)生第一個(gè)水擊壓力峰，水擊壓力峰維持到壓力波傳播到貯箱出口為止，壓力峰值衰減主要是由于壓力波振蕩過程中節(jié)流圈、管道等元件能量耗散造成的。氧化劑和燃料水擊仿真幅值和振蕩頻率與試驗(yàn)數(shù)據(jù)吻合較好，驗(yàn)證了仿真模型的正確性。

圖3 燃料路閥前水擊Fig.3 Pre—valve pressure of fuel supply line

圖4 氧化劑路閥前水擊Fig.4 Pre—valve pressure of oxidant supply line

2.2 關(guān)機(jī)水擊影響因素分析

為獲得發(fā)動(dòng)機(jī)關(guān)機(jī)水擊的主要影響因素及其規(guī)律，以燃料供應(yīng)管路為研究對(duì)象，針對(duì)推進(jìn)劑流量、節(jié)流元件位置、流量、管路長度、閥門作動(dòng)時(shí)間、推進(jìn)劑含氣率等設(shè)計(jì)變量進(jìn)行仿真分析，無量綱仿真結(jié)果如圖5～圖9所示。

圖5 推進(jìn)劑流量的影響Fig.5 Influence of propellant flow rate

圖6 管路長度影響Fig.6 Influence of pipe length

圖7 節(jié)流元件位置的影響Fig.7 Influence of restrictor position

圖8 閥門關(guān)閉時(shí)間影響Fig.8 Influence of valve actuation time

圖9 推進(jìn)劑含氣率影響Fig.9 Influence of gas rate in propellant

流量分別為設(shè)計(jì)工況的90%、100%、110%時(shí)，燃料閥前水擊壓強(qiáng)如圖5所示。各工況下水擊增量峰值分別為2.88、3.175和3.481 MPa。相比于設(shè)計(jì)工況，90%流量和110%流量下，燃料路水擊壓力分別增加了-9.29%和9.64%。水擊增量與流量成正比，水擊頻率和收斂特性與流量無關(guān)。因此，在發(fā)動(dòng)機(jī)設(shè)計(jì)過程中，為避免關(guān)機(jī)水擊過大，可采取低工況關(guān)機(jī)的措施，以降低關(guān)機(jī)風(fēng)險(xiǎn)。

由于發(fā)動(dòng)機(jī)設(shè)計(jì)過程中，流量一般由推力決定，改變流量意味著需要調(diào)整發(fā)動(dòng)機(jī)工況，難度較大。由流量公式Qm=ρuA可知，密度不變時(shí)，流量與流速和流通面積相關(guān)。由此，可知水擊增量可能與流速或管路截面積相關(guān)。為進(jìn)一步研究流速、管路直徑對(duì)水擊的影響，保持流量不變，流量分別為設(shè)計(jì)工況的90%、100%、110%時(shí)，燃料閥前水擊壓強(qiáng)如圖5所示。各工況下水擊增量峰值分別為2.888、3.175、3.509 MPa。相比于設(shè)計(jì)工況，90%流速和110%流速下，燃料路水擊壓力分別增加了—9.04%和10.52%。水擊增量與流速成正比，水擊頻率和收斂特性與流速無關(guān)。因此，當(dāng)發(fā)動(dòng)機(jī)工況不變時(shí)，可通過改變流速，有效降低發(fā)動(dòng)機(jī)關(guān)機(jī)水擊壓強(qiáng)。

為研究管路長度對(duì)發(fā)動(dòng)機(jī)關(guān)機(jī)水擊的影響，保持其他變量不變，管路長度分別為設(shè)計(jì)管長的100%、80%、30%時(shí)仿真結(jié)果如圖7所示。各工況下水擊增量峰值分別為3.175、3.175、2.469 MPa。相比于設(shè)計(jì)工況，80%管長和30%管長下，燃料路水擊壓力分別增加了0%和10.52%。造成這種現(xiàn)象的原因是當(dāng)管路長度足夠長，水擊壓強(qiáng)有足夠的發(fā)展時(shí)間時(shí)，水擊增量峰值等于最大水擊壓力；當(dāng)管路相對(duì)較短，導(dǎo)致水擊周期小于兩倍的閥門動(dòng)作時(shí)間時(shí)，水擊壓強(qiáng)來不及發(fā)展到峰值，已經(jīng)受到閥門打開的影響開始下降，其峰值小于最大水擊壓力。此外，管路長度變短時(shí)，水擊頻率增大。

為研究局部流阻位置對(duì)發(fā)動(dòng)機(jī)關(guān)機(jī)水擊的影響，保持其他變量不變，分別將節(jié)流圈置于貯箱出口和燃料閥前，關(guān)機(jī)水擊仿真結(jié)果如圖7所示。節(jié)流圈調(diào)整到閥前時(shí)，關(guān)機(jī)水擊振蕩頻率保持不變，壓強(qiáng)峰值增大了2.5%，振蕩收斂速度變慢，說明局部流阻元件位于貯箱出口有利于抑制關(guān)機(jī)水擊。

為研究閥門作動(dòng)時(shí)間對(duì)發(fā)動(dòng)機(jī)關(guān)機(jī)水擊的影響，保持其他變量不變，閥門關(guān)閉作動(dòng)時(shí)間12、24、36 ms時(shí)，關(guān)機(jī)水擊仿真結(jié)果如圖8所示。圖中可以看出，水擊周期為53.2 ms，當(dāng)閥門動(dòng)作時(shí)間小于半個(gè)水擊周期時(shí)，水擊為完全水擊，水擊壓強(qiáng)峰值保持不變；當(dāng)閥門動(dòng)作時(shí)間大于半個(gè)水擊周期時(shí)，水擊表現(xiàn)為間接水擊，水擊壓強(qiáng)峰值低于完全水擊壓強(qiáng)。同時(shí)可以看出，閥門動(dòng)作時(shí)間不影響水擊周期，但閥門動(dòng)作越慢，水擊振蕩相位越滯后。

3 結(jié)論

采用一維有限體積法建立了液體火箭發(fā)動(dòng)機(jī)試驗(yàn)臺(tái)系統(tǒng)關(guān)機(jī)水擊仿真模型，仿真與試驗(yàn)對(duì)比驗(yàn)證了仿真模型的正確性。在此基礎(chǔ)上，針對(duì)工程上常用的水擊抑制手段進(jìn)行仿真分析，主要結(jié)論如下：

1)減小系統(tǒng)流量可以降低閥前關(guān)機(jī)水擊壓強(qiáng)，水擊頻率和收斂性保持不變,發(fā)動(dòng)機(jī)設(shè)計(jì)過程中可采取降工況關(guān)機(jī)的措施，抑制水擊。

2)當(dāng)發(fā)動(dòng)機(jī)工況、流量無法調(diào)整時(shí)，可通過增大管路直徑、降低推進(jìn)劑流速的方法抑制關(guān)機(jī)水擊，且該方法操作簡單，在工程上較為可行。

3)改變管路長度，對(duì)水擊頻率有一定影響，當(dāng)管路足夠長時(shí)水擊壓強(qiáng)峰值保持不變，說明供應(yīng)管路長度與水擊壓強(qiáng)峰值不直接相關(guān)，但當(dāng)管路長度小于一定值，導(dǎo)致水擊周期小于兩倍的閥門動(dòng)作時(shí)間時(shí)，關(guān)機(jī)水擊來不及充分發(fā)展，管路越短，水擊峰值越小。

4)局部流阻元件位于貯箱出口，水擊峰值較低，頻率不變，收斂較快，是抑制關(guān)機(jī)水擊的簡便措施之一。

5)當(dāng)閥門動(dòng)作時(shí)間小于半個(gè)水擊周期時(shí)，水擊峰值保持不變，大于半個(gè)水擊周期時(shí)，水擊表現(xiàn)為間接水擊，水擊峰值低于完全水擊壓強(qiáng)。在條件允許情況下，應(yīng)使閥門動(dòng)作時(shí)間大于半個(gè)水擊周期。

6)推進(jìn)劑中注入一定量氣體，可以改變混合介質(zhì)的彈性模量，降低介質(zhì)聲速，顯著降低水擊壓強(qiáng)和頻率。