一种基于阶梯轴-柔性盘耦合系统的动力学建模方法与流程

本发明涉及一种基于阶梯轴-柔性盘耦合系统的动力学建模方法，属于机械动力学技术领域。

背景技术：

目前，现有的轴盘耦合系统动力学建模方法主要有两种方法：基于商用有限元分析软件和基于半解析方法的动力学建模，其中，基于商用有限元分析软件是将cad三维模型导入商用有限元分析软件或者直接在有限元软件中建立三维模型，选择合适的单元及合适的材料参数，对三维模型进行网格划分，建立有限元模型，设置合适的约束并选择合适的求解方法对系统的动力学特性进行分析。但利用现有的商用有限元分析软件对阶梯轴-柔性盘的系统进行动力学特性分析时，建模过程复杂且繁重，计算效率低。

基于半解析方法的动力学建模是将轮盘和转轴动能与势能分别相加，利用hamilton能量原理和galerkin方法得到系统的动力学方程，求解系统的固有频率。但半解析方法一般无法考虑多段变截面的转轴而且无法考虑转轴多个方向振动，无法全面描述转轴的振动。

目前基于阶梯轴-柔性盘耦合系统进行动力学建模的技术处于空白状态。

技术实现要素：

(一)要解决的技术问题

为了解决现有技术的上述问题，本发明提供一种阶梯轴-柔性盘耦合系统的动力学建模方法，以达到在保证系统主要振动模态的前提下，考虑系统旋转效应的影响，利用lagrange方程和galerkin方法得到系统的动力学方程。

(二)技术方案

为了达到上述目的，本发明采用的主要技术方案包括：

一种阶梯轴-柔性盘耦合系统的动力学建模方法，其包括以下步骤：

s1、构建阶梯轴-柔性盘动力学建模所需的三维坐标系，所述三维坐标系包括整体坐标系oxyz和轮盘局部坐标系oxyz；

s2、对阶梯轴-柔性盘系统的结构参数和材料参数进行测定，其中包括轮盘的内径ri，单位为m，轮盘的外径ro，单位为m，轮盘厚度hd，单位为m，与转轴左端轴向距离为sd，单位为m，弹性模量ed，单位为pa，泊松比vd，轮盘密度ρd，单位为kg/m³；转轴在弯曲方向约束弹簧的刚度kx，单位为n/m、ky，单位为n/m，约束弹簧的刚度cx，单位为(n·s/m)、cy，单位为n·s/m，轴承与转轴左端轴向距离为s1，单位为m和s2，单位为m，第i段半径为ri，单位为m，长度为li，单位为m，阶梯轴弹性模量es，单位为pa，泊松比vs，密度ρs，单位为kg/m³；

s3、通过轮盘上任意一点在整体坐标系oxyz中的位移向量r对时间的一阶偏导，得到该点的速度，再依据动能计算公式得到轮盘的动能；

s4、采用有限元方法，对转子系统进行建模，得到转轴的动能；

s5、将所述转轴的动能和轮盘的动能相加，得到系统的动能表达式；

s6、考虑轮盘在旋转过程中的离心刚化效应和旋转软化效应，得出轮盘的势能；

s7、考虑转轴的扭转和弯曲，得到转轴的势能；

s8、计算轴承的势能；

s9、将所述转轴的势能、轮盘的势能和轴承的势能相加，得到系统的势能表达式；

s10、将系统的动能表达式和势能表达式带入lagrange方程，引入瑞利阻尼，得到系统的运动微分方程；

s11、利用galerkin方法对轮盘的运动方程进行离散化处理；

s12、将系统的运动微分方程写成矩阵的形式；

s13、在轮盘轴向以及转轴径向施加力，计算转轴在外力下的幅频响应。

如上所述的动力学建模方法，优选地，在步骤s3中，所述轮盘的动能计算公式如公式(1)所示：

其中，t1为轮盘自身动能；t2为轮盘与转子耦合产生的动能；xs、ys、zs、θxd、θyd和θzd分别为盘所在位置处轴的x、y和z方向的平动及绕x、y和z方向转动位移；ud为圆盘上任意一点p在垂直于圆盘平面的位移；r和θ为点p在极坐标系下的坐标；ω为轮盘旋转角速度；符号(·)表示对时间的1阶偏导。

如上所述的动力学建模方法，优选地，步骤s4中，所述转轴的动能为如公式(2)所示：

其中，为转轴的速度向量，ms为转轴的质量矩阵，转轴采用两节点梁单元，每个节点有6个自由度，分别是x、y和z方向的平动及绕x、y和z轴方向的转动θx、θy和θz；在坐标系oxyz中，xa、ya、za、xb、yb和zb分别为a、b节点的x、y和z方向的位移，θax、θay、θaz、θbx、θby和θbz分别为a、b节点的x、y和z方向的转角；

梁单元的位移向量如公式(3)所示：

u^e＝[xayazaθaxθayθazxbybzbθbxθbyθbz]^t(3)

式中，上标e代表单元。

如上所述的动力学建模方法，优选地，步骤s6中，所述轮盘的势能如公式(4)所示：

其中，是拉普拉斯算子，表达式如公式(5)所示；为圆盘的抗弯刚度；nr和nθ为极坐标上正交合应力，表达式如公式(6)所示：

其中，

如上所述的动力学建模方法，优选地，在步骤s7中，所述转轴的势能的表达式如公式(7)所示：

如上所述的动力学建模方法，优选地，在步骤s8中，所述轴承的势能的表达式如公式(8)所示：

式中，kb为轴承刚度矩阵，kb＝diag[0…kb1…0…kb2…0]，kb1和kb2分别为轴承1和2的单元矩阵；

轴承采用线性弹簧-阻尼模型，轴承的单元刚度矩阵及单元阻尼矩阵的具体表达式如公式(9)所示：

如上所述的动力学建模方法，优选地，在步骤s10中，所述lagrange方程如公式(10)所示：

式中，l表达式为：l＝ttotal-vtotal，ttotal为系统的动能，vtotal为系统的势能。

如上所述的动力学建模方法，优选地，步骤s11中，所述galerkin方法的离散化过程引入正则坐标qξ(t)和qη(t)，将轮盘横向位移如公式(11)所示：

ud＝cosmθφd(r)qξ(t)+sinmθφd(r)qη(t)(11)

式中：m表示圆盘的节径数；φd为弹性盘假设模态组成的列向量，表达式如公式(12)所示：

其中，r1i为弹性盘的横向变形函数，(i＝1,2,…,nd)；nd为弹性盘的模态截断数，r1i的表达式如公式(13)所示：

r1i＝[sinβ(r-ri)-sinhβ(r-ri)]+α[cosβ(r-ri)-coshβ(r-ri)](13)

式中，cosβ(r-ri)chβ(r-ri)＝-1，

如上所述的动力学建模方法，优选地，在步骤s12中，所述阶梯轴-柔性盘耦合系统的运动微分方程如公式(14)所示：

式中，m为整个耦合系统的质量矩阵；g为整个耦合系统的阻尼矩阵，包括轴承的阻尼矩阵、圆盘和转轴的陀螺矩阵；d为整个耦合系统的黏性阻尼矩阵，采用瑞利阻尼；k为整个耦合系统的刚度矩阵；f为整个耦合系统的外部激振力向量。

其中，所述瑞利阻尼的表达式如公式(15)所示，

d＝αm+βk(15)

fn1和fn2表示系统的前两阶固有频率，ξ1和ξ2是对应的模态阻尼比。

如上所述的动力学建模方法，优选地，步骤s13中，所述系统的幅频响应计算过程为：首先给定正弦激励幅值，然后假设系统动力学方程的解为q＝q0e^iωt，并将其带入公式(14)中，由此可得系统动力学方程的解如公式(16)所示：

q0＝(-ω²m+iωc+k)^-1f(16)。

(三)有益效果

本发明的有益效果是：

本发明的方法节省了阶梯轴-柔性盘系统动力学实验所需的成本费用；只需修改结构尺寸和材料参数即可得到不同轴盘耦合系统的动力学模型，操作简便；其考虑了轮盘柔性和轴盘耦合效应，这更加接近轴盘真实的耦合系统；与借助传统的商用有限元软件来分析阶梯轴-柔性盘系统的动力学特性相比，本发明具有更高的计算效率；同时，本发明还能进行阶梯轴-柔性盘系统的碰摩响应分析，从而为含轴盘耦合的系统结构提供设计优化。

附图说明

图1为本发明中阶梯轴-柔性盘耦合系统动力学建模的流程图；

图2为阶梯轴-柔性盘耦合系统结构示意图；

图3为质量矩阵组集示意图；

图4为阶梯轴-柔性盘耦合系统的有限元建模示意图；

图5为阶梯轴-柔性盘系统实验台示意图；

图6为阶梯轴-柔性盘耦合系统与ansys振型对比；

图7为阶梯轴-柔性盘耦合系统与实验振型对比；

图8为阶梯轴-柔性盘耦合系统动频对比；

图9为阶梯轴-柔性盘耦合系统的幅频响应对比。

具体实施方式

为了更好的解释本发明，以便于理解，下面结合附图，通过具体实施方式，对本发明作详细描述。

实施例1

一种基于阶梯轴-柔性盘耦合系统动力学建模方法，如图1所示，包括如下步骤：

步骤1：构建阶梯轴-柔性盘动力学建模所需的三维坐标系，包括：整体坐标系oxyz，轮盘局部坐标系oxyz；

步骤2：对阶梯轴-柔性盘系统的结构参数和材料参数进行测定，其中包括轮盘的内径ri/m，轮盘的外径ro/m，轮盘厚度hd/m，与转轴左端轴向距离为sd/m，弹性模量ed/pa，泊松比vd，轮盘密度ρd/(kg/m³)；转轴在弯曲方向约束弹簧的刚度kx/(n/m)、ky/(n/m)，约束弹簧的刚度cx/(n·s/m)、cy/(n·s/m)，轴承与转轴左端轴向距离为s1/m和s2/m，第i段半径为ri/m，长度为li/m，阶梯轴弹性模量es/pa，泊松比vs，密度ρs/(kg/m³)；

步骤3：通过轮盘上任意一点在整体坐标系oxyz中的位移向量r对时间的一阶偏导，得到该点的速度，再依据动能计算公式得到轮盘的动能如公式(1)所示；

所述轮盘的动能计算公式如公式(1)所示：

其中，t1为轮盘自身动能；t2为轮盘与转子耦合产生的动能。xs、ys、zs、θxd、θyd和θzd分别为盘所在位置处轴的x、y和z方向的平动及绕x、y和z方向转动位移；ud为圆盘上任意一点p在垂直于圆盘平面的位移；r和θ为点p在极坐标系下的坐标；ω为轮盘旋转角速度；符号(·)表示对时间的1阶偏导；

步骤4：采用有限元方法，对转子系统进行建模，得到转轴的动能，得出转轴动能表达式如公式(2)所示：

其中，为转轴的速度向量，ms为转轴的质量矩阵，转轴采用两节点梁单元，每个节点有6个自由度，分别是x、y和z方向的平动及绕x、y和z轴方向的转动θx、θy和θz；在坐标系oxyz中，xa、ya、za、xb、yb和zb分别为a、b节点的x、y和z方向的位移，θax、θay、θaz、θbx、θby和θbz分别为a、b节点的x、y和z方向的转角；

梁单元的位移向量如公式(3)所示：

u^e＝[xayazaθaxθayθazxbybzbθbxθbyθbz]^t(3)

式中，上标e代表单元。

步骤5：考虑轮盘在旋转过程中的离心刚化效应和旋转软化效应，得出轮盘的整体势能，所述轮盘势能如公式(4)所示

其中，是拉普拉斯算子，表达式如公式(5)所示；为圆盘的抗弯刚度；ed和vd为圆盘的弹性模量和泊松比；nr和nθ为极坐标上正交合应力，表达式如公式(6)所示：

其中，

步骤6：考虑转轴的扭转和弯曲，得到转轴的势能，所述转轴势能的表达式如公式(7)

步骤7：轴承势能的表达式如公式(8)：

式中，kb为轴承刚度矩阵，kb＝diag[0…kb1…0…kb2…0]，kb1和kb2分别为轴承1和2的单元矩阵。

轴承采用线性弹簧-阻尼模型，轴承的单元刚度矩阵及单元阻尼矩阵的具体表达式如公式(9)：

步骤8：将系统的动能表达式和势能表达式带入lagrange方程，引入瑞利阻尼，得到系统的运动微分方程，所述lagrange方程为公式(10)：

式中，l表达式为：l＝ttotal-vtotal，ttotal为系统的动能，表达式为：

vtotal为系统的势能，表达式为：

nd为转轴上轮盘数目。

步骤9：利用galerkin方法对轮盘的运动方程进行离散化处理，所述galerkin方法的离散化过程为引入正则坐标qξ(t)和qη(t)，将轮盘横向位移写为：

ud＝cosmθφd(r)qξ(t)+sinmθφd(r)qη(t)(11)

式中：m表示圆盘的节径数；φd为弹性盘假设模态组成的列向量，表达式为：

式中：r1i为弹性盘的横向变形函数，(i＝1,2,…,nd)；nd为弹性盘的模态截断数。r1i的表达式为：

r1i＝[sinβ(r-ri)-sinhβ(r-ri)]+α[cosβ(r-ri)-coshβ(r-ri)](13)

式中，cosβ(r-ri)chβ(r-ri)＝-1，

步骤10：引入瑞利阻尼，得到阶梯轴-柔性盘耦合系统的运动微分方程：

式中，m为整个耦合系统的质量矩阵；g为整个耦合系统的阻尼矩阵，包括轴承的阻尼矩阵、圆盘和转轴的陀螺矩阵；d为整个耦合系统的黏性阻尼矩阵，采用瑞利阻尼；k为整个耦合系统的刚度矩阵；f为整个耦合系统的外部激振力向量。其中瑞利阻尼矩阵的表达式如下：

d＝αm+βk(15)

其中，fn1和fn2表示系统的前两阶固有频率，ξ1和ξ2是对应的模态阻尼比。

步骤11：设置外激励向量为零，确定系统不同转速下的固有频率。

步骤12：给定外加激励，计算系统在不同外加激励位置下的幅频响应。

实施例2

本实施例中阶梯轴-柔性盘耦合系统结构示意图如图2所示，其中，(a)为阶梯轴-盘结构示意图，(b)为阶梯轴结构尺寸图，(c)为轮盘结构尺寸图，单位为mm，动力学建模方法，包括以下步骤：

步骤1：获取阶梯轴-柔性盘耦合系统结构参数和材料参数，本发明假定轮盘和转轴是各向同性的线弹性材料，本构关系满足hooke定律，阶梯轴-柔性盘耦合系统的结构参数如图2所示，阶梯轴-柔性盘耦合系统的材料参数如图1所示：

表1阶梯轴-柔性盘耦合系统系统材料参数

步骤2：确定轮盘上任意一点在整体坐标系oxyz中的位移向量，并由此计算得到轮盘的动能表达式：

其中，t1为轮盘自身动能；t2为轮盘与转子耦合产生的动能；xs、ys、zs、θxd、θyd和θzd分别为盘所在位置处轴的x、y和z方向的平动及绕x、y和z方向转动位移；ud为圆盘上任意一点p在垂直于圆盘平面的位移；r和θ为点p在极坐标系下的坐标；ω为轮盘旋转角速度；符号(·)表示对时间的1阶偏导。

步骤3：采用有限元方法，对转子系统进行建模，得到转轴的动能，得出转轴动能表达式如公式(2)所示：

其中，为转轴的速度向量，ms为转轴的质量矩阵，转轴采用两节点梁单元，每个节点有6个自由度，分别是x、y和z方向的平动及绕x、y和z轴方向的转动θx、θy和θz。在坐标系oxyz中，xa、ya、za、xb、yb和zb分别为a、b节点的x、y和z方向的位移，θax、θay、θaz、θbx、θby和θbz分别为a、b节点的x、y和z方向的转角。

梁单元的位移向量为：

u^e＝[xayazaθaxθayθazxbybzbθbxθbyθbz]^t(3)

式中，上标e代表单元。

步骤4：考虑轮盘在旋转过程中的离心刚化效应和旋转软化效应，得出轮盘的整体势能，轮盘势能如公式(4)所示

其中，是拉普拉斯算子，表达式如公式(5)所示；为圆盘的抗弯刚度；ed和vd为圆盘的弹性模量和泊松比；nr和nθ为极坐标上正交合应力，表达式如公式(6)所示：

其中，

步骤5：考虑转轴的扭转和弯曲，得到转轴的势能，所述转轴势能的表达式为公式(7)所示：

步骤6：轴承势能的表达式为：

式中，kb为轴承刚度矩阵，kb＝diag[0…kb1…0…kb2…0]，kb1和kb2分别为轴承1和2的单元矩阵。

轴承采用线性弹簧-阻尼模型，轴承的单元刚度矩阵及单元阻尼矩阵的具体表达式如下：

轴承刚度kxb1、kyb1、kxb2和kyb2均为5.3×10⁸n/m，kzb1、kzb2为0n/m；轴承阻尼cxb1、cyb1、cxb2和cyb2均为2.1×10³n·s/m，czb1、czb2为0n·s/m。

步骤7：将系统的动能表达式和势能表达式带入lagrange方程，引入瑞利阻尼，得到系统的运动微分方程，所述lagrange方程为：

式中，l表达式为：l＝ttotal-vtotal，ttotal为系统的动能，表达式为：

vtotal为系统的势能，表达式为：

nd为转轴上轮盘数目。

步骤8：利用galerkin方法对轮盘的运动方程进行离散化处理，所述galerkin方法的离散化过程为引入正则坐标qξ(t)和qη(t)，将轮盘横向位移写为：

ud＝cosmθφd(r)qξ(t)+sinmθφd(r)qη(t)(11)

式中：m表示圆盘的节径数；φd为弹性盘假设模态组成的列向量，表达式为：

式中：r1i为弹性盘的横向变形函数，(i＝1,2,…,nd)；nd为弹性盘的模态截断数。r1i的表达式为：

r1i＝[sinβ(r-ri)-sinhβ(r-ri)]+α[cosβ(r-ri)-coshβ(r-ri)](13)

式中，cosβ(r-ri)chβ(r-ri)＝-1，

步骤9：引入瑞利阻尼，得到阶梯轴-柔性盘耦合系统的运动微分方程：

式中，m为整个耦合系统的质量矩阵；g为整个耦合系统的阻尼矩阵，包括轴承的阻尼矩阵、圆盘和转轴的陀螺矩阵；d为整个耦合系统的黏性阻尼矩阵，采用瑞利阻尼；k为整个耦合系统的刚度矩阵；f为整个耦合系统的外部激振力向量。

其中瑞利阻尼矩阵的表达式如下：

d＝αm+βk(15)

其中，fn1和fn2表示系统的前两阶固有频率，ξ1和ξ2是对应的模态阻尼比。

f为系统外激力向量，若外部激振力施加在轮盘上，系统外部激振力向量f为：

f＝f·rq

式中：f为外部激振力幅值。

若外部激振力施加在转轴弯曲方向上，系统外部激振力向量f为

f＝f·xs

步骤10：取galerkin方法的截断阶数为3以得到轮盘具体的质量矩阵、阻尼矩阵和刚度矩阵。在阶梯轴-柔性盘耦合系统的质量矩阵(矩阵组集示意图如图3所示)、阻尼矩阵和刚度矩阵(矩阵组集与图3类似，此处不再赘述)的值如下：

(1)转轴单元矩阵：

单元质量矩阵m^e：

式中：

单元刚度矩阵k^e：

式中：

单元陀螺矩阵g^e：

式中：es为材料杨氏模量；gs为剪切模量；ρs为材料密度；υs为材料泊松比；as为单元横截面积；

l为单元长度；js为极惯性矩；i为截面惯性矩；ω为旋转角速度；

另外，实心梁单元的剪切系数：

轮盘刚性部分采用集中质量模型，轮盘的单元质量矩阵及单元阻尼矩阵的具体表达式如下：

(2)轮盘矩阵：

(3)轴盘耦合矩阵：

步骤11：计算特征方程系数行列式的特征值λ，取其虚部的绝对值除以2π，并进行从小到大排序，获得一组固有频率ωk，其中，k表示系统模态的第k阶，k＝1,2,…。

步骤12：给定正弦激励幅值计算系统的幅频响应。假设系统动力学方程的解为q＝q0e^iωt，并将其带入式(14)中，由此可得系统动力学方程的解为：

q0＝(-ω²m+iωc+k)^-1f(16)

步骤13：通过有限元软件建立相同的有限元模型(如图4所示)，利用现有技术中有限元软件中的模态求解和谐响应求解功能计算系统的固有频率如表1、图6和图8所示，和谐响应如图9所示，其中图6中a代表为本发明方法，b代表有限元方法；图8为在刚性盘和柔性盘假设下系统固有频率随转速的变化，图8中为阶梯轴刚性盘采用本发明方法，为阶梯轴柔性盘采用本发明方法，为阶梯轴柔性盘采用有限元方法，为阶梯轴刚性盘采用有限元方法；图9为分别在转轴和轮盘处施加力计算系统的幅频响应，(a)为在转轴上施加力，(b)为在轮盘上施加力；图9中为阶梯轴刚性盘采用本发明方法，为阶梯轴柔性盘采用本发明方法，为阶梯轴柔性盘采用有限元方法，为阶梯轴刚性盘采用有限元方法。将有限元方法所获得的数据与通过本发明获得的数据进行对比，固有频率曲线与谐响应曲线均能很好的吻合，从而证明了本发明所提方法的正确性。

表1阶梯轴-柔性盘系统固有频率对比

步骤14：通过搭建相同尺寸的阶梯轴-柔性盘实验台，如图5(a)所示，其中1为电机，2为转轮盘1，3为转轮盘1，4为阶梯轴，5为轴承；(b1)为阶梯轴上锤击点以及采集点示意图，圆圈位置为信号采集位置，其余图上黑点为锤击点，锤击方向沿箭头方向，(b2)为盘1上锤击点以及采集点示意图，(b3)为盘2上锤击点以及采集点示意图；(c)为dh5956数据采集前端控制器；(d)为pcb086c01力锤；(e1)为磁块传感器，用于检测阶梯轴上的振动信号，(e2)和(e3)为磁块传感器，型号为ca-yd-182，用于检测轮盘上的振动信号，采用的型号分别为ca-yd-12573708和ca-yd-12573711。采用现有技术锤击法，测得静止状态下固有频率和振型(如表1和图7所示)，与本发明方法对比，吻合很好，再一次验证本发明提供方法的正确性。

步骤15：采用本发明方法和有限元方法建立相同尺寸的阶梯轴-刚性盘系统，求解固有特性和幅频响应，对比柔性轮盘的系统得出，轮盘柔性和轴盘耦合效应对于系统高阶振动影响很大。

通过本实施例的结果可以获得以下结论：

(1)轴盘耦合效应和轮盘柔性在高频振动下对系统固有特性影响很大，而在低频振动下，轴盘耦合效应和轮盘柔性对固有特性的影响不再明显；

(2)本发明所提方法与有限元软件ansys所得动频曲线和谐响应曲线均吻合较好，从而验证本发明所建立的动力学模型的正确性。其中，在计算固有频率过程中，本发明所提方法耗时0.28s，有限元软件耗时8.4s；在计算谐响应的过程中，本发明所提方法耗时11.8s，有限元软件耗时160.5s。说明，采用本发明所提方法相对于利用有限元软件能明显的提高计算效率。

以上所述，仅是本发明的较佳实施例而已，并非是对本发明做其它形式的限制，任何本领域技术人员可以利用上述公开的技术内容加以变更或改型为等同变化的等效实施例。但是凡是未脱离本发明技术方案内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与改型，仍属于本发明技术方案的保护范围。