旋转尾翼弹马格努斯效应数值模拟
2018-03-31王学德常思江刘士杰
吴 放,王学德,常思江,刘士杰
(1.南京理工大学 能源与动力工程学院,江苏 南京 210094;2.西北工业集团有限公司 科研二所,陕西 西安 710043)
尾翼稳定是弹箭飞行稳定方式的一种,主要是依靠尾翼在飞行过程中的受力作用产生可抑制攻角发散的俯仰力矩。在实际中,为了减小弹体气动外形不对称或质量偏心等对其弹道特性的影响,往往使尾翼稳定弹在飞行过程中具有一定的转速。但是当攻角和旋转同时存在时,会产生一个额外的侧向力,即马格努斯力。这是由流场附面层的畸变等因素造成弹体两侧压力分布不对称所导致的。在通常情况下,马格努斯力约为弹体所受法向力的1%~10%[1],马格努斯力矩矢量在攻角平面内垂直于弹轴。
尾翼式旋转弹箭的马格努斯力矩往往较大。对于弹体的马格努斯效应来说,其影响因素包括附面层位移厚度、气流切应力、压力梯度变化的不对称,附面层与涡的相互影响等。对于弹翼的马格努斯效应来说,其主要影响因素包括弹翼的附加攻角、附加速度、弹翼后缘压力的差异等。因此,旋转弹箭的马格努斯效应具有极大的复杂性和不确定性。研究旋转弹箭的马格努斯效应,减小其对弹箭飞行的影响,对弹箭设计有重大意义。
国内方面,我国对马格努斯效应的研究开始时间较晚,在20世纪80年代才开始有学者进行研究。苗瑞生等[2]在1987年总结了国内外对于旋转弹箭马格努斯效应研究的发展进程和研究内容,讨论了旋转弹箭空气动力设计的相关问题,为之后的研究指出了数个方向。王树军等[3]在2009年研究了横向喷流对旋转弹箭的影响,发现正攻角时力放大因子和旋转速度正相关,负攻角时负相关。李衡等[4]在2011年研究了自动驾驶仪在马格努斯效应下的工作状态,发现加速度反馈可抑制马格努斯效应对弹箭的影响。赵博博等[5]在2015年对扭曲尾翼弹箭进行了气动机理研究,结果表明扭曲尾翼可改善马格努斯效应的影响。刘周等[6]在2016年采用DDES方法对高速旋转弹箭进行了数值模拟,发现分离点会对马格努斯效应产生较大影响。王刚等[7]在2017年对旋转弹箭进行非定常模拟,证实了ARMA建模法对飞行轨迹预测的正确性和快速性。
国外方面,美国陆军弹道研究所(BRL)的Martin最先进行马格努斯效应的模拟计算[2],他在研究小攻角下旋转的无限长圆柱时,采用小扰动法和细长体理论求出马格努斯力和力矩。20世纪80年后对马格努斯效应的数值模拟进行了多种研究。Despirito[8]在2008年对比了RANS和RANS/LES模型亚声速时对ANF射弹的模拟结果,在低马赫数时RANS/LES偏离实验值较多。Klatt等[9]在2012年采用RANS模型对旋转弹箭进行流场模拟,并与风洞实验进行对比,研究了大攻角下马格努斯效应的变化规律。Benzi等[10]在2016年数值模拟了超声速下弹箭表面摩擦力和热传递过程,讨论了Knudsen数和马赫数的相互关系。
考虑到国内对于尾翼式旋转弹箭马格努斯效应的研究较少,对于马格努斯效应的实验测试也相对粗略,有必要对此开展深入的理论研究。本文针对某低速旋转尾翼稳定弹,利用计算流体力学方法,在不同工况条件下进行了数值模拟,由此探讨了该弹箭的马格努斯效应及其影响因素。研究结果对深入理解和应用旋转弹箭的马格努斯效应具有一定的指导作用。
1 数值计算方法
本文采用Navier-Stokes方程、Spalart-Allmaras湍流模型、滑移网格技术对尾翼旋转弹的流场进行数值模拟。通过定常计算,模拟绕流旋转弹箭的流场。计算结果收敛后将其作为初始条件,通过滑移网格技术进行非定常计算,模拟弹箭绕弹轴自旋的流场。滑移网格要求生成包围弹箭的运动区和外部静止区。数值计算时,弹箭和内部运动区域同时以给定角速度旋转以模拟弹箭旋转,外部区域静止。相较于旋转坐标系方程将非定常问题转化为定常问题,滑移网格方法的非定常计算结果更加精确可靠。
所有力和力矩相对于弹体坐标系进行测量,攻角平面为xOy平面,弹箭关于xOy平面和xOz平面对称。
本文中马格努斯力动态导数Cz,p和马格努斯力矩动态导数Cm,y,p的公式[11]为
(1)
(2)
式中:Cz为马格努斯力系数,Cm,y为马格努斯力矩系数。
2 几何模型
本文使用ANF弹箭模型。图1为弹箭几何模型细节。弹箭弹体直径d=31.75 mm。以弹径d为参考长度,弹箭由2.836d的头部、7.164d的圆柱体和4个对称的弹翼组成。弹头部为10°圆锥,弹箭顶端曲率半径r=0.002d。弹翼尺寸为0.94d×1d,具有尖锐前缘,后缘厚度为0.06d。ANF的重心位于从弹头起6.1d处,它的轴向惯性矩为0.001 329 kg/m2。
3 计算网格和边界条件
计算网格如图2所示。内部运动区域网格数为240万,周向节点100个,径向节点70个,轴向节点150个。外部静止区域网格数为160万,总网格数为300万。
内部运动区和外部静止区的交界面均采用滑移边界条件。外部静止区的外边界采用压力远场边界条件。将来自远场边界的第1层网格设置为吸收层,并且添加阻尼项以防止可能污染流动的数值波反射。弹箭壁面为黏性绝热壁面,其中垂直于壁面的初始网格节点间距为0.01 mm,满足边界层分辨率y+≤1的标准。
美国的阿诺德工程和开发中心采用阿诺德空军基地的冯卡门超声速风洞对ANF型弹箭进行过大量实验,具有较为完整的实验数据。本文所用实验数据均来自参考文献[11]。风洞中马赫数Ma、雷诺数Re、总压p0、总温T0如表1所示。
表1 风洞条件
计算工况:Ma=2.49,Ω=0.015,攻角α=20°,40°,60°,80°;Ma=2.49,Ω=0.025,攻角α=20°,40°,60°,80°。
4 计算结果及分析
4.1 流场分析
图3为弹箭在攻角α=40°,Ω=0及Ω=0,0.025时,x=9.2d处的马赫数云图。
从图3中可知,当弹箭无旋时,附面层不发生畸变,且关于攻角平面xOy对称分布。弹翼下方的空气在翼尖附近会翻向上方,使得翼尖附近上表面处的压力逐渐与下表面压力相等。当弹箭以Ω=0.025旋转时,附面层发生畸变,左侧变厚,低压区域减小;右侧变薄,低压区域增大。产生此种现象的原因:弹箭在正攻角条件下飞行时,由于空气的黏性作用,弹箭表面附近的空气也跟随弹箭一起旋转,产生了环流;弹箭左侧空气流动方向与环流方向相反,使气流速度和马赫数减小,附面层变厚;弹箭右侧空气流动方向与环流方向相同,使气流速度和马赫数增大,附面层变薄;气流压力与速度成反比,弹箭左右两侧产生压力差,从而形成了侧向力。
图4、图5是弹箭在α=20°,60°,Ω=0.015,0.025时的马赫数云图。
从图4和图5中可知,Ω=0.015时,弹翼附近流场畸变程度小于Ω=0.025时的畸变程度。Ω=0.015时,弹箭左右两侧的压力差小于Ω=0.025的压力差。因此,随着转速增大,弹箭受到的侧向力也会增大。
将图4(a)和图5(a)对比,图4(b)和图5(b)对比可知,当α=20°时,弹箭下表面和上表面压力差较小。然而当α=60°时,弹箭下表面压力远大于上表面压力。这是由于随着攻角的增大,弹箭与自由来流夹角增大。α=60°时弹箭左右两侧的低压区域大于α=20°时。α=20°时,在弹翼根部会产生涡,而α=60°时则不会。这是由于α=60°时,弹翼上表面气流速度较小,弹翼三维效应不明显,速度环量小。
图6是弹箭在α=40°,Ω=0.015时的马赫数云图。从图6可知,弹箭超声速飞行时,气流流经弹箭下表面时形成斜激波。随后在气流通过弹头与弹体结合处之后,气流经过膨胀后形成膨胀波。不同的是,在较大的正攻角条件下,气流流经弹箭头部上表面时,气流外折,形成膨胀波。在弹尾处,由于气流黏性的作用,弹翼上、下表面流速不同,气流分离,在弹底形成低压区,生成涡流,从而产生压差阻力。
4.2 气动特性分析
4.2.1 马格努斯力和力矩系数
图7~图9是在α=20°,40°,60°,Ω=0.015时,本文CFD模拟的马格努斯力系数和马格努斯力矩系数随滚转角Φ变化的曲线,并与文献[11]中数据进行比较。
从图中可知,当α=20°时,在旋转角度Φ=0°~40°区间内,模拟值与参考值符合情况良好;在Φ=40°时出现一定误差;Φ=40°之后模拟值稍小于参考值。α=40°时,在旋转角度Φ=20°~40°区间内,模拟值小于参考值,且存在可以接受的误差,其他角度区间符合情况良好。α=60°时,模拟值和参考值符合得十分良好。综上可认为本文CFD模拟值和文献中参考值符合较好,验证了本文算例的可靠性。
图10、图11是Ω=0.015,0.025时不同攻角下的弹箭马格努力系数和马格努斯力矩系数随滚转角Φ变化的曲线。
从图10、图11中可以看出,马格努斯力系数和力矩系数在第1个四分之一周期迅速减小,之后以稳定的周期性规律发生变化。由于弹箭上存在4个尾翼,系数变化通常在半个旋转周期内呈现2个峰和2个谷。攻角增大时,马格努斯力系数和力矩系数也会随之增大,在α=40°~60°之间取最大值。由于系数变化具有周期性,所以在定常计算收敛后,只需要计算旋转180°的数据,即可得出弹箭的飞行状态。
由图12可以看出,随着攻角增大,平均马格努斯力系数和平均马格努斯力矩系数的绝对值也逐渐增大,在α=40°~60°达到峰值,随后随攻角增大而减小。
图13是不同转速时平均马格努斯力系数和平均马格努斯力矩系数图。图中线段是对数据的最小二乘法线性拟合,因为马格努斯力系数和力矩系数在旋转速率Ω=0时为0,所以线性拟合在原点截断。
4.2.2 马格努斯动态导数
马格努斯动态导数由式(1)、式(2)确定。图14、图15为以攻角为自变量的动态导数,即马格努斯力动态导数Cz,p、马格努斯力矩动态导数Cm,y,p。
对于马格努斯力动态导数,CFD模拟结果在α=0°~40°范围内与实验数据符合良好,之后CFD模拟值高于实验值,在α=60°时误差约为18.9%。在α=0°~40°范围内,马格努斯力动态导数随攻角增大而增大,增大幅度逐渐减小,在α=40°~60°达到峰值,之后随着攻角增大而减小。对于马格努斯力矩动态导数,CFD模拟结果在α=0°~30°范围内与实验数据符合良好,在α=40°时误差较大,约为10.6%,之后CFD模拟结果与实验数据依然符合良好,但在α=80°时误差较大,约为12.9%。在α=0°~40°范围内,马格努斯力矩动态导数的绝对值随攻角增大而增大,增大幅度逐渐减小,在α=40°~60°达到峰值,之后随着攻角增大而减小。
雷诺数对马格努斯效应也有显著影响[12],有:
(3)
(4)
式中:l为全弹长。
高原环境下弹箭所处大气环境在铅垂方向变化剧烈,空气密度降低,雷诺数变化十分显著,对马格努斯动态导数产生巨大影响,使其产生非线性的变化。同时,弹箭所处大气的密度较低,使弹箭稳定性降低,在干扰力矩作用下,攻角的扰动逐渐增大,弹箭附面层分离加强,进一步影响马格努斯效应。
4.2.3 其他气动力系数
由于旋转并不会改变弹箭上、下表面的压力差,故升力不随旋转速率变化。同样,旋转速率改变也不会使阻力系数和俯仰力矩系数变化。
图16~图18是弹箭气动力系数数据图。图中Cl为升力系数,Cx为阻力系数,Cm,z为俯仰力矩系数。
由图可知,CFD模拟数据与实验数据、参考数据符合良好。攻角增大,弹箭迎风面和背风面压力差增加,使升力增大,达到峰值后减小。阻力也随着攻角增大而增大,达到峰值后减小。且压差阻力在总阻力中所占的比例增大。由图18可以看出,在α=0°~40°,俯仰力矩系数的绝对值随攻角增大而增大,之后随攻角增大而减小,α=70°后再次增大。
5 结论
本文通过计算流体力学方法数值模拟得到ANF射弹在Ma=2.49时,旋转速率Ω=0.015,0.025,攻角α=20°,40°,60°,80°时的气动力系数,研究了马格努斯效应及其影响因素,得到以下结论:
①在同一旋转速率Ω时,马格努斯力和力矩系数的绝对值随攻角的增大而增大。在α=40°~60°范围内达到峰值。在同一攻角时,马格努斯力和力矩系数具有自旋依赖性,会随着旋转速率Ω的改变而改变,其绝对值在较高旋转速率Ω时通常以一定幅度增加,呈线性关系。
②马格努斯力、力矩的动态导数的绝对值在α=0°~40°范围内随攻角增大而增大,且增大幅度逐渐减小,在α=40°~60°达到峰值,之后随着攻角增大而减小。
③在高原环境下,马格努斯效应对雷诺数十分敏感,研究其产生的影响,对改良现有弹种在高原环境下作战能力具有指导意义,在后续工作中会对此进行深入研究。
[1] GRAFF G Y,MOORE F G. Empirical method for predicting the Magnus characteristics of spinning shells[J]. AIAA Journal,2012,15(10):1 379-1 380.
[2] 苗瑞生,吴甲生. 旋转导弹的气动布局[J]. 力学进展,1987,17(4):479-488.
MIAO Ruisheng,WU Jiasheng. Aerodynamics of spinning projectiles[J]. Advances In Mechanics,1987,17(4):479-488. (in Chinese)
[3] 王树军,胡俊,吴甲生,等. 旋转导弹横向喷流/外流干扰的数值模拟研究[J]. 系统仿真学报,2009,21(11):3 472-3 475,3 478.
WANG Shujun,HU Jun,WU Jiasheng,et al. Numerical simulation of interaction between lateral jet and external flow over spinning missile[J]. Journal of System Simulation,2009,21(11):3 472-3 475,3 478. (in Chinese)
[4] 李衡,薛林,高庆丰. 马格努斯效应对旋转弹自动驾驶仪影响的研究[J]. 现代防御技术,2011,39(6):100-105.
LI Heng,XUE Lin,GAO Qingfeng. Research on the effect of Magnus effects on rolling missile autopilot[J]. Modern Defence Technology,2011,39(6):100-105. (in Chinese)
[5] 赵博博,刘荣忠,郭锐,等. 扭曲尾翼飞行器滚转特性[J]. 航空动力学报,2015,30(1):142-148.
ZHAO Bobo,LIU Rongzhong,GUO Rui,et al. Numerical prediction of the Magnus effect for twist fin swept flight projectile[J]. Journal of Solid Rocket Technology,2015,30(1):142-148. (in Chinese)
[6] 刘周,谢立军,杨云军,等. 弹丸旋转空气动力效应非定常数值模拟[J]. 航空学报,2016,52(5):1 401-1 410.
LIU Zhou,XIE Lijun,YANG Yunjun,et al. Unsteady numerical simulation of aerodynamics effect of a spinning projectile[J]. Acta Aeronautica et Astronautica Sinica,2016,52(5):1 401-1 410. (in Chinese)
[7] 王刚,邢宇,李亚楠. 旋转弹气动力建模与飞行轨迹仿真[J]. 航空学报,2017,53(1):108-117.
WANG Gang,XING Yu,LI Yanan. Aerodynamic modeling and flight trajectory simulation of spinning projectile[J]. Acta Aeronautica et Astronautica Sinica,2017,53(1):108-117. (in Chinese)
[8] DESPIRITO J. CFD prediction of Magnus effect in subsonic to supersonic flight[C]//The 46th AIAA Aerospace Sciences Meeting and Exhibit,Aerospace Sciences Meetings. Nevada:AIAA,2008:768-776.
[9] KLATT D,HRUSCHKA R,LEOPOLD F,et al. Numerical and experimental investigation of the Magnus effect in supersonic flows[C]//The 30th AIAA Applied Aerodynamics Conference,Fluid Dynamics and Co-located Conferences. Louisiana:AIAA,2012:326-341.
[10] BENZI J,GU Xiaojun,BARBER R W,et al. High-speed rarefied flow past a rotating cylinder:the inverse Magnus effect[J]. AIAA Journal,2016,54(5):1 670-1 680.
[11] VISHAL A B. Numerical prediction of roll damping and Magnus dynamic derivatives for finned projectiles at angle of attack[C]//The 30th AIAA Applied Aerodynamics Conference. Louisiana:AIAA,2012:683-695.
[12] 郭锡福. 火炮武器系统外弹道试验数据处理与分析[M]. 北京:国防工业出版社,2013:281.
GUO Xifu. Exterior ballistic test data processing and analysis for gun weapon systems[M]. Beijing:National Defense Industry Press,2013:281. (in Chinese)