The critical velocity for adhesion of a microparticle on impact with a substrate is formulated and internal dynamics that lead to critical points are graphically illustrated. The modelling involves energy balance; kinetic, plastic and heat. Analysis of interfacial high strain-rate, the plastic equivalent strain and duration of plastic deformation brings out the internal material dynamics that lead to the critical point for material jetting at the particle substrate interface. Through non-dimensional analysis, the critical velocity is a point where plastic ﬂow velocities, exceed acoustic ﬂow velocities. This illustration may underpin underlying thermal transport and vibrational signatures in microparticle impacts. Keywords: ballistic impact-adhesion, structural impact-excitation, heat transport, high strain rate, critical velocity, erosive wear, wave propagation, acoustics 1 Introduction When a microparticle hits a target surface, it may rebound, adhere or erode from the surface. Ballistic impacts of microparticles are a key manufacturing process in short peening (Al-Obaid, 1995) and additive methods (Wang et al., 2023). Supersonic impacts of dust storms involving microparticles must be miti- gated in aerospace vehicle structures to minimise changes in the dynamics and structural properties such as surface friction and surface roughness and eventual structural damage (Palmieri and Lambert, 2001); the need for mitigation of wear damage. Supersonic impact bonding only occurs at a certain velocity range which has been termed window of deposition (Schmidt et al., 2006). This velocity range is bounded by the critical velocity (minimum velocity) and the erosion velocity (maximum velocity) (Alkhimov et al., 1994; Moridi et al., 2013; Hassani-Gangaraj et al., 2018a; Li et al., 2006). Below the critical velocity the plastic deformation is too low to cause bonding whereas above the erosion velocity hydrodynamic penetration leads to strong erosion. Between these two characteristic velocities optimum conditions for deposition are reached (Schmidt et al., 2006). Both the critical and erosion velocity are characterized by spray material, powder properties, particle and substrate temperature (Wang et al., 2023; Chaban et al., 2021) and particle size (Dowding et al., 2020). Due to the small size of the particles, the interface physics maybe governed by microscopic heat transport, there is therefore a need for a careful study so as to inform the compu- tational modelling approach including material constitutive methods Rahmati and Ghaei (2014) and the experimental interrogation methods. 2 Methods 2.1 Solution procedures Analysis of the solution procedures were based on the implementation details in the Abaqus solver as detailed in the user manual in (Simulia, 2012). Computational models used are shown in Figure 1. The material parameters for this study are listed in Table 1. The details of the numerical modelling were ﬁrst given in a doctoral thesis (Moonga, 2019) by the ﬁrst author, some parts repeated here to make this paper 1 Figure 1: Computational model (a) 20 µm particle diameter copper, (b) 10 µm particle diameter copper on copper substrate with inﬁnity boundaries on the sides and bottom surfaces. self-contained. The general heat equation solved, equation 1, is an ∂ ∂x k∂T ∂x +∂ ∂y k∂T ∂y +∂ ∂z k∂T ∂z +Q=ρc∂T ∂t(1) energy conservation equation that applies Fourier’s law qx=−k∂T ∂x,qy=−k∂T ∂yandqz=−k∂T ∂zas the heat conduction components in the x,y, and zdirections and kis the thermal conductivity. Qis internal heat generation per unit volume, ρis the material density, cthe speciﬁc heat capacity of the material,∂T ∂t is the rate of change of temperature,∂T ∂zis the directional temperature gradient. In adiabatic analysis the increase in temperature is calculated directly at the material integration points based on adiabatic thermal energy increase caused by inelastic deformation (Simulia, 2012). Adiabatic anal- ysis assumes there is no heat conduction during the material deformation. This assumption can be assumed to hold in short durations of high speed manufacturing processes involving large inelastic strains. When both inelastic generation and heat conduction are required, a fully coupled temperature-displacement analysis must be performed using fully coupled thermal-stress analysis (Simulia, 2012). Plastic straining gives rise to a heat ﬂux per unit volume of Qpl=βσydεpl dt=βσ˙εpl(2) where Qplis the heat ﬂux (generation) that is added into the thermal energy balance, βis the user-speciﬁed material dependent inelastic heat generation fraction, σis the stress, and ˙εplis the rate of plastic straining. In an adiabatic analysis the heat equation solved at each integration point is therefore given as ρc∂T ∂t=βσydεpl dt(3) where ρis the material density and Cis the speciﬁc heat which could also be temperature dependent. Equation 3 is directly obtained from equation 1 where the heat conduction components qx, qyandqzare zero; and Qis deﬁned by equation 2. Whereas when both inelastic generation and heat conduction are required the heat equation solved will be given by equation 1. ∂ ∂x k∂T ∂x +∂ ∂y k∂T ∂y +∂ ∂z k∂T ∂z +βσydεpl dt=ρc∂T ∂t(4) 2.2 Mathematical modeling of critical points All the mathematical modeling were based on the analysis of energy balance with interest at the particle substrate interface around the critical and erosion velocity. 2 Material Cu Al Ni Ti Mild Steel Zn Density ρp,Kgm−38960 2710 8890 4510 7870 7140 6660 Thermal conductivity K 386 220 Speciﬁc Heat Cp,J(KgoC)−1383 920 456 528 481 Heat of Fusion J(kg)−1207 000 321 000 Melting Point K 1356 916 1726 1923 1793 Elastic Modulus GPa 124 68.9 207 116 200.1 108 Poisson’s ratio 0.34 0.3 0.31 0.34 0.30 0.25 Bulk modulus GPa 140 71.7 182 110 140 70 Bulk speed of sound Co,ms−13952.8 5143.7 4524.6 4938.6 4217.7 3131.1 3256.7 J C Plasticity A, B,MPa, 90, 292 148 n, C, m 0.31,0.025,1.09 Table 1: Properties of copper (Johnson and Cook, 1985), Aluminum (Gupta et al., 2006), Nickel (Razavipour and Jodoin, 2021), Titanium, Mild Steel (Bae et al., 2008) 3 Results and Discussion 3.1 Plastic strain rate in the substrate interfacial material In order to model this problem, let γsbe the energy transferred to the substrate and used in plastic defor- mation from initial energy of the impacting particle. Then from the kinetic energy of the particle the plastic deformation energy will be given by γs×1 2ρpVnpvp2=σysdεp dt.tc Vpls (5) Isolation of impact velocity yields equation 6 vp=2×σys×˙εs.tc×Vpls γs×ρpVnp1/2 (6) Since around the critical velocity, tc≈dp/vp, see for instance in Hassani-Gangaraj et al. (2018b) and Maev and Leshchynsky (2009), hence equation 6 can be written as vp≈2×σys×( ˙ε.dp)×Vpls γs×ρpVnpvp1/2 (7) Letµs= 2×Vpls/((Vnp)), then equation 7 can be given as vcr≈σys×( ˙εs.dp)µs) γsρp1/3 (8) deﬁning the critical velocity for adhesion. Previous studies have pointed out the relation of critical speed of adhesion to the bulk speed of sound in the material. At the critical velocity of adhesion, ˙εs.dpµsshould relate to the bulk speed of sound Cs, andvpis equal to vcr. It is important to take note that here the variables used are the nominal variables involving an undeformed particle volume. The variations in the density and the elasticity of the particle through Poisson’s ratio due to elasticity and deformation may determine the exact location. To further substantiate this claim, equation 8 is further analysed in Figure 2 using the non-dimensional values which gives a characterization of the nano scale structural dynamics. The lower bound is deﬁned by the response variable Co/vpand the upper bound for erosion by (Co/vp)2, severity of impact is given by (γsρpv2 p)/(σys) where γsis given as 1/2 as an approximation ( γsat critical was numerically observed to be higher than γp). The assertion that ˙εs.dpµsrelates to Csat the critical velocity, brings a new discourse to a consideration of the internal structural dynamics during this process. The results of which may inform experimental 3 Figure 2: Critical and erosion velocity equations for adhesion of microparticles on a substrate surface. interrogation methods for the transient impact dynamic state. On the plasticity of the material, since the plastic equivalent strain, can be estimated from ((γsρpvp2))/σys= 2×(εpl)(Vpls/Vnp), it seems for some materials the spall strength of the substrate could be estimated from σys/γsa similar argument can be made about the particle. The limitation is that σysused in this analysis is a quasi-static value of the material plasticity at room temperature. Thermal softening would lower the dynamic yielding to zero in the material model at melting, see the Johnson cook equation Johnson and Cook (1983). A comparison on these values can be made between various materials for instance between Cu and Al to see how this value, σys/γs, would relate to the spall strength. In a numerical analysis, data was obtained from the impact simulations of a 10 µm diameter copper particle on copper substrate with material properties listed in Table 1. For the impact velocity of 412.6ms−1, the maximum plastic equivalent strain was recorded, the plastic energy in the particle, and the plastic energy in the substrate. The parameters estimated in the computational analysis are shown in equation 9. vp=2×σys ρpεPl γsVpls Vnp 1/2 (9) From the energy balance γswas0.564andγpwas0.291, where γswas the ratio of plastic energy in the substrate to the kinetic energy of the particle and γpwas the ratio of plastic energy in the particle to the kinetic energy of the particle. The total plastic energy αpl=γs+γpis the ratio amount of kinetic energy used in the plastic straining of both the particle and the substrate. The volume ratio Vpls/Vnpwas1.463, plastic equivalent strain εplwas3.263plastic energy in substrate 2.25x10−7J, plastic energy in the particle 1.16x10−7J. This data was all based on numerical simulations. The plastic aﬀected volume in the substrate was calculated from the following equation: Epls=σys(εpl)Vpls, where Eplsis the numerically calculated plastic energy in the substrate, εplwas set to the maximum plastic equivalent strain. 3.2 Interfacial temperature increase in a single particle impact Now we move from the plastic deformation to the heat generation. For heat generation rate in the sub- strate interfacial material an additional parameter β= 0.9heat conversion parameter is added. Below the melting temperature, the change in temperature at the interface ∆Tin the substrate can be accounted for as ∆T=βεpsσysVpls ρcpsVths(10) 4 Figure 3: (a) shows the temporal evolution of temperature, stress and strain and (b) typical temperature values from the simulations at 500 ms−1, of a 20 µm Cu particle impacting a Cu substrate in a non-adiabatic analysis. The image in (c) gives the plastic equivalent strain and (d) typical temperature values from the sim- ulations at 412.6 ms−1impact of a 10 µm diameter Cu particle impacting a Cu substrate in a non-adiabatic analysis (all temperatures in degree Celsius, initial temperature was 25 degrees Celsius). Where Vplsis the plastic aﬀected volume, and Vthsthermal heating aﬀected volume and other parameters as previously deﬁned. The plastic strain εpsis velocity dependent and is also dependent of the assumption of adiabatic or non-adiabatic heat transfer. See equation 9 for the details on the plastic equivalent strain as detailed in the previous section. The yield stress used in equation 9 is obtained from the Johnson cook model material parameters, for instance 90MPa for copper, Table 1. With adiabatic modelling, the change in temperature at the critical velocity is much higher than for non-adiabatic analysis. In this section ex- amples of temperature evolution and plastic equivalent strains studied for the non-adiabatic solution are shown in Figure 3. From the simulations at 412.6ms−1impact of a 10µmdiameter Cu particle impacting a Cu substrate in a non-adiabatic analysis the maximum plastic equivalent strain attained was 3.263 at 20.5ns after impact; the maximum temperature was 159 degrees Celsius at 7.33nsafter impact. The other issue is that of thermal resistance at the interface, this value is normally not known in continuum computational studies. In Figure 3 interfacial jetting as seen induces tensile stresses, an analysis of these are reported in (Moonga and Jen, 2017). 3.2.1 Interfacial Melting and Erosion velocity in relation to the critical velocity A similar approach as from the previous sections is used here. Let γsbe the ratio of the kinetic energy used in plastic deformation of the substrate. Then from the initial kinetic energy of the particle β×γs×1/2ρpVnpv2 e=ρsVths[Cps(Tms−To) +Hfs] (11) Hereβis the tailor quinney constant which accounts for amount of plastic strain converted to heat, γsis the amount of plastic strain energy in the substrate with reference to the initial kinetic energy of the particle. Hfsis the heat of fusion, Tois the room temperature, Tmthe melting temperature, Cpthe speciﬁc heat capacity, Vthis the aﬀected volume in the substrate around the interface, ρpis the density of the particle material, ρsis the density of the substrate material, dpis the particle diameter, veis the impact velocity at onset of erosion; Vthis the thermal heating aﬀected volume in the substrate and can be written as Vth=AcdthhereAcis the contact area between the substrate and the deformed particle; dthis the thermal depth. From equation 11 ve=2ρsVths[Cps(Tms−To) +Hfs] γsβρpVnp1/2 (12) 5 For an approximate solution, here the contact area and thermal depth have previously been estimated in (Hassani-Gangaraj et al., 2018b) as Ac=πd2 p/2,dth= (ks/(ρsCs)×tc)1/2where ksis the thermal conductivity, tcis the time of contact between the substrate and the particle. For a numerical example, assuming Cu; with dp= 15 micron, Ac= 3.534×10−10m2,dth= 1.341µm; letγs= 0.5,β= 0.9, and Hfs= 207000 J/kg . The erosion velocity can now be calculated considering the nonzero heat of fusion as giving an erosion velocity veof854.2ms−1with similar results as obtained by (Hassani-Gangaraj et al., 2018b). In the computational studies we conducted in trying to understand the thermal model, strain value has been observed to depend on adiabatic or non-adiabatic conditions and initial temperatures. Equation 6 through to equation 12 and dynamic relations given in Figure 2 can be summarised, in the ﬁrst part, for onset of erosion to be governed by the equation 13 γseρpv2 e σysVnp Vths=bCo ve2Vpls Vths = 2×εplVpls Vths =ρs[Cps(Tms−To) +H(fs)] βσys(13) and in the second part, onset of adhesion to be governed by equation 14 where Tcrsis the substrate tem- perature at the critical velocity. γscrρpv2 cr σysVnp Vths =aCo vcrVpls Vths = (2×εpl)Vpls Vths =ρs[Cps(Tcrs−To)] βσys. (14) Equation 13 and equation 14 can be considered to be the nominal equations where the underlying dynamics that govern rebound, adhesion and erosion can be traced. It should be interesting to see how this result re- lates to higher speed shock velocities such as these (Stewart et al., 2023; Ma et al., 2020; Guo et al., 2024). The development of experimental in-situ observations of micro particle impact velocity and rebounding velocity reported in (Hassani-Gangaraj et al., 2018b; Reiser and Schuh, 2023) provided more insight on this problem. In-situ observation of interfacial temperature evolutions at these length scales and time scales has been diﬃcult (Chaban et al., 2021), as this would require ultrafast thermal measurement processes. An- other area that would make signiﬁcant contribution is that of the acoustics and vibrational measurements; but has not yet been explored for this problem. Plastic strain measurements in ultrafast experimental pro- cesses are normally used to determine the dynamic yield stress, and these normally use larger specimens, and the interest is normally not the transient dynamics of the deformation process leading to adhesion or rebound. Whilst the idea of interfacial jetting by hydrodynamic pressure had been proposed previously (Hassani-Gangaraj et al., 2018b), the relations given in the forgoing equations point to the fundamental physics. 4 Proposed future directions 4.1 Microscopic heat transfer The need for increased understanding of the energy transport mechanisms at small length scales gave rise to the ﬁeld of microscale heat transfer which involves study of individual heat carriers (Bejan and Kraus, 2003). Microscale heat transfer can be considered when the continuum model is nolonger valid at small length scales or small time scales compared to the intrinsic properties of the heat carriers. In the classical continuum model the heat transfer analysis applies the conservation of energy equation 3 coupled with Fourier’s law for thermal conduction in the solution procedures. The ﬁrst attempt at microscale heat trans- fer modeling is to modify the continuum model in such a way to account for microscale thermal transport relations. One common approach, though diﬃcult, applies the Boltzmann transport equation. When both these methods fail, molecular dynamics can be adopted. Microscopic heat transfer analysis requires an- alyzing microscopic energy carriers and the basic heat transfer mechanism (Bejan and Kraus, 2003). The basic thermal carriers in metals are the free electrons, while in semiconductors and insulators are the lattice vibrations (phonons) which travel like acoustic waves. Theoretical descriptions of electrons and phonons are based on the assumption that the material is of crystalline structure (Bejan and Kraus, 2003). In microscopic analysis, Catteneo’s equation 15 can be used to account for ﬁnite speed of thermal energy propagation (Bejan and Kraus, 2003). 1 τ∂q ∂t+q=−K ∆T(15) 6 where τis the relaxation time of the heat carrier. Substitution into the conservation equation gives a hy- perbolic diﬀerential equation. Equation 15 reduces to Fourier’s law when relaxation time is much less than the time scales of interest. At present it’s not yet clear whether during high speed impacts the interfacial temperatures are properly deﬁned using adiabatic analysis or using non adiabatic analysis where heat ﬂow through conduction is accounted using Fourier’s law. It could be that a better estimate of the heat trans- port could be made through equation 15. This equation would account for the time delay for the heat ﬂux to develop after temperature gradient has been applied. 4.2 Vibrations and acoustics Micro and nano particle vibrational signatures (Tang et al., 2023); and wave propagation mechanism (Kar- man and Duwez, 1950) in the particle will deﬁnitely inform on the underlying dynamics pointed to by the equations and methods deﬁned in this report, especially in relation to high strain rate deformation. 5 Conclusion There is an observed relationship between plastic strain rates and acoustics in the impacting particle by illustrative underlying dynamics that may lead to rebound, bonding or erosive wear. The window of depo- sition was thus deﬁned by the critical and erosion velocities deﬁned by continuum material parameters. Varying these parameters would therefore vary the location of the critical and erosion velocities, shifting the erosion or deposition window. It should be possible also to vary the width of the deposition window. At this point it was not conclusive of the heat transport across the particle substrate interface, this requires some further examination. The heat generation and thermal analysis at the interface should be able to provide a connection to plastic strains, plastic strain rates and energy distribution after impact. The levels of temperatures attained will deﬁnitely provide more details to the adhesion energy, and the strength of the bonding. This analysis should be able to give better connection to the underlying dynamics. Bonding is an interfacial process whilst rebound is based on the bulk of the particles recoverable energy in inter- action with the substrates recoverable energy. Adhesion delay has been observed experimentally; where the ﬁrst particles rebound however once the surface chemistry has changed bonding takes oﬀ, indicating the importance of interfacial physics in bonding. This study advances the understanding on the high strain deformation during micro particle impacts.