TC4钛合金是由美国在1954年研制成功的α-β型钛合金, 含有6%的α稳定元素Al和4%的β稳定元素V。由于具有优异的综合性能、良好的高温强度和蠕变抗力,TC4钛合金适合制造从-196℃到450℃温度范围内工作的各种零件。TC4钛合金主要用于制造飞机结构中的各种梁、框架、滑轨、起落架梁、航天火箭的壳体及压力容器、海洋工程中的近海石油钻井平台、电力工业中的燃气轮机叶片、造船工业中的船舶推进器等。研究表明, 在高温服役条件下, 蠕变断裂和拉伸断裂是金属构件的主要失效形式。因此, 研究TC4合金的黏塑性本构模型具有重要意义。
国内外诸多学者对TC4合金的本构模型进行了不同层面的研究和探讨。Jorge等[1]使用Johnson Cook、修正的Johnson Cook、Field Backofen、和修正的Arrhenius模型预测了TC4合金在高温下的拉伸行为, 研究表明修正的Arrhenius模型预测结果最好。Cheng等[2]提出了一种考虑损伤、塑性和应变率相关性的本构模型, 用于描述TC4合金拉伸行为。Paul等[3]发展了一个考虑初始晶粒尺寸的本构模型, 并通过该模型准确模拟了TC4合金在拉伸过程中的硬化和软化阶段。Wojciech等[4]使用Rusinek Klepaczko模型预测了TC4合金的准静态和动态的拉伸行为, 并发现该模型可以同时准确地预测TC4合金不同应变率下的拉伸行为。一些学者推导了本构模型的有限元数值实现方法, 并模拟了TC4合金的塑性成形过程[5-8]。此外, 损伤模型也被大量用在模拟金属材料的延性断裂问题中。Cao等[9]使用改进的Lemaitre损伤模型计算金属的塑性行为和低应力三轴度下的延性断裂问题。Soyarslan等[10]对Lemaitre损伤模型进行了修正, 并使用该模型模拟了金属材料的加工和断裂过程。
一些学者针对蠕变模型进行了深入地研究, Michele等[11]提出了一种基于多轴延性耗竭的连续损伤模型, 并对应力强度因子、C*积分等断裂参数进行了评价。周广磊等[12]建立了温度-应力耦合下的蠕变损伤模型, 并通过试验验证了该模型的准确性。还有学者[13-15]采用蠕变损伤本构模型研究了蠕变条件下裂纹的扩展行为, 并采用数值模拟的方法进行了预测。
目前针对TC4合金在高温下的本构模型研究主要集中在高温拉伸方面, 而可以将蠕变和高温拉伸同时计算的模型仍然较少, 因此本文提出一种不统一的考虑损伤的黏塑性本构模型用于计算TC4合金在高温下的拉伸和蠕变行为。
1 TC4试样蠕变及高温拉伸试验 1.1 TC4单轴蠕变试验 采用锻造TC4合金进行单轴蠕变试验, TC4钛合金的化学成分如表 1所示, 试验件几何尺寸如图 1所示。选取450℃下500和600 MPa 2种应力状态进行单轴蠕变试验。
表 1 TC4钛合金各元素的质量分数 %
Ti |
Al |
V |
Fe |
C |
H |
O |
rest |
5.5~6.75 |
3.5~4.5 |
≤0.30 |
≤0.08 |
≤0.015 |
≤0.20 |
TC4试样的蠕变曲线如图 2所示, 其在500和600 MPa下的蠕变断裂时间分别为17.3和4.2 h, 蠕变断裂应变分别为24.91%和23.40%。2种加载条件下试样的断裂应变相似, 但500 MPa下TC4试样的蠕变断裂时间为600 MPa下的4.12倍。
1.2 TC4高温拉伸试验 拉伸试验在25℃, 450℃和600℃下进行, 拉伸速率均为1 mm/min, 拉伸试样尺寸与蠕变试样尺寸相同。TC4试样的应力应变曲线如图 3所示, 随着温度的升高, 试样的应变硬化阶段减少, 应变软化阶段逐渐增加。对于600℃下的应力应变曲线, 应变软化阶段要远长于应变硬化阶段。TC4试样的拉伸性能如表 2所示。在25℃下, TC4试样弹性模量、屈服强度、抗拉强度和断裂应变分别为105 073.48, 846.42, 924.34 MPa和21.14%。与25℃下的拉伸性能相比, TC4试样在450℃和600℃下的弹性模量分别下降31.02%和59.43%;屈服强度分别下降39.40%和89.11%;450℃下试样的延性相比25℃下增加了12.86%, 而600℃下试样的高温延性很高, 变形量超出引伸计量程, 断面收缩率约为95%, 断裂应变远大于450℃下试样的断裂应变, 试样的断后形貌如图 4所示。
表 2 TC4试样拉伸性能
试验温度/℃ |
弹性模量/MPa |
屈服强度/MPa |
抗拉强度/MPa |
断裂应变/% |
25 |
105 073.48 |
846.42 |
924.34 |
21.14 |
450 |
72 481.88 |
513.01 |
605.02 |
24.26 |
600 |
42 632.66 |
92.20 |
241.39 |
|
2 本构模型与数值实现 2.1 统一的黏塑性本构模型 对于统一的黏塑性本构模型, 非弹性应变包括黏塑性应变和蠕变应变, 并通过一个统一的表达式进行计算。本文将Lemaitre损伤模型[16]加入到了Chaboche黏塑性本构模型[17]中用于预测试样的蠕变行为以及拉伸过程中的软化阶段。本构模型总结如下
|
(1) |
|
(2) |
|
(3) |
|
(4) |
s为应力偏张量, σeq为等效应力, σm为静水压力, αk为第k个背应力, 中的第三项用来描述时间硬化的回复作用。k0为初始屈服应力, Ck, γk, n, K, β, u, S, c为材料参数, p为累计塑性应变, D为损伤变量。
2.2 不统一的黏塑性本构模型 为比较统一与不统一模型的预测精度, 本文提出一种考虑蠕变损伤和随动硬化的不统一的黏塑性本构模型。将总应变率分解为两部分: 弹性应变率和非弹性应变率
|
(7) |
将非弹性应变率分解为蠕变应变率和塑性应变率之和
|
(8) |
采用多轴延性耗竭蠕变本构模型[18]如(9)~(14)式计算蠕变应变率。
|
(9) |
|
(10) |
|
(11) |
|
(12) |
|
(13) |
|
(14) |
式中: A, n, m为材料参数; t为蠕变时间; 为蠕变损伤率; 为蠕变应变率。εf*和εf分别为多轴蠕变失效应变和单轴蠕变失效应变。
塑性应变率可以通过屈服函数和一致性条件计算, 考虑损伤和随动硬化的屈服函数可以表示为
|
(15) |
|
(16) |
|
(17) |
|
(18) |
α为背应力, k0为初始屈服应力, Ck和γk为材料常数, p为累计塑性应变, D为损伤变量, 且损伤演化可以表示为
|
(19) |
式中:S, c为材料参数;pD为临界累计塑性应变。
考虑损伤的有效应力可以通过(20)式计算
|
(20) |
式中,De为弹性张量。
2.3 有限元实现 以不统一的模型为例推导其有限元实现方法。采用径向回退映射法更新应力和其他状态变量, 径向回退映射法是一种隐式积分方法[19]。首先假设在增量步[tn, tn+1]内连续体处于弹性状态, 定义此时的应变和应力分别为弹性试应变和试应力。
|
(21) |
|
(22) |
|
(23) |
式中: (·)trial表示变量的试状态; De为材料的弹性矩阵。
以上变量称为连续体的试状态, 如果试状态满足屈服函数即试状态处于弹性区以内, 则试状态为增量步[tn, tn+1]的真实解。如果试状态不满足屈服函数即
|
(24) |
则需要对连续体的试状态进行塑性矫正如(25)~(26)式所示, 径向回退映射法的示意图如图 5所示。
|
(25) |
|
(26) |
在完全隐式的向后欧拉方法中, 在增量步结束时需要更新每一个变量。例如从第n步到第n+1步, 径向回退映射法可以写成如下形式
|
(27) |
|
(28) |
|
(29) |
以上变量的残差可以表示为
|
(30) |
|
(31) |
|
(32) |
使用牛顿-拉夫逊方法求解以上非线性方程组, 第k次的迭代可以写为
|
(33) |
|
(34) |
|
(35) |
式中,。
通过(33)~(35)式可以得到各变量的增量
|
(36) |
|
(37) |
|
(38) |
式中
|
(39) |
|
(40) |
|
(41) |
|
(42) |
|
(43) |
|
(44) |
|
(45) |
为了获得全局迭代过程中的二阶收敛速率, 需要计算一致性切线模量。对(27)~(29)式进行微分得到
|
(46) |
|
(47) |
|
(48) |
将(46)~(48)式重新整理可得
|
(49) |
|
(50) |
|
(51) |
|
(52) |
式中
|
(53) |
|
(54) |
|
(55) |
|
(56) |
|
(57) |
|
(58) |
|
(59) |
将通过以上过程计算得到的应力用于计算蠕变损伤和蠕变应变, 计算流程图如图 6所示。
2.4 参数确定方法 本文采用遗传算法[20]计算本构模型的材料参数。遗传算法的主要步骤分为: 编码、随机产生初始种群、解码、计算适应度函数、选择、交叉、变异、收敛性判别和输出最优解, 其计算流程图如图 7所示。
本文适应度函数的确定方法为将2.1节提出的本构方程化简为一维形式, 分别对应不同温度下试样的单轴拉伸和蠕变试验, 得到一组微分方程, 将解码得到的材料参数初始值代入该方程并使用四阶龙格库塔方法求解, 得到应力的理论值。目标函数和适应度函数选择为
|
(60) |
|
(61) |
式中:σexp为实测应力值;σth为理论应力值。适应度函数f的最大值对应着本构方程参数的最优解。
2.5 计算结果 本文提出的黏塑性本构模型通过有限元软件Abaqus的UMAT子程序进行实现, 并使用三维八节点等参单元进行计算。统一黏塑性模型的材料参数如表 3~4所示, 不统一黏塑性模型的材料参数如表 5~6所示。TC4试样拉伸和蠕变曲线的计算结果如图 8所示。由计算结果可知, 统一的模型可以较为精确地预测TC4合金的拉伸行为, 但却不能准确模拟其蠕变行为, 计算结果的最大误差为43.7%。这是由于所使用的Lemaitre损伤模型不能同时准确预测拉伸与蠕变过程中的损伤。不统一模型计算所得的拉伸应力应变曲线在25℃和600℃下与试验结果吻合良好, 450℃下模拟曲线的屈服点比试验曲线的屈服点稍高。对于计算所得的蠕变曲线第一和第二阶段的蠕变速率与试验结果相吻合, 但第三阶段的蠕变速率要高于试验结果, 计算结果的最大误差为7.7%。与统一的黏塑性本构模型相比, 不统一的黏塑性本构模型虽然表达式比较复杂, 将非弹性应变分解为蠕变应变和塑性应变两部分, 但是可以针对蠕变应变和塑性应变采用不同的计算模型, 因此相比统一的黏塑性本构模型具有更高的计算精度。总的来说, 本文提出的不统一的黏塑性本构模型可以有效预测TC4合金蠕变曲线的3个阶段以及高温拉伸过程中的软化阶段。
表 3 TC4试样统一模型材料参数
温度/℃ |
C1 |
γ1 |
C2 |
γ2 |
n |
K |
k0/MPa |
25 |
1 023.731 |
2.63 |
387.22 |
26.13 |
4.00 |
74.32 |
846.42 |
450 |
3 968.71 |
22.34 |
2 533.72 |
6.49 |
1.09 |
508.11 |
479.17 |
600 |
11 215.31 |
170.29 |
7 696.48 |
424.21 |
8.00 |
43.21 |
92.20 |
表 4 TC4试样统一模型材料参数
温度/℃ |
β1 |
u1 |
β2 |
u2 |
S |
c |
25 |
7.24×10-4 |
1.39 |
3.49×10-4 |
1.20 |
9.60 |
3.0 |
450 |
5.17×10-4 |
2.53 |
8.064×10-4 |
2.69 |
1.98 |
3.3 |
600 |
5.51×10-4 |
1.52 |
1.18×10-4 |
1.29 |
8.17 |
1.0 |
表 5 TC4试样延性耗竭模型材料参数
E/MPa |
μ |
A/(MPa-n·h-1) |
n |
m |
εf/% |
72 481.88 |
0.3 |
1.37×10-19 |
6.25 |
-0.48 |
23 |
表 6 TC4试样不统一模型塑性参数
温度/℃ |
C1 |
γ1 |
C2 |
γ2 |
k0/MPa |
S |
c |
25 |
2 350.92 |
18.94 |
777.13 |
20.38 |
846.42 |
5.81 |
0.63 |
450 |
22 111.44 |
482.89 |
2 326.49 |
17.79 |
513.01 |
9.03 |
1.00 |
600 |
12 522.31 |
248.16 |
2 233.62 |
73.31 |
92.20 |
1.00 |
0.82 |
3 结论 1) TC4试样在450℃和500 MPa的下蠕变断裂时间是600 MPa下的4.12倍。与25℃下的拉伸性能相比, 450℃和600℃下TC4试样的弹性模量分别减低31.02%和59.43%;屈服强度分别降低39.40%和89.11%;450℃下TC4试样的断裂应变增加12.86%, 600℃下TC4试样的断裂应变远高于25℃和450℃下的断裂应变。
2) 不统一的黏塑性本构模型虽然表达式比较复杂, 将非弹性应变分解为蠕变应变和塑性应变两部分, 但是可以针对蠕变应变和塑性应变采用不同的计算模型, 因此相比统一的黏塑性本构模型具有更高的计算精度。
3) 本文提出的不统一的考虑损伤的黏塑性本构模型以及推导的数值实现方法可以有效预测TC4合金蠕变曲线的3个阶段以及高温拉伸过程中的软化阶段。
The viscoplastic constitutive model of TC4 alloy under high temperature