毛聰聰，阮文俊，孫雪明，賈昊楠

(1.南京理工大學(xué) 能源與動(dòng)力工程學(xué)院，江蘇 南京 210094；2.江西星火軍工工業(yè)有限公司,江西 南昌 331709；3.中國(guó)人民解放軍32181部隊(duì),河北 石家莊 050000)

火箭彈發(fā)射過程中，發(fā)射筒會(huì)經(jīng)受高溫、高壓、高速燃?xì)馍淞鞯臎_擊。在對(duì)火箭發(fā)射筒進(jìn)行力學(xué)分析時(shí)，也需將發(fā)射筒內(nèi)表面壓力數(shù)據(jù)作為邊界條件，為不同狀態(tài)下的發(fā)射筒有限元分析提供載荷數(shù)據(jù)[1-4]。尤其對(duì)帶導(dǎo)轉(zhuǎn)槽的火箭武器系統(tǒng)進(jìn)行發(fā)射過程分析時(shí)，需要考慮由于導(dǎo)轉(zhuǎn)槽存在而導(dǎo)致的不對(duì)稱力帶來的影響，因此有必要對(duì)帶導(dǎo)轉(zhuǎn)槽的火箭發(fā)射筒三維內(nèi)流場(chǎng)及發(fā)射裝置內(nèi)的不對(duì)稱力進(jìn)行研究。

國(guó)內(nèi)外許多科研人員對(duì)燃?xì)馍淞鞯臄?shù)值模擬做了大量研究。張磊等[5]通過大渦模擬研究了不同噴管口徑對(duì)燃?xì)馍淞髁鲌?chǎng)的影響。鄧春麗等[6]建立了導(dǎo)彈發(fā)射過程二維流場(chǎng)計(jì)算模型，研究了不同時(shí)刻發(fā)射筒的內(nèi)流場(chǎng)和外流場(chǎng)特征。張俊等[7]通過研究某固體火箭發(fā)動(dòng)機(jī)燃?xì)饬鲌?chǎng)的分布，得到了流場(chǎng)軸向與徑向的燃?xì)馍淞鲄?shù)分布規(guī)律。王華等[8]建立了多噴管流場(chǎng)計(jì)算模型，研究了燃?xì)馍淞髁鲌?chǎng)的速度與溫度分布。Sinha等[9]對(duì)燃?xì)馍淞鳑_擊發(fā)射筒過程進(jìn)行數(shù)值模擬，將不同網(wǎng)格數(shù)量和湍流模型得到的計(jì)算結(jié)果進(jìn)行了對(duì)比分析。

為了使數(shù)值模擬更加接近實(shí)際的發(fā)射過程，筆者采用動(dòng)網(wǎng)格技術(shù)對(duì)火箭彈發(fā)射時(shí)的運(yùn)動(dòng)過程進(jìn)行模擬[10-12]。通過建立火箭彈發(fā)射過程發(fā)射筒內(nèi)流場(chǎng)計(jì)算模型，使用SSTk-ω湍流模型對(duì)非穩(wěn)態(tài)流場(chǎng)進(jìn)行數(shù)值模擬，并對(duì)帶導(dǎo)轉(zhuǎn)槽的火箭發(fā)射筒內(nèi)流場(chǎng)和非對(duì)稱力進(jìn)行研究。

1 理論基礎(chǔ)

1.1 控制方程

質(zhì)量守恒方程：

(1)

動(dòng)量方程：

(2)

能量方程：

(3)

1.2 湍流模型

筆者采用的湍流模型為SSTk-ω湍流模型，與標(biāo)準(zhǔn)k-ω模型相比，該模型在湍流粘度定義中考慮了湍流剪切應(yīng)力的輸運(yùn)過程，是對(duì)標(biāo)準(zhǔn)k-ω模型的增強(qiáng)。表達(dá)式為

(4)

(5)

式中：GK表示湍流的動(dòng)能；Gω為ω的方程；Γk和Γω分別代表k與ω的發(fā)散項(xiàng)；Dω代表正交發(fā)散項(xiàng)。

2 計(jì)算模型及網(wǎng)格劃分

2.1 發(fā)射筒實(shí)體模型及計(jì)算模型

某單管火箭武器發(fā)射筒的實(shí)體模型如圖1所示。將發(fā)射筒尾部的凸臺(tái)結(jié)構(gòu)簡(jiǎn)化，并將導(dǎo)轉(zhuǎn)槽簡(jiǎn)化為連續(xù)的螺旋槽后，根據(jù)實(shí)體模型建立的三維仿真計(jì)算物理模型如圖2所示。

因?yàn)椴恍枰紤]火箭彈受到燃?xì)馍淞鞯臎_擊狀況，所以將火箭彈簡(jiǎn)化為噴管。在發(fā)射筒力學(xué)分析時(shí)由于需要考慮導(dǎo)轉(zhuǎn)槽處的應(yīng)力集中，因此可以通過三維帶導(dǎo)轉(zhuǎn)槽的發(fā)射筒內(nèi)流場(chǎng)計(jì)算，得到螺旋導(dǎo)轉(zhuǎn)槽處不同時(shí)刻表面的壓力分布。

2.2 初始及邊界條件

流場(chǎng)初始化壓力大小為101.325 kPa,溫度為300 K。因涉密原因，文中噴管的入口壓力和火箭彈運(yùn)動(dòng)速度都經(jīng)過無量綱化處理。將噴管進(jìn)口設(shè)置為壓力入口，入口燃?xì)鉁囟葹? 500 K。火箭彈發(fā)射時(shí)的入口壓力隨時(shí)間變化曲線由實(shí)驗(yàn)測(cè)定，選取實(shí)驗(yàn)曲線中壓力到達(dá)0.036 6的時(shí)刻作為計(jì)算的起始時(shí)刻，火箭彈運(yùn)動(dòng)出筒時(shí)刻作為計(jì)算的終止時(shí)刻，共歷時(shí)130 ms。從開始計(jì)算到火箭彈運(yùn)動(dòng)出筒的噴管入口壓力時(shí)間曲線如圖3所示。將發(fā)射筒壁面設(shè)定為絕熱無滑移壁面，發(fā)射筒尾端設(shè)置為壓力出口，壓力大小為101.325 kPa，溫度為300 K。根據(jù)噴管入口壓力隨時(shí)間變化曲線得到火箭彈運(yùn)動(dòng)速度隨時(shí)間變化的曲線，如圖4所示。

2.3 計(jì)算條件設(shè)定

流場(chǎng)計(jì)算選擇基于壓力的瞬態(tài)求解器，將燃?xì)猱?dāng)作理想氣體可壓縮流處理。不考慮火箭彈發(fā)射時(shí)的空氣阻力；不考慮燃?xì)饣瘜W(xué)反應(yīng)和固體顆粒的影響；不考慮燃?xì)馀c發(fā)射筒內(nèi)壁及周圍空氣的傳熱；不考慮火箭彈的旋轉(zhuǎn)，僅考慮火箭彈沿水平方向的運(yùn)動(dòng)。流體粘性系數(shù)采用Sutherland公式計(jì)算。

2.4 網(wǎng)格劃分

在劃分網(wǎng)格時(shí)，由于較難生成單塊的高質(zhì)量結(jié)構(gòu)化網(wǎng)格，因此采用分塊網(wǎng)格對(duì)計(jì)算域進(jìn)行劃分，4個(gè)網(wǎng)格劃分區(qū)域分別為噴管、尾部計(jì)算域、導(dǎo)轉(zhuǎn)槽和筒體。筆者采用Icem軟件劃分六面體結(jié)構(gòu)網(wǎng)格，對(duì)邊界區(qū)域添加了邊界層網(wǎng)格后,初始時(shí)刻的網(wǎng)格總數(shù)為60萬(wàn)個(gè)。初始時(shí)刻網(wǎng)格劃分如圖5所示。

2.5 動(dòng)網(wǎng)格設(shè)置

為了處理火箭彈發(fā)射時(shí)的網(wǎng)格變化，在流場(chǎng)計(jì)算時(shí)通過動(dòng)網(wǎng)格對(duì)火箭彈運(yùn)動(dòng)進(jìn)行模擬。噴管網(wǎng)格水平移動(dòng)，尾部計(jì)算域網(wǎng)格被拉長(zhǎng)，筒體與導(dǎo)轉(zhuǎn)槽網(wǎng)格固定不動(dòng)。動(dòng)網(wǎng)格更新方法包括彈性光順法、動(dòng)態(tài)層法和局部重構(gòu)法，根據(jù)本文的網(wǎng)格劃分方式，采用動(dòng)態(tài)層法模擬火箭彈運(yùn)動(dòng)。該方法的基本思想是：根據(jù)與運(yùn)動(dòng)邊界相鄰網(wǎng)格的高度，判斷網(wǎng)格需要分割還是合并，當(dāng)網(wǎng)格高度滿足式(6)時(shí)，網(wǎng)格將會(huì)被分割；而當(dāng)網(wǎng)格高度滿足式(7)時(shí)，兩層網(wǎng)格將會(huì)合并。

h≥(1+α)h0，

(6)

h≤βh0，

(7)

式中：α為網(wǎng)格層的分裂因子；β為網(wǎng)格層的合并因子；h0為理想單元網(wǎng)格高度。

筆者取h0=3 mm,α=0.4，β=0.2。通過Profile文件將火箭彈的運(yùn)動(dòng)速度賦給噴管及尾部計(jì)算域，噴管不斷向前運(yùn)動(dòng)，尾部計(jì)算域網(wǎng)格通過動(dòng)態(tài)層法進(jìn)行更新生成。初始時(shí)刻的運(yùn)動(dòng)域網(wǎng)格如圖6所示，火箭彈發(fā)射時(shí)運(yùn)動(dòng)域的網(wǎng)格變化如圖7所示。

3 數(shù)值模擬結(jié)果及分析

3.1 不同時(shí)刻流場(chǎng)縱截面壓力云圖

計(jì)算時(shí)刻分別為50、90和130 ms時(shí)的流場(chǎng)縱截面壓力云圖，如圖8所示。

因涉密原因，流場(chǎng)計(jì)算得到的壓力、速度及扭轉(zhuǎn)力都經(jīng)過了無量綱化處理。由于噴管內(nèi)壓力變化較大，因此為了更清楚顯示發(fā)射筒內(nèi)的壓力分布，將無量綱壓力在0.366以上的部分用相同顏色顯示。燃?xì)庠趪姽軆?nèi)由亞音速加速到超音速后,在噴管出口形成了超音速燃?xì)馍淞鳌０磁蛎洺潭葎澐?，燃?xì)馍淞鳛榍放蛎浬淞鳌Ｈ細(xì)馍淞鲝膰姽艹隹诹鞒鰰r(shí)的壓力高于周圍環(huán)境壓力，所以一出噴管首先會(huì)膨脹加速。燃?xì)饨?jīng)過噴管出口兩側(cè)的膨脹扇區(qū)后,繼續(xù)膨脹加速同時(shí)壓力不斷降低。當(dāng)噴管出口的膨脹波系到達(dá)射流邊界后，向內(nèi)反射為一系列的壓縮波，壓縮波足夠強(qiáng)時(shí)會(huì)匯合形成斜激波。燃?xì)饨?jīng)過激波后，壓力增加、溫度升高同時(shí)速度降低。當(dāng)燃?xì)鈮毫υ黾拥礁哂诃h(huán)境壓力后，燃?xì)馀蛎浖铀偻瑫r(shí)壓力降低。隨著這樣膨脹壓縮過程的反復(fù)循環(huán)，發(fā)射筒內(nèi)形成了一系列膨脹波與壓縮波結(jié)構(gòu)。

3.2 中心軸線上的速度與壓力分布

圖9和圖10分別為不同時(shí)刻流場(chǎng)中心軸線上的速度與壓力分布。從圖9和圖10中可以看出，燃?xì)馍淞鲝膰姽芰鞒龊?，噴管出口附近的膨脹波和壓縮波強(qiáng)度是最大的，因此燃?xì)庠趪姽芨浇乃俣扰c壓力波動(dòng)也最為明顯。在50 ms時(shí)，燃?xì)庠谥行妮S線上經(jīng)歷了1～2個(gè)膨脹-壓縮波；在90 ms時(shí)，燃?xì)庠谥行妮S線上經(jīng)歷了5～6個(gè)膨脹-壓縮波；在130 ms時(shí)，燃?xì)庠谥行妮S線上經(jīng)歷了12～13個(gè)膨脹-壓縮波。因發(fā)射筒長(zhǎng)度有限，所以當(dāng)燃?xì)獾竭_(dá)發(fā)射筒出口附近時(shí)，燃?xì)獾乃俣扰c壓力仍有較大波動(dòng)。由于粘性效應(yīng)和周圍空氣的作用產(chǎn)生了湍流剪切層，湍流剪切層內(nèi)的粘性阻尼使得膨脹波與壓縮波結(jié)構(gòu)不斷瓦解，因此中心軸線上的燃?xì)馑俣扰c壓力呈現(xiàn)振蕩衰減的趨勢(shì)。

3.3 不同時(shí)刻發(fā)射筒內(nèi)表面壓力分布

分別取計(jì)算時(shí)刻為50、90和130 ms時(shí)的發(fā)射筒內(nèi)表面壓力云圖，如圖11所示。燃?xì)鈴膰姽芰鞒龊?，位于噴管出口后的發(fā)射筒區(qū)域直接受到燃?xì)馍淞鞯膹?qiáng)烈沖擊，因此其表面壓力明顯高于噴管前部未直接受到燃?xì)馍淞鳑_擊部分的壓力。

當(dāng)燃?xì)饨?jīng)歷了在壁面附近的激波后壓力增加，其對(duì)發(fā)射筒壁面的壓力也隨之增加。50 ms時(shí)，燃?xì)鈱?duì)發(fā)射筒壁面的最大壓力為0.078；90 ms時(shí)，燃?xì)鈱?duì)發(fā)射筒壁面的最大壓力為0.065；130 ms時(shí)，燃?xì)鈱?duì)發(fā)射筒壁面的最大壓力為0.049。當(dāng)燃?xì)鈮毫υ黾拥礁哂谥車諝獾膲毫?，燃?xì)馀蛎浖铀俨⑶覊毫Σ粩嘟档?，其?duì)壁面的壓力也逐漸降低。燃?xì)庠诎l(fā)射筒內(nèi)經(jīng)歷了一系列膨脹壓縮波，所以燃?xì)鈱?duì)發(fā)射筒內(nèi)表面的壓力在整體上呈現(xiàn)振蕩衰減的趨勢(shì)。

3.4 不同時(shí)刻導(dǎo)轉(zhuǎn)槽內(nèi)的壓力及扭轉(zhuǎn)力分布

不同時(shí)刻導(dǎo)轉(zhuǎn)槽表面壓力隨距離變化曲線如圖12所示。

在發(fā)射筒右視圖中，位于導(dǎo)轉(zhuǎn)槽左側(cè)的壁面為本文的導(dǎo)轉(zhuǎn)槽左側(cè)面，位于導(dǎo)轉(zhuǎn)槽右側(cè)的壁面為本文的導(dǎo)轉(zhuǎn)槽右側(cè)面。50 ms時(shí)，燃?xì)鈱?duì)導(dǎo)轉(zhuǎn)槽頂面、右側(cè)面和左側(cè)面的最大壓力分別為0.072，0.070和0.065；90 ms時(shí)，燃?xì)鈱?duì)導(dǎo)轉(zhuǎn)槽頂面、右側(cè)面和左側(cè)面的最大壓力分別為0.065，0.060和0.055；130 ms時(shí)，燃?xì)鈱?duì)導(dǎo)轉(zhuǎn)槽頂面、右側(cè)面和左側(cè)面的最大壓力分別為0.045，0.044和0.042。因50 ms時(shí)噴管入口壓力是3個(gè)時(shí)刻中最大的，故燃?xì)鈱?duì)導(dǎo)轉(zhuǎn)槽3個(gè)表面的最大壓力也是3個(gè)時(shí)刻中最大的。在不同時(shí)刻，導(dǎo)轉(zhuǎn)槽頂面受到的最大壓力值始終是3個(gè)表面中最大的，而導(dǎo)轉(zhuǎn)槽右側(cè)面受到的最大壓力略大于導(dǎo)轉(zhuǎn)槽左側(cè)面。

從不同時(shí)刻導(dǎo)轉(zhuǎn)槽表面的壓力分布曲線中可以看出，導(dǎo)轉(zhuǎn)槽表面的壓力沿著中心軸線方向振蕩變化，但并未呈現(xiàn)振蕩衰減的趨勢(shì)。這是由于導(dǎo)轉(zhuǎn)槽表面會(huì)多次反射壓縮波，當(dāng)這些反射的壓縮波疊加后，燃?xì)庠俅谓?jīng)過這些壓縮波會(huì)增加對(duì)導(dǎo)轉(zhuǎn)槽表面的壓力，因此導(dǎo)轉(zhuǎn)槽表面的壓力并未沿著軸線方向振蕩衰減。

導(dǎo)轉(zhuǎn)槽內(nèi)扭轉(zhuǎn)力隨時(shí)間變化的曲線如圖13所示。導(dǎo)轉(zhuǎn)槽兩個(gè)側(cè)面受到燃?xì)馍淞髦苯幼饔妹娣e隨著火箭彈向前運(yùn)動(dòng)而逐漸增加，故導(dǎo)轉(zhuǎn)槽兩個(gè)側(cè)面受到的燃?xì)庾饔昧σ苍诓粩嘣黾印?dǎo)轉(zhuǎn)槽兩個(gè)側(cè)面上受到的燃?xì)馍淞髯饔昧Σ钪凳沟冒l(fā)射裝置內(nèi)產(chǎn)生了不對(duì)稱力即扭轉(zhuǎn)力，而右側(cè)面受到的燃?xì)馍淞髯饔昧κ冀K大于左側(cè)面，因此導(dǎo)轉(zhuǎn)槽內(nèi)扭轉(zhuǎn)力的方向始終垂直于右側(cè)面。扭轉(zhuǎn)力隨著時(shí)間呈現(xiàn)振蕩上升的趨勢(shì)，在火箭彈即將出筒時(shí)刻達(dá)到最大，最大的扭轉(zhuǎn)力值為769。

4 結(jié)論

筆者通過動(dòng)網(wǎng)格技術(shù)模擬了火箭彈從發(fā)射到出筒過程燃?xì)馍淞鲗?duì)發(fā)射筒的沖擊。通過流場(chǎng)計(jì)算得到不同時(shí)刻發(fā)射筒內(nèi)表面的壓力分布，將其加載到包含導(dǎo)轉(zhuǎn)槽結(jié)構(gòu)的三維發(fā)射筒有限元分析模型中進(jìn)行計(jì)算。可以得到發(fā)射筒在發(fā)射壓力作用下的應(yīng)力應(yīng)變情況，從而確定發(fā)射筒在發(fā)射壓力作用下的相對(duì)薄弱部位以及安全儲(chǔ)備。從本文的數(shù)值模擬中可以得到以下結(jié)論：

1)超音速燃?xì)馍淞髋c環(huán)境大氣之間的壓力差，導(dǎo)致發(fā)射筒內(nèi)產(chǎn)生了一系列膨脹壓縮波。由于流體粘性效應(yīng)和周圍大氣的作用，膨脹波與壓縮波結(jié)構(gòu)在發(fā)射筒內(nèi)逐漸崩解。

2)燃?xì)庠诎l(fā)射筒內(nèi)經(jīng)歷了一系列膨脹波和壓縮波，因此不同時(shí)刻中心軸線上的燃?xì)鈮毫?、燃?xì)馑俣群桶l(fā)射筒內(nèi)壁面上的壓力，在整體上呈現(xiàn)振蕩衰減的趨勢(shì)。

3)導(dǎo)轉(zhuǎn)槽內(nèi)壓縮波的反射疊加，使得不同時(shí)刻導(dǎo)轉(zhuǎn)槽表面上的壓力并未沿著中心軸線振蕩衰減。導(dǎo)轉(zhuǎn)槽的存在導(dǎo)致火箭發(fā)射筒內(nèi)產(chǎn)生了扭轉(zhuǎn)力，并且隨著火箭彈向前運(yùn)動(dòng)，扭轉(zhuǎn)力振蕩增加。