动叶可调增压风机转子的动态特性研究

图 1 增压风机三维实体模型 动叶可调增压风机转子的动态特性研究动叶可调增压风机转子的动态特性研究 成玫 (上海卫星工程研究所) 摘摘 要：要：应用有限元建立动叶可调增压风机转子的模型并研究其动力学特性。首先在 ABAQUS有限元软件平台下，建立了动叶可调增压风机转子的三维实体模型，通过对 VDLOAD 子程序的二次开发实现了非线性滚动轴承力载荷的有效模拟， 实现了动叶可调增压风机转子在滚动轴承支撑条件下的动力学三维仿真分析。其中，对于滚动轴承模型，充分考虑了轴承间隙、轴承滚珠与滚道的非线性赫兹接触力以及滚动轴承的实际接触角,并将改进后的滚动轴承模型与已有的结论作比较，达到了较好的一致性，说明用该模型来分析转子系统的动特性是可行的。继而以该模型为基础，在转子不平衡和滚动轴承力的情况下，研究转子系统的动力学特性。并观察轴心轨迹图、频谱图，发现了一些规律，为以后在现场中出现类似的故障信号作指导，节省维修时间和费用，提高工作效率。 关键词：关键词：转子，ABAQUS，滚动轴承，轴心轨迹 前言前言 有限元技术是工程技术领域进行科学计算的极为重要的方法之一， 利用有限元可以分析十分复杂的工程结构的机械性能， 也可以对各种工程事故进行技术分析， 它的实质就是把具有无限个自由度的连续系统， 理想化为只有有限个自由度的单元集合体， 使问题转化为适合于数值求解的结构型问题， 它的模拟是对真实情况的数值近似； 并且随着电子计算机技术和有限元数值仿真技术的发展， 对大型设备进行数值仿真的计算技术更加趋于完善， 计算结果的可信度也大大提高。现今的有限元分析软件很多，有 ABAQUS、ANSYS 等等，本文就利用ABAQUS 有限元软件对动叶可调增压风机的转子进行了动特性研究。 1.1. 动叶可调增压风机转子的计算模型动叶可调增压风机转子的计算模型 1.1 1.1 有限元模型的建立有限元模型的建立 对动叶可调增压风机转子采用了实体建模，模型中包括螺栓、倒角、轴肩、退刀槽等细节结构，以期计算模型尽可能地接近真实结构。在 ABAQUS 有限元软件中对于叶片这种曲面的绘制不理想，因此先利用 Catia 软件对增压风机进行了三维实体建模（如图 1） ，然后将其模型导入 ABAQUS软件中。 几何模型简化 图 2 简化后的三维实体模型 图 4 网格的划分 图 3 进一步简化后的 三维实体模型 数值仿真计算结果的正确与否以及精度的高低在很大程度上取决于有限元单元网格的质量。 本文需要计算的是动叶可调增压风机转子的振动特性， 所以简化的原则是不能改变结构的振动特性、刚度特性以及其他不影响计算的结构特性。在我们现阶段的研究中，并没有考虑由叶片产生的Alford力对转子轴系动力学特性的影响，而只是考虑了转子的偏心和轴承力这 2 个因素，同时这里的叶片是变截面弯扭叶片，难于划分规则的网格，从而会降低计算精度。在传统的计算中，一般会按照叶片和轮盘的质量及转动惯量，将其模化为一点，附加在转子模型上进行求解，但这样就无法反映出分布质量和分布转动惯量对临界转速的影响，以及由轮盘和叶片产生的离心力对轴的刚度的影响。大量的计算实例证明，有限元三维实体模型将更贴近实际1-3 有限元网格划分 。因此本文中没有将叶片作为附加质量建立在转子之上，而是采用了文献4中提出的等密度模化法， 忽略了实际的叶轮与转子之间的键联结关系， 并将增压风机转子与叶片构建为一个整体， 且保证模化前后的质量和转动惯量不变（如图 2） 。另外，对于动叶可调增压风机转子这种复杂的设备， 在进行有限元网格划分之前有必要对一些严重影响网格划分质量的部分（如倒角、轴肩、退刀槽等）进行适当的简化。在保证质量不变的情况下，本文将这些地方也进行了简化，改善了网格划分质量（如图 3） 。 有限元网格质量的好坏是进行有限元数值仿真的基础，对于动叶可调增压风机转子系统，本文选择了六面体缩减单元 C3D8R， 整个模型共有 2538 个单元，2992 个节点（如图 4） 。而且，由于缩减积分单元在显示动力学计算的过程中存在沙漏现象，即整个过程的动能和势能总和超过了总能量值，使得结果 精 度 变 差 ， 所 以 本 文 采 用 了“Second-accuracy”的沙漏控制方法，对于网格的划分及计算的选择，是需要针对各自特定的使用场合，通过不断的试用和比较，来选取最合适自己所分析的问题的方案。 1.2 1.2 边界条件处理边界条件处理 在本文中对于轴承处的边界条件处理， 分为了两种情况讨论， 分别为铰支和考虑滚动轴承力的支承。首先，由于六面体单元的节点只有 3 个自由度，没有旋转自由度，所以本文在轴上设置了参考点，通过对参考点施加恒定的转动速度来带动整个轴的旋转。 在铰支情况下我们采取的边界条件如下：全约束靠近涡轮盘一侧的与轴承配合的连接面，而另一侧与轴承配合的连接面处允许其轴向窜动。 对于轴承力， 综合考虑了轴承径向间隙、 非线性赫兹接触以及支承刚度的周期变化等因素后，建立的滚动轴承支承力5()()+=+=+=+=iNiiibNiiYbYiNiiibNiiXbXbbbbrYXCfFrYXCfFsinsincoscossincos123022112302212222如下： 式中：2bXF是2X方向的滚动轴承力分量，2bYF是2Y方向的滚动轴承力分量，bN为滚动体的个数，0r是滚动轴承的径向间隙，轴旋转的同时，角i只随时间而变化，均匀排列的 各 滚 珠 之 间 的 夹 角 为bN2， 在 轴 静 止 时 ， 第i个 滚 珠 的 角 位 置 为()()bbNiiN, 2, 112=，当轴以角速度旋转时，假设逆时针方向旋转为正，在时刻t角位置i的值为： ()tiNcagebi+=12 其中ioicageRRR+=，式中，oR、iR为外、内圈滚道半径；cage，分别为滚珠中心(即保持架中心)的角速度、转子的速度。 因此()()boiibiNitRRRiN, 2, 112=+= 滚动轴承的实际接触角与轴承间隙时有关的， 但在上式中并没有考虑滚动轴承的实际接触角，本文将滚动轴承的实际接触角也计入公式中。设滚动轴承的接触角为，则 )1arccos(0boidRRr= 因此轴承力可写为 ()()+=+=+=+=iNiiibNiiYbYiNiiibNiiXbXbbbbrYXCfFrYXCfFsincossincoscoscossincos123022112302212222 为了验证本文改进后的滚动轴承力是否正确， 与文献6中， Mevel建立的滚动轴承力进行了比较。转子的旋转速度为低转速 300 r/min时，可以观察到由于轴承内部刚度周期变化所引起的变柔度振动VC振动（varying compliance vibration） ，VC振动的频率为旋转频率的BN倍，BN是和轴承有关的系数，它取决于轴承尺寸。图 6 为采用本文模型计算得到的X、Y方向振动位移，从图 6 中可以清楚的看出，X、Y方向的运动是周期的，并十分明显的表 图 8 施加重力加速度 现出了滚珠的通过频率，当一个滚珠离开载荷区的时候，转子下降，并接触一个新的滚珠，使颤振迅速出现和消失。该计算结果与Mevel6的计算结果（如图 5）达到了较好的一致性，说明本文进行改进后的滚动轴承力的计算公式的正确性。同时，由图 7 的频谱上可以看出，轴承的振动表现为旋转频率的BN倍及其谐波。Fukata7的研究也表明，当转速在远离X和Y方向两个临界转速时，运动是周期的，表现出滚珠的转动频率和它的谐波，这再一次验证了本文改进后计算结果的正确性。 图 5 Mevel 计算结果绘制的转子响应 0100200300400500600700-0.4-0.3-0.2-0.100.10.20.3t /sAmplitude 0100200300400500600700-1.4-1.38-1.36-1.34-1.32-1.3-1.28t /sAmplitude 图 6 本文计算结果绘制的转子响应 0246810120123f/nPower Spectrum 图 7 功率谱图 1.3 1.3 载荷施加载荷施加 在铰支情况下，施加的载荷是：整个转子受重力的作用（如图 8） ，然后可以通过修改密度来达到实现偏心的目的。 对于滚动轴承力支承条件下，则将滚动轴承力施加在轴承处， 但 在 ABAQUS 中， 像这种包含位移和速度参数的载荷并不能直接加载，需要借助于ABAQUS 提供的二次开发平台 subroutine 来实现，具体来说就是选用其中的 VDLOAD 子程序来写入。 程序以 Fortran 软件为平台接入，VDLOAD 子程序的格式如下： subroutine vdload ( 图 9 ABAQUS 的滚动轴承力加载界面 C Read only (unmodifiable)variables - 1 nblock, ndim, stepTime, totalTime, 2 amplitude, curCoords, velocity, dirCos, jltyp, sname, C Write only (modifiable) variable - 1 value ) include vaba_param.inc dimension curCoords(nblock,ndim), velocity(nblock,ndim), 1 dirCos(nblock,ndim,ndim), value(nblock) character*80 sname do 100 k=1, nblock USER CODING TO DEFINE VALUEUSER CODING TO DEFINE VALUE 100 continue return end 其中 nblock， ndim， stepTime， totalTime， amplitude， curCoords， velocity， dirCos，jltyp，sname 为 subroutine 提供的软件接口，通过这些参数可以调用计算过程中模型某个或某些节点的坐标、速度等信息，还包含运算的时间步长、载荷的类型的信息。 在 ABAQUS 中，并不能对单个节点加载这种载荷，只能在杆单元、面单元和体单元上加载，因此我们借助杆单元来将滚动轴承力加载在转子与轴承接触的地方。具体的方法是在ABAQUS 中建立两段杆单元， 分别与转子 tie 在一起， tie 在一起的节点具有完全相同的位移及转速。将滚动轴承力加载在这两段杆单元上，如图 9。 2.2. 增压风机转子的动态响应分析增压风机转子的动态响应分析 2.1 2.1 转子不平衡转子不平衡 在 ABAQUS 中，通过修正一部分密度值，可以实现转子不平衡的设置。对于上述建立的增压风机转子系统的三维实体模型， 采用这个方法实现了转子的不平衡， 并进行了动特性分析，图 10-图 12 是在不同的转速下的 X 方向时域波形图、转子轨迹图及频谱图。综合比较，发现转子系统的振动信号的原始时间波形为简谐波， 其轴心轨迹呈圆形， 并且在振动信号的频谱图中，以转子系统的基频成分为主，在升降速的过程中，当转子系统的工作转速小于临界转速时，振动幅值随转速的增加而上升（如图 10-11） ，但在转子系统的工作转速大于临界转速之后，振动幅值随转速的增加法反而减小了（如图 11-12） 。 05010015020000.020.040.060.080.10.120.14f /HzAmplitude 图 10 转子转速为 30 rad/s 05010015020000.020.040.060.080.10.120.140.16f /HzAmplitude /m 图 11 转子转速为 350 rad/s 010020030040050060000.020.040.060.080.10.120.14f /HzAmplitude /m 图 12 转子转速为 400 rad/s 2.2 2.2 轴承间隙轴承间隙 轴承是机械系统中重要的支承部件， 其性能与工况的好坏将直接影响到与之相连的转子的性能，并且轴承的径向间隙是决定轴承-转子系统动态响应的一个重要参数，因此，本文对滚动轴承间隙的变化对转子系统动力学特性的影响也进行了研究。图 13 是在相同的转速下，不同轴承间隙的转子频谱图。通过比较，发现随着滚动轴承间隙的增加，频谱图中的振动幅值将逐渐增加，说明随着滚动轴承间隙的增加，滚动轴承的动态刚度减小了，从而导致转子系统的运动稳定性逐渐变差。 0501001502000 0.020.040.060.080.1 0.120.140.16f /HzAmplitude0501001502000 0.020.040.060.080.100.120.140.16f /HzAmplitude0501001502000 0.020.040.060.080.100.120.140.16f /HzAmplitude 图 13 转子在不同间隙下的频谱图 3.3. 结论结论 首先在 CATIA 软件中建立了动叶可调增压风机的转子系统的三维实体模型， 通过一定的简化原理， 对原三维实体模型进行简化并导入 ABAQUS 中， 然后通过对 ABAQUS 的加载单元进行了二次开发，利用 VDLOAD subroutine 加载了滚动轴承力，即而进行了在转子不平衡和不同的滚动轴承间隙下的有限元仿真计算， 最后对轴承力的改进以及对仿真计算结果进行分析得到如下结论： 1) 综合考虑了轴承径向间隙、 非线性赫兹接触、 支承刚度的周期变化以及滚动轴承的实际接触角等因素后， 建立的滚动轴承支承力与文献相比的结果较为符合， 说明改进后的滚动轴承力模型的有效性。 2) 在转子不平衡情况下， 转子系统的振动信号的原始时间波形为简谐波， 其轴心轨迹呈圆形，并且在振动信号的频谱图中，以转子系统的基频成分为主，在升降速的过程中，当转子系统的工作转速小于临界转速时， 振动幅值随转速的增加而上升， 但在转子系统的工作转速大于临界转速之后，振动幅值随转速的增加法反而减小。 3) 随着轴承间隙的增加，使得滚动轴承的动态刚度减小，转子系统的振动幅值增加，导致转子系统的运动稳定性逐渐变差。