一种基于晶体塑性理论的蠕变疲劳剩余寿命评定方法与流程

1.本发明涉及蠕变疲劳的寿命预测领域，尤其涉及一种基于晶体塑性理论的蠕变疲劳剩余寿命评定方法。


背景技术：

2.对于汽轮机、热交换器和核反应堆等关键限寿部件，其在运行过程中除了承受恒定荷载外，通常还会承受装置起停和温度波动等引起的交变载荷，因此在服役过程中可能伴随着严重的蠕变疲劳载荷交互作用，降低这些关键部件的服役寿命。此外，对低耗能和长寿命的需求不断增加，使得这些关键部件的运行条件更加苛刻，为了保证高温结构在蠕变疲劳下完整性和可靠性，建立精准的剩余寿命评估方法是非常重要的。
3.为了保证高温构件在服役过程中的安全性和可靠性，国内外学者已先后提出了不同类型的剩余寿命评定方法，可以分成两大类。第一类为宏观参数法，通过一些试验参数(如z参数、剩余强度等)建立唯象的本构方程去评估材料的剩余寿命，这类方法无法描述材料微观组织的演化规律；第二类为微观组织退化法，即通过一些代表性的微观组织参量在变形过程中的损伤演化规律来表征材料退化的过程，但这类方法无法定量地描述材料的损伤程度。因此如何精确地评估材料在服役过程的剩余寿命，尤其是剩余蠕变疲劳寿命，对这些关键限寿部件是非常重要的。
4.近些年，有限元软件的快速发展为更好地模拟复杂工况的力学响应提供了可靠的保障。其中，abaqus有限元软件允许用户通过二次开发接口补充abaqus前后处理模块中无法实现的功能。基于此，如何利用abaqus 有限元软件来评定材料的蠕变疲劳剩余寿命是亟需解决的问题。


技术实现要素：

5.针对上述现有技术中的不足，本发明提供一种基于晶体塑性理论的蠕变疲劳剩余寿命评定方法，能够描述材料微观组织的演化规律，并定量地描述材料的损伤程度，从而更好地实现蠕变疲劳损伤分析以及蠕变疲劳剩余寿命评估，具有直观、适用性强、精确度高的优点。
6.一种基于晶体塑性理论的蠕变疲劳剩余寿命评定方法，包括：
7.步骤s1，通过abaqus有限元软件的单元用户子程序uel编译基于位错密度的晶体塑性本构方程。
8.步骤s2，拟合蠕变疲劳迟滞回线以确定所述晶体塑性本构方程的材料参数。
9.步骤s3，根据所述晶体塑性本构方程以及所述晶体塑性本构方程的材料参数，对待测材料代表性体积单元进行晶体塑性有限元模拟，获得待测材料在不同蠕变疲劳工况下的力学响应。
10.步骤s4，根据待测材料在不同蠕变疲劳工况下的力学响应，提取待测材料代表性体积单元在每个循环周次下的累积能量耗散，确定待测材料的局部能量耗散体积占比以及
待测材料的容许蠕变损伤和容许疲劳损伤。
11.步骤s5，以所述容许蠕变损伤和所述容许疲劳损伤为基本坐标轴，以所述局部能量耗散体积占比为第三轴，绘制三维蠕变疲劳剩余寿命评定图，计算待测材料的蠕变疲劳剩余寿命。
12.进一步地，所述步骤s1包括：
13.步骤s11，建立晶体塑性本构方程中的流动准则方程如下：
[0014][0015]
式中，为参考塑性滑移率，k为玻尔兹曼常数，θ为绝对温度，f0为热激活自由能，τ
α
第α滑移系的分解剪切应力，表示为σ为应力值，s
α
表示α滑移系的滑移面法向，m
α
表示α滑移系的滑移方向，表示张量的并矢积，b
α
为第α滑移系的背应力，s
α
为第α滑移系的滑移阻力， p和q为材料常数。
[0016]
进一步地，步骤s12，根据所述流动准则方程，在所述晶体塑性本构方程中引入静态回复项的背应力方程如下：
[0017][0018]
式中，h
b
为背应力硬化常数，r
d
为滑移阻力相关动态回复系数，r
s
为背应力静态回复系数。
[0019]
进一步地，步骤s13，建立晶体塑性本构方程中的滑移阻力演化方程如下：
[0020][0021][0022][0023][0024]
式中，表示滑移阻力变化率，μ是剪切模量，b是柏氏矢量，h
αβ
表明滑移系α和β之间的硬化矩阵，和分别为滑移系β上的统计存储位错的刃型部分、螺型部分和几何必须位错部分，k
e
表示晶体在α滑移系上刃型位错产生的系数，k
s
表示晶体在α滑移系上螺型位错产生的系数，d
e
和d
s
分别表示刃型位错和螺型位错的湮灭系数，是滑移方向m
α
上的纯螺型密度分量，和分别是滑移方向t
α
和n
α
上的纯刃型密度分量。
[0025]
进一步地，步骤s14，建立无位错势垒的临界剪应力方程如下：
[0026]
[0027]
式中，为无位错势垒时的临界剪应力，c
τ
材料常数，n表示滑移系个数。
[0028]
进一步地，所述步骤s2中拟合蠕变疲劳迟滞回线的方法为试参法。
[0029]
进一步地，所述步骤s3包括：通过abaqus软件计算待测材料代表性体积单元的每个积分点的应力应变值，并将每个积分点的应力应变值平均化，获取待测材料在宏观尺度上的力学响应。
[0030]
进一步地，所述步骤s4中确定待测材料的局部能量耗散体积占比的方法包括：
[0031]
步骤s41，根据所述步骤s3中待测材料代表性体积单元的晶体塑性有限元模拟结果，提取待测材料在蠕变疲劳保载时间前后的累积能量耗散。
[0032]
步骤s42，根据所述步骤s41中的累积能量耗散，计算待测材料的局部能量耗散体积占比。
[0033]
所述步骤s41中待测材料在蠕变疲劳保载时间前后的累积能量耗散为：
[0034][0035]
式中，w为累积能量耗散，α表示α滑移系，τ
α
为第α滑移系的分解剪切应力，为参塑性滑移率，n为滑移系个数，t为时间。
[0036]
所述步骤s42中待测材料的局部能量耗散体积占比为：
[0037][0038][0039]
式中，v
led
为局部能量耗散体积，v
led，crit
为局部能量耗散体积的临界值， a和e分别是代表性体积单元中的局部累积能量耗散面积和相应单元个数。
[0040]
进一步地，所述步骤s5中的容许蠕变损伤和容许疲劳损伤分别为：
[0041][0042][0043]
式中，n
i
为蠕变疲劳裂纹萌生寿命，w
c，cyc
和w
f，cyc
分别为某一循环下的蠕变和疲劳指示因子，w
c，crit
和w
f，crit
分别为蠕变和疲劳指示因子的临界值。
[0044]
进一步地，所述步骤s6中待测材料的蠕变疲劳剩余寿命为：
[0045][0046][0047]
式中，为容许蠕变损伤，为容许疲劳损伤，n
r
为蠕变疲劳剩余寿命，n
i
为蠕变疲劳裂纹萌生寿命，表示归一化的局部能量耗散体积。
[0048]
本发明的一种基于晶体塑性理论的蠕变疲劳剩余寿命评定方法，能够得到材料在蠕变疲劳载荷下的应力应变行为以及损伤演化规律，其中基于位错密度的晶体塑性理论写
入材料用户自定义的子程序uel，并能够构建三维的蠕变疲劳剩余寿命评估图，具有直观、适用性强、精确度高的优点。
附图说明
[0049]
图1为根据本发明的一个实施例的基于晶体塑性理论的蠕变疲劳剩余寿命评定方法的流程图；
[0050]
图2为根据本发明的一个实施例的基于晶体塑性理论的蠕变疲劳剩余寿命评定方法的有限元模型的周期性边界条件和加载条件；
[0051]
图3为根据本发明的一个实施例的基于晶体塑性理论的蠕变疲劳剩余寿命评定方法的疲劳和蠕变疲劳模拟曲线的拟合图；
[0052]
图4为根据本发明的一个实施例的基于晶体塑性理论的蠕变疲劳剩余寿命评定方法的累积能量耗散随着循环周次线性增加的轨迹图；
[0053]
图5为根据本发明的一个实施例的基于晶体塑性理论的蠕变疲劳剩余寿命评定方法的三维蠕变疲劳剩余寿命评定图。
具体实施方式
[0054]
下面根据附图1
‑
5，给出本发明的较佳实施例，并予以详细描述，使能更好地理解本发明的功能和特点。
[0055]
本发明所使用的一种基于位错密度的晶体塑性有限元模型，需借助 abaqus软件提供的单元用户自定义子程序uel进行编译，最后通过调用编译后的uel实现蠕变疲劳复杂工况的晶体塑性有限元模拟。请参阅图1，本发明提供的一种基于晶体塑性理论的蠕变疲劳剩余寿命评定方法，包括以下步骤：
[0056]
步骤s1，通过abaqus有限元软件的单元用户子程序uel编译基于位错密度的晶体塑性本构方程。
[0057]
其中，步骤s1具体包括：
[0058]
步骤s11，建立晶体塑性本构方程中的流动准则方程如下：
[0059][0060]
式中，表示第α滑移系的塑性滑移率，为参考塑性滑移率，k为玻尔兹曼常数，θ为绝对温度，f0为热激活自由能，τ
α
为第α滑移系的分解剪切应力，表示为σ为应力值，s
α
表示α滑移系的滑移面法向，m
α
表示α滑移系的滑移方向，表示张量的并矢积，b
α
为第α滑移系的背应力，s
α
为第α滑移系的滑移阻力，p和q为材料常数；
[0061]
步骤s12，根据上述流动准则方程，在晶体塑性本构方程中引入静态回复项的背应力方程如下：
[0062]
[0063]
式中，表示第α滑移系的背应力变化率，h
b
为背应力硬化常数，r
d
为滑移阻力相关动态回复系数，r
s
为背应力静态回复系数。
[0064]
步骤s13，建立晶体塑性本构方程中的滑移阻力演化方程如下：
[0065][0066][0067][0068][0069]
其中，表示滑移阻力变化率，μ是剪切模量，b是柏氏矢量，h
αβ
表明滑移系α和β之间的硬化矩阵，和分别为滑移系β上的统计存储位错的刃型部分、螺型部分和几何必须位错部分，k
e
表示晶体在α滑移系上刃型位错产生的系数，k
s
表示晶体在α滑移系上螺型位错产生的系数，d
e
和d
s
分别表示刃型位错和螺型位错的湮灭系数，是滑移方向m
α
上的纯螺型密度分量，和分别是滑移方向t
α
和n
α
上的纯刃型密度分量。
[0070]
步骤s14，建立无位错势垒的临界剪应力方程如下：
[0071][0072]
其中，为无位错势垒时的临界剪应力，c
τ
材料常数，n表示滑移系个数。
[0073]
上述步骤s11
‑
s14中建立的公式(1)
‑
(7)共同构成晶体塑性本构方程
[0074]
步骤s2，通过试参法拟合蠕变疲劳迟滞回线，从而确定晶体塑性本构方程的材料参数。
[0075]
步骤s3，对待测材料代表性体积单元进行晶体塑性有限元模拟，获得待测材料在不同蠕变疲劳工况下的力学响应。具体包括：结合步骤s1中晶体塑性本构方程以及步骤s2中确定的材料参数，通过abaqus软件计算待测材料代表性体积单元的每个积分点的应力应变值，并将每个积分点的应力应变值平均化，获取待测材料在宏观尺度上的力学响应，通过与试验结果对比验证有限元模拟结果的准确性。
[0076]
步骤s4，提取待测材料代表性体积单元在每个循环周次下的累积能量耗散，确定待测材料的局部能量耗散体积占比以及待测材料的容许蠕变损伤和容许疲劳损伤。
[0077]
其中，确定待测材料的局部能量耗散体积占比的方法包括：
[0078]
步骤s41，提取待测材料在蠕变疲劳保载时间前后的累积能量耗散，分析保载时间前后累积能量耗散的变化规律。
[0079]
步骤s42，将步骤s41中的累积能量耗散看作蠕变疲劳的指示因子，按照下式计算待测材料的局部能量耗散体积占比：
[0080]
[0081][0082]
式中，w为累积能量耗散，v
led
为局部能量耗散体积，表示归一化的局部能量耗散体积，v
led，crit
为局部能量耗散体积的临界值，a和e分别是代表性体积单元中的局部累积能量耗散面积和相应单元个数，累积能量耗散 w公式为：
[0083][0084]
式中，w为累积能量耗散，α表示α滑移系，τ
α
为第α滑移系的分解剪切应力，为参塑性滑移率，n为滑移系个数，t为时间。
[0085]
步骤s5，根据待测材料在不同蠕变疲劳工况下的力学响应，提取待测材料代表性体积单元在每个循环周次下的累积能量耗散，确定待测材料的容许蠕变损伤和容许疲劳损伤。其中，容许蠕变损伤和容许疲劳损伤分别为：
[0086][0087][0088]
式中，n
i
为蠕变疲劳裂纹萌生寿命，w
c，cyc
和w
f，cyc
分别为某一循环下的蠕变和疲劳指示因子，w
c，crit
和w
f，crit
分别为蠕变和疲劳指示因子的临界值。
[0089]
步骤s6，以步骤s5计算出的容许蠕变损伤和容许疲劳损伤为基本坐标轴，以步骤s4计算出的局部能量耗散体积占比为第三轴，绘制三维的蠕变疲劳剩余寿命评定图，依据该评定图获得蠕变疲劳剩余寿命。其中，蠕变疲劳剩余寿命的计算公式为：
[0090][0091][0092]
式中，n
r
为蠕变疲劳剩余寿命。
[0093]
实验结果
[0094]
下面采用本发明提供的基于晶体塑性的蠕变疲劳寿命预测方法针对镍基高温合金inconel 718，并采用试验获得的裂纹萌生寿命信息来验证本发明的有效性。其中，蠕变疲劳试验在高温650℃下进行，棒状试样的加载方式为应变控制，应变比为
‑
1、0和
‑
∞，加载的应变范围为1.2％、1.6％、2.0％，保载时间施加在最大拉应变处，保载时间分别为60s、120s、300s。
[0095]
本发明的基于晶体塑性理论的蠕变疲劳剩余寿命评定方法需要对同样材质的试样进行不同工况下的高温下蠕变疲劳试验，获得的试验结果用于在本发明的蠕变疲劳剩余寿命评定方法的步骤s2中确定晶体塑性本构方程所需的材料参数。通过试参法调整蠕变疲劳试验的模拟结果，使之与蠕变疲劳试验的滞回曲线和应力松弛曲线相吻合，拟合结果如图2所示。其中，拟合结果为：f0＝295kj
·
mol
‑1，p＝0.96，q＝1.12，τ0＝485mpa，c
r
＝0.0027n/mm，ρ0＝1
×
8.18mpa，h
b
＝850mpa，r
d
＝ 8mpa，r
s
＝0.001/s，h0＝
360mpa，k
e
＝1240，d
e
＝24.8，d
s
＝124。
[0096]
图3给出了累积能量耗散随着蠕变疲劳循环周次变化的轨迹图。可以看出，随着保载时间的增加，累积能量耗散也相应增加。此外，累积能量耗散与循环周次呈线性关系。
[0097]
将提取的累积能量耗散代入公式(8)和公式(9)，可以获得局部累积能量耗散体积比，图4给出了局部累积能量耗散体积随着循环周次的变化图，可以看出，局部累积能量耗散体积随着循环周次线性增加。
[0098]
图5给出了三维的蠕变疲劳剩余寿命评定图，此图是在传统的蠕变疲劳损伤交互图基础上发展而来的，增加了以局部累积能量耗散体积占比为坐标的第三轴。可以看出，只需根据局部累积能量耗散体积占比，便可以很容易确定材料的剩余寿命。
[0099]
以上记载的，仅为本发明的较佳实施例，并非用以限定本发明的范围，本发明的上述实施例还可以做出各种变化。即凡是依据本发明申请的权利要求书及说明书内容所作的简单、等效变化与修饰，皆落入本发明专利的权利要求保护范围。