带扰流片旋转稳定弹动态稳定性

杨杰， 常思江， 魏伟

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

0 引言

通过在常规无控弹上加装弹道修正组件，弹道修正弹逐步发展成为一类新型弹药。弹道修正弹介于常规无控弹与导弹之间，其成本低、毁伤率高、效费比高，因而发展较为迅速广泛，其中在针对旋转弹的改造中，双旋弹道修正弹极具发展潜力[1-2]。

根据改造方式的不同，双旋弹道修正弹主要分为两类：一类是在旋转弹头部加装弹道修正引信，如精确制导组件(PGK)；另一类是在旋转弹尾部加装扰流片。安装在头部的修正执行机构为鸭舵，包括固定舵和可动舵，固定舵舵偏不可调，使得旋转弹修正能力有限[3-4]，可动舵舵偏可调但使得舵机结构复杂[5-6]。与鸭舵相比，在转向能力和阻力损失方面，安装在尾部的扰流片与其性能相当，但扰流片结构简单，修正能力更强，且攻角响应方向与扰流片提供的升力方向一致，具有显著优势[7]。

带扰流片旋转稳定弹的结构如图1所示，其由前体(即弹体)和后体两部分组成，扰流片安装在旋转弹后体船尾部分，前体与后体通过轴承连接。旋转弹未发射时，扰流片收缩在后体控制舱内，全弹外表与普通旋转弹相同。旋转弹发射后，扰流片暂不弹出，飞行过程中后体通过电机反转，与前体形成差动旋转，而前体保持较高转速以保证全弹陀螺稳定性。当旋转弹飞行到预定区域后，后体通过电机反转至稳定状态后相对空间不转，扰流片展开，提供所需的控制力和控制力矩，进行弹道修正。

图1 带扰流片旋转稳定弹结构示意图Fig.1 Structure diagram of spin-stabilized projectile with spoiler

带扰流片旋转稳定弹的气动分析计算表明[8]，在超声速条件下，扰流片的存在使旋转弹尾部局部产生剧烈的压力变化，扰流片迎风表面压力最大且靠近扰流片前方的区域压力也较大，扰流产生的马赫波系在弹体表面呈类似半椭圆形分布，随着马赫波系距离扰流片越来越远，相应位置的弹体表面压力逐渐降低，在亚、跨声速条件下，扰动向前传播的距离更远。在任何马赫数下，在扰流片两侧的压力分布对称，在旋转弹同一截面沿周向的压力分布在扰流片处最大，越靠近扰流片压力越大。由此可见，扰流片所造成的压力分布不对称产生了相应的控制力和控制力矩，进而影响带扰流片旋转稳定弹的运动特性。

目前国内外针对带扰流片旋转稳定弹的研究较少，Rogers[7]通过研究发现，扰流片在亚音速范围内具有一定的控制效率，随着马赫数的增加，控制效率呈指数级增加，该机构对高超声速旋转弹的控制研究具有重要意义。美国陆军研究实验室的Celmins等[9]基于40 mm旋转稳定弹提出几种带扰流片旋转稳定弹的结构设计方案，通过实验验证了该结构可靠的性能。随后Fresconi等[10]通过实验和数值计算得到空气动力数据，研究了扰流片的结构参数与修正能力之间的关系。Wey等[11]在旋转弹船尾安装滚转解耦的扰流片，通过空气动力分析和七自由度弹道仿真说明了该方案的作用原理和可靠的性能，研究表明，在不同射击条件下，155 mm旋转弹的弹道散布都可被修正，且扰流片只需在弹道末端作用。Fresconi等[12]将扰流片应用于旋转稳定弹的控制系统，以155 mm旋转弹为例，通过6自由度模型研究制导律并计算分析其修正能力，数值计算研究了旋转弹的角运动特性，仿真结果表明该结构具有较大的修正能力，能够有效减小散布且未出现飞行不稳定现象。除此之外，其他将扰流片应用于尾翼弹的研究[13-18]同样表明，扰流片具有较好的修正能力和可靠的性能。

在上述对带扰流片旋转稳定弹的研究中，多集中于对其修正原理及气动特性的分析，对旋转弹弹道特性的研究结果往往局限于数值仿真，缺乏对动态稳定性的深入研究，包括线性和非线性角运动特性等。为深入探索该类弹道修正弹的动态稳定性，本文建立了带扰流片旋转稳定弹的动力学模型，推导出带扰流片旋转稳定弹的线性角运动方程，并针对该类旋转弹的控制原理，研究了弹体攻角在扰流片作用下的瞬态响应和稳态响应，随后考虑了扰流片气动系数的非线性，分析了旋转弹结构参数、飞行参数和气动参数对旋转弹系统分岔特性及稳定性的影响，并对其修正能力进行了分析研究，以期为今后研制该类旋转弹提供相应的理论依据。

1 动力学建模

如图1所示，由于扰流片的质量和体积都很小，可忽略质量不对称并假设控制过程无延迟。扰流片展开后，旋转弹整体产生气动不对称，形成附加的阻力、升力和俯仰力矩。

动力学建模用到的坐标系有：地面坐标系、弹道坐标系、弹轴坐标系、前体坐标系和后体坐标系。前体坐标系与后体坐标系的定义如下：

前体坐标系Oxfyfzf与常规弹体坐标系Oxbybzb重合[19]，而后体坐标系Oxayaza则是在前体坐标系基础上将Oyf轴和Ozf轴绕Oxf轴转过角度Δγ得到，坐标系示意图如图2所示，Δγ为差动滚转角，有

图2 坐标系示意图Fig.2 Coordinate system

Δγ=γa-γf，

(1)

式中：γa、γf分别为后体和前体的滚转角，γf与弹体滚转角γ相同。其余坐标系的定义可参见文献[19]。

图3为零攻角时作用于后体的附加力和力矩示意图，设沿弹轴从弹尾向前看有一垂直于弹轴的横截面，假设扰流片位于后体坐标系Oya轴位置，ΔRx表示附加阻力，ΔRy表示附加升力，ΔMz表示附加力矩。定义Oya轴与弹轴坐标系Oη轴的夹角为后体扰流片的滚转角γs，从η轴算起，顺时针为正。

图3 零攻角时附加力和力矩示意图Fig.3 Schematic diagram of additional force and moment for AOA=0°

基于上述坐标系，可推导得到带扰流片旋转稳定弹的6自由度动力学模型，由于本文主要对角运动特性进行研究，因此不再对质心运动方程进行赘述，仅给出带扰流片旋转稳定弹的角运动方程：

(2)

式中：θe表示速度高低角；t表示时间；Fyt表示作用在旋转弹上的空气动力在弹道坐标系Oyt轴下的分量；v表示弹丸质心速度；m表示弹丸质量；ψa表示速度方向角，

(3)

Fzt表示作用在旋转弹上的空气动力在弹道坐标系Ozt轴下的分量。

(4)

式中：ωξ、ωη、ωζ表示弹丸绕质心转动的角速度；C表示极转动惯量；A表示赤道转动惯量；Mη表示作用在旋转弹上的空气动力矩在弹轴坐标系Oη轴上的分量[20]；φa表示弹轴方位角，

(5)

(6)

式中：Mζ表示作用在旋转弹上的空气动力矩在弹轴坐标系Oζ轴上的分量。

(7)

式中：φe表示弹轴高低角。

为研究带扰流片旋转稳定弹的动态稳定性，下面根据以上建立的动力学模型，基于线性假设对攻角的解析解进行推导分析，并结合非线性理论对其动力学分岔特性进行分析，最后进行修正能力分析。

2 攻角运动特性

考虑到带扰流片旋转稳定弹的角运动是在理想弹道附近进行的，在建立角运动方程时略去小量及小量乘积，有θ≈θe，则(2)式～(5)式可简化为

(8)

式中：ψc表示速度方向角分量；δe表示高低攻角；δa表示方向攻角；g表示重力加速度。

(9)

(10)

式中：φc表示弹丸姿态角分量。

(11)

根据复摆动角φ、复偏角ψ和复攻角Δ的定义及相互关系，将(8)式～(11)式写成复数形式，可得

(12)

(13)

联立(12)式和(13)式，忽略小量乘积，并将自变量从时间改为弧长s，得到攻角运动方程

(14)

考虑由扰流片引起的角运动，线性角运动方程(14)式可写为

(15)

下面对扰流片弹出后弹体攻角的瞬态响应和稳态响应进行分析。

2.1 扰流片弹出产生的瞬态解

旋转弹在飞行过程中，扰流片弹出对旋转弹进行控制，相当于旋转弹受到一个阶跃激励。当旋转弹受到阶跃激励时，角运动方程可写成如下形式：

(16)

式中：ε(s)为单位阶跃函数；si为阶跃激励开始作用时的弹道弧长。

在[0,si]弧长间隔内，旋转弹未受到激励作用，处于无控状态；在s>si弧长段，激励一直存在。根据杜哈梅积分[21]，可以将该激励视为一系列脉冲激励的迭加，即

(17)

式中：x(t)为系统响应；p(τ)为外部激励；h(t-τ)为τ时刻由单位脉冲引起的系统响应，有如下表达式，

(18)

于是得到阶跃激励时旋转弹角运动方程的解为

(19)

(19)式得到的攻角瞬态解与6自由度弹道仿真得到的攻角进行对比，图4为60 s启控、扰流片滚转角为90°时的对比曲线。

图4 60 s启控、90°滚转角时攻角对比曲线Fig.4 Comparison curves of angle of attack for t=60 s and γs=90°

由图4可知，由(19)式计算出的攻角瞬态响应与刚体弹道数值仿真得到的攻角是吻合的，证明了该解析解的正确性。结合kc和bc的表达式与(19)式分析可知，扰流片附加升力系数导数和附加力矩系数导数越大，攻角响应幅值越大，而计算流体力学模拟计算表明该系数与扰流片高度近似呈正比关系，因此扰流片高度越高，攻角响应幅值越大。

2.2 扰流片弹出产生的稳态解

Δ″+(H-iP)Δ′-(M+iPT)Δ=(kc+ibcP)eiγs，

(20)

式中：等号右边为定值，利用常数变易法求得特解为

(21)

(21)式的分母实数化得到

(22)

式中：

(23)

实部与虚部分离，可得

(24)

(25)

(24)式和(25)式即为扰流片弹出后产生的攻角稳态解。60 s时启控，扰流片在0°方位上产生阶跃激励时，数值积分和理论推导计算出的攻角响应的瞬态解和由(24)式和(25)式计算得到的稳态解如图5所示。

图5 稳态解与瞬态解对比曲线Fig.5 Curves of comparison between steady state solution and transient solution

由图5可知，稳态解反映了旋转弹攻角在阶跃激励下摆动结束后的平衡位置。相比之下，解析计算所得瞬态解的平衡位置与稳态解更接近，这是因为在理论推导过程中，稳态解和瞬态解的计算采取了相同的假设。

弹体攻角的稳态响应为一复攻角，由(24)式和(25)式可求得其幅角为

(26)

(23)式代入(26)式，并根据反正切相加减定理[22]，可得

(27)

由(27)式可知，Δst的幅角始终与扰流片滚转角γs相差一个大约为180°的相位角。

3 非线性分岔特性

随着现代外弹道理论的发展，与弹箭非线性运动相关的研究逐步开展[23]，而分岔理论是非线性分析的一个有力工具。分岔理论的应用有助于寻找引起弹箭飞行特性突变的弹箭气动参数和结构参数分岔点，为参数设计、避免非线性飞行突变不稳提供依据。为研究分析带扰流片旋转稳定弹的结构参数和气动力参数对弹箭角运动特性和飞行稳定性的影响，下面应用分岔理论对该旋转弹进行非线性分岔特性分析。

3.1 研究方法

分岔问题起源于力学失稳现象的研究，可以分为静态分岔和动态分岔。若任意小的参数变化会使结构不稳定系统的相轨迹拓扑结构发生突然变化，则称这种变化为分岔。因此可将分岔的定义叙述为：

对于含参数的系统：

(28)

式中：x为状态变量，x∈Rp；μ为分岔参数，μ∈Rq；p、q表示维数，Rp、Rq表示欧几里德空间。当参数μ连续变动时，若系统(28)式相轨迹的拓扑结构在μ=μ0处发生突然变化，则称系统(28)式在μ=μ0处出现分岔，μ0称为分岔值或临界值，(x,μ0)称为分岔点。在参量μ的空间Rq中，由分岔值构成的集合称为分岔集。在(x,μ)的空间Rp×Rq中，平衡点和极限环随参数μ变化的图形称为分岔图。

鞍结分岔是指控制参数变化过程中系统因形成鞍结点而出现的分岔，该分岔点也被称为极限点，属于静态分岔；霍普夫分岔是指系统参数变化经过临界值时，平衡点由稳定变为不稳定并从中生长出极限环，是一类动态分岔。

由于带扰流片旋转稳定弹结构非对称、气动非对称，其动力学系统属于高维非线性系统，求得其解析解比较困难，因此本文采用数值方法对其分岔特性进行研究分析。数值方法首先计算得到系统解随分岔参数变化的曲线即平衡点曲线，通过每一点的特征值判断其稳定性、分岔类型及分岔位置，最终得到系统分岔图。

3.2 非线性角运动方程

带扰流片旋转稳定弹的动力学系统是一个高维的非线性系统，状态变量较多，因此需要进行简化将其降维。而旋转弹的稳定性与其飞行过程中的角运动密切相关，因此根据动力学中的角运动方程构建低维模型。

带扰流片旋转稳定弹在飞行过程中，近似认为δe=φe-θe，δa=φa-ψa，根据(2)式～(7)式，简化得到其非线性角运动方程：

(29)

式中各力和力矩表达式见文献[20]，其中：c′y=cy,l+cy,n·δ2，δ表示总攻角，cy,l和cy,n分别为升力系数导数的线性项和立方项；m′z=mz,l+mz,n·δ2，mz,l和mz,n分别为静力矩系数导数的线性项和立方项。扰流片外露高度Hs是对扰流片气动特性影响最大的参数[7]，为便于研究结构参数Hs对旋转弹非线性动力学系统的影响，建立扰流片气动参数模型如下：

(30)

式中：Δcy,l和Δcy,n分别为单位高度扰流片引起的附加升力系数导数的线性项和立方项；Δmz,l和Δmz,n分别为单位高度扰流片引起的附加静力矩系数导数的线性项和立方项，该数值通过计算流体力学方法数值模拟得到。

基于简化得到的低维非线性模型(29)式，对其进行分岔特性分析与仿真，研究带扰流片旋转稳定弹在不同参数下的分岔特性，仿真时旋转弹速度取400 m/s，转速取1 300 rad/s，其他必需的旋转弹参数和气动参数见表1和表2，表1中：m″y表示马格努斯力矩系数2阶联合偏导数；m′zz表示赤道阻尼力矩系数导数。

表1 主要计算参数Tab.1 Main calculation parameters

表2 扰流片气动参数Tab.2 Aerodynamic parameters of spoiler

3.3 分岔特性分析

旋转弹处于有控状态时，扰流片弹出且转角固定，扰流片产生持续的附加力矩对弹丸进行控制。当扰流片方位角为0°时，计算得到(29)式的系统平衡点为(-0.16 rad,0.006 rad,-0.009 rad/s,-0.04 rad/s)T，分别与δe、δa、ωη、ωζ对应。

为研究速度、气动参数等因素对系统特性的影响，首先计算系统的分岔图，图6为攻角随扰流片高度变化的分岔曲线。考虑到研究对象实际的物理意义，扰流片高度不可能为负值，因此下文仅针对Hs>0 mm的情况进行分析。

图6 系统分岔图Fig.6 Bifurcation diagram of system

由图6可知，系统的平衡点数量和稳定性随着扰流片高度的变化而变化。当0 mm≤Hs≤3 mm时，系统有一个稳定的平衡点；当3 mm≤Hs≤12.5 mm时，系统有两个稳定平衡点和一个不稳定平衡点。从平衡点的分布可以看出：稳定平衡点1为期望的飞行状态，带扰流片旋转稳定弹具有良好的动态特性，旋转弹随着Hs的逐渐增大产生更大的攻角；而稳定平衡点2在较小的Hs下产生较大的攻角，不利于旋转弹飞行。对于本文算例而言，为保证带扰流片旋转稳定弹按照期望的平衡状态飞行，应使扰流片外露高度在0～3 mm的范围内。

3.4 不同因素对系统特性的影响

带扰流片旋转稳定弹的动力学系统分岔图可能随其飞行参数、气动参数等表现出不同的特性，分析这些参数的影响可以为旋转弹气动外形和弹体结构的设计提供参考依据。由于稳定平衡点1是期望的平衡点，因此下文仅对稳定平衡点1的变化曲线进行分析。

3.4.1 速度变化对系统特性的影响

旋转弹在不同条件下的初速以及飞行中的速度变化都可能存在差异，因此首先分析速度不同对系统特性的影响。图7为速度不同时攻角随扰流片高度变化的平衡点曲线。

图7 不同速度下系统平衡点曲线Fig.7 Equilibrium curves of system at different velocities

由图7可知：带扰流片旋转稳定弹速度的变化对系统的稳定域影响不大，飞行稳定情况下的最大扰流片外露高度变化不大；在γs=0°条件下，Hs不变时，速度越大，方向攻角越小，而高低攻角变化不大。

3.4.2 扰流片方位角变化对系统特性的影响

旋转弹在飞行过程中，通过不断改变扰流片方位角从而对旋转弹修正方向进行控制，因此需要分析不同方位角对系统特性的影响。图8为扰流片方位角不同时计算得到的系统平衡点曲线。

图8 不同扰流片方位角系统平衡点曲线Fig.8 Equilibrium curves of system at different γs

由图8可知：当γs=0°时，旋转弹产生与扰流片方位相反的高低攻角和近似为0°的方向攻角；当γs=90°时，旋转弹产生与扰流片方位相反的方向攻角和近似为0°的高低攻角。随着γs的变化，高低攻角和方向攻角随γs呈正弦变化。扰流片外露高度Hs越大，产生的攻角幅值越大。

3.4.3 气动系数变化对系统特性的影响

带扰流片旋转稳定弹在飞行过程中，角运动主要受到气动力矩的影响，而扰流片主要通过产生的附加力矩对旋转弹进行控制，因此各力矩系数可能会对系统特性产生影响。图9为弹体静力矩系数导数线性项不同时系统分岔图；图10为扰流片附加静力矩系数导数线性项不同时系统分岔图；图11为弹体静力矩系数导数立方项不同时系统分岔图；图12为扰流片附加静力矩系数导数立方项不同时系统分岔图。

图9 不同mz,l系统平衡点曲线Fig.9 Equilibrium curves of system at different mz,l

图10 不同Δmz,l下系统平衡点曲线Fig.10 Equilibrium curves of system at different Δmz,l

图11 不同mz,n下系统平衡点曲线Fig.11 Equilibrium curves of system at different mz,n

图12 不同Δmz,n下系统平衡点曲线Fig.12 Equilibrium curves of system at different Δmz,n

由图9可知：随着静力矩系数导数线性项的变化，系统能产生的最大攻角变化不明显；mz,l增大会使系统的稳定域增大，在Hs相同条件下，mz,l增大会使攻角的幅值减小。

由图10可知：随着附加力矩系数导数线性项的变化，系统能产生的最大攻角变化不明显；Δmz,l增大会使系统的稳定域减小，在Hs相同条件下，Δmz,l增大会使攻角的幅值增大。

由图11可知，随着静力矩系数导数立方项绝对值的增大，系统能产生的最大攻角幅值逐渐减小，且系统的稳定域逐渐减小。在Hs相同条件下，当Hs较小时，mz,n不同产生的攻角差异不大；当Hs较大时，mz,n绝对值越大，产生的攻角幅值越大。

由图12可知，随着附加静力矩系数导数立方项绝对值的增大，系统能产生的最大攻角幅值逐渐增大，且系统的稳定域逐渐增大。在Hs相同条件下，当Hs较小时，Δmz,n不同产生的攻角差异不大；当Hs较大时，Δmz,n绝对值越大，产生的攻角幅值越小。

4 修正能力分析

对带扰流片旋转稳定弹的攻角运动和动力学分岔特性进行研究之后，下面根据其动力学模型对其修正能力进行数值仿真分析。图13为旋转弹在弹道末端距离地面不同高度hc启控时的修正范围对比曲线；图14为旋转弹在不同扰流片高度Hs下的修正范围对比曲线。

图13 不同hc下修正范围对比曲线Fig.13 Comparison curves of correction range under different hc

图14 不同Hs下修正范围对比曲线Fig.14 Comparison curves of correction range under different Hs

由图13可知，距离地面高度hc=500 m时启控修正半径约为500 m，hc=1 000 m时启控修正半径达1 000 m，hc=1 500 m时启控修正半径达1 500 m，因此距离地面高度越高时启控，旋转弹修正范围越广，修正效率较高。

由图14可知，Hs=7.75 mm时修正半径约为1 000 m，Hs=11.5 mm时修正半径达1 600 m，由此可见扰流片高度越高，所能提供的修正能力越强。

5 结论

本文建立了带扰流片旋转稳定弹的动力学模型。针对带扰流片旋转稳定弹的工作模式，求解了旋转弹攻角的瞬态响应和稳态响应。通过分岔计算分析了不同参数对系统特性的影响，得到如下结论：

1) 扰流片产生阶跃激励控制信号时，理论推导得到的攻角瞬态响应解具有较高的精度，与数值解较吻合，稳态解则反映了攻角在阶跃激励下摆动结束的平衡位置，且扰流片高度越高，攻角响应的幅值越大,攻角响应方向与扰流片转角存在一个接近180°的相位差。

2) 有控状态下，带扰流片旋转稳定弹的动力学系统存在两个稳定的平衡点，其中只有一个平衡点为期望的平衡点，并反映了期望的扰流片高度范围。因此应选择合适的扰流片高度，既能产生所需攻角又不致弹体失稳。当扰流片转角固定时，速度对系统稳定域影响不大；扰流片转角变化对系统稳定域影响不大。

3) 当扰流片转角固定时，随着静力矩系数导数线性项增大，期望平衡点的稳定域逐渐增大，在Hs相同条件下所能产生最大攻角幅值减小；随着附加静力矩系数导数线性项增大，期望平衡点的稳定域逐渐减小，在Hs相同条件下所能产生最大攻角幅值增大。

4) 当扰流片转角固定时，随着静力矩系数导数立方项增大，期望平衡点的稳定域逐渐减小，系统能产生的最大攻角幅值逐渐增大；随着附加静力矩系数导数立方项增大，期望平衡点的稳定域逐渐增大，系统能产生的最大攻角幅值逐渐增大。