CN110532732A - 一种叶片-机匣碰摩关系的确定方法 - Google Patents

Abstract

本发明涉及一种叶片‑机匣碰摩关系的确定方法，包括：获取扭型板的运动状态信息，得到叶片在机匣中运动时叶片的运动状态信息；获取圆柱壳的运动状态信息，得到叶片在机匣运动时机匣的运动状态信息；获取叶片‑机匣系统的整体运动状态信息；并根据扭型板和圆柱壳的几何关系确定叶片与机匣的碰摩关系；基于所述叶片与机匣的碰摩关系，确定所述叶片‑机匣是否满足生产条件。本发明采用的叶片‑机匣碰摩关系的确定方法考虑了叶片的安装角和扭角使获取的叶片与机匣的碰摩关系具有更高的精度。同时本发明采用具有弹性支撑边界的圆柱壳来模拟柔性机匣，得到叶片在柔性机匣中的碰摩关系。

Description

一种叶片-机匣碰摩关系的确定方法

技术领域

本发明涉及一种叶片-机匣碰摩关系的确定方法。

背景技术

叶片-机匣之间小的间隙可以减少航空发动机气体泄露，提高发动机效率。然而，这将过小的间隙将增加叶片-机匣碰摩发生的可能性。碰摩会产生复杂的振动响应，缩短叶片和机匣的使用寿命。因此人们越来越重视由叶片-机匣碰摩导致的叶片与机匣的振动损伤了。假设旋转叶片为固定在刚性圆盘上的悬臂梁，许多学者研究了叶片，圆柱壳特性。

但是以往的研究中不能考虑叶片沿弦向的振动。三维有限元叶片模型的建模需要足够的网格密度和质量，这可能导致摩擦仿真计算资源较高且不能定性地考虑叶尖沿弦向振动对摩擦相应的影响。且以往的研究中，都是对刚性机匣的研究，如采用不同的环形阵列单元模型、柔性环半解析模型、弯曲梁有限元模型和三维实体有限元模型。对柔性机匣的研究还是空白。

发明内容

(一)要解决的技术问题

为了解决现有技术的上述问题，本发明提供一种叶片-机匣碰摩关系的确定方法。

(二)技术方案

为了达到上述目的，本发明提供一种叶片-机匣碰摩关系的确定方法方法，包括以下步骤：

采用扭型板的运动状态模拟叶片在机匣中的运动状态，采用圆柱壳的运动状态模拟叶片在机匣中运动时的机匣运动状态，所述扭型板的一端与圆盘的盘周面呈角度连接，所述圆盘设定所述圆柱壳内部，且所述圆盘的中心在所述圆柱壳中心轴线上，包括以下步骤：

S1、基于扭型板在运动过程中的能量信息获取扭型板的运动状态信息，并根据所述扭型板的运动状态信息得到叶片在机匣中运动时叶片的运动状态信息；

基于弹性约束下圆柱壳的能量信息获取圆柱壳的运动状态信息；并根据所述圆柱壳的运动状态信息得到叶片在机匣运动时机匣的运动状态信息；

其中，所述圆柱壳的两侧分别沿轴向、周向和径向设定平动弹簧以及圆柱壳的旋转方向设定扭转弹簧；

S2、基于叶片在机匣中运动时叶片的运动状态信息和机匣的运动状态信息，获取叶片-机匣系统的整体运动状态信息；

并根据扭型板和圆柱壳的几何关系确定叶片与机匣的碰摩关系；

S3、基于所述叶片与机匣的碰摩关系，确定所述叶片-机匣是否满足生产条件。

优选的，在步骤S1前还包括：

依据扭型板的静态参数和动态参数，建立扭型板对应的能量方程，其中所述能量方程包括：扭型板的动能、扭型板应变能、扭型板的离心势能、法向摩擦碰摩力和切向碰摩力以及气动力在弯曲方向所做的功；

相应的，所述步骤S1包括：

基于所述能量方程依据哈密顿原理得到扭型板的运动方程；采用伽辽金法对所述扭型板运动方程进行离散化和模态截断处理以及进行动力学分析获取扭型板的运动状态信息，并根据所述扭型板的运动状态信息得到叶片在机匣中运动时叶片的运动状态信息。

优选的，其特征在于，

建立扭型板坐标系系统包括：扭型板全局坐标系OXYZ、扭型板旋转坐标系o'x'y'z'、扭型板局部坐标系oxyz；

基于所述扭型板坐标系系统，建立所述扭型板对应的动能方程：

其中，b为扭型板宽度；h为扭型板厚度；R_d为圆盘半径；θ为扭型板旋转角度；扭角β(x)是从β₀到β₁线性变化的:

其中为扭型板沿径向预扭率，β₀和β₁分别是扭型板根处和扭型板尖处的交错角；β'是扭型板尖相对于扭型板根的扭角；

基于板壳理论，建立扭型板对应的的应变能V₁方程和扭型板离心势能V₂方程式：

其中：

ε_y＝v′_y；

E为扭型板的杨氏模量；

υ扭型板的泊松比；

其中，L为扭型板的长度；f_c1和f_c2为扭型板在两个方向的离心力；

其中，所述f_c1和f_c2的表达式分别为：

其中，扭型板单位面积上的气动力F_e在弯曲方向所做的功W为：

其中F_n为法向摩擦碰摩力，F_t为切向碰摩力；

其中单位面积的气动力F_e的表达式为：

其中，F_e1为气动力幅值，n_e表示障碍物数量；

优选的，所述步骤S1包括：

基于扭型板的动能、扭型板应变能、扭型板的离心势能、法向摩擦碰摩力和切向碰摩力以及气动力在弯曲方向所做的功，根据哈密顿原理得到叶片的运动方程表达式为：

其中，V＝V₁+V₂；

采用伽辽金方法对扭型板运动方程进行离散，并基于正则坐标U_(m-1)N+n(t)、V_(m-1)N+n(t)、W_(m-1)N+n(t)，采用组合梁函数法，获取叶片的径向位移、弯曲位移和摆动位移；

其中，

cosh(α_m)cos(α_m)＝-1；

cosh(ψ_n)cos(ψ_n)＝1；

其中，M和N分别表示叶片在坐标x，y方向上的模态截断数；则，根据哈密顿原理对所述叶片运动方程进行变分处理和伽辽金截断处理获取所述叶片的整体运动微分方程为：

其中M_b、G_b、D_b和K_b分别是叶片的质量矩阵、科氏力矩阵、阻尼矩阵和刚度矩阵；

其中，η是正则坐标向量；F_b是外力矢量，包括非摩擦力矢量F_non-rub和摩擦力矢量F_rub，F_b＝F_non-rub+F_rub；

其中，K_b＝K_be+K_bs+K_bc+K_bacc；

K_be为结构刚度矩阵，K_bs为旋软化矩阵，K_bc为离心刚化矩阵和K_bacc为加速度导致的刚度矩阵。

优选的，在步骤S1前还包括：

依据弹性约束的圆柱壳和所述圆柱壳上的平动弹簧的静态参数和动态参数以及壳层理论，建立圆柱壳和弹簧的能量方程，其中，所述能量方程包括：圆柱壳的动能、圆柱壳的应变能、弹簧的变形能；并依据所述能量方程获取圆柱壳的能量信息；

相应的，所述步骤S1包括：基于瑞立法原则对所述圆柱壳的能量信息进行处理，获取弹性约束下圆柱壳的运动状态信息，并根据所述圆柱壳的运动状态信息得到叶片在机匣运动时机匣的运动状态信息。

优选的，其特征在于，

建立所述圆柱壳正交坐标系Oxzθ；

基于所述圆柱壳坐标系系统，建立所述圆柱壳对应的动能方程和应变能方程包括：

其中，R_c、L_c和h_c分别表示圆柱壳的中曲面的半径、长度和厚度；ρ_c、c_c和E_c分别为圆柱壳的材料密度、泊松比和杨氏模量；圆柱壳的轴向、周向、径向位移分别为u、v和w；在圆柱壳的两侧分别沿轴向、周向和径向引入平动弹簧k^u、k^v和k^w，旋转方向扭转弹簧k^θ；无量纲长度ξ＝x/L_c；

其中

其中

建立用任意周向波数p表示的圆柱壳位移方程：

其中，p为任意周向波数，ω_c为圆柱壳的固有频率；U_p(ξ)、V_p(ξ)和W_p(ξ)是相应的三个变形方向的轴向模态函数；

其中：

其中，a_i、b_i和c_i是待定系数；和是施密特正交化的特征正交多项式；TR是实际计算中截断的项数；

其中，根据公式可以得到多项式的正交集：

其中

对正交多项式进行归一化处理得到

则，特征正交多项式集满足

其中

基于用任意周向波数p表示的圆柱壳位移方程，所述圆柱壳的动能可表示为：

其中是无量纲频率参数；

其中有量纲频率与无量纲频率之间的转换关系为

基于用任意周向波数p表示的圆柱壳位移方程表达式，所述圆柱壳的应变能可表示为：

基于所述圆柱壳坐标系系统，以及圆柱壳和人工弹簧的静态参数和动态参数，获取圆柱壳的自由边界处的人工弹簧中的储存的变形能表达式：

其中，为x＝0处圆柱壳的边缘轴向的平动弹簧刚度；

为x＝0处圆柱壳的边缘周向的平动弹簧刚度；

为x＝0处圆柱壳的边缘径向的平动弹簧刚度；

为x＝0处圆柱壳的边缘平动弹簧扭簧刚度；

为x＝l处圆柱壳的边缘轴向的平动弹簧刚度；

为x＝l处圆柱壳的边缘周向的平动弹簧刚度；

为x＝l处圆柱壳的边缘径向的平动弹簧刚度；

为x＝l处圆柱壳的边缘平动弹簧扭簧刚度；

基于用任意周向波数p表示的圆柱壳位移方程表达式，所述圆柱壳的自由边界处的人工弹簧中的储存的变形能表达式为：

其中无量纲弹簧刚度可表示为：

基于所述圆柱壳的动能、圆柱壳的应变能和圆柱壳上弹簧的变形能获取圆柱壳的能量函数：

Γ＝T_c-U_c-U_s。

优选的，所述步骤S1包括：

基于圆柱壳的能量函数，依据瑞立法原则对圆柱壳的能量函数求导使Γ最小；

其中

获取圆柱壳的运动学方程：

其中K_cs是约束弹簧产生刚度矩阵；

K_c和M_c分别为圆柱壳的刚度矩阵和质量矩阵；

a_i、b_i和c_i为待定系数构成的Ritz向量；

其中，a_i＝[a₁,a₂,…a_TR]^T；b_i＝[b₁,b₂,…b_TR]^T；c_i＝[c₁,c₂,…c_TR]^T。

优选的，所述步骤S1之前还包括：

根据代表叶片的扭型板和代表机匣的圆柱壳在运动中的几何位置关系，以及代表叶片的扭型板和代表机匣的圆柱壳的动态参数和静态参数，建立扭型板和圆柱壳共同坐标系；

将所述扭型板上的板尖等分点设定为n；

根据扭型板和代表机匣的圆柱壳在运动中的几何位置关系，建立扭型板和圆柱壳之间的间隙函数的表达式：

其中，r_g为静止时刻扭型板轨道半径，r_g＝L+R_d；g₀为同轴圆柱壳与盘的平均间隙，g₀＝R_c-h_c/2-R_d-L；cⁱ _min是第i个板尖点与机匣内壁的初始最小间隙，假设每个板尖点与对应机匣点的初始最小间隙相同；为初始相位角；Ω是扭型板的转速；uⁱ _b为第i个板尖点的径向位移；uⁱ _c为第i圆柱壳点的径向位移；

获取在每个载荷步下，第i个板尖点与其对应的圆柱壳点之间的动态间隙函数的表达式为

基于所述板尖点与其对应的圆柱壳点之间的间隙确定板尖点和相应的圆柱壳点之间的碰摩条件；

其中，在每个载荷步，第i个板尖点的渗透深度δⁱ，

则

当穿透深度δⁱ第i个板尖点大于0，法向和切向碰摩力加载到板尖点i上，同时大小相等方向相反的切向碰摩力量加载到相应的圆柱壳i上，没有发生碰摩时δi＝0；

其中Fⁱ _n和Fⁱ _t分别为第i个板尖点的法向碰摩力和切向碰摩力；

μ是板尖和圆柱壳之间的摩擦系数。

优选的，其特征在于，所述步骤S1包括：根据公式(A)获取的扭型板运动状态信息，并根据所述扭型板的运动状态信息表示叶片在机匣中运动时叶片的运动状态信息；

根据公式(B)获取圆柱壳的运动状态信息，并根据所述圆柱壳的运动状态信息表示叶片在机匣中运动时机匣的运动状态信息；

所述步骤S2包括：将表示叶片在机匣中运动时叶片的运动状态信息的公式(A)和表示叶片在机匣中运动时机匣的运动状态信息的公式(B)

组集得到公式(C)，并根据公式(C)获取叶片-机匣系统的整体运动状态信息；

其中，叶片-机匣系统的整体运动方程中的x_b表示扭型板运动公式(A)中的η，代表叶片的正则坐标向量，叶片-机匣系统的整体运动方程中的x_c表示弹性约束下圆柱壳运动公式(B)中的向量；

基于叶片机匣系统的整体运动方程根据所述叶尖点与其对应的机匣点之间的间隙确定叶尖点和相应的机匣点之间的碰摩条件确定叶片与机匣的关系。

优选的，其特征在于，

所述扭型板上的板尖等分点n设定为21；则

对应与板尖点i的圆柱壳点i的x坐标为：

对应与板尖点i的圆柱点i的θ坐标为：

(三)有益效果

本发明的有益效果是：本发明采用的叶片-机匣碰摩关系的确定方法考虑了叶片的安装角和扭角使获取的叶片与机匣的碰摩关系具有更高的精度。同时本发明采用具有弹性支撑边界的圆柱壳来模拟柔性机匣，得到叶片在柔性机匣中的碰摩关系。

附图说明

图1为本发明叶片-机匣碰摩关系确定方法的流程图；

图2为本发明扭型板叶片示意图；

图3为本发明含弹性约束的圆柱壳结构示意图；

图4为本发明扭型板与圆柱壳碰摩示意图；

图5为本发明第二实施例中叶片的前三阶固有频率；

图6(a)为本发明第二实施例中不同的弹簧刚度下的无量纲频率；

图6(b)为本发明第二实施例中不同的弹簧刚度下的无量纲频率；

图6(c)为本发明第二实施例中不同的弹簧刚度下的无量纲频率；

图6(d)为本发明第二实施例中不同的弹簧刚度下的无量纲频率；

图7(a)为本发明第二实施例中机匣对称支承边界模态截断数TR的敏感性；

图7(b)为本发明第二实施例中机匣非对称支承边界模态截断数TR的敏感性；

图8为第二实施例中机匣的振型其中(a1)为现有技术中对称边界条件下的有限元模型，(a2)本发明对称边界条件下的模型，(b1)现有技术中非对称边界条件下的有限元模型，(b2)本发明非对称边界条件下的模型；

图9为现有技术中旋转叶片与柔性机匣有限元模型；

图10为本发明第二实施例中叶片-机匣碰摩的解析法求解流程图；

图11(a)为本发明碰摩点1径向位移响应；

图11(b)为本发明非碰摩点1'径向位移响应；

图12(a)为本发明叶尖点1径向位移响应；

图12(b)为本发明机匣点1径向位移响应；

图13(a1)为本发明碰摩中各叶尖点的碰摩力分布图；

图13(b1)为本发明碰摩中各尖点的侵入深度图；

图13(a2)为有限元模型得到的碰摩中各叶尖点的碰摩力分布图；

图13(b2)为有限元模型得到碰摩中各尖点的侵入深度图；

图14(a1)为本发明得到的叶尖径向位移；

图14(b1)为本发明得到的叶尖弯曲位移；

图14(a2)为有限元模型得到的叶尖径向位移；

图14(b2)为有限元模型得到的叶尖弯曲位移；

图15(a1)为本发明得到的机匣碰摩点径向位移；

图15(b1)为本发明得到的机匣非碰摩点径向位移；

图15(a2)为有限元模型得到的机匣碰摩点径向位移；

图15(b2)为有限元模型得到的机匣非碰摩点径向位移；

图16(a1)为现有模型得到的叶尖点碰摩力图；

图16(b1)为现有模型得到的叶尖点侵入深入图；

图16(c1)为现有模型得到的叶尖点弯曲位移波形图；

图16(d1)为现有模型得到的叶尖点弯曲位移频谱；

图16(e1)为现有模型得到的机匣点径向位移波形图；

图16(f1)为现有模型得到的机匣点径向位移频谱；

图16(a2)为本发明得到的叶尖点碰摩力图；

图16(b2)为本发明得到的叶尖点侵入深入图；

图16(c2)为本发明得到的叶尖点弯曲位移波形图；

图16(d2)为本发明得到的叶尖点弯曲位移频谱；

图16(e2)为本发明得到的机匣点径向位移波形图；

图16(f2)为本发明得到的机匣点径向位移频谱；

图17为依据本发明得到的不同时刻下的机匣径向变形情况；

图18为依据本发明得到的不同时刻下的机匣变形和机匣的节径分量。

具体实施方式

为了更好的解释本发明，以便于理解，下面结合附图，通过具体实施方式，对本发明作详细描述。

实施例一

步骤一：参见图2所示，建立坐标系系统，其中，OXYZ,o'x'y'z'和oxyz分别表示扭型板的全局坐标系、旋转坐标系和局部坐标系。o”uvw是一个与扭角有关的变化坐标系。平面o”vw平行于平面oyz，v轴和y轴之间的夹角为β(x)-β₀。R_d为圆盘半径；θ是一个与时间有关的旋转角度。对任意长度‘x’而言，扭角β(x)是从β₀到β₁线性变化的:

这里是扭型板沿径向预扭率，β₀和β₁分别是扭型板根处和扭型板尖处的交错角。β'是扭型板尖相对于扭型板根的扭角。

旋转坐标系o'x'y'z'与全局坐标系OXYZ的转换关系为:

由于扭型板扭角的存在，扭型板上任意一点Q在全局坐标系OXYZ中的位移矢量为:

其中

点Q的速度可以表示为：

旋转扭型板的动能为：

其中b为扭型板宽度，h为扭型板厚度。

基于板壳理论，考虑扭转角的影响，扭型板的应变能V₁可表示为:

其中

其中E为扭型板的杨氏模量；υ扭型板的泊松比；

扭型板离心势能V₂为：

其中，L是扭型板的长度；f_c1和f_c2扭型板在两个方向的离心力。f_c1和f_c2可以被表示为：

法向摩擦碰摩力F_n、切向碰摩力F_t和扭型板单位面积上的气动力F_e在弯曲方向所做的功W可定义为:

其中，单位面积的气动力F_e表达式为

这里F_e1为气动力幅值，F_e1＝0.003MPa，n_e表示障碍物数量。

由哈密顿原理得到扭型板的运动方程，其表达式为:

其中V＝V₁+V₂。

采用伽辽金方法对扭型板运动方程进行离散。引入正则坐标U_(m-1)N+n(t)、V_(m-1)N+n(t)、W_(m-1)N+n(t)，采用组合梁函数法，得到径向位移、弯曲位移和摆动位移:

其中M和N分别表示板叶片在x、y方向上的模态截断数。

cosh(α_m)cos(α_m)＝-1；

cosh(ψ_n)cos(ψ_n)＝1；

扭型板的整体运动微分方程可表示为:

其中M_b、G_b、D_b和K_b分别是质量矩阵、科氏力矩阵、阻尼矩阵和刚度矩阵；η是正则坐标向量；F_b是外力矢量，包括非摩擦力矢量F_non-rub和摩擦力矢量F_rub，F_b＝F_non-rub+F_rub；K_b＝K_be+K_bs+K_bc+K_bacc,这里K_be，K_bs，K_bc和K_bacc为结构刚度矩阵、旋软化矩阵、离心刚化矩阵和加速度导致的刚度矩阵。

参见图3所示，其中R_c、L_c和h_c分别表示圆柱壳的中曲面的半径、长度和厚度。ρ_c、c_c和E_c分别为圆柱壳的材料密度、泊松比和杨氏模量。Oxzθ为正交坐标系。圆柱壳的位移分别用u、v和w表示。在圆柱壳的两侧分别沿轴向、周向和径向引入平动弹簧k^u、k^v和k^w，旋转方向引入扭转弹簧k^θ。

为方便计算,引入无量纲长度ξ＝x/L_c，圆柱壳的动能可以表示为

圆柱壳的应变能可表示为:

其中应变ε_x、ε_θ和ε_xθ为：

根据Sanders’壳层理论，中曲面的应变和及中曲面曲率和τ⁰定义为:

圆柱壳的应力-应变关系可以表示为:

圆柱壳位移可以用任意周向波数p表示:

ω_c是机匣的固有频率；U_p(ξ)、V_p(ξ)和W_p(ξ)是相应的三个变形方向的轴向模态函数。它们可以满足不同的边界支撑条件。

其中a_i、b_i和c_i是未知系数；和是施密特正交化的特征正交多项式。TR是实际计算中截断的项数，随着TR值的增大，计算结果的收敛性和精度将得到提高。

给定一个多项式根据下面的公式可以得到多项式的正交集：

其中

对上面的正交多项式进行归一化

然后，特征正交多项式集满足

在用该方法构造的正交多项式族中，只要初始多项式满足圆柱壳的边界条件，那么多项式族中的其他元素也可以满足圆柱壳的边界条件。自由边界的第一个多项式可以表示为

将用任意周向波数p表示的圆柱壳位移公式代入圆柱壳的动能公式，圆柱壳的动能可表示为:

其中是无量纲频率参数，而有量纲频率与无量纲频率之间的转换关系为

将用任意周向波数p表示的圆柱壳位移公式代入圆柱壳应变能公式，圆柱壳的应变能公式表示为:

参见附图3，考虑到机匣的弹性约束在实际中更为合理，因此在自由边界处采用了人工弹簧。人工弹簧中储存的变形能可以表示为：

其中，为x＝0处圆柱壳的边缘轴向的平动弹簧刚度；

为x＝0处圆柱壳的边缘周向的平动弹簧刚度；

为x＝0处圆柱壳的边缘径向的平动弹簧刚度；

为x＝0处圆柱壳的边缘平动弹簧扭簧刚度；

为x＝l处圆柱壳的边缘轴向的平动弹簧刚度；

为x＝l处圆柱壳的边缘周向的平动弹簧刚度；

为x＝l处圆柱壳的边缘径向的平动弹簧刚度；

为x＝l处圆柱壳的边缘平动弹簧扭簧刚度。

将用任意周向波数p表示的圆柱壳位移公式代入人工弹簧中储存的变形能公式，则变形能公式表示为：

其中无量纲弹簧刚度可表示为：

圆柱壳的能量函数为：

Γ＝T_c-U_c-U_s

基于瑞利法的原则，通过对每个系数求导使Γ最小：

其中

圆柱壳动力学方程方程为:

其中K_cs是约束弹簧产生刚度矩阵，K_c和M_c分别为圆柱壳的刚度矩阵和质量矩阵。a_i、b_i和c_i为待定系数构成的Ritz向量，可以表示为：

a_i＝[a₁,a₂,…a_TR]^T,b_i＝[b₁,b₂,…b_TR]^T,c_i＝[c₁,c₂,…c_TR]^T；

本实施例中，考虑a_i，b_i和c_i,非奇异解，弹性约束下圆柱壳的频率方程为：

步骤二：由于叶片的动力学微分方程和圆柱壳的动力学方程都是矩阵方程，所以将叶片的动力学微分方程和圆柱壳动力学方程进行组集可得到叶片-机匣的整体运动方程。

其中，叶片-机匣系统的整体运动方程中的x_b表示叶片运动微分方程中的η代表叶片的正则坐标向量，叶片-机匣系统的整体运动方程中的x_c表示弹性约束下圆柱壳的动力学方程中的向量。

参见附图4所示，根据薄壁圆柱壳理论，圆柱壳内壁点的位移可用相应的中曲面点的位移代替。因此，通过判断圆柱壳和扭型板尖点与对应的中曲线点之间的位置关系，可以判断摩擦是否发生。用扭型板模拟叶片，用圆柱壳模拟机匣。每个时刻对应n个叶尖点的n个圆柱壳点的运动状态不同，因此可能会在扭型板与圆柱壳之间发生点碰或局部碰摩。点碰是指一个扭型板尖点与一个圆柱壳点之间的碰摩，局部碰摩是指两个或多个扭型板尖点与对应的圆柱壳点之间的碰摩。如附图4所示。附图4中O_c和O分别为静止圆柱壳的中心和圆盘的中心；r_g为静止时刻扭型板尖轨道半径，r_g＝L+R_d；g₀为同轴圆柱壳与圆盘的平均间隙，g₀＝R_c-h_c/2-R_d-L。cⁱ _min是第i个扭型板尖点与圆柱壳内壁的初始最小间隙，假设每个扭型板尖点与对应圆柱壳点的初始最小间隙相同。为初始相位角；Ω是扭型板的转速；uⁱ _b为第i个扭型板尖点的径向位移；uⁱ _c为第i圆柱壳点的径向位移。

根据几何关系，方程可以写成：

扭型板尖点i和圆柱壳之间的间隙函数gⁱ _gap可以被表示为:

由于薄壁圆柱壳主要集中在径向变形上，假设在摩擦过程中，只在与扭型板尖点i对应的圆柱壳i处产生径向位移。在每个载荷步下，第i个扭型板尖点与其对应的圆柱壳点之间的动态间隙函数的表达式为：

在每个载荷步，需要判断每个扭型板尖点和相应圆柱壳点之间的碰摩，第i个扭型板尖点的渗透深度δⁱ可以表示为:

当穿透深度δⁱ第i个扭型板尖点大于0，法向和切向碰摩力加载到扭型板尖点i上，同时大小相等方向相反的切向碰摩力量加载到相应的机匣点i上，没有发生碰摩时δi＝0。

其中和别为第i个扭型板尖点的法向碰摩力和切向碰摩力。μ是扭型板尖和圆柱壳之间的摩擦系数。

本实施例中，考虑到计算效率和精度，在后续的分析中，扭型板尖上的等分点设为21。

扭型板尖点对应的圆柱壳点位置坐标推导如下：初始时刻扭型板尖点11与机匣点相对，圆柱壳点11坐标为(L_c/2,R_c,0°)。推导对应与其他扭型板尖点的圆柱壳点坐标(L_i,R_c,θ_i)(i＝1,2,3,…,21)。由于安装角和扭角的存在，扭型板尖和水平位置之间的角度是β₀+β'(见图3)。对应与扭型板尖点i的圆柱壳点i的x坐标为：

对应与扭型板尖点i的圆柱壳点i的θ坐标为：

并根据扭型板和圆柱壳的几何关系确定叶片与机匣的碰摩关系。

步骤三、基于所述叶片与机匣的碰摩关系，确定所述叶片-机匣是否满足生产条件。

本实施例中采用的叶片-机匣碰摩关系的确定方法考虑了叶片的安装角和扭角使获取的叶片与机匣的碰摩关系具有更高的精度。同时本实施例采用具有弹性支撑边界的圆柱壳来模拟柔性机匣，得到叶片在柔性机匣中的碰摩关系。

实施例二

本实施例以悬臂板理论为基础，采用半解析法，用含安装角的扭板对旋转叶片进行数值模拟。将模型的固有频率与有限元模型的固有频率进行比较，验证本实施例中叶片动力学模型的正确性。利用ANSYS软件中的shell181单元建立了悬臂叶片的有限元模型。这个有限元模型有20×20个元素和441个节点，每个节点有6个自由度(DOFs)。旋转叶片的参数如表1所示。叶片在Ω＝10000rev/min下的收敛性结果在表2中列出。通过比较叶片的前三阶模态，模态截断阶数从M＝N＝3到M＝N＝6被检查。收敛分析表明，M＝N＝5对扭型叶片是合适的，最大误差约为2.5％(见表2)。叶片的前三阶固有频率(f_n1、f_n2和f_n3)见附图5。f_n1、f_n2和f_n3分别表示一阶弯曲模态、一阶扭转模态和二阶弯曲模态的频率。可以看出，叶片的前三个固有频率会随着转速的增加而增加。提出的模型和有限元模型的计算结果如表3所示，提出的模型的固有频率相对于有限元模型的最大误差为2.6％(见表3)。

表1.旋转扭型叶片的参数

表2.在转速Ω＝10000rev/min时不同模态截断数下的叶片的前三阶固有频率

注:()为从有限元模型得到的结果。

表3叶片的固有频率

本实施例为了验证所提出的弹性约束圆柱壳动力学模型的正确性，对参考文献中的一个实例进行了分析。圆柱壳的物理参数如下：L_c＝0.2m、R_c＝0.1m、h_c＝0.247×10^-3m、ρ_c＝2796kg/m³、c_c＝0.3和E_c＝71.02×10⁹N/m²。假设机匣两端具有均匀的边界条件，即和一组弹簧的刚度值从0，1，10⁵变化，其他弹簧的刚度值设为1。在不同弹簧刚度下无量纲频率如图6所示。

可以看出，本文的固有频率与参考文献的固有频率吻合较好。利用圆柱壳理论，对具有弹性约束的机匣模型进行了验证。无量纲频率随周向波数p的增加而减小，然后增大。无量纲刚度与其它无量纲刚度相比，对圆柱壳的频率有较大的影响。由图6(a)可知，在时圆柱壳的最小频率对应的周向波数p为7，在最小频率对应的周向波数p为8。这一现象表明，周向波数p随最低频率的变化而变化。径向弹簧刚度和扭簧刚度对机匣的固有频率有轻微的影响，如图6(c)和图6(d)所示。

为了进一步验证所提模型的准确性，并为叶片-机匣摩擦响应分析提供一定的参考，对下一节所使用的机匣固有频率进行了分析。机匣的参数如表4所示。一般情况下，机匣前端的柔性要大于后端的，导致机匣两端边界不一致。本节采用对称边界条件和非对称边界条件对模型进行验证，假设对称边界机匣两端均为k^u＝k^v＝k^w＝1×10⁷N/m,k^θ＝1×10⁵N·m/rad；非对称边界机匣前端为k^u＝k^v＝k^w＝1×10⁶N/m,k^θ＝1×10⁶N·m/rad，机匣后端为k^u＝k^v＝k^w＝1×10⁷N/m,k^θ＝1×10⁵N·m/rad。首先分析了模态截断TR个数的灵敏度，如图7所示。由图可知，随着模态截断TR数目的增加，两种支护边界的固有频率趋于稳定，因此，在接下来的分析中，将TR设为7。

表4.机匣参数

利用ANSYS软件中的shell181单元建立了机匣的有限元模型，并对提出的模型进行了验证。在FE模型中，机匣的圆周方向被划分为360个部分，因此在机匣的两端都有360个节点。机匣两端各节点有4个弹簧元件，包括3个平动弹簧元件和1个扭簧元件，扭簧沿切向运动。该模型与有限元模型的机匣固有频率比较结果如表5和表6所示。具有对称边界条件和非对称边界条件的机匣振型如图8所示

由表5和表6可以看出，PM(本发明提出的方法)与FEM(有限元模型)之间的最大误差为1.6％，进一步验证了该弹性支承模型的准确性。对称边界条件下机匣的最低阶频率对应的模态为(1，3)，而非对称边界条件下机匣的最低阶频率对应的模态为(1，4)非对称边界条件。与对称边界条件相比，具有非对称边界条件的机匣的固有频率会减小(见表6)

表5.在对称边界条件下提出的模型和有限元模型固有频率比较

表6.在非对称边界条件下提出的模型和有限元模型固有频率比较

本实施例中，旋转叶片与柔性机匣的有限元模型如图9所示。Shell181单元用于构建的含扭角和安装角的叶片(β₀＝30°，β＝5°)。采用181单元建立了具有对称边界条件的柔性机匣，并结合Combin14弹簧单元建立了弹簧。具体建模过程如下：机匣圆周方向等分360份，机匣两端的轴向方向等分4份，机匣中间轴向缝份20份以实现网格细化，机匣的单元总数是10080个。机匣弹性支承，在机匣两端的各节点设置4根弹簧(包括3根平移弹簧，1根扭转弹簧)，弹簧的另一端施加固定约束。叶片长度方向等分20个部分。同时将叶片弦向平均分为20个部分，保证叶尖处分布21个节点。叶片单元总数为400个，叶根节点全约束。几个关键的机匣点位置定义如下：对应于叶片初始位置叶尖中心点11的机匣点11坐标为(L_c/2,R_c,0°)。初始最小间隙为叶尖点11与机匣点11之间，碰摩发生在最小间隙附近，因此将机匣点11称为机匣碰摩点。机匣点1和21沿叶片弦向排列，也是机匣的摩擦点。机匣点11'的坐标为(L_c/2,R_c,90°)，机匣点1′、11′、21′为机匣的非碰摩点。

采用半解析法的旋转叶片-机匣摩擦流程图如图10所示。该流程图详细描述了半解析法的数值求解过程。叶片的响应的求解使用了整体叶片质量矩阵M_b、科氏力矩阵G_b、阻尼矩阵D_b和结构刚度矩阵K_b。为了区分叶片与机匣的运动分量，下标3表示叶片的运动。值得注意的是，在法向摩擦力和切向摩擦力的作用下，机匣的模态转换坐标是不同的，因此分别计算了两种碰摩力作用下机匣的运动状态。下标1表示法向碰摩力作用下机匣的运动分量，下标2表示切向碰摩力作用下机匣的运动分量。本实施例将圆柱壳各模态下的质量矩阵M_c、刚度矩阵K_c和阻尼矩阵D_c是不同的分块矩阵，因此实施例采用分块法求解机匣响应。在每个载荷步下，在每个模态p下的机匣的运动状态是分别计算的。当p从1循环到P时，这个荷载步的计算才算完成。θ′是机匣提取点的角度，ωt是叶片的转动角度。

在图10中，和分别是机匣的径向力矢量和切向力矢量。机匣上的随动点的定义为为每一时刻叶尖点正对的那个机匣点。

机匣提取点的径向位移表达式为：

其中是一个TR×1维的轴向模态函数。

机匣上随动点的径向位移转换关系为：

本实施例假设在叶片旋转过程中，叶片表面受均匀气动力作用。在叶片旋转过程中考虑了离心刚化、旋转软化和科氏力的影响。由于安装误差的存在，轮盘中心与机匣中心存在静态不对中，叶尖与机匣内壁间隙不均匀。叶片旋转过程中会产生叶尖径向位移，因此在最小间隙附近叶尖与机匣之间可能发生局部摩擦。因此，有必要在该区域附近进行网格细化，以准确确定叶尖与机匣之间的碰摩相互作用。当叶片旋转1度后，计算叶尖点和对应的机匣点的径向位移，通过判断叶尖与机匣内壁的动态间隙函数gⁱ _gap来确定碰摩的发生。

在实际运行过程中，航空发动机叶片需要有一定的安装角和扭角才能满足设计要求。因此，有必要研究柔性机匣与扭型叶片间的碰摩响应过程。叶片参数如表1所示，机匣参数如表4所示，其他仿真参数如下：

在应用该方法之前，首先要对数值收敛性进行评估(见图10)，模态截断数P对机匣响应的收敛结果如图11所示。随着P值的增加，机匣碰摩点1和机匣非碰摩点1'的响应趋于稳定。当P＝80时，结果基本相同，因此在下面的分析中选择P值为80。将提出的机匣模型模态截断阶数设为TR×P＝7×80

计算了不同时间步长下叶尖点1和机匣碰摩点1处的位移响应。一个周期的总数如下：n_s＝90，180，360，720。时间步Δt＝60/(Ωn_s)分别为6.67×10^-5s，3.33×10^-5s，1.67×10^-5s，8.33×10^-6s。详细的对比结果如图12所示，结果表明该模型具有良好的收敛性。如图12所示，当n_s＝360时，叶片和机匣的响应与n_s＝720时非常接近。考虑到计算效率和精度，在接下来的分析中，时间步长设为1.67×10^-5s(n_s＝360)。

图13、图14和图15分别是提出的模型和有限元模型下的旋转叶片与柔性机匣之间的摩擦响应。这些数据表明，有限元计算结果与目前的计算结果吻合较好。如图13(a1)和13(a2)所示，当含扭角和安装角的叶片与机匣碰摩时，每个叶尖点的碰摩力分布呈现一种沿叶尖弦向两端高中间低的趋势。从图13可以看出，不是所有的叶尖点都同时参与碰摩，这里有一定的时间差约为1×10^-4s。当叶片旋转到初始最小间隙附近时，叶尖点1首先与对应的机匣点发生碰摩，其他叶尖点依次参与碰摩。三个叶尖点的侵入深度曲线如图13(b1)和图13(b2)所示。叶尖点1的侵入深度最大，因此叶尖点1的法向碰摩力最大。叶片径向和切向位移响应如图14所示。由于扭角和安装角的影响，叶尖两端处的径向位移大于中间处的径向位移。在叶片旋转过程中，叶尖1点和叶尖21点在叶片弯曲方向上的离心力分量大小相等，方向相反。两个叶尖点在弯曲方向上具有相反的稳态分量(图14(b1)和图14(b2))。从图15(a1)可以看出，机匣中心点11处的径向位移略大于点1和点21。对比图15(a1)和图15(b1)，机匣碰摩位置的径向位移比非碰摩位置的径向位移大一个数量级。

采用所提出的模型和参考文献模型的叶片和机匣碰摩振动响应如图16所示。参考文献采用ANSYS软件中模拟碰摩过程，叶片和机匣分别简化为板和集中质量点。在相同的工况下，比较板壳碰摩与板集中质量点碰摩引起的振动响应。所提出的模型不仅考虑了机匣柔性的影响，而且考虑了机匣的弹性约束，更加符合实际工况。叶尖和机匣间的最小间隙cⁱ _min＝1μm，和其他参数使用表4中数据。集中质量点机匣的质量为1.67kg，弹簧刚度为2.3×10⁶N/m，阻尼为1000N·s/m。

图16(a1)和图16(a2)分别为两种模型下的叶片碰摩力分布。采用参考文献模型仿真，碰摩力和三个叶尖点的侵入深度均有一定比例的减小。但是由于本文提出的模型的每个机匣点的运动状态不同，所以碰摩力分布不是线性的，而是与每个时刻下叶尖点和对应机匣点的运动状态有关的。对比图16(c1)和图16(c2)可以看出，参考文献的模型得到的弯曲位移的响应大于提出的模型得到的弯曲位移响应。考虑机匣柔性的影响，碰摩力和叶片振动响应会有所减小。叶片的频率分量是一致的。两种模型在3f_r和9f_r处均有幅值放大，这是因为3f_r接近叶片的一阶弯曲频率(486.6Hz)，9f_r接近叶片的一阶扭转频率(1516.9Hz)，其中f_r(166.7Hz)为叶片的转频。为了比较这两种模型，只提取了机匣的碰摩点的响应。参考文献模型在f_r处出现了幅值放大(见图16(f1))，但提出的模型在4f_r处出现了幅值放大(见图16(f2))。这是因为4f_r(666.8Hz)更接近圆柱壳的最低阶固有频率(720.2Hz)。

图17为机匣在不同时刻的变形情况。机匣的变形用剖面图显示，截面图描述了机匣轴向L_c/2处各机匣点的径向位移。这个图是将径向位移放大1×10⁴倍后沿机匣半径方向绘制的。图17(a)中，a点为机匣碰摩点，B点、C点、D点为机匣非碰摩点。在t＝0.0723s时，叶尖与机匣之间存在碰摩，其余三个时刻碰摩消失。从图17(a)可以看出，a点的局部变形明显大于其他点。集中质量点模型下的机匣变形为整体的变形。碰摩时集中质量点机匣向右移动一定的距离，并且变形量小于圆柱壳机匣模型。因此，叶尖与机匣之间的碰摩力和侵入深度小于集中质量点模型的(见图17)。同时可以清楚地观察到碰摩诱发的机匣的节径振动。

图18(a1-d1)和图18(a2-d2)为不同时刻的机匣变形和机匣的节径分量。将圆柱壳径向变形的数据延拓10份，再通过傅里叶变换得到机匣的节径振动。碰摩诱发的机匣的径向振动变形表明，9节径模态(833.8Hz)占主导作用。这是因为5f_r(833.5Hz)的频率与机匣的第9阶固有频率(833.8Hz)最接近，对应于模态(1，9)。因此，最终的机匣变形是以9节径为主导的多节径耦合的组合形式。对比结果表明，考虑机匣柔性时，可以观察到碰摩诱发的机匣的节径振动，并且所提出的模型与集中质量点机匣模型相比，摩擦力和叶片振动响应都会减小。

运用本发明的一种叶片-机匣碰摩关系的确定方法，得到以下结论：

为了研究叶片-机匣碰摩对机匣节径振动的影响规律，建立了旋转叶片-柔性机匣动力学模型。在该模型中，叶片由一个含安装角和扭型板建立，基于Sanders’壳理论建立含弹性约束的柔性机匣模型。通过与有限元模型的固有特性及碰摩振动响应的比较，验证了本文提出的模型的正确性。通过计算叶尖点和对应的机匣点的径向位移，确定每一步载荷步下叶尖点与对应机匣点之间的间隙，进而确定是否发生摩擦。分析了叶片和柔性机匣的振动响应，并与将机匣简化为集中质量点模型的仿真结果进行了比较。讨论了转速、叶片形状、机匣的支承边界条件、叶尖几何形状等因素对机匣碰摩振动响应的影响。总结如下:

当含安装角和扭板叶片与机匣碰摩时，各叶尖点的摩擦力分布沿叶片弦向呈现出叶尖两端高，中部低的趋势。不是所有的叶尖点都同时参与碰摩，而是存在约为1×10^-4s的时间差。叶片振动在3f_r和9f_r处出振幅放大现象这是因为3f_r接近叶片的一阶弯曲频率(486.6Hz)，9f_r接近叶片的一阶扭转频率(1516.9Hz)，其中f_r(166.7Hz)为叶片的转频。

圆柱壳机匣模型可以考虑机匣在碰摩过程中的局部变形，可以很明显地观察到碰摩诱发的机匣的节径振动。在本文研究中，碰摩诱发的机匣变形是以9节点直径为主的多节径耦合的形式。与集中质量点模型相比，摩擦力和叶尖弯曲位移响应均有所减小。集中质量点模型得到的机匣响应在f_r处出现幅值放大现象，而圆柱模型下的机匣在4f_r处会出现幅值放大现象。

需要理解的是，以上对本发明的具体实施例进行的描述只是为了说明本发明的技术路线和特点，其目的在于让本领域内的技术人员能够了解本发明的内容并据以实施，但本发明并不限于上述特定实施方式。凡是在本发明权利要求的范围内做出的各种变化或修饰，都应涵盖在本发明的保护范围内。

Claims (10)

1.一种叶片-机匣碰摩关系的确定方法，采用扭型板的运动状态模拟叶片在机匣中的运动状态，采用圆柱壳的运动状态模拟叶片在机匣中运动时的机匣运动状态，所述扭型板的一端与圆盘的盘周面呈角度连接，所述圆盘设定在所述圆柱壳内部，且所述圆盘的中心在所述圆柱壳中心轴线上，其特征在于，包括以下步骤：

S1、基于扭型板在运动过程中的能量信息获取扭型板的运动状态信息，并根据所述扭型板的运动状态信息得到叶片在机匣中运动时叶片的运动状态信息；

基于弹性约束下圆柱壳的能量信息获取圆柱壳的运动状态信息；并根据所述圆柱壳的运动状态信息得到叶片在机匣运动时机匣的运动状态信息；

其中，所述圆柱壳的两侧分别沿轴向、周向和径向设定平动弹簧以及圆柱壳的旋转方向设定扭转弹簧；

S2、基于叶片在机匣中运动时叶片的运动状态信息和机匣的运动状态信息，获取叶片-机匣系统的整体运动状态信息；

并根据扭型板和圆柱壳的几何关系确定叶片与机匣的碰摩关系；

S3、基于所述叶片与机匣的碰摩关系，确定所述叶片-机匣是否满足生产条件。

2.根据权利要求1所述的方法，其特征在于，在步骤S1前还包括：

依据扭型板的静态参数和动态参数，建立扭型板对应的能量方程，其中所述能量方程包括：扭型板的动能、扭型板应变能、扭型板的离心势能、法向摩擦碰摩力和切向碰摩力以及气动力在弯曲方向所做的功；

相应的，所述步骤S1包括：

基于所述能量方程依据哈密顿原理得到扭型板的运动方程；采用伽辽金法对所述扭型板运动方程进行离散化和模态截断处理以及进行动力学分析获取扭型板的运动状态信息，并根据所述扭型板的运动状态信息得到叶片在机匣中运动时叶片的运动状态信息。

3.根据权利要求2所述的方法，其特征在于，

建立扭型板坐标系系统包括：扭型板全局坐标系OXYZ、扭型板旋转坐标系o'x'y'z'、扭型板局部坐标系oxyz；

基于所述扭型板坐标系系统，建立所述扭型板对应的动能方程：

其中，b为扭型板宽度；h为扭型板厚度；R_d为圆盘半径；θ为扭型板旋转角度；扭角β(x)是从β₀到β₁线性变化的:

其中为扭型板沿径向预扭率，β₀和β₁分别是扭型板根处和扭型板尖处的交错角；β'是扭型板尖相对于扭型板根的扭角；

基于板壳理论，建立扭型板对应的的应变能V₁方程和扭型板离心势能V₂方程式：

其中：

ε_y＝v′_y；

E为扭型板的杨氏模量；

υ扭型板的泊松比；

其中，L为扭型板的长度；f_c1和f_c2为扭型板在两个方向的离心力；

其中，所述f_c1和f_c2的表达式分别为：

其中，扭型板单位面积上的气动力F_e在弯曲方向所做的功W为：

其中F_n为法向摩擦碰摩力，F_t为切向碰摩力；

其中单位面积的气动力F_e的表达式为：

其中，F_e1为气动力幅值，n_e表示障碍物数量；

4.根据权利要求3所述的方法，其特征在于，所述步骤S1包括：

基于扭型板的动能、扭型板应变能、扭型板的离心势能、法向摩擦碰摩力和切向碰摩力以及气动力在弯曲方向所做的功，根据哈密顿原理得到叶片的运动方程表达式为：

其中，V＝V₁+V₂；

采用伽辽金方法对扭型板运动方程进行离散，并基于正则坐标U_(m-1)N+n(t)、V_(m-1)N+n(t)、W_(m-1)N+n(t)，采用组合梁函数法，获取叶片的径向位移、弯曲位移和摆动位移；

其中，

cosh(α_m)cos(α_m)＝-1；

cosh(ψ_n)cos(ψ_n)＝1；

其中，M和N分别表示叶片在坐标x，y方向上的模态截断数；则，根据哈密顿原理对所述叶片运动方程进行变分处理和伽辽金截断处理获取所述叶片的整体运动微分方程为：

其中M_b、G_b、D_b和K_b分别是叶片的质量矩阵、科氏力矩阵、阻尼矩阵和刚度矩阵；

其中，η是正则坐标向量；F_b是外力矢量，包括非摩擦力矢量F_non-rub和摩擦力矢量F_rub，F_b＝F_non-rub+F_rub；

其中，K_b＝K_be+K_bs+K_bc+K_bacc；

K_be为结构刚度矩阵，K_bs为旋软化矩阵，K_bc为离心刚化矩阵和K_bacc为加速度导致的刚度矩阵。

5.根据权利要求1所述的方法，其特征在于，在步骤S1前还包括：

依据弹性约束的圆柱壳和所述圆柱壳上的平动弹簧的静态参数和动态参数以及壳层理论，建立圆柱壳和弹簧的能量方程，其中，所述能量方程包括：圆柱壳的动能、圆柱壳的应变能、弹簧的变形能；并依据所述能量方程获取圆柱壳的能量信息；

相应的，所述步骤S1包括：基于瑞立法原则对所述圆柱壳的能量信息进行处理，获取弹性约束下圆柱壳的运动状态信息，并根据所述圆柱壳的运动状态信息得到叶片在机匣运动时机匣的运动状态信息。

6.根据权利要求5所述的方法，其特征在于，

建立所述圆柱壳正交坐标系Oxzθ；

基于所述圆柱壳坐标系系统，建立所述圆柱壳对应的动能方程和应变能方程包括：

其中，R_c、L_c和h_c分别表示圆柱壳的中曲面的半径、长度和厚度；ρ_c、c_c和E_c分别为圆柱壳的材料密度、泊松比和杨氏模量；圆柱壳的轴向、周向、径向位移分别为u、v和w；在圆柱壳的两侧分别沿轴向、周向和径向引入平动弹簧k^u、k^v和k^w，旋转方向扭转弹簧k^θ；无量纲长度ξ＝x/L_c；

其中

其中

建立用任意周向波数p表示的圆柱壳位移方程：

其中，p为任意周向波数，ω_c为圆柱壳的固有频率；U_p(ξ)、V_p(ξ)和W_p(ξ)是相应的三个变形方向的轴向模态函数；

其中：

其中，a_i、b_i和c_i是待定系数；和是施密特正交化的特征正交多项式；TR是实际计算中截断的项数；

其中，根据公式可以得到多项式的正交集：

其中

对正交多项式进行归一化处理得到

则，特征正交多项式集满足

其中

基于用任意周向波数p表示的圆柱壳位移方程，所述圆柱壳的动能可表示为：

其中是无量纲频率参数；

其中有量纲频率与无量纲频率之间的转换关系为

基于用任意周向波数p表示的圆柱壳位移方程表达式，所述圆柱壳的应变能可表示为：

基于所述圆柱壳坐标系系统，以及圆柱壳和人工弹簧的静态参数和动态参数，获取圆柱壳的自由边界处的人工弹簧中的储存的变形能表达式：

其中，为x＝0处圆柱壳的边缘轴向的平动弹簧刚度；

为x＝0处圆柱壳的边缘周向的平动弹簧刚度；

为x＝0处圆柱壳的边缘径向的平动弹簧刚度；

为x＝0处圆柱壳的边缘平动弹簧扭簧刚度；

为x＝l处圆柱壳的边缘轴向的平动弹簧刚度；

为x＝l处圆柱壳的边缘周向的平动弹簧刚度；

为x＝l处圆柱壳的边缘径向的平动弹簧刚度；

为x＝l处圆柱壳的边缘平动弹簧扭簧刚度；

基于用任意周向波数p表示的圆柱壳位移方程表达式，所述圆柱壳的自由边界处的人工弹簧中的储存的变形能表达式为：

其中无量纲弹簧刚度可表示为：

基于所述圆柱壳的动能、圆柱壳的应变能和圆柱壳上弹簧的变形能获取圆柱壳的能量函数：

Γ＝T_c-U_c-U_s。

7.根据权利要求6所述的方法，其特征在于，所述步骤S1包括：

基于圆柱壳的能量函数，依据瑞立法原则对圆柱壳的能量函数求导使Γ最小；

其中

获取圆柱壳的运动学方程：

其中K_cs是约束弹簧产生刚度矩阵；

K_c和M_c分别为圆柱壳的刚度矩阵和质量矩阵；

a_i、b_i和c_i为待定系数构成的Ritz向量；

其中，a_i＝[a₁,a₂,…a_TR]^T；b_i＝[b₁,b₂,…b_TR]^T；c_i＝[c₁,c₂,…c_TR]^T。

8.根据权利要求1所述的方法，其特征在于，所述步骤S1之前还包括：

根据代表叶片的扭型板和代表机匣的圆柱壳在运动中的几何位置关系，以及代表叶片的扭型板和代表机匣的圆柱壳的动态参数和静态参数，建立扭型板和圆柱壳共同坐标系；

将所述扭型板上的板尖等分点设定为n；

根据扭型板和代表机匣的圆柱壳在运动中的几何位置关系，建立扭型板和圆柱壳之间的间隙函数的表达式：

其中，r_g为静止时刻扭型板轨道半径，r_g＝L+R_d；g₀为同轴圆柱壳与盘的平均间隙，g₀＝R_c-h_c/2-R_d-L；cⁱ _min是第i个板尖点与机匣内壁的初始最小间隙，假设每个板尖点与对应机匣点的初始最小间隙相同；为初始相位角；Ω是扭型板的转速；uⁱ _b为第i个板尖点的径向位移；uⁱ _c为第i圆柱壳点的径向位移；

获取在每个载荷步下，第i个板尖点与其对应的圆柱壳点之间的动态间隙函数的表达式为

基于所述板尖点与其对应的圆柱壳点之间的间隙确定板尖点和相应的圆柱壳点之间的碰摩条件；

其中，在每个载荷步，第i个板尖点的渗透深度δⁱ，

则

当穿透深度δⁱ第i个板尖点大于0，法向和切向碰摩力加载到板尖点i上，同时大小相等方向相反的切向碰摩力量加载到相应的圆柱壳i上，没有发生碰摩时δi＝0；

其中Fⁱ _n和Fⁱ _t分别为第i个板尖点的法向碰摩力和切向碰摩力；

μ是板尖和圆柱壳之间的摩擦系数。

9.根据权利要求4或7或8所述的任一方法，其特征在于，所述步骤S1包括：根据公式(A)获取的扭型板运动状态信息，并根据所述扭型板的运动状态信息表示叶片在机匣中运动时叶片的运动状态信息；

根据公式(B)获取圆柱壳的运动状态信息，并根据所述圆柱壳的运动状态信息表示叶片在机匣中运动时机匣的运动状态信息；

所述步骤S2包括：将表示叶片在机匣中运动时叶片的运动状态信息的公式(A)和表示叶片在机匣中运动时机匣的运动状态信息的公式(B)

组集得到公式(C)，并根据公式(C)获取叶片-机匣系统的整体运动状态信息；

其中，叶片-机匣系统的整体运动方程中的x_b表示扭型板运动公式(A)中的η，代表叶片的正则坐标向量，叶片-机匣系统的整体运动方程中的x_c表示弹性约束下圆柱壳运动公式(B)中的向量；

基于叶片机匣系统的整体运动方程根据所述叶尖点与其对应的机匣点之间的间隙确定叶尖点和相应的机匣点之间的碰摩条件确定叶片与机匣的关系。

10.根据权利要8所述的方法，其特征在于，

所述扭型板上的板尖等分点n设定为21；则

对应与板尖点i的圆柱壳点i的x坐标为：

对应与板尖点i的圆柱点i的θ坐标为：