基于修正GTN 模型的不锈钢管剪切过程韧性断裂准则研究

2021-03-22董建鹏王时龙

工程力学 2021年3期

董建鹏，王时龙，周杰，杨波，马驰

(重庆大学机械传动国家重点实验室，重庆 400044)

乏燃料组件是核反应堆卸载出来的不能再维持临界反应的组件，目前世界范围内核发达国家多采用“闭式燃料循环”处理乏燃料组件[1]，该方法需要将乏燃料组件剪切成小段，以方便后续化学提取残留的核元素[2]。不锈钢管是乏燃料组件的重要组成部分，由于其自身较高的延展性，剪切过程中易造成刀具崩刃甚至局部断裂[3]，因此研究不锈钢管的剪切机理对提高刀具的使用寿命有重要作用。然而，近年来对管材旋压[4]、弯曲[5]、矫直[6]、冲压[7]等加工方式的研究较多，但对剪切过程的探索仍比较少，尤其是涉及不锈钢管剪切领域的工作。目前现有文献多集中在管材断裂前的弹塑性变形阶段[8 − 9]，缺乏对后续韧性断裂过程的研究。此外，在模拟变形过程中多采用连续损伤力学方法[10]，力学模型上未考虑实际材料的不均匀性，由此无法准确反映材料的变形破坏情况。

作者在之前的工作中，采用修正Gurson-Tvergaard-Needleman 模型(GTN-J 模型)研究了SUS304 不锈钢管剪切的断裂机理，但是该模型在最大剪切力的预测上仍存在较大的误差，其原因在于随着刀具的移动，钢管的强度随着损伤的累积不断减小，因此钢管不同位置的临界损伤阈值应该是不同的，而GTN-J 模型中的断裂准则被设置为恒定的常数(即当损伤超过该临界值，单元断裂)。鉴于此，本文在该模型基础上，同时考虑不锈钢管中空的特殊结构，对GTN-J 模型的断裂准则进行改进，从而建立适宜于不锈钢管剪切的损伤模型。通过VUMAT 编制该模型应力-应变算法并将其嵌入到ABAQUS/Explicit 中实现其数值求解，随后对比仿真及试验得到的剪切载荷曲线，分析了不同阶段钢管的裂纹扩展情况，并对比了改进前后断裂准则的预测效果。

1 本构模型及数值模拟实现

1.1 修正GTN 本构模型

GTN 模型被广泛应用到研究金属材料韧性断裂过程中[11 − 13]，其具体形式如下：

式中： f 代表当前材料的孔洞体积分数； fc及 fF分别代表材料孔洞开始连通及完全断裂时的孔洞体积分数。

GTN 模型仅适应于高应力三轴度状态下的断裂预测，对于如剪切等低应力三轴度状态，原始GTN 模型无法模拟材料变形后的损伤累积。鉴于此，Wei 等[15]在现有剪切修正GTN 模型基础上引入新的损伤变量表征低应力三轴度下的损伤累积，修正后屈服方程为(以下简称GTN-J 模型)：

式中： Ds为材料的剪切损伤； Dc是材料剪切损伤加速累积的临界参数； DF表示材料完全断裂失效时的剪切损伤。

GTN-J 模型中采用两种独立的损伤变量表征材料复杂应力状态下的断裂机理：1)孔洞体积分数；2)剪切损伤。

1.1.1 孔洞体积分数演化

对于金属材料，仅在拉应力下出现孔洞萌生，由此孔洞萌生系数表示为：

式中： fN为孔洞形核系数； εN为孔洞萌生的平均等效塑性应变； SN为孔洞萌生标准差。

1.1.2 剪切损伤演化

GTN-J 模型中剪切损伤累积规律形式如下：

式中： kw为材料剪切常数； wσ为应力状态参数，用来区分材料受载过程中的不同应力状态，保证剪切损伤在复杂应力状态下能够正常累积，其函数形式如式(11)所示。

1.2 断裂准则修正

在GTN-J 模型中，当单元的孔洞体积分数超过临界值 fc或剪切损伤超过临界值 Dc时，判定临界值 fc或剪切损伤超过临界值 Dc时，判定该单元失效(即发生断裂)。整个加载过程中，两种损伤的临界值保持不变。

然而，对于不锈钢管来说，钢管不同部位的临界值是不同的。如图1 所示为刀具剪切钢管上/下半部分时剪切力 Fs的分解。 Fs可分解为切向剪切力 Fst和径向剪切力 Fsr，其中切向剪切力 Fst主导剪切过程中管材裂纹的扩展，当 Fst超过材料的许用应力极限时，裂纹沿管材切向扩展演化。

图1 钢管剪切过程剪切力分解Fig.1 Schematic diagram of shear force decomposition for tube shearing

剪切过程中刀刃与钢管从点接触发展到线接触，接触线上剪切力分布并不均匀。如图2 所示，沿着厚度方向将剪切区域分成n 个等长度的微元，各微元受力分别为Ft1,Ft2,···,Ftn，剪切力的切向分力可表示为：Fst=Ft1+Ft2+···Ftn。随着剪切刀向下移动，与刀具接触的部分会出现裂纹，一旦沿厚度方向某一处微元达到材料承载极限出现裂纹，裂纹会很快扩展至整个厚度方向，因此可以近似认为同一厚度方向的微元会同时达到材料承载极限，从而同时被撕裂，考虑到同种材料的承载极限是相同的，假设管材厚度方向单位长度的许用应力极限为 σs，则各微元在断裂时受力可表示为：Fti=tσs/n。由此，剪切力的切向分力可表示为：Fst=Ft1+Ft2···Ftn=n×t×σs/n=σst。

图2 剪切过程接触线上剪切力的不均匀分布Fig.2 Irregular distribution of shear force on contact line during shearing process

本文所使用管材厚度为1 mm，由此管材剪切力可由式(12)近似计算得到(剪切时剪切管材对称的两侧，因此需在公式中乘2)。此外随着损伤的累积，管材的许用应力极限随着刀具的位移不断减小，试验发现，当剪切管材下半部时(即x>R)，许用应力极限有明显地降低(在剪切下半部分时试验剪切力明显降低)。因此剪切过程剪切力估计值修正为式(13)所示。需要注意的是剪切力与临界损伤呈正相关，即剪切力越大，对应部位的临界损伤也应该越大。考虑到管材剪切过程中主要由剪切损伤主导，研究中仅对剪切损伤进行修正，根据剪切力的变化构造了一个阈函数式(14)来表征管材不同部位的剪切损伤临界值。图3 所示为剪切力与剪切损伤的变化趋势。

图3 剪切力及剪切损伤变化趋势Fig.3 Variation trends of shear force and shear damage

式中： β为剪切力 Fs与切向分力的夹角； R为管材外半径； x是刀具的位移(刀具刚接触管材上半部，x=0 )； σs1和 σs2分别为管材上/下半部分的许用应力极限；与此对应的Dc1和Dc2分别代表管材上/下半部分的临界剪切损伤值； kf为定义的损伤系数； a1表示刀具位移的修正系数用来消除由于管材本身加工、制造、安装等造成的误差。

1.3 数值模拟实现算法

由于ABAQUS 本身中未嵌入GTN-J 本构模型，本文采用完全隐式向Euler 应力更新算法Achouri 等[16 − 17])，通过编制GTN-J 模型的VUMAT用户子程序，最终实现该模型的数值运算。如图4所示为其数值算法流程，整个调用过程划分为三个阶段：弹性阶段、塑性阶段及断裂阶段。

2 剪切试验及有限元模拟

2.1 试验方案

图5(a)所示为后处理中需要剪切的组件芯管，其由不锈钢管及脆性芯块构成。剪切过程中，脆性芯块在前期工艺中易受压破裂，因此本文对其简化处理，剪切试验中采用在不锈钢管两端填充支撑块的方案。

图4 数值算法流程图Fig.4 Flow chart of numerical algorithm

剪切试验采用SUS304 不锈钢管，剪切装置包含剪切刀、压紧块、固定刀及支撑块四部分。试验过程中剪切刀及固定刀被夹紧在在材料疲劳试验机(Landmark 370.1)上下夹块处，剪切装置安装如图5(b)所示，剪切装置的重要结构尺寸如图5(c)及表1 所示。

剪切试验剪切两种规格的不锈钢管；1)外径16 mm，厚度1 mm；2)外径12 mm，厚度1.5 mm。剪切第二种规格时需要在钢管外增加一套筒，如图6 所示。剪切前调整刀具位置使其离不锈钢管上表面3 mm～5 mm，选取50 mm/s 及5 mm/s 两种剪切速度进行剪切，并在剪切过程中采集剪切刀具载荷-位移数据。

2.2 有限元模型

为了提高有限元计算效率，对剪切模型进行简化，建立如图5(d)所示模型。

考虑到剪切装置的刚度远大于不锈钢管的刚度，因此将其设置成刚体，且仅保留靠近剪切区域的结构。剪切装置各零件网格类型及大小设置如表2 所示。

图5 剪切试验装置及有限元模型Fig.5 Shear device and finite element model

表1 剪切设备结构尺寸Table 1 Structural dimensions of shearing device

图6 外径12 mm 的不锈钢管剪切图示Fig.6 Diagram of shearing stainless steel tube with outer diameter of 12 mm

对于不锈钢管，其左侧区域由于支撑块I 的存在剪切过程几乎不发生变形，由此该区域设置成刚性单元；右侧区域仅保留靠近剪切区域的部分，采用Mises 塑性本构(仅产生塑性变形)；中间区域应用本文改进后的GTN-J 本构模型。不锈钢管不同区域本构模型设置如图7 所示。

图7 SUS304 不锈钢管不同区域本构设置Fig.7 Constitutive models for different regions of SUS304 tube

表2 剪切装置网格划分Table 2 Mesh generation of shear device

2.3 本构模型参数确定

关于GTN 模型参数的确定已有较多研究[18 − 19]，而对于修正GTN 模型，其共需提前确定22 个参数的值，为解决数据唯一性问题，可将其分成5 类：

1) 5 个表征材料性能的参数( E , ν , σ0, K , n)

该类参数可以通过SUS304 标准拉伸试验确定[20]，具体如下：弹性模量 E=195.85 GPa，泊松比 ν=0.29，初始屈服强度 σ0=390.3 MPa，材料应变强化常数 K=1235 和 n=0.38；

2) 2 个GTN 模型参数( q1, q2)

针对金属材料，Tvergaard 建议 q1和 q2的值分别为1.5 和1[2]；

3) 6 个与孔洞体积分数相关的参数( f0, fc, fF, εN,SN, fN)

Zhang 和Cong[20]采用GTN 模型模拟SUS304板材液压胀形时得到了SUS304 不锈钢材料孔洞体积分数相关参数的值。

4) 6 个与剪切损伤相关的参数( D0, Dc, DF, kw,c1, k)

剪切损伤参数最早由Wei 等[15]引入到GTN模型中，不同于孔洞体积分数，该参数并没有实际的物理意义，无法通过试验得到，本文首先通过单个单元测试分析了不同参数对剪切效果的影响程度结合Wei 等[15]现有的研究确定剪切损伤相关参数的值。

如图8 所示，建立一立方单元，对其施加纯剪切边界，对比不同参数下对应力-应变演化过程的影响。由图8、图9 可明显看出 DF, c1, k的变化对应力-应变演化基本无影响，而 D0, Dc, kw相互影响，对同一种材料来说，仅需确定 D0, kw则 Dc参数的值即可唯一确定。

考虑到金属材料相近的特性，对于SUS304 不锈钢材料，采用Jiang 的推荐值，即 D0=0.01,DF=0.15, c1=0.5, k=0.2, kw=4，由此本文提出的Dc1及Dc2的值可通过有限元模拟方法唯一确定，即先设置2 个参数的初始值进行有限元模拟，然后纯剪切条件下DF/c1/k对应力-应变演化的影响调节各参数的值直到模拟的载荷-位移曲线与剪切试验数据差值最小。

5)与断裂准则相关的参数( kf, a1)

同理，2 个与断裂准则相关的参数也可通过有限元模拟方法得到，根据式(14) kf和 a1两个参数与临界损伤参数Dc1及Dc2相关联，当Dc1及Dc2确定， kf和 a1即可唯一确定。

图8 纯剪切条件下DF,c1,k对应力-应变演化的影响Fig.8 Influence of DF,c1,k on stress-strain evolution under pure shear condition

图9 纯剪切条件下Dc,kw对应力-应变演化的影响Fig.9 Influence of Dc,kw on stress-strain evolution under pure shear condition

表3 罗列了有限元最终采用的修正GTN-J 的模型参数。

表3 修正GTN-J 模型参数Table 3 Modified GTN-J model parameters

3 结果讨论

3.1 原始断裂准则与修正后预测结果对比

剪切过程中不锈钢管剪切力-位移模拟曲线与试验结果对比如图10 所示。整体上来看，根据裂纹沿管材周向的演化路径，可以将整个剪切过程分为5 个阶段：1)开始剪切阶段，此时刀具尚未接触不锈钢管；2)初始裂纹萌生阶段，不锈钢管与刀具接触部分被撕裂，产生最大剪切力；3)第I 阶段撕裂过程，裂纹沿(a−b1/a−b2)在不锈钢管周向向外扩展；4)第II 阶段撕裂过程，裂纹沿(b1−c1/b2−c2)在不锈钢管周向向内扩展；5)断裂，由于支撑块的存在，不锈钢管残留部分一同被撕裂分离。

图10 模拟的载荷-位移曲线与试验数据对比(D=16 mm, t=1 mm, V=50 mm/s)Fig.10 Comparison between simulated load-displacement curves and experimental results (D=16 mm,t=1 mm, V=50 mm/s)

为了对比两种断裂准则的优劣性，考虑到修正后的断裂准则对管材上/下半部分采用不用的临界损伤值，因此使用原始断裂准则模拟两次，分别设置剪切损伤临界值为Dc=0.038 和Dc=0.11。通过对比明显可以得到本文提出的断裂准则的模拟结果与试验结果误差更小，预测效果更好(表4)：相较于原始断裂准则，断裂时的位移预测误差由28.8%降低至4.9%；最大剪切力预测误差由27.5%降低至17.9%。

最大剪切力的预测误差相较试验值偏小的原因如下：不锈钢管由不锈钢薄板经多次辊轧弯曲成形，经过多次加工硬化，其力学性能相较初始的薄板有很大的提高。而模拟时不锈钢管的应力应变曲线由薄板拉伸试验得到(由加工不锈钢管的钢板截取)，未考虑其成形过程的加工硬化，因此模拟值相较试验值偏小。

表4 两种断裂准则模拟结果对比Table 4 Comparison of simulated results between original and modified fracture criteria

3.2 修正模型的广泛适用性

为了进一步验证模型的广泛适用性，本文增加了一组管材直径12 mm，厚度1.5 mm，剪切速度50 mm/s 的剪切试验。

由图11 的对比结果可以看出：从第I 阶段到第IV 阶段，修正模型的模拟结果与试验结果基本一致，仅第V 阶段模拟曲线与试验曲线有差异。但该差异是由于工装中支撑块与钢管配合不紧密造成剪切过程中支撑块后移所致。如图12 所示，在钢管即将断裂时，由于支撑块与钢管配合并不紧密，支撑块后移，钢管无法直接断裂，随着剪切刀下降，被剪切部分向下变形翻转，但此时裂纹并不扩展，直到刀具接触变形后的残留部分将钢管最终撕断，因此试验中剪切曲线在即将断裂时剪切力出现一个平稳台阶，随后小幅度上升后降为0。

图11 模拟的载荷-位移曲线与试验数据对比(D=12 mm, t=1.5 mm, V=50 mm/s)Fig.11 Comparison between simulated load-displacement curve and experimental results(D=12 mm, t=1.5 mm, V=50 mm/s)

图12 不锈钢管(外径12 mm 厚度1.5 mm)剪切过程支撑块后移现象Fig.12 Back movement of support block when shearing stainless-steel tube of D=12 mm and t=1.5 mm

若支撑块不后移，剪切力曲线会持续下降直至为0(同剪切外径16 mm 厚度1 mm 钢管现象一致)，试验曲线变化趋势如图11 中虚线所示，模拟曲线与试验曲线仍基本一致。综上，可以认定整个剪切过程中本文提出的模型的模拟结果与试验结果基本一致，提出的模型是具有广泛适用性的。

3.3 剪切速度的影响分析

如图13 所示为5 mm/s 速度下的剪切试验与修正模型模拟结果的对比(D=16 mm/t=1 mm)，可以得到修正模型可应用在不同速度工况下的模拟中，速度的变化对分析结果并无明显影响。

图13 模拟的载荷-位移曲线与试验数据对比(D=16 mm, t=1 mm, V=5 mm/s)Fig.13 Comparison between simulated load-displacement curves and experimental results (D=16 mm,t=1 mm, V=5 mm/s)