相变潜热对焊接数值模拟的影响∗

许波，周建平，许 燕，李建军，王恪典,2†

(1.新疆大学机械工程学院，新疆 乌鲁木齐 830047；2.西安交通大学机械工程学院，陕西 西安 710049)

0 引言

焊接过程中，焊件不仅温度发生变化，且在一定温度范围内焊缝区金属还会发生相变.热源靠近时，焊缝区附近母材金属熔化吸收相变潜热，热源离开时，母材金属凝固而释放相变潜热[1].在早期的焊接数值计算中，受研究手段和计算能力限制，研究人员通常简化研究对象，会忽略相变潜热；随着研究方法的不断完善和计算机水平的提高，计算过程中逐渐考虑相变潜热对焊接数值模拟的影响.相变潜热的几种常用处理方法有：温度回升法[2]、显热熔法[3]、等效热容法[4]、热焓法[5]和液体分数法[6].本文借助有限元软件ANSYS进行数值计算，以热焓法考虑相变潜热，通过是否考虑相变潜热两种情况进行数值计算，并进行实验验证模拟的准确性和可行性，探究相变潜热对焊接数值模拟的影响.

1 实验设备及材料

实验设备包括松下TA-1800型焊接机器人和YD-350GL型焊机.焊接过程用Thermovision A40M型红外热像仪进行拍摄，记录焊接过程中不同时刻温度场分布结果，并通过配套软件FLIR R ＆ D software处理拍摄到的温度场数据.实验方法为熔化极气体保护焊（GMAW），焊接电压22.9 V，焊接电流154 A，焊接速度为0.007 m/s，保护气为混合气（75%Ar+25%CO2）.实验材料选择X80管线钢，尺寸为150 mm×200 mm×10 mm.在钢板上开V型槽，进行单道焊接实验，焊接前对坡口及焊缝两侧用丙酮处理，坡口尺寸如图1所示.

图1 焊缝几何模型Fig 1 Geometric Model of Weld

2 有限元分析

2.1 温度场控制方程

焊接温度场的分析属于典型的非线性瞬态热传导问题.求解能量守恒方程可得整个区域的热焓分布H(T)，通过已知焓和温度的关系即可确定整个区域的温度场分布，进而确定相变界面的位置.能量守恒方程通用表达式如下[5]：

式中：ρ为密度，H(T)为热焓，c为比热容，u为流体速度，k为热导率，g(T)为内热源.

2.2 应力场

焊接过程中，焊缝金属在热循环下发生弹性变形、塑性变形、体积变化、蠕变等现象.焊接热弹塑性应力计算中，材料屈服服从Von Mises屈服准则，首先通过式（2）计算各节点位移增量，再根据位移增量计算出应力增量.材料处于弹性或塑性状态时应力应变关系和某一单元平衡方程如下[6]：

式中：dσ为应力增量；dε为应变增量；D为弹性或塑性矩阵；C为温度相关向量；dFe为单元节点上力的增量；dRe为单元初应变等效节点力增量；Ke为单元刚度矩阵；dδe为节点位移增量.

2.3 相变潜热

本文选用热焓法考虑焊接过程中的相变潜热.热焓法的优势是不需要实时跟踪两相界面的变化，简化了相变潜热问题的处理过程.热焓可表示为密度和比热容的乘积对温度的积分，定义式如下[7]：

式中：ρ(T)为随温度变化的密度；c(T)为随温度变化的比热容；T为温度.

2.4 热源模型

合适的热源模型很大程度决定模拟结果的准确性和可信度[8].考虑热源移动时前方的温度梯度大于后方温度梯度以及气体保护焊中电弧挺度对熔池深度方向的影响，因此选用Goldak等人提出的一种体积分布热源模型—双椭球体热源模型，该热源方程如下[9]：

前半椭球内部的热流密度分布如式（5）：

后半椭球内部的热流密度分布如式（6）：

式中：a，b，cf，cr为双椭球热源的椭球宽度、椭球深度、前半椭球长和后半椭球长，其中，cf=0.5b，cr=2b；ff是前半椭球的热输入比，fr是后半椭球的热输入比，η为电弧热效率，取0.85，U为电弧电压，I为电弧电流，x，y，z为直角坐标系坐标，v为焊接速度.

2.5 边界条件

焊接热传导分析中，热力学边界条件包括对流换热qconv和辐射换热qrad，焊缝中心以辐射换热为主，远离焊缝区域以对流换热为主，因此考虑总的换热系数h[10,11]：

式中：h为总的换热系数；T为焊件表面温度.

2.6 网格划分

本文借助ANSYS有限元软件进行数值模拟.考虑焊接对称性，只需建立一半模型.为提高计算效率，将焊缝中心网格划分细密，远离焊缝的网格划分逐渐粗大，网格划分模型如图2所示.温度场计算时，单元选择具有8节点、每个节点只有一个温度自由度的Solid70单元，将计算数据保存在结果文件中.应力场计算时，将单元转换为相应的结构单元Solid45单元，然后将温度数据作为载荷施加在焊件上，计算出焊接应力值.运用ANSYS自带的生死单元技术考虑焊缝逐步填充过程.由于此时焊缝前方热通量较小，模拟结果有更高的峰值温度并且更接近实测温度结果[12,13].模拟初始条件为：环境温度25℃.

图2 网格划分Fig 2 Mesh generation

3 结果与分析

3.1 相变潜热对温度场影响

图3(a)为实测温度结果，由于FLIR R ＆ D software显示最低温度为300℃，所以300℃以下皆为黑色.从图3(a)可以看出，随热源的移动，熔池表面呈现双椭圆状，熔池最高温度为2 216℃.图3(b)为粗晶区一点热循环曲线，从中可以看出，考虑相变潜热计算出的温度数据与实验温度在升温阶段和降温前半段比较吻合，模拟计算结果在后半段降温较快.综合而言，模拟计算结果与实际结果吻合，从而验证模型的准确性和可行性.不考虑相变潜热计算出的峰值温度比实验结果高100℃左右，而在冷却阶段温度变化率基本保持一致.

图3 温度场实测结果Fig 3 The measured results of the temperature field

图4(a)、图4(b)为同一时刻下两种情况的温度场计算结果，图中灰色区域为熔池（高于1 350℃），从中可以看出，考虑相变潜热的情况下，熔池最高温度为2 216℃，与实验温度值相吻合；不考虑热焓熔池最高温度为2 367℃，二者相差近150℃，其主要原因是母材金属熔化过程会吸收相变潜热，导致用于升高温度的热量减小，因此计算出的熔池峰值温度相对较低.

图4 温度场模拟结果Fig 4 Temperature field simulation results

焊接实验完成后借助线切割机切割焊件，制作金相试样.配置浓度为3%硝酸酒精溶液，通过擦拭金相试样表面进行金相腐蚀，用金相显微镜观察焊缝形貌，如图5所示，不同区域因显微组织不同而分界明显.从中心向边缘依次为焊缝区、热影响区和母材区.

将两种情况下计算出的焊缝熔池尺寸与实际焊缝尺寸对比，结果如表1所示.从中可以看出考虑相变潜热计算出的熔池尺寸偏小，主要是由于母材金属熔化时吸收相变潜热，导致熔池向周围传递的热量减少，高温加热区域变小，最终整个熔池尺寸偏小.考虑相变潜热计算出的熔池熔宽与实际熔宽吻合，而熔池深度相对于不考虑相变潜热时，与实际熔深相差较大.综合来看，考虑相变潜热时，计算出的熔池形貌尺寸与实际焊缝形貌尺寸更为接近.

图5 焊缝宏观形貌Fig 5 Macro-morphology of weld

表1 焊缝尺寸对比Tab 1 Weld Size Contrast

3.2 相变潜热对应力场影响

图6(a)至图6(d)分别对比了两种情况下横向残余应力、厚度方向残余应力、纵向残余应力和等效残余应力的模拟计算值.从图中可以看出，考虑相变潜热和不考虑相变潜热计算出的残余应力趋势基本保持一致.焊接残余应力峰值位于热影响区附近，且存在较大的应力梯度，应力状态也较复杂.因此，焊缝热影响区是整个焊件的薄弱区，容易出现应力缺陷.

从图6(a)可以看出，考虑相变潜热时的最大横向残余应力为255 MPa，而不考虑相变潜热时最大横向残余应力为262 MPa，且前者应力峰值更靠近焊缝中心；图6(b)为沿厚度方向残余应力对比情况，从中可以看出，受相变潜热影响，不考虑相变潜热所得最大拉应力相对较大，而最大压应力相对较小.横向残余应力和厚度方向残余应力峰值存在差异，与两种情况下焊件经历的不同热循环有关，不考虑相变潜热时焊件拥有更高的峰值温度，则其残余应力峰值越高.对于纵向残余应力，如图6(c)所示，是否考虑相变潜热对其影响不大，只在焊件边缘区域有略微影响.等效残余应力如图6(d)所示，应力峰值基本一致，在距焊缝中心15 mm处，不考虑相变潜热应力值比考虑相变潜热应力值低34 MPa，以及边缘地区有略微差异.总体而论，相变潜热对焊接残余应力的计算有影响，但影响效果不明显.

图6 焊接残余应力Fig 6 Welding residual stress

4 结论

采用数值模拟与实验相结合的方法，计算出有无相变潜热两种情况下的温度场及应力场分布情况，得出相变潜热对焊接数值模拟的影响情况，结论如下：

（1）建立了三维有限元数值分析模型，模拟结果与实验结果吻合较好，证明了模型的准确性和可行性.

（2）不考虑相变潜热所得熔池峰值温度比实际峰值温度高150℃，热影响区峰值温度比实际峰值温度高100℃，且考虑相变潜热计算出的熔池尺寸与实际焊缝形貌尺寸更为接近.

（3）相变潜热对横向焊缝残余应力和厚度方向残余应力峰值有略微影响，而对纵向残余应力峰值和等效残余应力峰值影响较小.