Multiphysics Simulation and Experimental Investigation of the Densification of Metals by Spark Plasma Sintering

Multivariant experimental investigations and multiphysics microstructural modeling of the spark plasma sintering process of metallic powders have been performed up to a relative density of approximately 80%. In comparison, the effect of sintering temperature, pressure, and particle size on the interparticle contact area growth and axial shrinkage of cylindrical specimens of copper and nickel particles is measured in laboratory scaled tests. Herein, for the first time all relevant for sintering phenomena are considered simultaneously: the fully coupled thermo‐electro‐mechanical modeling of the spark plasma sintering processes, additionally taking into account for lattice, grain boundary, surface diffusion, electromigration, and thermomigration, has been carried out. The computational analysis of various physical phenomena allows to identify dominant and insignificant mechanisms. The two‐level numerical simulation includes the modeling of the sintering setup at the macroscopic level and the neck formation process in particle chain systems at the microscopic level. The results of the numerical simulations show a very good agreement with the experimental data. Therefore, the impact of electrical and mechanical loads as well as of particle size on microscopic distribution of temperature, inelastic strain, and on densification has been studied by the finite element simulations.

Numerical methods are usually applied for the investigation of sintering processes within the sintering tool or for the whole machine comprising the powder sample(s).The computational technology mostly applied is based on the finite element method.[41][42][43][44] The more complex problem, taking into account both coupled diffusion and mechanics, in which the analysis of the sintering process is carried out at the microscopic level considering the inhomogeneity of the stress state, was investigated in refs.[45,46].As sintering process is a complicated issue, it is important to apply multiphysics modeling and simulations of powder systems for various technological processes. [47]he aim of the present work is to evaluate dominant mechanisms of sintering as well as to analyze an inhomogeneity of temperature, electric field, displacement, stress, plastic, and creep strains by finite element simulations.The validity of the simulations is shown by comparing results of numerical modeling of SPS with experimental results.
The article is structured as follows: after description of the experimental procedure in Section 2, the coupled thermoelectro-mechanical modeling and simulation accounting surface, grain boundary and volume diffusion, electromigration, and thermomigration are performed in Section 3 with detailed analysis of influence of pressure, temperature, and different sintering mechanisms on the neck radius evolution.The comparison between experiment and numerical simulation is presented in Section 4.

Experimental Section
For the experiments, air atomized copper powder provided by Ecka Granules and nickel powder from Hoeganaes have been used.The copper powder was sieved to obtain particle sizes of 20-25, 45-50, and 100-106 μm.The nickel powder had sizes of 15-35 μm.Hence, it was only sieved to obtain particle sizes of 20-25 μm.The powders were reduced by heating to 450 °C for 1 h in a hydrogen atmosphere.To retain only spherical particles, the 100-106 μm copper powder was put onto a slightly inclined glass plate.Particles rolling down were used for experiments.The other (nonspherical) particles were discarded.
The powders were subsequently stored inside the argon-filled glove box that was later employed to perform the experiments in an inert atmosphere.
The experiments were performed in the custom-made setup depicted schematically in Figure 1 (see also ref. [48]).Punches from a molybdenum alloy (TZM) of approximately 2 mm diameter were located at the top and bottom of the tool.A fused silica tube served as current insulating die and allowed for the measurement of the sample temperature by a Keller PZ 20 AF 3 pyrometer.An ATP K-2528 force sensor attached to the upper punch was used to monitor the applied force (and implicitly the pressure).Graphite foil discs were used as separator between the powder and the punches on the bottom and top of the specimen.Springs attached to the setup were employed to set the desired pressure prior to the experiment.It was not possible to adjust the pressure throughout the SPS process.A customized Knürr-Heinzinger PCY 18-400 power source attached to the punches was used to generate a pulsed current with 10 ms pulses separated by 5 ms or 140 ms pauses; in the following these pulse patterns will be denoted as 10 ms þ 5 ms and 10 ms þ 140 ms.Note that for 5 ms pauses the current actually never fell to 0 A. The sample was heated by the pulsed current to the desired process temperature and the effective current was adjusted to maintain the process temperature.The actual temperature and pressure were recorded throughout the experiment by a LabJack T-7 Pro USB A/D converter.A Pico Technology 2206B USB oscilloscope recorded the current and voltage across the sample with an attached Pico Technology TA 375 100 MHz current clamp and the second oscilloscope channel attached to the TZM punches.By monitoring the sample with an USB microscope, the setup worked as a simple dilatometer.
Figure 1b shows a photograph of the glowing sample during an experiment.If long pauses between the pulses were employed, a substantially higher current was necessary to maintain the sample temperature.Hence, the graphite foils glowed brightly during the pulses as shown in Figure 1c and confirmed their role as primary source of heat.
A typical progress of temperature (red curve), peak current (orange), effective current (blue), and minimal current (dark cyan) throughout an experiment with a sintering temperature of 600 °C, a pressure of 5 MPa, and 10 ms þ 5 ms pulses is shown in Figure 1d.The corresponding effective current density and pressure are shown in Figure 1e.
After the desired process time the current was switched off and the samples were removed from the fused silica die.To measure the progress of sintering, the samples were fractured in the center after immersion in liquid nitrogen.Scanning electron microscopy (SEM) was used to image the fractured specimens as shown in Figure 1f,g.
Results of the measurements of the neck radii in the SEM pictures for different loading conditions are given in Table 1.The range of variation in neck radii (contact zones) for different samples under the same loading conditions is shown in Table 1 using the symbol AE.

Diffusion Mechanisms in the Sintering Model
The main sintering mechanisms and the used equations for diffusion are given in the present section.SPS as densification process represents combination of various diffusion mechanisms known from pressure-less sintering [49] (involving, e.g., surface and volume (lattice) diffusion), and from a pressure-assisted sintering process (grain boundary diffusion, plasticity, and creep) as well as potentially other effects like electromigration and thermomigration.
In general, the flux vector j due to surface diffusion, grain boundary diffusion, or volume diffusion is defined by the scaled gradient of the differences of the chemical potentials of atoms μ a and of vacancies μ v [50] Note that the diffusion flux vector is given in atoms per unit area and unit time.D ¼ D 0 e À Q RT is the effective diffusion coefficient, where D 0 is the maximal diffusion coefficient (at infinite temperature) of the material; k = 1.38065Â 10 À23 J K À1 is the Boltzmann constant, R = 8.31446 JK À1 mol À1 is the ideal gas constant, T is the absolute temperature, and Q is the activation energy.The parameter n a represents the number of atoms per unit volume (n a ¼ 1=Ω, where Ω is the atomic volume) for volume diffusion, and is the number of atoms per area of interface (n a ¼ δ=Ω, where δ is the effective thickness of the interface) for free surface diffusion and grain boundary diffusion.
In the present research, multiple sintering mechanisms of material transport are taken into account: surface diffusion, grain boundary diffusion, volume (lattice) diffusion, electromigration, and thermomigration.The various types of diffusion have different mass transfer paths.Their respective driving forces-given by various chemical potentials-are listed below.
The resulting atomic flux density vector and displacement rates (velocity vector) are determined by the contributions of all relevant transport mechanisms (five driving forces and three mass transport routes) In Equation ( 2) and (3), the upper index of flux j denotes the driving force (surface curvature κ, mechanical stresses σ, gradient of the temperature T, electric field E, gradient of the concentration c) and the lower index for the transport path (surface s, grain boundary gb, volume (lattice) v).The used equations for the chemical potentials, flux density vectors, and displacement rates for surface diffusion are taken from refs.[50,51], for volume diffusion from refs.[45,52], for grain boundary diffusion from refs.[53,54], and summarized in Table 2.
Here, γ s is the specific surface free energy of a crystalline solid; κ ¼ 1=R I þ 1=R II is the mean curvature, where R I and R II are the two principal radii of surface curvature (the curvature is positive for a convex pore); p ¼ À1=3trσ is the hydrostatic pressure and μ a0 is the chemical potential of the atoms in a stress-free state.σ nn ¼ n ⋅ σ ⋅ n denotes the normal stress on the grain boundary with the normal vector n.For the stress tensor σ, it is assumed that positive values correspond to tensile stresses.The differential operator ∇ s ¼ ð1 À nnÞ ⋅ ∇ is the surface gradient, 1 denotes the unit tensor, and ∇ gb ¼ ð1 À nnÞ ⋅ ∇ is a surface gradient calculated along the grain boundary.In Equation (17), (19), and (20), f is the fraction of the vacancy volume to the atomic volume Ω.The normal velocity v σ gb is positive when the grains move apart.Q v is the enthalpy difference between vacancy movement and formation, Z is the effective valence of a migrating ion, and e is the elementary charge.
Surface diffusion causes changes in the shape of the particles through rearrangement of matter that leads to interparticle neck growth, without changing the distance between the centers of the particles (shrinkage).Grain boundary diffusion is accompanied by the removal of material from the boundary region that leads to shrinkage (densification) during sintering process.
The mass conservation law results in where j Σ is defined by Equation ( 2).The diffusion constants (see Table 3) for copper used in the numerical simulations are taken from the literature. [27,55,56]

Thermo-Electro-Mechanical Field Problem
In the present section, the balance equations for the coupled thermo-electro-mechanical field problem and the constitutive equations are given.The governing equation of heat conductionobtained on the basis of the energy balance and the Fourier's lawthe Maxwell electrostatics equations-simplified for a conductive material to the law of charge conservation-and the mechanical equilibrium equations are defined by the relations summarized in Table 4.
Here, T is the temperature, ρ is the material density, C p denotes the specific heat capacity (at constant pressure), λ is Only one specimen tested.
Table 2. Equations for chemical potentials, flux density, and velocity for different mechanisms of diffusion.

Chemical potentials Atomic flux density Velocity
Surface diffusion [50,51] Volume diffusion [45,52] Grain boundary diffusion [53,54] Volume thermomigration [64][65][66] Volume electromigration [64,67] the thermal conductivity, φ is the electric scalar potential, and ρ e indicates the electrical resistivity, providing the link between the electric field intensity vector to the current density vector E ¼ À∇φ ¼ ρ e J.The Joule heat per unit volume and unit time generated by the electric current is defined by In Equation ( 23), σ is the Cauchy stress tensor, which is related to the Kirchhoff stress tensor by the relation τ ¼ Jσ, with is the deformation gradient.r ¼ rðR, tÞ is the radius vector of a material point in the current configuration, R is the radius vector in the initial configuration, and u ¼ r À R is the displacement vector.The analysis of the SPS problem at the microscopic level shows that large local plastic strains and small elastic strains occur. [46]In the considered cases, the maximal elastic strain intensity is approx.200 times smaller than the maximal plastic strain intensity.In the case of small elastic strains, as a consequence of the multiplicative decomposition of the deformation gradient F ¼ F e ⋅ F p into elastic and plastic parts, an additive decomposition of the strain rate tensor d ¼ ð∇vÞ S ¼ d e þ d p can be obtained.When temperature and viscous deformations are taken into account, the additive decomposition of the strain rate tensor has the form For the simulation of the inelastic deformation processes at high temperatures, a viscoelastic/viscoplastic material model is used accounting for both plastic and viscous properties (above and below the yield stress). [57] To calculate elastic strain rate d e , the hypoelastic rate constitutive equations are used as ad hoc extension of the infinitesimal theory which is true for small elastic strain case where τ J ¼ τ : À w ⋅ τ À τ ⋅ w is the Jaumann rate of the Kirchhoff stress, w ¼ ðv∇Þ A is the spin tensor, the symbol "⋅" denotes the scalar product (dot product), and "⋅⋅" is a double contraction (double dot product).The fourth-order elastic tensor is given by where E is the Young's modulus and v is the Poisson's ratio.
The symbols ⊗ and ⊗, ⊗ denote the direct dyadic product ðA ⊗ BÞ ijkl ¼ A ij B kl and the indirect dyadic products The thermal strain rate tensor d T , introduced in the Equation ( 6), is defined by the relation where α is the tangential coefficient of thermal expansion.The plastic strain rate tensor d p is calculated by the normal flow rule where the plastic potential f corresponds to the von Mises criterion with linear isotropic hardening f ¼ ffiffiffiffiffiffiffiffiffiffiffi ffi The deviator of the Kirchhoff stress tensor is defined by the expression s ¼ τ À 1 3 1trτ: The Odquist parameter is determined by the history of plastic strains ε dt.The plastic multiplayer λ from Equation ( 10) satisfies the Kuhn-Tucker loading/unloading conditions λ ≥ 0, f ðτÞ ≤ 0, λf ðτÞ ¼ 0 along with the consistency requirement λf : ðτÞ ¼ 0. The viscous strain rate tensor d v is defined by a Norton powertype equation The right Hencky logarithmic strain tensor ε ¼ ln U is used for further analysis.Here, U is the right stretch tensor appearing in the polar decomposition F ¼ R ⋅ U.The logarithmic strain tensor ε is related with d by ε The temperature-dependent material parameters E, v, α, H, A 0 , Q, n, σ Y of the constitutive model ( 7)- (11) for polycrystalline copper are taken from the literature [27,55,56,[58][59][60][61] and used in the numerical simulations in Section 3.4, 4, and 5.
To model heat flow and electric current across the interfaces, special surface contact elements are used to all the matching surfaces.The surface contact resistance is defined by the electrical contact conductance, ECC [Ω À1 m À2 ] and the thermal contact conductance, TCC [W m À2 K À1 ], which are considered both temperature and pressure dependent.The electrical current density j n at each pair of contacting surfaces is defined by expression where V 2 À V 1 is the potential difference across the interfacial gap.The amount of heat flux q n between two contacting surfaces is defined by where ϕ is a fraction of electric dissipated energy converted into heat (Joule heating).The values of ECC and TCC used in computations are taken from the article. [48]

Finite Element Model and Boundary Conditions
The two-level simulations are used for the detailed investigation of the sintering processes and to describe the powder response during the experiments.At the macrolevel, the thermoelectric analysis of a part of the miniature SPS machine (Figure 2) has been carried out.At the microlevel, the thermo-electro-mechanical analysis of the particle chain model (Figure 3) is fulfilled.
The axisymmetric finite element model for the central part of the SPS machine, including the copper powder sample, the  surrounding quartz glass tube, the graphite paper, and the TZM punches, is shown in Figure 2a.The numerical simulations have been performed with the finite element program ANSYS using a mesh with 1262 isoparametric quadratic 8-node elements (PLANE223).The thermoelectrical boundary conditions, considering convection and radiation on the external boundary of central part of the tool setup and electric and thermal contacts (see Equation ( 12) and ( 13) on the graphite paper, are given in Figure 2b.The coefficient α was chosen to be equal 20 W (m 2 K) À1 based on estimates of heat transfer during free convection of a vertical cylinder.The ambient temperature T ∞ was 30 °C.The electric current I(t) is obtained by measurements during sintering experiments (see Figure 1d). [62]Within the first 50…70 s, the electric current I(t) heats up the sample from room temperature to a preset value of 600 °C, 700 °C, or 800 °C (see Table 1 and Figure 1d).Then the sample temperature remains approximately constant until t = 600, 1200, or 2400 s (see Table 1 and Figure 1d).The resulting temperature field distribution in the tool setup at 1200 s is shown in Figure 2c.A comparison of the calculated temperature with experimental measurements at the points A and B obtained by means of the pyrometer has demonstrated a very good agreement (difference less 1.5%).Figure 2d shows the detailed distribution of the temperature field in the sintered sample at t = 1200 s.The temperature difference within the sample is approx.3 K.The temperature evolution in the center of sample is specified as discrete boundary condition in the center of the particle for the particle chain model.
The particle chain model-e.g., consisting of three particles (Figure 3a)-is used to describe the powder response at the microlevel in the experimental setup.The aim of the calculations is to analyze the distribution of inhomogeneous fields within the particle and to compare the contributions of different mechanisms for the interparticle neck growth.The finite element mesh for the axisymmetric particle model is shown in Figure 3b.The finite element simulations have been performed using a mesh with 6880 nodes and 2353 isoparametric quadratic 8-node elements (PLANE223).The boundary conditions used for the solution of the coupled thermo-electro-mechanical (TEM) problem ( 21)-( 23) are given in Figure 3c.
particle is applied in the particle center, where p is the given constant pressure on the cylindrical sample of the powder.p ¼ p N particles πR 2 sample πR 2 particle is a variable mean pressure in the equatorial cross section of the particle.Simultaneously, due to periodicity of the arrangement, a symmetry condition is prescribed on the upper boundary (OA in Figure 3c); this means that this boundary remains straight and horizontal.
For the convergence of the iterative procedure of the unsteady geometrically and physically nonlinear coupled problem over a time interval of 2400 s, 24 000 load increments were required.

Results of Finite Element Simulations of Copper Particle Chain
The results of the finite element thermo-electro-mechanical simulations of the copper particle in the middle of the particle chain (see blue box in Figure 3a) are depicted in Figure 4 for the loading case corresponding to the experiment (see Table 1, p = 10 MPa, T max = 600 °C).The six different field distributions: 1) vertical displacement field u z (Figure 4a), 2) vertical stress σ zz (Figure 4b), 3) plastic strain von Mises intensity ε p i (Figure 4c), 4) creep strain von Mises intensity ε v i (Figure 4d), 5) electric current density j jj (Figure 4e), and 6) electric scalar potential φ (Figure 4f ) are shown for a particle radius of R = 50 μm under a pressure of p = 10 MPa and an electric current I = 0.14 A at three different times t = 100, 600, and 1200 s.With increasing time at constant pressure and temperature, the neck radius and the creep strain are increased, whereby the stress relaxation is observed.The temperature field distribution in the particle after sintering for 100, 600, and 1200 s is very close to being homogeneous with a mean temperature of approximately 600 °C.The difference of the temperature inside the particle in steady state at t = 1200 s is about 0.0015 °C (approximately 0.00025%).The temperature maximum is located at the contact zone.
Due to large level of inelastic strains (maxε ), the geometrically nonlinear formulation was used to solve the task.The neck radius x in nonlinear formulation is 23.3 μm, while x in linear formulation is 24.5 μm, which is a deviation of 5.2%.A comparison of geometrically linear and nonlinear solutions for the neck radius evolution under p = 10 MPa, T max = 600 °C is shown in Figure 5.
The maximal vertical displacement |u z | in nonlinear formulation is 5.3 μm, while |u z | in linear formulation is 6.1 μm, which is even 15% larger.

Comparison between Experiment and Simulation
The comparison of FE simulation results with experimental data for the neck radius x of a copper particle system (R = 50 μm) versus time under various pressures and electric currents is shown in Figure 6.The comparison of experimental and numerical results demonstrates a very good agreement.
Figure 7 gives experimental results for the neck radius of a particle-punch system for R = 10 μm and R = 50 μm.
The results of comparing the vertical displacement (shrinkage of the sample) computed in FE simulation with the experimental data show the presence of additional displacement associated with reordering of particles, which occurs mainly at the initial stage (see Figure 8).For the shown loading case of p = 10 MPa, the additional displacement corresponds to the À6.7 μm for both considered temperatures (marked in Figure 8   x, µm t, s  the calculated curves are shifted downward by 6.7 μm, they show a very good agreement with the experiment, indicating that particle rearrangement is completed rapidly.

Influence of Pressure
The FE simulation, as well as the experiments, was carried out for different levels of pressure from 0.5 to 20 MPa.The pressure level is limited by the capabilities of the miniaturized testing machine [48] and the use of the quartz glass tube.Table 5 provides a comparison of the results at two levels of compressive pressure (0.5 and 5 MPa).The obtained solution demonstrates significant pressure sensitivity.Thus, the radius of the neck x increases by 2.9 times, the vertical displacements by 6.8 times, the plastic strain increase by 3.3 times, and the creep strain by 3.4 times with an increase in pressure from 0.5 to 5 MPa.The character of the distribution of displacement fields, deformations, and stresses is the same in both cases.However, the part of the volume subject to inelastic deformation increases significantly (see Table 5).
The neck (contact) radius x is very sensitive to the applied pressure.With increasing pressure p, the neck radius nonlinearly increases as shown as a function of sintering time in Figure 9a.
The neck radius evolution is principally defined by creep and plastic strains (see Figure 9b).This graph shows the maximum values of inhomogeneous fields in the powder particle.The location of these maximums may change over time.The dependency of neck radius x and the vertical displacement u z from pressure for different times is shown in Figure 10a,b.
While it is possible to express the neck radius evolution with a power-law approximation the neck growth exponent cannot be used to identify the dominating mechanism of densification in pressure-assisted sintering.Also, it demonstrates linear dependency from pressure (see Figure 11).Extrapolation to zero pressure gives a value of n%9.Note that the values of the exponent n, obtained in the multiphysics simulations of the particle behavior, slightly exceed the simplified analytical estimates for the two-particle model, which offers predictions in the range of n = 2-7. [63]

Influence of Temperature
The numerical simulation, as well as the experimental investigations, was carried out for different levels of temperature varied from 600 to 800 °C.Table 6 provides a comparison of results at two levels of temperature (600 and 700 °C) at 1200 s (pressure p = 5 MPa in both cases).The obtained solutions demonstrate temperature sensitivity.Thus, the radius of the neck x increases by 33%, the vertical displacements to 2.1 times, the plastic deformation increase by 0.3%, and the creep deformation by 32% with an increase in temperature from 600 to 700 °C.The character of the distribution of displacement fields, deformations, and stresses is kept.However, the part of the volume subject to inelastic deformation increases (see Table 6).
Table 5. Influence of pressure.Results of FE simulations of a copper particle (R = 50 μm) under different pressures p = 0.5 MPa and p = 5 MPa at 1200 s (temperature T max = 600 °C).
Creep strain field p = 0.5 MPa 5.5 0.29 3. The neck (contact) radius x is also sensitive to the applied temperature.It linearly increases with increasing temperature as shown as a function of sintering time in Figure 12a.
The neck radius evolution is defined principally by creep and plastic strains (see Figure 12b).This graph shows the maximum values of inhomogeneous fields in the powder particle.The location of these maximums may change over time.
The dependency of the neck radius x and the vertical displacement u z from temperature for different times is shown in Figure 13a,b.
A simple power-law approximation can also be used to describe the neck radius evolution with a temperature-dependent neck-growth exponent where the exponent demonstrates linear dependency from temperature (see Figure 14).

Influence of Diffusion Mechanisms
The dominant sintering mechanisms have been studied at the microlevel for the particle chain model under thermomechanical loading of the center of the sample in the SPS-machine with the following approach: first, the solution of the nonlinear coupled thermo-electro-mechanical problem was obtained using the finite element tool ANSYS.Second, the finite difference solution of the surface and grain boundary diffusion problem ( 16), (18)  was realized in an in-house Cþþ code tool by means of staggered scheme. [46]Finally, the solutions for volumetric diffusion (17), electromigration (20), and thermomigration (19) coupled with mechanical stress are obtained using ANSYS.
The results for the considered conditions (T max = 700 °C, p = 5 MPa) in the time interval up to 1200 s showed the predominance of the growth of the interparticle neck and shrinkage due to creep and plasticity (see Figure 15), that is in agreement with Ashby's maps of sintering mechanisms. [27,55]Surface and grain boundary diffusion contribute to changes in the neck radius of less than 3%.Neck growth due to volume (lattice) diffusion, thermomigration, and electromigration does not exceed 0.0001%.
The importance of taking into account the temperature dependence of mechanical properties in computations is indicated by  Table 6.Influence of temperature.Results of FE simulations of a copper particle (R = 50 μm) under different temperatures T = 600 °C and T = 700 °C at 1200 s (pressure p = 5 MPa). x Creep strain field T = 600 °C 16. the fact that in the dominant mechanism of neck growth due to plasticity, only 48% (indicated by (p) in Figure 15a) is caused by the appearance of pressure, and a change of 62% (indicated by (T ) in Figure 15a) is due to the decrease of the yield strength and of the hardening modulus with increasing temperature.The maximum volume diffusion strain (see Table 7) is 10 5 times lower than the creep strain (see Table 6 and Figure 12b).
For deformations due to electromigration and thermomigration, the difference is even more pronounced: they are 10 6 and 10 10 times less compared to creep deformation.The maximum diffusion flux were: volumetric j jj%10 À9 mol (m 2 s) À1 , electromigration j jj%10 À10 mol (m 2 s) À1 ], and thermomigration j jj%10 À14 mol (m 2 s) À1 ] (see Table 7).The obtained results are in very good agreement with the estimations presented in. [64]t should be noted that, according to ref.
, [46] a decrease in the particle size and of the pressure would lead to an increasing influence of the diffusion mechanisms.Also, the effect of electromigration, as explained by ref., [64] will be only noticeable for a significant increase of the electric current.
It should also be noted that in all cases there is a noticeable inhomogeneity of all the considered fields within the particle volume (see Table 7 and Figure 4).Extreme values of strains and stresses, diffusion fluxes, and driving forces are observed at the outer region of the contact zone.

Concluding Remarks
In the present research, both experimental and numerical results of the neck radius growth within the SPS process of polycrystalline copper powder are given: the numerical results are obtained  by a multiphysics microstructural thermo-electro-mechanical model additionally accounting for lattice diffusion, grain boundary diffusion, surface diffusion, electromigration, and thermomigration.A very good agreement was obtained between the numerical and the experimental results for the neck radius evolution.This enables the simulation tool to intensively study of different densification behavior in SPS.The two-level numerical simulations are carried out for the detailed study processes which were observed in the experiments: At the macrolevel, the thermoelectric analysis of a part of the miniature SPS-machine has been performed, At the microlevel, the thermo-electro-mechanical analysis accompanied by additionally accounting for lattice, grain boundary, surface diffusion, electromigration, and thermomigration of the particle chain model is fulfilled.The remarkable inhomogeneity of strain and stress fields, electric current density, temperature gradient, and diffusion fluxes is observed in the 3D domain of the particle.
Even with an explicit consideration of different densification mechanism known from the particle chain model and additional contribution due to the electric current, the dominant mechanisms for neck growth and shrinkage under the considered conditions are creep and plasticity.The results of the computations have demonstrated the importance of multiphysics simulation based on the three-way coupling of thermal, electrical, and mechanical behavior and the relevance of taking into account geometric nonlinearity.Considered in the study, pulsed electric current patterns (with equal effective current) do not affect neck growth at least for copper and under the conditions realized.For pure copper, this might only change, if capacitor discharge equipment with extremely high electrical pulse power, i.e., higher voltage and current, is used.For less thermally and electrically conductive materials, or for thicker oxide layers on the surface of the particles, the situation can also differ.The latter is currently under investigation.
The study of the impact of electrical and mechanical loads, particle size, temperature, and material on inelastic strain distribution and on densification quantified the increase of the densification and neck growth with temperature and pressure.The neck growth exponent n decreases both with increasing temperature and pressure.The obtained results will enable the further improvement of the process to achieve shorter process times, less energy consumption, and even better properties than currently possible.

Figure 1 .
Figure 1.a) Schematic representation of the experimental setup; b) setup during experiment p = 5 MPa, T = 600 °C; c) specimen during experiment p = 5 MPa, T = 600 °C, pulse 10 ms þ 140 ms; d) characteristic dependencies of minimal, maximal, and effective electric current and temperature from time for experiment p = 5 MPa, T = 600 °C, pulse 10 ms þ 5 ms; e) characteristic dependency of pressure and current density from time for experiment p = 5 MPa, T = 600 °C; f ) sample fracture surface after testing p = 5 MPa, T = 600 °C; and g) microscopic view for the neck radius measurement.

Figure 2 .
Figure 2. Experimental setup for sintering: a) axisymmetric FE mesh; b) boundary conditions; c) temperature field distribution in the tool setup at t = 1200 s; and d) temperature field distribution in the sample at t = 1200 s.

Figure 3 .
Figure 3. Representative volume element for sintering simulation at the microlevel: a) geometry of three-particle chain model; b) axisymmetric finite element mesh; and c) boundary conditions.

Figure 4 .
Figure 4. Results of simulations of a copper particle (R = 50 μm) under pressure of p = 10 MPa at t = 100, 600, and 1200 s for geometrically nonlinear formulation: a) absolute value of vertical displacement field |u z | distribution; b) vertical stress σ zz ; c) plastic strain intensity ε p i ; d) creep strain intensity ε v i ; e) electric current density |j|; and f ) absolute value of electric scalar potential |φ|.

Figure 5 .
Figure 5.Comparison of geometrically linear and nonlinear solutions for the neck radius evolution (p = 10 MPa, T max = 600 °C).

Figure 9 .
Figure 9. Evolution of a) neck radius and b) plastic and creep strain for various pressures (T max = 600 °C).

Figure 10 .
Figure 10.Influence of pressure on a) neck radius and b) vertical displacement for various times.

Figure 11 .
Figure 11.Dependence of exponent n in the approximation x R À Á n ¼ BR Àm t

Figure 12 .Figure 13 .
Figure 12.Evolution of a) neck radius and b) plastic and creep strain for various temperatures (p = 5 MPa).

Figure 14 .
Figure 14.Dependence of exponent n in the approximation x R À Á n ¼ BR Àm t

(b) 15 .
proportion of mechanisms in the sintering process on a) neck radius and b) vertical displacement.S/GB and V/E/T are surface/grain boundary diffusion and volume diffusion/electromigration/thermomigration.Table7.Comparison of diffusion strains and diffusion fluxes for different diffusion mechanisms: pressure-driven volume (lattice) diffusion, electromigration, and thermomigration for a copper particle (R = 50 μm) at 1200 s (T = 600 °C, p = 5 MPa).Diffusion strain [%]Diffusion flux[mol (m 2 s) À1 ] Driving force Volume (lattice) diffusion j σ v ¼ ÀD v ∇c À D v kT ∇trσ Electromigration j E v ¼ ÀD v ∇c À D v Ze kT ∇φ Thermomigration j T v ¼ ÀD v ∇c À D v Q v kT ∇T T

Table 1 .
Averaged neck radiuses x of particles and its variation after different sintering experiments. a)

Table 4 .
Balance equations for heat conduction, electrostatics, and mechanical equilibrium.