CN109657397B - 基于频响函数的汽轮机叶片-转子系统稳定性的预测方法 - Google Patents
基于频响函数的汽轮机叶片-转子系统稳定性的预测方法 Download PDFInfo
- Publication number
- CN109657397B CN109657397B CN201811643941.1A CN201811643941A CN109657397B CN 109657397 B CN109657397 B CN 109657397B CN 201811643941 A CN201811643941 A CN 201811643941A CN 109657397 B CN109657397 B CN 109657397B
- Authority
- CN
- China
- Prior art keywords
- blade
- damping
- excitation
- shroud
- stability
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/17—Mechanical parametric or variational design
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q10/00—Administration; Management
- G06Q10/04—Forecasting or optimisation specially adapted for administrative or management purposes, e.g. linear programming or "cutting stock problem"
Abstract
本发明公开了基于频响函数的汽轮机叶片‑转子系统稳定性的预测方法,它解决了现有技术中转子系统无法准确稳定性预测的问题,具有能避免设计阶段无法准确预测机组稳定性造成的“先天不足”,将大幅提高设计效率降低设计成本的有益效果,其方案如下:基于频响函数的汽轮机叶片‑转子系统稳定性的预测方法,利用模态正交性求解大型多自由度汽轮机转子系统的频响函数,计算获得转子系统特征乘子,进而利用特征乘子与单位圆的关系判定系统稳定性。
Description
技术领域
本发明涉及非线性转子动力学领域,特别是涉及基于频响函数的汽轮机叶片-转子系统稳定性的预测方法。
背景技术
本部分的陈述仅仅是提供了与本公开相关的背景技术信息,不必然构成在先技术。
核电、热电联产等大型汽轮机组因绿色、高效等突出优势成为未来能源行业的主要发展趋势。转子系统一旦发生振动超标失稳,将对整个机组的安全“一票否决”,尤其大型汽轮机转子具有弹性长轴系、大柔性叶片等结构特点,且主要采用围带实现抑振增稳,在设计阶段准确预测转子系统稳定性对确保机组安全运行具有重要意义。
目前关于机组稳定性预测方法的报道多出现于论文中。文献中公开的大型机组稳定性预测方法主要分为以下两种:简化叶片的轴-轮盘稳定性预测方法和只考虑叶片的颤振预测方法。例如题为1000MW超超临界汽轮机汽流激振特性研究与应用、叶片-转子-轴承耦合系统的非线性动力学特性研究和高速透平机械稳定性评价与控制原理及方法研究等论文[1-3]中忽略叶片柔性及几何非线性,将其简化为集中质量或单摆模型,无法考虑流场在叶片表面产生的压力脉动诱发的围带摩擦阻尼及其对机组稳定性的影响。题为二维叶栅气固耦合颤振分析、航空发动机叶片非线性动力学分析和基于流固耦合的叶片颤振分析等论文[4-6]以及公开号为CN101908088A和CN101599104A等专利[7,8]中仅给出了叶片的颤振预测方法,忽略了轴系结构,偏离机组实际结构和运行工况,无法真实预测整个机组转子系统的稳定性。
弹性长轴系、大柔性叶片及围带等结构特点及其与流场间的流固耦合作用使大型机组转子系统在不同运行工况下的动刚度、阻尼(尤其围带摩擦阻尼)具有强激励(含转速)依赖非线性,过度简化则无法求解。现有的转子系统稳定性预测及优化设计方法,不能综合考虑机组与流场的流固耦合作用以及轴系与叶片间变形耦合的影响,尤其对于采用围带结构增稳减振的转子系统,在设计阶段无法准确预测稳定性难以避免只凭经验设计的“先天不足”。
因此,需要对基于频响函数的汽轮机叶片-转子系统稳定性的预测方法进行新的研究设计。
发明内容
为了克服现有技术的不足,本发明提供了基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,有效避免设计阶段无法准确预测机组稳定性造成的“先天不足”,将大幅提高设计效率降低设计成本。
基于频响函数的汽轮机叶片-转子系统稳定性的预测方法的具体方案如下:
基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,建立考虑叶片柔性和围带摩擦阻尼的汽轮机叶片-转子系统运动微分方程,利用其模态正交性进行解耦,求解叶片-转子系统频响函数,计算系统特征乘子,进而利用特征乘子与单位圆的关系判定系统稳定性。该预测方法,丰富非线性转子动力学稳定性预测研究理论,更为在设计阶段准确预测稳定性以实现大型汽轮机组结构优化设计提供方法,可有效避免设计阶段无法准确预测机组稳定性造成的“先天不足”,降低相关设计成本。
进一步地,上述预测方法,包括如下内容:
1)建立考虑弹性转轴与柔性叶片的变形耦合非线性和流体激励诱发的围带摩擦阻尼的激励依赖非线性的大型汽轮机叶片-转子系统运动微分方程;
2)获得不同运行工况下系统运动微分方程中的流体激励Ff(t),计算转子系统运动微分方程中的围带摩擦阻尼CS(F,t);
3)根据步骤2)获得的流体激励Ff(t)、围带摩擦阻尼CS(F,t),利用模态正交法解耦系统运动微分方程,求解各运行工况下综合考虑转轴系-叶片变形耦合和激励依赖围带阻尼激励依赖非线性、轴系叶片变形耦合影响的大型汽轮机叶片-转子系统频响函数;
4)通过求解叶片-转子系统频响函数计算系统特征乘子μF,并根据系统特征乘子μF与单位圆的关系判定系统在不同运行工况下的稳定性;
5)通过上述步骤2)-4)扫描计算所有运行工况,获得叶片-转子系统在各运行工况下的稳定性极限阈值,并将其组合绘制叶片-转子系统流致失稳三维极限曲面。
进一步地,所述步骤1)中所建立的系统运动微分方程是考虑围带摩擦阻尼的,运动微分方程表示为:
C(F,t)=CR+CO+CS(F,t),K=KR+KO (2)
式中,F(t)为系统的激励载荷矢量;Q(t)为系统的位移矢量;M、C和K分别为系统的质量、阻尼和刚度矩阵,其中,刚度矩阵中考虑了系统的结构刚度KR和支承刚度K0;阻尼项中考虑了系统的结构阻尼CR、支承阻尼C0及围带摩擦阻尼CS(F,t),且
F(t)=Ff(t)+Fm(t) (4)
其中,Ff(t)为湿蒸汽非平衡凝结流场产生的多频流体激励,Fm(t)为转子的转速激励、由不对中和不平衡等产生的机械激励。
进一步地,所述步骤1)所建立的系统运动微分方程中需辨识获得系统的质量M、结构阻尼CR、结构刚度KR、支承阻尼C0和支承刚度K0,具体过程如下:
1-1)根据汽轮机组转轴的长度和直径、支承跨度、叶片参数(安装角、扭转角、展弦比)以及围带参数(接触角、接触间隙),利用三维建模软件建立叶片围带接触角β、围带间初始间隙e0的围带-叶片-转子系统几何模型,利用有限元分析软件将系统离散为三维实体单元,得到式(1)-(2)中的系统质量M、结构阻尼CR和结构刚度矩阵KR以及单元形函数N;
1-2)将转子两端轴承支承近似为刚性支承,选择轴承参数(类型、载荷、润滑油型号和转速等)计算支承阻尼C0和支承刚度K0。
进一步地,所述步骤2)中利用数值方法获得系统运动微分方程中的流体激励Ff(t),具体如下:
2-1)在计算流体软件中建立叶-转子系统片流体计算域并划分网格,根据不同运行工况(进气速度V、进气角η)设置计算边界条件和转速Ω,选择适合的流体属性和湍流模型,进行非定常流场分析,获得叶片表面瞬态分布压力pi(x,y,z,t);
2-2)基于反距离加权插值法和虚功原理计算作用于叶片-转子系统三维实体单元的流体激励Ff(t):
Ff(t)=NTdF(t)T (5)
式中,pi(x,y,z,t)和di分别为流体计算域节点i(离散点)上的瞬态压力及其距叶片结构节点j(插值点)的距离;n为离散点的数量;N为叶片六面体单元型函数。
进一步地,所述步骤2)中围带摩擦阻尼CS(F,t)的计算步骤如下:
基于FFT多谐波平衡法,将系统运动微分方程中的复杂多频激励项F(t)转换为式(7)所示不同频率成分谐波激励分量的叠加,假设每个谐波分量均激起叶片-转子系统的简谐振动,则在第h阶谐波分量激励下,第i个叶片围带接触面上任一点在惯性坐标系(O-XYZ)下的振动位移响应qh i(t)可写为式(8):
δ=ln(A1/A2) (10)
式中,是系统模态质量;ωnr为忽略围带摩擦阻尼时系统的固有频率,δ为对数减幅率,A1和A2t和(t+T)时刻系统最大振动位移响应节点的振动位移的幅值,其中T为主振周期;将各阶谐波分量下的等效摩擦阻尼叠加即得到多频激励F(t)下的围带摩擦阻尼为
进一步地,所述步骤3)的具体过程如下:
根据步骤2)中求得的流体激励Ff(t)和围带摩擦阻尼CS(F,t)代入系统运动微分方程中,即得综合考虑激励依赖围带阻尼、转轴-叶片变形耦合影响的大型汽轮机叶片-转子系统的振动方程,其特征方程含有多(k)个特征根,即k个固有频率ωn1,ωn2,…ωnr…ωnk和相应的振型根据振型正交性,有:
据此,解耦后的系统振动方程可写为:
则第r阶模态频响函数H(ω)r用实部和虚部形式表示为:
H(ω)r=ReH(ω)r+ImH(ω)r (14)
其中
式中,ω为系统激励角频率,ωnr为系统第r阶模态固有频率,阻尼比ξr=cr/cfr=cr/2mrωnr为系统阻尼cr与第r阶临界阻尼系数cfr之比;定义λr=ω/ωnr为对应于第r阶模态的频率比,代入式(13)-(16)中得大型多自由度汽轮机叶片-转子系统第r阶模态频响函数(H(ω)r)的实部和虚部可分别写为:
进一步地,所述步骤4)中系统特征乘子μF求解过程如下:
当系统以接近第r阶模态固有频率ωnr的颤振频率ωcr振动时,系统运动微分方程式(13)的特征方程可写为:
det[I-diagΛrH(ωcr)]=0 (r=1,2,3,...,k) (19)
式中,I为k×k阶单位矩阵,diagΛr、H(ωcr)为k×k阶对角阵。则第r阶模态的特征根为:
Λr=[H(ωcr)]-1=[H(ω)r]-1 (20)
据矩阵运算法则,复数矩阵H(ω)r的逆矩阵可写为:
[H(ω)r]-1=(A+iB)-1
=(A+BA-1B)-1-iA-1B(A+BA-1B)-1 (21)
式中A和B代表频响函数的实部ReH(ω)r与虚部ImH(ω)r(式(18));将式(17)-(18)代入式(20)-(21)中,即可计算获得系统第r阶模态特征根;
系统特征乘子(μF)为:
μF=exp(ΛrT) (22)
其中,T为主振周期,通常情况下T=2πΩ/60。将第r阶模态特征根Λr代入(22)中即可获得机叶片-转子系统的特征乘子。
进一步地,所述步骤4)中系统稳定性的判定标准如下:
a)临界失稳:|μF|落在单位圆上;
b)Hopf分岔失稳:|μF|溢出单位圆的方式为共扼复数;
c)1周期分岔失稳:|μF|溢出单位圆的方式为正实数;
d)倍周期分岔失稳:|μF|溢出单位圆的方式为负实数;
若系统稳定,则返回步骤2),进行下一组运行工况下系统稳定性判定;若系统失稳,则在三维坐标系下标记当前工况点。
进一步地,所述下一组运行工况参数为进气速度Vi+ΔV、进气角ηi+Δη和转速Ωi+ΔΩ,其中,ΔV、Δη和ΔΩ分别为各当前工况参数Vi、ηi和Ωi相应的微增量。
与现有技术相比,本发明的有益效果是:
1)本发明通过预测方法的给出,解决了传统方法无法综合考虑大型机组转子系统因长轴系(弹性)、大柔性叶片及围带等结构特点及其与流场间的流固耦合作用导致的强激励(含转速)依赖非线性对机组稳定性的影响的缺陷,有效提高了多自由度非线性叶片-转子系统稳定性预测的精度。
2)本发明不仅丰富非线性转子动力学稳定性预测研究理论,更为在设计阶段准确预测稳定性以实现大型汽轮机组结构优化设计提供方法,可有效避免设计阶段无法准确预测机组稳定性造成的“先天不足”,将大幅提高设计效率降低设计成本。
附图说明
构成本申请的一部分的说明书附图用来提供对本申请的进一步理解,本申请的示意性实施例及其说明用于解释本申请,并不构成对本申请的不当限定。
图1为本发明的基于频响函数的汽轮机叶片-转子系统稳定性预测方法流程图;
图2为本发明的一个具体实施例的某1000MW汽轮机叶片-转子系统示意图;
图3为本发明的一个具体实施例的末级叶片流体计算域示意图;
图4为本发明的一个具体实施例的末级叶片流场表面分布压力示意图;
图5为本发明的一个具体实施例的三维稳定性极限曲面图。
具体实施方式
应该指出,以下详细说明都是例示性的,旨在对本申请提供进一步的说明。除非另有指明,本文使用的所有技术和科学术语具有与本申请所属技术领域的普通技术人员通常理解的相同含义。
需要注意的是,这里所使用的术语仅是为了描述具体实施方式,而非意图限制根据本申请的示例性实施方式。如在这里所使用的,除非上下文另外明确指出,否则单数形式也意图包括复数形式,此外,还应当理解的是,当在本说明书中使用术语“包含”和/或“包括”时,其指明存在特征、步骤、操作、器件、组件和/或它们的组合。
正如背景技术所介绍的,现有技术中存在的不足,为了解决如上的技术问题,本申请提出了基于频响函数的汽轮机叶片-转子系统稳定性的预测方法。
本申请的一种典型的实施方式中,如图2所示,基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,建立考虑叶片柔性和围带摩擦阻尼的汽轮机叶片-转子系统运动微分方程,利用其模态正交性进行解耦,求解系统频响函数,计算系统特征乘子,进而利用特征乘子与单位圆的关系判定系统稳定性。
上述方法,如图1所示,包括如下内容:
1)以某1000MW汽轮机组为具体实施例,如图2所示,建立考虑弹性转轴与柔性叶片的变形耦合非线性和流体激励诱发的围带摩擦阻尼激励依赖非线性的大型汽轮机叶片-转子系统运动微分方程:
C(F,t)=CR+CO+CS(F,t),K=KR+KO (2)
式中,F(t)为系统的激励载荷矢量;Q(t)为位移矢量,M、C、K分别为系统质量、阻尼和刚度矩阵;其中,刚度项中考虑了系统结构刚度KR和支承刚度K0;阻尼项中考虑了系统结构阻尼CR、支承阻尼C0及围带摩擦阻尼Cs(F,t),且
F(t)=Ff(t)+Fm(t) (4)
其中,Ff(t)为湿蒸汽非平衡凝结流场产生的流体激励,Fm(t)为转速激励和不对中、不平衡等转子故障产生的机械激励;
2)辨识系统运动微分方程(1)中的系统质量M、系统结构阻尼CR、系统结构刚度KR、系统支承阻尼C0、系统支承刚度K0,具体过程如下:
i)根据1000MW汽轮机转轴的长度和直径、支承跨度、叶片参数(安装角、扭转角、展弦比)以及围带参数(接触角、接触间隙),利用三维建模软件建立末级叶片围带接触角β、围带间初始间隙e0的围带-叶片-转子系统几何模型,利用有限元分析软件ANSYS Workbench中的网格划分功能,将系统离散为三维实体单元(Solid 185),得到式(1)-(2)中的系统质量、结构阻尼和结构刚度矩阵M、CR、KR以及单元形函数N;
3)利用数值方法获得系统运动微分方程(1)中的流体激励Ff(t),具体过程如下:
i)在计算流体软件CFX中建立末级叶片流体计算域(如图3)并划分网格,根据不同运行工况设置计算边界条件(进气速度V、进气角η和转速Ω),选择适合的流体属性和湍流模型,进行非定常流场分析,获得叶片表面瞬态分布压力pi(x,y,z,t),如图4所示;
ii)基于反距离加权插值法和虚功原理计算作用于叶片三维实体单元的流体激励Ff(t):
Ff(t)=NTdF(t)T (5)
式中,pi(x,y,z,t)和di分别为流体计算域节点i(离散点)上的瞬态压力及其距叶片结构节点j(插值点)的距离;n为离散点的数量,具体实施例中取n=4;N为叶片六面体单元形函数。
4)计算系统运动微分方程(1)中的围带摩擦阻尼CS(F,t),具体过程如下:
基于FFT多谐波平衡法,将式(1)中的复杂多频激励项F(t)转换为不同频率成分谐波激励分量的叠加(如式(7)),假设每个谐波分量均激起叶片的简谐振动,则在第h阶谐波激励下,第i个叶片围带l侧上任一点在惯性坐标系(O-XYZ)下的振动位移响应qh i(t)可写为式(8):
δ=ln(A1/A2) (10)
式中,是系统模态质量;ωnr为忽略围带摩擦阻尼时系统的固有频率,δ为对数减幅率,A1和A2分别为系统振动响应在初始时刻和截止时刻的幅值。将各阶谐波分量下的等效摩擦阻尼叠加即得到复杂多频激励F(t)下的围带摩擦阻尼为
工程中一般取前三阶谐波分量(n=3)即可满足要求。
步骤4)中,利用有限元分析软件ANSYS Workbench中的模态分析模块,定义支承刚度K0和支承阻尼Co,施加转速等效离心载荷激励,得到不同工况下系统前r阶模态质量mr和固有频率ωnr;利用ANSYS中的瞬态分析模块,通过在围带接触面添加接触的方式定义摩擦系数μ和接触刚度系数Kh;定义支承刚度K0和支承阻尼Co;施加转速和各阶谐波激励载荷,得到系统振动响应曲线;在系统的最大振动位移响应曲线读取出系统振动响应在初始时刻和截止时刻的幅值A1和A2,代入式(10)计算对数减幅率δ,再将δ、mr和ωnr代入式(9),即可求解出不同工况复杂多频激励诱发的围带摩擦阻尼Cs(F,t)。
5)求解综合考虑激励依赖围带阻尼、轴系叶片变形耦合影响的大型多自由度汽轮机叶片-转子系统频响函数,具体过程如下:
将步骤2)求得的M、CR、KR、C0、K0,步骤3)求得的流体激励Ff(t)和步骤4)求得的激励依赖围带摩擦阻尼CS(F,t)代入系统运动微分方程(式(1))中,即得综合考虑激励依赖围带阻尼、轴系叶片变形耦合影响的大型汽轮机叶片-转子系统的振动方程。其特征方程含有多(k)个特征根,即k个固有频率ωn1,ωn2,…ωnr…ωnk和相应的振型根据振型正交性,有:
据此,解耦后的系统振动方程可写为:
则第r阶模态的频响函数(H(ω)r)用实部和虚部表示为:
H(ω)r=ReH(ω)r+ImH(ω)r (14)
其中
式中,ω为系统激励角频率,ωnr为系统第r阶模态固有频率,阻尼比ξr=cr/cfr=cr/2mrωnr为系统阻尼cr与第r阶临界阻尼系数cfr之比;定义λr=ω/ωnr为对应于第r阶模态的频率比,代入式(13)-(16)中得大型多自由度汽轮机叶片-转子系统第r阶模态的频响函数(H(ω)r)的实部和虚部可写为:
6)求解大型多自由度汽轮机叶片-转子系统特征乘子μF,具体过程如下:
当系统以接近第r阶模态固有频率ωnr的颤振频率ωcr振动时,系统运动微分方程(式(13))的特征方程可写为:
det[Ι-diagΛrH(ωcr)]=0(r=1,2,3,...,k) (19)
式中,I为k×k阶单位矩阵,diagΛr、H(ωcr)为k×k阶对角阵。则第r阶模态的特征根为:
Λr=[H(ωcr)]-1=[H(ω)r]-1 (20)
据矩阵运算法则,复数矩阵(H(ω)r)的逆矩阵可写为:
[H(ω)r]-1=(A+iB)-1
=(A+BA-1B)-1-iA-1B(A+BA-1B)-1 (21)
式中,A、B代表频响函数的实部ReH(ω)r(式(17))与虚部ImH(ω)r(式(18))。将式(17)-(18)代入式(20)-(21)中,即可计算获得系统第r阶模态的特征根。
系统特征乘子(μF)为:
μF=exp(ΛrT) (22)
其中,T为主振周期,通常情况下T=2πΩ/60。将第r阶模态特征根Λr代入(22)中即可获得大型多自由度汽轮机叶片-转子系统特征乘子。
7)据系统特征乘子μF与单位圆的关系判定系统在任一运行工况(进气速度Vi、进气角ηi和转速Ωi)(i=1,2,3,…..n)下的稳定性,具体判定准则为:
a)临界失稳:|μF|落在单位圆上;
b)Hopf分岔失稳:|μF|溢出单位圆的方式为共扼复数;
c)1周期分岔失稳:|μF|溢出单位圆的方式为正实数;
d)倍周期分岔失稳:|μF|溢出单位圆的方式为负实数。
若系统稳定,则返回步骤(3),进行下一组运行工况(进气速度(Vi+ΔV)、进气角(ηi+Δη)、转速(Ωi+ΔΩ))(其中,ΔV、Δη和ΔΩ分别为各工况参数Vi、ηi和Ωi相应的微增量)下系统稳定性判定;若系统失稳,则在三维(进气速度、进气角和转速)坐标系下标记当前工况点(Vi、ηi和Ωi)。
8)重复步骤7),通过扫描所有运行工况(V、η、Ω),获得叶片-转子系统各运行工况的稳定性极限阈值,并将其组合绘制叶片-转子系统流致失稳三维(V-η-Ω)极限曲面,如图5所示。
以上所述仅为本申请的优选实施例而已,并不用于限制本申请,对于本领域的技术人员来说,本申请可以有各种更改和变化。凡在本申请的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本申请的保护范围之内。
Claims (8)
1.基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,建立考虑叶片柔性和围带摩擦阻尼的汽轮机叶片-转子系统运动微分方程,利用其模态正交性进行解耦,求解叶片-转子系统频响函数,计算系统特征乘子,进而利用特征乘子与单位圆的关系判定系统稳定性;
所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法包括如下步骤:
1)建立考虑弹性转轴与柔性叶片的变形耦合非线性和流体激励诱发的围带摩擦阻尼的激励依赖非线性的大型汽轮机叶片-转子系统运动微分方程;
2)获得不同运行工况下系统运动微分方程中的流体激励Ff(t),计算系统运动微分方程中的围带摩擦阻尼CS(F,t);
3)根据步骤2)获得的流体激励Ff(t)、围带摩擦阻尼CS(F,t),利用模态正交法解耦系统运动微分方程,求解各运行工况下综合考虑转轴-叶片变形耦合和围带阻尼激励依赖非线性影响的大型汽轮机叶片-转子系统频响函数;
4)通过求解叶片-转子系统频响函数计算系统特征乘子μF,并根据系统特征乘子μF与单位圆的关系判定系统在不同运行工况下的稳定性;
5)通过上述步骤2)-4)计算所有运行工况,获得叶片-转子系统在各运行工况下的稳定性极限阈值,并将其组合绘制叶片-转子系统流致失稳三维极限曲面;
所述步骤1)中所建立的系统运动微分方程是考虑围带摩擦阻尼的,运动微分方程表示为:
C(F,t)=CR+CO+CS(F,t),K=KR+KO (2)
式中,Q(t)为系统的位移矢量;M、C和K分别为系统的质量、阻尼和刚度矩阵,其中,刚度矩阵中考虑了系统的结构刚度KR和支承刚度KO;阻尼项中考虑了系统的结构阻尼CR、支承阻尼CO及围带摩擦阻尼CS(F,t),F(t)为系统的激励载荷矢量,且:
F(t)=Ff(t)+Fm(t) (4)
其中,Ff(t)为湿蒸汽非平衡凝结流场产生的多频流体激励,Fm(t)为转子的转速激励、由不对中和不平衡产生的机械激励。
2.根据权利要求1所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,所述步骤1)所建立的系统运动微分方程中需辨识获得系统的质量M、结构阻尼CR、结构刚度KR、支承阻尼CO和支承刚度KO,具体过程如下:
1-1)根据汽轮机组转轴的长度和直径、支承跨度、叶片参数以及围带参数,利用三维建模软件建立具有接触角β和初始间隙e0的围带的叶片-转子系统几何模型,利用有限元分析软件将系统离散为三维实体单元,得到式(1)-(2)中的系统质量M、结构阻尼CR和结构刚度矩阵KR;
1-2)将转子两端轴承支承近似为刚性支承,选择轴承参数计算支承阻尼CO和支承刚度KO;
叶片参数包括安装角、扭转角和展弦比,围带参数包括接触角和接触间隙,轴承参数包括类型、载荷、润滑油型号和转速。
3.根据权利要求1所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,所述步骤2)中利用数值方法获得系统运动微分方程中的流体激励Ff(t),具体如下:
2-1)在计算流体软件中建立叶片-转子系统流体计算域并划分网格,根据不同运行工况设置计算边界条件,运行工况包括进气速度V、进气角η和转速Ω,选择适合的流体属性和湍流模型,进行非定常流场分析,获得叶片表面瞬态分布压力pi(x,y,z,t);
2-2)基于反距离加权插值法和虚功原理计算作用于叶片-转子系统三维实体单元的流体激励Ff(t):
Ff(t)=NTdF(t)T (5)
式中,pi(x,y,z,t)和di分别为流体计算域节点i上的瞬态压力及其距叶片结构节点j的距离;i为离散点,j为插值点,n为离散点的数量;N为叶片六面体单元型函数。
4.根据权利要求1所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,所述步骤2)中围带摩擦阻尼CS(F,t)的计算步骤如下:
基于FFT多谐波平衡法,将系统运动微分方程中的复杂多频激励项F(t)转换为式(7)所示不同频率成分谐波激励分量的叠加,假设每个谐波分量均激起叶片-转子系统的简谐振动,则在第h阶谐波分量激励下,第i个叶片围带接触面上任一点在惯性坐标系(O-XYZ)下的振动位移响应qh i(t)可写为式(8):
δ=ln(A1/A2) (10)
式中,是系统模态质量;ωnr为忽略围带摩擦阻尼时系统的固有频率,δ为对数减幅率,A1和A2分别为t和(t+T)时刻系统最大振动位移响应节点的振动位移的幅值,其中T为主振周期;将各阶谐波分量下的等效摩擦阻尼叠加即得到多频激励F(t)下的围带摩擦阻尼为
5.根据权利要求1所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,所述步骤3)的具体过程如下:
根据步骤2)中求得的流体激励Ff(t)和围带摩擦阻尼Cs(F,t)代入系统运动微分方程中,即得综合考虑激励依赖围带阻尼、转轴-叶片变形耦合影响的大型汽轮机叶片-转子系统运动微分方程,其特征方程含有k个特征根,即k个固有频率ωn1,ωn2,…ωnr…ωnk和相应的振型根据振型正交性,有:
据此,解耦后的系统运动微分方程可写为:
则第r阶模态频响函数H(ω)r用实部和虚部形式表示为:
H(ω)r=ReH(ω)r+ImH(ω)r (14)
其中
式中,ω为系统激励角频率,ωnr为系统第r阶模态固有频率,阻尼比ξr=cr/cfr=cr/2mrωnr为系统阻尼cr与第r阶临界阻尼系数cfr之比;定义λr=ω/ωnr为对应于第r阶模态的频率比,代入式(13)-(16)中得大型多自由度汽轮机叶片-转子系统第r阶模态频响函数的实部和虚部可分别写为:
6.根据权利要求5所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,所述步骤4)中系统特征乘子μF求解过程如下:
当系统以接近第r阶模态固有频率ωnr的颤振频率ωcr振动时,系统运动微分方程式(13)的特征方程可写为:
det[I-diagΛrH(ωcr)]=0(r=1,2,3,...,k) (19)
式中,I为k×k阶单位矩阵,diagΛr、H(ωcr)为k×k阶对角阵,则第r阶模态的特征根为:
Λr=[H(ωcr)]-1=[H(ω)r]-1 (20)
据矩阵运算法则,复数矩阵H(ω)r的逆矩阵可写为:
[H(ω)r]-1=(A+iB)-1
=(A+BA-1B)-1-iA-1B(A+BA-1B)-1 (21)
式中A和B代表频响函数的实部ReH(ω)r与虚部ImH(ω)r;将式(17)-(18)代入式(20)-(21)中,即可计算获得系统第r阶模态特征根;
系统特征乘子(μF)为:
μF=exp(ΛrT) (22)
其中,T为主振周期,通常情况下T=2πΩ/60,将第r阶模态特征根Λr代入式(22)中即可获得叶片-转子系统的特征乘子。
7.根据权利要求1所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,所述步骤4)中系统稳定性的判定标准如下:
a)临界失稳:|μF|落在单位圆上;
b)Hopf分岔失稳:|μF|溢出单位圆的方式为共扼复数;
c)1周期分岔失稳:|μF|溢出单位圆的方式为正实数;
d)倍周期分岔失稳:|μF|溢出单位圆的方式为负实数;
若系统稳定,则返回步骤2),进行下一组运行工况下系统稳定性判定;若系统失稳,则在三维坐标系下标记当前工况点。
8.根据权利要求7所述的基于频响函数的汽轮机叶片-转子系统稳定性的预测方法,其特征在于,所述下一组运行工况参数为进气速度Vi+ΔV、进气角ηi+Δη和转速Ωi+ΔΩ,其中,ΔV、Δη和ΔΩ分别为各当前工况参数Vi、ηi和Ωi相应的微增量。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811643941.1A CN109657397B (zh) | 2018-12-29 | 2018-12-29 | 基于频响函数的汽轮机叶片-转子系统稳定性的预测方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811643941.1A CN109657397B (zh) | 2018-12-29 | 2018-12-29 | 基于频响函数的汽轮机叶片-转子系统稳定性的预测方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN109657397A CN109657397A (zh) | 2019-04-19 |
CN109657397B true CN109657397B (zh) | 2020-06-30 |
Family
ID=66118012
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811643941.1A Active CN109657397B (zh) | 2018-12-29 | 2018-12-29 | 基于频响函数的汽轮机叶片-转子系统稳定性的预测方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN109657397B (zh) |
Families Citing this family (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110929419B (zh) * | 2018-12-29 | 2021-08-13 | 山东大学 | 基于围带零阻尼的汽轮机转子系统失稳极限快速预测方法 |
CN110579312B (zh) * | 2019-10-17 | 2021-04-20 | 江苏方天电力技术有限公司 | 一种无试重旋转机械多轮盘轴系动平衡故障检测方法 |
CN111299668B (zh) * | 2019-12-11 | 2021-02-05 | 南京航空航天大学 | 一种不等齿距铣刀的齿间角确定方法 |
CN112179480B (zh) * | 2020-08-26 | 2021-11-19 | 西安交通大学 | 一种谐波叠加拟合叶尖定时振动参数识别的方法、系统及电子设备 |
CN112784379B (zh) * | 2021-03-09 | 2023-05-23 | 重庆邮电大学 | 一种叶轮转子系统的优化设计方法 |
CN114018480B (zh) * | 2021-11-24 | 2023-07-18 | 中国科学院重庆绿色智能技术研究院 | 一种大型旋转机械的转子不平衡故障的实时诊断方法 |
CN114580121B (zh) * | 2022-05-05 | 2022-08-16 | 西安航天动力研究所 | 基于有限元法的转子系统动特性计算方法、设备及介质 |
CN115306754B (zh) * | 2022-10-12 | 2023-02-17 | 中国航发四川燃气涡轮研究院 | 基于声阵列的轴流风扇气动失稳辨识方法 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102938003A (zh) * | 2012-10-17 | 2013-02-20 | 北京航空航天大学 | 一种叶轮机械计入错频的气动弹性稳定性数值预测方法 |
KR20160027454A (ko) * | 2014-08-29 | 2016-03-10 | 한국원자력연구원 | 습증기 침부식 모사장치 및 습증기 제조방법 |
-
2018
- 2018-12-29 CN CN201811643941.1A patent/CN109657397B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102938003A (zh) * | 2012-10-17 | 2013-02-20 | 北京航空航天大学 | 一种叶轮机械计入错频的气动弹性稳定性数值预测方法 |
KR20160027454A (ko) * | 2014-08-29 | 2016-03-10 | 한국원자력연구원 | 습증기 침부식 모사장치 및 습증기 제조방법 |
Non-Patent Citations (2)
Title |
---|
Wet Steam Nonequilibrium Condensation Flow-Induced Vibrations of a Nuclear Turbine Blade;Bing Guo等;《Journal of Engineering for Gas Turbines and Power》;20180731;第140卷(第11期);全文 * |
核电汽轮机长叶片与轴振动耦合特性研究;乔鹏等;《自动化仪表》;20170131;第38卷(第1期);全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN109657397A (zh) | 2019-04-19 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN109657397B (zh) | 基于频响函数的汽轮机叶片-转子系统稳定性的预测方法 | |
CN110929419B (zh) | 基于围带零阻尼的汽轮机转子系统失稳极限快速预测方法 | |
CN111062177B (zh) | 一种基于围带阻尼的汽轮机转子系统稳定性动态优化方法 | |
Petrov | Analysis of flutter-induced limit cycle oscillations in gas-turbine structures with friction, gap, and other nonlinear contact interfaces | |
Fuhrer et al. | On the impact of simulation approaches on the predicted aerodynamic damping of a low pressure steam turbine rotor | |
CN114154362B (zh) | 基于能量耗散的涡轮叶片缘板阻尼特性的仿真评估方法 | |
Pavlenko et al. | Using Computer Modeling and Artificial Neural Networks for Ensuring the Vibration Reliability of Rotors. | |
Jung et al. | Rotordynamic modelling and analysis of a radial inflow turbine rotor-bearing system | |
Mayorca et al. | Uncertainty of forced response numerical predictions of an industrial blisk: Comparison with experiments | |
Breńkacz | The experimental identification of the dynamic coefficients of two hydrodynamic journal bearings operating at constant rotational speed and under nonlinear conditions | |
Ning et al. | Blade forced response prediction for industrial gas turbines: Part 2—verification and application | |
Srinivasan | Flutter and resonant vibration characteristics of engine blades: An IGTI scholar paper | |
Kumar et al. | Aeroelastic aspects of axial compressor stage with self-recirculating casing treatment | |
Hosseini et al. | Effect of scaling of blade row sectors on the prediction of aerodynamic forcing in a highly-loaded transonic turbine stage | |
JP6807950B2 (ja) | 軸流ターボ機械のブレードを輪郭付けするための方法 | |
Willeke et al. | Reduced order modeling of mistuned bladed disks considering aerodynamic coupling and mode family interaction | |
Fruth et al. | Influence of the Blade Count Ratio on Aerodynamic Forcing: Part II—High Pressure Transonic Turbine | |
Repetski et al. | The sensitivity analysis for life estimation of turbine blades | |
Liao et al. | Finite element analysis and Lightweight design of hydro generator lower bracket | |
Ogbonnaya et al. | Analysis of gas turbine blade vibration due to random excitation | |
Filsinger et al. | Practical use of unsteady CFD and FEM forced response calculation in the design of axial turbocharger turbines | |
Savchenko et al. | Effect of the contact surfaces orientation in the shrouded flanges and level of vibration excitation in the rotor blades on their vibration stress state | |
CN106227909A (zh) | 一种一次消除汽轮发电机组转子动态挠曲的方法 | |
He et al. | Use of the harmonic balance method for flexible aircraft engine rotors with nonlinear squeeze film dampers | |
Kaneko et al. | Effect of material damping of steam turbine vane on flutter suppression |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |