火炮高平机参数辨识与灵敏度

林通，钱林方，付佳维，王明明

(南京理工大学 机械工程学院，江苏 南京 210094)

0 引言

高平机具备火炮高低机和平衡机的功能，结构紧凑、操作便利，广泛应用于现代火炮设计中。高平机在赋予火炮高低射角的同时也承受起落部分质量和部分发射载荷，其动力学特性与起落部分的俯仰刚度密切相关，对火炮射击精度等性能有重要影响。因此，研究高平机动力学模型，对火炮全炮发射动力学模型精度的提高及火炮射击精度的分析有重要意义。

火炮发射时起落部分相对耳轴的不平衡力矩，使得起落部分相对耳轴旋转抖动。此时高平机液压锁闭锁，高平机上行腔和下行腔油液不再流动，腔室内液压油承受火炮发射载荷，高平机的支撑力使得身管在发射过程中保持预定射角。目前，国内外众多学者对液压油压缩特性及油缸建模进行了研究。在油液压缩特性方面，冯斌[1]分析了闭式液压系统液压油中溶解及游离气体含量引起的含气油液弹性模量变化。Bureek等[2]考虑未溶解空气的影响，利用压缩法测定纯油液体积模量，并对含气油液压缩过程中压强随时间的变化规律进行了数学模拟。魏超等[3]在分析4种基于集中参数法的液压油有效体积模量的稳态模型基础上，推导得出考虑含气量时变特征的有效弹性模量动态模型。唐东林等[4]依据油液中气相成分随压强的变化过程，推导出含气油液有效体积弹性模量的理论模型，并通过试验数据对模型进行了验证。在油缸建模方面，李士军等[5]不考虑油液体积模量的变化，基于油液刚度和油缸缸体刚度提出了油缸等效力学模型。冯豪[6]建立了多因素影响下闭锁油缸刚度的数学模型，并给出了各因素产生的刚度占油缸总刚度的百分比。韩贺永等[7]基于含气油液的有效弹性模量，分析了油液含气量和油缸活塞位移对油缸等效刚度的影响。目前，在常规火炮发射动力学建模中，高平机通常简化成定刚度弹簧阻尼，然而在火炮点火发射后的数秒内，高平机油缸液体压强在0至几十兆帕内反复变化，油缸刚度存在明显非线性，采用估算等效刚度的计算结果难以较好地反映实际高平机的运动。因此，如何准确快速地建立高平机模型、提高火炮发射动力学模型的精度，从而更好地反映真实火炮发射过程的动力学响应，是研究的重点之一。

本文以某火炮高平机为研究对象，在分析高平机结构、受力的基础上，建立考虑含气油液刚度和油缸膨胀刚度的高平机理论模型，基于多体动力学仿真软件ADAMS完成火炮高平机参数化建模。基于某车载炮射击试验高平机测试数据，辨识获得高平机模型参数，验证高平机模型及参数辨识结果的准确性；通过对高平机参数的灵敏度分析，获得影响火炮起落部分俯仰运动的主要因素。

1 高平机动力学模型

1.1 高平机结构描述

某火炮高平机及其结构如图1所示，高平机分为高平机油缸筒与油缸杆两部分，油缸筒与摇架铰接，油缸杆与上架铰接，通过筒杆之间的伸缩运动实现火炮高低射角的调整。图1中，高平机油缸筒由高平机外筒和内筒组成，油缸杆为高平机中筒，Dm、Do和dm分别为中筒内径、外筒内径和中筒外径，A、B和C为外筒、中筒、内筒将高平机内部分成的3个腔室，即高平机的上行腔、下行腔和平衡腔。通过控制A腔、B腔的油液进出实现火炮高低调炮，C腔与蓄能气相连，提供平衡力矩。

图1 高平机结构Fig.1 Structure of elevating equilibrator

基于高平机油缸结构，计算得到A腔、B腔和C腔工作面积分别为

(1)

(2)

(3)

记A腔、B腔和C腔的初始压强分别为pA0、pB0和pC0，A腔、B腔初始体积分别为VA0、VB0.火炮发射时，受反后坐力及后坐部分往复运动影响，起落部分相对上架绕耳轴发生转动。此时高平机双向液压锁关闭，液压缸液压回路处于闭锁状态，闭锁的A腔、B腔起支撑和固定射角作用，其中油液受到反复拉伸压缩。油缸中油液受压时体积缩小、压强增大，同时缸筒受油液压强增大影响发生膨胀；油液受拉时体积增大、压强降低，同时缸筒收缩。

1.2 含气油液刚度

纯油液在一定压强作用下体积减小，表现出可压缩性，其规律可表示为

(4)

式中：p为压强；Vl为纯油液体积；El为纯油液体积模量，

(5)

纯油液的体积模量约为1.2～2.1 GPa，受压强和温度的影响很小，本文将其视为定值El=1.8×109Pa[1-2，8]。由(4)式得到纯油液在压强p作用下的体积Vl表达式为

(6)

式中：Vl0为初始油液体积；p0为大气压强。

由于实际工程中液压油往往会混入一定量气体，当油缸压缩时气体体积也会相应减小。气体压缩过程视为绝热过程，满足如下状态方程：

(7)

式中：Vg为压强为p时的气体体积；λ为气体绝热指数；Vg0为大气压下气体体积；Vgd为压缩过程溶解于油液的气体体积。

大气压下油缸油液中气体体积占总体积V0的体积分数为该含气油液的含气量，记为α，有

Vl0=(1-α)V0，

(8)

Vg0=αV0，

(9)

则压强为p时油液的总体积V为

(10)

由于空气在新的压强达到溶解平衡需要一定时间[1]，而火炮在发射后2 s左右已恢复平稳，忽略空气溶解项，则含气油液的正切体积模量可以表示为

(11)

由(11)式可知，含气油液的正切体积模量随压强变化而变化，与含气量也有直接关系。图2给出了不同含气量含气油液的正切体积随压强的变化。

图2 不同含气量油液正切体积模量随压强的变化Fig.2 Changing curves of tangent bulk modulus with pressure under different air contents

记A腔油液在压强变化dp作用下体积发生的变化为dVA，引起油缸活塞的位移为dδAl(压缩为正)。此时A腔含气油液体积为VA、等效刚度为KAl，则有如下关系式：

(12)

式中：EA为A腔含气油液的正切体积模量。由(11)式可知：

(13)

式中：αA为A腔含气油液的含气量；pA为A腔的压强。

整理(12)式，可得A腔液压油缸的等效刚度KAl为

(14)

式中：δAl油液体积变化引起的活塞压缩量；V′A0为A腔含气油液在大气压下的体积，可由A腔初始压pA0、含气量αA和初始体积VA0计算：

(15)

同理，可得B腔液压油缸的等效刚度KBl为

(16)

式中：EB为B腔内含气油液的正切体积模量；VB为B腔含气油液体积；δBl为油液体积变化引起的活塞压缩量；V′B0为B腔含气油液在大气压下的体积。

1.3 油缸膨胀刚度

闭锁液压油缸在油液压力作用下发生膨胀，特别对于初始体积较小的B腔，油缸膨胀体积的影响已无法忽略。对于无杆腔缸筒，在压强变化dp作用下，径向膨胀量dDc计算公式[2，6]为

(17)

式中：Dc为缸筒内径；Ec为缸筒材料的杨氏模量；dc为缸筒外径；ν为缸筒材料的泊松比。则缸筒膨胀体积dVc为

(18)

式中：L为油缸有效长度。

(19)

由此可得油缸膨胀产生的刚度Kc表达式为

(20)

从(20)式可以看出，对于油缸长度变化很小的闭锁油缸，缸筒膨胀刚度可看成定值。

1.4 高平机模型

在ADAMS软件中，将高平机A腔、B腔等效为油液等效弹簧和油缸等效弹簧串联(见图3)，串联弹簧总变形量由高平机筒与杆相对运动决定。图3中KAc、KBc分别为A腔、B腔的缸筒膨胀刚度。

图3 液压油缸等效模型Fig.3 Equivalent model of hydraulic cylinder

由(13)式、(14)式可知，油液等效弹簧为变刚度弹簧，其刚度受到油液压强的影响，积分可得A腔油液等效弹簧的作用力FA为

(21)

式中：x为A腔油液等效弹簧的压缩量。

由于EA随p的变化而变化，积分需要通过迭代得到。为了便于后续模型参数化以及不同含气量的模型计算，引入模量退化因子K，将(21)式重写为如下形式：

(22)

式中：K为大于0且小于等于1的无量纲数，K=1时，(22)式表示含气量αA=0(即EA=El)的纯油液压缩。在纯油液体积模量El已确定的情况下，因子K只与含气量α和体积变化率V/V0有关，通过数值计算可得到K随α、V/V0变化的曲线如图4所示。模型计算时，可根据αA、VA/V′A0插值得到K值。

图4 K随体积变化率和含气量变化曲线Fig.5 Changing curve of K with volume change rate and air content

同理，可得B腔产生的作用力FB为

(23)

由于C腔油液与蓄能器中的气室相连，与A腔、B腔相比，C腔在较小体积变化下的压强变化量可忽略，将火炮发射时C腔的作用力FC考虑为常数：

FC=pC0SC.

(24)

由以上推导可得高平机提供的支撑力F表达式为

(25)

基于上述分析，在ADAMS软件中创建11个设计变量，分别对应高平机理论模型中的αA、αB、VA0、VB0、pA0、pB0、SA、SB、KAc、KBc和D，以便参数辨识和灵敏度分析时调整参数。

2 高平机模型参数辨识

2.1 问题描述

在高平机模型中，A腔、B腔油液含气量αA、αB、和高平机阻尼D无法直接通过试验数据得到。由于A腔、B腔结构复杂，其缸筒膨胀刚度KAc、KBc只能粗略估算。因此，将上述参数确定为待辨识的参数，需要通过高平机的试验数据建立辨识模型反求得到。将待辨识的参数写成向量的形式：

I=(αA,αB,KAc,KBc,D)T.

(26)

将参数代入模型，获得载荷作用下的高平机A腔、B腔压强随时间的响应历程pA(I)和pB(I).通过高平机上安装的压力传感器，测得载荷作用下A腔、B腔压强随时间变化的数据A和B.定义辨识误差为

(27)

式中：n为数据点个数；Ai和pAi(I)分别为A腔试验和仿真压强第i个数据点的值；Bi和pBi(I)分别为B腔试验和仿真压强第i个数据点的值。

将高平机模型参数辨识问题看作一个优化问题，以模型仿真的输出值与真实系统的测量值最接近为目标，利用寻优方法，在待辨识参数的区间范围内寻找出最优的值。建立高平机模型参数辨识问题如下：

搜索:I，
minJ(I)，
s.t.Imin≤I≤Imax，

(28)

式中：Imax、Imin为待辨识参数的上界和下界。

2.2 全炮发射模型及测试

本文高平机的参数辨识基于某车载炮射击试验的高平机测试数据，车载炮发射状态如图5(a)所示。车载炮结构复杂，将其实际结构适当简化，可获得较理想的多刚体动力学模型，相应的拓扑结构如图5(b)所示。图5(b)中：h1、h2、h3为大架、座盘、千斤顶与土壤间的非线性接触力；h4为大架与车体间的旋转副和等效弹簧阻尼；h5、h6为座盘和千斤顶与车体间的滑移副和等效弹簧阻尼；h7为车体与上架间的旋转副和方向机力；h8为上架与摇架间的旋转副；h9为高平机杆与上架间的旋转副；h10为高平机筒与杆间的滑移副和高平机力；h11为高平机筒与摇架间的旋转副；h12为后坐部分与摇架间的滑移副和反后坐力。

图5 火炮发射动力学模型Fig.5 Gun launching dynamic model

模型中炮膛合力随时间变化的数据预先算出，沿着炮膛轴线施加在炮闩上，推动后坐部分后坐。制退机力、复进机力为后坐位移和后坐速度的函数，由对应公式算出并分别施加在制退机和复进机轴线上。高平机筒与杆间加入高平机模型，并输出高平机A腔、B腔压强随时间变化数据。

火炮发射时，后坐部分在炮膛合力、制退机力和复进机力作用下完成后坐复进运动。由于相对耳轴的不平衡力矩以及后坐部分质心位置变化，起落部分会相对耳轴发生转动，并在高平机的作用下来回旋转摆动，最终恢复平衡。在此过程中，高平机A腔和B腔的压强反复变化，在高平机A腔、B腔安装压力传感器，测量压强随时间变化数据。

2.3 算例分析

基于2.1节的辨识方法，对某高平机参数进行辨识。高平机基本参数如下：SA=4.22×10-3m2、SB=4.22×10-3m2、VA0=6.70×10-3m3、VB0=1.59×10-3m3，El=1.8×109Pa，Ff=1 100 N.以第1组实测A腔、B腔的初始压强pA0=1.52 MPa、pB0=3.43 MPa作为模型初始条件。考虑到本文优化问题参数较少、区间连续，引入收敛速度较快的粒子群优化算法进行参数辨识[9-13]。种群规模20，最大迭代次数50.待辨识参数取值范围如表1所示。

表1 待辨识参数区间范围Tab.1 Range of parameters to be identified

参数辨识结果如表2所示。将表2参数代入高平机模型，仿真获得A腔、B腔压强数据随时间变化的曲线。模型仿真值与实测值以及不考虑含气量的纯油液模型计算结果对比如图6和图7所示。

表2 参数辨识结果Tab.2 Identificated results of paramenters

图6 第1组高平机A腔压强对比Fig.6 Comparison of pressures in Chamber A:Group 1

图7 第1组高平机B腔压强对比Fig.7 Comparison of pressures in Chamber B:Group 1

从图6和图7中可以看出，由于火炮起落部分的俯仰运动，高平机反复压缩拉伸，高平机A腔、B腔压强在火炮发射过程中反复变化。A腔受压压强升高的同时B腔被拉伸压强下降接近真空，反之亦然。初始体积较小的B腔压强变化较大，两腔压强在来回两次较大波动后趋于平稳。从图6和图7中3条曲线对比结果可看出，试验测试值与本文模型计算结果能较好吻合，曲线幅值与相位基本一致。如果不考虑油液含气量，由于高平机模型刚度偏大，仿真曲线较试验曲线幅值偏大、相位超前。对比结果表明考虑含气油液刚度和油缸膨胀刚度的高平机模型，能更好地模拟高平机内压强变化。

在相同射击条件下，将第2组测试数据A腔、B腔初压pA0=1.06 MPa、pB0=2.53 MPa作为高平机模型初始条件，将表2参数代入高平机模型，得到仿真与试验曲线如图8和图9所示。由图8和图9可以看出，对于不同A腔、B腔初始压强的仿真与试验曲线均能较好地吻合，验证了模型的精度，同时表明采用粒子群优化算法辨识获得的高平机参数具有较高精度和适用性。

图8 第2组高平机A腔压强对比Fig.8 Comparison of pressures in Chamber A：Group 2

图9 第2组高平机B腔压强对比Fig.9 Comparison of pressures in Chamber B：Group 2

3 高平机对俯仰运动影响因素分析

3.1 问题描述

由高平机模型可知，高平机刚度受到结构参数、液压油含气量等因素影响。高平机刚度直接关系起落部分的俯仰运动，而俯仰运动间接反映了炮口扰动情况。因此，通过灵敏度分析获得高平机参数中对起落部分运动影响较大的关键参数，研究高平机参数名义值的调整对起落部分运动的影响规律，可为高平机设计提供理论支持。

3.1.1 设计变量

对表3中所列11个高平机模型参数进行灵敏度分析，其中：高平机A腔液压油含气量、B腔液压油含气量、A腔初压和B腔初压，基于试验数据给出其可调整的区间；其余高平机结构参数的范围由高平机结构设计时的可调整余量决定，基于初始值给予其上下浮动20%调整区间。参数具体范围如表3所示。

表3 高平机参数灵敏度分析采样范围Tab.3 Sampling range of elevating equilibrator paramenters for sensitivity analysis

3.1.2 目标函数

影响火炮射击精度的主要指标是起落部分的角位移和角速度，由于高平机影响起落部分的俯仰运动，将模型输出摇架在弹丸出炮口时刻的垂向角位移θ和垂向角速度ω作为目标函数。

3.2 灵敏度分析方法

运用Isight软件集成参数化的ADAMS动力学模型，对于每个参数，在表3所示取值范围内进行最优拉丁方采样，使获得的样本点尽可能均匀分布在设计空间内[14-15]。将样本点代入模型，获得对应的目标函数响应值。

为了公平地反映各个变量对响应的贡献度，Isight软件将变量区间归一化至[-1，1]。对于含有n个参数的模型，运用最小二乘法拟合构建基于样本点的二次回归模型：

(29)

式中：y为模型输出；xi和xj为模型参数；β0、βi、βii和βij分别为常数项、1次项、2次项和交互项系数。

(30)

(31)

(32)

通过对N的绝对值排序，可得到1阶项、2阶项和交互项中对响应影响较大的项。

3.3 结果分析

结合设计变量、目标函数、灵敏度百分比计算方法，通过300个样本点计算得到高平机参数对起落部分垂向角位移θ和垂向角速度ω的灵敏度百分比，图10和图11给出了影响最大的10项灵敏度百分比。

图10 高平机参数对θ的灵敏度Fig.10 Sensitivity of elevating equilibrator paramenters to θ

图11 高平机参数对ω的灵敏度Fig.11 Sensitivity of elevating equilibrator paramenters to ω

由上述灵敏度分析结果可知，对起落部分垂向角位移θ影响较大的有7个参数，灵敏度百分比从大到小依次为pB0、αB、pA0、αA、D、SB、VB0.对垂向角速度ω影响较大的同为上述7个参数，灵敏度百分比从大到小依次为αB、pB0、D、αA、SB、VB0、pA0.相对1阶项而言，平方项和交互项对目标函数作用均不明显，只有pB0和αB交互项灵敏度稍大。参数灵敏度正负值αA、αB、VB0对目标函数影响为正，pA0、pB0、D、SB对目标函数影响为负。即减小A腔、B腔油液含气量，缩小B腔油液初体积；同时增加A腔、B腔初压，高平机阻尼和B腔活塞工作面积有利于减小弹丸出炮口时起落部分的垂向角位移和角速度。

4 结论

本文描述了某火炮高平机的基本结构，基于含气油液压缩模型及油缸膨胀刚度建立了高平机动力学模型，并在ADAMS软件中实现了高平机模型参数化；提出了高平机模型参数辨识方法，基于某车载炮射击试验数据，引入粒子群优化算法对高平机参数进行了辨识。通过灵敏度分析，获得了对火炮起落部分运动影响较大的高平机参数。得出主要结论如下：

1)辨识后的高平机动力学模型仿真结果与试验测试结果吻合良好，表明本文提出的基于含气油液刚度和缸筒膨胀刚度的高平机模型能较好地反映高平机在火炮发射状态下的动力学特性。

2)从仿真结果与试验数据可看出，高平机油液压强在火炮发射时变化剧烈。A腔、B腔油液反复受压提供支撑力使火炮恢复平稳，初始体积较小的B腔压强波动较大。

3)灵敏度分析结果表明，高平机参数中pB0、αB、pA0、αA、D、SB、VB0对弹丸出炮口时火炮起落部分的运动影响较大。减小A腔、B腔油液含气量，增加初压；增加B腔工作面积，减小初体积；增加高平机阻尼都能有效减小弹丸出炮口时起落部分垂向角位移和角速度。