爆炸焊接界面波物质点法三维数值模拟＊

王宇新，李晓杰，王小红，闫鸿浩，孙 明

（大连理工大学工程力学系，辽宁 大连116024）

在近几十年中，对于爆炸焊接的理论研究以及工业应用得到了迅速发展［1－2］。在炸药爆轰载荷作用下，2金属板之间高速碰撞形成金属射流和周期性波纹状的复合界面，通常被称为界面波。金属复合板强度以及复合率与界面波状态密切相关，金属爆炸焊接只有在复合界面上形成大小适当的波纹才能保证复合板的结合强度。关于爆炸焊接复合界面波的形成机理［3－4］有4种经典理论：（1）复板侵彻刻入机理，复板的周期性刻入基板并在驻点前沿形成驼峰而形成波状界面的特征；（2）亥尔姆霍尔兹失稳模型，复板再入射流和基板之间存在速度差而形成界面波；（3）应力波模型，在复板与基板的碰撞点处的可压缩冲击波与在自由边界反射的稀疏波相互作用而形成界面波；（4）卡门涡街模型，在金属板碰撞点附近材料的卡尔曼涡街流动形成界面波。然而，这些界面波理论模型还不完善，迄今为止，爆炸焊接界面波的形成机理仍没有统一的认识。

由于界面波的形成涉及多种物理现象，单纯使用固体力学计算模型对此很难进行分析，而采用爆炸焊接实验又难以观察到界面波的动态形成过程，因此，本文中基于冲击动力学和物质点法［5］对爆炸焊接界面波的形成过程进行三维数值模拟，并开展相关的爆炸焊接实验。

1 物质点法基本理论

采用物质点法求解冲击动力学问题时，首先定义背景网格，然后在背景网格内将连续体剖分为具有集中质量的离散质点，再利用形函数实现质点和背景网格节点之间的2次映射，通过求解连续方程和运动方程以及材料本构方程计算所有质点的应变和应力值［5］，在计算过程中背景网格固定不变。与其他数值计算方法相同，物质点法满足的连续方程、运动方程和能量方程分别为：

式中：ρ为质量密度，a为加速度，b为单位体积力，v为速度，σ为柯西应力张量，e为质量热力学能，应变率

将连续体离散为Np个质点后，为了利用这些质点来描述连续体的变形过程，假定质点p（p＝1，2，…，NP）在t时刻位置为xtp，质量为Mp，密度为ρtp，速度为vtp，柯西应力张量为σtp，将每个质点看作集中质量点，并且每个离散质点的质量保持不变，那么连续方程自然满足。为了求解运动方程则需要采用弱形式［6］：

式中：σs为比应力张量，σs＝σ／ρ；Ω为当前连续体构形，∂Ωc为指定应力边界；w为试函数；假设边界上应力为c，且cs＝c／ρ，w在指定位移边界上为零。这样每个质点密度写成δ函数形式［6］：

背景网格节点包含了多个参量，坐标x（t）i、位移u（t）i、速度v（t）i和加速度a（t）i以及试函数w（t）i等，统一使用变量r（t）i（i＝1，2，…，Nn）表示，则有［7］：

利用形函数实现从背景网格节点到质点的映射，更新质点的位移、速度、应变等参量。由于物质点法采用显式积分算法，将网格节点插值得到的参数代入式（4）可得：

将运动方程（7）写成以下形式：

采用行求和方式获得对角质量阵，这样式（8）和式（10）分别变为：

从而式（11）变为：

基于材料本构模型对质点应力进行更新，计算应力增量为：

式中：E是切线刚度张量。计算质点当前时刻的应力张量为：

计算当前时刻质点内能为：

物质点法的理论模型以及数值计算方法都较成熟，对于爆炸冲击问题通常采用显式积分算法完成。关于物质点法的具体计算步骤和求解过程可以参考文献［8－12］。

2 爆炸焊接数值计算模型

2.1 爆炸焊接布置形式

实际爆炸焊接工程中，一般采用飞板与基板平行的方式进行爆炸焊接，2金属板之间保持一定间距，为了避免飞板上表面被炸药烧伤，还要加一薄缓冲层，如图1所示，爆炸焊接装置主要由雷管、炸药层、缓冲层、飞板、支撑以及基板等组成。

图1 爆炸焊接布置简图Fig.1 Distribution pattern of explosive welding

在爆炸焊接过程中，雷管起爆后，炸药会在起爆端不远处达到稳定的爆轰状态，爆轰波以vd的速度向炸药未起爆区域传播，如图2所示。

图2 爆炸焊接过程示意图Fig.2 The process of explosive welding

炸药爆炸产物驱动飞板向下运动，在预置间隙中经历加速过程，与基板发生碰撞之前获得较高的飞行速度vp，称为碰撞速度。随着炸药稳定爆轰的推进，碰撞点会以稳定速度vc向未焊合区域行进，该速度称为碰撞点移动速度，平行布置时vc＝vd。在碰撞点附近会产生高压，同时金属复合界面会受到剧烈的剪切作用，进而发生强烈的塑性变形，导致界面材料升温，甚至发生熔化流动，并形成向前的微喷射射流使界面完成自清理，最终实现爆炸焊接。由于应用物质点法对爆炸焊接界面波的模拟研究计算量很大，不可能对炸药爆轰、基板与飞板碰撞、焊接界面形成的全部过程进行模拟。同时，基于爆炸焊接的实际情况，炸药爆轰驱动飞板，飞板在预置间隙加速后，爆轰产物对飞板背面施加的压力已经下降为0.1～1GPa量级，飞板以一定速度与基板碰撞后的压力却高达10～100GPa量级，远大于炸药爆轰对飞板的压力，爆轰压力在2金属板接触后对爆炸焊接界面的影响已经很小。因此，在对爆炸焊接界面的数值模拟过程中可忽略炸药爆炸驱动过程，单独对高速碰撞过程进行建模，爆炸焊接的计算模型如图3所示。

图3 爆炸焊接数值计算模型示意图Fig.3 Calculation model of explosive welding

飞板向下以一定速度和倾斜角度β碰撞基板，选择飞板与基板的倾斜角度应该在爆炸焊接窗口上下限的最大与最小弯折角之间，焊接窗口上下限是指保证爆炸焊接可以复合的飞板碰撞基板的最大与最小速度，很多金属的爆炸焊接窗口上下限已经通过实验测定。在已知炸药爆速的前提下，由飞板碰撞点速度计算公式［13］：

可以获得本文计算模型飞板的初始倾斜角度，这样设计的计算模型与轻气炮加速、电磁加速碰撞焊接情况是完全相符的。当确定了飞板与基板材料的爆炸焊接窗口碰撞点速度的上下限以及炸药的爆轰速度后，由式（20）就能够计算出飞板与基板的碰撞角，飞板的初始速度、倾斜角度和垂直向下的初始速度也可以由此确定。

2.2 爆炸焊接材料模型

确定爆炸焊接的计算模型后，假定飞板和基板材料都为普通碳钢，普通碳钢的爆炸焊接窗口下限为272m／s，飞板的初始速度设计为300m／s，碰撞点速度处于焊接窗口中间。由于爆炸焊接整个过程历时较短，在微秒量级范围内，所以忽略热传导并假定为绝热过程。2块金属板碰撞过程中，复合界面材料发生塑性流动。采用物质点法计算与有限元法类似，为了描述爆炸焊接复合界面处的金属塑性流动变形，材料模型采用Johnson－Cook塑性模型［14］：

式中：σy为Von Mises流动应力，εpe为等效塑性应变，为相对等效塑性应变率，取1s－1，B、n为材料应变硬化特性参数，C为材料应变率敏感特性参数；式中第1个乘积因子表示当ε·＊＝1.0时流动应力同等效塑性应变之间的关系，第2个乘积因子表示应变率强化效应。

爆炸焊接过程还涉及到材料在高压条件下的力学行为，所以对其还要运用相应的状态方程进行描述，这里材料状态方程选用Mie－Grüneisen状态方程。当材料处于拉伸膨胀状态时，状态方程形式为：

当材料处于压缩状态时，状态方程形式为：

式中：cg、a、s1、s2和s3是与材料冲击压缩特性有关的常数，cg是冲击波速度（us）－质点速度（up）曲线的截距，s1、s2和s3是us－up曲线斜率的系数，a 是一阶体积修正系数，γ0是Grüneisen系数，E 为材料比内能，体积变化率μ＝ρ／ρ0－1。

本文中普通钢状态方程模型参数［15－16］为：ρ＝7 850kg／m3，a＝0.43，G＝124MPa，A＝400MPa，B＝220MPa，n＝0.08，C＝0.025，cg＝4 573m／s，s1＝1.33，s2＝0，s3＝0，γ0＝1.07。

3 爆炸焊接界面波形成数值模拟

假设爆炸焊接计算模型的尺寸板X方向板长100mm，Y方向板厚20mm，Z方向板宽16mm，在Z方向的板两侧面施加固定约束条件。利用物质点法对本文中设计的爆炸焊接计算模型进行三维数值模拟之前，需要划分背景网格和质点单元等，三维背景网格单元为六面体立方体单元，网格单元大小为2mm，在X、Y、Z3个方向划分的背景网格单元数量分别为100、100和50，每个背景网格单元内的质点数量划分为8个，飞板和基板的材料模型选择Johnson－Cook模型，材料在高压下的力学行为选用Mie－Grüneisen状态方程进行描述。爆炸焊接计算模型的三维前处理结果见图4。

图4 爆炸焊接数值计算模型前处理Fig.4 Preprocess of explosive welding

爆炸焊接数值模拟是利用自主研发的物质点法程序SPS2.0完成的，该程序是应用C＋＋编制实现的，并提供了基本的后处理功能，图5给出了爆炸焊接数值模拟的三维变形图，图6给出了爆炸焊接数值模拟的中间剖面变形图。爆炸焊接中间剖面的等效应力σe后处理云图，见图7。爆炸焊接中间剖面的压力p场后处理云图，见图8。

图5 爆炸焊接复合界面波的形成Fig.5 Interfacial wave formation of explosive welding

图6 爆炸焊接复合界面波x－y中间剖面Fig.6 Interfacial wave formation of explosive welding at x－ysection

图7 爆炸焊接x－y中间剖面等效应力Fig.7 Effective stress of explosive welding at x－ysection

图8 爆炸焊接x－y中间剖面压力Fig.8 Pressure of explosive welding at x－ysection

为了验证物质点法计算的准确性，本文中设计了普通钢板爆炸焊接实验：飞板和基板采用平行布置，基板和飞板的长度均为1 000mm，宽度均为400mm，厚度均为20mm；炸药采用铵油，爆速为2 400m／s；2金属板之间的支撑间距为10～12mm，保证飞板碰撞速度vp在爆炸焊接窗口内，约300m／s。爆炸焊接实验后，取复合板中间部分做为试样进行研磨，利用电子金相显微镜获得复合界面金相显微组织照片，如图9所示。

图9 爆炸焊接复合界面波状态金相照片Fig.9 Metallograph of interfacial wave formation

将爆炸焊接复合界面的金相组织以及界面波纹形状（图9）与数值模拟结果（图6）进行比较可知，两者的界面波形状基本一致，由此说明采用物质点法对爆炸焊接界面波形成问题的模拟结果具有一定的准确性和可靠性。

分析复合界面的金相照片，可以发现界面处的金属材料有熔化现象，这与数值模拟结果（图7）是相符的。在飞板与基板碰撞点位置附近计算的等效应力等于零，说明飞板与基板在冲击碰撞位置的材料几乎处于流体状态。为什么在碰撞点附近区域出现熔化流体现象，其中重要原因是该处的压力最高。数值计算获得的压力大小在12GPa以上，远高于材料强度。通过数值模拟得到的爆炸焊接压力分布图（图8）可以直观地说明这一点。

4 结 论

（1）爆炸焊接界面波的形成机理是由于飞板与基板的高速碰撞，在碰撞点处压力最高，材料发生熔化形成流体涡旋变形并且迅速冷却而复合为一体，这个过程随着碰撞点的移动不断地向前推进，使得2金属板复合形成周期性波状界面。

（2）物质点法的三维数值模拟结果与爆炸焊接实验结果比较符合，由此说明物质点法用于模拟爆炸焊接界面波的形成过程具有实际工程应用价值。

（3）此外，基于本文数值模拟结果也证明了物质点法用于爆炸焊接和冲击动力学问题计算具有突出的优点：一是避免了有限元法在复合界面处难以克服的网格畸变问题；二是与其他无网格法相比，物质点法计算效率更高，本文数值模拟质点数量10万多，台式机计算时间约1h，而利用光滑粒子流体动力学法则需要几天时间。因此，物质点法在爆炸焊接、高速冲击、塑性大变形等力学问题计算中具有广泛的应用前景。

［1］邵丙璜，张凯.爆炸焊接原理及工业应用［M］.大连：大连工学院出版社，1987：1－5.

［2］Grignon F，Benson D，Vecchio K S.Explosive welding of aluminum to aluminum：Analysis，computations and experiments［J］.International Journal of Impact Engineering，2004，30（10）：1333－1351.

［3］Kahraman N，Gülenc B，Findik F.Joining of titanium／stainless steel by explosive welding and effect on interface［J］.Journal of Materials Processing Technology，2005，169（2）：127－133.

［4］张登霞，李国豪，周之洪，等.材料强度在爆炸焊接界面波形成过程中的作用［J］.力学学报，1984，16（1）：73－80.Zhang Deng－xia，Li Guo－hao，Zhou Zhi－hong，et al.Effect of material strength on forming process of explosive welding interface wave［J］.Acta Mechanica Sinica，1984，16（1）：73－80.

［5］Sulsky D，Chen Z，Schreyer H L.A particle method for history－dependent materials［J］.Computer Methods in Applied Mechanics and Engineering，1994，118（1／2）：179－186.

［6］Chen Z，Hu W，Shen L.An evaluation of the MPM for simulating dynamic failure with damage diffusion［J］.Engineering Fracture Mechanics，2002，69（17）：1873－1890.

［7］Wang Y X，Beom H G，Sun M，et al.Numerical simulation of explosive welding using the material point method［J］.International Journal of Impact Engineering，2011，38：51－60.

［8］Wiezckowski Z.The material point method in large strain engineering problems［J］.Computer Methods in Applied Mechanics and Engineering，2004，193（39／40／41）：4417－4438.

［9］Hu W，Chen Z.Model－based simulation of the synergistic effects of blast and fragmentation on a concrete wall using the MPM［J］.International Journal of Impact Engineering，2006，32（12）：2066－2096.

［10］黄鹏，张雄，马上.基于 OpenMP的三维显式物质点法并行化研究［J］.计算力学学报，2010，27（1）：21－27.Huang Peng，Zhang Xiong，Ma Shang.Parallelization of 3Dexplicit material point method using OpenMP［J］.Chinese Journal of Computational Mechanics，2010，27（1）：21－27.

［11］Zhang D Z，Zou Q S，VanderHeyden W B.Material point method applied to multiphase flows［J］.Journal of Computational Physics，2008，227（6）：3159－3173.

［12］Sulsky D，Schreyer L.MPM simulation of dynamic material failure with a decohesion constitutive model［J］.European Journal of Mechanics：A／Solids，2004，23（3）：423－445.

［13］Li X J，Mo F，Wang X H.Numerical study on mechanism of explosive welding［J］.Science and Technology of Welding and Joining，2012，17（1）：36－41.

［14］Daridon L，Oussouaddi O，Ahzi S.Influence of the material constitutive models on the adiabatic shear band spacing：MTS，power law and Johnson－Cook models［J］.International Journal of Solids and Structures，2004，41（11／12）：3109－3124.

［15］Steinberg D J.Equation of state and strength properties of selected materials［R］.Livermore：Livermore National Laboratory，1996.

［16］Johnson G R，Cook W H.Fracture characteristics of three metals subjected to various strains，strain rates，temperatures and pressures［J］.Engineering Fracture Mechanics，1985，21（1）：31－48.