基于ABAQUS二次开发的轧制过程组织性能预测

付晓杰，杜晓钟，黄庆学

(太原科技大学 重型机械教育部工程研究中心，太原 0300240)

摘要：考虑到实际板带轧制过程中组织性能预测的精度缺陷，基于低合金钢轧制过程动态再结晶型的真实应力-应变数学模型及微观组织模型，通过VUMAT及USDFLD数据接口，依托于ABAQUS有限元软件，对软件材料子程序及其微观组织模型进行了二次开发研究，最后对Q345B板带的热轧制过程进行了模拟验证。研究结果表明：模型的预测精度较高，板带横截面处各部位的静态软化率和奥氏体晶粒直径差别较大，板带中心处的奥氏体晶粒分布比较均匀，板带侧面及表面奥氏体晶粒细化程度不大，要靠后续的精轧过程继续细化。

关键词：材料本构模型；ABAQUS二次开发；组织性能预测；轧制过程模拟；动态再结晶

收稿日期：2015-02-22

基金项目：973计划前期研究专项(2012CB722801；山西省自然科学基金(20121060)

作者简介：付晓杰(1989-)，女，硕士研究生，研究方向为高精度板带轧制技术。

中图分类号：TG335.5文献标志码：A

在轧制过程中，金属变形抗力值的大小是计算轧制过程的力和功，制定合理的轧制工艺规程，设计和校核轧制设备的重要依据。衡量轧制产品质量的标准既取决于产品最终形状和尺寸的外部质量又取决于产品性能的内部质量。从宏观方面进行控制就要制定合理的轧制规程，这就要求对力能参数进行精确的计算。准确的本构关系模型是轧制过程中材料的动态响应与热力参数的媒介，是用数值分析方法对轧制过程进行数值模拟的前提条件。从微观方面进行控制就要控制内部微观组织，在化学成分确定的情况下，微观组织晶粒度是衡量内部微观组织的重要标准。在钢材料的本构关系研究中，学术界提出了不同类型的钢本构关系模型。Maheshwari等[1]提出的修正Johnson-Cook模型，此模型是基于隐式算法需要进行平衡迭代，计算慢，容易出现不收敛的问题。王培起等[2]在研究大型H型钢热轧工艺过程有限元仿真系统开发中，对修正的Johnson-Cook模型进行了ABAQUS的二次开发，但计算量大，耗时长。基于神经网络和模糊理论等的模型[3]，这些模型较为复杂，耦合到有限元软件中有着较大的困难。管克智、周纪华[4]建立不同钢种的变形抗力模型，大多数属于动态回复型，模型系数较多，计算起来比较繁琐。许勇顺、柳建韬等[5]建立的动态再结晶型的应力应变曲线数学模型，公式简单易懂，适合工程实际应用，此模型得到广泛应用。陈庆军等[6]在商业有限元软件MSC.Mrac中对文献[5]中的动态再结晶模型进行二次开发。马博等[7]以Visual C和OpenGL编程语言为开发平台，选用文献[5]中的动态再结晶模型进行了基于条元法多参数耦合模拟仿真系统的开发，实现了热力微观组织耦合。

首先引入这种适合于工程应用的动态再结晶型的应力-应变曲线数学模型，利用ABAQUS中的VUMAT的子程序接口，用 FORTRAN 语言编译代码，进行了材料本构模型的二次开发。然后利用ABAQUS的USDFLD的子程序接口，完成微观组织模型的二次开发。最后利用开发的子程序对Q345B材料的板带轧制过程进行模拟，实现热力微观组织耦合，提取模拟结果与实际对比，表明二次开发子程序的有效性，完成组织性能预测，可为技术研究及现场实践提供技术支持。

1本构关系数学模型

ABAQUS 为了方便用户开发自己研究的本构模型，允许用户应用FORTRAN 语言编写用户子程序代码来扩展主程序功能，VUMAT是ABAQUS提供给用户定义特殊材料属性的子程序接口，用户可以根据需要通过二次开发添加所需的本构模型[8-9]。

本文应用弹塑性应力—应变增量理论，采用混合强化模型，利用切向预测径向返回的方法[10-11]更新应力，实现VUMAT用户子程序开发。

1.1切向预测径向返回算法

在弹塑性有限元计算中，采用切向预测径向返回算法更新应力，取代隐式积分迭代计算弹塑性张量，提高计算效率。当材料屈服时，应变增量分为弹性应变和塑性应变，如图(1)所示。

图1　显式差分算法

在试探应力的基础上减去应用广义胡克定律多算的塑性应变对应的应力增量，最终得到弹塑性应变对应的应力增量，即:

(1)

1.2混合强化模型

混合强化模型同时考虑各项同性强化和随动强化两种法则，适应材料一般硬化特性要求。

由混合硬化的屈服条件推导出VonMises屈服准则，屈服函数可表示为：

(2)

采用弹塑性增量法分析变形过程中的应力变化，计算时首先假定应力符合线弹性公式，由Hook定律计算试探应力为：

(3)

式中：λ和μ为材料的Lame常数，

有效塑性应变张量形式[12-13]为:

(4)

式中：εp为有效塑性应变，εp为塑性应变，dγ是贯穿于整个塑性加载历史的非负标量函数。

采用混合硬化法则[14]，背应力写成张量的形式为:

(5)

式中：Ep为单轴应力-应变曲线斜率，Q为当前应力点法向方向的张量形式。

对于混合硬化模型，Mises屈服面的法线定义为:

(6)

若试探应力值未超过屈服应力值，真实应力就是试探应力；若试探应力值超过屈服应力值，材料发生塑性变形，此时的应力状态、背应力状态、等效塑性应变根据式(3)-式(6)更新为：

σnew=σold+dα

(7)

(8)

采用VUMAT的特定的子程序接口，实现本构模型的步骤如图2所示：

(1)从ABAQUS主程序中载入材料参数、应变增量和应力张量等已知量，作为VUMAT子程序的初始化中的变量值。

图2　VUMAT子程序流程图

(2)由ABAQUS主程序传入的已知量，然后通过式(3)计算试探应力。

(3)计算应力差张量ξ:

(9)

(4)计算式(9)，若f>0，则通过式(11)计算有效塑性应变增量，在子程序中用自定义符号deqps表示；若f>0，则deqps为0.

(5)计算塑性应变，在ABAQUS的自定义用户状态变量中，用SDV7来存储塑性应变值，stateold(i,7)表示可以从SDV7中读出的应变值。

(10)

(6)通过式(7)进行背应力更新和应力更新。

2奥氏体再结晶子模型

轧件内部微观组织晶粒度大小是影响产品最终力学性能和使用性能的重要因素之一。通过Fortran语言编写奥氏体再结晶函数，利用call函数在VUMAT子程序中调用，图(3)和图(4)示出了奥氏体再结晶函数，并采用式(11)-式(12)通过编写USDFLD用户子程序，完成奥氏体晶粒长大模型的二次开发，实现对组织性能的预测。

图3　小于动态再结晶临界应变

图4　大于动态再结晶临界应变

(11)

(12)

3Q345B板带轧制模拟

本文模拟了钢种Q345B在其实测热物性参数[15]条件下的轧制过程。利用ABAQUS的前处理进行建模，精轧模型如图(5)，在User Material中输入模型参数，本构模型参数见表1.在EDIT JOB菜单中的GENERAL子菜单的USER SUBROUTINE FILE对话框中选择用户子程序.for文件即可。

某钢厂的原料规格尺寸，原料尺寸选取为180×1050×2380 mm，轧制产品尺寸为6×2800×2400 mm，经计算粗轧7道次，精轧8道次，其轧制规程如表2.

表1　算例的本构模型参数

表2　寸为的压下规程

图5　精轧模型

4轧制过程分析结果

4.1轧制力与温度场分析

通过对Q345B钢进行上述规程的板带轧制模拟，图(6)给出了Q345B钢轧制道次轧件表面温度模拟结果与实测值的对比分析。

分析表明模拟的粗轧道次和精轧道次轧件表面温度与实测的轧件表面温度温差在左右，模拟的精度较高。图(7)给出了Q345B钢轧制力模拟结果与实测值的对比分析，分析表明有限元模拟的计算结果和实测的结果很接近，其最大偏差为7.71%。通过以上分析，验证了VUMAT子程序的可行性和准确性，模拟结果具有比较理想的精度。

4.2组织性能预测

首先在新开发的本构模型基础上，完成了15个

图6　厚6 mm Q345B 的宽薄板轧制过程

道次的板带轧制模拟，然后再利用子程序USDFLD开发的微观组织模型， 完成精轧一个道次的板带轧

图7　Q345B厚6 mm宽薄板轧制过程获得的

制模拟，图(8)-(13)示出了轧制间隙的横截面静态软化率和奥氏体晶粒尺寸变化情况。通过结果分析，由于精轧过程温度较低，只有板带中心发生了少量的动态再结晶。

图8　横截面的静态软化率云图

图9　横截面宽度方向的静态软化率

从静态软化率的分布可以看出，由于轧制过程中，板带表面和侧面温度降低的较快，在板带的宽度方向从中心到两侧面、厚度方向从中心到表面，静态软化率逐渐减小，两侧面静态软化率减小程度较大。

图10　横截面厚度方向的静态软化率

图11　横截面的平均晶粒尺寸云图

图12　横截面宽度方向的平均晶粒尺寸

从晶粒尺寸分布情况看，由于静态软化程度和温度的影响，静态软化程度越小，温度越低，初始晶粒尺寸占主要部分，在板带的宽度方向从中心到侧面、厚度方向从中心到表面平均晶粒尺寸增大。

图13　横截面厚度方向的平均晶粒尺寸

5结论

(1)新建了一种适用于动态再结晶型低合金钢的ABAQUS热轧制工艺模拟的变形抗力及微观组织模型，经分析验证模型的预测精度较高。

(2)通过对板带轧制过程的数值分析可得板带横截面处各部位的静态软化率和奥氏体晶粒直径差别较大，板带中心处的奥氏体晶粒平均直径为65 μm，分布比较均匀；板带侧面及表面奥氏体晶粒平均直径在90~120 μm，晶粒细化程度不大，要靠后续的精轧过程继续细化。

参考文献：

[1]MAHWAHWAR A K，PATHAK K K，RAMAKRISHANA N，et al.Modified Johson-Cook material flow model for hot defomation processing[J].Journal of Material Science，2010，45(4)：859-864.

[2]王培起.大型H型钢热轧工艺过程有限元仿真系统开发[D].山东济南：山东大学，2013.

[3]NIU Y，LUO J，LI M Q.An adaptive constitutive model in the isothemal compression of Ti600 alloy [J].Materails Scieence andEngineering，2010，527(75)：5924-5929.

[4]赵志业.金属塑性变形与轧制理论[M].2版.北京：冶金工业出版社，2012.

[5]许勇顺，柳建韬，聂明，等.金属热变形应力-应变曲线数学模型的研究与应用[J].应用科学学报，1997，15(4)：379-384.

[6]陈庆军.高强度宽薄板轧制过程有限元模拟及再结晶行为研究[D].北京：北京科技大学，2006.

[7]马博.低合金钢板带轧制过程微宏观多参数耦合建模[D].河北秦皇岛：燕山大学，2011.

[8]陈建康，黄筑平.单向应力条件下松弛时间率相关的非线性粘弹性本构模型[J].高分子学报，2003，3(4)：414-419.

[9]范庆来，栾茂田，杨庆.修正剑桥模型的隐式积分算法在ABAQUS 中的数值实施[J].岩土力学，2008，29(1)：269-273.

[10]王勖成.有限单元法.[M].北京：清华大学出版社，2003.

[11]LOU J，LI M，LI X，et al.Constitutive model for high temperature deformation of ti-tanium alloys using internal state variables[J].Mechanics of Materials，2010，(42)：157-165.

[12]马晓峰.ABAQUS 6.11有限元分析从入门到精通[M].清华大学出版社，2013.

[13]陈惠发，A.F.萨里普.弹性与塑性力学[M].北京：北京建筑工业出版社，2004.

[14]HAMADA M，TANAKA M.A Numerical Method Considering the B-auschinger Effect for Large Deflection Analysi-s of Elastic-Plastic Circular Plates.[J].Bulletin of JSME，2008，103(15)：1029-1040.

[15]曾正明.机械工程材料手册：金属材料[M].北京：北京机械工业出版社，2003.

Performance Prediction of Rolling Process based on

ABAQUS Secondary Development

FU Xiao-jie，DU Xiao-zhong，HUANG Qing-xue

(University of Science and TechnologyTaiyuan，The heavy mechanical engineering research center

of Ministry of Education ，Taiyuan 030024，China)

Abstract：Considering precision defects of organization performance prediction in actual strip rolling process，based on real stress-strain mathematical model of the dynamic recrystallization model and microstructure model in rolling process of low alloy steel，and by VUMAT and USDFLD data interface，based on the finite element software ABAQUS，

there is a secondary development research for material subroutine and its microstructure model.Finally，it is verification of simulation by Q345B strip hot rolling process.The results show that the prediction precision is higher，and static softening rate and the austenitic grain diameter difference is bigger from each part of the plate with cross section.In strip center the austenitic grain size distribution is more uniform ，and degree of grain refinement for austenitic grain is not big from strip profile and surface，so it is necessary to rely on the follow-up finishing process.

Key words：material constitutive model，secondary development of ABAQUS，organization performance prediction，rolling process simulation，dynamic recrystallization