張曉東，張培林，傅建平，王 成，楊玉棟

（軍械工程學(xué)院一系，河北 石家莊 050003）

火炮會在射擊瞬間產(chǎn)生巨大的后坐沖擊力，制退機是消耗后坐能量、控制平穩(wěn)后坐的關(guān)鍵部件。后坐沖擊時間僅為0.1～0.2s，后坐速度可達到約10m／s［1］。在此工況下，制退機內(nèi)部液體呈三維非定常湍流流動，并存在運動界面和真空中的射流［2］，流動現(xiàn)象復(fù)雜，研究其內(nèi)部流場規(guī)律非常困難。

隨著計算機技術(shù)、數(shù)值計算技術(shù)和湍流理論的快速發(fā)展，數(shù)值模擬手段已逐漸被應(yīng)用于制退機內(nèi)部的復(fù)雜流動研究。從制退機局部結(jié)構(gòu)到簡化的二維模型，從層流模型到非定常湍流模型的計算［3－7］都取得了很多有益的結(jié)論。但到目前為止，數(shù)值模擬結(jié)果與實驗結(jié)果仍存在較大的誤差，主要是計算的建模誤差所致［8］。究其原因，一是模型過于簡化，二是受計算理論和計算條件所限。

基于Boussinesq假設(shè)的k－ε雙方程湍流模型是當(dāng)前流體機械內(nèi)部流場計算最常的用湍流模型，主要有標(biāo)準(zhǔn)k－ε 模型［9］，以及改進的 RNGk－ε 模型［10］、Realizable k－ε 模型［11］等。高樂南［4］、陳剛［5］利用標(biāo)準(zhǔn)k－ε模型對制退機簡化模型的定常湍流流動進行了數(shù)值模擬。吳利民［6］利用RNGk－ε模型對制退機簡化模型的內(nèi)部湍流流場進行了數(shù)值模擬。以上采用的部分簡化模型如圖1所示。

圖1 2種簡化的制退機模型Fig.1 Two simplified models of recoil brake

本文中基于制退機內(nèi)部流場三維模型，利用動網(wǎng)格和滑移網(wǎng)格技術(shù)，對后坐沖擊過程中制退機內(nèi)部湍流流場進行數(shù)值模擬。為檢驗數(shù)值模擬的準(zhǔn)確性，分別采用3種不同的k－ε湍流模型進行計算，并將數(shù)值模擬結(jié)果與實驗結(jié)果進行對比，以期為制退機內(nèi)部流場計算的湍流模型選擇提供參考。

1 基本控制方程與湍流模型

假設(shè)制退機內(nèi)部液體進行不可壓縮三維非定常湍流流動，在直角坐標(biāo)系下建立矢量形式的流體連續(xù)性方程和動量方程分別為

式中：ρ為流體密度，t為時間，u為速度矢量，p 為壓強，xi、xj、xk為坐標(biāo)分量，Si、Sj、Sk分別為廣義源項分量，μeff為有效粘性系數(shù)。μeff＝μ＋μt，μ為分子粘性系數(shù)，μt為湍動粘度。

在渦粘模型方法中，不直接處理應(yīng)力項，而是把湍流應(yīng)力表示成湍動粘度的函數(shù)其中Cμ為計算常數(shù)。對以上控制方程進行時間平均，引入Reynolds應(yīng)力項，就可以基于k－ε湍流模型建立湍動能k與耗散率ε的輸運方程。再聯(lián)立μt與湍流模型，進而可對控制方程進行求解。

1.1 標(biāo)準(zhǔn)k－ε模型

標(biāo)準(zhǔn)k－ε模型［9］是目前使用最廣的湍流模型，具有適用范圍廣、計算精度合理的優(yōu)點，是一種針對高Re數(shù)的湍流計算模型。當(dāng)遇到彎曲壁面流動、強旋流和逆壓梯度較大的問題時，計算精度會降低。當(dāng)流動不可壓，且不考慮源項時，定義k與ε的輸運方程為

式中：Gk為湍動能k的產(chǎn)生項，各常數(shù)C1ε＝1.44，C2ε＝1.92，Cμ＝0.09，σk＝1.0，σε＝1.3。

1.2 RNGk－ε模型

RNGk－ε模型［10］通過重正化群方法，修正了湍動粘度，考慮了平均流動中的旋轉(zhuǎn)及旋流流動情況，在ε方程中出現(xiàn)了考慮非平衡應(yīng)變率影響的附加耗散生成項，反映了主流的時均應(yīng)變率Eij，可以更好地處理具有強曲率影響和壁面約束的湍流分離流動。RNGk－ε模型定義k與ε的輸運方程為

1.3 Realizable k－ε 模型

針對標(biāo)準(zhǔn)模型對時均應(yīng)變率很大時，可能導(dǎo)致負的正應(yīng)力的情況，Realizable k－ε模型［11］取消了湍動粘度計算中Cμ的常數(shù)取值，引入相關(guān)的旋轉(zhuǎn)和曲率，使其與應(yīng)變率相關(guān)聯(lián)。Realizable k－ε模型定義k與ε的輸運方程為

2 建模、網(wǎng)格劃分與數(shù)值求解

根據(jù)某型火炮制退機的實際結(jié)構(gòu)，建立了三維計算模型，模型參數(shù)采用真實尺寸，如圖2所示。

圖2 制退機三維計算模型Fig.2 3－D computation model of recoil brake

2.1 網(wǎng)格劃分

因制退機內(nèi)部結(jié)構(gòu)復(fù)雜，既有規(guī)則的圓柱形內(nèi)腔，又有傾斜的流液孔和變截面的環(huán)形縫隙。因此，合理的網(wǎng)格劃分將是提高流場計算準(zhǔn)確度的重要前提。本文中采用分塊對接的網(wǎng)格劃分方法，將制退機內(nèi)部分為若干區(qū)域，對形狀簡單規(guī)則、會產(chǎn)生大變形的區(qū)域，如工作腔、非工作腔、復(fù)進節(jié)制腔等運用結(jié)構(gòu)網(wǎng)格進行劃分。對形狀復(fù)雜、不涉及變形的區(qū)域，如制退桿活塞、調(diào)速筒等運用適應(yīng)性更強的非結(jié)構(gòu)網(wǎng)格進行劃分。通過滑移網(wǎng)格界面可實現(xiàn)運動區(qū)域與靜止區(qū)域間的數(shù)據(jù)實時傳遞。考慮到制退機內(nèi)部結(jié)構(gòu)為中心軸對稱，為方便計算，計算模型選為整體結(jié)構(gòu)的1／2，計算網(wǎng)格數(shù)為365 842。

2.2 網(wǎng)格運動方法

為真實模擬后坐沖擊過程中制退機內(nèi)部液體流動，本文利用動態(tài)層變網(wǎng)格技術(shù)，通過非結(jié)構(gòu)網(wǎng)格區(qū)域整體運動與結(jié)構(gòu)網(wǎng)格區(qū)域內(nèi)網(wǎng)格的增減，實現(xiàn)了在給定速度下的制退機內(nèi)流場運動。結(jié)構(gòu)網(wǎng)格更新采用動態(tài)層變法［12］，該方法根據(jù)與運動物面相鄰的網(wǎng)格層高度來決定增加或減少網(wǎng)格層數(shù)，動網(wǎng)格模型指定一個理想層高，鄰近動邊界的網(wǎng)格單元根據(jù)層高度來分裂出新單元層或與鄰近層合并。以復(fù)進節(jié)制腔區(qū)域為例，后坐過程中處于被拉伸狀態(tài)，則動邊界相鄰的單元層被擴展至hmin，hmin為單元層最小高度，滿足hmin＞ （1 ＋αs）hi，αs為分裂因子，hi為理想單元層高。若條件滿足，該層將根據(jù)指定條件（恒定高度或恒定比例）進行分裂，即在層j將分裂為2層：1層高為hi，1層高為（hmin–hi）。反之，被壓縮的單元層將與鄰近的單元層合并成1個新層。

2.3 數(shù)值求解方法

利用有限體積法將控制方程離散為計算模型網(wǎng)格單元上的代數(shù)方程，擴散項、壓力項采用二階中心差分格式離散，動量項采用一階迎風(fēng)格式離散，湍動能、湍流耗散率采用二階迎風(fēng)格式離散，近壁區(qū)域采用增強型壁面函數(shù)方法處理，壓力速度耦合方程采用SIMPLEC算法求解。

3 計算結(jié)果與分析

采用標(biāo)準(zhǔn)k－ε模型、RNGk－ε模型、Realizable k－ε模型分別對制退機內(nèi)部流場進行數(shù)值計算，得到工作腔、復(fù)進節(jié)制腔的壓力數(shù)據(jù)，并與實驗測試的壓力結(jié)果進行了對比，如圖3所示。

從圖3可以看出，后坐過程中制退機內(nèi)的工作腔和復(fù)進節(jié)制腔壓力的計算結(jié)果均較好地反映出實際壓力的變化規(guī)律，采用標(biāo)準(zhǔn)k－ε模型的計算結(jié)果與實驗結(jié)果具有更好的總體符合度，其最高壓力點也比其余2種模型的計算結(jié)果更接近實驗值，而RNGk－ε模型和Realizable k－ε模型的計算結(jié)果基本一致，但與實驗結(jié)果的誤差大于標(biāo)準(zhǔn)k－ε模型的，并且相對誤差隨壓力的增減而相應(yīng)地增減。

圖3 壓力的計算值與實驗值的比較Fig.3 Pressures calculated by different k－εmodels compared with experimental results

后坐過程中，在0～13.5ms是彈丸膛內(nèi)運動期，即后坐加速運動期，3種模型計算結(jié)果與實驗值符合較好；在13.5～116.2ms的后效期及慣性后坐期，即后坐減速運動期，標(biāo)準(zhǔn)k－ε模型的計算結(jié)果最接近實驗值，其他2種模型的計算誤差明顯增大；在116.2～172.0ms的后坐結(jié)束期，速度變化較為平穩(wěn)，3種模型的計算值與實驗值基本吻合，誤差都較小。表1給出了工作腔、復(fù)進節(jié)制腔壓力計算結(jié)果與實驗結(jié)果的誤差。經(jīng)對比分析可知，采用標(biāo)準(zhǔn)k－ε模型的數(shù)值計算結(jié)果與實驗結(jié)果吻合最好。

表1 不同湍流模型下的制退機內(nèi)部壓力計算誤差Table1 Pressure computation errors by different turbulence models of recoil brake

4 結(jié) 論

（1）利用動網(wǎng)格與滑移網(wǎng)格技術(shù)，對后坐沖擊過程中，分別采用標(biāo)準(zhǔn)k－ε模型、RNGk－ε模型和Realizable k－ε模型對真實結(jié)構(gòu)的制退機內(nèi)部三維湍流運動流場進行了數(shù)值模擬，結(jié)果與實驗曲線的變化趨勢基本一致，較好地反映了制退機內(nèi)部壓力的變化規(guī)律。

（2）對比分析工作腔、復(fù)進節(jié)制腔壓力計算結(jié)果與實驗結(jié)果，得出在彈丸膛內(nèi)運動期與后坐結(jié)束期，3種模型的壓力計算值均與實驗值較好吻合，誤差較小，而在后效期與慣性后坐期，標(biāo)準(zhǔn)k－ε模型計算值最接近實驗值，誤差不超過5%，小于其他2種模型。

［1］二Ο二研究所.火炮射擊參數(shù)的常規(guī)測試［M］.北京：第五機械工業(yè)部，1980：147－149.

［2］高樹滋，陳運生，張月林，等.火炮反后坐裝置設(shè)計［M］.北京：兵器工業(yè)出版社，1995：54－63.

［3］Chen C J，Neshat H N，Li P.The Finite Analytic Method［M］.Iowa Institute of Hydraulic Research，The University of Iowa，1983.

［4］高樂南.火炮駐退機湍流流場的預(yù)測分析［D］.南京：華東工學(xué)院，1990.

［5］陳剛.駐退機內(nèi)不可壓湍流流動研究［D］.南京：華東工學(xué)院，1991.

［6］吳利民.駐退機湍流場檢測技術(shù)研究及其數(shù)值模擬［D］.南京：南京理工大學(xué)，1999.

［7］鄭建國.火炮制退機流場的數(shù)值模擬［J］.力學(xué)與實踐，2001，23（2）：30－32.ZHENG Jian－guo.Numerical simulation of flow field in a gun recoil mechanism［J］.Mechanics in Engineering，2001，23（2）：30－32.

［8］張楠，沈泓萃，姚惠之.阻力和流場的CFD不確定度分析探討［J］.船舶力學(xué)，2008，12（2）：211－224.ZHANG Nan，SHEN Hong－cui，YAO Hui－zhi.Uncertainty analysis in CFD for resistance and flow field［J］.Journal of Ship Mechanics，2008，12（2）：211－224.

［9］Launder B E，Spalding D B.Lectures in mathematical models of turbulence［M］.London：Academic Press，1972.

［10］Yakhot V，Orzag S A.Renormalization group analysis of turbulence：Basic theory［J］.Journal of scientific Computing，1986，1（1）：3－11.

［11］Shih T H，Liou W W，Shabbir A，et al.A new k－εeddy viscosity model for high reynolds number turbulent flows－model development and validation［J］.Computers ＆ Fluids，1995，24（3）：227－238.

［12］張志榮，冉景煜，張力，等.內(nèi)燃機缸內(nèi)CFD瞬態(tài)分析中動態(tài)網(wǎng)格劃分技術(shù)［J］.重慶大學(xué)學(xué)報：自然科學(xué)版，2005，28（11）：97－100.ZHANG Zhi－rong，RAN Jing－yu，ZHANG Li，et al.Techniques of dynamic mesh in the CFD transient analysis of an engine［J］.Journal of Chongqing University：Natural Science Edition，2005，28（11）：97－100.