多发弹气动参数联合辨识方法研究

刘洋， 常思江， 魏伟

(1.南京理工大学 能源与动力工程学院, 江苏 南京 210094; 2.瞬态冲击技术重点实验室, 北京 102202;3.北京理工大学 机电学院, 北京 100081)

0 引言

气动辨识是飞行器系统辨识中的关键问题，是炮弹、火箭弹、导弹、飞机等获取自身气动力参数的重要手段，在各类弹箭的气动外形设计、系统参数确定、设计定型、射表编制等过程中具有不可替代的作用。气动辨识是以一定的弹箭运动模型和辨识算法为基础，从野外靶场或靶道实测弹箭运动数据中提取出弹箭气动力参数的过程。理论上讲，气动辨识本质上就是利用输入、输出信息确定系统结构参数。从工程应用角度看，与导弹、飞机等有控飞行器相比，常规炮弹、火箭弹等无控飞行器的气动辨识难度更大，其原因主要在于：1)描述无控飞行器的动力学模型往往具有较强的非线性结构[1]，数学处理难度较大；2)常规炮弹、火箭弹发射时的瞬时状态往往具有随机性且难以准确测量，并且由于稳定飞行过程中不受控，无法为系统提供有效的输入；3)某些弹丸(如枪弹)的气动力参数非线性特性明显，对弹道影响较大，而非线性气动力研究本身就具有很大难度；4)常规炮弹、火箭弹等开展飞行试验时，由于体积空间、试验成本、发射条件、飞行环境等各方面的限制，往往造成测量数据的种类、精度及采样频率等，难以同导弹、飞机等相提并论，从而给气动辨识带来困难。

由于常规炮弹、火箭弹等无控飞行器的射程、密集度、稳定性等与气动力参数密切相关，尽可能准确地获取气动力参数，对性能评估、方案改进优化等具有重要作用，故有必要深入开展气动辨识方法与技术研究。从现有文献看，相关气动辨识研究的重点主要是辨识算法，例如Chapman-Kirk方法(C-K法)[2]、Levenberg-Marquardt方法(L-M法)[3-5]、极大似然法[6-9]以及某些全局优化算法[3,10-11]等。实际工程中，为了确定弹箭气动力参数，一般都是在不同初速、射角等条件下，成组发射弹箭并进行测量。对于1组射弹，由于气象变化、起始扰动、弹箭结构参数等具有不同程度的随机性，即便测量设备具有相同误差，同一组试验的每发弹也不可能测得完全相同的弹道数据。已有实践表明，同一组试验中每发弹的气动力参数辨识结果往往有一定差异，对于某些动态气动力参数，如俯仰阻尼力矩系数、马格努斯力矩系数等，发与发之间的差异甚至很大。目前，绝大多数相关研究关注的都是单发弹丸的气动辨识过程，主要涉及模型和辨识算法，并未考虑如何优化处理多发弹丸辨识结果的差异。例如，在对气动力参数精度要求甚高的射表试验数据处理中[12-13]，就是对多发弹测量数据辨识出的气动力参数取平均值，作为该弹的气动力参数以供后续使用。然而实际科研中发现，采用平均值气动力参数去重构弹道，经常出现计算值与单发测量值差异较大的现象。这表明，尽管基于单发测量数据的独立辨识结果具有最优性，但平均值气动参数对于多发测量数据来说并不是全局最优的。目前常用的C-K法、L-M法、极大似然法等都属于局部优化方法，在用于单发测量数据独立辨识时具有精度较高、收敛性好等优点，但用于多发弹气动参数联合辨识却具有较大局限性，辨识过程极不稳定。近年来全局优化算法在气动辨识中的应用研究逐渐增多，但大多还是解决单发数据气动辨识的状态初值问题、收敛性问题等[3,10]。在考虑多发弹气动参数联合辨识方面，国外有一些初步研究。如文献[3]针对仿真数据，利用L-M法和遗传算法分别对多发弹的气动参数进行了联合辨识，然而，L-M法收敛于局部极值，而遗传算法的结果残差比L-M法更大，这说明两种算法都未能在多发弹气动参数联合辨识中取得较好的全局最优解。文献[4]中给出了L-M法和有限差分法的气动辨识结果，比文献[3]有较大改进，但也未能解决局部极值问题。文献[14]提出采用多组试验数据联合辨识某无控空间探测器的气动参数，但其重点在于将气动参数表示成攻角和马赫数的非线性函数，未考虑可能存在的局部解问题，并未给出具体的全局策略和算法。

本文针对枪弹、常规炮弹、火箭弹等无控弹箭的气动参数联合辨识问题开展研究，提出一个可同时处理多发弹测量数据并给出唯一具体气动力参数值的全局优化策略。以单发弹独立辨识为基础、以目前国际优化领域新兴的差分进化算法为工具，建立了利用多发弹测量数据联合辨识气动力参数的全局优化流程，并通过对某大口径炮弹攻角纸靶试验数据的处理，验证了本文所提策略的正确性和有效性。

1 多发弹测量数据联合辨识策略

1.1 辨识流程设计

弹丸气动辨识流程主要由以下5个部分组成[5]：测量数据预处理、数学模型、辨识算法、准则函数和辨识结果后处理。测量数据预处理和辨识结果后处理主要为剔除野值与曲线的平滑滤波，数学模型则根据辨识问题的不同可采用4自由度、5自由度、6自由度外弹道方程或攻角运动方程等，准则函数的选取与辨识算法有关。因此，多发弹测量数据联合辨识的流程设计重点是辨识算法。

辨识算法按寻优点个数可分为局部优化算法和全局优化算法，按搜索方向可分为梯度法和非梯度法，在气动辨识中，一般选用基于梯度的局部优化算法或非梯度的全局优化算法。

文献[4]中选取L-M法作为多发弹测量数据联合辨识的算法，该算法是局部优化算法中的一种，具有计算时间较短、计算效率较高的优点，在对单发弹的辨识中得到了广泛应用。但是，利用局部优化算法对多发弹测量数据进行辨识，存在以下2个缺点：1)多发弹联合辨识中求得的矩阵可能达到几十甚至上百阶[3]，若矩阵接近奇异矩阵，计算可能发散；2)局部优化算法寻得全局最优的概率较低[3-4,15]，当预先给出的待辨识参数与实际值相差较大时，容易陷入局部最优解。此外，由于实测数据的误差种类较多，且有些误差难以定量甚至难以定性(如纸靶试验中弹丸穿靶瞬间对其自身的干扰)，可能使局部最优解个数增多，从而导致局部优化算法搜索到全局最优解的概率更低。

因此，局部优化算法不适用于多发弹测量数据联合辨识。故考虑使用全局优化算法作为多发弹测量数据联合辨识的算法，该类算法计算稳定性更好，搜索到全局最优解的概率更高，但搜索维度较多，计算效率较低。

为改善全局优化算法的性能，本文提出一种新的辨识流程：先利用局部优化算法分别对各发弹进行快速辨识(以下简称独立辨识)，利用辨识结果生成全局优化算法中的搜索空间；再使用全局优化算法对多发弹进行联合辨识(以下简称联合辨识)。具体辨识流程如图1所示。图1中n表示试验弹丸的发数，1发弹对应1套试验数据、气象数据和结构参数，n发弹共有n套数据。待辨识参数为状态初值和气动参数。状态初值为弹道方程中的积分初值或攻角方程的初始幅值、初始相位等，由于每发弹的发射条件、起始扰动等都有可能不同，导致每发弹的状态初值差异较大，而某些状态初值，如弹丸角运动的初始幅值与初始相位等难以直接测量，故状态初值应作为待辨识参数与气动参数同时进行辨识。每组状态初值中包含N1个待辨识参数，每组气动参数中包含N2种待辨识气动参数，如阻力系数、升力系数导数等。

图1 辨识流程Fig.1 Identification procedure

图1中的独立辨识中仅对单发弹测量数据辨识，可根据实际情况选用C-K法、极大似然法等辨识算法，辨识结果为联合辨识提供搜索空间，该搜索空间用于联合辨识中待辨识参数的初始化，使其位于较小的约束范围内，进而减少迭代次数，提升计算效率与稳定性；联合辨识中对多发弹测量数据同时辨识，辨识算法选用全局优化算法。

因此，多发弹测量数据联合辨识的具体步骤为：

1)对多发弹测量数据进行预处理，剔除野值、平滑测量曲线等，并根据测量数据种类及待辨识气动参数构建数学模型；

2)使用局部优化算法进行独立辨识，得到n组状态初值与n组气动参数；

3)利用独立辨识结果生成联合辨识的搜索空间，使用全局优化算法对多发弹相关气动参数进行联合辨识，得到n组状态初值与1组气动参数；

4)对辨识结果进行平滑、插值等后处理。

1.2 准则函数与最优解定义

在整个辨识过程中，独立辨识与联合辨识的准则函数应在形式上保持一致，如：不能同时使用最小二乘准则和极大似然准则。本文使用最小二乘准则作为辨识准则。

最小二乘准则可用(1)式表示：

(1)

在辨识、优化等问题中，J值越小，表明全局优化性越好[16]。在约束域内：若J是最小值，则对应的解是全局最优解；若J是极小值且非最小值，则对应的解是局部最优解。将此概念扩展到多发弹测量数据联合辨识的问题中，利用辨识出的1套气动参数和n套状态初值重构出n条弹道，若对应的J值是约束域内的最小值，则该气动参数与状态初值共同构成全局最优解，反之则为局部最优解。以上概念也作为后续评价全局性与最优解的依据。

2 差分进化算法的应用

2.1 基本算法

差分进化算法是近年来国际优化领域新兴的一种全局优化算法，相比于遗传算法、粒子群算法等，差分进化算法具有更好的收敛速度、鲁棒性和全局寻优能力[16]，适合计算量较大的联合辨识。

差分进化算法主要是由“变异”、“交叉”和“选择”3个部分构成。但是在联合辨识中，种群的搜索空间与初始化也较为重要。

差分进化算法的“变异”遵循(2)式：

(2)

式中：k为迭代次数；θ是待变异的个体，θr1和θr2是随机的2个个体；a是[0,1]之间的1个随机数。与粒子群算法等具有一定方向性算法不同，差分进化算法的变异方向是完全随机的，这意味着算法具有更好的全局搜索能力。

“交叉”是指2个个体间随机交换若干参数生成新的个体，也可以进行算数重组或者连续重组。“选择”可以选用模拟退火的选择方式，防止算法因为早熟而产生局部最优解[17]。

2.2 算法应用

在全局寻优过程中，较大的搜索空间不仅增加迭代次数、降低辨识效率，对于某些弹道方程，如果气动参数与状态初值的迭代初值不在合理区间内，可能导致计算发散。对于带约束的优化问题，若能将搜索空间预先确定在约束范围内，可大幅减少无效迭代。因此，在联合辨识中，获得合理搜索空间是非常重要的环节。

差分进化算法应用于图1中的联合辨识部分，具体流程如图2所示。

图2 差分进化算法的应用流程Fig.2 Flowchart of differential evolution algorithm

如图2所示，联合辨识的搜索空间由独立辨识结果确定。由于每发弹状态初值不同，故以状态初值的辨识结果为中心按一定比例扩展生成新的搜索空间。经反复调试发现扩展比例与状态初值类型有关，当状态初值为可测量(如速度、位置等)，扩展比例可取±2%～±5%；当状态初值为不可测量(如角运动初始幅值和初始相位)，扩展比例可取±20%～±30%. 需说明的是，在调试过程中当采用±50%、±30%和±20%这3组不同扩展比例进行辨识，得到了相同的辨识结果，但迭代次数分别为 26 541次、9 538次和4 076次。这表明，对于本文所提策略，扩展比例大小对计算速度影响较大，但对计算精度影响较小。扩展比例根据具体问题调试选取，调试的一般原则是在保证多次测试均能成功获得最优解的前提下，尽可能地减小其数值[16]。气动参数的搜索空间可由独立辨识结果的最大值和最小值确定，为增大全局最优解位于搜索空间的概率，该搜索空间可进行适当放大。因此，状态初值的搜索空间共有n×N1维，气动系数的搜索空间共有N2维。

由于差分进化算法具有一定的随机性，为确保搜索到全局最优解，当搜索空间维度较多时，应增大个体数量N，使种群在初始化搜索空间内的密度更大、分布更均匀，增大搜索到全局最优解的概率[16]。此外，由于独立辨识采用的是局部优化算法，其给出的初始搜索空间有可能未包含全局最优解，为解决这一问题，本文的策略是利用同一程序进行多次辨识，前几次辨识不对搜索空间施加边界约束，可使迭代向初始化搜索空间外进行。若前几次辨识结果差异很小且位于搜索空间内，该结果大概率为全局最优解，则后续几次辨识中可对搜索空间施加“出界反射”的边界约束，即当“变异”后的个体在某个维度上超出搜索空间时，则在该维度上将其“反射”回搜索空间内，以提高计算效率[16]；若前几次辨识结果差异较大或不在搜索空间内，后续若干次辨识仍然不对搜索空间施加边界约束，以增加搜索到全局最优解的概率。当多次辨识完成后，若辨识结果差异很小(如小于0.01%)，可认为辨识结果是全局最优解；若辨识结果差异较大，应选取物理意义正确且准则函数更小的解作为辨识结果，但该解是否为真正意义上的全局最优解，还需结合气动工程计算、计算流体力学数值计算等进行验证。

辨识过程中，局部优化算法和全局优化算法的收敛条件均为

Jk

(3)

式中：Jk为准则函数的k次迭代值；b1=10-3. 然而，实际中往往达不到这样的条件，所以可采用以下条件中的任意一种：

|Jk+1-Jk|

(4)

|θk+1-θk|

(5)

式中：Jk+1为准则函数的k+1次迭代值；b2=10-5，b3=10-7. (4)式代表准则函数已经收敛，而(5)式代表待辨识参数已经收敛。在局部优化算法中，达到任意一种收敛条件都说明迭代已经收敛，即使该解可能为局部最优解。

对于差分进化算法，个体数量不止1个，所以无法使用(5)式，且由于可能存在局部最优解，使得相差较大的θ能计算出很接近的J值，所以(4)式也不能单独使用。在联合辨识中，考虑到整个种群应在所有维度上收敛，即有

(6)

式中：D()表示取方差；E()表示均值；i表示维度；b4=10-3.(6)式代表种群中所有个体在第i个维度达到收敛，当n×N1+N2个维度都满足(6)式时，可认为种群在所有维度上收敛。当可能存在较多局部最优解时，可以考虑同时使用多个种群进行辨识。当辨识结果满足(5)式而不满足(6)式时，说明辨识结果可能是局部最优的；只有同时满足(4)式和(6)式，才能说明辨识结果是全局最优的。

综上分析，在独立辨识中，可以同时使用(3)式、(4)式和(5)式，满足任意1个条件即认为计算收敛；在联合辨识的差分进化算法中，需要同时满足(4)式和(6)式，才能认为收敛到全局最优解。b1、b2、b3和b4可视具体情况而定，基本原则是：当弹丸发数较少或测量点较少时，为保证充分收敛，b1、b2、b3和b4的取值应更小(如比上述所给数值再小1个数量级)；当弹丸发数较多或测量点较多时，为优化计算时间，b1、b2、b3和b4的取值可适当放大，但一般不能大于上述所给数值。

3 一个典型气动力矩参数辨识问题

为了验证本文建立的多发弹测量数据联合辨识策略，本节考虑一个典型的气动力矩参数辨识问题[17]，即利用纸靶试验(测攻角)辨识弹丸静力矩系数导数、俯仰阻尼力矩系数导数及马格努斯力矩系数导数，将上述辨识流程应用于实测纸靶数据处理。

3.1 气动辨识用数学模型

由于纸靶试验往往采取平射，单发弹丸飞行时间较短，可认为全弹道马赫数基本不变、转速无衰减；同组试验的多发弹之间初速相差很小，认为具有相同马赫数。可采用复攻角运动方程作为气动力矩参数辨识用数学模型。

弹丸复攻角运动方程的齐次形式[18]如下：

(7)

式中：Δ为复攻角；s为弹道弧长；i为虚数单位；H为角运动阻尼作用项；P为转速作用项；M为静力矩作用项；T为升力和马格努斯力矩耦合作用项；α为纵向攻角；β为横向攻角；Clα为弹丸升力系数导数；CD为阻力系数；ρ为大气密度；d为弹径；m为弹丸质量；l为特征长度；Iy为赤道转动惯量；Ix为极转动惯量；p为弹丸转速；v为弹丸运动速度；CMα为静力矩系数导数；CMqα为俯仰阻尼力矩系数导数；CMpα为马格努斯力矩系数导数。

该方程的解析解[18]为

(8)

式中：KS0和KF0分别为慢圆运动和快圆运动的初始幅值;φS0和φF0分别为慢圆运动和快圆运动的初始相位;λS和λF分别为慢圆运动和快圆运动的阻尼指数;φ′S和φ′F分别为慢圆运动和快圆运动的角频率。

旋转稳定弹的陀螺稳定因子Sg为

(9)

工程上一般取Sg>1.3为陀螺稳定性条件，只有满足该条件，才能保证弹丸形成周期运动[18]。因此，陀螺稳定性条件即可作为气动辨识的约束条件：

P2-5.2M>0，

(10)

在整个辨识过程中，P和M的取值需始终满足约束条件(10)式。

因此，计算对应s处的α和β时，先利用(7)式计算出H、P、T和M，并检验其是否满足约束条件(10)式：若满足约束，则将其代入(8)式中进行计算，若不满足约束，则对CMα、CMqα和CMpα重新取值，直到其满足约束条件。

3.2 试验条件及测量数据

在某大口径榴弹的纸靶试验中，共射击1组6发弹丸。试验前，对每发弹进行了静态参数测量，且每发弹的纸靶布置方式相同：第1张纸靶布置在距离炮口25 m处，测量段s为25～125 m，每间隔5 m布置1张纸靶，弹丸平均飞行时间t≈0.11 s，平均马赫数Ma=2.679 2，弹丸初速的最大值与最小值之差为3.1 m/s，符合3.1节中描述的假设。试验完成后读取纸靶数据并对其进行曲线平滑和野值剔除，处理后每发弹有15～17个测量点。每发弹的m、Ix、Iy、和v均采用实测值；由于没有直接测量转速p，故将p作为待辨识参数。弹道测量数据为α、β和s，辨识中p和v作为常数处理。

表1中给出了独立辨识中的48个待辨识参数。

表1 独立辨识的待辨识参数

而联合辨识中的待辨识参数为

共有33个待辨识参数。

3.3 辨识的收敛性检验

利用数值仿真MATLAB软件编程并对多发弹测量数据进行辨识，为充分地利用计算资源，使用软件中的Parallel Pool模块进行多线程并行计算。在独立辨识中，使用不同线程对多发弹并行计算；在联合辨识中，使用不同线程对多个个体并行计算。

在独立辨识中，使用文献[19]中方法进行辨识，6发弹各自的准则函数收敛曲线如图3所示。

图3 独立辨识收敛性Fig.3 Convergence of independent identification

由图3可以看出，6发弹的准则函数都不满足(3)式，而是满足(4)式或(5)式。由于预先给定的初始相位、幅值和气动参数与实际值(真值)相差较大，而迭代增量由导数确定，故图3中各条曲线的截距不同，而斜率相同。

搜索空间维度D=33，根据2.2节中的辨识策略，选取个体数量N=2 000，并使用同一程序连续进行10次辨识。前5次辨识没有进行边界约束，其辨识结果均位于初始化搜索空间内且相对误差小于0.01%，可认为该结果很可能为全局最优解，后5次辨识施加“出界反射”的边界约束以提高计算效率。10次辨识结果的相对误差小于0.01%，故认为该辨识结果即为全局最优解。

准则函数最大值Jmax、最小值Jmin、平均值Jmean和方差D(J)能够反应整个种群在迭代过程中的收敛情况。Jmin所对应个体是潜在的最优解，而Jmax与Jmean则能代表种群的初始分布情况与收敛速度。图4给出了最后1次辨识中J和D(J)随迭代次数变化的曲线。

图4 联合辨识收敛性Fig.4 Convergence of combined identification

由图4可以看出：Jmax在迭代之初的值很大，但是下降速度非常快，k=1 000时几乎与Jmean和Jmin重合，且D(J)<10-3，准则函数基本收敛；k=2 000之后3条曲线就完全重合；k=3 000时，D(J)<10-5，整个种群的准则函数已经完全收敛。计算完成后，Jmin对应的气动参数与状态初值即为多发弹联合辨识的全局最优解。

除了J需要收敛外，种群中的所有个体θ也需收敛。种群位置的关系与辨识的收敛性相关，种群位置越集中，收敛性越好。图5给出了不同k时种群在气动参数搜索空间中的位置，为方便观察种群位置变化，在图5(b)～图5(f)中将其投影到CMα、CMpα和CMqα相互正交形成的3个面上。

图5 联合辨识气动参数收敛性Fig.5 Convergence of aerodynamic parameters for combined identification

从图5(a)～图5(f)可以看出各气动参数的迭代过程：当k=1时，种群的初始分布较为均匀，能充分地填充整个搜索空间；当k=179时，种群已向某一区域收敛；当k=1 000时，相比于初始搜索空间，整个种群几乎收敛为一个点；随着k的增大，种群的收敛性越来越高，但是在k=2 000之后，收敛方向发生改变。当搜索空间内仅有1个最优解时，在达到一定迭代次数后，种群的收敛方向应仅向内部收缩，而种群位置几乎不会移动[16]，如果种群位置发生移动，说明J的梯度方向发生改变，这种改变可能陷入局部最优解。

从图5(b)和图5(c)的投影可以看出，从k=179到k=1 000的收敛过程即为整个种群向内部收缩的过程；从图5(d)～图5(f)的投影可以看出，虽然从k=2 000到k=4 000同样是整个种群向内部收缩，但种群整体的位置也发生了改变，说明该位置附近可能存在局部最优解。在多次辨识中均出现了这种情况，说明当多发弹同时辨识时，若使用局部优化算法或搜索能力不足的全局优化算法时，很容易陷入局部最优解[3]。

在联合辨识中，由于无需雅克比矩阵，且搜索空间的选取较为合理，辨识过程中始终满足(10)式的约束，未出现计算发散；而独立辨识中却出现了发散，人为减小步长后才能收敛。这表明全局优化算法的稳定性要优于局部优化算法。

3.4 辨识结果与分析

为便于讨论，以单发弹独立辨识结果重构的弹道称为独立辨识弹道，以联合辨识结果重构的弹道称为联合辨识弹道。理论上讲，1种弹形对应唯一的气动参数，实际中由于多发弹辨识出的气动参数在数值上并不完全一致，在对辨识结果的处理中，往往取气动参数的平均值或者加权平均值，作为该弹的气动参数[12-13]。将独立辨识结果中的气动参数取平均值，利用该气动参数与辨识结果中的状态初值重构的弹道称为平均值计算弹道，其准则函数计算方法与联合辨识弹道相同。

由于单发弹的独立辨识不具有全局性，故以平均值计算弹道和联合辨识弹道为主要研究对象，独立辨识弹道作为参考，Jt为单发弹辨识弹道准则函数，Jm为平均值计算弹道准则函数，Jc为联合辨识弹道准则函数，独立辨识结果与平均值计算结果如表2～表4所示。

表2 状态初值的独立辨识结果

表3 气动参数的独立辨识结果

表4 气动参数平均值计算结果

从表3和表4中可以看出：CMα的均值与极差分别为3.428与0.46，极差与均值之比为13.4%，CMqα的极差与均值之比为121%，CMpα的极差与均值之比为130%；平均值计算弹道的准则函数Jm=4.996 5，远大于独立辨识弹道的准则函数之和Jt=1.685 1，这说明独立辨识弹道的准则函数虽然较小，但多发弹辨识结果相差较大，故平均值气动参数计算弹道与测量值相差较大。

利用独立辨识结果生成初始搜索空间进行多发弹联合辨识，经反复调试，选取状态初值和气动参数的扩展比例均为±20%. 结果如表5和表6所示。

表5 状态初值的联合辨识结果

表6 气动参数的联合辨识结果

对比表3、表4和表6可以看出，尽管两种方法辨识出的状态初值和气动参数相差不大，但联合辨识弹道的准则函数值Jc=2.089 8却远小于平均值计算弹道的准则函数Jm. 需要说明的是，这6发弹独立辨识的准则函数值并不相同，表3中给出的是6发弹的准则函数之和。为便于研究，这里选择独立辨识结果中准则函数最小的第4发弹作为研究对象，研究攻角变化规律与辨识误差。

将第4发弹的独立辨识弹道、联合辨识弹道、平均值计算弹道与测量值进行对比，结果如图6所示，图6中δ为攻角。测量段s取25～125 m，为便于研究，2条辨识弹道和平均值计算弹道范围取s为0～225 m，由于弹丸初速和转速较高，该段弹道上的速度变化约为1%，转速变化约为0.5%，工程上可作为常数处理。

图6 辨识弹道与平均值计算弹道Fig.6 Identified and calculated trajectories

从图6中可以看出：独立辨识弹道和联合辨识弹道与测量值较为接近，而平均值计算弹道与测量值相差较大；3条曲线在第1个攻角周期时相差不大，峰值位置和高度也很接近，随着距离的增加，差异逐渐变大，这说明在待辨识参数中，KS0、KS0、φS0和φF0的误差对弹道影响较小，而CMα、CMqα、CMpα和p的误差对弹道影响很大。

由于测量值间隔较大，且在攻角峰值和零点处缺少测量值，而该发弹的独立辨识弹道与测量值最接近，故以图6中独立辨识弹道为基准，研究第4发弹联合辨识弹道与平均值计算弹道的运动规律。

在测量段内，联合辨识弹道每个周期的运动距离与独立辨识弹道相差0.9%，前3个周期攻角峰高度值分别相差0.118°、0.016°和0.064°，峰值位置分别相差0.35 m、0.06 m和0.49 m；平均值计算弹道每个周期的运动距离与独立辨识弹道相差2.3%，并从第2个周期开始，攻角峰值高度、位置开始出现偏离，这种偏离随着距离增加，前3个周期攻角峰值高度分别相差0.066°、0.174°和0.256°，峰值位置分别相差0.73 m、1.78 m和2.83 m. 在测量段之后，由于误差逐渐累积，在第5个运动周期时，平均值计算弹道对应的峰值位置将提前4.93 m，对于本次试验，这意味着将提前1张纸靶出现攻角峰值。平均值计算弹道第4、第5个峰值高度与独立辨识弹道的峰值高度之差达到15%和19%，而联合辨识弹道仅有6%和8.8%.

以攻角计算值δ与对应弹道弧长处的攻角测量值δe之差作为绝对误差Δδ，以绝对误差与测量值之比Δδ/δe作为相对误差，具体如表7所示。

表7 攻角误差

从表7中可以看出，独立辨识弹道的误差最小，联合辨识弹道误差次之，平均值计算弹道误差最大。相对误差最大的点出现在δ=0°附近，因而该点的相对误差较大，绝对误差较小。独立辨识弹道和联合辨识弹道的最大绝对误差分别为0.213°和0.209°，而平均值计算弹道的最大绝对误差达到了0.631°，由于攻角最大测量值仅为2.9°，故该误差不容忽视。攻角相对误差和绝对误差随距离的变化如图7所示。

图7 攻角误差Fig.7 Errors of angle of attack

从图7中可以看出，独立辨识弹道、联合辨识弹道和平均值计算弹道相对误差的2个峰值分别出现在s=45 m和s=90 m处，这两处均在零点附近，除此之外，最大相对误差分别为14.7%、14%和40.6%. 最大绝对误差分别占攻角最大测量值的7.34%、7.21%和21.76%，结合表7中的误差大小可知，平均值计算弹道的误差约为其他2条弹道的2～3倍。

对于第4发弹丸，独立辨识弹道和联合辨识弹道与测量值的误差较小，而平均值气动参数计算的弹道误差则相对较大。其他5发弹与第4发弹类似，只是误差数值大小不同。对于全部弹丸，Jm比Jc大139.1%，因此，联合辨识所得气动参数更接近实际，全局性更好。

值得注意的是，在气动辨识中，只有同时计算出CMα、CMqα、CMpα和转速p，才能得到准确的弹丸角运动周期和阻尼指数。而KS0、KS0、φS0和φF0仅与试验条件有关，每发弹丸都可能有所不同，故只有同时辨识出以上所有参数，才能得到全局最优解。因此可以认为，只有联合辨识的结果为全局最优，平均值气动参数不是全局最优。

由于加工制造误差及发射条件不完全相同，每发弹丸在实际飞行过程中的真实气动特性不可能完全相同，且由于飞行环境差异、测量误差等影响，每发弹丸的辨识结果仅是对该发弹实际条件下的最优估计，不一定是真值，而本文对于多发弹联合辨识结果，实际上是对这一型或这一批次弹丸气动特性的最优估计，对于未进行飞行试验的同型弹丸来说，该辨识结果是对其气动特性的最可靠预测。

4 结论

本文提出了一种利用全局优化算法对多发弹测量数据进行联合辨识的策略，设计了对应的辨识流程。利用某纸靶试验中多发弹测量数据对该辨识策略进行了验证，并与现有的气动辨识方法进行对比，得出以下主要结论：

1)应用本文辨识策略不仅获得了弹丸气动参数与状态初值的全局最优解，还通过利用局部优化算法构建初始搜索空间，进一步提高了计算效率和稳定性，有效地解决了多发弹测量数据联合辨识气动参数的局部最优问题，所得辨识结果是对该型弹丸气动特性的最可靠预测。

2)与现有方法相比，应用本文辨识策略针对多发弹测量数据辨识，所得准则函数更小，重构的弹道与测量值更接近，能更为准确反映弹丸运动规律。

3)差分进化算法应用于多发弹测量数据联合辨识气动参数，具有较好的计算稳定性与收敛性，且算法本身无需求导和大矩阵求逆，进一步增强了计算稳定性，适于数据量较大的气动辨识问题，具有较好的实用性和通用性，为多发弹联合辨识气动力参数提供了新的思路。