随着现代火炮向着高初速、高射频、高射击精度和高机动性方向发展，研究人员对于火炮射击时的稳定性和安全性提出了更高的要求[1]。身管振动引起弹丸出炮口姿态的扰动，最终影响射弹的散布规律，一直是学者们关注的热点。筋桶摇架结构作为承载身管的主要构件，其在冲击作用下产生的振动也会影响身管振动特性。以往的研究工作较少考虑火炮其他构件的运动对身管振动的影响。因此，综合考虑身管和筋桶结构之间的耦合作用，对于研究身管振动特性是十分必要的。
关于身管振动的研究，前人研究的重点多在于身管振动分析上。文献[2]对火炮进行模态分析，然后将炮身的后坐运动简化为简谐运动，对炮身部分进行谐响应分析，得到激励频率对炮口振动的影响规律。文献[3]针对弹丸膛内高速运动冲击下的身管振动问题，建立了移动弹丸作用下身管的横向振动模型，对不同弹速激励下的身管振动响应进行了数值仿真。文献[4]利用有限元仿真对某新型火炮身管组进行分析，分别考虑单炮箍及双炮箍对身管组模态频率的影响，得到身管组在2种情况下不同阶数的模态频率及振型图。综上所述，身管的振动一直是研究的重点，但是目前多数研究是针对身管自身的研究，未能够对身管和筋筒结构的耦合系统进行分析。本文利用多体系统传递矩阵法[5]，针对链式炮的特点，建立了链式炮系统的总传递方程和链式炮发射动力学仿真系统，得到了振动特性和起始扰动以及火炮的动力响应等的仿真结果。为火炮的结构设计和优化工作提供参考依据，对于提高火炮射击精度具有应用价值。
1 链式炮理论模型
1.1 理论模型
应用多体系统传递矩阵法，根据身管后坐式链式炮拓扑模型示意图1中彼此间有相对运动部件的几何结构、装配关系，以及各个部件间的联接方式。将三维实体模型转化为由体和铰组成的多体系统动力学模型，根据各个部件的自然属性，将它们分别视为刚体、弹性体、弹性元件等力学元件，并依次编号。
图1 身管后座式链式炮拓扑模型示意图
图2为基于传递矩阵法建立的模拟链式炮动力学模型。其中：安装体(0)、炮箱(2)、身管延伸部(8)和炮口制退器(14)简化为刚体；身管(9～13)、筋桶结构(3～5)简化为弹性体；安装体接触(1)、缓冲器(6)、支点(7)简化为弹性元件。
图2 模拟链式炮动力学模型示意图
根据动力学模型及各元件的编号，定义边界点和联接点的状态向量为：
Z2，1、Z2，3、Z3，4、Z4，5、Z5，7、Z2，6、Z8，6、Z8，9、Z9，10、Z10，11、Z11，12、Z12，13、Z13，14、Z14，0的形式均与Z0，1相同。
1.2 单元传递矩阵
1) 安装体到炮箱的传递关系为：
Z2，1=U1Z0，1
(1)
式中：
Kx、Ky和Kz为线刚度；Kx′、Ky′和Kz′为角刚度; I为单位矩阵; O为零矩阵。
2) 炮箱的传递关系为：
Z2，1=U2ZO，2
(2)
式中：
O1为缓冲器输出点；O2为筋筒结构输出点; m为刚体质量；ω为固有频率； J1为相对输入点I的惯量矩阵；
(b1，l，b2，l，b3，l)和(cc1，cc2，cc3)分别为输出点Ol和质心C在以输入点I为原点的坐标系中的坐标。
3) 筋桶结构的传递关系为：
Z3，4=U3Z2，3
(3)
Z4，5=U4Z3，4
(4)
Z5，7=U5Z4，5
(5)
式中，U3为空间振动Euler-Bernoulli梁的传递矩阵，具体形式可以查阅文献[5]得到，U4、U5的形式和U3相同。
4) 炮箱到身管延伸部的传递关系为：
Z8，6=U6Z2，6
(6)
式中，U6的形式和U1相同。
5) 身管延伸部的传递关系为：
Z8，9=U8Z8，6
(7)
式中：
J1为相对输入点I的惯量矩阵；
(b1，b2，b3)和(cc1，cc2，cc3)分别为输出点O和质心C在以输入点I为原点的坐标系中的坐标。
6) 身管的传递关系为：
Z9，10=U9Z8，9
(8)
Z11，10=U10Z9，10+E5U7Z5，7
(9)
Z12，11=U11Z11，10
(10)
Z13，12=U12Z12，11
(11)
Z14，13=U13Z13，12
(12)
式中，U9、U10、U11、U12和U13的形式和U3相同；U7的形式和U1相同；
7) 炮口制退器的传递关系为：
Z14，0=U14Z13，14
(13)
式中，U14的形式和U8相同。
根据式(1)～式(13)，得到链式炮的总传递方程。
(14)
式中：
E6=U14U13U12U11
E7=U10U9U8U6E3+E5U7U5U4U3E4
安装体(0)固定，位移和角位移恒等于零；炮口制退器(14)为自由边界，力和力矩恒等于零。将上述边界条件代入总传递矩阵，整理得：
(15)
式中，
ua，b为矩阵Uall中第a行第b列的元素。
为使式存在非零解，必须满足以下条件：
(16)
求解式(16)，得到动力学模型的固有频率ωk(k=1，2，…)。再将ωk代入式(15)，得到与ωk相对应的特征向量
进而由传递关系计算各联接点的状态向量，即得到与固有频率ωk相对应的振型。
2 发射动力学
实际过程中阻尼对系统的振动特性影响很小，所以忽略阻尼的影响，固模拟链式炮系统的动力学方程为：
Mvtt+Kv=f
(17)
式中，M为质量参数矩阵; K为刚度参数矩阵； v为位移列阵； f为外力列阵，i=(2，3，4，5，8，9，10，11，12，13，14)。
对于炮箱、身管延伸部和炮口制退器i=(2，8，14)，有：
Fi为外力主矢，mi为外力矩主，
为外力对第一个输入点矩的坐标。
对于弹性体i=(3，4，5，9，10，11，12，13)，有：
令：
(18)
式中， V为特征向量，
且满足正交性：
MVk·Vp=δk，pMp
(19)
(20)
式中，
将式(19)和式(20)代入式(17)，得：
(21)
针对每一个固有频率ωk(k=1，2，3,…)，求解式(21)得：
(22)
式中， fp(t)=f·Vp。
再将式(22)代入式(18)得到系统的动态响应。
(23)
3 链式炮发射动力学仿真
通过内弹道计算得到某链式炮的弹底压力，内弹道计算结果由程序计算所得。根据内弹道参数计算出流场初始化条件，利用Fluent软件对炮口流场进行非定常仿真计算得到炮口制退器产生制退力。将弹底压力和炮口制退力分别加载得到如图3的系统炮膛合力。
图3 系统炮膛合力曲线
依据链式炮元件和系统传递方程和传递矩阵，编制了链式炮振动特性仿真系统，并对链式炮振动特性进行了数值仿真。仿真固有频率如表1。
表1 仿真固有频率
模态阶数123456固有频率ω/(rad·s-1)3.77.810.513.116.223.1
建立了基于多体系统传递矩阵法的链式炮发射动力学数值仿真系统。数值仿真了某链式炮在射角和方向角均为0°时发射动力学过程，得到了起始扰动以及火炮的动力响应。并与相关的实验进行了对比分析，验证了仿真结果的正确性。图4为炮口扰动计算结果。
图4 炮口扰动曲线
根据设计要求，设置初始条件为：立靶距离200 m，弹丸初速为800 m/s，链式炮的射速为120发/min。将每次火炮发射时炮口扰动的位移值代入外弹道计算软件中，可以得到图5的10连发立靶坐标值。
图5 仿真立靶密集度效果图
图6为实际立靶密集度效果图。
图6 实际立靶密集度效果图
立靶密集度：水平EZ=0.15，竖直EY=0.19。
立靶密集度：水平EZ=0.19，竖直EY=0.18。仿真立靶密集度与试验立靶密集度结果相比差距小，而且试验链式炮的极限后座位移40 mm与系统仿真结果相符，表明了基于多体系统传递矩阵法的发射动力学仿真系统的正确性，为提高链式炮射击精度等动态性能的动态设计提供了重要的理论依据。
当前筋桶结构的长度为855 mm，利用链式炮发射动力学仿真系统做了筋桶长度与模拟仿真密集度的仿真结果图。由图7可以看出，在一定范围内增加筋桶结构的长度会改善炮口扰动，这是由于筋桶结构为身管提供了一个支点，使得炮口振动减少，有利于提高立靶密集度；但过于长的筋筒结构又由于自身的振动反过来影响身管的振动，使得炮口振动增大，不利于提高立靶密集度。
图7 筋桶长度与立靶密集度关系曲线
4 结论
利用多体系统传递矩阵法，针对链式炮的特点，建立了链式炮系统的总传递方程和链式炮发射动力学仿真系统，数值仿真了某链式炮在射角和方向角均为0°时发射动力学过程，得到了起始扰动以及火炮的动力响应。得出结论：将筋桶结构的长度增加6 mm可以进一步提高立靶密集度。
[1] 张相炎，郑建国，袁仁枢.火炮设计理论[M].北京：北京理工大学出版社，2014.
[2] 李猛.某火炮动态特性与射击频率的匹配性分析[J].兵器装备工程学报，2017，38(01):44-46..
[3] 史跃东，罗忠，钱超.运动激励下身管动力学特性研究[J].四川兵工学报，2013，34(04):1-3.
[4] 黄晨烨，戴劲松，尹强，等.炮箍个数对身管组模态频率及振型的影响分析[J].兵工自动化，2020，39(05):49-53.
[5] 芮筱亭，贠来峰，陆毓琪，等.多体系统传递矩阵法及其应用[M].北京：科学出版社，2008.
[6] 谢锋，魏立新，朱延飞，等.转管火炮炮口扰动对射击密集度的影响分析[J].舰船科学技术，2018，40(21):153-155.
[7] LI Bing，QUAN Shuanglu，JIN Wei，et al.Identification of Clearance and Contact Stiffness in a Simplified Barrel-Cradle Structure of Artillery System[J].2015，7(2):1490-1503.
[8] WANG Pingxin，RUI Xiaoting，YU Hailong，et al[J]Transverse vibration of the upper track of a tracked vehicle with tubular busing in the track pin.2019，140:504-519.
[9] 陈东阳，肖清，谢俊超，等.基于多体系统传递矩阵法的舵系统振动特性分析[J].哈尔滨工程大学学报，2020，41(02):184-193..