0 引言
制导炮弹具有远程穿透、高速飞行、精准打击的优点,在现代化战争中扮演着至关重要的角色[1-3]。远程制导炮弹飞行过程中转速高,控制难度大,如何降低转速,提高测量精度、高效控制是当前研究的热门问题[4]。
北京理工大学的陈庆森等[5]利用多刚体动力学的方法研究了封控弹弹伞系统末弹道的动态特性,使封控弹末弹道满足设计要求,减旋叶片能使封控弹减旋到理想值。南京理工大学的逯麒等[6]提出了一种适用于电磁发射超高速制导炮弹的滚转运动稳定区域和减旋策略,主动起旋与减旋止摆控制策略为电磁发射超高速制导炮弹飞行控制提供新颖可行的方法。王中原等[7]针对阻力环叠加减旋翼片的二维弹道修正技术方案,进行横向弹道修正能力及飞行稳定性分析,并通过炮射试验验证该方法的准确性和有效性。邱荷等[8]根据典型末敏弹的结构特点及减速减旋段的运动特性,提出了一种利用Matlab和VR技术对其进行可视化仿真研究的方法。中北大学的刘一鸣[9]、祝敬德[10]等设计了一种适用于弹载环境的减旋控制系统,采用半捷联模拟控制减旋电路,有效提高了测量精度。
本文中描述一种利用船尾装置进行弹体减旋的制导炮弹外形,这种减旋船尾采用轴承与船尾和前端非硬连接的装配方式,船尾前端的轴承连接装置,能够保证弹体出膛后达到减旋效果。从而降低弹体的转速,达到可控的目的。
1 物理模型
首先对具有减旋船尾的制导炮弹气动外形进行设计,采用三维建模软件进行结构件设计及装配,然后将建好的炮弹模型导入GAMBIT中。本文中研究的制导炮弹由弹头、战斗部、发动机、密封圈、船尾、闭气盖等部分组成,制导炮弹整体气动外形结构如图1所示。
2 数值方法
2.1 工程计算方法
在空气动力学理论的基础上进行工程计算,并对炮弹结构和流场进行适当的假设和简化,最后进行求解运算。外弹道研究方法必须满足计算速度快,计算简捷,计算精度高3个特点。所以,工程计算方法被广泛应用在外弹道研究中。相对于外形而言,工程计算方法可以使阻力、升力、翻转力矩等系数的计算结果和正常值的误差控制在10%以内,甚至更小。所以,工程计算方法比较适用于弹体结构变化对气动特性的影响,并且通过计算可以获得炮弹的气动参数等信息。
图1 炮弹三维模型
Fig.1 Three-dimensional model of missile
根据文献[11]中所推导,炮弹组件的阻力系数表达式为:
CxBW=1.1(CxOBW+CxiBW)
(1)
升力系数为:
(2)
其中:CxOBW为弹体的零升阻力系数;CxiBW为尾翼的诱导阻力系数;CyB为弹体的升力系数;CyW为尾翼的升力系数;SW为尾翼的表面积;SM为弹体的表面积;ψ为尾翼对数系数(本文中ψ=0.75);
为弹体与尾翼的干扰因子修正系数。
2.2 控制方程
连续性方程和动量方程[11]为:
(3)
(4)
通常使用Boussinesq假设来求解能量方程,一般认为雷诺应力和平均速度梯度成正相关,其表达式为[12]:
(5)
(6)
2.3 湍流模型
湍流作为流体的一种流动形式,通常是无序的、多角度、不规则的流动。一般而言,湍流具有良好的扩散特性。在空气动力学中,湍流一般是指某特定区域、较小时间内的风速波动。湍流的物理模型由不同大小的漩涡组成,并且漩涡大小、流向等特征都是无序的,无法用数学公式表示出来。
在空气动力学的研究中Spalart-Allmaras模型、k-ε模型、k-ω模型和雷诺应力模型经常使用,本文中采用Spalart-Allmaras模型。Spalart-Allmaras模型是用于航空领域的较为简便的单方程模型,用于解决墙壁束缚(wall-bounded)流动方面的问题,在解决逆压梯度的边界层问题展现出很好的效果。在空气动力学中,一般使用该模型解决飞行器等绕流流场计算问题。
(7)
其中:
为湍流运动黏度;Gv为湍流黏度的增加项;Yv为湍流黏度的减少项;Vv为分子运动的黏度;
为用户自定义源项。
2.4 网格划分及边界条件
为了精确获得阻力系数、升力系数和力矩系数等的计算值,在GAMBIT划分网格时,需要在弹体附近布置较密的网格。针对弹体计算域进行建模时,也需要格外注意,一般选取计算域为大圆柱体和小圆柱体2个。对大圆柱体而言,长度为弹体长度的7.5倍,直径为弹体直径的40倍;小圆柱体即为内部加密区,长度为弹体的2.5倍,直径为弹体直径的5倍,如图2所示。在网格划分时越靠近弹体的网格越密,远离弹体的网格可以相对稀疏一些,其中弹头部分的网格越密,计算结果越准确。为了加快计算机运算速度,一般特地将计算域分割成两个部分,只对一半进行仿真计算。
图2 计算域和网格划分示意图
Fig.2 Computational domain and grid partition diagram
模型网格划分模型如图3所示,整个计算域采用了网格加密的办法,使得计算更加精确。整个计算域网格数在100万以上。将划分好网格的炮弹流场域导入Fluent中进行计算。
图3 弹体周围网格划分示意图
Fig.3 Diagram of meshing around projectile
设置边界条件的时候选用密度基耦合显式求解器进行求解,对来流和物面分别选用远场边界条件和黏性边界条件:
(8)
其中:ui=(i=1,2,3)为坐标方向的速度分量; F为压力; ρ为密度;n为法线方向。
3 计算结果分析
使用Fluent软件对减旋弹体进行气动仿真分析时需要计算不同攻角和不同马赫下的工况值,分别选取马赫数为1 Ma,1.5Ma,2Ma,2.5Ma,3Ma,3.5 Ma,4Ma;攻角α取1°,2°,3°,4°;共28种状态,对应升力系数、阻力系数、力矩系数、X方向速度分布云图、全弹表面压力分布云图等计算仿真分析如下。
3.1 动态监控窗口分析
以4°攻角,2.5Ma为例,经过Fluent 计算,得到以下阻力系数,升力系数和力矩系数等。通过观察可知,经过3 000次迭代之后发现,各系数均变化不大,说明此时处于收敛状态。
500次迭代完成之后, 4°攻角、1.5Ma时阻力系数的收敛值是0.174 0;4°攻角、2.0Ma时阻力系数收敛值是0.255 1;4°攻角、2.5Ma时阻力系数收敛值是0.324 6; 4°攻角、3.0Ma时阻力系数收敛值是0.466 5。通过3 000次迭代发现,阻力系数的变化曲线图在200次迭代,基本呈现直线状态,即此时在经过迭代阻力系数的值基本不发生变化,说明阻力系数的计算值达到收敛状态。将4°攻角、1.5Ma,4°攻角、2 Ma,4°攻角、2.5 Ma,4°攻角、3 Ma这4种不同马赫的阻力系数值在Excel表格中对比发现,在4°攻角一定的情况下,随着马赫数的升高,阻力系数值也在不断升高。阻力系数随迭代过程变化曲线如图4所示。
500次迭代完成之后, 4°攻角、1.5Ma时升力系数的收敛值是0.066 99;4°攻角、2.0Ma时升力系数收敛值是0.136 2;4°攻角、2.5Ma时升力系数收敛值是0.230 3;3.0Ma时升力系数收敛值是0.408 8。通过3 000次迭代发现,升力系数的变化曲线图在200次迭代,基本呈现直线状态,即此时在经过迭代升力系数的值基本不发生变化,说明升力系数的计算值达到收敛状态。将4°攻角、1.5Ma,4°攻角、2Ma,4°攻角、2.5Ma,4°攻角、3Ma这4种不同马赫的升力系数值在Excel表格中对比发现,在4°攻角一定的情况下,随着马赫数的升高,升力系数值也在不断升高。升力系数随迭代过程变化曲线如图5所示。
图4 阻力系数随迭代过程变化曲线
Fig.4 Curve of resistance coefficient with iterative process
图5 升力系数随迭代过程变化曲线
Fig.5 The variation curve of lift coefficient with iteration process
500次迭代完成之后, 4°攻角、1.5Ma时力矩系数的收敛值是0.086 68;4°攻角、2.0Ma时力矩系数收敛值是0.200 3;4°攻角、2.5Ma时力矩系数收敛值是0.355 0;4°攻角、3.0Ma时力矩系数收敛值是0.646 0。通过3 000次迭代发现,力矩系数的变化曲线图在200次迭代,基本呈现直线状态,即此时在经过迭代力矩系数的值基本不发生变化,说明力矩系数的计算值达到收敛状态。将4°攻角、1.5Ma,4°攻角、2 Ma,4°攻角、2.5 Ma,4°攻角、3 Ma这4种不同马赫的力矩系数值在Excel表格中对比发现,在4°攻角一定的情况下,随着马赫数的升高,力矩系数值也在不断升高。力矩系数随迭代过程变化曲线如图6所示。
图6 力矩系数随迭代过程变化曲线
Fig.6 Curve of torque coefficient versus iterative process
3.2 空气绕炮弹流动情况窗口分析
图7(a)—图7(d)四幅图分别为4°攻角下1.5 Ma、2 Ma、2.5 Ma、3 Ma时X方向速度分布云图。
当马赫数Ma>1时,即速度V高于音速时,高动态下的弹体同时受到摩擦作用、涡阻作用以及激波作用。此时,随着周围环境压强Pa,密度ρ,温度T的升高,弹体表面的流速会瞬间下降。在Fluent求解运算中一般选择满足气体状态方程的理想气体(ideal-gas),在数学理论中满足气体状态方程的理想气体的激波作用没有厚度,是不连续的,但在高动态弹体实际飞行中,还应考虑密度和黏性的作用,实际的热传导速率也不容忽视。
虽然激波作用干扰高动态飞行的弹体,但实际的影响非常小,所以在数值计算中,得到的激波值较小,并且随着马赫数的升高,激波值逐渐减小。高马赫飞行的炮弹会受到气流压缩等影响,计算机计算的参数也会发生突变现象。仔细观察图像发现,炮弹的头部与尾部均可以看到头部激波和尾部激波,并且马赫数越高,产生的激波效果越明显。
图7 X方向速度分布云图
Fig.7 Velocity distribution nephogram in X direction
图8(a)—图8(d)四幅图分别为4°攻角下1.5 Ma、2 Ma、2.5 Ma、3.0 Ma时全弹表面压力分布云图。炮弹表面的压力分布是炮弹设计的重要参考指标。由图8可以看出,炮弹在飞行中弹头所受的压力最大,1.5 Ma时弹头所受压力值为2.27e+0.5 Pa,2 Ma时弹头所受压力值为5.55e+0.5 Pa,2.5Ma时弹头所受压力值为7.47e+0.5 Pa,3.0Ma时弹头所受压力值为9.82e+0.5 Pa。船尾所受压力较小,分别为1.03e+0.5、1.71e+0.5、6.02e+0.4、2.42e+0.3 Pa。
图8 全弹表面的压力分布云图
Fig.8 Cloud picture of pressure distribution on the projectile surface
由此可知,随着马赫数的增大,弹体表面压力值逐渐增大,其中弹头所承受的压力值最大,在材料处理方面,也要对弹头的材料进行特殊处理,保证弹体飞行过程中的稳定性。
3.3 炮弹模型受力分析
通过Fluent软件仿真,在cd-history、cl-history和cm-history 中读取可以得到最终收敛的阻力系数、升力系数和力矩系数值。为了节省计算机运行时间和降低网格划分难度,只对半个计算域进行计算,这样避免了流场中压力波等值边界的反射引起数值计算的不稳定。所以,真实的阻力系数(Cd)、升力系数(Cl)和力矩系数(Cm)应为计算机仿真得到计算值的2倍。
从图9中可知,攻角一定时,随着马赫数的增大,阻力值(Cd)也在不断上升,且上升程度由迅速逐渐变为缓慢,但最终呈上升趋势;马赫数一定时,随着攻角的增大,阻力值先下降后上升到峰值,然后下降,总体呈上升趋势。
图9 阻力随马赫数和攻角的变化曲线
Fig.9 The variation curve of resistance with Mach number and attack angle
从图10中可知,攻角一定时,随着马赫数的增大,升力值(Cl)也在不断上升,虽然在0°攻角时,上升率很小,但是升力系数在其余攻角时随着马赫数的增加而不断增大,且攻角越大,升力系数的增长趋势越快。马赫数一定时,随着攻角的增大,升力值也在不断上升,且马赫越大,升力系数的增长趋势越快。
图10 升力随马赫数和攻角的变化曲线
Fig.10 The variation curve of lift with Mach number and attack angle
从图11中可知,攻角一定时,随着马赫数的增大,力矩系数值(Cm)也在不断上升,且随着攻角的增大,力矩系数的增量差值也逐渐变大,攻角越大,升力系数的增长趋势越快。马赫数一定时,随着攻角的增大,力矩系数逐渐呈上升趋势,且上升程度逐渐缓慢。
图11 力矩随马赫数和攻角的变化曲线
Fig.11 The variation curve of torque with Mach number and attack angle
4 结论
为解决高转速导致远程制导炮弹控制难度大的问题,设计了一种利用船尾装置进行弹体减旋的气动外形结构,对气动参数进行了仿真和分析,根据仿真结果得出以下结论:
1) 具有减旋船尾结构的制导炮弹气动外形设计合理,气动参数满足稳定飞行条件,可为远程尾控制导炮弹设计提供参考依据。
2) 远程尾控制导炮弹弹头所承受的压力值较大,在材料处理方面,需要对弹头的材料进行特殊处理,以便于更好地适应飞行环境。
[1]袁荣亮.大国竞争下美国导弹防御系统发展动向及趋势分析[J].现代防御技术,2022,50(2):39-44.
YUAN Rongliang.Analysis on the development and trend of American missile defense system under great power competition[J].Modern Defense Techniques,2022,50(2):39-44.
[2]杨志豪,傅健,王良明.高超声速制导炮弹复合制导方法研究[J].兵器装备工程学报,2019,40(9):65-69.
YANG Zhihao,FU Jian,WANG Liangming.Research on compound guidance method of hypersonic guided projectile[J].Journal of Ordnance Equipment Engineering,2019,40(9):65-69.
[3]刘海涛,王旭刚,柳文安.超远程制导炮弹气动特性及仿真分析[J].兵器装备工程学报,2022,43(5):132-137.
LIU Haitao,WANG Xugang,LIU Wenan.Aerodynamic characteristics and simulation analysis of extended range guided munition[J].Journal of Ordnance and Equipment Engineering,2022,43(5):132.-137.
[4]李述清,张志利,魏玉淼,等.导弹精度评估方法仿真分析[J].系统仿真学报,2021,33(8):1825-1831.
LI Shuqing,ZHANG Zhili,WEI Yumiao,et al.Simulation analysis of assessment method for missile accuracy[J].Journal of System Simulation,2021,33(8):1825-1831.
[5]陈庆森,宋荣昌,关平,等.封控弹弹伞系统末弹道动态特性仿真[J].兵工学报,2015,36(S1):289-294.
CHEN Qingsen,SONG Rongchang,GUAN Ping,et al.Dynamic characteristics simulation for terminal trajectory of bomb parachute system in blockade and control ammunition[J].Journal of Ordnance Engineering,2015,36(S1):289-294.
[6]逯麒,陈坤,栗保明.电磁发射超高速制导弹药急速滚转稳定性仿真研究[J/OL].火炮发射与控制学报:1-7[2022-07-04].
LU Qi,CHEN Kun,LI Baoming.Simulation research on rapid rolling stability of electromagnetically launched hyper velocity projectile[J/OL].Journal of Artillery Launch and Control:1-7[2022-07-04].
[7]王中原,常思江.横向弹道修正的一种快速计算方法[J].兵工学报,2014,35(6):940-944.
WANG Zhongyuan,CHANG Sijiang.An approximated calculation method for lateral trajectory correction[J].Journal of Ordnance Engineering,2014,35(6):940-944.
[8]邱荷,刘荣忠,郭锐,等.末敏弹的可视化仿真研究[J].计算机仿真,2013,30(1):137-140.
QIU He,LIU Rongzhong,GUO Rui,et al.Study on visual simulation of terminal-sensitive projectile[J].Computer Simulation,2013,30(1):137-140.
[9]刘一鸣,李杰,王一焕,等.旋转弹用轴向减旋平台控制系统设计[J].电子器件,2017,40(4):968-972.
LIU Yiming,LI Jie,WANG Yihuan,et al.Design of the control system for axial anti-rotation platform used in rotating ammunition[J].Electronic Devices,2017,40(4):968-972.
[10]祝敬德,李杰,张松,等.半捷联模拟控制减旋电路设计[J].实验室研究与探索,2015,34(6):46-49,62.
ZHU Jingde,LI Jie,ZHANG Song,et al.Design of rotaryreduction analog control circuit for semi-strapdown system[J].Laboratory Research and Exploration,2015,34(6):46-49,62.
[11]沈仲书,刘亚飞.弹丸空气动力学[M].北京:国防工业出版社,1984.
SHEN Zhongshu,LIU Yafei.Projectiles aerodynamics[M].Beijing:National Defense Industry Press,1984.
[12]唐玉发,张合,徐国泰.基于MATLAB与ANSYS的火箭弹气动特性工程计算与数值仿真[J].南京理工大学学报,2013,37(5):682-686.
TANG Yufa,ZHANG He,XU Guotai.Engineering calculation and numerical simulation of rocket aerodynamic characteristics based on MATLAB and ANSYS[J].Journal of Nanjing University of Science and Technology,2013,37(5):682-686.