一种高维非线性系统全局形态分析方法

技术领域

本发明涉及机械设备动力学分析技术领域，具体地是指一种高维非线性系统全局形态分析方法。

背景技术

非线性能量阱由于其具有附加质量小、能量定向传递、吸振频带宽等优点，近二十年来一直是被动控制技术领域的研究热点，并相继在机械工程、土木工程、航天工程、声学控制、振动能量采集等诸多领域得到了应用。非线性能量阱系统属于高维强非线性范畴，其中存在丰富的非线性动力学行为。但是，目前大部分关于非线性能量阱系统的研究工作，都是在确定性系统理论的范围内进行的，而实际上非线性能量阱应用于干扰源众多的物理系统时，初始条件具有很大的不确定性，极易引起系统产生非线性振荡、跃迁至大振幅状态甚至崩溃。因此，研究非线性能量阱系统的全局性态是使其安全可靠运行的重要保证。

近年来，对低维非线性动力系统的研究已取得一系列重要成果，但研究高维非线性系统的全局形态分析方法还不是很多，国际上处于发展阶段，国内尚处于起步阶段。针对开展高维非线性能量阱系统全局分析时，存在内存消耗大以及计算效率低的问题。

发明内容

本发明的目的在于克服现有技术中的不足之处，而提出一种高维非线性系统全局形态分析方法，用于定性分析周期、准周期或混沌等不同响应机制，能够将交叉映射序列合并，避免元胞被重复计算。

为实现上述目的，本发明所设计的一种高维非线性系统全局形态分析方法，其特殊之处在于，所述方法包括如下步骤：

1)建立待减振机械设备动力学系统模型，并改写为空间状态方程组形式；

2)将动力学系统连续空间离散化为胞化空间，选定分析平面，建立元胞矩阵Z；

3)进行参数初始化处理，建立op向量用于标识矩阵Zs的空行，并将映射序列轨迹长度初始化为Ls＝0；

4)根据行标索引r递增顺序依次从元胞矩阵Z中筛选尚未处理的元胞用于填充标识矩阵Zs的空行；通过扫描矩阵pc，若返回值B＝1，则判定元胞为已处理胞，并赋予相同的全局性态标识数I，若B＝0，则判定元胞为未处理胞，并填充至标识矩阵Zs中；当所有待处理胞的返回值B＝0时，表明矩阵Zs填充完毕，并更新矩阵pc的周期数号和步号；

5)将填充完毕的标识矩阵Zs并行数值计算得到胞映射的像胞，并重新赋值给标识矩阵Zs；

6)对经过像胞处理的标识矩阵Zs进行扫描判定；若像胞为新的未处理胞，将其添加到矩阵pc中，并形成新的周期序列；若像胞为已处理胞，其演化过程对应三种情形：

(i)像胞对应矩阵pc中一条已完成映射序列，则赋予该像胞相同周期号与步号；

(ii)像胞与原胞为同一个周期胞，则将其添加到矩阵pc中，形成一条新的周期序列，并赋予新的周期组号；

(iii)像胞与另一条正在处理映射序列相交叉，则将其添加至交叉映射序列中，更新交叉映射序列元胞个数Ks，清空像胞所在行，并从此行开始一条新的映射序列；

7)映射序列后处理；检查映射序列轨迹长度并识别运动状态，当映射序列轨迹大于最长步数max时，则终止映射，将其视为混沌运动；当所有元胞被遍历计算时，则依据每个元胞的周期号组以及步数号组标识数，得到完整的吸引子以及吸引域信息。

进一步地，所述步骤1)的具体步骤包括：

11)建立待减振机械设备的吸振系统动力学模型：吸振系统由待减振机械设备m

12)对吸振系统作假设包括：(1)待减振机械设备以及非线性能量阱形状规则、质量分布均匀，外界激励力作用于机械设备质心，只产生垂向振动；(2)忽略弹簧、阻尼等元件质量的影响；(3)不考虑机械设备和非线性能量阱重力以及弹簧预应力对静平衡位置的影响；(4)所承载基座为质量无限大的刚体；

23)将系统运动学方程改写为空间状态方程组形式

式中，t＝Ω

优选地，所述步骤5)的并行计算的扩展控制方程为：

其中，

优选地，所述步骤5)中胞映射的过程为：从原相空间任意一点出发，沿相轨迹数值积分一段时长后的终点坐标

与现有技术相比，本发明设计的一种高维非线性系统全局形态分析方法，具有如下优点：

1、本发明提出了一种并行多自由度胞映射方法，并利用该方法分析了非线性能量阱系统的全局性态；该方法在多自由度胞映射方法的基础上，引入并行重构算法，对胞映射和胞处理过程进行改进，能够同时计算多个胞映射步骤。

2、本发明用于定性分析周期、准周期或混沌等不同响应机制，应用于非线性能量阱系统的全局性态分析，在离散空间内通过观察离散流的周期性来判断动力学系统的运动状态，刻画出典型参数范围内多稳定吸引子的全局性态，揭示系统吸引子与吸引域之间的确定性关系。

3、本发明能够确定吸引子数目、大小以及对应的吸引域，刻画非线性能量阱系统振动响应完整样貌。

4、本发明在映射过程中将交叉映射序列合并，避免元胞被重复计算。

附图说明

图1为本发明高维非线性系统全局形态分析方法的流程图。

图2为本发明研究对象的机械设备耦合非线性能量阱的吸振系统动力学模型。

图3为图2的实物图。

图4为机械设备响应随激励力幅值变化的全局性态图。

图5为对应的最大Lyapunov指数图谱。

图6为不同激励力幅值下机械设备响应类型。

图7为MDCM法与PMDCM法计算过程示意图。

图8为实施例一的系统响应不同吸引子相轨迹示意图。

图9为实施例一的系统响应共存吸引子及吸引域示意图。

图10为实施例二的系统响应不同吸引子相轨迹示意图。

图11为实施例二的系统响应共存吸引子及吸引域示意图。

具体实施方式

以下结合附图和具体实施例对本发明作进一步的详细描述。

本发明提出的高维非线性系统全局形态分析方法(PMDCM)中定义的主要变量以及含义如下：

Z：矩阵用于存储胞化空间子集S中的所有元胞，Z∈M(M,N)

Z

pc：有序数组用于存储所有已处理胞，pc∈M(N

op：向量用于标识Z

g：数组用于标识pc中已处理胞的周期数，

inds：数组用于标识pc正在处理胞序列中胞的步号，inds∈M(N

M：子集S中元胞的个数，

N：动力学系统维数，

N

N

N

N

Ks：每条映射序列元胞个数，

Ls：每条映射序列最长轨迹长度，

max：自定义的映射序列最长步数，

r：矩阵Z

ri：矩阵Z中未处理胞的最小行标

如图1所示，本发明提出一种高维非线性系统全局形态分析方法，包括如下步骤：

1)建立待减振机械设备动力学系统模型，并改写为空间状态方程组形式；

2)将动力学系统连续空间离散化为胞化空间，选定分析平面，建立元胞矩阵Z；将动力学系统连续相空间离散化为胞化空间，并在选定分析平面内得到元胞矩阵Z。

3)进行参数初始化处理，建立op向量用于标识矩阵Zs的空行，并将映射序列轨迹长度初始化为Ls＝0，为后续填充矩阵Zs做准备。

4)元胞前处理及填充矩阵Zs。根据行标索引r递增顺序依次从元胞矩阵Z中筛选尚未处理的元胞用于填充标识矩阵Zs的空行；通过扫描矩阵pc，若返回值B＝1，则判定元胞为已处理胞，并赋予相同的全局性态标识数I，若B＝0，则判定元胞为未处理胞，并填充至标识矩阵Zs中；当所有待处理胞的返回值B＝0时，表明矩阵Zs填充完毕，并更新矩阵pc的周期数号和步号；

5)并行数值积分计算。将填充完毕的标识矩阵Zs并行数值计算得到胞映射的像胞，并重新赋值给标识矩阵Zs；在并行计算过程中需引入如下扩展控制方程：

其中，

6)像胞扫描及处理过程。对经过像胞处理的标识矩阵Zs进行扫描判定；若像胞为新的未处理胞，将其添加到矩阵pc中，并形成新的周期序列；若像胞为已处理胞，其演化过程对应三种情形：

(i)像胞对应矩阵pc中一条已完成映射序列，则赋予该像胞相同周期号与步号；

(ii)像胞与原胞为同一个周期胞，则将其添加到矩阵pc中，形成一条新的周期序列，并赋予新的周期组号；

(iii)像胞与另一条正在处理映射序列相交叉，则将其添加至交叉映射序列中，更新交叉映射序列元胞个数Ks，清空像胞所在行，并从此行开始一条新的映射序列；

7)映射序列后处理。检查映射序列轨迹长度并识别运动状态，当映射序列轨迹大于最长步数max时，则终止映射，将其视为混沌运动；当所有元胞被遍历计算时，则依据每个元胞的周期号组以及步数号组标识数，得到完整的吸引子以及吸引域信息。

本发明的提出基于Poincaré映射，Poincaré映射是一种定性分析非线性系统动力学特性以及拓扑性质的有效手段，其基本思想是基于微分同胚概念对n维连续系动力学系统进行离散系降维，通过构造一个n-1维离散时间系统的闪频映射取代原连续时间系统的流，从而在离散空间内通过观察离散流的周期性来判断动力学系统的运动状态。

考察如下任意n维非自治动力学系统

其中，F(x,t)为具有周期性的向量函数，存在F(x,t+T

则经过一个周期T

其中，

Poincaré映射关键之处在于构造合理的Poincaré截面从而获得映射点集信息，该方法与系统是否存在小的摄动量无关。

本实施例以待减振机械设备动力学系统为研究对象，其机械设备耦合非线性能量阱的吸振系统动力学模型如图2所示，实物图如图3所示，主系统由m

为简化分析，作如下几点假设：

(1)待减振机械设备以及非线性能量阱形状规则、质量分布均匀，外界激励力作用于机械设备质心，只产生垂向振动；

(2)忽略弹簧、阻尼等元件质量的影响；

(3)不考虑机械设备和非线性能量阱重力以及弹簧预应力对静平衡位置的影响；

(4)所承载基座为质量无限大的刚体。

考虑图2所示的系统动力学模型，此时激励形式满足f

首先引入长度量纲l

t＝Ω

由于在两个时间尺度条件下，存在微分运算关系

将式(6)和式(7)代入到式(5)中，可得到简谐激励条件下系统无量纲化的运动学方程

其中，

对式(8)而言，机械主系统质量可视为单位1，ξ

重新考察系统(8)，将其改写为空间状态方程组形式

其中，x＝(z

按照上述原理，若将时间序列对应的所有映射点集位移分量作为纵坐标，控制参数作为横坐标，便构成了二维平面内系统的全局性态图。

对系统进行全局性态分析。固定参数ξ

图4为机械设备响应随激励力幅值变化的全局性态图，图5为对应的最大Lyapunov指数图谱。从图中可知，无论是向前延拓或向后延拓，激励力幅值在大参数范围内变化时系统均呈现非常丰富的动力学特性，周期运动(周期1、周期2、周期3等)、准周期以及混沌运动均存在于系统解枝中，其中稳定周期运动占主体地位，而准周期运动、混沌运动主要位于某些特定激励力幅值区间内。以向前延拓计算结果为例，不同激励力幅值下机械设备响应类型如图6所示，其中蓝色曲线代表相轨迹，红色标记代表Poincaré截点。

另外，在图4中还可观察到倍周期分岔、逆分岔以及混沌窗等复杂动力学现象，而促使本为线性特征的主系统出现准周期、混沌等复杂行为的诱因主要是非线性能量阱的调控作用，导致系统响应在特定激励力幅值处出现了分岔。

对比向前与向后延拓计算结果不难看出，非线性能量阱系统在多个控制参数区域内出现了多稳定吸引子共存现象，典型区域包括：

①共存区域1：0.19≤f≤1.11，两个周期1吸引子共存；

②共存区域2：10.5≤f≤11.86，周期1吸引子与准周期吸引子共存；

③共存区域3：11.87≤f≤16.96，两个周期1吸引子共存；

④共存区域4：19.01≤f≤19.26，两个周期3吸引子共存；

⑤共存区域5：19.34≤f≤19.62，一个周期3吸引子与混沌吸引子共存；

⑥共存区域6：25.19≤f≤26.01，一个周期2吸引子与混沌吸引子共存。

多稳定吸引子共存是由于非线性系统多值性本质决定的，每个吸引子有其对应吸引域，吸引域代表汇聚到同一吸引子的所有初始条件集合，而初始条件决定了系统响应的发展过程和最终性态。因此，进一步确定吸引子数目、大小以及对应的吸引域对于刻画非线性能量阱系统振动响应完整样貌异常重要。

本发明在多自由度胞映射法(Multi-Degrees-of-Freedom Cell Mapping，MDCM)的基础上，引入并行计算扩展控制方程对MDCM法胞处理过程进行改进，拓展了一种并行多自由度胞映射法(Parallelized Multi-Degrees-of-Freedom Cell Mapping，PMDCM)。

将n维非自治动力学系统的式(2)写为分量形式

胞映射法的基本思想是通过引入胞的概念将连续相空间离散化，从而将胞状态用胞化空间内一系列有间隔整数唯一标识，记元胞z＝[z

c

那么，每一个元胞边界所包含的原相空间状态向量为

(z

若h

可用二维胞化空间重构，元胞子集S满足

其中，

经过上述过程，从原相空间任意一点出发，沿相轨迹数值积分一段时长后的终点坐标

MDCM法在具体实施过程中，将子集S中的原胞分为原始胞、正在处理胞、已处理胞和再次被处理的已处理胞四类，其实质是从所选分析平面中某一胞中心点出发，在胞化空间中不断向前映射直至映射轨迹再次经过已处理胞时映射过程终止，并依次从新的胞中心点重新开始映射，直至所有元胞计算完毕。每条完整映射序列可记为

z→ζ(z)→ζ

其中，m为映射步数。通过逐一分析每条映射序列的演化过程，获取每个元胞的周期信息，最终可得到整个胞化空间的全部吸引子和对应吸引域。PMDCM法则是对MDCM法的胞映射过程和胞处理过程进行改进，现以二维分析平面中3×3大小的子集S为例，通过示意图来直观描述两种计算方法的区别，如图7所示。

在图7中，通过本发明PMDCM法并行计算数N

一是在MDCM法中，每次只能处理一个元胞并向前持续映射，直至映射轨迹收敛于吸引子时终止，且映射过程与扫描处理过程是交替进行的；而在PMDCM法中，多条映射序列可并行向前映射，一旦当前映射轨迹经过已处理胞时映射过程随即终止，不再向前映射，且映射过程与扫描过程是同时进行的。

二是PMDCM法可能存在两条或两条以上映射轨迹相互交叉的情形，因此需要增加新的胞扫描判定过程，其功能在于将交叉映射序列合并，避免元胞被重复计算。

虽然两种胞映射方法最终计算结果一致，但显然前者存在被重复处理的元胞，而后者能保证每个元胞只被计算一次，从而能够极大的节省计算时间。

实施例一

选取f＝0.68，此时系统响应存在两个稳定周期1吸引子，不同吸引子相轨迹如图8所示，其中A

实施例二

选取f＝11.3，此时系统响应存在一个稳定的周期1吸引子和一个稳定的准周期吸引子，不同吸引子相轨迹如图10所示。若选择分析平面-5≤z

最后需要说明的是，以上具体实施方式仅用以说明本专利技术方案而非限制，本领域的普通技术人员应当理解，可以对本专利的技术方案进行修改或者等同替换，而不脱离本专利技术方案的精神和范围，其均应涵盖在本专利的权利要求范围当中。