一种基于分段权函数的涡轮盘裂纹扩展模拟方法与流程

本发明是一种针对航空发动机涡轮盘的裂纹扩展模拟方法，它是一种考虑涡轮盘裂纹扩展平面应力分布情况基于分段权函数求解裂纹强度应力从而模拟涡轮盘裂纹扩展的方法，属于航空航天发动机技术领域。

背景技术：

根据断裂力学可知，应力强度因子与结构形状、尺寸、载荷及裂纹长度密切相关，为保证裂纹扩展分析的计算效率，工程上通常采用简化模型代替实际结构，以远端均布应力代替实际载荷，从而建立应力强度因子的经验公式，如现有文献"newmanjc,rajuis.stressintensityfactorequationsforcracksinthree-dimensionalfinitebodiessubjectedtotensionandbendingloads[r].hampton,virginia,usa:nasalangleyresearchcenter,1984."针对三维矩形板裂纹模型，建立了多种类型裂纹在承受拉伸和弯曲载荷下的应力强度因子近似表达式。

对于复杂结构而言，简化模型端面应力往往呈现非均匀应力分布，导致经验公式中的均布载荷假设不成立，从而使裂纹扩展预测结果存在一定偏差。因此，如何在非均布应力场下准确求解应力强度因子，是决定裂纹扩展模拟精度、效率的关键一环。现有文献"buecknerhf.anovelprincipleforthecomputationofstressintensityfactors[j].journalofappliedmathematicsandmechanics,1970,50(10):529-546."提出权函数方法，用以计算复杂应力分布下的应力强度因子。现有文献"kiciaka,glinkag,emanm,etal.weightfunctionsandstressintensityfactorsforcornerquarter-ellipticalcrackinfinitethicknessplatesubjectedtoin-planeloading[j].engineeringfracturemechanics,1998,60(2):221-238."建立了复杂非均匀应力分布下三维矩形板裂纹体的应力强度因子通用权函数计算方法。然而该方法只适用于应力沿板厚度方向分布的情况，而对于航空发动机涡轮盘等部位，其裂纹面应力分布主要是沿着径向分布的周向应力，而将此裂纹进行简化成三维矩形板裂纹时，裂纹面上应力分布则沿矩形板宽度方向分布，因此该方法不太适用。

技术实现要素：

本发明技术解决方案：克服现有涡轮盘裂纹扩展模拟方法的不足，提供一种任意非线性应力分布裂纹应力强度因子的权函数表达式方法，能够充分考虑涡轮盘裂纹扩展平面应力分布情况采用分段权函数模拟涡轮盘偏心孔处裂纹扩展，提高涡轮盘裂纹扩展模拟预测精度。

本发明技术解决方案：一种基于分段权函数的涡轮盘裂纹扩展模拟方法，结合涡轮盘静强度分析结果，确定危险部位及相应的应力分布，根据应力强度因子模拟结果，确定任意非线性应力分布权函数，最后结合涡轮盘危险部位平面应力分布采用分段拟合权函数开展涡轮盘裂纹扩展模拟，实现步骤如下：

(1)基于涡轮盘静强度分析，确定危险部位。为提高有限元计算效率、减少计算量，取涡轮盘一个扇区区域作为分析模型，施加周期对称边界条件；为便于网格划分将涡轮盘切分为榫槽、偏心孔、盘心和盘身多个子模块，采用六面体单元划分网格，提高计算精度；查阅高温合金手册获得材料参数，包括弹性模量、泊松比和密度，覆盖涡轮盘全部服役温度；载荷为涡轮盘设计工作状态的转速和叶片离心力(转换为榫槽接触面拉应力)，温度场由涡轮盘典型位置温度通过拟合获得；位移边界条件为约束安装边的轴向和周向位移；进行静强度分析，确定涡轮盘危险部位。

(2)基于应力强度因子模拟结果，确定任意非线性应力分布权函数。

对于孔边椭圆角裂纹，如图2(a)所示，记孔外表面裂纹前缘为a点(对应裂纹长度a)、孔内表面裂纹前缘为b点(对应裂纹长度b)。

根据权函数法定义，a、b两点的应力强度因子可表示为：

其中，σ(x)为a、b两点的应力，ma(x,a)、mb(x,a)分别为a、b两点的权函数，可表示为：

其中，a为裂纹长度，m1a、m2a、m3a、m1b、m2b、m3b为代表结构特征的待定系数，通过参考载荷状态下的应力强度因子有限元解加以确定。

分别选取沿简化模型x向呈均布、线性、二次分布的裂纹面正应力作为参考载荷状态，如图4所示，用以确定上述待定系数。其中，σ0为常数，a为裂纹长度，为几何修正因子，分别对应i＝0,1,2分别表示参考载荷为均布、线性分布和二次分布，q为第2类完全椭圆积分近似表达式的平方，q＝e(k)²，时，时，

将图4中的应力分布及式(3)、式(4)带入式(1)、式(2)，积分得：

其中，kix(i＝0,1,2,x＝a,b)为a、b两点分别在参考载荷为均布、线性分布和二次分布时的应力强度因子，fix(i＝0,1,2,x＝a,b)为几何修正因子。

将m1a、m2a、m3a、m1b、m2b、m3b表示为f0a、f1a、f2a、f0b、f1b、f2b的形式：

影响应力强度因子的几何参数主要有：裂纹长度a和b、矩形板厚度t，采用有限元奇异单元法计算孔边椭圆角裂纹分别在裂纹面承受均布、线性分布和二次分布应力时的应力强度因子以建立应力强度因子表达式中的几何修正因子。基于参考应力强度因子解，采用多项式函数分布拟合获得了形式较为统一的、高精度的参考应力强度因子表达式。

取fia、fib(i＝0,1,2)为(b/t)的三次多项式，即：

其中，αij、βij(i＝0,1,2、j＝0,l,3)为几何修正因子fix(i＝0,1,2,x＝a,b)的拟合系数，可表示为(b/a)的二次多项式，即：

其中，为αij、βij的拟合系数。

franc3d是目前唯一可以计算一般各向异性材料中三种断裂模式应力强度因子k的断裂力学软件，且其理论模型适用性强、软件参数设置灵活、计算精度较高，因此本发明选取franc3d结合有限元软件ansys进行断裂力学计算。

通过有限元建立有限厚度(t＝7.8mm)平板角裂纹模型，在σ0＝100mpa下分别计算图4所示均布、线性、二次分布载荷下a、b两点的应力强度因子、几何因子及参数αij、βij。采用最小二乘法拟合确定pk^ij取值。

同理，对于孔边穿透裂纹，在权函数法中，应力强度因子、几何因子、拟合参数之间的关系可表示为：

其中，k为应力强度因子，σ(x)为应力，m(x,a)为权函数，可表示为：

其中，m1、m2、m3为代表结构特征的待定系数，通过参考载荷状态下的应力强度因子有限元解加以确定。

其中，fi(i＝0,1,2)为几何修正因子。

将m1、m2、m3表示为f0、f1、f2的形式：

取fi(i＝0,1,2)为(b/t)的二次多项式，即：

其中，αik(i＝0,1,2)为fi的拟合系数。

通过有限元建立有限厚度(t＝7.8mm)平板穿透裂纹模型，在σ0＝100mpa下分别计算图4所示均布、线性、二次分布载荷下裂纹尖端应力强度因子、几何因子及拟合参数αik取值。

以上，计算得到孔边椭圆角裂纹分别在裂纹面承受均布、线性分布和二次分布应力时的应力强度因子解及经验关系式，将其分别带入式(5)～式(8)，可获得a，b两点权函数表达式的待定系数，将待定系数带入式(3)和式(4)，最终获得a，b两点的权函数；同理得到孔边穿透裂纹分别在裂纹面承受均布、线性分布和二次分布应力时的应力强度因子解及经验关系式，将其分别带入式(14)～式(16)，获得权函数表达式的待定系数，将待定系数带入式(13)，最终获得孔边穿透裂纹的权函数。

(3)基于步骤(2)建立的权函数，开展涡轮盘裂纹扩展模拟。分析步骤(1)确定的危险部位应力特点，确定裂纹扩展的方向(径向或者周向，若最大周向应力大于最大径向应力，径向裂纹扩展相比于周向裂纹扩展出现的可能性更高)并绘制相应的裂纹扩展平面分布。依据裂纹扩展应力强度因子的权函数方法求解过程，对于不同的裂纹长度a，根据裂纹扩展平面应力分布规律将平面应力分布表示为分段函数形式，且分段函数为(1-x/a)的多项式，即：

其中，s0、s1、s2为给定裂纹长度a下，依据应力分布数据确定的拟合参数。当拟合参数确定后，分别求解均布、线性、二次应力分布下的裂纹扩展应力强度因子，并根据叠加原理计算三者共同作用时的应力强度因子。

基于上述发明方法，计算涡轮盘危险部位的裂纹扩展应力强度因子。初始1/4圆形角裂纹半径选定1～3mm，每次计算间隔给定裂纹扩展增量0.05～0.2mm，得到涡轮盘危险部位表面裂纹长度(a)与内表面裂纹长度(b)随循环数n的变化关系，从而模拟得到涡轮盘裂纹扩展规律。

本发明与现有技术相比的优点在于：本发明基于分段权函数的涡轮盘裂纹扩展模拟方法，与传统的涡轮盘裂纹扩展模拟方法相比，利用franc3d+ansys的有限元方法计算应力强度因子解并此基础上获得了任意非线性应力分布角裂纹应力强度因子的权函数，采用均布、线性、二次分布的裂纹面正应力作为参考载荷状态建模，可以进行任意载荷状态下的裂纹扩展模拟分析，应用范围广；进行权函数参数拟合时采用分段拟合，精度更高。

附图说明

图1为本发明一种基于分段权函数的涡轮盘裂纹扩展模拟方法流程图；

图2为孔边裂纹简化模型；(a)孔边椭圆角裂纹，(b)孔边穿透裂纹；

图3为基于分段权函数法和有限元法的涡轮盘裂纹扩展模拟结果。

图4为参考载荷及参考应力强度因子图。

具体实施方式

下面结合附图，对本发明一种基于分段权函数的涡轮盘裂纹扩展模拟方法的技术方案做进一步说明。

考虑涡轮盘裂纹扩展平面应力分布情况，本发明提出的基于分段权函数的涡轮盘裂纹扩展模拟方法，流程图如图1所示。

(1)基于涡轮盘静强度分析，确定危险部位。为提高有限元计算效率、减少计算量，取涡轮盘1/10扇区区域(36°/360°)作为分析模型，施加周期对称边界条件；为便于网格划分将涡轮盘切分为榫槽、偏心孔和盘心等子模块，采用六面体单元划分网格，提高计算精度；查阅高温合金手册获得gh4149镍基高温合金材料参数，包括弹性模量、泊松比和密度，覆盖涡轮盘全部服役温度；载荷为涡轮盘设计工作状态的转速(20900r/min)和叶片离心力(转换为榫槽接触面拉应力495mpa)，温度场由涡轮盘典型位置温度通过拟合获得；位移边界条件为约束安装边的轴向和周向位移；进行静强度分析，确定涡轮盘危险部位。

(2)基于应力强度因子模拟结果，确定任意非线性应力分布权函数。

对于孔边椭圆角裂纹，如图2中(a)、(b)所示，记孔外表面裂纹前缘为a点(对应裂纹长度a)、孔内表面裂纹前缘为b点(对应裂纹长度b)。

根据权函数法定义，a、b两点的应力强度因子可表示为：

其中，σ(x)为a、b两点的应力，ma(x,a)、mb(x,a)分别为a、b两点的权函数，可表示为：

其中，m1a、m2a、m3a、m1b、m2b、m3b为代表结构特征的待定系数，通过参考载荷状态下的应力强度因子有限元解加以确定。

分别选取沿简化模型x向呈均布、线性、二次分布的裂纹面正应力作为参考载荷状态，如图4所示，用以确定上述待定系数。

其中，σ0为常数，a为裂纹长度，为几何修正因子，分别对应i＝0,1,2分别表示参考载荷为均布、线性分布和二次分布，q为第2类完全椭圆积分近似表达式的平方，q＝e(k)²，时，时，

将图4中的应力分布及式(3)、式(4)带入式(1)、式(2)，积分得：

其中，kix(i＝0,1,2,x＝a,b)为a、b两点分别在参考载荷为均布、线性分布和二次分布时的应力强度因子，fix(i＝0,1,2,x＝a,b)为几何修正因子。

将m1a、m2a、m3a、m1b、m2b、m3b表示为f0a、f1a、f2a、f0b、f1b、f2b的形式：

影响应力强度因子的几何参数主要有：裂纹长度a和b、矩形板厚度t，采用有限元奇异单元法计算孔边椭圆角裂纹分别在裂纹面承受均布、线性分布和二次分布应力时的应力强度因子以建立应力强度因子表达式中的几何修正因子。矩形板厚度和裂纹尺寸配置方案为：b/t＝0.2，0.4，0.6，0.8；b/a＝0.2，0.4，0.6，0.8，1，2，3，4。基于参考应力强度因子解，采用多项式函数分布拟合获得了形式较为统一的、高精度的参考应力强度因子表达式。

取fia、fib(i＝0,1,2)为(b/t)的三次多项式，即：

其中，αij、βij(i＝0,1,2、j＝0,l,3)为几何修正因子fix(i＝0,1,2,x＝a,b)的拟合系数，可表示为(b/a)的二次多项式，即：

其中，pk^ij(k＝1,2,l,6)为αij、βij的拟合系数。

franc3d是目前唯一可以计算一般各向异性材料中三种断裂模式应力强度因子k的断裂力学软件，且其理论模型适用性强、软件参数设置灵活、计算精度较高，因此本发明选取franc3d结合有限元软件ansys进行断裂力学计算。

通过有限元建立有限厚度(t＝7.8mm)平板角裂纹模型，在σ0＝100mpa下分别计算图4所示均布、线性、二次分布载荷下a、b两点的应力强度因子、几何因子及参数αij、βij。采用最小二乘法拟合确定pk^ij取值。

同理，对于孔边穿透裂纹，在权函数法中，应力强度因子、几何因子、拟合参数之间的关系可表示为

其中，k为应力强度因子，σ(x)为应力，m(x,a)为权函数，可表示为：

其中，m1、m2、m3为代表结构特征的待定系数，通过参考载荷状态下的应力强度因子有限元解加以确定。

其中，fi(i＝0,1,2)为几何修正因子。

将m1、m2、m3表示为f0、f1、f2的形式：

取fi(i＝0,1,2)为(b/t)的二次多项式，即：

其中，αik(i＝0,1,2)为fi的拟合系数。

通过有限元建立有限厚度(t＝7.8mm)平板穿透裂纹模型，在σ0＝100mpa下分别计算图4所示均布、线性、二次分布载荷下裂纹尖端应力强度因子、几何因子及拟合参数αik取值。

以上计算得到孔边椭圆角裂纹分别在裂纹面承受均布、线性分布和二次分布应力时的应力强度因子解及经验关系式，将其分别带入式(5)～式(8)，可获得a，b两点权函数表达式的待定系数，将待定系数带入式(3)和式(4)，最终获得a，b两点的权函数；同理得到孔边穿透裂纹分别在裂纹面承受均布、线性分布和二次分布应力时的应力强度因子解及经验关系式，将其分别带入式(14)～式(16)，获得权函数表达式的待定系数，将待定系数带入式(13)，最终获得孔边穿透裂纹的权函数。

(3)基于步骤(2)建立的权函数，开展涡轮盘裂纹扩展模拟。考查涡轮盘危险部位偏心孔处，该处最大周向应力大于最大径向应力，径向裂纹扩展相比于周向裂纹扩展出现的可能性更高，以涡轮盘偏心孔处的径向裂纹扩展为例。

基于涡轮盘静强度分析结果，绘制涡轮盘辐板处径向裂纹扩展平面应力分布。依据裂纹扩展应力强度因子的权函数方法求解过程，对于不同的裂纹长度a，根据裂纹扩展平面应力分布规律将平面应力分布表示为分段函数形式，且分段函数为(1-x/a)的多项式，即：

其中，s0、s1、s2为给定裂纹长度a下，依据应力分布数据确定的拟合参数。当拟合参数确定后，分别求解均布、线性、二次应力分布下的裂纹扩展应力强度因子，并根据叠加原理计算三者共同作用时的应力强度因子。

基于上述发明方法，计算涡轮盘偏心孔处径向裂纹的裂纹扩展应力强度因子，初始1/4圆形角裂纹半径2mm，每次计算间隔给定裂纹扩展增量0.1mm，得到辐板表面裂纹长度(a)与偏心孔内表面裂纹长度(b)随循环数n的变化关系，并与有限元法进行比较，如图3所示，本发明提出的方法与有限元法预测结果偏差较小且计算效率更高。

提供以上实施例仅仅是为了描述本发明的目的，而并非要限制本发明的范围。本发明的范围由所附权利要求限定。不脱离本发明的精神和原理而做出的各种等同替换和修改，均应涵盖在本发明的范围之内。