一类弹性碰撞振动系统周期倍化分岔预测及其神经网络控制

1、一类弹性碰撞振动系统周期倍化分岔预测及其神经网络控制摘要：针对一类单自由度含间隙和预紧弹簧的弹性碰撞振动系统的分岔控制问题，提出了一种基于Lyapunov指 数及径向基函数神经网络的分岔预测及控制方法$首先建立了系统的PoincarO映射，推导了弹性碰撞振动系统周 期运动存在的条件，研究了在主要分岔参数平面中的动力学分布；其次利用Lyapunov指数分析了系统的稳定性, 提出利用追踪Lyapunov指数谱分岔点来预测周期倍化分岔发生的方法；最后基于径向基函数神经网络设计了参 数反馈分岔控制器、基于周期倍化分岔点处的最大Lyapunov指数构造适应度函数，并利用Lyapun。v指数判断是 否实现

2、了分岔控制，以引导自适应混合引力搜索算法对控制器的参数进行优选，从而实现周期倍化分岔控制。关键词：非线性振动；非光滑系统;周期倍化分岔; PoincarO映射; Lyapun。v指数谱; RBF神经网络引言分岔是非线性系统所具备的独特现象且已经成 为非线性动力学不可或缺的组成部分，分岔理论的 研究不仅揭示了系统的各种运动状态之间的相互联 系和转化，而且与混沌密切相关$对于非线性系统， 分岔现象可能产生有害的动力学行为，需要避免或 抑制；又或为了使系统产生人们所需要的分岔行为， 需要设计适当的控制器以改变非线性振动的分岔特 性$因此，对非线性系统的分岔分析和控制的研究 具有重要的科学意义和广阔的

3、应用前景）%*$研究分 岔控制，可以有效地避免、延缓和消除分岔所导致的 不良后果，对提高系统的稳定性和可靠性具有理论 指导意义$非光滑动力系统在机械、电路等领域十分普遍， 它会导致类似于复杂非线性系统具有的分岔和混沌 运动，然而很多分岔特性及机理又与普通光滑非线 性系统完全不同）%$* $非光滑系统中含有的间隙、 预紧、干摩擦等非光滑因素使其PoincarO映射在控 制目标附近不可微，故基于局部线性化映射的各种 控制策略及其推广形式不能胜任这类系统的混沌运 动控制）%*$周期倍化分岔过程是一条通向混沌的典型道 路，因此消除或延迟周期倍化分岔的发生是控制混 沌发生的一个有效方法$ Abed等提出

4、了周期倍 化分岔的局部镇定问题，并采用反馈控制延迟倍周 期分岔的发生$唐驾时等或研究了 Logistic模型 的倍周期分岔的控制问题，设计了各种线性控制器， 使倍周期分岔延迟或提前出现甚至消失$罗晓曙 等利用系统的状态反馈和参数调节的方法，有效 地实现了离散非线性动力系统的倍周期分岔的延迟 控制和混沌吸引子中不稳定周期轨道的控制$王学 梅等逐*根据一般迭代映射的倍周期分岔定理，从数 学上论证了电压型不连续导电模式（DCM）Boost和 Buck变换器中倍周期分岔现象产生的条件，由此 揭示了 DC；DC变换器中倍周期分岔现象发生的机 理$姜海波等基于Floquet理论揭示了 Logistic 映

5、射周期解的分岔机理$卫晓娟等应用基于 RBF神经网络的智能优化控制方法研究一类含间 隙碰撞振动系统混沌运动的控制，将混沌运动控制 为预期的规则运动$文18提出了一种不依赖被 控系统数学模型的无模型自适应参数反馈混沌控制 方法，研究了一类单自由度非光滑系统的混沌控制 问题$文19通过分析对称性破缺分岔机制，采用 了一个直接的、有效的线性控制器，精确控制了一类三次方对称离散混沌系统发生对称性破缺分岔和倍 周期分岔时分岔点的位置！ Go*n- 等2o#在碰撞瞬 时，引进了具有动态变量的超越映射，计算了碰撞系 统的李雅普诺夫指数谱！金俐等21对3维刚性约 束和分段光滑非光滑动力系统引进局部映射，利用

6、映射分析方法得出了非光滑系统Lyapunov指数谱 的通用计算方法。本文针对一类单自由度含间隙和预紧弹簧的弹 性碰撞振动系统的分岔控制问题，应用智能控制方 法，提出了一种基于Lyapunov指数及径向基函数 神经网络的分岔预测及控制方法。本文内容安排如 下：首先建立了系统的 Poincare 映射，推导了弹性 碰撞振动系统周期运动存在的条件。利用Lyapunov 指数分析了系统的稳定性；其次提出追踪 Lyapunov指数谱分岔点来预测周期倍化分岔发生 的方法；最后基于径向基函数神经网络设计了参数 反馈分岔控制器，并基于周期倍化分岔点处的最大 Lyapunov指数构造适应度函数(即利用Lyapu

7、nov 指数预测周期倍化分岔的发生，以及判断是否实现 了分岔的控制)，以引导自适应混合引力搜索算法对 控制器的参数进行优选，从而实现了周期倍化分岔 的控制！1系统模型及其运动方程1.1力学模型图1为一个含间隙及预紧弹簧的弹性碰撞振动 系统模型，左边是质量为M的物块由刚度为K1的 线性弹簧和阻尼系数为41的线性阻尼器连接于支 承，并受到简谐激励sin(6T + ；)的作用(F.为 简谐激励力幅值，6为简谐激励力频率为初始相 位)右边是一个带有预压缩弹簧K2和阻尼系数 为42的碰撞面(无质量)，用来缓冲碰撞。取物块静 平衡位置为坐标原点，弹簧K2被预压缩来缓冲振 动，设质量块的位移为X,间隙为，弹

8、簧X的预压 缩量为D !当X V2时，这是一个简单的谐波振子，方程 如下MX +41X + K1X = F. cos (6T + ;)(1)当X时物块与右边墙面接触发生，因为在 接触表面没有质量，方程如下MX + (41 +42)X+ (K1 +K2X +K2 (D ) = F.cos (6T + ；)(2)取无量纲量：(a)系统模型图(a) Schematic diagram of the system(a)系统模型图(b)刚度模型(b) The stif&iess model图系统模型Fig. 1 System model41,KiDW =则系统的无量纲运动微分方程为+ 2)+ )=cos

9、(/ +(41,KiDW =则系统的无量纲运动微分方程为)+ 2#(1+&广)+ (1 + &k)=(3)cos(/ + ；) &k(W !)()B !)1.2方程的解及周期运动存在条件物块与碰撞面发生接触阶段(即)B !时)，由 初始条件)/ to) =)o ,)/ to) = %o , ； /to) = ；o(4)当#V匕顼四时，方程(3)的通解为 1+#1(/;/o,%o) = e+3(?cos ( $d2 ( / /o) +但si,(皿2 / /$)# +cos ($ + ；o) +Bzsin (/ +；o)1 +&k(5)式中3 = #(1+&r), % = 3 (1 +&k) (1

10、 +&r)2(5)式中A = 1 + &k $ R =2(1 + &k $2)2 + 4#2 (1 + &广)2$2 2(1 + &k (1 + &k $2 )2 + 4#2 (% + 任)2$2。物块与碰撞面没有发生接触阶段(即)V !时)，由初始条件)(1 ) =)1 ,)(1 ) = %1 , ； (1 ) = $ 1 + ；0 = C$(6)当1时，方程(3)的通解为:= e+2aiCosCdi / /%) +Ssin(sd (/ /%) +A1 cos (+co)+D1Sin (/ + co)(7)式中式(5)和式中式(5)和(7)中(1 & 1,2)分别为积分常数和振幅常数文中用q

11、 = in/en来表示碰撞振动系统的周期 运动，其中 $3 为(=*() 8 R2 X S,)= !,)B 0 Poincare截面的点，表示物块与右边碰撞面的碰 撞次数 % e3 为(3 = * (),) 8 R2 X S1 = 0mod(/ = ) Poincare 截面上的点，表示系统的激励力$周期数。可以通过坐标变换将初始位移)0变换到！ 处,q = 1/3周期运动表示如果取物块与碰撞面碰 撞后瞬时为时间坐标原点(即(=0 )，经过(1时间, 物块离开碰撞面；之后经过(2时间，物块又一次和 碰撞面碰撞，那么时间(1和(2必满足(1 + (2 =(T = 23&/ 3(8)则系统q =

12、1/3周期运动的边界条件、周期条件及 运动衔接条件为)(0) = !, )(23&/$) = !,)(0) =)(23&$) = %0, )(t1 ) =)(t1+ ),)(t1 ) =)(t1+ )(9)式中(1表示物块在阶段以负向速度离开碰撞 面的瞬时时刻；(1+表示物块进入)V!阶段的瞬时时 刻。为了便于系统周期解的推导，引入下列符号：1 = 6 / , 2 = 6 3 , Cd1 = cos(3d/) ,Sd1 = sin $1)，Cd2 =cos $d2t)，&2 = sin $d2()，C3 = cos ( ()，s3 = sin ()， d11 = 2Cd2，d12 = 2Sd2

13、 ,d，21 = 2 (C d2 + $d，2 Sd2 )，d-222 ( Sd2 + 3dj Cd2 ) , 11 = A2C3 + B2S3 ,12 = ：2S$B2C3 , 21 = B23C$ A2SS3 , 22=B23S3 ：2$C $，=1Cd1，612=1Sd1，621 = 1(Cd1=B23S3 ：2$C $，=1 ( Sd1 + 3dCd1 )，C11 = A1 C$ + B1 S$， C12 =一 A1 S$ + B1 C$，21 =12，22 = 11，#11 = #22 =cos ($11)，#12 = sin $Z1)，#21 = sin (3Z1)。则系统 的q

14、= 1/3周期解可以表示为3 (t)=()j+E ()3 (t)=(80，( (10)(8 (1,D (Oj+E OFC， 式中 3 = )，)T，D () = (dy- )2x2，E ()= (6ij)2X2，j = 2，!*2 T， E() = (0,gii + gi2 = 0. 010864 丰 0Cl. 2 gii + 屈一I 虹 | = 0. 008943 0条件C1. 1，C1.2满足，所以在参数(I)时，系统存在 q = 1/1周期运动。那么在 (p Poincare 截面上存在 周期一不动点。在此不动点附近系统的Jacobian 矩阵D#(p-) = .#p-) -=-* ,p

15、=p*(16)dX 的一个特征值P%从一1处穿越单位圆，其他特征值 在单位圆内，系统发生周期倍化分岔。其中式(16) 中P表示分岔参数，-表示 (p Poincare 截面上的不 动点坐标同理，在基准参数(I)下，变化3,当3* =0. 01797时同样满足条件C1. 1,C1. 2,此处不再 赘述。接下来数值计算系统其他参数取(I),在$3) 参数平面内的运动分布，如图2所示。图中用不同 颜色标识出了不同周期运动状态所处参数区域。在 $ 8 2.4,3. 2&,3 8 0,0.25参数区间内，周期倍 化分岔序列存在于图中左下角，即间隙和激励力频 率均较小的情况下图2 $ ,3)参数平面运动分

16、布Fig. 2 ($ ,3) parameter plane motion distribution(b)随S变化(b) Against with 8图3系统分岔图和Lyapunov指数谱图Fig. 3 The bifurcation d i agram and Lyapunov exponent spectra of the system图3(a)和(b)为在(pPo incare截面上物块的速 度z随激励频率$和间隙3变化的周期倍化分岔序 列图及Lyapunov指数谱图，即在图2中沿$ = 2. 6 及3= 0. 1两条虚线变化时的分岔图和Lyapunov 指数谱图随着$和3的分别增加，q

17、=1/1周期运 动经历倍化分岔序列，演变为q = 2/2 ,q = 4/4,q = 8/8 ,(周期运动直至混沌。随后混沌吸引子与其 吸引域边界发生边界激变,造成混沌运动的消失,系 统演变为q = 2/3周期运动0根据Lyapunov指数 判据：当系统所有Lyapunov指数兀 0 ,则系统渐 进稳定；当最大Lyapunov指数人 = 0，其他V 0 ,则系统发生分岔；若至少有一个兀A 0 ,则系统进 入混沌运动状态0图中PDl(i = 1,2 ,()表示周期 倍化分岔点，CP$(i = 1,2 ,()表示Lyapunov指数 突变点各点具体数值如表1和2所示。表1周期倍化分岔点Tab. 1

18、Period ic doubling bifurcation pointsPD1Pd2PD3PD#表2 Lyapunov指数突变点Tab. 2 Lyapunov exponential mutation pointsCP1CP2CP3CP4$2.418762.466522.481442.4854830.007980.041490.050400.05268通过图3可以发现在发生周期倍化分岔前,系 统的Lyapunov指数谱已经发生分岔,而周期倍化 分岔是非线性系统通向混沌的典型路径。若能够在 系统发生周期倍化分岔前得到预警,那么在对系统 进行分岔或混沌控制时就可以提前动作，从而降低 控制成本和难

19、度,实现对系统的预测控制。定义 Lyapunov指数谱差值随分岔参数族p变化的函 数为dA(p) = 1(p) +X2(p)(17)式中 1为在参数P时刻系统的最大Lyapunov指 数,为另外一个Lyapunov指数。可以用来确 定Lyapunov指数的突变点。本文在计算Lya； punov指数谱时选取的初始扰动(z；01，技2. 433002. 475562. 484542. 48634 017970. 046830. 052170. 05319)= (0. 001,0. 001)。两个Lyapunov指数未突变或分岔前，dA(A)9 0. 008959 4 0. 01 ,如图 4(a)和

20、(b)的 局部放大图(al )和(bl )所示$当d(G) 0. 008959 4 0. 01时，两个Lyapunov指数开始分 离，各自向相反的方向逐渐增大$随后AX(p)的值 连续增加，当到达周期倍化分岔点时，即当1100 时，达到局部最大值=0. 423$随后又逐渐 减小至2. 433002. 475562. 484542. 48634 017970. 046830. 052170. 05319定义：突变阈值dACp为系统的Lyapunov指数 谱发生分岔时的值；分岔阈值dAPD为系统发生 周期倍化分岔时的值$故本文所述系统dACp = 0. 01，d = 0. 423 $ 另外令.巳+

21、1 巳1m尸q _尸qLc3+23+1(18)计算得g。4 3. 447,g. 4 3. 969$其中g。为系 统随参数3增加时Lyapunov指数相邻两突变点之 间差值的比率& g.为系统随参数.增加时Lya- punov(18)通过分析可以得出，对于本文中的非光滑系统,随着分岔参数的变化，当dA通过分析可以得出，对于本文中的非光滑系统,定的周期运动状态；随着系统参数的变化，当逐 渐增大到=Apd ,且1 = 0时系统开始周期 倍化分岔序列，Lyapunov指数突变点的变化倍率为 g。4 3. 693,g. 4 3. 908 $虽然本文中dACP的阈值 为随分岔参数变化时Lyapunov指数

22、的观察值，但 随参数族p变化曲线在Lyapunov指数突变点 CP,($ = 1,2,)邻域内的斜率为1,也就是说其突变点的位置不会随的阈值的变化而变化，所以 Lyapunov指数相邻两突变点之间差值的比率g不 会随的阈值变化而变化。但的阈值及g。和 g.的普适性需要进一步研究$3 基于RBF神经网络的非光滑系统 分岔控制器设计及优化非光滑系统中某些主要分岔参数的变化会引起 周期倍化分岔序列的发生，在倍周期分岔点处系统 的最大Lyapunov指数为零，所以本文基于参数反 馈分岔控制法的原理，利用系统的Lyapunov指数 的变化规律来控制倍周期分岔现象$本文基于径向 基函数神经网络(22) (

23、Radial Basis Function Neural Network,简称RBFNN)设计参数反馈分岔控制器，图4系统Lyapun。v指数谱差值随。和.变化趋势Fig. 4 The variation trend of system Lyapunov exponent spectrum difference with。and .以输出一个小扰动施加于系统的可控参数，通过对 系统可控参数进行动态微幅调整，从而将周期倍化 分岔后产生的L = 2/2,q = 4/4,运动控制到q = 1/1运动状态，即实现了周期倍化分岔的控制$根据上节周期倍化分岔预测所得结论，当cUAcUcp时 系统即将发生周

24、期倍化分岔序列；当dA = dApD ,且 |妇1 = 0时系统进入周期倍化分岔序列#所以可以 在dAw dA图4系统Lyapun。v指数谱差值随。和.变化趋势Fig. 4 The variation trend of system Lyapunov exponent spectrum difference with。and .3.1控制器设计本文设计的基于RBFNN的参数反馈分岔控制 器为具有单稳层的三层前馈网络，包括输入层、隐含 层和输出层#通过输入层空间到隐含层空间的非线 性变换以及隐含层空间到输出层空间的线性变换， 实现输入层空间到输出层空间的映射#根据倍周期 分岔控制目标，将灸次迭代

25、后稳定的q = l/l周期运 动的状态x =际,丁作为控制器的输入,控制器的输出定义为系统激励频率或间隙的微幅调 整量#由此，确定RBF神经网络输入层为个节 点、输出层为l个节点，控制器结构如图5所示#图5 RBF神经网络结构图Fig. 5 The RBF neural network structure diagram控制系统结构框图如图6所示图6控制系统结构框图Fig. 6 The block diagram of the control system图6中Z图5中Xi和X为输入向量，假设有.个隐节 点，RBF控制系统结构框图如图6所示图6控制系统结构框图Fig. 6 The block

26、diagram of the control system图5中Xi和X为输入向量，假设有.个隐节 点，RBF隐层的径向基函数向量为h = &h,T ,h,为 隐层的第R个节点的高斯函数值，且h, = exp （ N ）（19）1llCl.）式中c = &C$,2 =:表示第/个输入Cllc.的第,个神经元的高斯函数坐标值的中心，$ = l, ，3, = l。向量 b = &!l ,.T 表 示第，个节点高斯函数的宽度。牝为控制器输出的扰 动量（即施加于系统可控参数上的微幅调整量）,为 保持控制的有效性,设定最大扰动量为Fmax,则 Umax F Fmax :=叫，T为连接隐层 和输出层的权矩

27、阵。3.2基于AHGSA的控制器参数优化采用文&l7提出的自适应混合引力搜索算法 （简称AHGSA算法）进行控制器参数优化,以避免 控制器参数选择的主观性和盲目性# AHGSA算法 是GSA算法经过改进而提出的,其优点是克服了 GSA算法全局探索能力强而局部开发能力弱的特 点,提高了算法的整体寻优效率#3.2. l 适应度函数的建立假设粒子群规模为广,粒子$的位置用5维向量 *,表示,s = nX m+2m （ n为分岔控制器输入变量 的个数,m为RBF神经网络隐层节点数）,则粒子群 的位置可以用矩阵*X表示。本文选定系统的最大Lyapunov指数作为判断依据，建立控制器参数优化选择时应满足的

28、适应度函数为0#(*$) = .l,k = l.2 , | U(k)k = l( 20 )式中Z*为（p Poincare( 20 )的最大Lyapunov指数（即图3中的PD$ （/ = l , 2 , ）;Amax为（p Poincare截面上倍化分岔后的周期 运动的最大Lyapunov指数（即未受控系统的最大 Lyapunov指数）；Amax为受控系统的最大Lyapunov 指数；e为一个小正数；.l , .2为相对权重; u（1为k时刻控制器的输出;L为输入/输出数据序列长度3.2.2基于AHGSA的控制器参数优化流程采用AHGSA算法使式（20）取得最小值的全 局最优位置Pg即是RBF神经网络参数反馈分岔控 制器的最优参数_ !和c,其具体优化流程如下:1）随机初始化粒子的位置和速度；2）按式（20）计算每个粒子的适应值；3）对每个粒子，将其适应值与个体历史最优、 以及与全局最优进行比较，以更新当前的个体最优 位置与群体最优位置；4）更新每个粒子的速度和位置；5）如未达到预先设定的停止准则，则返回步骤 2）,若达到则停止计算&8仿真研究在确定神经网络隐层节点数时，在满足控制系 统性能要求的前提下，以取尽可能紧凑的网络结构 为原则，将网络隐层节点选定为5个& AHGSA算 法的参数设置为：种群规模为30,最大迭代次数为