上海交通大学学报, 2025, 59(8): 1181-1191 doi: 10.16183/j.cnki.jsjtu.2024.260

机械与动力工程

淬火分配钢成形极限实验与晶体塑性模拟

杨浩1, 唐伟琴,2

1.合肥通用机械研究院有限公司, 合肥 230061

2.上海交通大学 机械与动力工程学院, 上海 200240

Experiment and Crystal Plasticity Simulation of Forming Limit of Quenching Partitioning Steel

YANG Hao1, TANG Weiqin,2

1. Hefei General Machinery Research Institute Co., Ltd., Hefei 230061, China

2. School of Mechanical Engineering, Shanghai Jiao Tong University, Shanghai 200240, China

通讯作者: 唐伟琴,助理研究员;E-mail:weiqint@sjtu.edu.cn.

责任编辑: 王历历

收稿日期: 2024-07-2   修回日期: 2024-08-19   接受日期: 2024-09-4  

基金资助: 合肥通用机械研究院有限公司博士科技基金项目(2023010792)

Received: 2024-07-2   Revised: 2024-08-19   Accepted: 2024-09-4  

作者简介 About authors

杨浩(1992—),工程师,从事金属材料塑性力学和氢脆研究.

摘要

淬火分配(QP)钢因塑性变形过程会诱发马氏体相变,因而具有高强度和高延伸率的特点,但其成形性能目前尚未明确.通过Nakajima实验获得了QP1180钢板在不同应变路径下的极限应变,采用考虑相变的晶体塑性有限元模型耦合Marciniak-Kuczynski理论(CPFEM-PT-MK),分析了织构演化和相变对QP1180钢板成形极限的影响.结果表明:当应变路径ζ=0.1时,QP1180钢板的极限主应变最低,CPFEM-PT-MK模型能够较好地预测QP1180钢板的成形极限;不同应变路径下,QP1180钢各相的织构演化存在明显差异,考虑织构演化时,其成形极限更高;当不发生相变时,极限主应变最低点在应变路径ζ=0位置,与发生相变时明显不同.此外,相变并不总是提高QP1180的成形极限,其效果与应变路径有关.

关键词: 晶体塑性有限元; 马氏体相变; Marciniak-Kuczynski理论; 织构演化; 成形极限

Abstract

Quenching-partitioning (QP) steel combines ultrahigh strength with good ductility due to the martensitic transformation during plastic deformation. However, the formability of the QP1180 steel remains unclear. In this paper, the ultimate strains of the QP1180 steel under different strain paths are obtained through Nakajima experiment. The effects of the texture evolution and phase transformation on the forming limit of QP1180 steel are analyzed by using a crystal plastic finite element model coupled with the Marciniak-Kuczynski theory (CPFEM-PT-MK). The results show that the ultimate principal strain of QP1180 steel is the lowest under the strain path ζ=0.1, and the established CPFEM-PT-MK model successfully predicts the forming limit of the QP1180 steel sheet. The texture evolutions of the constituent phases in the QP1180 steel are different under various strain paths. According to the simulation, the texture evolutions enhance the forming limit of the QP1180 steel under various strain paths. Without phase transformation, the minimum limited major strain of the QP1180 steel is located at the strain path of ζ=0, which is significantly different from that when phase transformation occurs. Furthermore, the phase transition, related to the specific strain path, does not always enhance the forming limit of the QP1180 steel.

Keywords: crystal plasticity finite element; martensitic transformation; Marciniak-Kuczynski (MK) theory; texture evolution; forming limit

PDF (12960KB) 元数据 多维度评价 相关文章 导出 EndNote| Ris| Bibtex  收藏本文

本文引用格式

杨浩, 唐伟琴. 淬火分配钢成形极限实验与晶体塑性模拟[J]. 上海交通大学学报, 2025, 59(8): 1181-1191 doi:10.16183/j.cnki.jsjtu.2024.260

YANG Hao, TANG Weiqin. Experiment and Crystal Plasticity Simulation of Forming Limit of Quenching Partitioning Steel[J]. Journal of Shanghai Jiaotong University, 2025, 59(8): 1181-1191 doi:10.16183/j.cnki.jsjtu.2024.260

淬火分配(quenching-partitioning, QP)钢因其在保证高强度的同时具有良好塑性,被认为是极具前景的车身零部件用材.QP钢初始包含铁素体、马氏体和残余奥氏体多相微观组织,变形过程需要各相相互协调,且塑性变形过程中残余奥氏体会转变为马氏体.由于其多相微观结构及变形诱发马氏体相变[1-2],变形过程中的力学响应十分复杂,导致成形极限难以预测[3].为了探究QP钢在不同应变路径下的极限应变及影响因素,国内外学者从实验和理论的角度开展了研究.

在实验方面,主要采用平底凸模(Marciniak方法)或球形凸模(Nakajima方法)对不同宽度试样进行刚模胀形实验,以获取不同应变路径加载过程.结合不同应变测量技术,如数字图像相关(digital image correlation,DIC),测量失效位置的主应变和次应变从而得到相应应变路径下的极限应变,最后将所有应变路径的极限应变汇总,得到成形极限曲线(forming limit curve,FLC).例如,Chen等[4]利用Marciniak方法并结合DIC技术研究QP980钢的成形极限,并与同等强度级别的商用DP980钢比较,发现QP980的成形极限明显高于DP980.QP980在接近平面应变路径加载下的成形极限(FLD0)约为0.18,而DP980只有0.12.Nakajima实验[5-7]获得的QP980钢FLD0结果约为0.2,与Marciniak实验结果接近.尽管针对QP980钢成形极限已开展不少实验研究,但强度级别更高的QP钢如QP1180钢的成形极限目前尚不清楚.

由实验方法建立完整的FLD往往耗时且困难,为了能在理论上解释颈缩现象并预测FLD,学者们提出了多种预测方法.这些方法分为3类:第1类以Keeler[8]为代表,通过总结大量实验数据得到计算FLC的经验公式;第2类以Hill[9]、Stören和Rice[10]为代表,以连续介质力学为基础并假设材料没有任何缺陷,通过变形系统失稳的物理现象来理解成形极限,即“分叉”理论;第3类以Marciniak和Kuczyński[11]为代表,通过引入材料缺陷建立模型,预测成形极限,即MK模型.MK模型因形式简单且预测精度高而得到广泛应用.通过将MK模型与材料本构模型结合,可以预测材料在不同应变路径下的成形极限[12-13].例如,余海燕等[14]将MK理论与现象学本构结合,成功拟合了相变诱发塑性(transformation induced plasticity,TRIP)钢的成形极限,并给出相变提高TRIP钢成形极限的结论.然而,现象学本构模型难以反映材料塑性变形过程中的微观本质,计算结果难以具备预测性.

本文以第3代先进高强钢QP1180钢为研究对象,首先通过Nakajima实验获取QP1180钢板的FLD;然后,将相变晶体塑性有限元模型(crystal plasticity finite element method considering phase transformation,CPFEM-PT)与MK理论相结合,建立CPFEM-PT-MK集成计算方法,仅通过拟合平面应变路径下的极限应变,就能准确预测其它应变路径下的极限应变,从而获得完整的FLD;最后,基于CPFEM-PT-MK方法,探讨相变和织构演化对QP1180钢板成形性能的影响,发现相变并非总是提高TRIP钢的成形极限.

1 FLD实验研究

1.1 实验方案

实验材料为宝钢提供的1.4 mm厚QP1180冷轧钢板,采用ISO 12004-2-2008标准中的Nakajima测试方法,如图1所示.凸模直径为100 mm,凹模直径为105 mm,圆角半径R为10 mm,凸模与试样中间涂有润滑层以减小摩擦对实验结果的影响,本文采用的润滑层为0.1 mm厚聚四氟乙烯和少量二硫化钼润滑脂.试样几何尺寸如图1(b)所示,直径d为180 mm,圆角半径R为25 mm,宽度W分别为20、40、60、80、100、120、140、180 mm,宽度方向与轧制方向(rolling direction,RD)一致.利用电火花切割方法将试样从钢板中切下,然后对切割面用180~1 200 目的砂纸依次从粗到细进行打磨,去除毛边,以防实验过程中试样从边缘缺陷处开裂导致实验失败.

图1

图1   Nakajima实验示意图(mm)

Fig.1   Schematic diagram of Nakajima test (mm)


实验前,首先利用电化学腐蚀方法在试样表面印刷尺寸为2.5 mm的网格,以便后续应变分析.实验过程的压边力为1 MN,凸模速度为1.5 mm/min,当成形力下降值达到300 N时实验停止,认为试样已发生破裂.每种宽度尺寸试样重复3次实验,如图2所示.实验完成后,用相机从两个不同的角度拍摄裂纹附近的网格,利用网格应变自动分析系统分析试样颈缩或断口附近位置的应变,图3给出宽度为20 mm试样颈缩附近区域的主应变云图结果.

图2

图2   QP1180成形极限实验试样

Fig.2   Samples of QP1180 of forming limit tests


图3

图3   宽度20 mm试样颈缩附近区域的主应变分布

Fig.3   Distribution of major strain near the necking of a 20 mm-wide specimen


1.2 实验结果

分析不同宽度尺寸试样破裂附近位置的应变,得到QP1180钢板的FLD,如图4所示.不同应变路径下,QP1180的极限主应变有所不同,主应变最低的位置对应为接近应变路径ζ=0靠右侧位置,这是由于非线性加载过程造成的.根据Aretz[15]和Butcher等[16]的研究,Nakajima实验过程存在两个阶段:第1阶段球形凸模与板料的接触面较小,此时试样变形量较小,中心区域的应变路径接近等双向拉伸应变状态;第2阶段球形凸模与板料的接触面较大,此时试样变形量较大,中心区域近似为线性应变加载状态.由于加载路径发生变化,所以导致极限主应变最低点处于平面应变靠右位置.在FLD右侧,随着应变路径ζ=εminor/εmajor的增加,QP1180的极限主应变先降到最低点,然后不断增大,在等双向拉伸应变路径(ζ=1)下,极限主应变约为0.31,εmajorεminor分别为主应变和次应变;在FLD左侧,随着ζ减小,极限主应变不断增大,在单向拉伸应变路径ζ=-0.5 下,极限主应变约为0.28;在平面应变路径ζ=0下,极限主应变约为0.15.

图4

图4   QP1180成形极限实验结果

Fig.4   Experimental FLD results of QP1180


2 CPFEM-PT 模型与MK理论耦合方法

详细的CPFEM-PT本构模型及计算方法在文献[17]中已有详细描述,本文仅给出CPFEM-PT与MK理论的耦合方法.将MK理论与CPFEM-PT本构模型结合,形成CPFEM-PT-MK集成计算方法.在沟槽内外分别建立代表体积元(representative volume element, RVE)模型;根据沟槽内外区域的应力平衡和变形协调,对RVE模型分别施加不同的边界条件,实时更新缺陷演化并比较沟槽内外应变率差异;根据破坏准则最终确定不同应变路径下的极限应变,最后将不同应变路径下的极限应变汇总,从而得到FLD.

2.1 MK理论

MK理论假设由于材料的不均匀性导致板材在变形前存在缺陷,此缺陷可以采用几何沟槽代替.沟槽部分属于薄弱区域,随着变形进行,沟槽内应变率逐渐大于沟槽外;当沟槽内外应变率差别达到一定程度时,材料在沟槽处发生失效,此时沟槽外对应应变为极限应变.根据MK理论,将板材区域分为两部分:无缺陷的区域A即槽外区域和有缺陷的区域B即槽内区域,如图5所示.

图5

图5   MK理论示意图

Fig.5   Schematic diagram of MK theory


区域A的初始厚度为${t}_{0}^{\mathrm{A}}$,区域B的初始厚度为${t}_{0}^{\mathrm{B}}$,定义初始缺陷f0=${t}_{0}^{\mathrm{B}}$/${t}_{0}^{\mathrm{A}}$.区域B的局部坐标系(n,t,z)与区域A的全局坐标系(X,Y,Z)间的关系为:zZ重合,nX之间夹角为θ,0°θ≤20°[18].区域A承受面内等比例主应变和次应变加载,主应变的方向与X方向重合,次应变的方向与Y方向重合,主应变增量(d${{\epsilon }^{\mathrm{A}}}_{XX}$)与次应变增量(d${{\epsilon }^{\mathrm{A}}}_{YY}$)之间在整个加载过程中保持恒定比例:

ζ=$\frac{\mathrm{d}{{\epsilon }^{\mathrm{A}}}_{YY}}{\mathrm{d}{{\epsilon }^{\mathrm{A}}}_{XX}}$, -0.5≤ζ≤1

根据力平衡条件和变形协调要求,区域A、B在交界面处存在如下关系[19]:

$\left.\begin{array}{l}{{F}^{\mathrm{A}}}_{nn}={{F}^{\mathrm{B}}}_{nn}\\ {{F}^{\mathrm{A}}}_{nt}={{F}^{\mathrm{B}}}_{nt}\\ \mathrm{d}{{\epsilon }^{\mathrm{A}}}_{tt}=\mathrm{d}{{\epsilon }^{\mathrm{B}}}_{tt}\end{array}\right\}$

式中:FAnnFBnn分别为区域A、B在交界面处沿着n方向的载荷;FAntFBnt分别为区域A、B在交界面处沿着t方向的载荷;d${{\epsilon }^{\mathrm{A}}}_{tt}$、d${{\epsilon }^{\mathrm{B}}}_{tt}$分别为区域A、B沿着交界面t方向的应变增量.随着变形的逐渐增加,沟槽角度θ会发生变化,根据Signorelli等[18]的研究,θ遵循如下演化方程:

$\tan \theta=\exp \left(\varepsilon_{X X}^{A}-\varepsilon_{Y Y}^{A}\right) \tan \theta_{0}$

式中:θ0为沟槽初始角度.当区域B沿厚度方向的应变率与区域A沿厚度方向的应变率比值满足${{\stackrel{·}{\epsilon }}^{\mathrm{B}}}_{zz}$/${{\stackrel{·}{\epsilon }}^{\mathrm{A}}}_{ZZ}$≥10,此时区域A的主应变和次应变即为极限应变[19].由于实际板料的缺陷是随机分布的,所以对于每种应变路径,从0°至20°每隔5°计算一个初始角θ0,最终取较小的极限应变作为此应变路径下的极限应变[18].

2.2 CPFEM-PT-MK计算方法

分别建立区域A、B的RVE模型,并施加周期性边界条件,如图6所示.区域A、B的初始织构关系与沟槽角度有关,例如沟槽角度为20° 时,区域B各相的织构为区域A各相织构绕Z轴旋转20° 得到.

图6

图6   RVE模型节点位置

Fig.6   Joint location in RVE model


首先,计算区域A的变形,根据应变路径ζ对区域A施加相应的边界条件:

$\left.\begin{array}{l}{U}_{XYZ}^{000}=0\\ {U}_{XY}^{001}=0\\ {U}_{XZ}^{010}=0\\ {U}_{YZ}^{100}=0\end{array}\right\}$, $\left.\begin{array}{l}{U}_{X}^{100}=U\\ {U}_{Y}^{010}=L\left[{\left(1+\frac{U}{L}\right)}^{\zeta }-1\right]\end{array}\right\}$

式中:U为位移量,上标表示节点,下标表示方向,如${U}_{XYZ}^{000}$表示节点000沿XYZ 3个方向的位移;L为RVE模型的尺寸.先将区域A加载到一个变形量足够大的程度,然后提取A区域RVE变形过程的应变和应力历程,每个历程下应变和应力分量分别为所有单元应力和应变对应分量的加权平均:

$\left.\begin{array}{l}{{\epsilon }^{\mathrm{A}}}_{ij,g}=\frac{1}{N}\stackrel{N}{\sum _{k=1}}{\epsilon }_{ij,k}\\ {{\sigma }^{\mathrm{A}}}_{ij,g}=\frac{1}{N}\stackrel{N}{\sum _{k=1}}{\sigma }_{ij}^{k}\end{array}\right\}$

式中:ij代表坐标分量XYZ,分别对应1、2、3;N为RVE模型的单元数量;下标g表示全局坐标系.将全局坐标系下的应力和应变张量旋转到局部坐标系(l)下:

$\left.\begin{array}{l}{{\epsilon }^{\mathrm{A}}}_{l}=R{{\epsilon }^{\mathrm{A}}}_{g}{R}^{\mathrm{T}}\\ {{\sigma }^{\mathrm{A}}}_{l}=R{{\sigma }^{\mathrm{A}}}_{g}{R}^{\mathrm{T}}\end{array}\right\}$

式中:R为全局坐标系和局部坐标系之间的旋转矩阵,有

R=$\left[\begin{array}{lll}\mathrm{c}\mathrm{o}\mathrm{s}\theta & \mathrm{s}\mathrm{i}\mathrm{n}\theta & 0\\ -\mathrm{s}\mathrm{i}\mathrm{n}\theta & \mathrm{c}\mathrm{o}\mathrm{s}\theta & 0\\ 0& 0& 1\end{array}\right]$

由局部坐标系下的应变和应力可计算区域A在交界面处的FAnnFAntUAtt

$\left.\begin{array}{l}{{F}^{\mathrm{A}}}_{nn}={L}^{2}{{\sigma }^{\mathrm{A}}}_{11,l}\mathrm{e}\mathrm{x}\mathrm{p}\left({{\epsilon }^{\mathrm{A}}}_{22,l}\right)\mathrm{e}\mathrm{x}\mathrm{p}\left({{\epsilon }^{\mathrm{A}}}_{33,l}\right)\\ {{F}^{\mathrm{A}}}_{nt}={L}^{2}{{\sigma }^{\mathrm{A}}}_{12,l}\mathrm{e}\mathrm{x}\mathrm{p}\left({{\epsilon }^{\mathrm{A}}}_{22,l}\right)\mathrm{e}\mathrm{x}\mathrm{p}\left({{\epsilon }^{\mathrm{A}}}_{33,l}\right)\\ {{U}^{\mathrm{A}}}_{tt}=L\left(\mathrm{e}\mathrm{x}\mathrm{p}\right({{\epsilon }^{\mathrm{A}}}_{22,l})-1)\end{array}\right\}$

根据式(5)及初始缺陷f0可求得初始时刻区域B在交界面处沿着n方向的载荷FBnn,沿着t方向的载荷FBnt和沿着t方向的位移UBtt.公式如下:

$\left.\begin{array}{l}{{F}^{\mathrm{B}}}_{nn}=\frac{{{F}^{\mathrm{A}}}_{nn}}{{f}_{0}}\\ {{F}^{\mathrm{B}}}_{nt}=\frac{{{F}^{\mathrm{A}}}_{nt}}{{f}_{0}}\\ {{U}^{\mathrm{B}}}_{tt}={{U}^{\mathrm{A}}}_{tt}\end{array}\right\}$

随着变形的增加,缺陷f按以下公式演化:

$f=f_{0} \exp \left(\varepsilon_{z z}^{\mathrm{B}}-\varepsilon_{Z Z}^{\mathrm{A}}\right)$

f0通过拟合ζ=0应变路径下的极限应变得到,f0=0.991.对于区域B的RVE模型,边界条件设置如下:

$\left.\begin{array}{l}{U}_{XYZ}^{000}=0\\ {U}_{XY}^{001}=0\\ {U}_{XZ}^{010}=0\\ {U}_{Z}^{100}=0\end{array}\right\}$, $\left.\begin{array}{l}{F}_{X}^{100}={{F}^{\mathrm{B}}}_{nn}\\ {F}_{Y}^{100}={{F}^{\mathrm{B}}}_{nt}\\ {U}_{Y}^{010}={{U}^{\mathrm{B}}}_{tt}\end{array}\right\}$

计算出区域B的变形后,判断${{\stackrel{·}{\epsilon }}^{\mathrm{B}}}_{zz}$/${{\stackrel{·}{\epsilon }}^{\mathrm{A}}}_{ZZ}$是否达到临界值(此处临界值取10).若该比值达到临界值,则此时区域A对应的主应变和次应变即为当前应变路径及初始沟槽角度θ0下的极限应变;反之,则需更新缺陷和区域B的边界条件后继续计算.对于给定的应变比ζ和初始沟槽角度θ0,计算流程如图7所示.通过遍历所有角度并取最小的极限应变作为该应变路径下的极限应变,再遍历所有应变路径从而得到FLD.需要注意的是,每一增量步下的f均需要区域B上一步的计算结果${{\epsilon }^{\mathrm{B}}}_{zz}$,因此需要不断对区域B的计算结果进行后处理.该过程通过Python搭建自动化前处理、计算与后处理平台得以实现.

图7

图7   CPFEM-PT-MK方法计算流程

Fig.7   Calculation procedure of CPFEM-PT-MK method


2.3 失效判据

对于一个含有沟槽的板料,在承受面内拉伸变形时,由于沟槽内区域厚度较薄,所以为了维持力平衡,槽内区域必定需要更大的内应力.而在内应力变大的条件下,沟槽内部区域必然承受比沟槽外部区域更大的变形,这又将使得沟槽内的厚度减薄速度比沟槽外部更快,从而使得缺陷比初始更大;缺陷程度增大.为了保持力平衡,则需要更大的变形量,以达到更高的应力强化水平,因此沟槽内部与沟槽外部在缺陷与流动应力增长速率的差异,导致变形量差异越来越大;最终,材料将在某一时刻在沟槽内部发生破坏.图8给出在ζ=1应变路径下,沟槽内外沿厚度方向应变及两者应变率比值随变形的演化,变形程度采用主应变进行度量.图中:ε33为厚度方向应变;ε33Aε33B分别为区域A、B厚度方向应变.由图可见,在变形开始阶段,沟槽内外的厚向应变随变形的增加比较接近;当主应变超过0.07以后,槽内区域的厚向应变逐渐高于槽外,表明槽内的减薄速度更快,且两者的差异随变形增加越来越大,但此阶段沟槽内外厚向应变速率比值仍缓慢上升;当材料接近破坏时,槽内的厚向应变迅速增大,沟槽内外厚向应变速率比值迅速增大,表明材料发生集中性失稳,此时主应变即为ζ=1应变路径下的极限主应变,通过应变比可以得到极限次应变.对于其他应变路径,计算与分析方法类似.

图8

图8   ζ=1应变路径下沟槽内外厚向应变及应变率比值随主应变的演化

Fig.8   Thickness strain inside and outside grooves and ratio of strain rate versus major strain under strain path of ζ=1


3 结果与讨论

3.1 CPFEM-PT-MK预测结果

基于CPFEM-PT-MK方法预测得到的FLD结果如图9所示,其中,ζ=0应变路径下的极限应变为拟合结果,用于确定初始缺陷f0,其余应变路径下的结果均为预测结果.由图可见,预测结果与实验结果基本吻合,只需拟合FLD0,即可通过CPFEM-PT-MK方法得到整条FLD曲线.从计算结果来看,当应变比接近ζ=0.1时,极限主应变最低,根据Wu等[20]的研究,这是由于考虑了弹性变形的影响.应变比ζ在0.1~0.9范围内,随着ζ的增加,极限主应变逐渐增大,但与ζ之间并非线性关系,而近似为抛物线.

图9

图9   QP1180钢板成形极限的预测结果

Fig.9   Predicted FLD of QP1180 blank


3.2 应力应变曲线及变形织构分析

不同应变路径下的等效应力等效塑性应变曲线和残余奥氏体(retained austenite,RA)体积分数演化如图10所示,这里仅给出ζ=-0.5,ζ=0和ζ=1的情况,以便于比较,下同.从相变规律来看,3种应变路径下的相变速度由快到慢依次为:ζ=1,ζ=0,ζ=-0.5;从应力应变曲线来看,ζ=1和ζ=-0.5 应变路径下的应力要高于ζ=0应变路径下的应力,这与相变速度不完全一致,说明QP1180的宏观应力应变曲线不仅受到相变影响,还与应变路径相关.

图10

图10   不同应变路径下应力应变曲线和RA体积分数演化

Fig.10   Stress strain curves and volumetric evolutions of RA under different strain paths


此外,不同应变路径下的应力应变曲线均存在一个明显的过渡转变区,即图10中虚线框部分,这是由于相变接近饱和时,RA体积分数演化规律发生改变导致的.宏观应力为各相应力的平均:

σeq=${\sigma }_{\mathrm{e}\mathrm{q}}^{\mathrm{R}\mathrm{A}}$wRA+${\sigma }_{\mathrm{e}\mathrm{q}}^{\mathrm{p}\mathrm{B}\mathrm{C}\mathrm{C}}$wpBCC+${\sigma }_{\mathrm{e}\mathrm{q}}^{\mathrm{N}\mathrm{M}}$wNM

式中:${\sigma }_{\mathrm{e}\mathrm{q}}^{\mathrm{R}\mathrm{A}}$${\sigma }_{\mathrm{e}\mathrm{q}}^{\mathrm{p}\mathrm{B}\mathrm{C}\mathrm{C}}$${\sigma }_{\mathrm{e}\mathrm{q}}^{\mathrm{N}\mathrm{M}}$分别为RA、初始体心立方(body centered cubic,BCC)和新生马氏体(newly formed martensite,NM)相的应力;wRAwpBCCwNM分别为RA、初始BCC和NM相的体积分数.

不同应变路径下RA、初始BCC和NM的von_Mises等效应力-等效塑性应变的演化如图11所示,各相的应力为该相所有晶粒的平均值,应变为所有相的平均值,这里的平均值是指首先求出每个应力分量的平均值,再根据平均应力分量求出von_Mises等效应力.

图11

图11   各相应力随宏观等效塑性应变的演化

Fig.11   Stress evolution of each phase with respect to macroscopic equivalent plastic strain


可以看出,与宏观应力相同的是,各相在不同应变路径下的应力也存在差异.对于RA和初始BCC相,ζ=0应变路径下的应力水平明显低于ζ=-0.5 和ζ=1应变路径;而对于NM相,则是ζ=-0.5 应变路径下的应力水平较高,ζ=0和ζ=1应变路径下的应力水平较低.这表明各相在不同应变路径下的应力水平相对大小关系存在差异.与宏观应力不同,各相应力并没有图10中宏观应力的过渡转变区特性.因此,宏观应力中的过渡转变区不是由于所含相应力变化造成的.根据式(12),宏观应力还与各相体积分数有关.其中,wpBCC不发生变化,wRA随变形增加不断减少,而wNM随变形增加而增加.wRAwNM最终均趋于饱和,且wRAwNM的总和始终为0.145.由图10可知,在宏观应变εeq为0~0.1时,wRA近似线性减少.随后,不同RA单元的相变速度相继达到饱和,导致RVE中的RA单元相变速度迅速降低.这使得相变带来的强化作用显著减弱,从而使宏观应力出现一段过渡转变区.对比不同应变路径下应力-应变曲线的过渡转变区可知,在ζ=-0.5 应变路径下,进入过渡转变区时的变形量比ζ=1应变路径下进入过渡转变区时对应的变形量更大.这表明相变速度越慢,则过渡转变区越晚出现.

不同应变路径下RA和BCC的变形织构存在显著差异,如图12所示.令φ1Φφ2分别为Bunge欧拉角的3个参量——章动角、进动角和自转角.对于BCC相,在ζ=-0.5 应变路径下,变形织构主要集中在3条线上:① (φ1=0°,Φ=45°,φ2=45°)与(φ1=45°,Φ=0°,φ2=0°)的连线附近;② (φ1=90°,Φ=0°,φ2=45°)与(φ1=45°,Φ=0°,φ2=90°)的连线附近;③ (φ1=45°,Φ=90°,φ2=90°)与(φ1=90°,Φ=45°,φ2=0°)的连线附近.ζ=0应变路径下,变形织构主要集中在以下区域:① (φ1=0°,Φ=45°,φ2=45°)和(φ1=45°,Φ=60°,φ2=60°)的连线附近;② (φ1=45°,Φ=60°,φ2=60°)和(φ1=90°,Φ=60°,φ2=45°)的连线附近.ζ=1应变路径下,变形织构主要在(Φ=55°,φ2=55°)线附近.对于RA相,在ζ=-0.5 应变路径下,变形织构主要集中在以下区域:① (φ1=55°,Φ=45°,φ2=90°)和(φ1=90°,Φ=30°,φ2=45°)的连线附近;② (φ1=45°,Φ=90°,φ2=45°)和(φ1=45°,Φ=45°,φ2=0°)的连线附近.ζ=0应变路径下,变形织构主要集中在以下区域:① (φ1=30°,Φ=45°,φ2=90°)和(φ1=90°,Φ=30°,φ2=45°)的连线附近;② (φ1=30°,Φ=55°,φ2=0°)和(φ1=55°,Φ=90°,φ2=45°)的连线附近.ζ=1应变路径下,变形织构主要在以下区域:① (φ1=45°,Φ=45°,φ2=90°/0°)和(φ1=55°,Φ=90°,φ2=45°)附近.

图12

图12   εeq=0.5时BCC相和RA相在不同应变路径下的变形织构

Fig.12   Deformation texture of BCC and RA phase under different strain paths at εeq=0.5


3.3 织构演化和相变对FLD的影响

材料在变形过程中会出现择优取向,且择优取向与变形方式有关.为了研究织构演化对FLD的影响,计算模型中不考虑晶粒取向更新,保持每个晶粒的取向在变形过程中保持不变,并计算FLD,结果如图13所示.可以看出,织构不发生演化时,FLD整体低于织构发生演化的情况,这表明织构演化提高了QP1180的成形极限.不同应变路径下织构演化对FLD的影响有所差异,织构演化对接近ζ=0应变路径的成形极限几乎没有影响.

图13

图13   织构演化对FLD的影响

Fig.13   Effect of texture evolution on FLD


相变对QP1180钢板FLD的影响如图14所示.不发生相变时,极限主应变最低点在ζ=0应变路径,而不发生相变时在ζ=0.1应变路径,这可能是两者织构演化存在一定差异造成的[21].在FLD左侧ζ=-0.5 应变路径下,发生相变时的极限主应变更低;在FLD右侧,当应变路径ζ<0.2时,发生相变时的极限主应变要高于不发生相变时的极限主应变,而当ζ>0.2时,发生相变时的极限主应变比不发生相变时的极限主应变略低.这说明并非所有应变路径下相变均有利于提高成形极限,这与目前文献中的结论[14,22-23]有所不同.这些文献的结论为相变始终提高材料的成形极限,表现为发生相变时的FLD始终比不发生相变的FLD更高.另外,由于相变在不同应变路径下的速度有所差异,所以导致两者的差距不均匀.文献得出该结论的依据为:发生相变时材料的加工硬化率n更高,而n值更高的材料,通常发生颈缩时对应的极限应变应该越大.

图14

图14   相变对FLD的影响

Fig.14   Effect of phase transformation on FLD


对于发生相变的材料,是否发生相变时材料的n值始终高于不发生相变时的n值,这一点值得讨论.图15对比了3种应变路径下有无相变的von_Mises等效应力及n值随等效塑性应变的变化,并将每种情况下达到极限应变时沟槽内外对应在曲线上的位置进行标记.图中:w.o.表示不发生相变,w.表示发生相变.可以看出,尽管发生相变时的应力始终高于不发生相变时的情况,但两者的n值存在交叉点,位于交叉点变形量之前,发生相变情况的n值更高,而当变形量超过交叉点时,则是不发生相变情况的n值更高.由图10可知,当发生相变且相变接近饱和时,宏观应力存在过渡转变区,此阶段n值迅速下降,从而与不发生相变时的n值发生交叉.在ζ=0应变路径下,极限应变在n值交叉点以前变形阶段出现,因此发生相变时的成形极限更高,而在ζ=-0.5 和ζ=1应变路径下,极限应变出现在n值交叉点以后变形阶段,因此发生相变情况的成形极限更低.

图15

图15   相变对应力应变曲线和n值的影响

Fig.15   Effect of transformation on stress strain curves and n value


4 结论

本文通过Nakajima实验获得QP1180钢板的FLD;并结合MK理论与CPFEM-PT晶体塑性本构模型建立了CPFEM-PT-MK方法,该方法仅通过拟合FLD0确定初始缺陷f0,就可以准确地预测QP1180钢板的FLD;利用模型分析了各相在不同应变路径下的应力应变曲线和织构演化,最后探讨了织构演化和相变对QP1180FLD的影响.具体结论如下:

(1) CPFEM-PT-MK方法能够准确预测QP1180钢板的FLD.根据模型的计算结果,极限主应变最低点不在ζ=0应变路径,而位于近似ζ=0.1应变路径,极限主应变最小值约为0.1;在ζ=-0.5 应变路径下,极限主应变约为0.3;在ζ=1应变路径下,极限主应变约为0.32.

(2) QP1180钢板宏观应力应变曲线存在一个过渡转变区,而RA、初始BCC和NM相的应力应变曲线没有过渡转变区.该过渡转变区与相变接近饱和有关.当相变接近饱和时,由于RA体积分数的变化速率逐渐降低,相变带来的强化作用迅速减小,所以导致过渡转变区的加工硬化率迅速降低.在不同应变路径下,QP1180宏观应力应变曲线及各相应力应变曲线存在差异,在ζ=0应变路径下最低.不同应变路径下RA和BCC相的织构演化存在明显差异.不考虑织构演化时,QP1180成形极限更低.

(3) 当不发生相变时,极限主应变最低点在ζ=0应变路径.相变并非总是提高QP1180的成形极限,而是与应变路径有关.在ζ=0应变路径附近,相变会提高极限应变;而在ζ=-0.5 和ζ=1应变路径下,相变反而降低了极限应变.根据应力应变曲线和n值的对比,发生相变时QP1180的宏观应力存在一个过渡转变区,此阶段n值迅速下降,与不发生相变情况下的n值出现交叉.当极限应变位于n值交叉点以前的变形阶段时,相变能够显著提高成形极限;而当极限应变位于n值交叉点以后的变形阶段时,相变可能会降低成形极限.

参考文献

徐祖耀.

马氏体相变研究的进展(二)

[J]. 上海金属, 2003, 25(3): 1-10.

[本文引用: 1]

XU Zuyao.

Progress in martensitic transformations (II)

[J]. Shanghai Metals, 2003, 25(3): 1-10.

[本文引用: 1]

徐祖耀.

马氏体相变研究的进展(一)

[J]. 上海金属, 2003, 25(3): 1-8.

[本文引用: 1]

XU Zuyao.

Progress in martensitic transformations(Ⅰ)

[J]. Shanghai Metals, 2003, 25(3): 1-8.

[本文引用: 1]

YI H L, SUN L, XIONG X C.

Challenges in the formability of the next generation of automotive steel sheets

[J]. Materials Science & Technology, 2018, 34(9): 1112-1117.

[本文引用: 1]

CHEN X P, NIU C, LIAN C W, et al.

The evaluation of formability of the 3rd generation advanced high strength steels QP980 based on digital image correlation method

[J]. Procedia Engineering, 2017, 207: 556-561.

[本文引用: 1]

MOHAMMED B, PARK T, POURBOGHRAT F, et al.

Multiscale crystal plasticity modeling of multiphase advanced high strength steel

[J]. International Journal of Solids & Structures, 2018, 151: 57-75.

[本文引用: 1]

HU X.

Mechanical behavior and formability of quench and partition steel sheets

[J]. Strength of Materials, 2020, 52(1): 32-39.

GAO X L, MIN J Y, ZHANG L, et al.

Prediction and experimental validation of forming limit curve of a quenched and partitioned steel

[J]. Journal of Iron & Steel Research, International, 2016, 23(6): 580-585.

[本文引用: 1]

KEELER S P. Plastic instability and fracture in sheets stretched over rigid punches[D]. Cambridge, USA: Massachusetts Institute of Technology, 1961.

[本文引用: 1]

HILL R.

On discontinuous plastic states, with special reference to localized necking in thin sheets

[J]. Journal of the Mechanics & Physics of Solids, 1952, 1(1): 19-30.

[本文引用: 1]

STÖREN S, RICE J R.

Localized necking in thin sheets

[J]. Journal of the Mechanics & Physics of Solids, 1975, 23(6): 421-441.

[本文引用: 1]

MARCINIAK Z, KUCZYŃSKI K.

Limit strains in the processes of stretch-forming sheet metal

[J]. International Journal of Mechanical Sciences, 1967, 9(9): 609-620.

[本文引用: 1]

FYLLINGEN Ø, HOPPERSTAD O S, LADEMO O G, et al.

Estimation of forming limit diagrams by the use of the finite element method and Monte Carlo simulation

[J]. Computers & Structures, 2009, 87(1/2): 128-139.

[本文引用: 1]

WANG H, WU P D, BOYLE K P, et al.

On crystal plasticity formability analysis for magnesium alloy sheets

[J]. International Journal of Solids & Structures, 2011, 48(6): 1000-1010.

[本文引用: 1]

余海燕, 高云凯.

基于M-K模型的相变诱发塑性钢板的成形极限研究

[J]. 中国机械工程, 2007, 18(1): 109-113.

[本文引用: 2]

YU Haiyan, GAO Yunkai.

Study on forming limit diagram for transformation-induced plasticity sheet steel based on M-K model

[J]. China Mechanical Engineering, 2007, 18(1): 109-113.

[本文引用: 2]

ARETZ H.

An extension of Hill’s localized necking model

[J]. International Journal of Engineering Science, 2010, 48(3): 312-331.

[本文引用: 1]

BUTCHER C, KHAMENEH F, ABEDINI A, et al.

On the experimental characterization of sheet metal formability and the consistent calibration of the MK model for biaxial stretching in plane stress

[J]. Journal of Materials Processing Technology, 2021, 287: 116887.

[本文引用: 1]

杨浩, 唐伟琴, 李大永.

考虑相变的晶体塑性有限元模型开发及验证

[J/OL]. 塑性工程学报. https://link.cnki.net/urlid/11.3449.TG.20240923.0925.002.

URL     [本文引用: 1]

YANG Hao, TANG Weiqin, LI Dayong.

Development and verification of crystal plastic finite element model considering phase transition

[J]. Journal of Plasticity Engineering. https://link.cnki.net/urlid/11.3449.TG.20240923.0925.002.

URL     [本文引用: 1]

SIGNORELLI J W, BERTINETTI M A, TURNER P A.

Predictions of forming limit diagrams using a rate-dependent polycrystal self-consistent plasticity model

[J]. International Journal of Plasticity, 2009, 25(1): 1-25.

[本文引用: 3]

BONG H J, LEE J, HU X H, et al.

Predicting forming limit diagrams for magnesium alloys using crystal plasticity finite elements

[J]. International Journal of Plasticity, 2020, 126: 102630.

[本文引用: 2]

WU P D, NEALE K W, VAN DER GIESSEN E.

On crystal plasticity FLD analysis

[J]. Proceedings of the Royal Society of London Series A: Mathematical, Physical & Engineering Sciences, 1997, 453(1964): 1831-1848.

[本文引用: 1]

RATCHEV P, VAN HOUTTE P, VERLINDEN B, et al.

Prediction of forming limit diagrams of Al-Mg rolled sheets taking texture into account

[J]. Texture, Stress, & Microstructure, 1994, 22(4): 219-231.

[本文引用: 1]

PARK T, HECTOR L G, HU X H, et al.

Crystal plasticity modeling of 3rd generation multi-phase AHSS with martensitic transformation

[J]. International Journal of Plasticity, 2019, 120: 1-46.

[本文引用: 1]

CONNOLLY D S, KOHAR C P, MISHRA R K, et al.

A new coupled thermomechanical framework for modeling formability in transformation induced plasticity steels

[J]. International Journal of Plasticity, 2018, 103: 39-66.

[本文引用: 1]

/