0 引　言

在研究非线性系统特性时，经常会出现状态空间的一个感兴趣的区域中出现多个吸引子共存，称之为多稳态现象^[1-2]。多稳态现象的危害在实际工程中，高速飞行器的颤振问题、高速列车的蛇形运动失稳、海洋结构中流固耦合失稳等问题，是由于系统对初始条件敏感的原因，即多稳态现象的存在，比如双极限环^[3]、平衡点和极限环共存^[4-5]。当多稳态现象在电力系统中出现时，将导致系统运行参数持续无规则的振荡，严重危害系统的运行安全。在非线性系统中，随机因素广泛存在，比如机翼遇到的高空湍流、列车遇到的来自轨道的激励等^[6]。初始条件使得系统的最终去向很不确定，往往具有很复杂的行为，导致系统的运行状态存在不可预料因素，很可能与设计者的最初预想完全不一致，对工程应用也将造成巨大的威胁，因此需要透彻地了解此类非正常现象。

超空泡是一种奇特的物理现象，空泡将整个水下航行体包裹起来，最大限度地避免水的黏性阻力，实现高速航行^[7-8]。通常用空化数来表征相似的空泡状态，它也是空泡流中最重要的相似参数，空化数的定义为^[3]：σ=2(p_∞−p_c)/ρV²，P_∞是无扰动流体中的压力，P_c为空泡内压，ρ为水的密度，V为航行体在理想流体中的运动总速度。然而，超空泡航行体在实测中往往受初始运动参数(初始深度、初始垂直速度、初始俯仰角、初始俯仰角速度)的影响而运动失稳，也具有初始条件敏感性。而对于多稳态现象的研究则主要存在于洛伦兹系统以及电力电子等非线性系统中，用以探究系统容易受内在参数变化影响的原因^[9-10]，有利于促进其效率变换特性的改善，对于进一步的稳定设计和混沌利用具有重要的指导意义。超空泡航行体也是一种强非线性系统，其中存在多种非线性现象，为了设计稳定可靠的超空泡航行体控制系统，需要对其运动状态的不可预料因素进行详细研究，从而避免产生一些非理想的破坏性结果。因此完全有必要对其多稳态现象进行分析和讨论。

1 超空泡航行体动力学模型

在推导动力学方程组之前为简化模型给出几个工程上常用的假设条件：

1）航行体为刚体，其外形关于自身的纵平面对称；

2）不考虑波浪和洋流对航行体运动的影响，认为运动处于理想流体中；

3）航行体运动过程中，其质量不受燃料消耗的影响，认为它的质量恒定；

4）航行体是水中由发射架水平发射；

5）由于整个航行体几乎都处于超空泡内，可认为其重心与浮心重合。

如图1所示，长度为L的航行体由两段构成：后端长度2/3L和半径R的圆柱段，以及前端长度1/3L的圆锥段。重心距离头部为17/28L，作用在航行体上的力主要有：空化器上的升力F_cavitator、尾翼上的升力F_fins、航行体质心位置的重力F_gravity，以及滑行力F_planing^[8]，V为纵平面内航行体头部空化器的合速度，w为垂直速度，方向垂直于航行体的中心轴线指向下，θ是航行体的俯仰角，q为俯仰角速度，航行体深度为z。

当超空泡航行体在无控状态下运动时，航行体的重力迫使它在空泡内滑行是不稳定的，为保证航行体的稳定运动，设置了反馈控制器，其控制输入为尾翼偏转角δ_e和空化器偏转角δ_c。Dzielski等^[8]提出了经典控制律δ_e=0，δ_c=15z^−30θ−0.3q，但由于尾翼偏转角为0，往往会使航行体缺乏尾翼提供的支持力，当空化器提供的升力不能平衡航行体的重力时，将导致航行体在重力的作用下浸入空泡而失稳。因此，本文令控制律为δ_e1=kq，δ_c1=15z−30θ−0.3q，k为航行体俯仰角速度q的反馈控制增益。为了研究超空泡航行体的动力学特性，以航行体头部圆盘形空化器顶端面的圆心为体坐标系原点，根据表1，通过分析计算航行体所受的流体动力^[8]，得到以z，w，θ，q为状态变量，以k和σ为可变参数的四维动力学模型为：

$\begin{split}{l}\dot{z}=&w-V\theta \text{，}\\ \dot{w}=&-15{b}_{22}z+{a}_{22}w+30{b}_{22}\theta +(0.3{b}_{22}+{a}_{24}+k{b}_{21})q+\\ &G+{d}_{2}{F}_{planing}，\\ \dot{\theta }=&q\text{，}\\ \dot{q}=&-15{b}_{42}z+{a}_{42}w+30{b}_{42}\theta +(0.3{b}_{42}+{a}_{42}+k{b}_{41})q+\\ &{d}_{4}{F}_{planing}。\end{split}$

(1)

式中：C=635.6588(1+σ)R_n², a₂₂=−4.2079CV+0.4250, a₂₄=−6.1540CV+14.7101V, a₄₂=1.5154CV, a₄₄=−1.3046CV,b₂₁=0.6167CV², b₂₂= −2.8054 CV², b₄₁=−0.7249CV², b₄₂=2.2406CV², d₂=−1.2266,d₄=1.4492。

超空泡航行体模型参数值如表1所示^[8]。

2 动力学行为分布图

基于动力学模型(1)，随机选取其状态变量的z，w，θ，q的初始值，依照Lyapunov稳定性理论^[11]将模型的稳定解、周期解、混沌解分别用R，G，Y在图2中表示出来，绘制出分岔参数空化数σ和尾翼偏转角反馈控制增益k所定义的动力学地图，刻画了系统动力学行为对σ和k的依赖关系。图中，当σ，k的取值对应于R区域的点时，航行体能够稳定运动；在G区域中任选一点(σ, k)，航行体周期振荡；当σ，k在Y区域内取值时，航行体则会发生剧烈的振动与冲击，进而倾覆；B区域表示系统发散。

该图相对完整地反映了当参数σ，k同时变化时，系统在动力学行为上所处的不同状态。水平切面是系统随σ变化的分岔图，竖直切面是系统随k变化的分岔图。利用动力学行为分布图能够确定航行体稳定运动的参数取值范围。当空化数σ一定时，在稳定运动和周期运动的范围内调节尾翼偏转角反馈控制增益k的取值，能够有效实现超空泡航行体的稳定航行，对航行体的稳定控制具有指导意义。观察图2可以发现：

1）航行体具有R区域稳定、G区域周期、Y区域混沌3种稳态运动，分别在3个区域中任意各取一点，随机选取其状态变量的z，w，θ，q的初始值，其相轨在平面上的投影如图3所示。区域R中所取的参数为σ=0.02356，k=2.05，如图3(a)所示，垂直速度w和俯仰角θ在反馈控制律的作用下被吸引到平衡点上，此时航行体对应的运动状态如图3(b)所示，尾翼一直处于穿过空泡浸入水中的状态，航行体在空泡内的位置和姿态固定，由头部的空化器和尾部的滑行力维持系统的平衡，航行体处于小攻角斜向运动的稳定航行状态。在区域G中随机选取参数σ=0.02869，k=5.36，如图4(a)所示，系统映射形成了闭合的极限环，极限环的出现说明系统发生了周期性振荡。此时航行体对应的运动状态如图4(b)所示，航行体的垂直速度w沿极限环不断变化，航行体尾翼相应地不时摆动与空泡碰触，引起航行体滑行力周期性的变化，航行体周期振荡。在区域Y中随机选取参数σ=0.03261，k=0.58，图5中混沌吸引子的出现表明航行体的运动具有复杂的非线性动力学行为，航行体尾翼非周期地摆动，不断地与空泡碰触，穿过空泡浸入水中的浸没深度超过了尾翼的长度，导致航行体倾覆。

在实际工程应用中，航行体的稳定运动和周期振荡对应的都是稳定的弹道，而混沌振荡的出现回导致航行体的倾覆，必须进行有效的控制以避免这种情况发生。

2）σϵ [0.01980, 0.03255]，k_qϵ [8.18, 24.53]范围内，图中显示R区域的稳定点密集地散布在G区域的周期区域中。航行体发射时初始参数的微小变化会引起航行体运动状态的改变，表明该范围可能出现稳定平衡点与周期吸引子共存的多稳态现象。

在该范围中选取σ=0.02670，k=−9.7，代入式(1)，经计算可得此处的平衡点(z^*, w^*, θ^*, q^*)= (0.0186，0.6510，0.0081，0)。图6为不同初始运动参数下的相轨迹图。如图6(a)所示，当初始运动参数为α₁(z₀，w₀，θ₀，q₀)= (−0.1426，−1.8572，−0.2654，−1.3698)时，在反馈控制律的作用下系统状态变量收敛于R区域的稳定点；航行体处于小攻角斜向运动的稳定航行状态。如图6(b)所示，当初始运动参数为α₂ (z₀，w₀，θ₀，q₀)= (0.5377，1.8339，−2.2588，0.8622)时，系统状态变量收敛于B区域的周期吸引子。航行体尾部不时地摆动与空泡碰触，引起滑行力周期性的变化，航行体周期振荡运动。图6(c)是两者共存的情形在z-θ-w平面上的投影，R区域的点表示稳定平衡点，环1是周期吸引子。上述讨论表明，系统表现出对初值敏感的多吸引子共存现象，航行体被发射的初始深度，垂直速度，俯仰角以及俯仰角速度取值不同时，其轨线可能趋于2种类型的吸引子，即稳定平衡点与周期吸引子，2种吸引子表明航行体运动状态的不同。

3）σϵ [0.03258, 0.03563]，kϵ [−15.12, 26.75]范围内，图中显示G区域的周期区域中穿插着Y区域的混沌点。初始运动参数设置不当会导致航行体的倾覆，表明该范围可能出现周期吸引子与混沌吸引子共存的多稳态现象。

在该范围中选取σ=0.03458，k=−6.54，将其代入式(1)，经计算可得此处的平衡点(z^*, w^*, θ^*, q^*)= (0.0015，−0.0245，−0.0003，0)。图7为不同初始运动参数下的相轨迹图。如图7(a)所示，当初始运动参数为β₁(z₀，w₀，θ₀，q₀)= (−1.0891，0.0326，0.5525，1.1006)时，系统方程的解收敛于一个R区域的周期吸引子，航行体处于周期振荡的运动状态；如图7(b)所示，当初始运动参数为β₂ (z₀，w₀，θ₀，q₀)= (−1.2141，−1.1135，−0.0068，1.5326)时，系统的解收敛于B区域的混沌吸引子。混沌吸引子的出现说明航行体已运动失稳。图7(c)是两者共存的情形在z-θ-w平面上的投影，环1表示周期吸引子，环2是混沌吸引子。上述讨论表明，当超空泡航行体发射的初始条件不同时，可能会导致运动状态的不同，甚至可能引起航行体的倾覆，存在极限环与混沌吸引子共存的现象。

3 超空泡航行体的多稳态运动

为了进一步研究航行体的运动状态，验证多稳态现象。

在σ=0.02670，k=−9.7处，航行体运动特性的仿真结果如图8所示。图中实线表示初始运动参数α₁时，系统4个状态变量z, w, θ, q，随时间变化的运动情况，虚线则表示初始运动参数α₂时，系统各个变量的时域响应情况。

初始运动参数为α₁时，在控制律的作用下，状态变量z, w, θ, q迅速稳定在平衡点(0.0772, 3.0198, 0.0365, 0)上，航行体的尾部穿过空泡浸入水中，产生的滑行力与空化器和尾翼提供的升力共同平衡航行体重力，使航行体稳定运动。

初始运动参数为α₂时，4个状态变量均围绕平衡点(−0.0282, −1.2007, −0.0145, 0)周期振荡；航行体尾部不断与空泡壁面发生碰撞。

在σ=0.03458，k=−6.54处，航行体运动特性的仿真结果如图9所示。实线表示初始运动参数为β₁时，系统4个状态变量z, w, θ, q，随时间变化的运动情况；虚线为初始运动参数为β₂时系统各个变量的时域响应情况。

当初始运动参数为β₁时，4个状态变量z，w，θ，q分别以平衡点(0.0019，−0.0138，−0.0002，0)为中心周期振荡，在重力的作用下，航行体的尾部与空泡壁面发生碰撞，穿透空泡浸入水中，产生滑行力。随着浸没深度的增加，尾部在滑行力的作用下又被迅速弹回空泡内，继而与空泡另一壁面碰撞，如此往复，表明航行体在空泡内周期振动。

当初始运动参数为β₂时，航行体状态变量随时间的变化在平衡点附近发生非周期振荡，航行体尾部与空泡壁面发生轻微的碰撞，表明航行体此时在空泡内处于非周期的混沌状态。由于系统各个状态变量处于非周期性不稳定状态，航行体的运动也极容易受到外界干扰，从而导致失稳。

4 结　语

本文基于超空泡航行体的动力学模型，通过动力学行为分布图探讨了航行体动力学行为的一般特性，基于相轨迹图确认了多种吸引子共存，并利用时域仿真验证了超空泡航行体的多稳态现象，得到以下结论：

1）在超空泡航行体系统中，动力学行为分布图可以展示可能存在共存吸引子的参数区域，为发射稳定运动的超空泡航行体提供了参数设置的参考依据。

2）超空泡航行体具有稳定、周期、混沌3种稳态运动，在适当的参数条件下，存在多种多稳态运动，其中包括稳定平衡点与周期吸引子共存、周期吸引子与混沌吸引子共存的情况。当系统参数取定时，航行体在不同初始状态下可能表现出完全不同的运动状态。

3）在实际工程应用中，可根据动力学行为分布图避免选择混沌吸引子和发散区域对应的参数，即避免使尾翼偏转角控制增益k>32，或空化数σ>0.03258。