朱光濤,馬永忠,馬志偉,羅忠新
(武警工程大學,西安 710086)
基于Matlab的某型橡皮動能彈外彈道仿真研究*
朱光濤,馬永忠,馬志偉,羅忠新
(武警工程大學,西安 710086)
為了研究某型橡皮動能彈彈丸的運動特性,文中利用Matlab仿真軟件對彈丸的運動過程進行仿真,通過對仿真結果的分析,得到了彈丸速度隨時間和距離變化的規律以及彈丸在不同距離上對有生目標的作用效果,解決了彈丸的飛行動能受空氣阻力和重力影響,不同射擊角度彈丸飛行動能衰減也不一樣的問題,驗證了彈丸在初速度為100 m/s時能夠滿足15~50 m的有效作用距離指標,為下一步優化改進提供理論基礎。
Matlab;動能彈;外彈道仿真
動能彈是利用彈丸的飛行動能打擊有生目標,使之致傷致痛,從而失去抵抗能力或行動受到抑制的一種警用非致命性彈種。動能彈的致傷物質——彈丸,由于彈藥本身非致命的要求,多采用密度較小的柔性材料制成,如橡膠、塑料、木材等。在研究某型橡皮動能彈時,彈丸的飛行動能會受到空氣阻力和重力影響,不同射擊角度彈丸飛行動能衰減也不一樣[1],為解決這些問題,確保某型橡皮動能彈的非致命效果,文中對該型橡皮動能彈進行外彈道仿真。
使用Matlab仿真軟件對彈丸的運動過程進行仿真,對彈丸進行受力分析是十分必要的。彈丸在飛行過程中受到的力主要是空氣阻力和重力,重力加速度為一個常量即g=9.8 kg/m2,所以下面主要對空氣阻力進行分析[2-3]。
1)空氣阻力與彈丸特性、空氣特性以及彈丸和空氣之間相對運動特性三方面有密切的關系。根據有關試驗研究可知,空氣阻力一般計算公式為:
(1)
式中:Rx為空氣阻力;S=πd2/4為彈丸特征面積(m2),d為彈丸的最大直徑;ρ為當地的空氣密度(kg/m3);Cx0(v/c)為攻角δ=0時彈丸的阻力系數(阻力系數Cx0(v/c)的值可通過查表獲得)。
2)當空氣對彈丸的阻力R不與彈軸也不與速度v矢量共線反向時,彈軸與速度方向存在的夾角δ即為攻角,攻角δ使得阻力R分為作用效果截然不同的兩個分力,即切向阻力Rx和升力Ry。
①切向阻力Rx的作用效果是使得速度v減小,其計算公式為:
(2)
當攻角δ較小且不跨音速時,攻角不為零時的阻力系數Cx與攻角為零時的阻力系數Cx0之間存在以下關系:
Cx=Cx0(1+Kδ2)
(3)
同時容易證明阻力的攻角系數K為:
(4)

綜上所述,可以將切向阻力Rx的計算公式寫為:
(5)
式中:Kx(M,δ)為切向阻力特征數。
②升力Ry作用效果是使得速度v改變方向,其計算公式為:
(6)

兩者之間的關系式如下:
(7)
當攻角較小時:
(8)
(9)
式中:Ky(M)為升力特征數;l為全彈長。
文中研究的橡皮塊彈速度較低,故忽略空氣阻力中的渦流阻力和波動阻力,只考慮摩擦阻力,同時由于彈丸形狀的不規則性,空氣阻力對彈丸產生扭轉力矩,飛行過程中必然會發生翻轉,如果考慮飛行中的攻角δ變化,則會導致計算和分析變得非常復雜,為了研究彈丸各參數變化對彈丸速度衰減和運動軌跡的影響,將前述模型進行適當簡化:以單個子彈丸為研究對象,不考慮彈丸的翻轉,視其運動為質點運動,選取扇形面為迎風面,并忽略紙質包覆體撕裂對彈丸飛行的影響,各彈丸運動特性完全相同。
2.1 彈丸運動方程的建立
彈丸運動過程中,主要受空氣阻力和重力的作用,假設彈丸初速為v0,飛行方向與水平面的夾角為α,彈丸質量為m,以槍口為坐標原點建立x-y平面坐標系(如圖1所示),具體分析如下:

圖1 平面坐標系
在y軸方向上,彈丸受到重力和空氣阻力的作用,由牛頓第二定律可知:
ma=(-1)n·mg-F阻
(10)

(11)
彈丸的運動過程可以分為兩種情況[4]:
1)當彈丸偏向上飛行時,即當y軸存在正軸方向的速度分量時,飛行角度0≤α≤π/2,重力與空氣阻力同向,取n=1。
2)當彈丸偏向下飛行時,即當y軸存在負軸方向的速度分量時,飛行角度-π/2≤α≤0,重力與空氣阻力反向,取n=2。
由于彈丸在x軸方向上只受到空氣阻力的作用,則有:
(12)
式中:vx為在x軸方向上的速度分量;vy為在y軸方向上的速度分量;S′為彈丸的迎風面積;Cx0為彈丸的阻力系數;ρa為當地空氣密度;g為重力加速度。

對于彈丸運動過程的情況1)而言,vy先減小為零,之后又在y軸負方向上做自由落體運動,因此對其運動過程分上升和下降兩種情況考慮:
上升:
(13)

其上升到最高點時,vy1=0,對應,即0
(14)
同理,對于彈丸運動過程的情況2)而言:
(15)

x軸方向上水平速度分量的通解為:
(16)

基于速度與位移之間的關系,則由水平速度分量vx和垂直速度分量vy,即可求得彈丸在x軸和y軸上任意時刻的位移,所以:
對于彈丸運動過程的情況1)而言,其運動軌跡方程為:
(17)
對于彈丸運動過程的情況2)而言,其運動軌跡方程為:
(18)
通過方程(13)~(18),可以計算出彈丸出槍口后任何時刻的速度及其運動軌跡。
2.2 計算方法的選取
選用Matlab軟件中的內置程序ode45函數作為求解算法[5-6],ode45函數又稱Runge-Kutta-Fehlberg法的自適應步長算法。該方法使用4階和5階Runge-Kutta公式來改進解并且監測計算精度,同時通過自適應改變的步長h減少其求解的總步數,提高了計算效率,完全可以滿足彈丸運動過程中仿真計算的需要。
根據西亞切阻力定律可知:
當v=50~150 m/s時,Cx0=0.255;當v=150~250 m/s時,Cx0=0.259。
彈丸初速設定為100 m/s,所以選取阻力系數Cx0=0.255。選取4個不同的發射角度,即彈丸的飛行方向與水平面的夾角α取4個不同值,對彈丸的速度衰減及運動軌跡的規律進行研究,其計算參數如表1所示。

表1 計算參數
將參數代入計算模型,對彈丸在2 s內速度衰減以及在y軸方向1.5 m內的運動軌跡進行比較分析,結果如圖2所示。

圖2 發射角不同時彈丸的速度衰減和運動軌跡
由圖2(a)可知,發射角α=45°時,彈丸速度衰減最慢,同時發射角α=30°的速度衰減規律與α=45°的速度衰減規律相似,兩條曲線基本重合,發射角α=0°時,速度衰減最快,當t=0.8 s左右時,發射角α=0°的速度衰減變為最慢,發射角α=15°的彈丸速度衰減變為最快,t=2 s時,發射角α=0°的速度保持最好v=45 m/s,其余發射角速度保持相似v=40 m/s。
由圖2(b)可知,取落高h=1.5 m,發射角α=15°時,彈丸運動軌跡先向上運動大約60 m后才開始下落的射程最遠超過100 m;發射角α=30°時,彈丸運動軌跡先向上運動大約50 m后才開始下落的射程超過90 m;發射角α=45°時,彈丸運動軌跡先向上運動大約40 m后才開始下落的射程為80 m,發射角α=0°時,彈丸運動軌跡直接向下,射程最近為48 m左右。
由此可知,射擊姿勢確定后,平射彈丸的射程最近,但是落地彈丸的速度還很高,浪費了相當多的能量;以小角度射擊時,彈丸的射程最遠,可以達到100 m以上,且速度保持較好。應根據處置情況的不同,確定最佳的射擊角度。
由前述彈丸運動方程可以得到不同發射角對應的速度隨距離變化曲線,如圖3所示。

圖3 不同發射角對應的彈丸速度隨距離變化曲線
由圖3可知,發射角α=45°時,彈丸速度隨距離衰減最快,隨著射角的減小,彈丸速度隨距離衰減的速度減慢;發射角α取0°、15°和30°,且水平運動距離小于50 m時,彈丸的速度隨距離變化曲線基本重合,50 m之后速度差異逐漸增大。
分析可知,發射角α≤30°時,且水平運動距離D≤50 m,彈丸速度隨距離衰減規律相似,可以使用統一標準確定彈丸對目標打擊的非致命性。
根據GJBZ 20262—95防暴動能彈威力標準中的規定,防暴動能榴(霰)彈的非致命性標準,即動能5 J≤Ed≤120 J和比動能0.5 J/ cm2≤ed≤12 J/cm2。在已知的初始條件下:質量m=5.5 g,扇形面S1=πr2/4=3.14×1.8×1.8÷4=2.54 cm2,側面S2=l×d=1.8×1.8=3.24 cm2,弧面S3=πrd/2=3.14×1.8×0.9=5.09 cm2。則可求出動能、比動能對應的速度要求。
在動能標準下的速度要求:42.5 m/s≤v≤208.9 m/s
在比動能標準下速度要求:扇形面為21.5 m/s≤v1≤105.2 m/s;側面為24.3 m/s≤v2≤118.9 m/s;弧面為30.4 m/s≤v3≤149.0 m/s。
當射角α取0°、15°和30°時,15 m處的速度為90 m/s,完全滿足v1≤105.2 m/s的小于中度損傷的標準,50 m處的速度為70 m/s,完全滿足v3≥42.5 m/s的大于無損傷的標準。當射角α=45°時,15 m處的速度為85 m/s,完全滿足v1≤105.2 m/s的小于中度損傷的標準,50 m處的速度為65 m/s,完全滿足v3≥42.5 m/s的大于無損傷的標準。所以取彈丸初速v0=100 m/s能夠滿足防暴彈動能打擊的非致命性要求。
文中選用Matlab軟件中的內置程序ode45函數作為求解算法,對橡皮動能彈彈丸的運動過程進行了仿真,進而對其運動特性進行了計算分析,得出了不同射角條件下彈丸速度衰減的規律,以及彈丸速度隨距離變化的規律,并結合動能彈打擊的非致命性要求,驗證了彈丸在初速度為100 m/s,能夠滿足15~50 m的有效作用距離指標,確保了彈丸打擊的非致命性。
[1] 韓子鵬. 彈箭外彈道學 [M]. 北京: 北京理工大學出版社, 2008: 20-38.
[2] 劉琨. 子母彈子彈拋撒散布動力學分析 [D]. 南京: 南京理工大學, 2007.
[3] 浦發. 彈道系數與空氣阻力定律問題的探討 [J]. 兵工學報彈箭分冊, 1981(3): 1-12.
[4] 肖亞, 郭三學, 歐陽的華, 等. 防暴彈破片速度衰減規律 [J]. 遼寧工程技術大學學報(自然科學版), 2013, 32(5): 619-622.
[5] 王能超. 計算方法簡明教程 [M]. 北京: 高等教育出版社, 2004: 80-98.
[6] 馬利兵, 林都. 基于MATLAB的外彈道模型仿真研究 [J]. 中北大學學報, 2006, 27(5): 412-415.
Simulation of Certain Type of Rubber Kinetic Energy Projectile Trajectory Based on Matlab
ZHU Guangtao,MA Yongzhong,MA Zhiwei,LUO Zhongxin
(Engineering University of China Armed Police Force, Xi’an 710086, China)
In order to study motion characteristics of a type of rubber bullet kinetic projectile, Matlab was used to simulate projectile motion, according to analysis of simulation results, both law of projectile velocity changing with time and distance and effects on animate target at different distance were obtained, these solve kinetic energy of projectile when flight affected by gravity and air resistance, and different projectile kinetic energy attenuation at different shooting angle, and it is verified the initial velocity 100 m/s of the projectile meets 15~50 m effective distance indicator, this provides a theoretical basis for further optimization improvement.
Matlab; kinetic energy projectile; trajectory simulation
2015-12-24
朱光濤(1991-),男,江西吉安人,碩士研究生,研究方向:非致命性武器。
TJ012.3
A