CN113177273A - 一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法 - Google Patents

Abstract

本发明公开了一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，该方法以不同运行负荷下的接触强度分析为基础，分别建立点—面接触、线—面接触或面—面接触的摩擦阻尼模型，并结合线性化等效原理和能量耗散原则，获得不同摩擦阻尼器结构的接触刚度、阻尼代理模型，并突破传统方法的局限性，考虑高阶谐波分量的影响，计算接触刚度、阻尼对不同谐波项的等效分量，构建了多工况无量纲等效系数数据库，在数值模型中调用数据库获取相应等效系数，进行接触刚度和阻尼分析计算，为旋转机械叶片摩擦阻尼器的工业减振分析提供代理模型数据库及快速等效求解方法。

Description

一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解 方法

技术领域

本发明属于旋转机械叶片振动特性分析方法，特别涉及一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法。

背景技术

在高转速、大功率旋转机械中，其关键部件的减振一直是需要重点关注的问题，例如燃气轮机、压气机和透平叶片等在近共振点处的剧烈振动会产生过度应力而导致结构失效。摩擦阻尼器由于其结构简单、减振效果良好而被广泛应用于工业旋转机械叶片的减振设计中。摩擦阻尼器本质上是由一块金属组成的，常见的形式有阻尼围带、凸台拉金、松拉金、阻尼块等，其通过摩擦作用来耗散振动能量，能够有效抑制部件的振动水平，提高结构的安全可靠性。

在对摩擦阻尼器进行工程设计时，通常认为接触面达到理想的面-面接触状态，并基于基本负荷下受到的正压力载荷，将设计工况下的接触刚度和阻尼转化为接触面的等效刚度系数和等效阻尼系数，代入强非线性动力学方程求解，获得理想接触状态下***的摩擦阻尼特性。但为适应灵活调峰运行要求，机组需完成快速、深度变负荷运行，此时的工况严重偏离设计工况，并产生交变载荷，使接触面从理想的面-面接触变为点-面接触、线-面接触状态，且不同激振频率下阻尼接触面消耗的能量不同，设计工况下获得的等效刚度和等效阻尼不能适用在整个激振频率范围；而且不同激振力和正压力下，燃气轮机、压气机和透平叶片摩擦阻尼器、叶片叶根与轮缘之间接触界面的刚度和阻尼效果也会随着接触状态的改变发生变化；另一方面，目前的刚度阻尼计算方法主要考虑激励的一阶谐波项，而忽略非线性振动产生的超谐波项以及次谐波项，但当运动周期中滑移和粘滞状态并存且粘滞状态的权重较大时，就必须需要考虑高阶谐波分量对***振动响应的影响，这些项在共振频带内往往会对振动响应造成很大影响，造成较大的共振频带以及响应估计偏差。

发明内容

本发明的目的是提出一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，以不同运行负荷下的接触强度分析为基础，分别建立点-面接触、线-面接触或面-面接触的摩擦阻尼模型，并结合线性化等效原理和能量耗散原则，获得不同摩擦阻尼器结构的接触刚度、阻尼代理模型，并突破传统方法的局限性，考虑高阶谐波分量的影响，计算接触刚度、阻尼对不同谐波项的等效分量，增强等效方法估计对不同工况的鲁棒性，为旋转机械叶片-摩擦阻尼器(如燃气轮机、压气机和透平叶片摩擦阻尼器、叶片叶根与轮缘之间的接触摩擦等)的工业减振分析提供代理模型数据库及快速等效求解方法。

为达到上述目的，本发明采用以下技术方案：

一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，包括以下步骤：

步骤1，根据摩擦阻尼器材料获得相关参数：接触面间摩擦系数μ₀、泊松比ν、弹性模量E和剪切模量G，并通过分析获得不同转速工况下叶片摩擦接触面间的正压力载荷Q，得到接触面等效接触半径R和切向接触刚度k_c；

步骤2，根据运行负荷进行接触强度分析，根据接触面间的应力分布选择相应的摩擦模型来描述接触面间的摩擦阻尼特性，获得微动滑移、宏观滑移状态下的迟滞回线，建立接触面的摩擦阻尼函数模型f-s；

步骤3，采用一阶谐波平衡法，对相对位移和摩擦力进行傅里叶展开，根据线性化等效原则及能量耗散原理，计算得到摩擦接触面间的等效刚度系数K_e和等效阻尼系数C_e；

步骤4，重复步骤1～3对多种工况下的不同材料、接触特征以及压力载荷的摩擦阻尼器接触刚度和阻尼进行转化，并进行无量纲化处理，构建不同摩擦阻尼器结构的无量纲等效刚度系数

和无量纲等效阻尼系数

数据库；

步骤5，考虑高阶谐波因素，采用多谐波平衡法，根据线性化等效原理及步骤2的摩擦阻尼函数模型计算得到摩擦接触面间的等效刚度系数K_e和等效阻尼系数C_e；

步骤6，对多种工况下的不同材料、接触特征以及压力载荷的摩擦阻尼器接触刚度和阻尼进行转化，并进行无量纲化处理，构建高阶谐波下不同摩擦阻尼器结构的无量纲等效刚度系数

和无量纲等效阻尼系数

数据库；

步骤7，在数值分析模型中建立弹簧阻尼单元，根据模型不同材料、接触特征以及压力载荷，从数据库中调用相应的等效刚度系数和等效阻尼系数平均到弹簧阻尼单元各个方向上，进行接触刚度和阻尼分析计算。

本发明进一步的改进在于，步骤2中的摩擦阻尼函数模型包括：点-面摩擦接触模型、线-面摩擦接触模型和面-面摩擦接触模型。

本发明进一步的改进在于，步骤2中在快速升降负荷中摩擦接触面恰好接触，接触面应力分布成集中点状分布时，建立点-面摩擦接触模型，其描述函数取决于相对位移幅值大小：当接触面间的相对位移幅值D≤1.5D₀时，摩擦阻尼器处于微动滑移状态，此时，接触面间的摩擦迟滞回线为：

曲线ABC：

曲线CEA：

当D＞1.5D₀时，阻尼结构处于微动滑移宏观滑移共存阶段，此时，接触面间的摩擦迟滞回线为：

曲线A’B’C’：

曲线C’C”：

曲线C”E”A”：

曲线A”A’：f＝1(π+θ₀＜ωt≤2π)

其中

当接触面达到设计工况下的锁紧状态，接触面的应力分布不再发生变化时，建立面-面摩擦接触模型，其描述函数取决于相对位移幅值大小：D≤D₀时对应的无量纲摩擦迟滞回线为：

曲线ABCDE：

曲线EFGA：

D＞D₀时对应的无量纲摩擦迟滞回线为：

曲线DEF：

曲线FG：

曲线GABC：

曲线CD：

当接触面应力为线性分布时，建立线-面摩擦接触模型，根据力学分析过程，得到摩擦力f和阻尼器最右端位移D随滑移长度δ变化的表达式。

本发明进一步的改进在于，步骤3包括：

针对点-面摩擦接触模型和面-面摩擦接触模型：

步骤3.1，采用一阶谐波平衡法分别对相对位移s和摩擦力f进行傅里叶展开：

s＝D cosωt (1)

f＝f_k(D)cosωt+f_c(D)sinωt (2)

步骤3.2，根据傅立叶级数系数的运算，得到f_k(D)和f_c(D)：

步骤3.3，由线性化等效原理得到：

步骤3.4，将式(1)代入式(4)得到摩擦力的表达式为：

步骤3.5，将式(5)与式(2)相比，得到等效刚度系数和等效阻尼系数分别为：

步骤3.6，将式(3)代入式(6)，得到等效刚度系数和等效阻尼系数的表达式为：

针对线-面摩擦接触模型：

步骤3.1，根据能量耗散原理，得阻尼器在一个周期内的耗散能计算式：

步骤3.2，由力学过程分析，得到等效刚度系数和等效阻尼系数分别为：

本发明进一步的改进在于，步骤3.3中采用线性化等效原则，摩擦接触面简化为一个无质量弹簧阻尼***，接触面间的非线性摩擦力等效成弹性力和阻尼力的线性叠加。

本发明进一步的改进在于，步骤3.1中采用能量耗散原理，通过计算迟滞回线包围的面积来得到阻尼模型在一个振动周期内消耗的能量。

本发明进一步的改进在于，步骤4中将接触面的刚度和阻尼特性与阻尼器的几何参数分离，获得不同摩擦阻尼函数模型的等效系数无量纲形式；其中，

点-面摩擦接触模型的无量纲等效刚度系数

和无量纲等效阻尼系数

面-面摩擦接触模型的无量纲等效刚度系数

和无量纲等效阻尼系数

线—面摩擦接触模型的无量纲等效刚度系数

和无量纲等效阻尼系数

本发明进一步的改进在于，步骤5包括：

步骤5.1，给定谐波数N_h，采用多谐波平衡法分别对相对位移s和摩擦力f(s)进行N_h阶傅里叶级数展开：

步骤5.2，根据摩擦阻尼函数模型得到f(s)，并根据傅立叶变换得到式(9)右端各项系数：

步骤5.3，由线性化等效原理，将式(8)、(9)代入式(4)，通过比较对应项系数相等，计算得到各阶等效刚度系数k_en和等效阻尼系数c_en：

步骤5.4，将相对位移、摩擦力表示为矩阵形式，应用Newton-Raphson方法迭代求解s：

式(11)中，k表示当前迭代步数；R＝(ME²+CE+KI)s+f-F为方程残差项，当R＝0时收敛；

步骤5.5，将式(4)写成矩阵形式得到：

f＝(K_e+C_eE)s (12)

由式(10)得到等效刚度系数矩阵为：

等效阻尼系数矩阵为：

本发明进一步的改进在于，步骤6中对多谐波平衡法得到的等效系数无量纲化，得到各阶无量纲等效刚度系数

和无量纲等效阻尼系数

的表达式为：

式(13)中各项系数为：

无量纲等效刚度系数矩阵表示为：

无量纲等效阻尼系数矩阵表示为：

本发明至少具有如下有益的技术效果：

本发明提供的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，用参数化理论分析求解代替工程经验设计，获得机组不同运行负荷下的等效刚度和等效阻尼，并考虑了不同材料、不同正压力载荷和接触状况的影响；

进一步，提供了点—面接触、线—面接触和面—面接触三种摩擦阻尼模型，涵盖机组负荷变化过程中的基本接触特征，分析时可根据叶片—摩擦阻尼器的实际接触强度分析结果选择相应的摩擦阻尼模型，从而准确获得各状态下的摩擦迟滞回线及其描述函数。

进一步，在一阶谐波平衡法分析过程中采用线性化等效原则，将摩擦接触面简化为一个无质量弹簧阻尼***，由受力平衡分析将接触面间的非线性摩擦力等效成弹性力和阻尼力的线性叠加，从而将复杂的非线性力学问题转化为线性求解，获得点—面接触和面—面接触模型的等效刚度和等效系数，既能简化计算量又满足结果精度要求。

进一步，在一阶谐波平衡法分析过程中采用能量耗散原理，通过计算迟滞回线包围的面积可直接得到阻尼模型在一个振动周期内消耗的能量，并由力学分析得到线—面接触模型的等效刚度和等效系数。

进一步，将接触面的刚度和阻尼特性与阻尼器的几何参数分离，获得不同摩擦阻尼函数模型的等效系数无量纲形式，方便归纳比较，在无量纲形式的基础上乘以不同的材料参数即可得到对应材料的等效刚度系数和等效阻尼系数的表达式。

进一步，在燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法中考虑了高阶谐波分量的影响，计算接触刚度、阻尼对不同谐波项的等效分量，可提高摩擦阻尼数值模型的分析精度，增强等效方法估计对不同工况的鲁棒性。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面对实施例或现有技术描述中所需要使用的附图做简单的介绍；显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来说，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明实施例的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法的流程图；

图2是本发明实施例中的点—面接触摩擦模型：(a)为点—面接触摩擦模型粘滞滑移接触示意图；(b)为切向摩擦力时域曲线；(c)微动滑移阶段摩擦力随位移变化迟滞回线；(d)共存状态下摩擦力随位移变化迟滞回线；(e)为点—面接触摩擦模型的无量纲等效刚度系数曲线；(f)为点—面接触摩擦模型的无量纲等效阻尼系数曲线。

图3是本发明实施例中的微动滑移状态下面—面摩擦接触模型：(a)为面—面摩擦模型粘滞滑移接触示意图；(b)为切向摩擦时域曲线；(c)微动滑移阶段摩擦力随位移变化迟滞回线。

图4是本发明实施例中的微动滑移宏观滑移共存状态下面—面摩擦接触模型：(a)为面—面摩擦模型粘滞滑移接触示意图；(b)为切向摩擦时域曲线；(c)共存状态下摩擦力随位移变化迟滞回线；(d)为面—面摩擦模型的无量纲等效刚度系数曲线；(e)为面—面摩擦模型的无量纲等效阻尼系数曲线。

图5是本发明实施例中的微动滑移宏观滑移共存状态下线—面摩擦接触模型：(a)为线—面摩擦接触模型示意图；(b)为初始加载状态下滑动区域和滞止区域分布；(c)为切向力从最大值减小过程中滑动区域和滞止区域分布；(d)为切向力从最小值增大过程滑动区域和滞止区域分布；(e)为线—面摩擦接触模型迟滞回线。

图6是多谐波平衡法求解等效刚度和等效阻尼系数的流程示意图。

图中：1——球面摩擦副；2——平面摩擦副；3——粘滞区；4——滑移区；5——摩擦应力分布。

具体实施方式

为使本发明实施例的目的、技术效果及技术方案更加清楚，下面结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述；显然，所描述的实施例是本发明一部分实施例。基于本发明公开的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的其它实施例，都应属于本发明保护的范围。

请参阅图1所示，本发明提供的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其主要步骤为：

步骤1，根据摩擦阻尼器材料获得相关参数：接触面间摩擦系数μ₀、泊松比v、弹性模量E和剪切模量G，并通过分析获得不同转速工况下摩擦接触面间的正压力载荷Q，得到接触面等效接触半径R和切向接触刚度k_c；

步骤2，根据运行负荷进行接触强度分析，根据接触面间的应力分布选择相应的摩擦模型来描述接触面间的摩擦阻尼特性，获得微动滑移、宏观滑移状态下的迟滞回线，建立接触面的摩擦阻尼函数模型f-s；

步骤3，采用一阶谐波平衡法，对相对位移和摩擦力进行傅里叶展开，根据线性化等效原则及能量耗散原理，计算得到摩擦接触面间的等效刚度系数K_e和等效阻尼系数C_e；

步骤4，重复步骤1～3对多种工况下的不同材料、接触特征以及压力载荷的摩擦阻尼器接触刚度和阻尼进行转化，并进行无量纲化处理，构建不同摩擦阻尼器结构的无量纲等效刚度系数

和无量纲等效阻尼系数

数据库；

步骤5，考虑高阶谐波因素，采用多谐波平衡法，根据线性化等效原理及步骤2的摩擦阻尼函数模型计算得到摩擦接触面间的等效刚度系数K_e和等效阻尼系数C_e；

步骤6，对多种工况下的不同材料、接触特征以及压力载荷的摩擦阻尼器接触刚度和阻尼进行转化，并进行无量纲化处理，构建高阶谐波下不同摩擦阻尼器结构的无量纲等效刚度系数

和无量纲等效阻尼系数

数据库；

步骤7，在数值分析模型中建立弹簧阻尼单元，根据模型不同材料、接触特征以及压力载荷，从数据库中调用相应的等效刚度系数和等效阻尼系数平均到弹簧阻尼单元各个方向上，进行接触刚度和阻尼分析计算。

请参阅图2所示，在快速升降负荷中摩擦接触面恰好接触，接触面应力分布成集中点状分布时，建立点—面摩擦接触模型，其模型示意图如图2(a)所示，图2(b)为切向摩擦力时域曲线。其描述函数取决于相对位移幅值大小：规定当相对运动位移大于1.5倍的临界运动位移D₀时，接触面完全滑移，摩擦力达到最大值。当相对运动位移幅值D≤1.5D₀时，接触面不发生整体滑动，初始加载过程为曲线OAF，达到稳定状态后，曲线CEA代表激振力的加载过程，而曲线ABC则代表卸载过程，如图2(c)中所示。当相对运动位移幅值D>1.5D₀时，接触面发生整体滑动，初始加载过程为曲线OF，达到稳定状态后，曲线C”E”A”A’代表加载过程，A’B’C’C”代表卸载过程，如图2(d)中所示。

由弹性接触理论，对于两种材料相同的点—面接触对，接触半宽R可由式(2-1)计算：

点—面摩擦接触对在切向力作用下的切向接触刚度：

(2-2)式(2-2)中：r为球面半径，N为正压力，ν为泊松比，E为弹性模量，G为剪切模量，且对于各向同性材料，

定义点—面接触模型的临界相对位移幅值为D₀＝μ₀Q/k_c，无量纲相对运动位移幅值

无量纲摩擦力

无量纲微动滑移阶段摩擦力幅值

当接触面间的相对位移幅值D≤1.5D₀时，摩擦阻尼器处于微动滑移状态，此时，接触面间的无量纲摩擦力和相对位移表达为：

曲线ABC：

曲线CEA：

当D＞1.5D₀时，阻尼结构处于微动滑移宏观滑移共存阶段，此时对应的具有宏观滑移情况下无量纲形式的迟滞曲线表达式为：

曲线A’B’C’：

曲线C’C”：

曲线C”E”A”：

曲线A”A’：

其中

根据一阶谐波平衡法，当***受到简谐激励时，***的响应也是简谐变化的，且频率与激振力的频率相同，定义相对运动位移及其无量纲形式为：

接着，对摩擦力f进行一阶傅立叶展开可得：

f＝f_k(D)cosωt+f_c(D)sinωt (2-9)

根据傅立叶级数系数的运算，可以得到f_k(D)和f_c(D)：

根据线性化等效原理，将摩擦接触面简化为一个无质量弹簧阻尼***，且接触面间的非线性摩擦力可以由弹性力和阻尼力线性叠加来计算，即：

其中：K_e为摩擦接触面间等效刚度系数/N·m^-1，C_e为摩擦接触面间等效阻尼系数/N·s·m^-1；。

得到摩擦力的表达式为：

将式(2-9)与式(2-11)相比，可到等效刚度系数和等效阻尼系数分别为：

对式(2-12)进行无量纲化处理分别得到：

将摩擦迟滞曲线(2-3)～(2-8)代入式(2-13)进行积分运算即可得到无量纲等效刚度系数

和无量纲等效阻尼系数

分别如图2(e)和图2(f)所示。

请参阅图3所示，当接触面达到设计工况下的锁紧状态，接触面的应力分布不再发生变化时，建立面—面摩擦接触模型：当阻尼器处于仅有微动滑移状态时，面—面摩擦接触模型的滑移区域分布示意图和切向摩擦力变化曲线分别如图3(a)和图3(b)所示，单个周期内接触面间的迟滞曲线如图3(c)中所示。

由弹性接触理论，接触半宽R由式(3-1)计算得到：

其中L为圆柱长度，Q₀为接触中心的正压力。，定义面—面接触模型的临界相对位移幅值为

其描述函数取决于相对位移幅值大小：当相对运动位移幅值D≤D₀时，接触面不发生整体滑动，初始加载过程为曲线OA，达到稳定状态后，曲线ABCDE代表激振力的加载过程，而曲线EFGA则代表卸载过程。设接触区域的宽度为2a，则D≤D₀时对应的无量纲摩擦迟滞回线表达式如图3(c)中所示为：

曲线ABCDE：

曲线EFGA：

由一阶谐波平衡法计算得到无量纲的等效阻尼系数及等效刚度系数

分别为：

请参阅图4所示，当阻尼器处于微动滑移和宏观滑移两者共存状态时，面—面模型的滑移区域分布示意图和切向摩擦力变化曲线分别如图4(a)和图4(b)所示：当相对运动位移幅值D>D₀时，接触面发生整体滑动，初始加载过程为曲线OH，达到稳定状态后，曲线GABCD代表加载过程，DEFG代表卸载过程。则D>D₀时对应的无量纲摩擦迟滞回线表达式如图4(c)中所示为：

曲线DEF：

曲线FG：

曲线GABC：

曲线CD：

其中：

由一阶谐波平衡法计算得到无量纲的等效阻尼系数及等效刚度系数

分别为：

综合上述

和

两种滑移情况，得到面—面接触模型的无量纲等效刚度系数和等效阻尼系数曲线如图4(d)和图4(e)所示，表达式分别如下：

请参阅图5所示，当接触面应力为线性分布时，建立线—面摩擦接触模型并获得其摩擦迟滞回线。在线—面接触模型中，阻尼结构被简化成一个矩形板压于刚性平面之上，矩形板弹性模量为E；横截面积为A；长度为l；接触面之间的摩擦系数为μ，不论粘滞状态还是滑移状态，其大小均为定值。具体结构如图5(a)所示，垂直于矩形板表面上作用有分布的正压力q，并假设作用于矩形板表面上的正压力沿长度方向呈二次函数分布，其大小由下式确定：

此时，作用在阻尼器上的正压力Q可以表达为：

假设作用在矩形板上的切向力为F_a，呈正弦函数变化，其表达式为：

F_a=F sinωt

其中F_a为切向力幅值，ω为角速度，t为时间。

在线—面摩擦接触模型中摩擦力和相对位移的变化如下：初始阶段如图5(b)所示，切向力从0增大到最大值过程中，矩形板最右端从静止状态变为拉伸状态，并且随着切向力的增大，拉伸段的长度也在增长，当切向力达到最大值时，拉伸段长度达到最大值δ_a；如图5(c)所示，当切向力从最大值开始减小时，C区因外力减小而向x方向滑移，摩擦力方向与x方向相反；B区保持前一时刻的正应变状态，摩擦力方向和x方向一致；A区处于粘滞、无滑移、零应变状态，无摩擦力存在。矩形板最右端开始从拉伸状态变为压缩状态，阻尼器内部仍然处于拉伸状态，随着切向力进一步减小，压缩区域逐渐增大，直至滑移区域全部由拉伸状态变为压缩状态，此时，切向力也减小到最小值，压缩段长度达到最大值δ_d；如图5(d)所示，当力从最小值开始增大时，C区因外力增大而被拉伸，摩擦力方向和x方向一致，B区保持上一时刻的负应变状态.摩擦力方向和x方向相反；A区处于粘滞、无滑移、零应变状态，无摩擦力存在。矩形板最右端由压缩状态变为拉伸状态，并且随着切向力的增大，拉伸区域逐渐扩展，直至切向力达到最大时，整个滑移区域均为拉伸状态，其长度为δ_l。此时，切向力完成了一次循环。

由于阻尼结构在一个循环周期内消耗的能量和迟滞回线包围的面积相等，因此可以通过计算迟滞回线包围的面积来得到阻尼模型在一个振动周期内消耗的能量。设矩形板最右端的位移为u，拉升状态下滑移区域的长度为δ_i，压缩状态下滑移区域的长度为δ_d。根据力学分析过程，如图(5)所示，阻尼器在一个周期内的耗散能计算式如下式所示：

当F_a<μ₀Q时，阻尼器处于微动滑移阶段，即D≤1。在微动滑移单杆模型中，等效阻尼系数可以表达为：

等效刚度系数可以表达为：

此时，无量纲等效阻尼系数如下式所示：

无量纲等效刚度系数如下式所示：

当F_a>μ₀Q时，阻尼结构处于微动滑移宏观滑移共存状态，即D>1。此时阻尼结构在单个振动周期内消耗的能量由微动滑移和宏观滑移两部分组成，如下式所示：

W=W₁+W₂

阻尼结构由于纯微动滑移耗散的能量由下式计算获得：

W₁=W(δ_a)

阻尼结构由于纯宏观滑移消耗的能量由下式计算获得：

W₂=2(B-D)F_a

此时，在微动滑移宏观滑移共存阶段，线—面接触模型的无量纲等效阻尼系数表达为：

在微动滑移宏观滑移并存阶段，线—面接触模型等效刚度系数表达为：

其中

综合上述

和

两种滑移情况，得到线—面接触模型的无量纲等效刚度系数和等效阻尼系数分别如下：

结合图(6)所示，当运动周期中滑移和粘滞状态并存且粘滞状态的权重较大时，需要考虑高阶谐波分量对***振动响应的影响。将叶片摩擦阻尼器抽象为一个机械非线性振动***，其受到外力呈正弦函数变化，得到非线性动力学方程在时域上的一般表达式为：

将自变量s(t)用傅里叶级数展开至N_h阶如下式：

位移s(t)的一阶和二阶导数可表示如下：

接着，对摩擦力f(s)进行傅立叶展开可得：

式(5-5)右端各项系数可以通过f(s)的傅里叶变换得到：

将位移s(t)、摩擦力f(s)，激励F_a写成矩阵乘积的形式如下：

s(t)＝Φ^Ts，f(s)＝Φ^Tf，F_a＝Φ^TF (6-6)

其中，

Φ＝[1 cosωt sinωt…cos(N_hωt)sin(N_hωt)]^T

F＝[0 0 F…0 0]^T

位移s(t)的一阶和二阶导数可表示如下：

其中：

然后将式(5-6)代入式(5-1)中，并写成矩阵形式如下：

MΦ^TE²s+CΦ^TEs+KΦ^Ts+Φ^Tf＝Φ^TF (6-7)

由于式(5-7)对任意时间t都成立，所以各阶谐波系数完全相等，可以得到如下非线性代数方程组：

(ME²+CE+KI)s+f＝F

定义方程残差项为：

R＝(ME²+CE+KI)s+f-F

应用Newton-Raphson迭代求解式：

式中，k表示当前迭代步数。

与式(2-10)相同，根据线性化等效原理，将摩擦接触面简化为一个无质量弹簧阻尼***，且接触面间的非线性摩擦力可以由弹性力和阻尼力线性叠加来计算，即：

将式(5-2)、式(5-3)和式(5-5)代入式(5-8)得

通过比较对应项系数相等，计算得到各阶等效刚度系数k_en和等效阻尼系数c_en：

将式(5-8)写成矩阵形式为：

Φ^Tf＝K_eΦ^Ts+C_eΦ^TEs (6-9)

由于式(5-9)对任意时间t都成立，所以各阶谐波系数完全相等，可以得到如下非线性代数方程组：

f＝(K_e+C_eE)s (6-10)

式(5-10)中，等效刚度系数矩阵可表示为：

等效阻尼系数矩阵为：

同样，可得到各阶无量纲等效刚度系数

和无量纲等效阻尼系数

的表达式为：

式(5-11)中各项系数为：

无量纲等效刚度系数矩阵可表示为：

无量纲等效阻尼系数矩阵可表示为：

Claims (9)

1.一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，包括以下步骤：

步骤1，根据摩擦阻尼器材料获得相关参数：接触面间摩擦系数μ₀、泊松比v、弹性模量E和剪切模量G，并通过分析获得不同转速工况下叶片摩擦接触面间的正压力载荷Q，得到接触面等效接触半径R和切向接触刚度k_c；

步骤2，根据运行负荷进行接触强度分析，根据接触面间的应力分布选择相应的摩擦模型来描述接触面间的摩擦阻尼特性，获得微动滑移、宏观滑移状态下的迟滞回线，建立接触面的摩擦阻尼函数模型f-s；

步骤3，采用一阶谐波平衡法，对相对位移和摩擦力进行傅里叶展开，根据线性化等效原则及能量耗散原理，计算得到摩擦接触面间的等效刚度系数K_e和等效阻尼系数C_e；

步骤4，重复步骤1～3对多种工况下的不同材料、接触特征以及压力载荷的摩擦阻尼器接触刚度和阻尼进行转化，并进行无量纲化处理，构建不同摩擦阻尼器结构的无量纲等效刚度系数

和无量纲等效阻尼系数

数据库；

步骤5，考虑高阶谐波因素，采用多谐波平衡法，根据线性化等效原理及步骤2的摩擦阻尼函数模型计算得到摩擦接触面间的等效刚度系数K_e和等效阻尼系数C_e；

步骤6，对多种工况下的不同材料、接触特征以及压力载荷的摩擦阻尼器接触刚度和阻尼进行转化，并进行无量纲化处理，构建高阶谐波下不同摩擦阻尼器结构的无量纲等效刚度系数

和无量纲等效阻尼系数

数据库；

步骤7，在数值分析模型中建立弹簧阻尼单元，根据模型不同材料、接触特征以及压力载荷，从数据库中调用相应的等效刚度系数和等效阻尼系数平均到弹簧阻尼单元各个方向上，进行接触刚度和阻尼分析计算。

2.根据权利要求1所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤2中的摩擦阻尼函数模型包括：点—面摩擦接触模型、线—面摩擦接触模型和面—面摩擦接触模型。

3.根据权利要求2所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤2中在快速升降负荷中摩擦接触面恰好接触，接触面应力分布成集中点状分布时，建立点—面摩擦接触模型，其描述函数取决于相对位移幅值大小：当接触面间的相对位移幅值D≤1.5D₀时，摩擦阻尼器处于微动滑移状态，此时，接触面间的摩擦迟滞回线为：

曲线ABC：

曲线CEA：

当D＞1.5D₀时，阻尼结构处于微动滑移宏观滑移共存阶段，此时，接触面间的摩擦迟滞回线为：

曲线A’B’C’：

曲线C’C”：

曲线C”E”A”：

曲线A”A’：

其中

当接触面达到设计工况下的锁紧状态，接触面的应力分布不再发生变化时，建立面—面摩擦接触模型，其描述函数取决于相对位移幅值大小：D≤D₀时对应的无量纲摩擦迟滞回线为：

曲线ABCDE：

曲线EFGA：

D>D₀时对应的无量纲摩擦迟滞回线为：

曲线DEF：

曲线FG：

曲线GABC：

曲线CD：

当接触面应力为线性分布时，建立线—面摩擦接触模型，根据力学分析过程，得到摩擦力f和阻尼器最右端位移D随滑移长度δ变化的表达式。

4.根据权利要求3所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤3包括：

针对点—面摩擦接触模型和面—面摩擦接触模型：

步骤3.1，采用一阶谐波平衡法分别对相对位移s和摩擦力f进行傅里叶展开：

s＝D cosωt (1)

f＝f_k(D)cosωt+f_c(D)sinωt (2)

步骤3.2，根据傅立叶级数系数的运算，得到f_k(D)和f_c(D)：

步骤3.3，由线性化等效原理得到：

步骤3.4，将式(1)代入式(4)得到摩擦力的表达式为：

步骤3.5，将式(5)与式(2)相比，得到等效刚度系数和等效阻尼系数分别为：

步骤3.6，将式(3)代入式(6)，得到等效刚度系数和等效阻尼系数的表达式为：

针对线—面摩擦接触模型：

步骤3.1，根据能量耗散原理，得阻尼器在一个周期内的耗散能计算式：

步骤3.2，由力学过程分析，得到等效刚度系数和等效阻尼系数分别为：

5.根据权利要求4所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤3.3中采用线性化等效原则，摩擦接触面简化为一个无质量弹簧阻尼***，接触面间的非线性摩擦力等效成弹性力和阻尼力的线性叠加。

6.根据权利要求4所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤3.1中采用能量耗散原理，通过计算迟滞回线包围的面积来得到阻尼模型在一个振动周期内消耗的能量。

7.根据权利要求4所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤4中将接触面的刚度和阻尼特性与阻尼器的几何参数分离，获得不同摩擦阻尼函数模型的等效系数无量纲形式；其中，

点—面摩擦接触模型的无量纲等效刚度系数

和无量纲等效阻尼系数

面—面摩擦接触模型的无量纲等效刚度系数

和无量纲等效阻尼系数

线—面摩擦接触模型的无量纲等效刚度系数

和无量纲等效阻尼系数

8.根据权利要求7所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤5包括：

步骤5.1，给定谐波数N_h，采用多谐波平衡法分别对相对位移s和摩擦力f(s)进行N_h阶傅里叶级数展开：

步骤5.2，根据摩擦阻尼函数模型得到f(s)，并根据傅立叶变换得到式(9)右端各项系数：

步骤5.3，由线性化等效原理，将式(8)、(9)代入式(4)，通过比较对应项系数相等，计算得到各阶等效刚度系数k_en和等效阻尼系数c_en：

步骤5.4，将相对位移、摩擦力表示为矩阵形式，应用Newton-Raphson方法迭代求解s：

式(11)中，k表示当前迭代步数；R＝(ME²+CE+KI)s+f-F为方程残差项，当R＝0时收敛；

步骤5.5，将式(4)写成矩阵形式得到：

f＝(K_e+C_eE)s (12)

由式(10)得到等效刚度系数矩阵为：

等效阻尼系数矩阵为：

9.根据权利要求8所述的一种燃气轮机压气机和透平叶片接触刚度和阻尼分析求解方法，其特征在于，步骤6中对多谐波平衡法得到的等效系数无量纲化，得到各阶无量纲等效刚度系数

和无量纲等效阻尼系数

的表达式为：

式(13)中各项系数为：

无量纲等效刚度系数矩阵表示为：

无量纲等效阻尼系数矩阵表示为：

Images