火炮反后坐装置安装在炮身和架体之间，用于提供弹性力和制动力控制后坐部分在火炮射击时的后坐运动，并使之复位^[1-2]。在此过程中，反后坐装置减小了火炮架体在射击时的受力，还把射击时的全炮后坐运动变为可控的炮身后坐运动，并能自动复位。反后坐装置中的参数包括基本物理参数、结构参数和经验系数。一些参数基于实验或经验的系数不易准确选定，因而对这些参数进行辨识具有一定的实用价值。杜中华等^[3]采用数值仿真的方法对某型反后坐装置的液压阻力以及密封件工作压强进行了分析；狄长春等^[4]对反后坐装置后坐复进过程中复杂流场性质及结构场动态特性进行了数值仿真。此外，火炮反后坐随着机械系统呈现出高速化、大型化和复杂化的方向发展，常规研究手段已经无法准确描述系统中存在的大量不确定性参数。对于不确定参数的优化方法可分为3类：概率不确定性优化方法、非概率不确定性优化方法和概率-非概率混合不确定性优化方法。工程实际中，往往只能得到有限的样本数据，很难获得不确定参数的精确概率分布，所以概率不确定性优化方法在工程实际中存在局限性。相对于概率不确定性优化方法，非概率不确定性优化方法具有对数据要求低，计算简单等优点，越来越多的科研人员对此进行了研究。其中采用区间来描述事物的本质和特征的区间数理论方法被应用于非概率不确定性优化。文献[5-6]分别提出了基于区间分析的方法进行参数估计和结构辨识；姜潮^[7]提出了基于区间的不确定性优化理论与算法；Jiang等^[8]提出一种改进的区间可能度计算方法提高计算精度；赵抢抢等^[9]利用区间序关系转换模型将区间不确定参数辨识问题转换成确定性优化问题，并结合差分进化算法实现了弹药协调器的参数辨识；王敏容等^[10]以可靠性指标为约束条件，提出了一种基于区间模型的结构非概率稳健优化设计方法。当用区间数描述系统的不确定参数时，由于目标函数相对于每一个设计变量的取值构成一个区间，所以辨识问题无法通过传统的确定性辨识方法求解。针对存在区间不确定参数的优化问题，可以通过非线性区间优化的数学模型转化为确定性问题进行求解。

本文为解决火炮后坐复进过程中反后坐装置的参数辨识问题，建立了火炮后坐复进过程的解析模型。为了提高辨识效率，同时建立了神经网络模型，并根据不确定参数的变化特性对相关性能参数进行了分类以及用区间数进行描述。针对在一次动作过程中可通过参数辨识的方法得到其确定值的第1类参数，利用区间可能度转换模型，将区间不确定性问题转换为确定性问题进行求解。以反映辨识结果与测试结果的相似程度的时间序列数据相似度作为优化目标，利用粒子群算法实现参数辨识。针对只能确定其变化区间的参数，利用粒子群算法结合区间优化方法优化其区间。

1 火炮反后坐装置运动过程数学建模 1.1 火炮复进机和制退机数学方程

火炮在后坐运动结束后，后坐部分在复进机力作用下回到待发射位置。取后坐部分为受力体, 后坐过程中的运动微分方程为^[1]：

$ m \mathrm{d} u / \mathrm{d} t=m\left(\mathrm{~d}^{2} x\right) /\left(\mathrm{d} t^{2}\right)=F_{p t}-F_{R} $

(1)

$ F_{p t}=F_{\phi h}+F_{f}+F+F_{T}-m g \sin \varphi $

(2)

式中：F_pt为炮膛合力；F_R为后坐阻力；$F_{\phi h} $为制退机力；F_f为复进机力；F为反后坐装置密封装置的摩擦力；F_T为导轨的摩擦力；u为后坐部分后坐速度；x为后坐部分后坐位移；m为后坐部分质量；φ为火炮射角。其中，$F_{f} 、F_{\phi h} $为：

$ F_{f}=A_{f} P_{f 0}\left(V_{0} /\left(V_{0}-A_{f} x\right)\right)^{n} $

(3)

$ \begin{array}{c} F_{\phi h}=\left( K_{1} \rho\left(_0^A -A_{p}\right) 3\right) /\left(2 a_{x}{}^{2}\right) u^{2}+ \\ \left(K_{2} \rho A_{f j}{} ^{3}\right) /\left(2 A_{1}{}^{2}\right) u^{2} \end{array} $

(4)

火炮复进过程中的运动微分方程为：

$ \begin{array}{c} m \mathrm{d} u / \mathrm{d} t=m\left(\mathrm{~d}^{2} x\right) /\left(\mathrm{d} t^{2}\right)= \\ F_{f}-F_{\phi h}-\left(F+F_{T}+m g \sin \varphi\right) \end{array} $

(5)

$ F_{f}=A_{f} P_{f 0}\left(V_{0} /\left(V_{0}-A_{f} x\right)\right)^{n} $

(6)

式中: A_f为复进机活塞工作面积；P_f0为复进机气体的初压力；V₀为复进机气体的初体积；n为复进机气体非线性指数。后坐时，由于制退杆从制退机内抽出而加大了制退机内的空间体积，使之形成制退机非工作腔内的真空。复进时，在制退机非工作腔真空消失以前，总的复进液压阻力中不含有制退机复进液压阻力；当制退机非工作腔的真空消失之后，总的液压阻力包含制退机复进液压阻力。

当$ 0 \leqslant \bar{x}<l_{p}$时，

$ F_{\phi f}=\left(K_{2} \rho A_{f j}^{3}\right) /\left(2 a_{f}^{2}\right) \bar{u}^{2} $

(7)

当$l_{p} \leqslant \bar{x} \leqslant l_{\lambda} $时，

$ \begin{array}{c} F_{\phi f}=\left(K_{2} \rho A_{f j}^{3}\right) /\left(2 a_{f}^{2}\right) \bar{u}^{2}+ \\ \left(K_{1} \rho\left(A_{o f}-a_{c}\right)^{2}\right) /\left(2\left(a_{x}+a_{c}\right)^{2}\right) \bar{u}^{2} \end{array} $

(8)

式中：l_p为制退机非工作腔真空消失点；l_λ为最大后坐行程；$\overline x $为以最大后坐点为起点的复进行程；$\overline u $为复进速度；ρ为液体密度；K₁为主流液压阻力系数；K₂为支流液压阻力系数；A_fj复进节制器活塞工作面积；A_of为复进时制退活塞工作面积；a_f为复进节制器流液孔面积；a_x为制退机流液孔面积，是$\overline x $的函数；a_c为2条纵向沟槽的截面积；文中的目标函数函数是通过数学模型所得的P₁腔室的仿真压力曲线与测试压力曲线之间的相似度。在发射试验中，后坐部分进行后坐复进运动，同时测得的P₁腔室的压力曲线。图 1为制退机结构及压力测试腔示意图。P₁为：

$ {P_1} = ({K_1}\rho {({A_0} - {A_p})^2})/(2a_x^2){u^2} $

(9)

式中：P₁为制退机P₁腔室的压力；A_o为制退杆活塞的工作面积；A_p为节制环孔面积

1.2 火炮反后坐装置待辨识参数分类

由1.1节可知，火炮反后坐装置系统包含参数众多。待辨识的参数有3个，分别是复进机气体初压P_f0、复进机气体非线性指数n和制退机主流液压阻力系数K₁。根据3个参数在发射过程中的变化特性分为2类。第1类参数是难以通过测量的方法获得，在动作过程中可近似看作确定的，可通过辨识的方法获得，即复进机气体初压P_f0。P_f0是复进机的基本参数，由于密封装置的泄露，结构尺寸误差等各方面原因，使得在发射过程中很难测得其准确值。本文根据在经验值的基础上适当增大范围的原则确定P_f0的初始区间。

第2类参数是区间不确定参数，无法确定其准确的分布规律，只能确定其变化的区间。对于火炮反后坐装置而言，2个区间不确定的参数分别为：复进机气体非线性指数n和制退机主流液压阻力系数K₁。复进机气体非线性指数n取决于复进机的散热条件和活塞运动速度。制退机液压阻力系数是反映制退机内液体真实流动因素的符合系数，这些因素有液体流动的沿程损失、流动的局部损失、流液孔处液流截面的收缩、液体可压缩性和流动的非定常流动等。本文在经验区间的基础上适当增大区间范围的原则确定K₁和n的初始区间，如表 1所示。

1.3 神经网络模型的建立

火炮后坐复进过程是复杂的运动过程，而辨识过程需要大量的循环计算，为了提高辨识效率，利用前述数学模型，建立神经网络模型。反向传播神经网络(BP神经网络)是一种多层前向反馈神经神经网络，具有强大的自适应、泛化和非线性逼近能力。BP神经网络通常采用梯度下降搜索算法，不断调整网络权值和阈值，使得网络的输出值与实际输出值的均方误差最小。文章在3个待辨识参数的初始区间进行拉丁超立法抽样1 000组，代入数学模型，得到1 000组P₁腔压力曲线，并计算1 000组仿真曲线与测试压力曲线的曲线相似度。以3个变量为输入，曲线相似度为输出，建立神经网络的1 000组样本，神经网络模型结构如图 2所示。

2 参数辨识问题的描述及求解方法 2.1 参数辨识问题描述

整个参数辨识过程分为2个步骤。1)对第1类参数的辨识可利用寻优的方法求解，即在主流液压阻力系数K₁和复进机气体非线性指数n在初始区间范围内变化的情况下，寻找出最优的复进机气体初压P_f0，使得辨识结果与测试结果的时间序列数据相似度区间的可能度最大。其数学模型为：

$ \left\{ {\begin{array}{*{20}{l}} {\mathop {\max }\limits_X P\left( {{f^I}(X,D) \ge {V^I}} \right)}\\ {D = \left( {{d_1},{d_2}, \cdots ,{d_l}} \right),{d_k} = \left[ {d_k^L,d_k^R} \right]} \end{array}} \right. $

(10)

式中：X=P_f0为设计变量，D=(K₁, n) 为不确定变量，由于D的取值为一区间数，因此对于任何一个确定的设计变量X，相似度f^I(X_i, D) 的取值构成一个区间数，且均在区间[0, 1]。将相似度区间转换成区间可能度$ P\left(f^{I}\left(X_{i}, D\right) \geqslant V^{I}\right), V^{I}$为性能区间，求解式后，可以得到一个最优设计向量X^*，使得目标函数的可能度取最大值$P_{\max }=P\left(f^{I}\left(\boldsymbol{X}^{*}, D\right) \geqslant V^{I}\right) $)。传统的确定性优化方法中，决策的判定都是基于目标函数在各个设计向量处的具体数值。针对此类目标函数为区间的参数辨识问题，本文提出了一种基于区间可能度数学转换模型将区间不确定优化问题转换为确定性问题的方法，再结合粒子群算法进行辨识。

2) 在辨识出第1类参数的基础上，对第2类参数的区间进行优化。每一组(K₁, n)对应一个相似度值f(K₁, n)，如图 3所示。

需要找到待辨识参数制退机的主流液压阻力系数K₁和复进机气体非线性指数n的最大区间，使得该区间内所有的相似度值均大于μ。其数学模型为：

$ \left\{ {\begin{array}{*{20}{l}} {\max {L_{ij \to rk}}(x,y) = {{\left( {{x_r} - {x_i}} \right)}^2} + {{\left( {{y_k} - {y_j}} \right)}^2}}\\ {{F_{ij \to rk}} = \min \left( {\mathit{\boldsymbol{A}}_{ij \to rk}^I(x,y)} \right) \ge \mu }\\ {x = \left( {{x_1},{x_2}, \cdots {x_m}} \right),{x_r} \ge {x_i}}\\ {y = \left( {{y_1},{y_2}, \cdots ,{y_l}} \right),{y_k} \ge {y_j}} \end{array}} \right. $

(11)

式中：$ x=k ; y=n$；$ \boldsymbol{A}_{i j \rightarrow r k}$为相似度矩阵；$ \boldsymbol{A}_{i j \rightarrow r k}^{I}$为相似度矩阵的区间；$ F_{i j \rightarrow r k}$为相似度矩阵区间的最小值$ L_{i j \rightarrow r k}$为相似度矩阵的长度。

$ \boldsymbol{A}_{i j \rightarrow r k}=\left[\begin{array}{cccc} f_{i j} & f_{(i+1) j} & \cdots & f_{r j} \\ f_{i(j+1)} & f_{(i+1)(j+1)} & \cdots & f_{r(j+1)} \\ \vdots & \vdots & & \vdots \\ f_{i k} & f_{(i+1) k} & \cdots & f_{r k} \end{array}\right] $

(12)

$ \boldsymbol{A}_{i j \rightarrow r k}^{I}=\left[\min \left(\boldsymbol{A}_{i j \rightarrow r k}\right) \max \left(\boldsymbol{A}_{i j \rightarrow r k}\right)\right] $

(13)

针对此类区间优化问题，本文采用式(11)的数学模型结合粒子群算法进行优化。

2.2 时间序列相似度计算

本文采用离散序列的一致性度量方法，即动态时间规整(DTW)计算仿真曲线与测试曲线的相似度来评价辨识结果。试验所得的压力数据序列长度与仿真所得的压力数据序列长度不一致，并且2组数据序列的时间轴无法完全对齐，所以一般的欧氏距离法存在不足。DTW方法用满足一定条件的时间规整函数描述两者之间的时间对应关系，它允许时间轴上的漂移。将DTW距离进行归一化处理，所得的值即为2组时间序列数据之间的相似度^[11-13]。动态时间规划的具体步骤为：

1) 假设2组时间序列分别为n维和m维，构造一个n×m的矩阵，用于存放两序列点对点之间的距离(一般可采用欧氏距离)，距离越小表示两点之间的相似度越高。

2) 把矩阵看成一个网格，算法的目的为找到一条通过此矩阵网格的最有路径，改路径通过的网格点即为两组时间序列经过对齐后的点对；

3) DTW算法定义一个归整路径距离，找到最优路径后，将所有相似点之间距离和进行归一化处理，来衡量2组时间序列之间的相似性。

2.3 区间可能度转换模型

根据区间数学，区间数被定义为一对有序的实数^[14]：

$ A^{I}=\left[A^{L}, A^{R}\right]=\left\{x \mid A^{L} \leqslant x \leqslant A^{R}, x \in \bf{R}\right\} $

(14)

式中上标I、L、R分别表示区间、区间下界和区间上界。区间的宽度定义为：

$ L\left(A^{I}\right)=A^{R}-A^{L} $

(15)

区间可能度方法是基于模糊集来构造区间可能度，是一种定量比较区间之间的优劣程度^[15]。将区间B视为文章辨识的性能区间，则区间A与区间B的位置关系如图 4所示。

记P(A^I≥B^I)为区间A^I≥B^I的可能度。本文采用了2种不同区间可能度的构造方法:

构造方法(Ⅰ):

$ \begin{array}{l} P\left(A^{I} \geqslant B^{I}\right)=\max \left[0, L\left(A^{I}\right)+L\left(B^{I}\right)-\right. \\ \ \ \ \ \ \ \ \ \left.\max \left(0, B^{R}-A^{L}\right)\right]\left({ }^{L}\left(A^{I}\right)+L\left(B^{I}\right)\right)-1 \end{array} $

(16)

构造方法(Ⅱ)：

$ P\left(A^{I} \geqslant B^{I}\right)=\left\{\begin{array}{ll} 1, & A^{L} \geqslant B^{L} \\ \frac{1}{2} \frac{\left(A^{R}-B^{L}\right)^{2}}{L\left(A^{I}\right) L\left(B^{I}\right)}, & A^{L}<B^{L}<A^{R} \\ 0, & A^{R} \leqslant B^{L} \end{array}\right. $

(17)

2.4 粒子群算法

本文对2类参数的辨识均采用了粒子群算法^[16-18]，该算法是智能优化算法的一种，有着精度高、收敛快的特点。粒子群算法优化的原理是通过不断更新粒子速度和位置，并计算其适应度值，直到适应度值符合要求为止。粒子群算法的流程如图 5所示。

根据2类参数的不同辨识要求，设置粒子群算法的参数，包括种群规模，最大迭代次数、学习因子、适应度函数以及粒子位置更新准则。

3 火炮发射试验与辨识结果分析 3.1 火炮发射试验

文中目标函数是通过数学解析模型所得的仿真曲线与试验曲线之间的相似度。测试曲线是在火炮发射试验过程中，测得的P₁腔室的压力曲线。图 6为试验过程示意图，图 7为测得的压力曲线。

3.2 辨识结果分析 3.2.1 第1类参数辨识过程及结果

设置粒子群智能算法的相关参数：种群规模为100;最大迭代次数为50；适应度函数即为区间可能度；学习因子c₁=0.5, c₂=0.5；更新权重w=0.8。粒子速度及位置更新准则为：

$ \begin{array}{c} v_{i d}(k+1)=w v_{i d}(k)+c_{1} r_{1}\left(p_{i d}-x_{i d}(k)\right)+ \\ c_{2} r_{2}\left(p_{g d}-x_{i d}(k)\right) \end{array} $

(18)

$ x_{i d}(k+1)=x_{i d}(k+1)+v_{i d}(k+1) $

(19)

式中: k是迭代次数；r₁、r₂均为[-1, 1]的均匀随机数；p_id、p_gd分别为个体最优值和全局最优值。分别采用2种区间可能度的构造方法对第1类参数进行辨识，性能区间V^I为[0.8, 1.0]。当采用第1种可能度构造方法时，可得复进机气体初压P_f0的最优值为7.588 2 MPa, 所用时间为6.708 s；当采用第2种区间可能度构造方法时，可得复进机气体初压P_f0的最优值为7.342 0 MPa，所用时间为6.563 s。然后，K₁和n在初始区间采样，分别同2种可能度构造法下得到的复进机气体初压P_f0最优值代入模型，可得P₁腔室压力曲线簇与测试曲线的对比，即如图 8和图 9所示。分别采用2种区间可能度的构造方法辨识出复进机气体初压值。由图 8和图 9可得，采用第2种可能度构造方法得到的复进机气体初压值代入神经网络模型，得到的P₁腔室压力曲线簇更加稳定，且与测试曲线的整体相似度更高。因此，复进机气体初压P_f0的最优值取7.342 0 MPa。以此为基础，进行第2类参数的区间优化。图 8中，该结果是在采用第1种可能度构造法的情况下获得。黑色曲线为测试曲线，其他不同颜色曲线为K₁和n在不同采样样本下得到的辨识曲线。图 9中，该结果是在采用第2种可能度构造法的情况下获得。黑色曲线为测试曲线，其他不同颜色曲线为K₁和n在不同采样样本下得到的辨识曲线。

3.2.2 第2类参数辨识过程及结果

设置粒子群智能算法的相关参数：种群规模为100；最大迭代次数为200；c₁=1.2, c₂=1.3；适应度函数为$L_{i j \rightarrow r k}(x, y)=\left(x_{r}-x_{i}\right)^{2}+\left(y_{k}-y_{j}\right)^{2} $；速度位置更新准则如式(18)和式(19)。μ分别取0.50、0.80、0.95进行优化，得到如表 2结果。

为验证该算法的优化精度，当μ取0.50时，参数K₁和n在优化后所得区间内随机采样1 000组，代入数学模型中，得到P₁腔室压力曲线簇，并计算相似度，得到相似度在优化后相似度区间范围内的概率为97.4%。同理，当μ分别取0.80、0.95时，参数K₁和n在优化后所得区间内随机采样1 000组，带入数学模型中，并计算P₁腔室压力曲线簇的相似度，可得相似度值在优化后相似度区间内的概率分别为97.2%和98.0%。3种情况相似度均较高，可得本文的算法具有较高的优化精度。

4 结论

1) 为了辨识火炮反后坐装置不确定性参数，文章建立了火炮发射过程中，反后坐装置的数学解析模型，并为提高辨识效率，建立神经网络模型。为了得到P₁腔室的压力曲线作为辨识过程中的目标曲线，进行了火炮发射试验。

2) 火炮复进机的气体初压属于第1类参数，针对此类不确定参数，利用区间可能度数学转换模型，将区间不确定优化问题转化为确定性优化问题进行求解。以测试曲线和数值仿真曲线的时间序列数据相似度作为优化目标，利用粒子群算法实现参数辨识。

3) 火炮制退机的主流液压阻力系数和复进机的非线性指数属于第2类参数。针对此类无法确定其准确的分布规律，只能确定其变化的区间，文章提出了一种基于粒子群算法结合区间优化算法优化其区间范围。辨识结果表明文章提出的同时针对2种不同类型参数的辨识方法具有可行性和有效性。