镍基高温合金IC10的本构模型研究

畅一鹏，张宏建，卢孔汉，温卫东，崔海涛

(南京航空航天大学 航空发动机热环境与热结构工业和信息化部重点实验室，江苏 南京 210016)

0 引言

随着新一代航空发动机的发展，现役的合金材料已逐步开始采用镍基高温合金。IC10作为国产镍基定向凝固高温合金的代表之一，目前已被广泛应用于诸如涡轮发动机等先进动力推进系统的热端部件中，因而开展其本构关系的建模预测研究有助于提高材料的工程应用价值。

目前国内外对镍基高温合金均已不同程度地开展了本构关系研究。岳珠峰等[1]对DD3合金开展了不同温度下的拉伸、蠕变本构关系研究，采用弹塑性晶体滑移理论建立了相应的本构关系。蔚夺魁[2]基于晶体塑性理论并对材料非线性运动硬化采用变量描述，在滑移剪切率演化方程中引入背应力项，建立了GH4169合金的高温循环应变硬化晶体塑性本构模型。上述研究均是基于晶体塑性框架并结合经典弹塑性本构理论而建立，而如果晶体塑性滑移理论从微观的晶体变形出发，通过描述晶体滑移系的开动以及位错等内变量的演化，则能够更加全面、准确地预测在不同载荷下的材料力学行为。

目前针对镍基高温合金IC10合金研究主要集中于单轴静强度实验、高温下蠕变机理以及动态回复方面研究[3-4]，而针对IC10合金的疲劳力学行为以及相应的循环塑性本构关系却鲜有报道。本文在晶体塑性理论的框架下，应用率相关的硬化方程，编写UMAT用户子程序，针对IC10合金在600℃下单轴拉伸以及循环力学行为开展数值模拟研究，讨论了晶体塑性滑移模型用于IC10合金高温疲劳力学行为的合理性。

1 晶体塑性本构关系的基本公式

材料变形应变梯度F可以由一个塑性部分Fp和一个刚体转动的弹性部分Fe相乘，得[3]：

F=Fe·Fp

(1)

图1为晶体变形的图像。其中晶体变形主要分为塑性与弹性两部分。

图1 晶体变形图

令L表示速度梯度：

(2)

各滑移系中由滑移引起的剪切应变与整体塑性变形之间的关系为：

(3)

速度梯度可以分解为一个对称部分和一个反对称部分：

(4)

上标“e”、“p”分别表示弹性速度梯度、塑性速度梯度。

(5)

(6)

(7)

(8)

以上公式为晶体塑性理论的基本公式，建立了滑移剪切率和宏观变形率之间的关系。

式(9)将应力率、变形率以及滑移剪切率联系到一起：

(9)

其中：

B(α)=W(α)σ-σW(α)

(10)

(11)

式中：EMT为瞬时弹性模量张量；D为变形率张量。

2 晶体塑性率相关本构关系

在上述晶体塑性理论中，求解应力-应变关系的关键就是解决滑移系剪切率的计算问题。

目前广泛应用的是率相关和率无关两种形式的晶体塑性滑移理论。其中率相关模型的滑移剪切率是唯一确定的，这在实际计算中，带来很多方便之处。

基于各滑移系中切应力τα提出率相关的幂函数硬化方程：

(12)

式(12)中，滑移系的分切应力与宏观应力之间的关系可表示为：

τα=σ/Pα

(13)

式中：τα为各个滑移的分切应力；Pα为取向因子；σ为晶轴系下的应力张量。

从式(12)中可以看出gα的演化在描述材料的本构模型非常重要。本文中采用gα的硬化函数公式为：

(14)

式中hαβ称为硬化系数，它决定了滑移系β中的滑移剪切率对滑移系α所造成的硬化。硬化系数是变形历史、变形温度和速度的函数，本文采用的形式为：

(15)

式中：qαβ是潜在硬化的极值；h0为单滑移系硬化率初值；τs为饱和参考剪切应力；β为硬化指数；k为温度影响因子[4]。

3 本构模型的有限元实施

本文采用商用有限元软件ABAQUS进行有限元分析，通过用户子程序UMAT接口，编制了晶体塑性模型完整算法，计算流程图如图2所示。

图2 UMAT流程图

IC10合金是一种典型的多相L12型材料，由于本文中只考虑在单轴[001]方向的循环行为，因此在模拟时，可认为IC10是正交各向异性材料，且循环行为中只有八面体滑移系上的滑移运动。在本文中结合单轴拉伸试验数据，参考温志勋[5]、周杰[6]改进的模型参数计算方法来获取本文模型参数。

本文基于晶体塑性理论编制了改进的本构模型程序，通过用户子程序UMAT嵌入有限元软件ABAQUS中，采用表1的材料参数对IC10合金高温环境下在不同应变下的非弹性响应力学行为开展数值模拟，其中IC10合金元胞建模如图3所示，在此基础上分别加载单轴以及循环载荷。本文中采用假设均匀化模型，因此在施加周期性边界条件后，可认为各单元受力均匀一致。

表1 IC10合金600 ℃下的参数赋值

图3 元胞模型示意图

3.1 单轴拉伸数值模拟分析

材料在600℃下受到单向拉伸载荷，计算模型采用的是均匀单胞模型，采用的是周期性边界条件，并且在[001]方向上施加拉伸载荷。根据实验状况，本次施加的拉伸载荷为位移载荷，u=0.06a。其中a为单胞的边长。

通过对几何模型施加载荷，并进行应力分析，绘制应力-应变曲线。模拟结果与试验结果对比如图4所示。可以看出，在单向拉伸条件下，模拟曲线与试验曲线吻合度较高，说明了模型的有效性以及数值仿真方法的可靠性。

图4 单轴拉伸下模拟值与实验值对比

3.2 单轴循环数值模拟分析

在循环模拟中依旧采用图3给出的均匀单胞模型，在[001]方向上施加循环位移载荷。根据疲劳试验的实际加载波形，本文在模拟过程中采用三角波形，单个循环加载幅值曲线如图5所示。对几何模型分别施加0.75%、0.85%、1%的循环载荷，绘制循环应力-应变曲线，并与文献[7]实验曲线比较结果如图6所示。从图中可以看出，采用的晶体塑性理论可以合理地模拟IC10高温合金循环加载下的力学行为，特别在材料弹性卸载与加载阶段，数值模拟曲线与实验曲线吻合度非常高，但是在屈服阶段还是存在一定的误差。误差的主要来源在于高温下的循环加载下，材料在弹性卸载、加载的过程中会有塑性的影响。此外，滑移阻力、滑移变形等微观变量在循环加载下的变化具有不确定性，使得预测结果与试验结果存在部分误差。因此，本文建立的本构方法有待进一步改进，在率相关的硬化方程中考虑更多的循环运动变量影响因子，或在本构模型中添加与运动路径相关的变形函数，会更好地反映材料的变形过程。

图5 单个循环周期加载幅值曲线图

图6 IC10合金600℃下不同应变水平数值模拟结果

4 结语

本文基于晶体塑性理论建立了相应的IC10合金本构模型，采用率相关的硬化方程，编制了本构模型的UMAT用户材料子程序，利用ABAQUS软件开展了IC10合金在600℃、不同载荷下的循环应力-应变响应曲线的数值模拟研究。结果表明：采用晶体塑性理论能较好地描述IC10合金在特定温度下的复杂力学行为；同时也进一步验证了模型算法的完整性与计算程序的正确性。