基于θ-映射法的ZSGH4169合金缺口试样蠕变变形模拟

李 迪1，魏大盛

(1.中国航发商用航空发动机有限责任公司，设计研发中心，上海 201108；2.北京航空航天大学能源与动力工程学院，北京 100191)

0 引 言

蠕变是指在一定温度下金属受持续应力作用而产生缓慢塑性变形的现象，是温度、应力和时间共同作用的结果[1]。蠕变研究一般从两方面着手：一是从微观角度出发，研究蠕变机理及冶金因素对蠕变特性的影响，以提高金属的蠕变抗力；二是从唯象角度出发，以试验为基础，根据观察到的宏观蠕变现象建立蠕变方程，并寻求工程上的应用。从蠕变速率来看，蠕变曲线包括速率减小、速率恒定、速率增加3个阶段。研究人员根据材料在蠕变试验中表现出来的不同特征，建立了多种蠕变模型，如Norton率模型[2]以及在其形式上发展而来的时间硬化/应变硬化模型、Omega模型[3-4]、θ-映射法模型[5-6]等。这些模型特点鲜明且具有特定的适用范围：Norton率模型只能描述第一阶段和第二阶段的蠕变变形；Omega模型是针对特征明显的第三阶段蠕变曲线建立的；θ-映射法模型可以描述完整的蠕变过程，且参数组合较为灵活，在工程上应用较多[7]。

对于实际工程结构，如航空燃气轮发动机，其热端部件(涡轮盘、涡轮叶片等)在高温、高载荷作用下通常会发生蠕变变形，成为结构失效的隐患。因此，航空燃气轮发动机设计时需要对蠕变行为进行深入分析[8]，并针对蠕变可能导致的失效模式给出相应的设计准则。然而，螺栓孔边、圆角、榫头/榫槽连接位置等应力集中部位的载荷通常相对较高，易产生较大的蠕变变形和应力松弛，从而对传统疲劳寿命评价方法产生影响。准确计算这些部位的蠕变变形及松弛应力，是开展工程结构件持久寿命[9-10]、蠕变疲劳交互影响[11-12]等研究的前提和基础。同时，由于缺口处的蠕变变形无法准确测量，采用数值分析方法给出缺口蠕变变形特征也是一项重要的基础研究。为此，作者基于ZSGH4169合金的蠕变数据，采用16参数的θ-映射法模型对不同温度、不同载荷下的蠕变曲线进行分析，在获得方程参数的同时考核模型的计算精度，并设计了具有不同应力集中系数的缺口试样，以净截面应力作为参考标准，采用θ-映射法模型对其蠕变变形及应力松弛行为进行了模拟。

1 试验材料

试验材料为经时效处理的GH4169合金锻件，牌号记为ZSGH4169，其化学成分如表1所示[13]。

表1 ZSGH4169合金的化学成分(质量分数)

试样外轮廓尺寸为80 mm×18 mm×2 mm，工作段宽度为6 mm，与两端销钉夹持段以半径为12 mm的圆角过渡，中部开单边V型或双边半圆型缺口(图1)，缺口尺寸如表2所示。

缺口试样应力参数如表3所示。其中：σ0为计算应力集中系数时的远端应力；σn为缺口净截面应力，单边缺口sn计算公式为s0/(1-DW),双边缺口sn计算公式为s0(1-2D/W),W为试样宽度，D为缺口深度；σmax为缺口根部最大应力，通过有限元计算获得；Kt为应力集中系数，为最大应力与静截面应力的比值。

图1 ZSGH4169缺口试样的形状与尺寸Fig.1 Shape and size of ZSGH4169 notched specimen

表2 缺口试样的缺口类型与尺寸

表3 缺口试样的应力参数

2 蠕变模型的建立

根据蠕变曲线特征，θ-映射法模型将蠕变过程假设为强化和弱化过程的叠加，蠕变速率恒定的第二阶段可以看作是强化和弱化平衡的过程。基于此，能够描述3阶段蠕变曲线的方程为

εc=θ1[1-exp(-θ2t)]+θ3[exp(θ4t)-1]

(1)

式中：εc为蠕变应变；t为蠕变时间；θi(i=1，2，3，4)为材料参数。

第一项描述了材料的硬化；第二项描述了材料的弱化；参数θ1和θ3描述了应变的变化量；θ2和θ4描述了蠕变速率。

θi由应力σ和温度T决定，关系式为

lgθi=ai+biσ+ciT+diσT

(2)

式中：ai，bi，ci，di(i=1，2，3，4)为材料参数。

在模型建立过程中，应力以材料在相应温度下的屈服强度作归一化，温度则以材料的熔点温度作归一化，这样优化得到的方程参数在数量级上不会相差过大，且计算过程中的蠕变变形数值更加稳定。以ZSGH4169合金在650，700 ℃时不同应力水平下的蠕变试验数据[13]为基础，将参数的获取过程视为多参数优化问题(共计16个参数，与温度、载荷以及时间相关)，选择Levenberg-Marquardt优化算法，采用分步优化策略[14]得到方程参数，如表4所示。参数获取过程是一个多参数优化过程，不同优化策略将获得不同的结果，因此存在多组参数满足要求，这里选择1 组计算精度相对较高的参数进行后续计算。

分别采用Fortran程序及MSC.MARC有限元程序对ZSGH4169合金蠕变数据进行模拟计算。Fortran程序的计算目的是实现θ-映射法，而有限元程序则是为工程结构的蠕变分析服务。计算过程中，两程序均采用表4中的参数。由图2可以看出：采用Fortran程序的计算结果(空心点)与有限元程序结合子程序(MSC.MARC)的计算结果(实线)完全吻合，说明程序编写准确；θ-映射法采用上述参数可较好地模拟出ZSGH4169合金在不同温度、不种应力下的蠕变变形行为，且能够完整描述出蠕变各个阶段的特征，特别是加速损伤的第三阶段。

表4 ZSGH4169合金的θ-映射法模型参数

图2 不同温度及应力下ZSGH4169合金蠕变试验结果与θ-映射法模型计算结果Fig.2 Creep test and θ-projection method calculation results of ZSGH4169 alloy under different temperatures and stresses

3 缺口试样的蠕变变形模拟

3.1 蠕变变形计算方法

使用表4中θ-映射法模型参数计算各缺口试样的蠕变变形，力边界条件见表5。不同缺口试样均以净截面应力作为对照标准，应力水平分别取300，375，450 MPa，结合5种应力集中系数(包含Kt=1.00的情况，即无缺口试样)，进行15种不同情况下的蠕变变形计算，蠕变计算时间取1 500 h。

表5 缺口试样在不同净截面应力下的远端应力

3.2 V型缺口试样的蠕变变形

对缺口有限元模型进行网格划分(图3)，重点分析沿缺口平分线向根部逐步深入的5个节点(编号为2328，2510，2509，2508，2507)的蠕变变化规律。由图4可以看出，越靠近试样缺口表面，V型缺口试样(Kt=3.13)的蠕变应变越大，应力松弛越显著。由于缺口位置为多轴状态，计算中给出的应力均指等效应力。拉伸方向(y方向)的应力松弛在大小上略有不同。

图3 V型缺口试样有限元模型沿缺口根部平分线上的节点分布示意Fig.3 Schematic of nodes distribution of finite element model of V-notched specimen along symmetric path at notch root

图4 650 ℃、净截面应力300 MPa下V型缺口试样不同节点的蠕变曲线及应力松弛曲线模拟结果Fig.4 Simulation results of creep (a) and stress relaxation (b) curves of different nodes of V-notched specimen undernet section stress of 300 MPa at 650 ℃

由图5可以看出：随着净截面应力的提高，V型缺口试样表面节点(编号为2328)的蠕变应变增加；由于应力提高使得高应力区面积增加，试样应力松弛的时间较低应力下的略有滞后。

3.3 不同缺口试样的蠕变变形

由图6可知，随着应力集中系数的增大，试样蠕变变形程度增加，应力松弛现象显著。无缺口试样(Kt=1.00)的蠕变应变最小，该条件下应力是恒定的，因此应力松弛曲线中不包含Kt=1.00的情况。

图5 650 ℃、不同净截面应力下V型缺口试样表面节点的蠕变曲线及应力松弛曲线模拟结果Fig.5 Simulation results of creep (a) and stress relaxation (b) curves of surface node of V-notched specimen underdifferent net section stresses at 650 ℃

图6 650 ℃、净截面应力300 MPa下不同缺口试样表面节点的蠕变曲线及应力松弛曲线模拟结果Fig.6 Simulation results of creep (a) and stress relaxation (b) curves of surface nodes of different notched specimens undernet section stress of 300 MPa at 650 ℃

图7 650 ℃、净截面应力450 MPa下不同试样表面节点的蠕变曲线及应力松弛曲线模拟结果Fig.7 Simulation results of creep (a) and stress relaxation (b) curves of surface nodes of different notched specimens undernet section stress of 450 MPa at 650 ℃

由图7可以看出，650 ℃、净截面应力为450 MPa时，缺口试样蠕变变形计算结果的变化趋势与净截面应力300 MPa下的类似。应力集中系数为2.08的缺口试样的蠕变曲线计算结果与无缺口试样的十分接近，这是由于其应力集中程度及缺口应力梯度不高，缺口根部曲率半径与其他试样的相比相对较大。上述分析说明根据不同载荷及试样几何特征，可以尝试采用无缺口试样来预测缺口试样的蠕变变形。

对无缺口试样和应力集中系数为2.08的缺口试样蠕变行为进行进一步分析。650 ℃、净截面应力为300 MPa下，t=0时，缺口试样表面节点的等效应力为606.4 MPa，拉伸方向应力为624.8 MPa；t=1 500 s时，表面节点等效应力为401.5MPa，拉伸方向应力为418.6 MPa。以t=1 500 s时的应力作为无缺口试样(Kt=1.00)的拉伸应力，以考察持续变载(应力松弛)对蠕变变形的影响。由图8可以看出，应力松弛对蠕变变形计算结果产生了较大影响，相同恒定应力条件下的蠕变应变有所不同。

图8 300 MPa净截面应力下缺口试样及不同拉伸应力下无缺口试样的蠕变曲线计算结果Fig.8 Calculation results of creep curves of notched specimen under net section stress of 300 MPa and non-notched specimens under different tensile stresses

试验所建立的模型可以通过修正θ-映射公式来反映上述影响。缺口试样的蠕变变形同应力集中系数(与缺口净截面应力、表面最大应力相关)和应力梯度(与缺口处几何形式相关)密切相关，因此缺口处蠕变计算公式可表示为

εc=f(σ，T，t)·f(Kt)·f(dσ/dx)

(3)

式中：f(σ，T，t)为传统蠕变模型；f(Kt)为针对应力集中系数(或缺口根部最大应力)的修正因子；f(dσ/dx)为针对应力梯度的修正因子。

以上仅给出了描述缺口蠕变可能的表达式，对于具体表达式和可行的试验验证方法尚未进行深入研究，且未对是否存在更好的优化方法及流程，使计算结果的精度更高进行探讨，这将是重要的后续研究工作。

4 结 论

(1)θ-映射法模型可通过获得的一组参数描述不同温度、不同应力下ZSGH4169合金的蠕变行为，计算精度较高，在提高计算效率的基础上，能够用于实际工程结构的蠕变分析。

(2) 应力集中程度越高，应力梯度越大，ZSGH4169合金缺口根部的蠕变变形越明显，应力松弛程度越高；应力集中系数较小、应力梯度较为缓和时，缺口根部的蠕变变形呈现与无缺口试样相接近的趋势。

(3) 引入应力集中系数(或最大应力)、应力梯度影响因子对传统蠕变模型进行修正，给出了缺口试样蠕变变形计算的一般表达式。