基于激波场张量捕捉的超声速弹丸跟踪*

顾国华，杨 溢，王积忠，陈 浩

（陆军军官学院，合肥 230031）

0 引言

空气动力学研究表明，当旋成体（如弹丸）以超声速在空中飞行时将产生激波，不同物体、不同飞行条件产生的激波不同，激波状况体现了飞行物的气动特征。为研究旋成体气动特性，国内外学者和团队从不同角度和方向在理论和实验上进行了较为深入的研究［1-4］。

从公开报道看，Munk［5］较早提出了旋成体气动特性，提出基于细长体位势流理论并进行数值模拟解算，可求解 0～5°迎角的气动问题。Jorgense［6］将横流理论气动分解为势流项和横流粘性项，并对相关阻力进行修正和分解。后来逐渐发展了激波捕捉方法，如Harten1983年提出的TVD差分格式，该方法是一种具有二阶精度的总变差不增算法；之后发展为对称型TVD。近代高分辨率的激波捕捉算法得以长足发展，大体可归纳为以下4种：

①TVD差分格式的升级，ENO格式和WENO格式。②Godunov间断分解方法。③FCT差分格式。④Roe差分格式。

为探索基于激波场张量的弹丸跟踪，本文提出建立三维Euler方程求解来确立激波捕捉张量，通过激波捕捉阵列多声信号融合处理，实现对高精度反演激波声场位置解算和跟踪。

1 弹丸激波场

由空气动力学原理可知，当弹丸以超声速在空中飞行时，势必造成对弹丸头部周边空气的挤压，而使得弹尾部空气稀疏。因而在头部形成压缩波（高压波），在尾部形成稀疏波。通常将其压缩波称为弹丸激波。是一种非线性大气声学现象，如图1所示。由激波形成的某一区域称为激波场。

在二维空间分析，激波的高压区传播速度大于声速，而低压区则低于声速。因而，通常形成起伏状波形，形似英文字母“N”，所以俗称为“N”波。

弹丸“N”波的形成过程大体上可分为4个阶段［7］如图 2所示。

2 弹丸激波捕捉

2.1 激波张量分析

二维Euler方程是分析空气动力学的基本模型。为了抑制振荡和过冲现象，将通量项分解为对流通量项和压力项，在文献［8-9］的基础上拓展，采用理想气体模型建立三维Euler方程（1）。采用张量捕捉方法的优势在于信号更全面、更精确，利于分析处理。

其中，

其中，ρ为压力，u，v和 w 分别为 x，y，z方向上的速度分量，e，E为单位体积内的内能和总能。

方程（1）可转化为守恒型形式：

由式（1）、式（2）可得

解得非线性雅可比系数特征值矩阵ΛA如下，同理可解得 ΛB、ΛC。

解得A特征矢量矩阵RA如下，同理可解得B、C 的特征矢量矩阵 RB、RC。

2.2 激波捕捉累积多声信号融合算法

根据文献［10］建立的二维守恒型Euler方程，在任一点（x，y）处，将体积采用有限体积法可得：

其中，Rx，y表示流出控制体的净通量，在某一定点（i，j）处，净通量以通量形式可表示为

在α角度上，令信号空间U1与V1张成，噪声子空间由U2与V2张成，将其进行进行奇异值分解［11］

以i方向的通量为例，将通量计算公式改为：

对于由M*N矩阵组成的声矢量传感器阵列

在i-1/2界面处的压强p可表示为

进行奇异值分解可得

式中：∑'1为非零奇异值矩阵。U'1和V'1张成信号子空间，U'2和V'2张成与信号子空间正交的噪声子空间。

考虑到实际接收的数据矩阵长度有限，因此，对循环相关矩阵取最大似然估计

其中，L为快拍数。

3 弹丸运动跟踪

超声速弹丸空中飞行时，运动轨迹影响因素多、不易捕捉、具有一定随机性。试验与反演以某弹丸对空射击背景，航路参数为：等速直线水平飞行1 000 m，高度500 m，捷径400 m，目标航速200 m/s，航向由北向南。通过测试与反演弹目偏差解算跟踪弹丸，检验弹丸运动跟踪滤波效果。

经激波声场3 500个采样点后提取弹丸运动轨迹如图3所示。采用ENO格式处理后的的轨迹如图4所示，采用Godunov法处理后的轨迹如图5所示，采用本文提出的激波张量捕捉法处理后的波形如下页图6所示。

仿真实验中，通过比较分析，在3 500点采样数据中，相对于初始数据，采用3种方法处理后的残差曲线统计，采用ENO格式法的跟踪误差为18.5%，采用Godunov法的跟踪误差为15.6%，采用本文提出的激波张量捕捉法的跟踪误差减小到8.9%。

4 结论

当前，采用弹丸激波场对弹丸进行跟踪的研究还不够深入。本文探讨了基于激波场张量的弹丸跟踪仿真，实验和反演表明，提出的激波张量捕捉法收敛效果好，精度较高，有效解决了复杂激波声场中弹丸跟踪精度和收敛速度的矛盾，为超声速运动目标定位和跟踪提供了一种新方法。

［1］URSZULA L，KRYSTIAN S.Wavelet based shock wave and muzzle blast classification for different supersonic projectiles［J］.Expert Systems with Applications，2014，41（11）：5097-5104.

［2］KNIGHT D，LONGO J，DRIKAKIS D，et al.Assessment of CFD capability for prediction of hypersonic shock interactions［J］.Progress in Aerospace sciences，2012（48）：8-26.

［3］GNANI F，LO K H，ZARE B H，et al.Shock wave diffraction in the presence of a supersonic co-flow jet［J］.Shock Waves，2016（3）：1-10.

［4］赵占龙.弹箭线性与非线性气动特性计算［D］.南京:南京理工大学，2012.

［5］MUNK M M.The aerodynamic force on airship hulls［R］.NACA-TR-184，1923.

［6］JORGENSE L H.Prediction of static aerodynamic characteristics for slender bodies alone and with lifting surfaces to very high angles of attack［R］.NASA TR-R474，1977.

［7］张飞猛.对空射击声学靶脱靶量测试系统的精度分析［J］.兵工学报，2000，21（1）:23-26.

［8］张德良.计算流体力学教程［M］.北京:高等教育出版社，2010：259-260，289.

［9］WAGNER B，SCHMIDT W.Theoretical investigations on shock wave-boundary layer interaction in cryogenic nitrogen［J］.Springer International Publishing，2015，365（3）:341-350.

［10］贾勇，孙刚，刘苏.一种改进混合迎风格式在翼型流场激波捕捉中的应用［J］.力学季刊，2007，28（2）:223-227.

［11］郑彤.复杂噪声背景下基于声矢量传感器阵列的汽车噪声源方位估计方法研究［D］.长春：吉林大学，2012.

［12］王勇，王学德，许国龙.超声速弹丸高空流场数值模拟与稀薄气体效应研究［J］.兵器装备工程学报，2016，37（10）:13-18.

［13］李树军，李开端.基于自适应α-β滤波的目标跟踪算法研究［J］.科学技术与工程，2012（12）：9367-9369.

［14］顾国华，陈栋，王彩，等.基于声学靶传感器的弹着点测试研究与实现［J］.电子测量技术，2007，30（2）:154-156.