电信号扰动下的轴系扭振摘要本文用一种改进的Riccati扭转传递矩阵结合Newmark-β方法研究非线性轴系的扭转振动响应。

首先，该系统被模化成一系列由弹簧和集中质量点组成的系统，从而建立一个由多段集中质量组成的模型。

第二，通过这种新发展起来的程序可以从系统的固有频率和扭振响应中消除累计误差。

这种增量矩阵法，联合结合了Newmark-β法改进的Riccati扭转传递矩阵法，进一步应用于解决非线性轴系扭转振动的动力学方程。

最后，将一种汽轮发电机组作为一个阐述的例子，另外仿真分析已被应用于分析典型电网扰动下的轴系扭振瞬时响应，比如三相短路，两相短路和异步并置。

实验结果验证了本方法的正确性并用于指导涡轮发电机轴的设计。

关键词：传递矩阵法；Newmark-β法；汽轮发电机轴；电学干扰；扭转振动1.引言转子动力学在很多工程领域起着很重要的作用，例如燃气轮机，蒸汽轮机，往复离心式压气机，机床主轴等。

由于对高功率转子系统需求的持续增长，计算临界转速和动态响应对于系统设计，识别，诊断和控制变得必不可少。

由于1970年和1971年发生于南加州Edison’sMohave电站的透平转子事故，业界的注意力集中在由传动行为导致的透平发电机组内的轴的扭转振动。

当代的大型透平发电机组单元轴系系统是一种高速共轴回转体。

它是由弹性联轴器连接，由透平转子，发电机和励磁机组成。

电力系统故障或操作条件的变化引起的机电暂态过程可能导致轴的扭转振动，而轴的扭转振动对于设计来说是非常重要的。

对于透平发电机轴系扭振的研究，如发生次同步谐振和高速重合，基本的是对固有频率和振动响应的计算的研究。

当前，有限元法和传递矩阵法是最流行的两种分析轴系扭振的方法。

有限元法（FEM）通过二阶微分方程构造出转子系统直接用于控制设计和评估，而传递矩阵法（TMM）解决频域内的动态问题。

TMM使用了一种匹配过程，即从系统一侧的边界条件开始沿着结构体连续的匹配到系统的另一端。

所有边界点边界条件的匹配满意度为解决方案的定位提供了基础，而且在一个特定点转子系统的状态通过转移矩阵在连续点之间进行转换。

这种方法尤其适合“链连接”结构，例如转子系统，而且更主要的是它不需要对大型数组的存储和操作[2]。

图1.多段集中质量模型TMM首先由Prohl提出。

直到现在，在扭转振动特性上已经做出了许多贡献[4–9]。

Pestel 和Leckie[10]为应用TMM法决定振动系统的固有频率和阵型提供了详尽的参考。

Rao [11]应用TMM去分析自由振动，瞬态响应，临界转速和转子扭振系统的稳定性。

In Ref. [12],Lee et al.通过扩倍状态变量的数量改进了连续系统的TMM法来适应线性转子轴承系统的同步椭圆轨道。

同时，他们的研究也考虑了转动惯量和回转剪切力的影响。

更进一步，连续系统的TMM法拓展到了不平衡轴[13]和非对称转子轴承系统[14]。

Aleyaasin et al. [15]应用TMM法对安装在油润滑轴承上的转子进行了弯曲振动分析。

Chen [16]研究了变截面且带有附加质量的气缸的扭转振动。

Koser and Pasin [17]通过分析法讨论了驱动轴和机械装置的扭振。

另外，TMM法和拉格朗日动力学法经发展用来研究在[18,19]中轴系的扭振和横振。

通常，有两种直接时间积分法去获得动态系统的瞬态响应。

第一种是显式法例如龙格-库塔法，第二种是隐式法例如Newmark法。

在隐式法中，Newmark法是建立在平均加速度的基础上来保证数值的稳定性同时拥有二阶精度，这种方法也是应用最为广泛的一种方法[20]。

在引用[21]中，Newmark-β法连同Newton-Raphson法被用来解决非线性差分方程的迭代问题。

在本文工作中推荐一种模化的Riccati传递矩阵法。

与Newmark-β法结合后被用来研究在电网故障和扰动下透平发电机组轴的扭振特征频率和响应。

在第二部分，转子系统被模化成由弹性弹簧相连接的连续集中质量系统，之后用模化后的Riccati传递矩阵法来驱动特征频率方程。

在第三部分，模化后的矩阵法联合Newmark-β法被进一步用来解决轴系的扭转振动响应。

在第四部分，通过建立一个实验转子系统来验证我们方法的准确性。

最后，对我们的结果进行讨论和总结。

2. riccati 扭转传递矩阵的模化2.1转子系统模型根据轴的结构特性，轴系可以被分为n 个扭转振动单元，连接所有单元带有大的惯性的中心线组成了系统的中心。

从图1可以看出，轴等价于一个多段集中质量模型。

它包含一个长度为L 扭转刚度为Ki 的统一多自由度的扭杆，带有n 段转动惯量为Ji (i = 1, 2, . . . n).的质量段。

感兴趣的是它建立了扭转振动系统的n 个特征频率的关系。

图形1中展示的系统的运动公式可以写成以下几种典型的形式：图2.单元的动态特性()()tJ t K t T θθ+= (1)公式中t T 表示作用在轴上的外力矩，θ和θ表示旋转角位移和角加速度，J 和K 分别代表惯性矩和刚度矩阵。

当外力矩t T =0时，系统的扭转传递矩阵自振频率和振型可以通过下面的方程得到：()()0J t K t θθ+= (2)2.2．Riccati 扭转传递矩阵在图2中，系统被进一步被分为一系列具有动态特征的单元体。

端点的状态向量代表了系统的扭振状态，传递矩阵代表了相邻两单元之间状态向量之间的关系。

图2中第i 个圆盘的运动方程可以写作：R L i i i i J T T θ=- (3)向量i θ表示第i 个扭转角位移，i T 表示第i 个内部扭转力矩，Ri T 表示右侧的扭转力矩，L i T 表示左边的扭转力矩。

两侧变量之间的关系可以用如下形式来表示：2101i RLT T J i iθθ-ω⎡⎤⎡⎤⎡⎤= ⎣⎦⎣⎦⎣⎦ (4)其中Li T =1R i T - (5)111L R i i i i T K θθ---=+(6)将方程(5)和(6)代入(4)得1221111i i i i L K T p J p J iK RT i θθ-- --⎡⎤⎡⎤⎡⎤⎢⎥= ⎣⎦⎣⎦⎢⎥⎣⎦ （7）P 表示角频率。

当我们引入Riccati 传递矩阵[]{}[]{}:RRi i ii f S S e =，[]{}Re if iZ =，{}{},i i ii T f e θ==时，方程（7）可以被下式代替：111221221RRU U ee U U ffi ii +⎡⎤⎡⎤⎡⎤ =⎣⎦⎣⎦⎣⎦(8)[]{}{}[][]21221111112R ii i i i R i i i i f U U S S U U e S ++++==+ (9)方程(9)可以作为用来计算轴系扭振自振频率的递推公式。

在特定步长下，如果两次迭代差值符号不同，那么在这两个差值之间必然有一个零点此时即为实验频率。

也就是说在两个实验频率值之间必然存在一个自振频率的准确值，可以通过微分计算来得到。

当自振频率确定下来之后其振型也就通过方程(7)随之确定下来。

2.3. 改进的riccati 扭转传递矩阵方法过程中在应用riccati 方法过程中，由于除法运算频率方程通常会导致无穷奇点，这些奇点会影响自振频率计算的有效性。

为了排除这些奇点，方程(9)被改进为下式：[][][]11121111120()Nii i i i i i i U U S R S abs U U S +=+==+∏(10)其中N 是特征单元总数量（见图1和图2），“Π”表示连乘，abs （x ）表示取x 的绝对值。

改进后得到的方程(10)不仅能保持之前经典传递矩阵的优点，而且还能避免丢失固有频率或者由于无穷奇点得到不存在的频率解。

3.增量传递矩阵法3.1. 增量方程转子系统是一个高度复杂的非线性振动系统。

对于这样一个系统，传统振动系统和微扰方法已经不再适用。

目前有两种方法来计算转子的瞬态振动响应，例如TMM 法和模态综合法，两者均靠数值积分辅助计算。

但是它们不能直接应用到转子系统的扭振研究中去。

在这一小节，我们结合Newmark-β法和TMM 法来建立转子系统的增量方程。

Newmark-β法是一种积分方法，可以有效地求解非线性方程。

这种方法中的迭代运算通常可以用加速度和速度来表示成下面形式：++2111()(1)2t t t t t t t qq q q qt t βββ∆∆=----∆∆ (11)++()t t t t t t t qq q q q t γγββ∆∆=+---∆ 1(1)2t qβ- (12)其中γ和β为常数，q 是广义坐标。

由式(11)和(12)，增量方程可以被写作：+++2111()()2t t t t t t t t t qq q q q q t t βββ∆∆∆∆=-=∆--∆∆ (13)+++2()(1)2t t t t t t t t t qq q q q q t tγγγβββ∆∆∆∆=-=∆---∆∆ (14)现在，方程(13)-(14)给出了在t ∆时间内速度和加速度的变化量，而且可以被用来解决非线性方程问题。

3.2 扭振增量传递方程由图1，典型的轴系可以被模化成多圆盘组装在一起的系统。

每个圆盘都在惯性力矩，阻力矩，两端截面产生的力矩和外力矩共同作用下。

根据力矩平衡条件，增量方程可以表达成下式：R i T ∆(t+Δt) =l i T ∆(t +Δt) +i J θ∆ (t + Δt) +i C θ∆ (t +Δt)-ΔTti (t+Δt) (15)其中R Li i θθ∆=∆,i C 是阻力系数，ti T ∆表示惯性矩。

将方程(13)—(14)代入(15)得：R L L i i i i i T T A B θ∆=∆+⋅∆+ （16）其中2111,22i i i i i i i i i i J C A B J C t t t t γγγθθθθββββββ⎡⎤⎛⎫⎛⎫⋅+⋅-=+=--⋅+⋅∆ ⎪⎢ ⎪⎥⋅∆⋅∆⋅∆⎝⎭⎝⎭⎣⎦()ti T t t -∆+∆ (17)对于无质量的弹性轴端 ，递增关系可以表示为下式：+1+1==+R L i i RL R i i i i T T T K θθ∆∆∆∆∆⎧⎨⎩ （18）至于方程（8），对第i 个圆盘，如果我们令i i f T =∆，i i e θ=∆那么传递矩阵可以被写作：{}111221221f eLLF U U ffU U F eei ii +⎡⎤ =+⎡⎤⎡⎤⎣⎦⎣⎦⎣⎦ （19）其中111U =，12i U A =，211i U K =，221ii AU K =+，,/fj i ei i i F B F B K == 将Riccati 传递关系L Lii i i f S e P =+代入方程（19）得到： ‘[]11112i i i i S U S U +=+[]2122i i i U S U ⋅+[]111121i i fi i i i i ei U P F P S U PF +++=-+⎡⎤⎣⎦ （20）[][][]1121221212221L L i i i i i i i i i i ei e U S U e U S U U P F --+=+-++那么,Li i i S P e 和可以有方程（20）推导出来。