专利名称：汽轮发电机组轴系扭振的分析方法
技术领域：
本发明属于电力和大型回转机械制造技术领域，尤其涉及一种汽轮发电机组轴系扭振的分析方法。
背景技术：
随着我国电力事业的发展，近年来国家兴建的火电厂装机容量不断增大，机组轴系长度增加，使得轴系横截面积相对下降，汽轮发电机组轴系不能再视为转动刚体，而应作为一个复杂的弹性连续体进行研究。与此同时，“坑口”电站日益增多，电厂与负荷中心的距离越来越远，输电网容量也不断扩大，电网的控制系统也使系统变得复杂，这使得诱发机组轴系扭振的因素日益增加。因此，扭振成为影响机组安全运行的主要问题之一。由于因扭振而造成的轴系损伤具有很大的隐蔽性。一旦由于扭振带来的轴系疲劳积累到一定程度，就可能导致轴系断裂、崩溃，其后果是灾难性的，损失也将极为惨重。可见，对汽轮发电机组轴系扭振进行分析，保证机组的安全经济运行，是非常必要和重要的。在扭振分析中，建模与动态响应计算是很重要的环节，轴系模化结果以及计算结果的精确程度直接影响分析结果，从而影响对机组的运行决策。
自1971年美国Mohave电站发生事故以来，美国、德国、日本等工业发达国家相继展开对汽轮发电机组轴系扭振的研究，已经取得了显著成效。国内在这方面的研究起步较晚，自20世纪80年代以来，发生了一系列由机组轴系扭振引起的部件损坏事故。自此之后，国内进行了大量的扭振研究，研究的主要方向是轴系扭振产生的机理、轴系扭振故有频率计算、扭振监测以及寿命损耗计算等。
常用的扭振力学模型有两种，一种是连续质量模型，又称分布质量模型，是利用有限元法将轴系离散表示成有限结构元素进行连续性分析，此方法有较高的精度，但是占用大量的计算机内存及运算时间，所以并没有被广泛采用；另一种是集中质量模型，即通常所说的质量-弹簧模型，是由一连串具有转动惯量的集中质量块和无质量但具有一定刚度的弹簧连接的系统。按照划分质量块的多少，集中质量模型又可以分为简单集中质量模型和多段集中质量模型，简单集中质量模型精度不够，只能用来做定性分析，而多段集中质量模型有较高的精度，计算过程也很简便，所以被广泛采用。无论是何种模型，对转子轮盘上叶片的模化是必不可少的。目前国内对叶片转动惯量的计算并没有一个精确的方法，如果计算不准确，对轴系扭振特性分析有很大影响。
汽轮发电机组扭振分析可分为稳定性分析和动态响应分析两种。稳定性分析的方法包括特征根分析法、等值阻抗法、复力矩系数法和传递矩阵法。扭振的响应计算是扭振分析的一个重要环节，包括计算各个瞬时轴系上各个结点的扭角、角速度和角加速度响应等。动态响应分析主要使用时域仿真法，计算中必须考虑系统的强非线性，与数值积分法或逐步积分法结合求解。
发明内容
本发明的目的在于，提出一种汽轮发电机组轴系扭振的分析方法，解决目前使用的汽轮发电机组轴系扭振分析方法中，叶片转动惯量计算不准确和动态响应分析计算中非线性因素考虑不全面的问题。
本发明的技术方案是，一种汽轮发电机组轴系扭振的分析方法，其特征是所述分析方法包括下列步骤 步骤1进行轴系本体分段，确定每一轴段外径、内径和轴段长度； 步骤2计算轴系当量参数，包括每一轴段的转动惯量和刚度； 步骤3采用有限积分法计算叶片的转动惯量； 步骤4将分段后的轴系，模化成为集中质量块-弹簧模型； 步骤5计算轴系扭振模型外力矩； 步骤6采用Newmark-β法和Riccati增量传递矩阵法计算轴系扭振响应。
所述转动惯量计算公式为其中R为半径；ρ为材料密度；r1、r2分别为轴段的内外半径。
所述刚度计算公式为k＝HJp/L，其中G为材料剪切弹性模量；L为轴段长度；Jp为截面极惯性矩，D1、D2分别为轴的内外径。
所述步骤3中，叶片的转动惯量计算公式为
其中，A1、A2分别为叶片顶部和根部的截面面积；r1、r2分别为叶片顶部和根部的半径，Zd为叶片排动叶数目。
所述步骤4还包括下列步骤 步骤41根据分段结果，将轴系模化为n个单元轴段，每一个单元轴段包含一个有转动惯量而无弹性变形的集中质量块和与之相邻的有弹性而无转动惯量的弹簧； 步骤42以具有较大转动惯量的部件中心线作为集中质量块的中心点； 步骤43轴系上安装的部件作为附加转动惯量，加到相应的集中质量块上，即I＝I1+I2，I1为集中质量块转动惯量，I2为叶片转动惯量，I为考虑叶片后集中质量块的转动惯量； 步骤44两相邻集中质量块中心点间连接轴段的刚度作为该两集中质量之间无质量弹簧的刚度值，即K为无质量弹簧的扭转刚度，Ki、Ki+1为两相邻集中质量块的扭转刚度。
所述步骤5中，计算轴系扭振模型外力矩包括计算作用在汽轮机转子上的蒸汽力矩和作用在发电机转子上的电磁力矩。
所述计算作用在汽轮机转子上的蒸汽力矩的公式为其中T0、T0i分别为第i级的额定进汽温度和实测进汽温度，p2i为第i级的实测排气压力，p0i为第i级的实测进汽压力。
所述计算作用在发电机转子上的电磁力矩的公式为Me＝(-xdid+xafif+xaDiD)iq-(xqiq+xaqiq+xaGiG)id，其中xd为定子绕组d轴自感系数，xq为q轴自感系数，xad为d轴互感系数，xaq为q轴互感系数，If为励磁电流，id、iq分别为定子绕组d轴、q轴电流，ID、IQ为阻尼绕组在d轴和q轴的电流；并且有xaf＝xaD＝xad，xaQ＝xaq。
所述步骤6中，Newmark-β法用于计算各轴段瞬时的扭转角速度和角加速度。
所述步骤6中，Riccati增量传递矩阵法用于计算各轴段瞬时的扭角增量。
本发明的效果在于，采用有限积分法计算叶片的转动惯量，获得的叶片转动惯量更加准确，同时，采用Newmark-β法和Riccati增量传递矩阵法结合的方法计算扭振动态响应，解决了扭振形态随时间和空间变化的问题和振动中的非线性问题。
图1是本发明提供的汽轮发电机组轴系扭振的分析方法流程图。
具体实施例方式 下面结合附图，对优选实施例作详细说明。应该强调的是，下述说明仅仅是示例性的，而不是为了限制本发明的范围及其应用。
图1是本发明提供的汽轮发电机组轴系扭振的分析方法流程图。如图1所示，本发明的具体实现方式包括 步骤1进行轴系本体分段，确定每一轴段外径、内径和轴段长度。
多段集中质量模型是将轴系分解成有转动惯量而无弹性变形的集中质量块和只有弹性而无转动惯量的弹性轴段。为保证精确度，首先需要将轴系详细分段，可以根据汽轮发电机组的设计图纸，凡是轴横截面有突变的地方及存在集中惯量的位置，都取作分段点；长度较长的等截面轴段(如发电机转子绕组轴段)划分为若干个等长度的轴段。
轴系本体分段的目的是得到每一轴段的内径、外径、轴段长度的几何参数，以便进行轴系当量参数(轴段的转动惯量、扭转刚度)计算。
步骤2计算轴系当量参数，包括每一轴段的转动惯量和刚度。
在多段集中质量模型基础上，轴系的扭转振动方程可描述为 Iθ″+Dθ′+Kθ＝M(2-1) 式中θ为扭转角位移向量； I为转动惯量矩阵，I＝[I1I2......In]T，n为轴系分段数； D为阻尼矩阵，D＝[D1D2......Dn]T，n为轴系分段数； K为刚度矩阵，M＝[M1M2......Mn]T，n为轴系分段数； M为外力矩向量，包括施加在各轴段蒸汽扭矩和电磁扭矩，M＝[M1M2......Mn]T，n为轴系分段数。
由(2-1)式可以看出，当转动惯量矩阵I、阻尼矩阵D、刚度矩阵K以及外力矩矩阵M已知，则方程可解。轴系当量参数包括转动惯量I以及刚度矩阵K，它们由轴系结构及轴系材料特性决定。
转动惯量计算公式为 式中r为半径；ρ为材料密度；r1、r2分别为轴段的内外半径。
刚度的计算公式为 k＝GJp/L (2-3) 式中G为材料剪切弹性模量；L为轴段长度；Jp为截面极惯性矩，D1、D2分别为轴的内外径。
每一轴段的r1、r2、L均可由步骤1的分段后得到，剪切弹性模量G以及材料密度ρ可跟据材料类型查表得到。
步骤3采用有限积分法计算叶片的转动惯量。
在扭振特性的计算过程中，安装有叶片的轴段的转动惯量是由两部分组成，一部分是轴段本体的转动惯量，另一部分是安装在轴段上叶片的转动惯量。将叶片作为附加转动惯量，并通过相应算法求出其转动惯量的大小后，与相应的轴段本体转动惯量相加得到轴段与叶片整体的转动惯量。
汽轮机叶片是一个形状复杂的几何体，通常无法列出其结构的曲线方程，因此，只能使用近似的方法来计算其转动惯量。采用线性处理的有限积分法，将叶片分成许多小的简单几何体，然后求出这些有限小的简单几何体对同一转轴的转动惯量之和。计算时，可将叶片不同半径ri处的截面积fi，用一个宽度为bi、半径为ri，表面积与叶片在该处的截面积及转动惯量相同的圆柱面来代替，即 式中，Zd为叶片排动叶数目。
因此，bi曲线所包络平面的转动惯量就等于叶片排的转动惯量WR02。
式中，ρ为叶片材料密度。
为了简化计算，假定叶片截面积随半径成线性变化，截面积的表达式为 A(r)＝kr+c(3-3) 式中，k，c为待定参数。则 式中，A1、A2分别为叶片顶部和根部的截面面积；r1、r2分别为叶片顶部和根部的半径。
在转动惯量的计算过程中，由于叶片和叶顶的截面积无法直接获取，因此需对其进行近似估算。在现场对其进行实测的过程中，采用软丝近似得到了叶片顶部和底部的截面积，并将获取的截面图形通过软件转化，最终计算出截面积大小。但通过计算我们发现，这一测量结果并不理想，因此还需通过叶片的质量对转动惯量进一步的修正，叶片质量为m，密度为ρ则 将软丝截取的截面图形面积代入(3-4)式求得k，将此k值代入(3-5)式求得c，我们将它记做c1，再将c1代入(3-5)式求得k1。
将软丝截取的截面图形面积代入(3-4)式求得c，将此c值代入(3-5)式求得k，我们将它记做k2，再将k2代入(3-5)式求得c2。令则 A(r)＝k′r+c′(3-6) (3-6)式就是经叶片质量修正后的叶片截面积计算公式，用此式计算后的结果与现场实测结果非常接近。
叶片截面积转化为轮盘面积，
式中，i为动叶的级数。
转动惯量计算公式I＝∫r2dm，其中叶片dm＝A轮盘·ρ·dr，ρ为叶片密度。
则叶片对轴的转动惯量I为
步骤4将分段后的轴系，模化成为集中质量块-弹簧模型。
将分段后的轴系，利用步骤2所计算的当量参数K、I，模化成为集中质量块-弹簧模型。轴系模化的具体方法是 步骤41根据轴系分段结果，将轴系模化为n个单元轴段，每一个单元轴段包含一个有转动惯量而无弹性变形的集中质量块和与之相邻的有弹性而无转动惯量的弹簧。
步骤42以具有较大转动惯量的部件中心线作为集中质量块的中心点。
步骤43轴系上安装的部件作为附加转动惯量，加到相应的集中质量块上，即I＝I1+I2，I1为集中质量块转动惯量，I2为叶片转动惯量，I为考虑叶片后集中质量块的转动惯量。
步骤44两相邻集中质量块中心点间连接轴段的刚度作为该两集中质量之间无质量弹簧的刚度值，即K为无质量弹簧的扭转刚度，Ki、Ki+1为两相邻集中质量块的扭转刚度。
经过以上方法模化，轴系变为由一系列单元轴段组成的多段集中质量模型(集中质量块-弹簧模型)。
以某电厂500MW机组为研究对象，对机组进行了模化，在不考虑叶片附加转动惯量情况下，计算出了轴系扭振的固有频率，并得到了扭振各阶固有振型。然后通过以上对叶片转动惯量的计算方法，对该电厂机组叶片进行现场测量并计算其附加转动惯量，并计算得到考虑附加转动惯量后的轴系扭振固有特性。
对上述两计算结果进行比较发现，在考虑了叶片附加转动惯量以后，固有频率的计算结果有所降低，而且振型也发生了变化，其计算结果更加接近实测的轴系扭振固有频率。由此可以得出结论考虑叶片附加转动惯量的影响及其惯量的准确计算，对于提高轴系扭振响应计算精度有很大帮助，从而使汽轮发电机组轴系扭振在线监测所得到的数据更加可靠，具有很好的工程实践意义。
步骤5计算轴系扭振模型外力矩。
由式(2-1)知，分析轴系的扭转振动必须得到在轴系上施加的外力矩M，才能根据步骤1所建立的多段集中质量模型求解轴系扭振响应。轴系上作用的外力矩包括两部分作用在汽轮机转子上的蒸汽力矩和作用在发电机转子上的电磁力矩。
(1)作用在汽轮机转子上的蒸汽力矩的计算 蒸汽在动叶片安装处集中质量块上施加的扭矩为 式中n为轴系转速；pi为级的内功率。
级的内功率可由级的有效比焓降和蒸汽流量求得 式中D为级的进气量，单位是t/h；Gi为蒸汽流量，单位是kg/s。
在临界工况下，级内压力与流量的关系如下 式中p0、p0i分别为第i级的额定进汽压力和实测进汽压力，MP；T0、T0i分别为第i级的额定进汽温度和实测进汽温度，K；Gc为第i级的临界流量，kg/s。
在亚临界工况下，根据弗留格尔公式 式中p2、p2i分别为第i级的额定排气压力和实测排气压力，MP；G为第i的额定流量；设计工况下和变工况下级的理想比焓降可分别表示为 则 以上额定数据Gc、p0、p2、Δht和T0可以通过热力计算得到，也可由汽轮机制造厂家提供，其余实际数据可以通过现场实测得到。
(2)电磁力矩计算 电磁转矩以发电机参数和定、转子电流表示为 Me＝ψdiq-ψqid (5-8) 即 Me＝(-xdid+xafif+xaDiD)iq-(xqiq+xaqiq+xaGiG)id (5-9) 式中xd为定子绕组d轴自感系数； xq为q轴自感系数； xad为d轴互感系数； xaq为q轴互感系数； If为励磁电流； id、iq分别为定子绕组d轴、q轴电流； ID、IQ为阻尼绕组在d轴和q轴的电流； 并且有xaf＝xaD＝xad，xaQ＝xaq。
计算得到蒸汽施加在叶片上的力矩后，将蒸汽力矩施加到动叶片所在的集中质量块上。计算得到施加在发电机转子绕组轴段的电磁力矩后，根据对绕组轴段的分段(步骤1)将电磁力矩平均施加到每个集中质量块上。至此，我们便得到了轴系模型的外力矩向量M，以便进行轴系响应计算。
步骤6采用Newmark-β法和Riccati增量传递矩阵法计算轴系扭振响应。
本发明采用Newmark-β法和Riccati增量传递矩阵法结合的算法来计算扭振响应。逐步积分法解决扭振形态随时间变化情形，增量传递矩阵法解决扭振形态的空间变化情形。传递矩阵法计算速度快，而且能够保持一定的精度，不必求解固有频率、主振型、模态刚度矩阵和模态质量矩阵等，而且传递矩阵的维数不随系统自由度的增加而增加。本发明采用的增量传递矩阵法不同于以往的传递矩阵法，此方法将扭矩、扭角以增量的方式列出传递矩阵，可解决振动中的非线性问题。
使用增量传递矩阵法求解时，经过步骤1、2、3，将轴系模化为由单元轴段组成的集中质量块-弹簧模型，每一个单元轴段包含一个有转动惯量而无弹性变形的集中质量块和与之相邻的有弹性而无转动惯量的弹簧。计算精度则随着划分质量块的数量的增加而逐渐提高。以下各公式中，K为扭转刚度，I为集中质量块的转动惯量，均由步骤2、3计算得知；TL(t)为施加在集中质量块的外扭矩，由步骤4计算得知；TiR、TiL为相邻质量快传递的扭矩，角标L表示质量快左截面，R表示质量快右截面；θ为轴段扭角，i表示第i轴段，i+1表示第i+1轴段。
首先，进行集中质量块的扭转受力分析。
分析集中质量块受力可知，作用在质量块上的扭矩包括了左右截面扭力矩、阻力矩以及外力矩，由平衡条件得知 改写成增量形式为 其中， 式中Ii为转动惯量，Ci为阻尼系数，θ为该微元段的扭角，并且 其次，进行无质量弹簧的扭转受力分析。
由于无质量弹簧的转动惯量忽略不计，因而有 其增量表达式为 式中Ki，i+1为扭转刚度， 接下来，计算单元轴段传递矩阵。
单元轴段传递矩阵是单元轴段间扭矩、扭角的传递关系。设f＝ΔT，e＝Δθ，则单元轴段传递矩阵可写为 式中U11＝1，U12＝Ai，Ff＝Ti，然后，对Newmark-β法进行变换。
根据Newmark-β法，有 式中γ，β是Newmark-β法的精度控制参数，可根据计算要求进行选择，当γ，β满足条件γ≥1/2，β≥0.25(0.5+γ)2时，它为无条件稳定自起步格式；当γ＝1/2时，其计算精度为二阶；否则，由于一些与“数值阻尼”及“周期延长”相关的误差，其精度降为一阶。γ＝1/2，β＝1/4时，退化为常平均加速度法；γ＝1/2，β＝1/6时，退化为线性加速度法；γ＝1/2，β＝1/8时，加速度为台阶形变化。q为转子广义坐标。上式说明，如果已知瞬时t某结点的位移qt、速度
加速度
且还知道t+Δt时刻该结点的速度
和加速度
可将上式改为增量表达式 最后，应用Riccati法求解轴系响应。
引入Riccati传递矩阵法，设Riccati变换为 fi＝Siei+Pi (6-14) 带入式(6-13)，则得到Si和Pi的递推公式及ei的递推公式 Ri+1＝[U11·P+Ff]i-Si+1[U21·P+Fe]i (6-16) 由边界条件(fi＝0，ei≠0且为任意值)，可得S1＝0，P1＝0，利用(6-15)式，则可顺次得到S2，P2，S3，P3......，SN+1，PN+1。
对于右端截面N+1，有 fN+1＝SN+1·eN+1+PN+1 (6-18) fN+1＝0 (6-20) 利用式(6-17)式，可从右到左算出各ei(i＝N，N-1，N-2，......，1)，相应也可算出各截面fi，所求结果即为各截面在t+Δt时刻各状态矢量增量。也就是说，从起始时刻t0的初始条件θ10开始(由轴端扭角传感器测量计算可得)，在已知轴系各物理参数I1，I2，…IN，C1，C2…CN，K1，2，K2，3…KN-1，N和外力矩TL1，TL2…TLN条件下，应用上述Riccati法，求出t0+Δt瞬时各轴段的扭角增量(需要注意，外力矩施加的轴段位置)，然后用Newmark-β法，可求出t0+Δt瞬时各轴段的扭转角速度和角加速度，进而求得t0+2Δt瞬时各结点的扭角等，如此反复，可得到各个瞬时轴系上各个结点的扭角、角速度和角加速度响应。
以上所述，仅为本发明较佳的具体实施方式
，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。
权利要求
1.一种汽轮发电机组轴系扭振的分析方法，其特征是所述分析方法包括下列步骤
步骤1进行轴系本体分段，确定每一轴段外径、内径和轴段长度；
步骤2计算轴系当量参数，包括每一轴段的转动惯量和刚度；
步骤3采用有限积分法计算叶片的转动惯量；
步骤4将分段后的轴系，模化成为集中质量块-弹簧模型；
步骤5计算轴系扭振模型外力矩；
步骤6采用Newmark-β法和Riccati增量传递矩阵法计算轴系扭振响应。
2.根据权利要求1所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述转动惯量的计算公式为其中R为半径；ρ为材料密度；r1、r2分别为轴段的内外半径。
3.根据权利要求1所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述刚度的计算公式为k＝GJp/L，其中G为材料剪切弹性模量；L为轴段长度；Jp为截面极惯性矩，D1、D2分别为轴的内外径。
4.根据权利要求1所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述步骤3中，叶片的转动惯量计算公式为
其中，A1、A2分别为叶片顶部和根部的截面面积；r1、r2分别为叶片顶部和根部的半径，Zd为叶片排动叶数目。
5.根据权利要求1所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述步骤4还包括下列步骤
步骤41根据轴系本体分段结果，将轴系模化为n个单元轴段，每一个单元轴段包含一个有转动惯量而无弹性变形的集中质量块和与之相邻的有弹性而无转动惯量的弹簧；
步骤42以具有较大转动惯量的部件中心线作为集中质量块的中心点；
步骤43轴系上安装的部件作为附加转动惯量，加到相应的集中质量块上，即I＝I1+I2，I1为集中质量块转动惯量，I2为叶片转动惯量，I为考虑叶片后集中质量块的转动惯量；
步骤44两相邻集中质量块中心点间连接轴段的刚度作为该两集中质量之间无质量弹簧的刚度值，即K为无质量弹簧的扭转刚度，Ki、Ki+1为两相邻集中质量块的扭转刚度。
6.根据权利要求1所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述步骤5中，计算轴系扭振模型外力矩包括计算作用在汽轮机转子上的蒸汽力矩和作用在发电机转子上的电磁力矩。
7.根据权利要求6所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述计算作用在汽轮机转子上的蒸汽力矩，其计算公式为其中T0、T0i分别为第i级的额定进汽温度和实测进汽温度，p2i为第i级的实测排气压力，p0i为第i级的实测进汽压力。
8.根据权利要求6所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述计算作用在发电机转子上的电磁力矩，其计算公式为Me＝(-xdid+xafif+xaDiD)iq-(xqiq+xaqiq+xaGiG)id，其中xd为定子绕组d轴自感系数，xq为q轴自感系数，xad为d轴互感系数，xaq为q轴互感系数，If为励磁电流，id、iq分别为定子绕组d轴、q轴电流，ID、IQ为阻尼绕组在d轴和q轴的电流；并且有xaf＝xaD＝xad，xaQ＝xaq。
9.根据权利要求1所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述步骤6中，Newmark-β法用于计算各轴段瞬时的扭转角速度和角加速度。
10.根据权利要求1所述的一种汽轮发电机组轴系扭振的分析方法，其特征是所述步骤6中，Riccati增量传递矩阵法用于计算各轴段瞬时的扭角增量。
全文摘要
本发明公开了电力和大型回转机械制造技术领域中的一种汽轮发电机组轴系扭振的分析方法。所述方法包括现将轴系本体分段，确定每一轴段外径、内径和轴段长度；再计算轴系当量参数，包括每一轴段的转动惯量和刚度；然后采用有限积分法计算叶片的转动惯量；之后将轴系本体分段后的轴系，模化成为集中质量块-弹簧模型；接下来计算轴系扭振模型外力矩；最后采用Newmark-β法和Riccati增量传递矩阵法计算轴系扭振响应。本发明将非线性因素考虑到轴系扭振的分析中，解决了扭振形态随时间和空间变化对分析结果产生的影响。
文档编号G01H1/00GK101603855SQ200910088780
公开日2009年12月16日 申请日期2009年7月16日 优先权日2009年7月16日
发明者顾煜炯, 何成兵, 尹绪雨 申请人:华北电力大学