随着社会的发展，具有节能、绿色等特点的高强钢板，在航空、航天、矿山机械等领域的广泛应用。在高强钢板生产工序中，矫直是一个重要精整环节，直接决定了产品的板型质量和残余应力。由于高强钢板的高强韧性的特点，矫直工艺一直以来是行业内的技术难题。残余应力作为矫直质量的主要指标，直接影响着板材的疲劳强度等参数，进而缩短产品的使用寿命^[1]。关于矫直与残余应力之间关系，国内外学者已进行了大量的研究。管奔等^[2-3]研究了矫直应力残留与反弯特性的耦合关系，以及矫直过程中的应力累积行为。赵军等^[4]通过理论验证了曲梁过弯矫直与直梁纯弯曲之间的等价原理，通过反弯实验证实了其正确性。胡军峰等^[5]针对船用高强钢板残余应力在矫直过程变化的模拟研究，发现矫直过程中上下表面发生对称的表层塑形变形，从而引起残余应力的变化。周存龙等^[6]得出横向残余应力是引起板材侧弯的诱因之一。Zhang等^[7]得出在板材弯曲过程中，残余应力和应变在初始加载位置沿厚度方向减小。詹光曹^[8]研究发现通过增加矫直道次可以降低板材内部残余应力。Lee 等^[9]通过2D数值建模的方法，建立了应力与弯矩之间的关系；通过实验和模拟进行模型修正，研究了板厚方向残余应力的变化。Carvalho等^[10]通过静态模拟分析的方法研究了钢轨残余应力矫直过程中的变化，以及影响残余应力的主要参数。日比野文雄^[11]通过对材料进行纯弯曲假设和理想弹塑性假设，理论计算了矫直后残余应力。Kubo等^[12]进行板材反复弯曲过程中应力应变的关系实验，推导出冷轧钢板的残余应力分布模型。以上研究主要以矫直实验和无残余应力矫直为主，没有系统性研究初始残余应力在矫直过程中的演化规律。

基于矫直理论分析，建立了淬火高强板矫直过程的残余应力计算模型，并采用高强板Q960E为研究材料，对计算模型进行了实验验证及有限元模拟验证，分析残余应力在淬火矫直过程的变化规律。

1 矫直过程残余应力计算模型 1.1 淬火后残余应力计算

假设在淬火过程中，温度关于中性层面对称，过程不发生弯曲，采用温度热传导计算如下。

板材投入液体中冷却，其厚度方向温度上为无内热源,设定热传导为1维热传导，采用1维非定常热传导方程式^[13]如下：

$\left\{\!\!\!\! \begin{array}{l} k\dfrac{{{\partial ^2}T}}{{\partial {y^2}}} = \rho c\dfrac{{\partial T}}{{\partial t}} ; \\ T = {T_2},\;y = \pm \dfrac{h}{2}; \\ \dfrac{{\partial T}}{{\partial t}} = 0,\;y = 0; \\ T = {T_1},\;t = 0 \end{array} \right.$

(1)

式中， $y $ 为厚度方向坐标轴， $h $ 为相对厚度， $t $ 为热传导时间。

计算得到温度T的变化式：

$\begin{aligned}[b] T =& {T_2} + \frac{{4 \Delta T}}{\text{π} }\bigg( {{{\rm{e}}^{ - {p_1}t}}\cos \frac{{\text{π} y}}{h}} \bigg. - \\ &\bigg. {\frac{{{{\rm{e}}^{ - {p_3}t}}}}{3}\cos \frac{{3\text{π} y}}{h} + \frac{{{{\rm{e}}^{ - {p_5}t}}}}{5}\cos \frac{{5\text{π} y}}{h} + \cdots } \bigg) \end{aligned} $

(2)

式中： $\Delta T = {T_1} - {T_2}$ ， ${T_1}$ 为初始温度， ${T_2}$ 为最终温度； ${{{p}}_1} = {\text{π} ^2}k/(\rho c{h^2})$ 为常数， ${{p}}_3=9 \times {{p}_1} $ ， ${{p}_5}=25 \times {{p}_1} $ ，k为导热系数，ρ为材料密度，c为比热容。钢板淬火过程中，通过高温软化进入平衡状态，在刚接触冷水时，表面快速冷却，造成平板截面上应力场的失衡；四周处于无约束状态，冷却的部分将产生弹性变形，达到新的平衡，应变、应力计算分别为式（3）和（4）：

${\varepsilon _x} = {\varepsilon _{\textit{z}}} = {\varepsilon _0} + \frac{y}{\rho } - \alpha T$

(3)

${\sigma _x} = {\sigma _{\textit{z}}} = \frac{E}{{1 - {\nu ^2}}}\left( {{\varepsilon _{xx}} + \nu {\varepsilon _{\textit{zz}}}} \right)$

(4)

板料四周为自由边界条件，处于静力平衡状态：

$\left\{\!\!\!\! {\begin{array}{*{20}{l}} {\displaystyle\int_{ - \frac{h}{2}}^{\frac{h}{2}} {{\sigma _x}} {\rm{d}}y = 0}, \\ {\displaystyle\int_{ - \frac{h}{2}}^{\frac{h}{2}} {{\sigma _x}} y{\rm{d}}y = 0} \end{array}} \right.$

(5)

对于式（3）、（4）、（5）联合求解得：

${\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\begin{aligned}[b] {\sigma _x} =& {\sigma _{\textit{z}}} = - \frac{{\alpha E}}{{1 - \nu }} + \frac{{\alpha E}}{{2\left( {1 - \nu } \right)h}}\int_{ - \frac{h}{2}}^{\frac{h}{2}} T {\rm{d}}y + \\ & \frac{{2y\alpha E}}{{2{h^3}\left( {1 - \nu } \right)}}.\int_{ - \frac{h}{2}}^{\frac{h}{2}} T y{\rm{d}}y \end{aligned}} $

(6)

淬火过程中，假定温度关于中性层面对称，降到室温 ${T_2}$ 的时间为 ${{t}_2}$ ，则式（6）右边第3项为0。化简后如式（7）所示，为淬火后残余应力分布方程式：

${\sigma _x} = {\sigma _{\textit{z}}} = - \frac{{4\alpha E\Delta T{{\rm{e}}^{ - {p_1}{t_2}}}}}{{\text{π} \left( {1 - \nu } \right)}}\left( {\frac{2}{\text{π} } - \cos \frac{{\text{π} y}}{h}} \right)$

(7)

式（3）～（7）中， ${\varepsilon _0}$ 为中心层应变， ${\varepsilon _{{x}}}$ 、 ${\varepsilon _{{\textit{z}}}}$ 为上、下表面的横向应力与纵向应力， $\nu $ 为泊松比， $E$ 为弹性模量， $\alpha $ 为热膨胀系数。

1.2 矫直过程计算模型 1.2.1 研究单元体划分

如图1所示，在矫直弯曲区间中，沿长（x）、宽（z）、厚（y）方向分割成n、m、t份。并对其做出如下假设：

1）弯曲前后板料的横截面依旧保持平面，并且仍然垂直变形后的中轴线。2）矫直过程中，拉伸与压缩区域的应力应变对应关系一致。3）材料为均质且连续的材质，符合Hooke定律，应力中心层与应变中心层相互重合。4）矫直过程塑性变形前后体积保持不变。5）变性材料符合Von−Mises屈服条件。

1.2.2 矫直过程中表面受力分析

矫直过程是弯曲的过程，选取初始状态如图2（a）所示，在经过反向弯曲曲率角的作用下形成图2（b）的状态。此时，中间层在图2（b）的状态下：上表面x方向进行压缩变形，下表面进行拉伸变形。上表面在z方向受到挤压，塑性变形向z方向拉长；下表面在z方向受拉伸，塑性变形向z方向挤压。

由矫直原理^[14]可知，截面在弯曲和卸载后x方向应变分布方程为：

$ {\varepsilon _{wx}} = \left\{\!\!\!\! {\begin{array}{*{20}{l}} {1 + \lambda \left( {\dfrac{{2yC}}{H} - 1} \right){\mkern 1mu} {\kern 1pt} ,\;\dfrac{{{H_{\rm{t}}}}}{2} \le y \le \dfrac{H}{2}};\\ {\dfrac{{2Cy}}{H} ,\;\; - \dfrac{{{H_{\rm{t}}}}}{2} \le y \le \dfrac{{{H_{\rm{t}}}}}{2}};\\ { - 1 - \lambda \left( {\dfrac{{2yC}}{H} + 1} \right),\; - \dfrac{H}{2} \le y \le - \dfrac{{{H_{\rm{t}}}}}{2}} \end{array}} \right. $

(8)

式中，H为板材厚度， ${H_{{\rm{t}}}}$ 为弹性区厚度，h为截面上任意位置厚度， $\lambda $ 为强化系数，C为曲率比。

在板材表面取微元体，第n辊的残余应力状态为 ${\sigma^n _{ij}}$ ，经过第n+1辊形成弯曲应力 ${\sigma^{n + 1} _{ij}}$ 。通过米塞斯屈服准则中应力偏张量的第2不变量 ${J_2'}$ 达到某一值C时，判断其发生塑性变形：

$f\left( {\sigma _{ij}^\prime } \right) = J_2^\prime = C$

(9)

假设弯曲过程中发生的为纯单向拉伸与压缩，根据3维状态纯单向拉伸压缩的体积不变原则^[15]，可得微元体在n+1矫直区间三向应变 ${\varepsilon _x}$ 、 ${\varepsilon _{{y}}} $ 、 ${\varepsilon _{\textit{z}}} $ 关系为：

${\varepsilon _x^{n + 1}} + {\varepsilon _y^{n + 1}} + {\varepsilon _{\textit{z}}^{n + 1}} = 0$

(10)

钢板在矫直过程中，受到矫直辊压力的作用，原始残余应力分布不再满足自身平衡条件。卸载矫直力，弹复后，钢板四周重新处于自由无约束，达到一个新的平衡，此时残余应力关于中性层应满足其自身自平衡条件如下。

力平衡条件：

$\left\{\!\!\!\! {\begin{array}{*{20}{l}} {\displaystyle\int_{ - h/2}^{h/2} {\int_{ - b}^b {{\sigma _x}{\rm{d}}{\textit{z}}{\rm{d}}y = 0} } }, \\ {\displaystyle\int_{ - h/2}^{h/2} {\int_0^l {{\sigma _{\textit{z}}}{\rm{d}}x{\rm{d}}y = 0} } } \end{array}} \right.$

(11)

弯矩平衡条件：

$ \left\{\!\!\!\! {\begin{array}{*{20}{l}} {\displaystyle\int_{ - h/2}^{h/2} {\int_{ - b}^b {{\sigma _x}y{\rm{d}}{\textit{z}}{\rm{d}}y = 0} } },\\ {\displaystyle\int_{ - h/2}^{h/2} {\int_0^l {{\sigma _{\textit{z}}}y{\rm{d}}x{\rm{d}}y = 0} } } \end{array}} \right. $

(12)

根据平衡条件式（11）和（12）判断是否发生中心层偏移，若发生偏移，重新计算应变。

在进行曲率比计算过程中，忽略中性层偏移与曲率积分求解迭代的耦合因素，计算得出矫直过程中各区间的曲率变化。由各个区间的计算结果，进行多个单元独立变形计算，根据体积不变原则和残余应力自平衡特性，获得截面应力最终结果，具体计算流程如图3所示。

2 计算结果分析

计算采用的材料为Q960E，其力学性能参数如表1所示。

利用淬火后残余应力分布模型计算，结果如图4所示。从图4中可以发现，在板材表面形成–276.3 MPa的压应力，板材内部呈现与相反的应力，并且以三角函数形式分布，反映了在冷却过程中温度的变化非线性，符合板材实际淬火应力分布。

平行辊矫直的压下方案中，平行压下相比于采用理想小变形调整方案矫直精度高^[16]，故所示采用平行辊矫直方案。矫直机为11辊辊式矫直机，辊径95 mm，边辊距为130 mm，辊距为100 mm，压下量为1 mm，入口为7 mm，出口7 mm。矫直过程中残余应力计算结果如图5所示。矫直过程的曲率变化，通过曲率积分迭代求解的方法，计算出2～10辊弯曲曲率为0.69、–1.38、1.60、–1.75、1.89、–1.85、1.84、–1.56、0.95。

从图5可知：表面残余应力通过第2辊矫直后变化不明显；矫直进行到第6辊，残余应力变化最大。上表面x方向残余应力降低到82.1 MPa，为拉应力；z方向降到83.8 MPa，为压应力。下表面x方向降到116.7 MPa，为压应力；z方向38.6 MPa，为拉应力。上下表面残余应力的结果不同，原因在于，上下表面矫直过程历程不同。此外，在矫直的全过程，心部的残余应力几乎没有改变，这是因为在矫直过程中，表面发生塑性变形，心部发生为弹性变形，在残余应力自平衡调节过程中，有些许改变，但是变化不大。这也验证了塑性变形是影响残余应力变化的主要因素。板材矫直后表面残余应力的大小，影响板材作为结构材料的使用，是衡量板材质量重要的技术指标。故采取合理的矫直方案，可以达到控制表面残余应力的可能，从而得到更优质量的产品。

3 实验及模拟验证 3.1 验证方法及流程

本文中计算模型的验证方法由实验和模拟两部分组成，具体进行流程如图6所示。

实验部分，采用Q960E随电阻加热炉加热到900 ℃，进行水淬处理，得到淬火应力的高强版，并通过实验室十一辊辊式矫直机进行矫直。全过程残余应力的采集，采用取点测试方法，通过4行5列点分布，X射线衍射法^[17]采集，取平均值方法。

模拟部分，规定在此坐标系中，横向应力为x轴方向应力，纵向应力为z轴方向应力，厚度方向应力为y轴方向应力。通过使用MSC.MARC非线性有限元软件，采用8节点单元，共27000个；热传导采用瞬态温度分析算法，对高强钢Q960E温度场模拟，对所得数据进行热应力顺序耦合分析^[18]，获得板材热处理模拟后的残余应力状态。将此残余应力状态作为矫直仿真的初始条件，进行矫直模拟，单元个数与淬火模拟相同，分析方法为弹塑性模型和大变形更新的拉格朗日算法。

3.2 实验设备

矫直设备采用实验室全液压式十一辊辊式矫直机，上排为5个矫直辊，下排为6个矫直辊，辊颈95 mm，边辊距130 mm，辊距为100 mm，压下量与计算一致。

残余应力检测方法为X射线检测法，型号为加拿大Proto公司生产的iXRD−300便携式残余应力测试仪，选取Cr靶，功率为20 kV，4 mA。测试开始之前，材料表面需要用砂纸打磨光洁，去除氧化皮层及表面锈层；然后，用腐蚀液（10%HNO₃+90%CH₃CH₂OH）进行腐蚀并且抛光，再由丙酮擦拭干净。

3.3 模拟结果

根据以上淬火条件进行实验和模拟仿真，分别得到各自的应力状态结果。将淬火后的Q960E板材表面进行X射线残余应力测量。两者结果及其误差如表2所示，实验与计算所得276.3 MPa相比较，误差为2.8%～9.1%，基本吻合。

图7为矫直过程中稳定矫直阶段，第6辊位置板材所受应力状态。在矫直过程中，与矫直辊接触的表面处于压应力状态，相背面的处于拉应力状态，形成了上表面为压应力–拉应力–压应力，下表面处于拉应力–压应力–拉应力的状态，并且经过每一个矫直区间，皆是如此，该应力状态也反映了矫直过程就是反复拉压的过程。

由图8矫直后表面残余应力分布看出：在板材的上表面，x方向应力减少到42～79 MPa，受拉应力状态；z方向应力减少到73～105 MPa，受压应力状态。在板材的下表面，x方向应力减少到56～115 MPa，受压应力状态；z方向应力减少到23～58 MPa，为拉应力状态。这是因为上下表面矫直过程中受力路径不完全相同，上表面经历5道次压弯和4道次拉伸状态,下表面经历4道次压弯和5道次拉伸状态。矫直结果显示，板材前后两端残余应力仍然比较大，是由于板材头和尾部属于空矫区域。利用X射线对矫直后的上下表面进行残余应力测量发现，测量结果与模拟结果相差不大，在允许误差范围之内。

此外，在板材前后尾部及其宽度边缘，形成较大应力。主要原因是，首尾矫直时形成空矫区，矫直效果不理想，宽度边缘属自由端面，矫直中应力对其影响较大。

3.4 计算、模拟与实验对比分析

图9为残余应力计算结果与模拟结果对比数据。

由图9可知，在矫直过程中，中心位置残余应力几乎没有改变。通过计算发现，只在表面发生塑性变形，体现了塑性变形是引起残余应力改变的重要因素。

上下表面的残余应力都降低到了120 MPa以下，应力分布与模拟结果吻合较好。其中，产生差异地方的原因是，在模拟过程中网格划分稀疏，并且在计算中忽略张力、球张量、鲍辛格效应等因素。

表3为矫直实验后沿着板材长方向上，上下表面横向与纵向应力测量结果，从数据可以看出，矫直计算模型与矫直实验所得结果，误差范围较小，属于合理的范围之内，验证了与计算模型的一致性。此外计算结果、实验结果都与矫直有限元模型结果相接近，也验证了连续耦合矫直模拟的可行性。

4 结　论

建立了淬火高强板矫直过程的残余应力计算模型，对计算模型进行了实验验证及有限元模拟验证，并分析了高强度钢板材的矫直过程及矫平过程中对板材残余应力的影响，得出以下结论：

1）本文计算模型与实验误差在15%以下，在可忽略的误差范围以内，实验和模拟也验证了计算模型的实用性。

2）计算、实验、模拟结果表明，辊式矫直过程是降低残余应力的有效方法。矫直前上下表面x、z方向都处于压应力–276 MPa左右。矫直后上表面残余应力x方向减少到42～79 MPa，z方向减少到–89～–122 MPa；下表面x方向减少到–56～–115 MPa，z方向减少到23～58 MPa。