采用电磁轴承控制柔性转子临界转速分布的建模和仿真分析

方 鹏，蒋科坚

（浙江理工大学信息学院，杭州310018）

采用电磁轴承控制柔性转子临界转速分布的建模和仿真分析

方 鹏，蒋科坚

（浙江理工大学信息学院，杭州310018）

提出了采用电磁轴承支承柔性转子，通过支承特性的调节，改变柔性转子的各阶临界转速的分布位置，使柔性转子顺利实现超临界运行。将电磁轴承支承特性建模融入到经典转子动力学的柔性转子有限元建模理论中去，构建了电磁轴承支承的柔性转子系统模型进行仿真。仿真结果表明：通过调节电磁轴承的等效刚度，可以明显改变柔性转子平动和锥动两个刚体临界转速，而对一阶弯曲和二阶弯曲两个弯曲临界转速影响不大。调节等效阻尼可明显减小转子过临界时的振动。

电磁轴承；柔性转子；磁悬浮；有限元；临界转速

0 引 言

旋转机械是工程中常见的一类机械结构，在石油、化工、能源、电力、航空航天等诸多领域都是至关重要的关键组成设备。高速、高旋转精度和轻量化是旋转机械设计追求的目标，致使转子系统由刚性转子发展为柔性转子设计。在转子动力学中，一般认为在考虑转子动力学设计时，如果转子旋转而造成的自身形变不能被忽略，该转子可认为是柔性转子。另一种定义，当转子的运行速度超过其一阶弯曲临界转速的70%，可认为是柔性转子。由于柔性转子的自身强度相对较低，运行时的形变和振动十分突出。

主动电磁轴承（active magnetic bearing，AMB，下文称电磁轴承）用电磁力把转子悬浮在固定位置，以得到与传统机械轴承相类似的支承效果。电磁轴承不仅能提供转子无接触的悬浮支承力，更重要的是，电磁轴承可通过改变控制系统参数来调节其支承特性，从而实现对转子振动的主动控制。在目前的转子动力学研究中，有关机械轴承弹性支承的柔性转子系统的研究已经颇具规模，基本上形成了完整的理论体系和分析方法。电磁轴承作为一种新的轴承技术，为转子支承和改善转子振动提出了新的研究内容。

在经典转子动力学中，对柔性转子系统动力学特性的分析方法主要有传递矩阵法和有限元法。传递矩阵法中阶数不随转子系统自由度增大而增大，因而编程简单，运算速度快，但传递矩阵法不易解决转子系统的非线性问题。有限元法是被公认为比传递矩阵法更为精确的建模方法，但是，其计算矩阵的阶数按转子节点数量的增加而剧烈增长，计算复杂度较大。随着计算机技术的迅猛发展，有限元法的计算瓶颈得到很大缓解，在近几十年来已经成为柔性转子系统动力学常用的分析工具。柔性转子在不同转速下表现出不同的模态振型，各个模态振型的转速对应于柔性转子的各阶临界转速。转子运行在临界转速时，会产生较大的转子共振，严重影响了转子高速运行的稳定性，甚至造成事故。因此，如何使转子安全地跨越临界以及超临界运行，是柔性转子动力学研究的主要问题之一。

柔性转子作为转子动力学的重要内容已有较长的研究历史，国内外文献较多。但是，有关电磁轴承支承的柔性转子研究相对较少［1］。Arias等［2］为双盘柔性转子系统建立了有限元模型，每个单元视为欧拉梁结构，进行不平衡振动的计算和控制。Tseng等［3］提出以影响系数法计算不平衡量，以卡尔曼滤波器实现测量噪声消除的转子实时动平衡方案。李红伟等［4］采用影响系数法和振型平衡法分别对转子的刚性模态和前两阶挠性模态进行了本机动平衡。Kang等［5］采用基于影响系数法的自动平衡方法，其中不平衡质量修正用神经PD控制。Dyer［6］指出用影响系数法对柔性转子振动进行控制，不可能通过在转子上有限数量的位置（节点）施加电磁力的方法彻底消除整个转子（连续的无限数量的节点）的振动。汪希平等［7］和周朝暾等［8］通过刚性转子实验发现利用电磁轴承刚度变化可以敏感地调节转子的刚体临界转速。谢振宇等［9］通过柔性转子实验发现转子的一阶弯曲固有频率值受电磁轴承控制参数变化影响很小。要使转子的一阶弯曲固有频率避开工作转速，只能通过改变转子的本身尺寸来实现。Xu（徐旸）等［10］开展了基于电磁轴承的柔性转子过弯曲临界的实验，完成了10 MW高温气冷堆发电机磁轴承实验台架的设计工作，转子重3.5 t，长4.6 m。

本文提出采用电磁轴承支承柔性转子，通过调节支承的刚度阻尼，改变柔性转子各阶临界转速的分布位置，从而使转子的运行转速始终避开临界转速，最终达到抑制振动的目的。为此，笔者将经典转子动力学［11-12］中柔性转子有限元建模理论与电磁轴承支承特性建模相结合，构建了电磁轴承支承的柔性转子系统模型。通过仿真，分析电磁轴承等效刚度和等效阻尼的调整对柔性转子前四阶临界转速变化的影响规律。

1 电磁轴承电磁力建模

电磁轴承是用可控的电磁力把转子悬浮在固定位置，如图1所示。转子位移分别由两个垂直方向的位移传感器检测。如果转子在受到扰动后偏离其期望位置，此时由位移传感器检测出转子的偏移位移，控制器根据相应的控制策略产生控制信号，然后经功率放大器产生驱动定子绕组的控制电流。电磁轴承定子绕组可以产生主动电磁力，使得转子返回到其原来位置，这样便可以使转子动态地稳定悬浮在期望位置了。

图1 电磁轴承原理

电磁轴承的电磁力经线性化后可描述为位移刚度系数和电流刚度系数的线性函数［7］，如下：

式（1）中，ki，kx是电磁轴承电磁力线性化后的电流刚度系数和位移刚度系数。μ0=4π×10-7（H/m）为真空磁导率，A为磁路有效横截面积，N为线圈匝数，i0为线圈偏置电流，Δi为控制电流，x0为电磁轴承的设计气隙，Δx为转子在坐标方向的偏移位移。

电磁轴承的主动控制能力体现在电磁力可根据转子位移情况实时的调节，其控制系统传递函数的频域方程形式为：

式（3）中，Re｛G（jω）｝和Im｛G（jω）｝分别表示控制系统传递函数G（jω）的实部和虚部。可见，只有在可以获得电磁轴承控制传递函数的前提下，才能理论计算电磁轴承的等效刚度和等效阻尼。如采用PID控制的电磁轴承，其理想的传递函数为Kp+Ki/s+ Kds，其中Kp、Ki及Kd分别是控制器的P、I及D参数，那么其等效刚度和等效阻尼可表示为：

式（2）中，G（jω）为控制系统整体传递函数的频域关系，一般包括位移传感器、AD变换、控制器算法、DA变换、功率放大器等各控制环节传递函数的总和。

如果把电磁轴承的支承特性等效为一个机械轴承的刚度和阻尼，称之为等效刚度 Ke和等效阻尼［7］，其表达式分别为：当然，式（4）只是理想的PID控制，若考虑时滞，滤波降噪等实际情况，那么表达式要复杂得多。

2 电磁轴承支承的柔性转子有限元建模

2.1 柔性转子有限元建模简介

实际的转子系统通常是由一些轴段、圆盘和轴承座等部件组成，是一个质量连续分布的弹性系统，具有无穷多个自由度。有限元法的基本思想是将整个转子按节点划分为有限个单元，各单元彼此在节点处联结，从而使一个质量连续分布转子的运动，离散为有限个自由度的运动。转子划分的单元数目越多，其模型越精确，但计算量也大为增加。

因为柔性转子有限元建模作为转子动力学分析的常用工具，在大多数转子动力学教材中［11-12］都有十分详细的叙述。本文限于篇幅，仅对其建模思想做简要叙述，以不破坏本文阅读的可理解性。本文只对电磁轴承对柔性转子的支承建模做详细论述。

柔性转子相邻两个节点之间的轴段作为有限元分析的基本单元，称之为弹性轴段单元，如图2，其中左端面为第n节点，右端面为第n+1节点。每个轴段的空间坐标可以用两个端面的位移向量表示，称为轴段的广义坐标。端面的位移向量可由端面轴心坐标x及y和偏转角θy及θx来表示，x和y方向分开表示，

图2 弹性轴段单元

有限元法的建模过程，先建立每个轴段单元的动力学方程。在不考虑剪切和扭转变形影响的条件下，一个弹性轴段单元在x和y两个方向上的运动微分方程分别为：

式（6）中，ω为转子转速，｛fx｝及｛fy｝分别为弹性轴段两个端面节点处在x及y两个方向上所受的作用力和力矩，如不平衡力等外力，称为广义外力矢量。［mx］、［dx］、［gx］、［kx］和［my］、［dy］、［gy］、［ky］分别为转子在x和y两个方向上的质量矩阵、阻尼矩阵、陀螺矩阵和刚度矩阵。对于单个轴段在一个坐标方向内的运动方程而言，这些矩阵都是4×4的矩阵，其元素由轴段的线密度、尺寸、弯曲刚度等参数决定。当转子为轴对称时，x和y方向的矩阵相同。具体表达式由于篇幅较大，本文不详细列出，可参见有关转子动力学教材［11］。

根据有限元法，在列出每个轴段的运动方程后，依据相邻轴段间的连接关系，各个轴段方程可整合为一个转子整体的大矩阵运动方程，有：

式（7）中，｛F｝为广义外力矢量，［M］、［K］、［D］和［G］分别为整合后的转子的质量矩阵、刚度矩阵、阻尼矩阵和陀螺矩阵，形式如下：

式（8）至式（10）中，小方块表示对应单个轴段的4×4的［m］、［g］及［k］矩阵。每个轴段的［m］、［g］及［k］矩阵按图所示以对角线排列，并且在相邻矩阵的2×2部分重叠相加。［M］、［K］和［G］矩阵的上半部分为x平面的运动方程，下半部分为y平面的运动方程。可见，［M］、［K］和［G］矩阵是一个4N× 4N矩阵，即有限元法中的矩阵阶次将以节点数4倍的2次方量级增长。由于电磁轴承的转子悬浮无接触，内阻尼为0，因此阻尼矩阵［D］为零矩阵。电磁轴承的等效阻尼以外阻尼的形式体现在电磁力中，下文将详细叙述。在只有弹性轴段的情况下，［M］、［K］、［D］和［G］矩阵都是4N×4N阶的对称稀疏矩阵。

至此，柔性转子有限元建模的基本框架已经完成。对于柔性转子的其他组成部件，如刚性圆盘，集中质量点、轴承支承、不平衡力以及外力作用都是以在［M］、［K］、［D］和［G］稀疏矩阵中添加相应元素项的形式来实现建模，因篇幅所限不一一叙述。本文下面详细叙述柔性转子的电磁轴承支承的建模。

2.2 柔性转子的电磁轴承支承建模

最直观的方法是可以把电磁轴承等效为机械轴承，即用等效刚度和等效阻尼表示，然后用传统柔性转子有限元建模中的轴承单元表示方法，根据轴承所在的节点，在刚度矩阵［K］和阻尼矩阵［D］中添加相应的元素项。然而这并不是一个好方法，因为等效刚度和等效阻尼的计算需要知道控制传递函数，而考虑实际因素的电磁轴承控制传递函数会非常复杂，难以计算。甚至，诸如模糊控制等智能控制策略根本无法归纳出一个传递函数表达式，所以无法计算等效刚度和等效阻尼。此法局限性太大。

因此，本文采取直接用电磁力建模的方法，基本思想是认为电磁轴承支承就是转子在相应电磁轴承节点位置所受的外力。例如，当节点a和节点b处各有一个电磁轴承的情况，在原转子运动方程式（7）中，广义外力矢量｛F｝需按下式叠加4个轴承电磁力项。

式（11）中，ki，ax，ki，ay，ki，bx，ki，by及kx，ax，kx，ay，kx，bx，kx，by分别表示在a及b节点上的两个电磁轴承各个径向自由度的电流刚度系数和位移刚度系数。Iax，Ibx，Iay及Iby表示在a节点及b节点上电磁轴承在4个自由度的控制电流。

式（11）采用的是线性电磁力的表示方法，也可根据模型精度需要，采用非线性电磁力的描述，也可叠加诸如振动控制等附加力。此建模方法运用较为灵活，不限制控制策略的表示形式，只要控制策略已知，就可以在仿真模型中实现。

3 算例介绍

目前，由于国内对于完善的电磁轴承支承的柔性转子实验平台有待建设，本文将以英国Bath大学机械工程系Keogh教授实验室的电磁轴承-柔性转子实验平台为仿 真算 例［13-17］。该 电磁 轴承 实验台经过多年的动力学实验和研究，其各项动力学参数已被精确测定。本文的仿真数据和仿真结果可以从该实验平台的实际数据中找到依据，以确保本文有关柔性转子的理论分析和仿真结论的可信度。

如图3为该实验平台的基本结构。它由一根均质等直径轴和4个刚性圆盘组成。轴长度为1 930 mm，轴粗50 mm。4个刚性圆盘分别分布在轴的两端和中间。圆盘外径250 mm，圆盘宽30 mm，每个盘质量10 kg。在转子有限元分析中，把电磁轴承处的转子轴颈部分也视为刚性圆盘，外径200 mm，宽度80 mm。整个轴的总质量约100 kg。仿真假设，转子的不平衡质量只集中在4个刚性圆盘上，而忽略其余部分的不平衡质量。刚性圆盘的不平衡质量可调节，将4个刚性圆盘的不平衡质量的大小都设置为1 g·mm。两个电磁轴承的气隙均为1.2 mm，保护轴承气隙0.75 mm。电磁轴承能提供最大1 500 N的提升力。电磁轴承工作在4.3 A的偏置电流下，位移刚度系数为2×106N/M，电流刚度系数为487 N/A。

图3 柔性转子算例的结构和尺寸（8个节点）

4 等效刚度对柔性转子系统各阶临界分布的影响

由等效刚度计算公式，选取不同的Kp参数，可分别得到0.2×106、0.4×106、0.8×106、1.0×106、2×106、3×106、5×106N/m 7个不同的等效刚度值。Kd参数统一取10，即对应等效阻尼为4 800 N·s/m。转子转速范围为0～600 rad/s。通过仿真，可得到在不同等效刚度条件下，转子在右侧电磁轴承y方向的振动响应曲线，共有7组。为了方便分析转子振动响应随支承刚度变化的规律，把7组振动响应曲线集中绘入图4中。可见，在0～600 rad/s的转速范围内，转子的振动响应曲线上都出现了3或4个峰值。这些峰值所对应的转速就是该刚度下转子的各阶临界转速的位置。

图4 不同等效刚度下转子的振动响应曲线的仿真结果

根据转子动力学理论，柔性转子的模态振型和对应的临界转速由转子结构和支承特性所确定。为了弄清图4中各个转子振动峰值是否与转子在此等效刚度下的临界转速值相吻合，以及各临界转速所对应的模态振型，利用本项目组开发的柔性转子分析软件RotFE计算了该柔性转子在同样七组等效刚度条件下的前四阶临界转速值，并绘制了其振型的三维图。因限于篇幅，只绘出4组为例，如表1所示。

表1 RotFE得出的临界转速和其振型三维图

由表1可知，柔性转子的前四阶临界转速的振型分别为锥形涡动、平行涡动、一阶弯曲和二阶弯曲。其中，锥形涡动和平行涡动为刚体临界振型，一阶和二阶弯曲为柔性弯曲临界振型。

在等效刚度较低时，转子刚体临界转速均低于弯曲临界转速。然而，转子的刚体临界转速对等效刚度的变化十分敏感。随着等效刚度的增加，锥动和平动两个临界转速的位置明显地往高频方向移动。相比较，转子的弯曲临界转速受支承刚度变化影响较小。

当等效刚度逐渐升高，转子的锥形和平行临界转速的位置也随之升高。在等效刚度为2×106N/ m时，转子的平行涡动临界转速超越了转子的一阶弯曲临界转速，此时转子自身有明显的弯曲形变又带有刚体平行涡动振型的特点，称其振型为“平行+弯曲”的综合振型。在等效刚度为5×106N/m时，转子的锥形涡动临界转速也超越了转子的一阶弯曲临界转速。同样，此时转子自身已有明显的弯曲形变又带有刚体锥形涡动振型的特点，称其振型为“锥形+弯曲”的综合振型。

为了比较图4的振动峰值数据和表1临界转速数据，把两组数据一起绘制成图5形式进行比较。可见，仿真所得的各阶临界转速值和RotFE计算结果有很好的吻合，只是在0.2×106、0.4×106、0.8× 106、1.0×106N/m4个较小等效刚度的仿真测试中，观测不到转子平动临界转速的振动峰值。分析其原因，可能是锥形涡动临界和平行涡动临界靠得太近，另外，在电磁轴承支承特性中的等效阻尼作用下，可能使振动峰值不明显。

图5 仿真结果与RotFE计算结果的比较

5 等效阻尼对柔性转子系统各阶临界分布的影响

由图4可知，在高刚度5×106N/m条件下，转子在前四阶临界转速的振动比低刚度时更强烈，峰值更清晰明显。因此，在等效阻尼对转子临界转速影响的仿真中，等效刚度值统一选取为5×106N/m。微分参数Kd分别取5、10、50、100，即对应等效阻尼分别为2 400、4 800、24 000、48 000 N·s/m。可得到转子在不同等效阻尼条件下，右侧电磁轴承y方向的振动响应曲线，并将4组振动图线集中绘入图6中。

图6 不同等效阻尼对转子振动响应的影响

结果表明，等效阻尼的强弱不会明显改变柔性转子各阶临界转速的位置，但增加等效阻尼，可以很好地抑制转子振动的峰值。当阻尼较小时，各阶临界转速位置的振动峰值明显。随着阻尼的增大，临界转速的振动峰值逐渐减小，一些临界转速位置的振动峰值变得难以分辨。如当等效阻尼De=48 000 N·s/m时，除了一阶弯曲临界转速有明显峰值外，其余的临界转速峰值均被抑制。这也说明一阶弯曲临界受等效阻尼的影响最小。但是，提高等效阻尼意味着控制系统的微分作用增强，会带来放大系统噪声的负面影响，反而不利于系统稳定。在实际电磁轴承应用中，等效阻尼的调节范围有限。

6 结 论

本文提出采用电磁轴承支承柔性转子，通过支承特性的调节，改变柔性转子的各阶临界转速的分布位置，使柔性转子工作转速始终避开临界转速，实现超临界运行的目的。本文对电磁轴承-柔性转子系统的前四阶临界转速随等效刚度和等效阻尼变化的规律进行了理论仿真，仿真结果表明：

a）转子的刚体临界转速分布对支承等效刚度的变化较为敏感，随支承等效刚度的增大，锥形临界和平行临界转速明显向高频方向移动。转子一阶弯曲临界转速受支承刚度变化影响较小，几乎不随支承刚度的变化而变化。转子二阶临界转速随支承刚度的调节有一定的变化，但无论是相对变化量还是绝对变化量都比刚性临界小得多。

b）等效阻尼的变化不会明显改变柔性转子系统各阶临界转速的位置。增加等效阻尼可以很好地抑制转子振动的峰值，对抑制振动有益，但同时也会放大电磁轴承控制系统的噪声，对系统稳定不利。

因此，可以通过调节转子系统的等效刚度和等效阻尼，来改变转子系统刚体临界转速，但对一阶弯曲临界转速位置的调节效果十分有限。对于亚临界或超临界工作的磁悬浮轴承转子系统，只能通过改变转子本身的结构，尽量使转子工作转速远离一阶弯曲临界转速。

［1］Arredondo I，Jugo J，Etxebarria V.Modeling and control of a flexible rotor system with AMB-based sustentation［J］.ISA Transactions，2008，47（1）：101-112.

［2］Arias M M，Silva N G.Finite element modeling and unbalance compensation for a two disks asymmetrical rotor system［C］//Electrical Engineering，Computing Science and Automatic Control，2008：CCE 2008 5th International Conference on.IEEE，2008：386-391.

［3］Tseng C Y，Shih T W，Lin J T.A kalman filter-based automatic rotor dynamic balancing scheme for electric motor mass production［C］//Materials Science Forum. 2006，505：997-1002.

［4］李红伟，徐 旸，谷会东，等.电磁轴承-挠性转子系统的本机动平衡方法［J］.中国机械工程，2008，19（12）：1419-1422.

［5］Kang Y，Lin T W，Chu M H，et al.Design and Simulation of a Neural-PD Controller for Automatic Balancing of Rotor［M］//Advances in Neural Networks-ISNN 2006.Springer Berlin Heidelberg，2006：1104-1109.

［6］Dyer S W，Adaptive Optimal Control of Active Balancing Systems for High-Speed Rotating Machinery［D］. Michigan：University of Michigan，1999.

［7］汪希平，朱礼进，于 良，等.主动磁轴承转子系统动力学特性的研究［J］.机械工程学报，2001，37（11）：7-12.

［8］周朝暾，汪希平，严慧艳.主动电磁轴承转子系统振动模态的分析研究［J］.机械科学与技术，2004，10（23）：1163-1166.

［9］谢振宇，徐龙祥，李 迎，等.控制参数对磁悬浮轴承转子系统动态特性的影响［J］.航空动力学报，2004，19（2）：174-178.

［10］Xu Y，Yu H，Zhao L，et al.Development of a magnetic bearing rotor dynamics analysis software package［C］.Proc.of the 12th Int Symp on Magnetic Bearings.Wuhan：［s.n.］，2010：53-70.

［11］钟一谔，何衍宗，王 正，等.转子动力学［M］.北京：清华大学出版社，1987：143-194.

［12］徐龙祥.高速旋转机械轴系动力学设计［M］.北京：国防工业出版社，1994：131-182.

［13］Cade IS，Keogh PS，Sahinkaya M N.Fault identification in rotor/magnetic bearing systems using discrete time wavelet coefficients［J］.Mechatronics，IEEE/ ASME Transactions on，2005，10（6）：648-657.

［14］Schlotter M，Keogh P S.Synchronous position recovery control for flexible rotors in contact with auxiliary bearings［J］.Journal of Vibration and Acoustics，2007，129（5）：550-558.

［15］Sahinkaya M N，Abulrub A H G，Keogh P S，et al. Multiple sliding and rolling contact dynamics for a flexible rotor/magnetic bearing system［J］.Mechatronics，IEEE/ASME Transactions on，2007，12（2）：179-189.

［16］Keogh P S.Transient rotor/active magnetic bearing control using sampled wavelet coefficients［J］.Journal of Engineering for Gas Turbines and Power，2007，129：549.

［17］Cole M O T，Keogh P S，Sahinkaya M N，et al.Towards fault-tolerant active control of rotor-magnetic bearing systems［J］.Control Engineering Practice，2004，12（4）：491-501.

Modeling and Simulation Analysis of Controlling Distribution of Critical Speed of Flexible Rotor with Electromagnetic Bearing

FANGPeng，JIANG Ke-jian
（The School of Information Science and Technology，Zhejiang Sci-Tech University，Hangzhou 310018，China）

This paper proposes to use electromagnetic bearing to support flexible rotor and change the distribution location of critical speed of each order of flexible rotor through the adjustment of supporting characteristics，thus making flexible rotor realize supercritical operation smoothly，incorporates the modeling of supporting characteristics of electromagnetic bearing into the modeling theory of finite element of flexible rotor of classic rotor dynamics and establishes flexible rotor system model of electromagnetic bearing support for simulation.The simulation result shows that the adjustment of equivalent stiffness of electromagnetic bearing can significantly change translational and conical clinical speed of flexible rotor and has little influence on first-order and second-order bending critical speed.The adjustment of equivalent damping can significantly reduce the vibration of rotor when it passes the critical point.

electromagnetic bearing；flexible rotor；magnetic levitation；finite element；critical speed

TH133

A

（责任编辑：陈和榜）

1673-3851（2014）03-0221-07

2013-11-04

国家自然科学基金（11272288），中国博士后科学基金（2013M531452），浙江省自然科学基金（LY12E05027），浙江理工大学科研基金（1104827-Y）

方 鹏（1989-），男，硕士研究生，主要从事数字信号处理，磁悬浮控制方面的研究。

蒋科坚，电子邮箱：jkjofzju@163.com