引言

由于索结构线形美观、柔度大、质量轻以及强度高等特点，被广泛运用于各类工程中.众所周知，各类索结构极易受到外界环境因素影响^[1]，例如在太阳辐射等影响下，温度场会不断变化，从而引起结构热胀冷缩.国内外学者对索结构在温度场中的静、动力特性开展了大量研究^[2-6].此外，受支座位移运动引发的参数共振，将导致拉索展现出复杂的非线性动力学行为^[7，8].倘若进一步考虑温度变化影响，系统的振动特性会产生明显的定性和定量变化，进而影响索结构的安全性和稳定性.

近年来随着全局分岔理论的不断完善，诸如混沌之类的复杂非线性动力现象被人们进一步了解和认知.张伟等^[9]全面系统总结了高维非线性系统中多脉冲混沌理论，指出了全局分岔理论在工程应用中的现状和展望.尤其对于悬索这类典型的、同时包含平方和立方非线性的系统，全局分岔理论可以很全面地描述其非线性动力学行为.Zhang和Tang^[10]基于全局摄动方法研究了悬索在参数激励和外界环境激励耦合作用下的全局分岔现象和混沌行为.Chen和Xu^[11]基于能量相位法研究了斜拉索的Shilnikov同宿轨道以及混沌动力学行为.Zhang等^[12]借助高维Melnikov方法证明了张紧弦在1∶2内共振情况下存在Shilnikov同宿轨道.

除悬索、斜拉索和弦等张力结构外，有学者通过广义Melnikov方法和能量相位法研究了板、壳和天线等基本结构的混沌动力学和多脉冲轨道.An和Chen^[13]研究了气动热弹性功能梯度材料（FGM）截锥壳在复杂载荷作用下1∶2内共振和主参数共振情况下全局分岔和多脉冲轨道.高美娟等^[14]依据三阶剪切变形层合板理论，研究压电复合材料层合板在1∶2∶4内共振情况下的全局动力学.孙莹等^[15]通过能量相位法研究热应力下环形桁架天线混沌运动.

针对悬索这类典型的工程柔性结构，Zhao等^[5，6]基于局部分岔分析发现:温度变化会引起动力系统参数的微小改变，进而导致其整体非线性振动特性产生明显定性和定量的改变.倘若进一步采用全局分岔分析，可以从能量轨道的角度更好解释以及探究动力系统复杂的运动行为^[9].由此可见:一方面温度变化对支座位移运动引发的参数激励下悬索动力学行为影响有待进一步研究; 另一方面基于全局分岔理论描述温度变化对系统动力学行为的影响更鲜有报道.因此本文将局部分岔分析^[6]拓展到全局分岔分析，采用Galerkin法和多尺度法并基于规范型理论和能量相位法，系统探究温度变化对悬索全局动力行为尤其是多脉冲轨道的影响，揭示温度变化下动力系统可能出现的大幅振动机理.

1 运动方程及摄动分析

1.1 数学模型

如图1所示，悬索水平悬挂于O和B之间，以O为原点，建立直角坐标系O-xy，L和b分别代表悬索的跨度和初始垂度.本文考虑支座位移引发的共振响应，图中Δ表示支座位移最大值，水平和垂直方向最大位移分量分别为Δ_x和Δ_y，支座运动Δ_A表示为:Δ_A（t）=ΔcosΩt，其中Ω表示支座运动的频率.

当环境温度发生整体改变时，基于增量热场理论，悬索将产生热应力构形^[2，3]，此时张力改变，垂跨比随之变化，悬索的轴向和竖直方向位移分别用u（x，t）和v（x，t）表示.当不考虑温度变化影响时，利用Hamilton变分原理，可得悬索面内运动微分方程^[16]:

其中:点表示对t求导，撇表示对x求导，EA表示抗拉刚度，H表示初始水平张力，c_v表示阻尼系数，悬索的静态构形用抛物线表示:y=4bx（L-x）/L²。

在支座运动引发的参数激励作用下，悬索的总加速度${\ddot{v}}_t$由两部分组成:静态位移产生的加速度${\ddot{v}}_0$以及动态应变产生的加速度${\ddot{v}}_{}$.由于垂度较小，分析时可以假设${\ddot{v}}_0$是线性分布，此时${\ddot{v}}_0$可表示为: ${\ddot{v}}_0（x，t）=-（1-x）{\mathrm{\Omega }}^2{\mathrm{\Delta }}_y\mathrm{c}\mathrm{o}\mathrm{s}\mathrm{\Omega }t$.

为描述简便，引入以下无量纲参数:

其中:ω表示固有频率，α表示静态拉伸系数，ζ_v示悬索垂向阻尼比.为表示简便，忽略上画线，可得支座运动激励下悬索面内运动方程:

已有研究表明^[5]:温度变化会改变拉索张力以及垂跨比，其影响可通过引入张力改变系数χ²来表示.因此在无量纲方程的基础上，可得考虑温度变化时悬索受参数激励的运动微分方程为^[5，6]:

采用Galerkin法对方程进行离散:

其中:考虑温度影响下的模态函数φ_n（x）参考文献^[5].

本文考虑仅发生内共振的两个模态，令N=2，将方程（5）带入方程（4）中，方程两边同时分别乘模态函数φ_n（x），并从0到1积分，可得两自由度的运动微分方程组:

其中:相关系数如附录所示.

1.2 摄动分析

采用多尺度法对方程（6）、（7）进行求解，引入摄动参数ε以及三个时间尺度: $T_0={\epsilon }^{1}t，T_1={\epsilon }^{2}t，T_2={\epsilon }^{3}t$，将参数激励项和阻尼项摄动到二阶: $K_n={\epsilon }^2{K}_n，P_n={\epsilon }^2{P}_n，{\mu }_{\mathrm{1，2}}={\epsilon }^2{\mu }_{\mathrm{1，2}}$.外激励和固有频率间关系通过引入调谐参数σ₁和σ₂来表示:Ω=2ω₁+εσ₁，Ω=ω₂+εσ₂.

方程的解可以假设为$q_i（t，\epsilon ）=\sum _{j=1}^3{\epsilon }^j{q}_{i，j}$，将q_i（t，ε）带入式（6）、（7）中，并令两边εⁱ的系数相等.令q_i（t，ε）两边摄动参数系数相等，求解可得q_n，1的表达式: $q_{n，i}=A_n\left(T_1，T_2\right){\mathrm{e}}^{\mathrm{i}{\omega }_n{T}_0}+{\stackrel{-}{A}}_n\left(T_1，\right.\left.T_2\right){\mathrm{e}}^{-\mathrm{i}{\omega }_n{T}_0}（n=\mathrm{1，2}）$，其中A和${\stackrel{-}{A}}_{}$为一对共轭复数.将q_n，1带入q_i（t，ε）中进行求解，通过消除长期项，可以得到方程（6）、（7）的平均方程:

引入直角坐标: $A_m（t）=\left[p_m（t）-\mathrm{i}{q}_m（t）\right]{\mathrm{e}}^{i\beta m}/2$，带入方程（8）、（9）中，分离虚部和实部，得直角坐标形式下的平均方程:

其中:$f=-P_2/2{\omega }_2，k=K_1/4{\omega }_1$，其余系数如附录所示.

方程（10）~（13）确定了参数激励下悬索面内振动的相位和振幅，因此通过该方程组研究悬索在温度场中多脉冲混沌动力学行为.为了得到归一化方程，需对平均方程（10）~（13）进行化简，使用Maple程序进行计算^[17]，引入${\mu }_i=\epsilon {\mu }_i（i=\mathrm{1，2}），f=\epsilon f$以及线性变换: $p_1=\left(1-k+{\sigma }_1/2\right)u，p_2=\sqrt{2I}\mathrm{c}\mathrm{o}\mathrm{s}\phi ，q_1={\mu }_{1}u+v，q_2=\sqrt{2I}\mathrm{s}\mathrm{i}\mathrm{n}\phi$.

方程（10）~（13）可表示为可积的Hamilton系统:

其中:${\stackrel{-}{\mu }}_1={\mu }_1^2+\left(k+{\sigma }_1/2\right)\left(k+{\sigma }_1/2-1\right)$.

因此系统（14）~（17）的哈密顿函数表示为:

2 无扰动系统动力学

假设系统（14）~（17）中的扰动参数ε=0，那么余下系统中（14）、（15）不含φ，将I视作常数，因此可将系统对于（u，v）平面和（I，φ）平面进行解耦，得:

此时系统（20）的哈密顿函数简化为:

其中:$\eta ={\stackrel{-}{\mu }}_1-2I{\gamma }_2$.

基于方程（21），绘制系统（20）在三维相空间中的流形结构图.如图2所示，当η＞0，系统（20）将出现异宿分岔.当η＜0，系统的唯一解为零解（u，v）=（0，0），该零解为鞍点.在曲线${\stackrel{-}{\mu }}_1-2I{\gamma }_2$上，即I=${\stackrel{-}{\mu }}_1$/2γ₂的零解通过叉型分岔由一个变为三个，不妨将三个解分别记为q₀（I）=（0，0），q_±（I）=（D，0），由方程（20）的零解条件可以得到D:

观察图2（b）可知q₀（I）为轨道中心，q_±（I）为轨道鞍点.由于异宿轨道可能产生混沌运动，故本文主要研究图2（b）这类情况.对于所有的I∈[0，I₁]，方程（20）有一个鞍点q_±（I）和它自身的异宿轨道相连.

四维空间（u，v，I，φ）由下式定义:

系统（20）异宿轨道的表达式:

其中:$\delta =（\sqrt{2}/2）T_1$.

将无扰动系统限制在M₀上再进行分析，对于任意I∈[0，I₁]，满足D_IH[q₀（I），I]=0的值I被称为共振值，记作I_r，当I=I_r时可以得到此时连接异宿轨道的相位漂移角:

3 多脉冲混沌轨道

由于不变流形M₀在微小的非零扰动ε影响下是不动的，所以M₀在小扰动的持续作用下变为M_ε:

将系统（14）~（17）限制在M_ε上后，对$\sqrt{\epsilon }$进行Taylor展开，同时引入变换$I=I_r+\sqrt{\epsilon }h$及τ=$\sqrt{\epsilon }$T₁，得:

其中: $\dot{h}$，$\dot{\phi }$分别表示对τ求导.当扰动参数ε趋向于零时，令方程（27）、（28）中$\dot{h}$=0，$\dot{\phi }$=0，得到两个不动点:

其中:不动点p是中心，q是鞍点.经过鞍点q存在一个同宿轨道与自身相连，同时不动点p周围有许多周期轨道.当扰动参数ε足够小时，中心不动点p会下沉变为p_ε，连接q_ε的同宿轨道破裂，q_ε不稳定流形逐渐趋近于p_ε.

方程（27）、（28）的哈密顿函数数值保持不变，可得φ_n与φ_s的关系式:

基于Haller和Wiggins给出的表达式^[18]，耗散情况下n阶能量差函数表示为:

耗散情况下n阶能量差分函数转化为:

方程（27）、（28）存在非退化平衡点p≡（h，φ_c）=（0，φ_c），通过中心点p的能量差函数横截零点可得:

假设μ₁=μ₂=μ，当Δφ满足条件Δφ≠2mπ/N，m∈Z，同时系统满足下列非退化条件，据此可判断系统是否在N-脉冲轨道上发生跳跃:

$\begin{array}{c}\sqrt{2I_r}\left\{\mathrm{c}\mathrm{o}\mathrm{s}\left[-\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{s}\mathrm{i}\mathrm{n}\left(\mu \sqrt{2I_r}/f\right)+n\mathrm{\Delta }\phi \right]-\right.\\ \left.\mathrm{c}\mathrm{o}\mathrm{s}\left[-\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{s}\mathrm{i}\mathrm{n}\left(\mu \sqrt{2I_r}/f\right)\right]\right\}-\end{array}$

最后，需要确定从慢流形上起跳的任意一个N-脉冲轨道的降落点必须落在焦点的吸引域A_ε之内.考虑区间[-π，π]内的一点，它的降落点与φ_c+NΔφ的相位差约2kπ，降落点可以表示为:

其中:φ_s和φ_c如公式（29）、（30）中所表示.如果${\phi }_f^N$＞φ_s，将${\phi }_f^N$减去2kπ后重新定义${\phi }_f^N$，验证降落点平移2kπ后是否还在鞍点φ_s附近，如果${\phi }_f^N$的能量大于鞍点φ_s的能量，则:

当降落点落在吸引域之内，系统可能存在混沌运动.

4 数值算例和分析

本文的数值算例中，悬索的参数选取^[19]: L=850m，ρ=8376kg/m³，E=1.794×10¹¹Pa，A=0.1159m²，μ₁= μ₂=0.005，α=1.2×10^-5/℃，结合实际工程情况，温度变化ΔT=±40°C.图3给出了温度变化以及Irvine参数（logλ²）对悬索的前六阶固有频率的影响.图中放大处，二阶正对称模态和一阶反对称模态频率之间呈现两倍关系，系统可能发生2:1内共振，温度变化下系统的各参数和系数大小如表1所示.

采用四阶龙格库塔法计算方程组（10）~（13），并计算该四维系统的Lyapunov指数.假设支座横向振动幅值Δ_x=0.0001，当调谐参数σ₂=-0.10时，选取脉冲轨道数N=3.此时三组γ₁均大于0，为了使系统在该环境下产生异宿轨道，必须让方程（22）中的D²为正数，即$\eta ={\stackrel{-}{\mu }}_1-2I{\gamma }_2>0$.

图4表示不同温度场中悬索产生异宿轨道的条件.当ΔT=-40℃或0°C时，η＞0的解集为闭集，而当ΔT=+40℃时，η＞0的解集为开集.选取-40℃这组参数，当σ₁在[0.30775，1.25072]范围内，系统（20）产生异宿轨道.当σ₁=0.483时，带入方程（34）得到满足能量差函数横截零点的激励f=1.4724×10^-4.此时共振值I=I_r=3.4671×10^-8，相位漂移值Δφ=2.0862.通过求解方程（36），可得: ${\phi }_f^N$=6.2409，φ_s=3.1595，cos（${\phi }_f^N$）-cos（φ_s）=1.9989，$\mu \sqrt{2I_r}\left({\phi }_f^N-{\phi }_s\right)/f=0.0276$.系统参数满足方程（37）、（38），即可能出现混沌响应.

图5分别表示ΔT=0℃，-40℃，40℃时，该动力系统的时程曲线、相位图和庞加莱截面.如图所示，温度变化ΔT为0℃以及40℃，图5（a）和5（c）时程曲线在一段时间运动后将逐渐趋于稳定，最大无量纲振幅减小约为0.0002.二者相位图呈现为闭合圆环，庞加莱截面均表现为杂乱但有限的点集.由表2可知，这两种温度变化条件下，其Lyapunov指数均为（零，负，负，负），即可判定，当0℃及40℃时，系统吸引子类型均为周期型.

然而，当ΔT=-40℃时，系数满足方程（37）-（38）条件，即图5（b）所示的运动存在多脉冲Shilnikov型轨道.对比图5中三组时程曲线，当ΔT=-40℃时，很长一段时间后运动仍没有趋于稳定，最大振幅约为0.001.此时图5（b）中相位图和庞加莱截面图均表现出系统存在无法预知的运动行为.同理计算得到的-40℃的Lyapunov指数为（正，零，负，负），如表2所示，其中最大Lyapunov指数大于零，即判定该运动存在混沌吸引子，并且该混沌为多脉冲Shilnikov型轨道混沌.

由此可见，当ΔT=-40℃时，系统存在混沌吸引子，且伴随无序的大幅振动.此时动力系统的振动幅值难以准确预测，由此导致的大幅振动可能会影响结构的安全.当不考虑温度变化（0℃）或者变化升高（40℃）时，系统参数不再满足能量差函数的零解条件，Shilnikov型多脉冲轨道消失，混沌运动也随之消失，最大振幅明显下降.由此可见温度变化有可能明显改变动力系统的吸引子类型以及系统的振动幅值.

5 结论

温度变化会直接改变悬索静态构形，导致其固有频率及非线性系数发生改变，影响频率间的公倍关系，进而引发系统全局动力学行为产生变化.当悬索存在Shilnikov型多脉冲轨道并产生混沌运动时，温度变化会影响悬索系统中Shilnikov型多脉冲轨道的产生条件.在异宿轨道参数限制下，相同的激励条件，温度变化可能会引起产生混沌响应的多脉冲轨道消失，导致动力系统展现出截然不同的动力学行为，改变系统的共振响应幅值，影响结构安全.

附录

$\begin{array}{c}{K}_i=-\alpha {\int }_0^1{\mathrm{\Delta }}_x{\phi }_1^{\mathrm{\text{'}}\mathrm{\text{'}}}{\phi }_{i}dx,(i=\mathrm{1,\; 2})\\ {\omega }_i^2=-{\chi }^2{\int }_0^1{\phi }_i^{\mathrm{\text{'}}\mathrm{\text{'}}}{\phi }_{i}dx-\left(\alpha /{\chi }^4\right){\int }_0^1{y}^{\mathrm{\text{'}}\mathrm{\text{'}}}\cdot {\int }_0^1{y}^{\mathrm{\text{'}}}{\phi }_i^{\mathrm{\text{'}}}dx{\phi }_{i}dx,(i=\mathrm{1,\; 2})\\ P_i=-\left(\alpha /{\chi }^2\right)y^{\mathrm{\text{'}}\mathrm{\text{'}}}{\mathrm{\Delta }}_x{\int }_0^1{\phi }_{i}dx+{\mathrm{\Omega }}^2{\mathrm{\Delta }}_z{\int }_0^1(1-x){\phi }_{i}dx,(i=\mathrm{1,\; 2})\\ {\mathrm{\Lambda }}_{11}=-\left(\alpha /2{\chi }^2\right)\left({\int }_0^1{y}^{\mathrm{\text{'}}\mathrm{\text{'}}}{\int }_0^1{\phi }_1^{\mathrm{\text{'}}2}\mathrm{d}x{\phi }_1\mathrm{d}x-2{\int }_0^1{\phi }_1^{\mathrm{\text{'}}\mathrm{\text{'}}}{\int }_0^1{y}^{\mathrm{\text{'}}}{\phi }^{\mathrm{\text{'}}}{}_1\mathrm{d}x{\phi }_1\mathrm{d}x\right)\\ {\mathrm{\Lambda }}_{13}=-\left(\alpha /2{\chi }^2\right){\int }_0^1{y}^{\mathrm{\text{'}}\mathrm{\text{'}}}{\int }_0^1{\phi }_2^{\mathrm{\text{'}}2}dx{\phi }_{1}dx\\ {\mathrm{\Lambda }}_{22}=-\left(\alpha /{\chi }^2\right){\int }_0^1{\phi }_2^{\mathrm{\text{'}}\mathrm{\text{'}}}{\int }_0^1{y}^{\mathrm{\text{'}}}{\phi }_1^{\mathrm{\text{'}}}dx{\phi }_{2}dx\\ {\mathrm{\Gamma }}_{ij}=-(\alpha /2){\int }_0^1{\phi }_i^{\mathrm{\text{'}}\mathrm{\text{'}}}{\int }_0^1{\phi }_j^{\mathrm{\text{'}}2}dx{\phi }_{i}dx,(i,j=\mathrm{1,\; 2})\\ {\gamma }_1=-3{\mathrm{\Gamma }}_{11}/8{\omega }_1+5{\mathrm{\Lambda }}_{11}^2/12{\omega }_1^3\\ {\gamma }_2={\mathrm{\Lambda }}_{11}{\mathrm{\Lambda }}_{13}/2{\omega }_1^3-{\mathrm{\Gamma }}_{12}/4{\omega }_1-{\mathrm{\Lambda }}_{13}{\mathrm{\Lambda }}_{22}/\left(2{\omega }_1^3-8{\omega }_1{\omega }_2^2\right)\\ {\gamma }_3=-{\mathrm{\Gamma }}_{21}/4{\omega }_2+{\mathrm{\Lambda }}_{11}{\mathrm{\Lambda }}_{22}/4{\omega }_2{\omega }_1^2-{\mathrm{\Lambda }}_{22}^2/\left(4{\omega }_1^2{\omega }_2-16{\omega }_2^3\right)\\ {\gamma }_4=-3{\mathrm{\Gamma }}_{24}/8{\omega }_2-{\mathrm{\Lambda }}_{13}{\mathrm{\Lambda }}_{22}\left(8{\omega }_2^2-3{\omega }_1^2\right)/8{\omega }_1^2{\omega }_2\left({\omega }_1^2-4{\omega }_2^2\right)\end{array}$

参考文献