小口徑火炮彈炮耦合動態響應有限元時程分析

蘇忠亭，徐 達，李曉偉，韓振飛

(1.裝甲兵工程學院 兵器工程系，北京 100072;2.北京京航計算通訊研究所，北京 100074)

火炮身管的柔性對射擊精度影響很大，彈丸在膛內運動時，膛線賦予彈丸旋轉運動，以保證彈丸飛行穩定。彈丸膛內運動時，彈帶和前定心部與炮膛內壁的接觸碰撞會導致身管激烈振顫，炮口振動引起彈丸跳角，直接影響射擊命中率。文獻[1]建立了彈丸膛內運動引起炮管振動的方程，沒有考慮彈丸與身管的碰撞－接觸;文獻[2]基于多體接觸理論，考慮了彈丸彈帶、前定心部與膛壁的碰摩，但其有限元法本質上仍將彈性身管簡單轉化為剛體接觸;文獻[3]采用等參六面體縮減積分單元，考慮了彈丸彈帶和前定心部與膛壁的碰撞及膛線對彈丸的扭轉作用，對彈丸在膛內的運動過程進行了仿真計算，但膛線身管結構復雜，很難建立準確的三維模型;文獻[4]研究了彈管間隙和加速移動質量對身管振動的影響，采用懸臂梁的振型作為身管的近似模態，無法得到身管的真實模態。由于小口徑火炮身管短、膛壁薄，適合采用有限單元法編程求解。本文利用Dubowsky模型，基于有限單元法，將彈帶定為彈塑性體，前定心部視為剛體，身管為柔體，分別研究了彈帶凹槽部、凸起部和前定心部對身管的瞬態沖擊作用，根據彈帶行進過程確定身管各單元在各時間步長內所受的載荷譜，最終建立了彈丸膛內運動引起的身管動態響應模型。

1 身管有限元模型

1.1 坐標系的建立

根據結構動力學可知，身管屬于內中空懸臂梁的空間軸對稱問題，采用圓柱坐標系比較方便。建立坐標系ozρθ，o為以炮膛軸線與膛線起始部端面的交點，oz軸為未變形的身管軸線，oρ軸為身管徑向，以o垂直向上為θ角零點，向右為起始方向，如圖1所示。

圖1 圓柱坐標系下的身管結構Fig.1 Tube structure in cylindrical coordinate system

身管有n條膛線，膛線在oρθ平面內展開曲線為f(z);身管內壁半徑為 r1，外壁半徑為r2，膛線纏度為η，深度為d，陽線寬為a，彈丸彈帶寬為b，前定心部距彈帶前端距離為c，射角為α;在t時刻，彈丸質量M在膛內軸向運動速度為v(t)，沿膛線旋轉速度為ω(t)，彈丸膛內運動距離為s(t)=，均由內彈道計算得出。

圖2 有限單元劃分Fig.2 Division of finite element

1.2 有限單元劃分

身管單元包括兩部分，陽線部分和其余部分。如圖2(a)所示。首先取陽線進行網格劃分，軸向以彈帶寬b為單元長度的六面體網格，如圖2(b)所示。在ozρθ坐標系中，將陽線去掉，對身管進行離散，取軸向長度b，夾角θ，徑向dr的六面體網格單元，如圖2(c)所示。兩種單元通過公用面的節點進行載荷傳遞。它分為兩種情況:陰線內壁的六面體網格單元直接承受彈帶凸起部載荷作用;與陽線單元直接連接及不與彈帶接觸的身管內部結構單元，其節點載荷為由陽線節點傳遞的載荷或陰線內壁六面體網格單元節點傳遞的載荷。該六面體網格單元在圓柱坐標系中與三維六面體等參單元的局部坐標系中規則的母單元等效，采用八節點L型插值函數N=[N1N2N3N4N5N6N7N8]作為單元位移函數，單元節點位移向量為={u1，…，u8}，其中ui為三個元素的列向量。則單元撓度為:

1.3 單元編碼

對于變速彈丸依托膛線旋轉前進的身管結構，單元的編碼非常重要。結合圓柱坐標系和時間軸，建立四維向量為單元編碼。設某單元編碼為 ei，j，k，t，其中，i，j，k分別代表整數表示的周向、徑向和軸向的遞進關系，t為實數，表征在t時刻與彈丸發生觸碰的單元。

對于i，由于彈丸彈帶的凹槽部和凸起部在內膛圓周對各單元進行觸碰，首先對彈帶各部編號，i={0，1，2，…，2n－1}，以θ角零點為0，順時針遞進。彈丸沿膛線運動時，與彈帶各編號對應的陽線和陰線在i上不變;

2 彈炮耦合過程中身管動態響應模型

2.1 彈丸作用下身管的徑向動力學方程

將身管簡化為膛線起始部一端固定的變截面彈性懸臂梁[5]，根據 Euler-Bernoulli梁理論，假設橫向剪切應變和橫向線應變為零，身管彎曲方程為[6]:

其中:f(z，t)為身管所受的接觸力，EI(z)為身管(梁)的抗彎剛度，h為身管速度衰減系數，由身管的結構和材料性能決定。u(z，t)為身管z位置處t時刻的撓度，Me為身管在單元e的質量分布。

初始條件為:

邊界條件為:

將彈丸看作移動質量－彈簧－阻尼，彈簧和阻尼分別通過彈帶凸起部和前定心部與身管發生觸碰。剛度系數和阻尼系數根據接觸體的材料彈性模量和幾何形狀合適選取。如圖3所示，彈丸質量M以速度v(t)在膛內運動，通過彈簧K、阻尼C與身管(梁)相互作用，陰線面所受彈丸彈帶或前定心部載荷為以身管軸線即懸臂梁中性軸為起點的徑向力，大小為:

其中:δ為Dirac delta函數;s(t)為t時刻彈丸膛內運動距離，s(t)=v(t)dt，由內彈道計算得出;ρ(t)為彈丸徑向位移為彈簧剛度系數為阻尼系數，塑性變形后的彈帶與內膛的接觸剛度根據銅－鋼的材料彈性模量和幾何形狀合適選取，前定心部與內膛的接觸剛度根據鋼－鋼的材料彈性模量和幾何形狀合適選取，阻尼系數采取工程上的經驗值。

圖3 彈丸移動彈簧－阻尼－質量系統作用下身管懸臂梁模型Fig.3 Tube cantilever beam model under the moving spring，damping＆ mass system of projectile

彈丸(彈帶或前定心部)徑向運動方程為:

整理得:

由于離散體系的平衡方程是在節點上進行分析的，因此對身管劃分的三維單元與平面微分方程之間主要依靠單元的等效節點載荷來銜接，根據靜力等效原則，將作用在身管(梁)上的非節點面力載荷轉化為等效節點載荷，三維單元的面力等效節點載荷為:

根據彈丸與身管相互作用的不同，將彈丸膛內運動過程分為三個部分:階段一為彈丸彈帶由坡膛擠進膛線過程;階段二為直膛段彈丸彈帶和前定心部與身管共同作用過程;階段三為直膛段彈丸彈帶單獨作用過程。如圖4所示。

圖4 移動彈丸各階段膛內運動過程Fig.4 Moving process in chamber of projectile

2.2 階段一身管內壁單元受力分析

在階段一，彈丸彈帶前端面由0到2b。當彈帶完全嵌入膛線時，彈帶首先為彈性變形，逐步達到屈服極限后過渡到塑性變形，當到達2b點時，彈帶塑性變形結束[7]。因本文主要對彈炮系統的徑向應力應變進行分析，彈帶擠進階段彈丸對膛線的沖擊主要體現在彈性變形階段[8]，在此階段，可采用彈性變形階段的彈簧－阻尼－質量系統模擬彈丸與身管的相互作用。

將曲面內壁簡化為平面，當彈丸運動到單元e時，聯立式(2)、(3)、(6)，此單元的外載荷向量為:

式中:Me為單元質量。

式中:Δl為彈帶半徑強制量。

彈丸沒有直接作用的身管單元外載荷為自重:

2.3 階段二身管內壁單元受力分析

在階段二，彈丸彈帶前端由2b到2ηr1－c。彈帶凸起部與陰線膛壁為彈性接觸，可采用彈簧－阻尼－質量系統進行模擬，將膛線定義為主面，彈帶定義為從面，運用單純主從搜索算法確定接觸點位置;彈帶凹槽部與陽線的徑向受力為正壓力。彈帶在軸向受到彈底壓力作用、彈前空氣阻力作用及彈丸和身管單元的摩擦阻力，而本文的彈炮系統模型描述的是圓柱坐標系下身管徑向的受力及運動特征，故不考慮彈底壓力和摩擦阻力作用。由于彈管間隙的存在，彈丸前定心部與陽線產生彈性碰撞，之后分離與另一部分陽線產生碰撞，如此反復。

圖5 陽線受力分析Fig.5 Force analysis of rifling

陽線受力情況如圖5所示，陽線導轉側壓力使身管產生左旋回轉力矩，對徑向變形無影響，故不考慮。因彈帶材料采用理想彈塑性線性強化模型，在塑性變形后應力應變仍呈斜率很小的線性關系。為簡化計算，取屈服應力σs為彈帶凹槽部等效應力。

故徑向正壓力為:

彈帶凸起部與陰線單元e的觸碰力模型與式(5)彈帶凸起部相同。

彈帶凹槽部與單元e的觸碰力模型為:

彈丸前定心部與單元e的觸碰力模型與式(5)相同，接觸分為兩部分:初始接觸采用單純主從搜索算法，膛線定義為主面，前定心部定義為從面，搜索初始接觸單元;前定心部與膛壁某單元碰撞后彈回，經過兩個單側間隙B、軸向位移和彈丸旋轉角度，與下一單元碰撞[9]。設t時刻與彈丸前定心部碰撞的陽線單元為ei'，j，k'，t，則下一次與彈丸前定心部碰撞的陽線單元為ei'+1，j，k'+1，t，其中，

式中:B為前定心部與陽線單側間隙。

2.4 階段三受力分析

在階段三，彈帶與身管接觸模型與階段二相同。彈丸前定心部脫離炮口，與陽線不再接觸碰撞。故階段三彈炮耦合數值模型為式(7)與式(9)。

2.5 動力學系統建模及求解

為便于數值求解，綜合以上三個階段，編寫MCK矩陣[10]，建立系統動力學方程:

因矩陣[M]、[C]和[K]中有部分是時變的，需要以分塊矩陣形式表示，故:

式(15)中，分塊矩陣中的上標e即表示了在時間t參與接觸的時變單元，附加的剛度、阻尼只與移動荷載直接作用的單元位移有關。

彈丸膛內運動屬于移動荷載問題，MATLAB軟件內部函數中有專門針對移動荷載問題的ODE求解器[11]。在MCK矩陣基礎上，基于四階五次Ruge-Kutta微分方程進行二次開發，對彈炮耦合系統進行數值求解。

MATLAB提供的ODE45的調用格式為:

其中:odefun為指示微分方程等號右側的函數句柄，在迭代中計算MCK矩陣中的時變分塊。變量tspan為計算范圍，變量y0為系統初始狀態變量的值。options可以定義函數運行時的參數。T為時間點的列向量，Y為與時間T對應的求解數組。

3 數值計算算例

取 θ=11.25°，M=84.2 kg，m=0.36 kg，r1=0.015 m，r2=0.025 m，d=0.000 3 m，a=0.002 m，b=0.001 5 m，c=0.08 m，η =30，B=0.000 5 m，E=300 GPa[12]，h=0.2[12]，dr=0.002 m，對于彈帶，對于前定心部，根據內彈道計算得到的彈丸膛內運動速度和旋轉速度，采用2.5節中數值求解方法，應用MATLAB生成單元信息并進行彈炮耦合動態響應計算，生成的單元個數為105605，節點個數為115 360，其中陽線單元為9 605。取炮口外壁周向θ=0°單元為研究對象，求解其徑向(即垂直向)位移。

實彈射擊時，選用技術狀態良好的武器系統進行炮口振動測試，高低、方位射角均為0°，彈種為曳光穿甲彈，單發射擊，在常規靶場進行。在炮口端正上方安裝3 000 g加速度傳感器，應用NI公司高分辨率數據采集卡4472 B，在PXI－8105平臺上采用LabVIEW語言編寫測試框圖程序，測試與研究對象位置基本對應的單元的垂直向振動位移，總共獲得了12組數據。12發彈的地面散布較小，12組振動數據的差異也非常小，由于彈丸出炮口時間較短(約6 ms)，該自動炮身管仍處于自由后坐階段，未壓縮緩沖器彈簧，從而與搖架、炮塔、車體等的相互作用并不明顯，對測試曲線中前7 s的幅值較低部分取均值并與仿真位移曲線對比如圖6所示。

由圖6可以看出，兩者非常接近，對兩者進行相關性分析，corrcoef=0.895，彈炮耦合數學模型較為準確。

分別對炮口內壁周向 θ=0°、θ=90°、θ=180°、θ=270°的單元進行分析，得到彈丸膛內運動時，身管的撓度曲線如圖7所示。

由圖7可以看出，彈丸繞膛線在膛內運動一周，對內壁各方向觸碰規律相似，但由于前定心部的碰撞，各方向又有所差別[13]。改動數學模型，不考慮階段二中前定心部的碰撞規律，取炮口內壁周向θ=0°為研究對象，得出振動曲線如圖8所示。

圖6 仿真模型的實測曲線對比驗證Fig.6 Comparison of simulated model and tested model

圖7 炮口內壁周向各單元撓度曲線Fig.7 Deflection curve of elements inside the muzzle

由圖8可以看出，不考慮前定心部的碰撞時，身管振動幅度和頻率都較小，可知前定心部的碰撞對彈炮耦合影響很大。

取炮口端面θ=0°處內外壁單元為對象，時間歷程內炮口振動曲線如圖9所示。

圖9 炮口內壁與外壁單元振動曲線對比Fig.9 Comparison of elements vibration inside the muzzle and outside the muzzle

由圖9可以看出，與炮口內壁相比，外壁單元的振動頻率和幅值都較低，表明與彈丸直接作用的炮口內壁對彈丸影響更大。

4 結論

考慮彈丸彈帶擠進過程和剛體前定心部對內膛的瞬態碰撞對身管振動的影響，將彈帶凹槽部視為剛塑性體，凸起部應用彈簧－阻尼－質量系統，身管簡化為柔性等截面懸臂梁，對從彈帶擠進膛線開始到彈帶出炮口結束的彈丸膛內運動對身管陽線和陰線單元的觸碰作用進行了分析。由于計算量大，比較適合小口徑火炮。以某型火炮為例進行了數值計算，在實驗驗證的基礎上對炮口振動進行了分析，得出以下結論:

(1)彈炮耦合過程對炮口橫向振動和縱向振動影響規律相當。在一發彈的大約6 ms的膛內運動過程中，炮口各周向單元的徑向位移變化規律相同，相對于沖擊碰撞，彈丸重力因素影響較小;

(2)彈丸前定心部對炮口振動影響很大，使振動幅值和頻率都增大。由于彈丸膛內高速運動的章動性，使彈丸前定心部趨于不碰撞，而本文的數學模型和圖8證明了由于入膛時隨機碰撞引起的一系列規律碰撞是存在的。

本文建立了適合小口徑火炮的較為精確的彈炮耦合數學模型，為建立發射動力學模型和射擊精度分析提供了彈炮方面的理論基礎，并對其他位置的碰撞分析具有一定的借鑒意義。

[1]周 叮，謝玉樹.彈丸膛內運動引起炮管振動的小參數解法[J].振動與沖擊，1999，18(1):76 －81.ZHOU Ding，XIE Yu-shu.Solution of barrel vibration caused by motion of bullets in bore by use of small parameter method[J].Journal of Vibration and Shock，1999，18(1):76 －81.

[2]劉 雷，陳運生，楊國來.基于接觸模型的彈炮耦合問題研究[J].兵工學報，2006，27(6):985 －987.LIU Lei，CHEN Yun-sheng，YANG Guo-lai.A study on the projectile-barrel coupling based on contact model[J].Acta Armamentarii，2006，27(6):985 －987.

[3]葛建立，楊國來，陳運生，等.基于彈塑性接觸/碰撞模型的彈炮耦合問題研究[J].彈道學報，2008，20(3):103－106.GE Jian-li，YANG Guo-lai，CHEN Yun-sheng，et al.A study on projectile-barrel coupling problem based on elastoplastic contact/impact Model[J]. Journal of Ballistics，2008，20(3):103－106.

[4]劉 寧，楊國來.彈管橫向碰撞對身管動力響應的影響[J].彈道學報，2010，22(2):67 －70.LIU Ning，YANG Guo-lai.Effect of lateral impact between projectile and barrel on dynamic response of tube[J].Journal of Ballistics，2010，22(2):67 －70.

[5]陳 強，楊國來，王曉鋒，等.變速移動彈簧阻尼質量系統作用下梁的動態響應[J].動力學與控制學報，2010，8(3):254－257.CHEN Qiang， YANG Guo-lai， WANG Xiao-feng，et al.Dynamic responses of beam subjected to a moving oscillator with non-uniform velocity[J].Journal of Dynamics and Control，2010，8(3):254 －257.

[6]Nakalswarny K K.Finite element analysis of a projectile during gun launch[D]. BachelorofEngineering in Mechanical Engineering B.I.E.T，Kuvempu University，Karnataka，2002.

[7]孫河洋，馬吉勝，張高明，等.彈帶擠進過程中膛線的應力分析[J].軍械工程學院學報，2008，20(4):43 －46.SUN He-yang， MA Ji-sheng， ZHANG Gao-ming， et al.Analysis of the land stress during the engraving process[J].Journal of Ordnance Engineering College，2008，20(4):43－46.

[8]彭 濤，王學軍，黃善文.基于ANSYS的彈帶擠進變形及應力分析[J].艦船電子工程，2009，29(11):156 －159.PENG Tao，WANG Xue-jun，HUANG Shan-wen.Analysis on material deformation and stress in the process of driving band engraving based on ANSYS[J].Ship Electronic Engineering，2009，29(11):156 －159.

[9]吳 宏.彈丸膛內外運動及其對射擊精度影響研究[D].南京:南京理工大學，2001.

[10]王少欽，夏 禾，郭薇薇，等.變速移動荷載作用下簡支梁橋的動力響應及共振分析[J].振動與沖擊，2010，29(2):26－28.WANG Shao-qin，XIA He，GUO Wei-wei，et al.Dynamic response and resonance analysis for a simply-supported bridge under speed-varying loads[J].Journal of Vibration and Shock，2010，29(2):26 －28.

[11]Liu K，Reynders E，Roeck G D，et al.Experimental and numerical analysis of a composite bride for hig-speed trains[J].Journal of Sound and Vibration，2009，320:201 － 220.

[12]史躍東，王德石.考慮慣性效應的移動彈丸作用下身管振動特性[J].兵工學報，2011，32(1):414 －419.SHIYue-dong， WANG De-shi. Study on vibration characteristicsofbarrelsubjected to moving projectile considering inertia effect[J].Acta Armamentarii，2011，32(1):414－419.

[13]黃建亮，陳樹輝.縱向與橫向振動耦合作用下軸向運動梁的非線性振動研究[J].振動與沖擊，2011，30(8):24－27.HUANG Jian-liang，CHEN Shu-hui. Study on nonlinear vibration of an axially moving beam with coupled transverse and longitudinal motions[J].Journal of Vibration and Shock，2011，30(8):24 －27.