[发明专利]用于机械非线性振动分析的方程求解方法

权利要求书

1.一种用于机械非线性振动分析的非线性动力学方程求解方法，其特征在于，包括以下步骤：

步骤1，根据待分析的机械非线性振动系统，构建获得系统时域上非线性有限元动力学方程，表达式为：

式(1)中，k、c、m分别表示系统的刚度、粘滞阻尼及质量矩阵；x表示系统内所有自由度的位移向量；f_n(x)表示非线性力向量；f_l(t)表示外部激励向量；

步骤2，通过多谐波平衡法，将式(1)中的x以及f_l(t)进行广义傅里叶变换，获得频域上非线性有限元动力学方程；

步骤3，利用自动微分方法求解步骤2获得的频域上非线性有限元动力学方程的雅可比矩阵；

步骤4，根据步骤3获得的雅可比矩阵，采用弧长法求解步骤2获得的频域上非线性有限元动力学方程，获得频域上非线性有限元动力学方程的解；

步骤5，根据步骤4获得的频域上非线性有限元动力学方程的解，通过广义傅里叶逆变换，获得时域上非线性有限元动力学方程的解；根据时域上非线性有限元动力学方程的解，获得频响曲线，进行机械非线性振动分析；

其中，步骤2具体包括：

步骤2.1，将式(1)中的x进行傅里叶展开，获得的表达式为：

式(2)中：

D＝[I,Icos(ωt),Isin(ωt),...,Icos(N_kωt),Isin(N_kωt)]；

ω为激振频率；

N_k为选定谐波分量数；

步骤2.2，将式(1)中的f_l(t)进行傅里叶展开，获得表达式为：

式(3)中：

D＝[I,Icos(ωt),Isin(ωt),...,Icos(N_kωt),Isin(N_kωt)]；

ω为激振频率；

N_k为选定谐波分量数；

步骤2.3，将傅里叶展开后的x以及f_l(t)代回方程式(1)，获得频域上非线性有限元动力学方程表达式为：

R(X)＝Z(ω)X+F_n(X)-F_l＝0 (4)

式(4)中，X为各自由度方向位移谐波分量幅值矢量；Z(ω)为频域上的动态刚度矩阵；F_n(X)为频域上的非线性力向量；F_l为频域上的激振力谐波分量幅值矢量；

其中，Z(ω)表示为：

步骤3具体包括：

应用牛顿法求解式(4)；迭代过程表示为：

式中，k表示当前迭代步数；

将式(5)与式(4)联立获得表达式为：

式(6)中，采用自动微分方法求解具体步骤包括：

对各自由度方向位移谐波分量幅值矢量X做逆向傅里叶变换获得x，表达式为

根据位移-非线性力关系，获得非线性力向量f_n(x)；对f_n(x)做傅里叶变换得到频域上的非线性力向量F_n(X)，

获得X→x→f_n(x)→F_n(X)的计算图；

通过正向遍历计算图，求出每个节点的值；

通过反向遍历计算图，计算出每个节点的偏导，获得

步骤4具体包括以下步骤：

步骤4.1，在式(4)中加入λ，得到含λ的方程组表达式为：

Ψ(X,λ)＝Z(ω)X+F_n(X)-λF_l＝0； (7)

构造一个N+1维向量：ν(X,λ)＝[ν₁,ν₂,…,ν_N,ν_N+1]^T；

其中，式中，符号∧表示删除该列；

ν(X,λ)是N+1维空间中解曲线X(s)和λ(s)的切向向量；

解曲线的单位切向向量表达式为：式中，||ν||代表ν的模；

对第m个弧长步Δs，约束方程表达式为：

式(7)和式(8)结合，获得表达式为：

式(9)的雅可比矩阵表达式为：

步骤4.2，第一次迭代时，由牛顿法获得式(7)在弱非线性区上的一个初始解X^*和λ^*；

步骤4.3，预估步，其中，下标m表示第m个弧长步；

步骤4.4，校正步，迭代过程表达式为式(11)；在每一迭代步中，将近似解增量投影至τ_m方向得到弧长增量；