冷喷涂(Cold spraying,CS)作为新型涂层制备方法，最早是由前苏联科学院西伯利亚分院的理论与应用力学研究所的Alkhimov和Papyrin于1990年提出的概念[1]。在冷喷涂过程中，喷涂粒子被加速到较高的速度(200~1 000 m/s)，在固态下碰撞基体(或已沉积涂层)，发生强烈的塑性变形，沉积形成低氧或无氧化涂层[2]。因此，传统热喷涂中常见的问题，例如高温氧化、蒸发、分解、残余拉应力等可得到减小或消除[3]。这使得冷喷涂可以制备各种各样的涂层[4,5,6,7,8]，包括绝大部分金属和合金 (Al,Cu,Ti,Mg,Zn,Fe,Ni,Ta,Bronzes,Ti-6Al-4V,不锈钢,MCrAlY等)、 陶瓷(TiO2,Zr2O3等)、复合涂层，例如金属-金属涂层(Al/Ti,Fe/Al,W/Cu,Al/Ni,Zn/Al,Cu/Cr,Al/Cu等)、金属-陶瓷涂层(Al/Al2O3,Al/SiC,Al/TiN,Ti/SiC,Ni/Al2O3,Ni/B4C,Cu/Al2O3等)和 金属-金属间化合物涂层(Al/Mg17Al12,Al/FeSiBNbCu,Mo/Terfenol - D,Fe/SmFe2等)，甚至还有金属陶瓷涂层(Cr3C2-NiCr,WC-Co等)、金属玻璃涂层(NiTiZrSiSn等)和纳米结构材料涂层(纳米晶Ni,Cu/MWCNTs,Fe(Al)/Al2O3等)。

虽然已通过试验和数值模拟的方法对冷喷涂粒子碰撞过程进行了研究[9,10,11,12,13]，但冷喷涂粒子的实际结合机理尚不清楚。目前，被普遍认可的观点是粒子局部严重的塑性变形导致粒子间的结合，如图 1所示，喷涂粒子和基体在碰撞过程中经历了严重的局部塑性变形，导致表面氧化物破裂，露出大量新鲜金属，使得粒子和基体金属间充分接触。洁净表面的接触以及大的接触压力是粒子/基体结合的必要条件。研究还表明，粒子沉积需要高于一定的“临界速度”，在局部剧烈变形区形成的金属射流有助于清除破碎的氧化膜[10, 13, 14, 15]。另外，Yin 等人[16]的计算结果表明，虽然粘性的金属射流有助于清除破碎的氧化膜，但由于金属射流只存在于接触区域的外围，导致会有部分氧化物残留，并在粒子沉积后在接触区域的中间积聚，基体材料硬度、粒子速度和喷涂角度均对残余氧化膜的最终形态有影响。此外，有关结合机理的一些重要问题仍尚不明确，例如，为什么沉积粒子的形状是现在这样的？在碰撞过程中界面处是否发生了反应？喷涂粒子的结合过程中塑性变形起什么作用？因此，有不少研究者采用数值模拟方法研究了冷喷涂粒子的高速碰撞过程，以期获得有用的信息。文中总结了对冷喷涂粒子的碰撞行为数值模拟研究以及临界速度的预测。

商用有限元模拟软件ABAQUS和LS-DYNA已广泛应用于冷喷涂粒子碰撞行为的研究。以前建立的模型，主要采用拉格朗日法[9, 11, 12, 13, 14, 20, 21, 22]。Assadi等人[9]首先用ABAQUS/Explicit显式算法，建立了二维轴对称和三维的碰撞模型，但其三维模型过于简化。Grujicic等人[13]采用ABAQUS软件，分析了单个粒子的碰撞过程。Yin等人[22]采用LS-DYNA同时研究了单粒子和多粒子的高速碰撞过程。以上均将碰撞过程假设为绝热过程，而没有考虑热传导的影响。此外，Bae等人[20]也应用ABAQUS软件模拟了高速粒子的碰撞过程，值得注意的是，Bae等进行了完全热力耦合分析[20]。

但是，在拉格朗日算法模型中，接触区域网格会发生严重的变形，导致计算结果在很大程度上取决于网格尺寸。Li等人[23]最早注意到网格尺寸的影响并提出采用任意拉格朗日欧拉法(Arbitrary Lagrangian Eulerian,ALE)来解决这一问题。Li等人[23]分别用拉格朗日法和任意拉格朗日欧拉法，系统地研究了Cu粒子碰撞Cu基体的变形行为。结果发现，拉格朗日模型的网格尺寸对界面区域的局部剪切失稳影响显著，而任意拉格朗日欧拉法则不受网格尺寸的影响，但该模型同样没有考虑热传导的影响。

Li等人[14]也采用拉格朗日法建立了三维模型来研究粒子的碰撞过程，但首次采用了材料损伤行为这一概念来解决网格畸变问题。由于正碰撞过程的对称性，采用1/4模型以缩短计算时间，如图 2所示。采用Johnson-Cook失效模型，界面接触设置为包括内部单元表面，以确保在单元破坏删除后新暴露单元仍能有效接触。值得注意的是，基于材料损伤模型的计算，会影响变形粒子的形状和物理量输出结果。

在一些早期的研究中[10,15]，采用了美国桑迪亚国家实验室开发的耦合热-流体力学(CTH)代码。该算法为两步欧拉法，建立的模型为二维轴对称几何模型，主要用于模拟多维、大变形及剧烈冲击等复杂问题，但同样被看作绝热过程。

ABAQUS和LS-DYNA同样有欧拉算法模块。Li等人[25]用ABAQUS/Explicit(Ver6.8)模拟了Cu粒子喷涂Cu基体的碰撞过程。由于该软件不支持二维空间的欧拉算法，因而将模型扩展为厚度为1 μm的欧拉空间，并约束了厚度方向的位移。模拟中假设粒子以垂直方向撞击基体，以便应用对称模型以缩短计算时间，如图 3所示，网格尺寸0.2 μm(网格分辨率1/100 dp)。该过程同样假定为绝热过程。新版本的ABAQUS(Ver6.12)更新了欧拉单元算法，考虑到热传导的影响，Li等人[26]也对其模型进行了更新。

除了欧拉法，无网格法中的光滑粒子法(Smoothed particle hydrodynamics,SPH)也是模拟冷喷涂碰撞过程的另一有效方法。光滑粒子法通过采用大量的粒子来模拟连续的物质，相比较欧拉法，光滑粒子法在预测复杂表面运动方面有强大的功能。Li等人[17,27]应用LS-DYNA软件中的光滑粒子法模块，研究了入射角对粒子变形行为的影响，所建立模型如图 4(a)所示。研究表明，纯光滑粒子模型和拉格朗日-光滑粒子耦合模型(如图 4(b)所示)均能合理的预测粒子变形行为。正如拉格朗日模型中网格尺寸影响输出结果，光滑粒子法模型中单个粒子的重量也对输出结果有一定影响。此外，光滑粒子法也可以模拟冷喷涂中多粒子的碰撞过程[27]。但是，纯光滑粒子法模型的计算量太大，这也促进了拉格朗日法和光滑粒子法的耦合，即：在大变形区采用光滑粒子法，其它区域采用拉格朗日法。

2.1 粒子变形形貌

Dykhuizen等人[15]最早采用CTH代码模拟了Cu粒子撞击不锈钢基体的过程，研究发现，高的碰撞速度造成了更深的冲击坑。Grujicic等人[10]也采用CTH代码研究了Cu、Al之间的碰撞行为，解释了为什么Cu粒子在Al基体上的沉积效率大于Al粒子在Cu基体上的 沉积效率。如图 5所示，这一差异可以归结为Cu的密度和强度比Al大，因而比Al有更大的动能，更大的结合面积。

为了定量表征粒子的变形程度，通常用扁平率作为指标，其定义为，扁平粒子直径D与粒子原始直径dp的比值(D/dp)。计算结果表明[9, 15, 28]，粒子的扁平率随粒子速度的增加而增大，这与试验结果基本吻合。Li等人[23]也将扁平率的计算值与试验值进行了对比，二者吻合良好。但是，在高速碰撞下显著的金属射流导致无法准确测定扁平粒子的直径，致使扁平率的计算变得困难。因此，他们首先提出了采用压缩率来表征粒子的变形程度[23]，即：

其中，hp是扁平粒子在碰撞方向上的高度。研究发现，压缩率比扁平率更容易估计，且不受网格尺寸和模拟方法的影响，因而比扁平率更适于表征粒子变形程度[23]。

2.2 碰撞界面温度

在冷喷涂粒子的高速碰撞过程中，必然导致界面处产生温升，这对界面结合至关重要[9, 15, 23, 28, 29]。Assadi等人[9]研究发现，整个接触界面上的应变和温度的分布不均匀，说明仅在部分界面处发生了结合。Dykhuizen等人[15]用CTH软件的计算结果也高估了界面的温度，因为没有考虑热传导。Li等人[23]用拉格朗日法模拟了Cu粒子碰撞Cu基体的过程，观察到狭窄的界面区域经历了剧烈的变形，在局部变形区温度明显上升。但是，如果考虑到块状粒子和基体间的热传导，温度峰值将会显著降低[24]。

Yin等人[22]研究了热软化Cu粒子在Cu基体上的沉积行为，基于模拟结果，提出了热软化区的新定义，即：在温度超过0.5Tm的区域，可能发生热软化现象。研究发现，对于未预热或低温预热粒子，热软化区主要位于狭窄的界面区；如果预热温度显著提高，则热软化区会扩展至整个变形粒子，导致整个粒子都可能经历热软化。此外，其研究也表明，随着粒子初始温度提高，局部最高温度逐渐上升；对于高温粒子，能观察到明显的金属射流，这也与试验观察结果相符[22]。

2.3 能量转化

众所周知，粒子沉积过程是一个能量转化过程，粒子的动能大部分转化为内能[30]。因而可以认为，由于Al的密度较低，使得Al粒子的动能较小，最终导致变形程度小[30]。粒子的初始动能EU转化为四部分能量：塑性耗散能EP、粘性耗散能EV、摩擦耗散能EF以及弹性应变能ER。粒子动能的耗散可表示为[31]：

在Bae等人[20]的研究中，认为塑性耗散能EP、粘性耗散能EV和摩擦耗散能EF是结合能，而弹性应变能ER是回弹能。由于EV和EF相对于EP很小，可以忽略不计，Bae等人[20]用塑性耗散能EP和弹性应变能ER相互比较，以估计粒子和基体间是结合或反弹。

巫湘坤等人[32]用LS-DYNA研究了高速粒子与基体碰撞时的能量转化过程，以及粒子动能在粒子与基体材料之间的分配情况。研究表明，碰撞过程中粒子的能量不断转移到基体材料中，粒子能量不断减少，而基体材料能量不断增加。定义碰撞之后分配到粒子的能量Eparticle和分配到基体的能量Esubstrate的比值为能量分配系数K[32]，即：

K越大，说明碰撞过程中粒子中分配的能量越高，基体获取的能量则越低，说明粒子塑性变形消耗的动能比例大，变形程度剧烈。对于同种材料粒子/基体的组合，随着粒子动能增加，变形加剧，K值不断降低，有利于喷涂粒子结合，粒子初始动能及粒子/基体材料的变形能力均影响K值大小，K值高低表征粒子与基体变形程度差异，临界速度碰撞时K值在0.4~0.6之间。 3 临界速度预测 3.1 拉格朗日法下的临界速度预测

临界速度是冷喷涂中的重要变量。在早期的试验中已进行了研究[1, 2, 15]。但实际上影响临界速度的因素很复杂，主要包括材料的力学性能、粒子初始温度、粒子形状、尺寸及氧化状态[33]。

因此，采用数值模拟法来研究冷喷涂粒子的临界速度。Assadi等人[9]首先提出把发生绝热剪切失稳作为粒子是否成功结合的判据，Assadi等[9]模拟了Cu粒子碰撞Cu基体的过程，如图 6所示，接触界面处应变和温度剧烈变化，在550~580 m/s这一狭窄的速度范围内，金属射流从塑性流转变为粘性流，表明发生了绝热剪切失稳。依据材料的物理和力学性能，总结出了几种金属材料的临界速度方程。其后，Grujicic等人[13]同样研究了接触界面处的绝热剪切带的演化过程，分析了喷涂粒子下表面和基体上表面的某些单元格中等效塑性应变率、等效塑性应变、温度和等效正应力随时间的变化情况。在接触界面处产生剪切带的最小速度如表 1所示，作为对比，Assadi等人[9]报道的临界速度也加添加到表中。总的来说，两数值间联系紧密，这也表明在粒子沉积过程中，剪切带确实起到至关重要的作用。

基于Assadi等人[9]的研究结果，Li等[24]用ABAQUS软件计算了不同网格尺寸下单元格畸变的初始速度，结果如图 7所示。作为对比，利用LS-DYNA软件[23]针对Cu的模拟结果也绘制到了图中。如果把网格单元畸变的初始速度作为临界速度，则对同一材料而言，临界速度随网格尺寸的减小而降低。研究同时发现，材料的临界速度大都在300~400 m/s之间，即使是高强度的In718合金，也远低于所报道的临界速度[13,34]。根据先前 的研究[35]，Cu的临界速度大约为300 m/s，且其临界速度与粒子的氧化程度紧密相关。Kang等人[39]随后报道了关于Al的相似的结论，并证实了Li的模型。因此，Assadi[9]和Grujicic[13]关于临界速度的计算结果与其试验结果相符只是一种巧合，如图 8所示[37]，作为对比，文献中关于Cu[38]和316L不锈钢[34]的相关数据也绘制到了图中。因此，预测的临界速度应当是没有考虑氧化膜的情况下的值。上述结果也表明，在拉格朗日法下，将网格尺寸外推至接近零，可能会得到临界速度。

采用单元删除法，Li等[14]计算了不同速度下不同尺寸的粒子的临界速度。假定粒子表面第一层单元开始被删除时的速度为临界速度，结果如图 9所示。作为对比，文献报道的试验值[34,35]也绘制到了图中。可以看出，随着粒子尺寸大幅度的增加(20 μm~20 mm)，临界速度缓慢减小，这表明粒子尺寸对临界速度的影响很小。小尺寸Cu粒子预测的临界速度的变化趋势与试验结果[34]相符的很好，但比Schmidt等人[34]报道的试验结果低。根据先前的研究[35]，含微量氧的Cu粒子的临界速度大约为300 m/s，Schmidt等人[34]的研究表明直径20 mm的Cu球的结合速度大约为250 m/s。正如先前的研究所提出的[35]，粒子氧化状态是导致不同研究者试验结果不同的主要原因，因此，该模型计算的临界速度是合理的，因为同样未考虑氧化膜的作用。

Schmidt等[34]研究了粒子尺寸对粒子沉积行为的影响，试验测得的临界速度如图 10所示。结果表明，对Cu和316L不锈钢而言，其临界速度均随粒子尺寸的增大而减小，但Li[23]的研究表明，这实际上是由于粒子被氧化所致，在速度为500 m/s的情况下，20 μm和80 μm的喷涂粒子等效塑性应变(PEEQ)和能量变化结果表明不同尺寸粒子的局部变形情况相似[23]。这些结果进一步证实了粒子尺寸对临界速度无明显影响。

Li等建立了欧拉法模型，并提出了相应的预测临界速度的方法[25, 39, 40]。他们采用欧拉法模拟了直径20 mm的Cu球撞击低碳钢基体的碰撞过程，结果与Schmidt等人[34]的研究结果相符。因此，已有的欧拉模型可用于研究冷喷涂粒子的碰撞行为。

Yu等人[39]计算了绝热条件下，喷涂粒子以不同速度撞击基体的过程，等效塑性应变云图如图 11所示。可以看出，随着碰撞速度提高，粒子变形加剧，坑深和扁平率均增大。在相对低的碰撞速度下(200~280 m/s)，最大塑性应变区被粒子和基体材料包围，如图 11(a~c)；当碰撞速度上升至大约290 m/s时，高应变区被挤出，形成射流，如图 11(d)；随后，射流随碰撞速度的增大而伸长，图 11(e~i)；但是，碰撞速度进一步提高，射流不再连续，而是在射流前端材料发生了断裂，如图 11(j~l)。

碰撞速度对PEEQ最大稳定值的影响[39]如图 12所示。在速度为290~400 m/s时，PEEQ最大稳定值几乎是水平的，这表明塑性应变强化和热软化之间达到了动态平衡。 因此，Yu等人[39]把转变速度(Transition velocity)当作临界速度，在该模型中，即为290 m/s。这一速度与图 11中的高应变区开始以塑性状态挤出形成射流的速度相符。基于这些研究结果，提出了利用欧拉法模拟的PEEQ最大稳定值预测临界速度的方法。

如前所述，欧拉法[27,39]忽略了热传导的影响，而这会降低计算结果的精确度。因此，Wang等人[40]将热传导考虑在内，建立了耦合热力学-欧拉模型(CTM-Eulerian)，对比了绝热和非绝热条件下粒子的碰撞行为。Wang等[40]模拟了直径20 μm的Cu粒子碰撞Cu基体的过程，等效塑性应变和温度的分布云图如图 13所示。可以看出，热传导强烈的影响着温度分布，但对粒子/基体间的变形则无明显影响。

Wang等人[40]同时研究了在上述条件下最高温度和PEEQ最大值的变化情况，如图 14所示。可以看出，绝热和非绝热条件下的等效塑性应变(PEEQ)变化趋势相似，但绝热条件下的最高温度变化的更快，温度也更高。这是由于在绝热条件下忽略了接触界面处热传导所造成的热损耗。这也表明，热传导会对临界速度的预测产生显著影响[41]。

Cu粒子以不同速度碰撞Cu基体时PEEQ最大值随接触时间的演变情况见文献[42]。在不同粒子速度下，PEEQ最大值均快速达到各自稳定值，表明PEEQ最大值与粒子碰撞速度密切相关。

图 15为绝热[39]和非绝热[40]条件下碰撞速度对PEEQ最大值的影响。在非绝热条件下，速度为300 m/s时PEEQ最大值达到峰值，表明塑性应变强化和热软化间达到动态平衡。但值得注意的是，非绝热条件下的临界速度比绝热条件下的(290 m/s)要稍大，这表明热传导削弱了接触界面上的热软化作用，故需要更大的速度以产生绝热剪切失稳。因此，由耦合热力学-欧拉模型计算得出的300 m/s比绝热欧拉模型得出的结果更合理。

基于以上可靠性分析，耦合热力学-欧拉模型(CTM-Eulerian)已用于模拟其它常见的喷涂材料，如Fe、Ni、SS304不锈钢、Al、In718高温合金和TC4。图 16为不同材料的PEEQ最大值随速度的变化趋势。类似于Cu，所有的材料都出现了等效塑性应变平台，这是达到临界速度的重要标志。显然，由于材料性能的差异，转变速度各不相同。预测的Fe、Ni、SS304不锈钢、Al、In718高温合金和TC4的临界速度分别为350、380、395、410、490和500 m/s，一些值似乎与试验值[34]比较吻合。

上述研究表明，对冷喷涂粒子碰撞行为的研究还在不断发展。不仅模拟了碰撞过程和接触界面处的反应，也发展了预测临界速度的更加有效的方法。但是仍有一些重要问题有待通过数值模拟进行详细研究，主要包括：

(1) 对临界速度预测结果验证的物理方法及其基础理论；

(2) 氧化膜的影响机理；

(3) 在原子尺度范围结合的本质；

(4) 涂层的组织演变和再结晶；

(5) 残余应力；

(6) 涂层(尤其是复合涂层)的建立过程；

(7) 涂层性能的预测。