技术领域
本发明涉及蠕变疲劳的寿命预测领域，尤其涉及一种基于晶体塑性理论的蠕变疲劳寿命预测方法。
背景技术
现代航空发动机、工业燃气机、超临界发电机等装置的工作环境日趋复杂，其中的关键限寿部件除了承受装置稳态运行的恒定荷载之外，通常还承受装置起停和温度波动等引起的交变载荷的作用，其服役过程伴随着严重的蠕变疲劳载荷交互作用，降低了这些关键限寿部件的服役寿命。因此，为了保证关键限寿部件的安全性和可靠性，构建准确的蠕变疲劳寿命预测模型是非常重要的。
为了准确地评估高温结构在蠕变疲劳载荷下的耐久性问题，自20世纪50年代以来，国内外学者已先后提出了几十种蠕变疲劳寿命预测模型。近些年，有限元软件的发展为了更好地理解复杂工况下的应力应变行为以及实现该工况下的寿命预测。ABAQUS软件允许用户通过二次开发接口补充ABAQUS蠕变前后处理模块中无法实现的功能。其中，基于FORTRAN语言的材料用户自定义的子程序UMAT扩展了ABAQUS在本构方程方面上的应用。用户可以基于UMAT根据自己的材料以及需求，来编译所需的本构模型并通过ABAQUS软件实现模拟。
现今针对蠕变疲劳分析及寿命预测的主要能够描述包含裂纹萌生和扩展阶段的连续损伤力学理论，其通过损伤变量引入统一粘塑性理论的方式描述材料在循环载荷下损伤累积直至断裂的过程。这类方法侧重于宏观层面的类构件或者构件的损伤分析，无法描述微观层面的应力应变响应以及损伤机理演化。由于蠕变疲劳裂纹萌生通常是在微观层面上，使用这类方法很好精确地判定裂纹萌生以及揭示损伤机理。因此，如何有效地评估裂纹萌生寿命以及分析蠕变疲劳损伤机理，对关键限寿部件是非常重要的。
发明内容
针对上述现有技术中的不足，本发明提供一种基于晶体塑性的蠕变疲劳寿命预测方法，能够更好地实现蠕变疲劳损伤分析以及蠕变疲劳剩余寿命预测，具有直观、适用性强、精确度高的优点。
为了实现上述目的，本发明提供一种基于晶体塑性的蠕变疲劳寿命预测方法，其特征在于，包括以下步骤：
S1：建立基于材料的电子背散射衍射信息的ABAQUS的代表性单元模型；
S2：将背应力模型修正为含有静态回复项的背应力模型并写入晶体塑性本构方程，该晶体塑性本构方程通过用户自定义的子程序UMAT定义，将该晶体塑性本构方程用于对所述S1中的ABAQUS代表性单元模型进行有限元计算，以此模拟得到所述材料在给定载荷条件下的蠕变疲劳迟滞回线；
S3：通过试参法拟合蠕变疲劳迟滞回线，来获取所述晶体塑性本构方程的材料参数；
S4：通过ABAQUS软件计算每个积分点的应力应变值，并将应力应变值平均化，获取材料的蠕变疲劳迟滞回线和后处理云图；
S5：从后处理云图中提取最大的塑性滑移和能量耗散，将这两个参数视为材料在循环加载过程中损伤驱动力，分析它们随着循环周次的变化规律，从而提出一种蠕变指示因子和疲劳指示因子；
S6：结合所述S5的蠕变指示因子和疲劳指示因子和线性累积损伤准则
进一步地，所述S1包括：
S11：通过matlab语言将材料的电子背散射衍射信息转化为ABAQUS可读取的input文件中的节点和单元以及相应的晶粒取向文件；
S12：通过matlab语言读取所述S11中的晶粒取向文件，从而编译每个晶粒的多个材料常数；
S13：通过编写ABAQUS的脚本程序获得代表性体积单元的周期性边界条件和蠕变疲劳加载条件；
S14：将S11中的节点和单元、S12中的材料属性、S13中的周期性边界条件和蠕变疲劳加载条件手动组成ABAQUS可读取的input文件，从而获得ABAQUS代表性体积单元模型。
优选地，在所述S11中，所述晶粒的取向参数采用三个欧拉角
进一步地，在所述S2中，晶体塑性本构方程包括主控方程，流动准则方程，滑移阻力演化方程以及引入静态回复项的背应力方程。。
进一步地，所述S2包括：
S21：建立描述蠕变疲劳行为的晶体塑性本构方程中变形速率梯度L
其中，L
S22：建立晶体塑性本构方程中的流动准则方程，流动准则方程为：
其中，
S23：建立晶体塑性本构方程中的滑移阻力演化方程，滑移阻力演化方程为：
h
其中，h
S24：建立晶体塑性本构方程中引入静态回复项的背应力方程，引入静态回复项的背应力方程为：
其中，h
进一步地，在所述S3中，通过试参法拟合蠕变疲劳迟滞回线是通过将S2得到的材料在给定载荷条件下的蠕变疲劳迟滞回线与该材料的蠕变疲劳试验的数据对比来实现的。
进一步地，在在所述S4中，结合所述S1中建立的ABAQUS的代表性单元模型、S2中构建的晶体塑性本构方程、S3中确定的晶体塑性本构方程的材料参数，并通过ABAQUS软件来计算得到每个积分点的应力应变值。
进一步地，所述步骤S5包括：
S51：从ABAQUS后处理云图中提取蠕变疲劳保载时间前后的最大塑性滑移和能量耗散，分析保载时间前后的塑性滑移和能量耗散的变化规律；
S52：将所述步骤S51中的塑性滑移和/或能量耗散视为蠕变指示因子和疲劳指示因子。
进一步地，所述步骤S6中，预测不同保载时间下的蠕变疲劳裂纹萌生寿命包括：构建蠕变疲劳裂纹萌生寿命预测模型和确定蠕变疲劳裂纹萌生寿命。
进一步地，所述S6包括：
S61：构建蠕变疲劳裂纹萌生寿命预测模型，所述蠕变疲劳裂纹萌生寿命预测模型为：
其中，
塑性滑移p和能量耗散W的公式为：
其中，
S62：确定蠕变疲劳裂纹萌生寿命，所述蠕变疲劳裂纹萌生寿命为：
其中，X
进一步地，所述基于晶体塑性的蠕变疲劳寿命预测方法还包括步骤S7：基于累积塑性滑移占比构建蠕变疲劳剩余寿命评定图。
本发明的基于晶体塑性的蠕变疲劳寿命预测方法利用塑性滑移和能量耗散作为疲劳指示因子和蠕变指示因子，可以将微观尺度的参数用做损伤参量来进行寿命的评估，这种方法的好处在于微观尺度上的预测，可以有效地获得蠕变疲劳裂纹萌生位置和精确地预测蠕变疲劳裂纹萌生寿命。本发明不仅可以预测寿命，还能够获得每个循环周次中的塑性滑移损伤积累图，由此预测裂纹萌生的位置。
本发明的基于晶体塑性的蠕变疲劳寿命预测方法能够得到材料在蠕变疲劳载荷作用下的应力应变行为以及损伤演化规律，其中晶体塑性理论耦合了修正的背应力模型并写入材料用户自定义的子程序UMAT。
附图说明
图1为根据本发明的一个实施例的基于晶体塑性的蠕变疲劳寿命预测方法的流程图；
图2为根据本发明的一个实施例的基于晶体塑性的蠕变疲劳寿命预测方法的ABAQUS代表性体积单元模型图；
图3为根据本发明的一个实施例的基于晶体塑性的蠕变疲劳寿命预测方法的有限元模型的周期性边界条件和加载条件；
图4为根据本发明的一个实施例的基于晶体塑性的蠕变疲劳寿命预测方法的蠕变疲劳模拟曲线的拟合结果图；
图5为根据本发明的一个实施例的基于晶体塑性的蠕变疲劳寿命预测方法的累积非弹性滑移和能量耗散随着循环周次线性增加的轨迹图；
图6给出了累积非弹性滑移和能量耗散随着蠕变疲劳循环周次变化的轨迹图，其中，图6(a)示出了累积非弹性滑移随着蠕变疲劳循环周次变化的轨迹，图6(b)示出了能量耗散随着蠕变疲劳循环周次变化的轨迹；
图7为根据本发明的一个实施例的基于晶体塑性的蠕变疲劳寿命预测方法的试验寿命和预测寿命的对比图，其中，图7(a)示出了基于非弹性滑移的预测结果；图7(b)示出了基于能量耗散的预测结果。
具体实施方式
下面根据附图，给出本发明的较佳实施例，并予以详细描述，使能更好地理解本发明的功能和特点。
请参阅图1为本发明公开的一种基于晶体塑性的蠕变疲劳寿命预测方法，其包括以下步骤：
S1：通过背散射电子衍射技术获得原始材料的EBSD(电子背散射衍射)信息，EBSD信息主要包括晶体取向和晶体结构信息，基于材料的EBSD信息建立ABAQUS代表性单元模型；
其中，S1具体包括：
S11：通过matlab语言将材料的EBSD信息转化为ABAQUS可读取的input文件中的节点和单元以及相应的晶粒取向文件，其中，晶粒的取向参数采用三个欧拉角
S12：通过matlab语言读取所述步骤S11中的晶粒取向文件，从而编译每个晶粒的多个材料常数；在本实施例中，每个晶粒的材料常数的数量为15个，这些材料常数用于在ABAQUS计算过程中供UMAT调用，用于控制UMAT中背应力、滑移阻力和滑移系等的开启，并用于反映不同晶粒的不同取向；
S13：通过编写ABAQUS的脚本程序获得代表性体积单元的周期性边界条件和蠕变疲劳加载条件；
S14：将S11中的节点和单元、S12中的材料常数、S13中的周期性边界条件和蠕变疲劳加载条件手动组成ABAQUS可读取的input文件，从而获得ABAQUS代表性体积单元模型，用于蠕变疲劳有限元模拟。
S2：将背应力模型修正为含有静态回复项的背应力模型并写入晶体塑性本构方程，晶体塑性本构方程通过ABAQUS有限元软件提供的材料用户自定义子程序UMAT编译，由此，修正后的晶体塑性本构方程包括主控方程，流动准则方程，滑移阻力演化方程以及引入静态回复项的背应力方程，该晶体塑性本构方程用来描述材料的蠕变疲劳迟滞曲线；最后，通过ABAQUS软件调用UMAT子程序，将该晶体塑性本构方程用于对步骤S1中的ABAQUS代表性单元模型进行有限元计算，一方面得到材料在给定载荷条件下的蠕变疲劳迟滞回线，另一方面获得材料的塑性滑移和能量耗散；
其中，步骤S2进一步包括步骤：
S21：建立描述蠕变疲劳行为的晶体塑性本构方程中变形速率梯度L
其中，L
S22：建立晶体塑性本构方程中的流动准则方程：
其中，
S23：建立晶体塑性本构方程中的滑移阻力演化方程：
h
其中，h
S24：由于考虑蠕变回复效应，建立引入静态回复项的背应力方程：
其中，h
S3：通过将S2得到的材料在给定载荷条件下的蠕变疲劳迟滞回线与该材料的蠕变疲劳试验的数据对比以通过试参法拟合蠕变疲劳迟滞回线，从而通过试参法拟合蠕变疲劳迟滞回线，来获取晶体塑性本构方程所需的材料参数
S4：通过ABAQUS软件计算每个积分点的应力应变值，并将应力应变值平均化，获取材料的蠕变疲劳迟滞回线和后处理云图。由此，得到了材料的宏观的力学响应。
在所述S4中，结合所述S1中建立的ABAQUS的代表性单元模型、S2中构建的晶体塑性本构方程、S3中确定的晶体塑性本构方程的材料参数，并通过ABAQUS软件来计算得到每个积分点的应力应变值。
S5：从后处理云图中提取最大的塑性滑移和能量耗散，将这两个参数视为材料在循环加载过程中损伤驱动力，分析它们随着循环周次的变化规律，从而提出一种蠕变指示因子和疲劳指示因子。
其中，所述塑性滑移和能量耗散作为ABAQUS的状态变量通过后处理输出。
所述步骤S5包括：
S51：从ABAQUS后处理云图中提取蠕变疲劳保载时间前后的最大的塑性滑移和能量耗散，分析保载时间前后的塑性滑移和能量耗散的变化规律；
S52：将所述步骤S51中的塑性滑移和/或能量耗散视为蠕变指示因子和疲劳指示因子(即蠕变疲劳损伤失效的指示因子)。
S6：基于所述S5的蠕变指示因子和疲劳指示因子和线性累积损伤准则，预测不同保载时间下的蠕变疲劳裂纹萌生寿命；
其中，线性累积损伤准则就是把累积疲劳损伤和蠕变损伤进行线性叠加，即
其中，步骤S6进一步包括：
S61：构建蠕变疲劳裂纹萌生寿命预测模型，所述蠕变疲劳裂纹萌生寿命预测模型为：
其中，
塑性滑移p和能量耗散W的公式为：
塑性滑移p和能量耗散W的计算公式编写在用户子程序UMAT中，通过ABAQUS有限元软件计算，提取蠕变疲劳每个周次下以及保载时间前后的塑性滑移p和能量耗散W，用于获得蠕变指示因子的临界值X
S62：结合S61所述的蠕变疲劳裂纹萌生寿命预测模型，确定蠕变疲劳裂纹萌生寿命，蠕变疲劳裂纹萌生寿命为：
其中，X
此外，还可以包括S7：基于累积塑性滑移占比构建蠕变疲劳剩余寿命评定图。
实验结果
下面采用本发明提供的基于晶体塑性的蠕变疲劳寿命预测方法针对镍基高温合金GH4169，并采用实验的裂纹萌生寿命信息来验证本发明的有效性。其中，蠕变疲劳试验在高温650℃下进行，棒状试样的加载方式为应变控制，应变比为-1，加载的应变范围为1％、1.4％、1.6％、2％，保载时间施加在最大拉应变处，保载时间分别为120s、300s、1800s。
本发明的基于晶体塑性的蠕变疲劳寿命预测方法需要对同样材质的试样进行高温下的不同应变幅的单轴蠕变疲劳试验，获得的试验结果用于在本发明的蠕变疲劳寿命预测方法的步骤S3中确定晶体塑性本构方程所需的材料参数。首先基于EBSD信息构建如图2所示的ABAQUS代表性体积单元模型，代表性体积单元模型的约束条件和加载方式如图3所示。通过试参法调整蠕变疲劳试验的模拟结果，使之与蠕变疲劳试验的滞回曲线和应力松弛曲线相吻合，拟合结果如图4和图5所示。其中，拟合结果为：
通过采用本发明的基于晶体塑性的蠕变疲劳寿命预测方法，可以根据蠕变疲劳损伤的演化规律，确定蠕变指示因子的临界值X
图6给出了累积非弹性滑移(即累积塑性滑移)和能量耗散随着蠕变疲劳循环周次变化的轨迹图，其中，其中图6(a)示出了累积非弹性滑移随着蠕变疲劳循环周次变化的轨迹，图6(b)示出了能量耗散随着蠕变疲劳循环周次变化的轨迹。可以看出，随着保载时间的增加，累积非弹性滑移和能量耗散也相应增加。此外，累积非弹性滑移和能量耗散均与循环周次呈线性关系，因此可以用步骤S6所述公式(12)计算不同工况下的蠕变疲劳萌生寿命。
图7给出了不同蠕变疲劳工况下的裂纹萌生寿命预测图，其中，图7(a)示出了基于非弹性滑移的预测结果；图7(b)示出了基于能量耗散的预测结果。可以看出，本发明的蠕变疲劳寿命预测方法通过提出的蠕变指示因子和疲劳指示因子来预测萌生寿命的结果，与实验的结果较为接近，在2倍误差带范围内，此外，基于能量耗散的寿命预测结果比基于非弹性滑移的预测结果更加精确。因此，证明了此数值模拟方法的具有较高的可靠性。
以上记载的，仅为本发明的较佳实施例，并非用以限定本发明的范围，本发明的上述实施例还可以做出各种变化。即凡是依据本发明申请的权利要求书及说明书内容所作的简单、等效变化与修饰，皆落入本发明专利的权利要求保护范围。
机译： 一种可模制的热塑性热塑性聚氨酯，其制备方法也是如此，其制品具有高模量，低蠕变和高疲劳寿命
机译： 基于威布尔函数和最大疲劳变形的混凝土疲劳寿命预测方法及装置
机译： 基于威布尔函数和残差变形的混凝土疲劳寿命预测方法及装置
