使用非接触测量和动态响应重构技术的涡轮机叶片疲劳寿命分析

摘要

一种用于对涡轮机叶片的应力和应变场进行动态重构的方法，包括以下步骤：提供来自涡轮机叶片上的至少一个位置的一组响应测量；基于上频率极限和下频率极限对所述一组响应测量进行带通滤波（32）；由所述一组响应测量的局部最小值和局部最大值来确定（33）所述一组响应测量的上包络和下包络；由所述一组响应测量的所述上包络和所述下包络计算（34）候选固有模态函数；当所述候选固有模态函数是实际固有模态函数时，为所述涡轮机叶片提供（37）N×N振型矩阵，其中N是所述涡轮机叶片的自由度数；以及由所述实际固有模态函数和所述振型矩阵中的振型来计算（38）所述涡轮机叶片上的另一位置的响应。

说明书

相关美国申请的交叉引用

本申请要求2012年10月5日提交的Guan等人的美国临时申请号61/710,179“System And Method For Turbine Blade Fatigue Life Analysis Using Non-Contact Measurement And Dynamical Response Reconstruction Techniques”的优先权，其内容通过引用整体并入本文。

技术领域

本申请涉及使用非接触测量和动态响应重构技术的涡轮机叶片疲劳寿命分析。

背景技术

疲劳是当材料受到周期载荷时所发生的渐进且局部结构损坏。传统地，大部分注意力曾经集中在这种情况，即在应力低并且变形主要是弹性的情况下需要超过10⁴次周期以达到失效，这种情况被公知为高周疲劳情况。在这种情况下，材料性能通常由S-N曲线表征，所述S-N曲线描绘了周期应力大小（S）与失效周期的对数尺度（N）的曲线图。S-N曲线由对将要被表征的材料的样品的测试获得，其中在所述样品中，通过测试机施加有规律的正弦曲线应力，所述测试机还可以对达到失效的周期数进行计数。

涡轮机叶片是涡轮机的关键部件。由于一个叶片的失效可能会损坏其它叶片并且导致整个系统的停机这一事实，因此叶片失效可能损失巨大；因此提高叶片的耐久性可以有助于减少涡轮机系统的寿命周期成本。涡轮机叶片的疲劳寿命是用于评估耐久性的重要量。在涡轮机叶片疲劳寿命分析中，使用S-N曲线或基于da/dN的方法对于应力、应力和应变场的准确和可靠的评价是很重要的，其表征作为周期的函数的裂缝长度增长。涡轮机的设计阶段和与操作阶段两者均要求对涡轮机叶片的疲劳寿命进行评价。在设计阶段，使用有限元模型在不同的操作条件执行数值实验以基于模型应力结果定位作为受到高应力的热点。在操作阶段，使用接触或非接触传感器进行测量以间接估计那些热点的应力场。传统分析可能具有监测系统可能仅监测几个表面位置的潜在局限性。传统的接触或非接触测量可能具有整个叶片的空间分散测量，并且仅有那些具有测量数据的位置可以用于进一步分析。

图1中显示了常规的涡轮机叶片非接触测量和疲劳分析系统和方法的工作流程。现在参照图1，在方框10处获得尖端定时测量。用于这些测量的探针为叶片每转一圈提供一个挠度测量。采样速率较低且与发动机旋转同步。通过知道在测量位置处的旋转速度和半径，可以将每一个叶片的到达时间转变成挠度。存在可以被分析的两种普遍的振动类型：（1）由方框11表示的同步振动和（2）由方框12表示的非同步振动。同步振动由相对于发动机壳体固定并通常位于流动场中的驱动器激励。这些振动的频率在相位方面关于发动机速度被锁定。通过诸如LSMF（最小二乘法模型拟合）和SDOF（单自由度）分析的基于最小二乘法建模技术来对同步振动加以分析。非同步振动由通常为非固定的并因此不与发动机速度相关联的驱动器激励。进入系统的能量使转子叶片沿其固有频率模式振动。非同步振动可以使用诸如TRAV（前进波）和IND（单独叶片）技术的基于FFT的分析技术，但是频率通常被混淆。与有限元模型（FEM）13结合的这些分析的输出产生振型14，所述振型14被表示为N×N的矩阵，其中N是自由度数。可以由FEM得到的振型可以用于计算应力-挠度比曲线15，所述应力-挠度比曲线15可用于计算测量位置16处的应力并计算用于高周期疲劳（HCF）寿命预测的S-N曲线17。

常规的系统和方法具有一些缺陷：

（1）由于硬件的限制，所述系统和方法只能测量几个点。整个叶片的应力或应变场难以推断。

（2）应力和应变场分析依赖于应力-挠度比，所述应力-挠度比并不是应力和应变的物理动态响应。

（3）必须认真设计测量位置或传感器布局以获得疲劳寿命的准确估算。由于仅可以测量几个位置并且必须优化所述位置，因此这引入了对传感器布局和有限元分析的组合优化挑战。

（4）整个叶片的应力或应变场的可视化可能是有挑战的，尤其是叶片表面下方的空间位置。例如，在常规的系统中，不能获得基于用户选择的叶片的截面的可视化。

发明内容

本文所述的本发明的示例性实施例总体上包括用于使用分散表面位置测量数据重构整个叶片应力和应变场的系统和方法。在重构之后，可以计算叶片几何形状的任意空间位置的疲劳寿命以用于设计优化和基于条件的维护。

根据本发明的一方面，提供了一种用于对涡轮机叶片的应力和应变场进行动态重构的方法，包括以下步骤：提供来自涡轮机叶片上的至少一个位置的一组响应测量；基于上频率极限和下频率极限对所述一组响应测量进行带通滤波；由所述一组响应测量的局部最小值和局部最大值来确定所述一组响应测量的上包络和下包络；由所述一组响应测量的所述上包络和所述下包络计算候选固有模态函数；当所述候选固有模态函数是实际固有模态函数时，为所述涡轮机叶片提供N×N振型矩阵，其中N是所述涡轮机叶片的自由度数；以及由所述实际固有模态函数和所述振型矩阵中的振型来计算所述涡轮机叶片上的另一位置的响应。

根据本发明的另一个方面，所述一组响应测量的所述上频率极限和所述下频率极限通过对所述一组响应测量进行傅里叶变换并由峰值拾取对主频分量进行识别来加以确定。

根据本发明的另一个方面，在所述涡轮机叶片上的所述至少一个位置是尖端位置。

根据本发明的另一个方面，所述一组响应测量的所述上包络和所述下包络由连接所述一组响应测量的所述局部最小值和所述局部最大值的三次样条函数来确定。

根据本发明的另一个方面，所述由所述一组响应测量的所述上包络和所述下包络计算候选固有模态函数的步骤包括：计算所述一组响应测量的所述上包络和所述下包络的平均包络；和通过从所述一组响应测量减去所述平均包络计算所述候选固有模态函数。

根据本发明的另一个方面，其中，如果所述候选固有模态函数不是实际固有模态函数，则对所述候选固有模态函数执行确定上包络和下包络的以及所述计算所述上包络和所述下包络的平均包络的步骤，并且通过从前一个候选固有模态函数减去所述前一个候选固有模态函数的所述平均包络来计算新的候选固有模态函数。

根据本发明的另一个方面，所述方法包括重复所述确定所述一组响应测量的上包络和下包络的步骤和所述由所述一组响应测量的所述上包络和所述下包络计算候选固有模态函数的步骤，直到不能从所述一组响应测量计算出进一步的固有模态函数为止。

根据本发明的另一个方面，通过求解（M^-1K）的特征系统计算N×N振型矩阵，其中，M是表示所述涡轮机叶片的质量的矩阵，而K是表示所述涡轮机叶片的刚度的矩阵。

根据本发明的另一个方面，计算所述涡轮机叶片上的另一位置的响应包括计算                                                ，其中下标e表示能够被测量的自由度，u表示测量所得不到的自由度，并且i表示第i个模态，φ_ik表示在第i个模态下第k个自由度的振型，而δ_ie表示在对于第e个自由度的第i个固有模态函数。

根据本发明的另一个方面，所述方法包括由计算出的响应计算基于S-N曲线的疲劳寿命估计。

根据本发明的另一个方面，提供了一种用于对涡轮机叶片的应力和应变场进行动态重构的方法，包括以下步骤：提供来自涡轮机叶片上的至少一个位置的一组响应测量；使用经验模态分解由所述一组响应测量来确定固有模态函数；以及通过计算来计算所述涡轮机叶片上的另一位置的响应，其中下标e表示能够被测量的自由度，u表示测量所得不到的自由度，并且i表示第i个模态，φ_ik表示在第i个模态下第k个自由度的振型，而δ_ie表示对于第e个自由度的第i个固有模态函数，其中所述振型是通过求解（M^-1K）的特征向量而计算出的N×N振型矩阵的分量，其中M是表示所述涡轮机叶片的质量的矩阵，而K是表示所述涡轮机叶片的刚度的矩阵。

根据本发明的另一个方面，使用经验模态分解来确定固有模态函数的步骤包括：由所述一组响应测量的局部最小值和局部最大值确定所述一组响应测量的上包络和下包络；计算所述一组响应测量的所述上包络和所述下包络的平均包络；以及通过从所述一组响应测量减去所述平均包络来计算所述候选固有模态函数。

根据本发明的另一个方面，所述方法包括通过对所述一组响应测量进行傅里叶变换并由峰值拾取对主频分量进行识别来确定所述一组响应测量的上频率极限和下频率极限；以及在使用经验模态分解来确定固有模态函数之前，基于所述上频率极限和所述下频率极限对所述一组响应测量进行带通滤波。

根据本发明的另一个方面，提供一种由计算机可读的非暂时程序存储设备，其能够有形地具体化由所述计算机运行的指令的程序以执行用于对涡轮机叶片的应力和应变场进行动态重构的方法步骤。

附图说明

图1示出了根据本公开的实施例的常规的系统和方法的工作流程。

图2示出了根据本公开的实施例的系统和方法的工作流程。

图3示出了根据本公开的实施例的动态响应重构方法的流程图；

图4示出了根据本公开的实施例的叶片FEM模型和测量位置。

图5示出了根据本公开的实施例的作为时间的函数的所施加的力。

图6显示了根据本公开的实施例的重构结果和理论解之间的对比的曲线图。

图7示出了根据本公开的实施例的当由重构结果和理论结果之间的校正测量时的统计性能。

图8是根据本公开的实施例的用于动态响应重构的系统的框图。

具体实施方式

本文所述的本发明的示例性实施例总体上提供了用于动态响应重构的系统和方法。虽然实施例易于具有各种修改和替换形式，但是在附图中以示例的方式显示了所述实施例的具体实施例，并且本文将详细描述所述具体实施例。然而，要理解的是不旨在将本发明受限于所公开的具体形式，而是相反，本发明涵盖落入本发明的精神和范围内的所有修改、等同物和替换形式。

根据本公开的实施例的系统和方法可以通过使用动态响应重构方法解决常规的系统和方法的缺陷。基本思想是基于有限元模型、测量响应和经验模态分解技术重构整个叶片的应力和应变场。图2中示出了根据本公开的实施例的系统和方法的工作流程。图2中的方框20到24对应于图1中的方框10到14，因此将不会进一步描述该方框20到24。在方框25中处，振型和FE模型与根据本公开的实施例的重构方法一起使用以基于尖端测量对整个涡轮机叶片计算应力。以下将参照以下所述的图3提供根据本公开的实施例的重构方法的细节。一旦已经对整个叶片计算应力，则可以在方框26处看到整个叶片的应力剖面图，并且可以在方框27处对高周疲劳（HCF）寿命预测计算S-N曲线。

根据本公开的实施例的动态响应重构方法基于具有间歇判据的经验模态分解（EMD）技术和结构模型，并且包括两个步骤。第一步骤是使用具有间歇判据的EMD将所测量的信号分解成被公知为固有模态函数（IMF）的模态响应。在第二步骤中使用产生的IMF来重构在传感器不可测的关键位置处的动态响应。

EMD的思想是将给定信号分解成对原始信号形成完全或近正交基的函数。具有零平均值并在被公知为IMF的零交叉之间仅具有一个极值的所产生的函数足以描述信号的基本动态。

固有模态函数（IMF）是满足以下所述的两个条件的函数：（1）在整个数据集中，极值的数量和零交叉的数量相同或最多相差一个；以及（2）在任意点处，由局部最大值限定的包络和由局部最小值限定的包络的平均值是零。通过这种限定，在每一个周期中由零交叉限定的IMF仅涉及一个振荡模态，且不会有复杂的骑行波（riding wave）。

通过信号分解的方式保证了EMD方法的充分性和完整性。使用IMF表示信号允许对时变频率进行保存，这在傅立叶或小波域中被隐藏。

为了从给定信号滤波出IMF，EMD方法使用了筛选处理。为了获得与系统的模态响应相对应的IMF，可以将由ω_int表示的间歇频率施加到筛选处理以确保IMF中的每一个都仅包含一个频率分量。被公知为间歇判据的这种思想将从IMF剔除低于或大于ω_int的所有频率分量，并且这可以使用带通滤波器在筛选处理期间完成。

根据本公开的实施例，图3是根据本公开的实施例的用于使用A处的测量来对B处的动态响应进行重构的动态响应重构方法的流程图。对于被表示为y(t)的时域信号，根据本公开的实施例的筛选处理如下所述。现在参考附图，在步骤31处，使用快速傅里叶变换（FFT）或FEM计算从测量数据中识别出主频分量。这些用于表示ω_int的上限ω_iH和下限ω_iL，ω_iL<ω_int<ω_iH。在FFT之后，对傅氏谱使用峰值拾取算法识别主频分量。接下来，在步骤32处，将ω_iL<ω_int<ω_iH频率范围内的带通滤波施加到原始时域数据y(t)，以获得与那些被识别的模态相关联的响应。接着，在步骤33处，由信号的局部最大值和最小值确定上包络和下包络。用于计算上包络和下包络的示例性非限制技术使用三次样条插值以连接局部最大值并连接局部最小值。设m₁表示上包络和下包络的平均包络。在步骤34处，第一个候选分量被计算为h₁＝y(t)－m₁，并且在步骤35处确定h₁是否为IMF。如果不是，重复步骤33和34，这次使用h₁作为新信号，代替y(t)。这里，新的平均包络被表示为m₁₁表示，并且新的候选分量是h₁₁＝h₁－m₁₁，其中m_1k是k次迭代的上包络和下包络的平均值。从步骤35重复步骤33和34直到产生的结果h_1k（即h_1k＝h_1(k-1)－m_1k）为IMF为止。来自数据的第一IMF分量可以被表示为f₁(t)=h_1k。其包含信号的最短周期分量，并且是模态坐标下的响应。

如上所述，根据本公开的实施例的动态响应重构方法仅基于特征时间尺度从数据分离出最好的局部模态。然而，根据本公开的实施例的筛选处理有两个作用：（a）消除骑行波；和（b）使不均匀振幅平滑。然而，第二作用可能会完全除去有意义的振幅波动。为了确保IMF分量保留振幅和频率调制两者的足够的物理意义，应该存在停止判据。这可以通过限制由两次连续筛选结果计算的标准差SD的大小来完成。基于通过从来自同一数据的1024个点移出五个点而计算的用于两个傅氏谱的SD，示例性非限制SD值可以被设定到0.2与0.3之间。

以上步骤33到35可以从步骤36开始重复以从剩余部分信号r₁(t)＝y(t)－f₁(t)获得第二IMF f₂(t)。可以重复这些步骤直到获得第n个IMF。根据本公开的实施例的筛选处理可以通过以下任意判据加以停止：要么当IMF分量f_n或剩余部分r_n变小使得其小于预定值时；要么当剩余部分r_n变成从其中不能提取更多IMF的单调函数时。即使对于具有零平均值的数据，最终的剩余部分仍然能够不等同于零；对于具有趋势的数据，则最终的剩余部分应该为所述趋势。

通过将根据本公开的实施例的筛选处理和EMD方法应用到传感器测量数据，可以在不需要留下时域的情况下提取所有所需的模态响应。如方程（1）所示，原始信号y(t)可以被表示为n个IMF和剩余项的总和：

            （1）

其中f_i(t)是第i个IMF，而r(t)是剩余部分。通过对不同的固有频率以不同的频率范围重复，可以获得所有模态响应。使用根据本公开的实施例的具有间歇判据的筛选处理，原始信号表达式可以被写成方程（2）

      （2）

其中x_i(t)是第i个模态响应（也是IMF）。项f_i(t), i=1,…,n×m是其它IMF但不是模态响应。

这些IMF具有多个特性：

（1）每一个IMF都包含信号的固有特性；

（2）一旦获得IMF，下一个IMF将在同一时刻具有不同的频率；以及

（3）每一个IMF系列的第一IMF被认为是模态响应的近似值。

以上获得的模态响然后用于使用变换方程对识别的模态中的每一个的期望位置处的模态响应进行外推。使用具有间歇判据的EMD方法已经将时间域中的传感器测量分解成多个IMF。可以利用变换方程对除了传感器位置之外的位置处的动态响应进行外推。这些变换方程可以从目标结构的有限元模型获得。

再次参照图3，根据本公开的实施例的外推过程如下。

考虑n-DOF动态系统，其质量、刚度和阻尼矩阵被分别表示为M、K和C。矩阵形式的运动方程可以被表示为：

   （3）

其中X是响应可变矢量，而F是载荷矢量。根据结构动态，在步骤37处，可以通过求解（M^-1K）的特征系统从涡轮机叶片的FEM获得该结构动态的N×N振型矩阵，其中N是自由度数。产生的振型矩阵φ被表示为方程（4）

                （4）

方程（4）的矩阵的每一列对应于一个特征向量且代表一种振型。例如，矩阵的第一列是第一振型，矩阵的第二列是第二振型等。在自由振动下，相对应的特征值与使用质量矩阵M和刚度矩阵K表示的结构的理论固有频率相关。相关性为f_i＝sqrt(D_i)/(2π)，其中D_i是第i个特征值，而f_i是第i个结构固有频率。

基于模态分析，变换方程可以被写成：

         （5）

其中下标e表示可以由传感器测量的DOF，u表示使用传感器的测量来说不可获得的DOF，而i表示第i个模态。φ_ie表示在第i个模态下第e个DOF的模态信息，而δ_ie表示在第i个模态下第e个DOF的模态响应。想到使用图3中的步骤来对测量响应进行滤波，并且滤波的输出为IMF。因此，由图3中得到的IMF是方程（5）中的δ_ie，其中方程（5）中的下标“i”表示第i个模态。

方程（5）构建了两个坐标之间在结构上的关系。在步骤38处，使用该方程，传感器不可获得的位置处的模态响应δ_iu可以由测量位置处的模态响应δ_ie计算。变换方程由结构的有限元分析获得，并且其简单性确保了计算是非常高效的。

给定振型矩阵φ和IMF，可以由方程（5）容易地求得期望位置处的模态响应。在已经计算出所有模态响应之后，可以应用模态叠加方法以获得时域中的动态响应。最后，进行模态叠加以将模态响应转换成物理响应。

作为示例，假设测量位置是第6个DOF（自由度），并且期望重构与第2个DOF对应的位置处的物理响应。因此下标e=6且下标u=2。假定通过求解（M^-1K）的特征向量得到矩阵φ。对于第一模态（i=l）：使用如图3所示的根据本公开的实施例的方法处理测量数据。果是由δ₁₆表示的IMF。从φ矩阵中获得φ₁₆和φ₁₂的值。由方程（5）计算δ₁₂。对于第二模态（i=2）：使用φ₂₆、φ₂₂和δ₂₆进行与第一模态相同的操作，以获得δ₂₂。对于第三模态（i=3）：使用φ₃₆、φ₃₂和δ₃₆进行与第一模态相同的操作以获得δ₃₂。重复以上步骤直到获得所有重要模态。想到重要模态是可以从测量数据的傅氏谱被识别的重要模态。以下第2个DOF重构响应简单地为δ₁₂+δ₂₂+δ₃₂+…的总和。

一旦应力和应变响应被重构，可以使用雨流计数法获得基于S-N曲线的疲劳寿命估计的雨流矩阵。如果存在类似于裂缝的瑕疵，则还可以应用基于da/dN的方法。

与常规的系统和方法比较，根据本公开的实施例的系统和方法具有以下特征：

（1）根据本公开的实施例的系统和方法需要一个测量位置的最小值来重构整个叶片的应力或应变场。

（2）通过根据本公开的实施例的方法获得的整个叶片的应力和应变场是物理响应，所述物理响应可以与其它分析（例如，使用物理响应的傅氏谱的频域分析）直接有关。

（3）测量的位置理论上不会影响根据本公开的实施例的重构方法的性能，而且不需要特殊优化来为测量确定最好的位置。

（4）因为根据本公开的实施例的系统和方法可以重构整个叶片的应力和应变场，因此可以容易地获得应力和应变分布概况、热点和最大应力的位置、以及其它有用的量。如图8所示，可以便利地实施整个叶片的3D可视化和给定切割平面的截面图。

（5）根据本公开的实施例的系统和方法不需要额外的硬件和建模。诸如有限元模型和测量数据的当前输入可以直接用在根据本公开的实施例的系统和方法中。此外，数值模拟可以与根据本公开的实施例的系统和方法一起使用以优化传感器布局。例如，基于用户给定或随机载荷剖面，可以通过数值模拟找到应力和应变比其它位置大的多的热点、或具有大挠度的点，并且那些位置可以是测量位置的合理候选。

本文以下给出了说明本公开的实施例的整体思想的说明性示例。示例性模型是受到冲击载荷的叶片。通过求解有限元模型的运动方程（EoM）来获得时域响应。选择叶片尖端周围任意位置的响应以代表尖端定时测量数据。使用根据本公开的实施例的方法，采用代表性测量数据重构叶片的其它位置的动态响应。结果与使用FEM的那些位置的直接解进行比较。图4（a）表示叶片FEM模型，而图4（b）示出了尖端定时测量位置。对x-y平面上的节点（即具有z=0的节点），限制了所有自由度（DOF）。整个叶片FEM模型有5820个活跃的DOF。叶片的材料是铝合金。图5是作为时间的函数的施加力的曲线图。图6示出了使用具有10%的噪点分量的测量数据在重构挠度响应与所选择的自由度（DOF）的理论解之间的比较。重构位置被任意选择。

相对于理论结果，可以使用重构结果的相关系数对根据本公开的实施例的方法的性能进行定量。选择位置的完美重构结果使得相关系数为1。由于测量噪点和建模不确定性，所述系数应该接近于1。图7示出了使用了具有10%噪点分量的测量数据的整体性能。图7（a）示出了所有DOF的相关系数，而图7（b）示出了系数的标准化直方图。如果相关系数大于0.9被认为是准确的重构，在本示例中所有DOF（总共5820个DOF）中的大约97%可以被准确地重构。

要理解的是本公开的实施例可以以各种形式的硬件、软件、固件、专用处理或其组合被实施。在一个实施例中，本公开可以在诸如被有形地具体化在计算机可读程序存储设备上的应用程序的软件上实施。应用程序可以被上传到包括任何适当体系结构的机器或由所述机器运行。

图8是根据本公开的实施例的用于动态响应重构的系统的框图。现在参照图8，用于实施本发明的计算机系统81尤其可以包括中央处理单元（CPU）82、存储器83和输入/输出（I/O）接口84。计算机系统81总体上通过I/O接口84联接到显示器85和诸如鼠标和键盘的各种输入设备86。支持电路可以包括诸如电高速缓存、电源、时钟电路和通信总线的电路。存储器83可以包括随机存取存储器（RAM）、只读存储器（ROM）、硬盘驱动器、磁盘驱动器等或其组合。本发明可以作为存储在存储器83中并由CPU 82运行的程序87来实施以处理来自信号源88的信号。因此，计算机系统81是通用计算机系统，当运行本发明的程序87时，所述通用计算机系统变成专用计算机系统。

计算机系统81还包括操作系统和微型指令代码。本文描述的各种处理和功能可以是微指令代码的一部分或经由操作系统运行的应用程序的一部分（或其组合）。此外，各种其它外围设备可以连接到诸如另外的数据存储设备和打印设备的计算机平台。

进一步要理解的是由于附图中所示的组成系统部件和方法步骤中的一些可以在软件中实施，因此系统部件（或处理步骤）之间的实际连接可能会基于对本发明进行编程的方式而不同。给出本文提供的本发明的教导，本领域的普通技术人员将能够构思本发明的这些或类似实施方式或构造。

虽然已经参照示例性实施例详细地说明了本发明，但是本领域技术人员将理解在不偏离本发明的如所述权利要求中所述的精神和范围的情况下可以对本发明进行各种修改和替换。