Multiscale Microstructures and Microstructural Effects on the Reliability of Microbumps in Three-Dimensional Integration

The dimensions of microbumps in three-dimensional integration reach microscopic scales and thus necessitate a study of the multiscale microstructures in microbumps. Here, we present simulated mesoscale and atomic-scale microstructures of microbumps using phase field and phase field crystal models. Coupled microstructure, mechanical stress, and electromigration modeling was performed to highlight the microstructural effects on the reliability of microbumps. The results suggest that the size and geometry of microbumps can influence both the mesoscale and atomic-scale microstructural formation during solidification. An external stress imposed on the microbump can cause ordered phase growth along the boundaries of the microbump. Mesoscale microstructures formed in the microbumps from solidification, solid state phase separation, and coarsening processes suggest that the microstructures in smaller microbumps are more heterogeneous. Due to the differences in microstructures, the von Mises stress distributions in microbumps of different sizes and geometries vary. In addition, a combined effect resulting from the connectivity of the phase morphology and the amount of interface present in the mesoscale microstructure can influence the electromigration reliability of microbumps.


Introduction
Three-dimensional (3D) integration is an emerging technology in which multiple materials, technologies, and functional components are vertically stacked and interconnected to form highly integrated micro-nano systems [1]. Three dimensional integration technology has the potential to increase system performance in many applications, from mobile to high-performance computing and high density memories. Techniques such as through-silicon-via (TSV) and fine-pitch bumping allow for a dense integration of heterogeneous technologies. As complementary metal oxide semiconductor (CMOS) technology barely provides the expected performance improvement, this dense integration of heterogeneous technologies is critical in maintaining Moore's momentum [2][3][4][5][6][7].
Due to technical challenges for stacking up wafers with TSVs for electrical feed through without any bumps, 3D integration with microbumps and thin chips has been favored since 2000, e.g., memoery-chips stacking, active interposers, and passive interposers [5]. Interconnects distribute the clock and other signals, provide ground, and supply power for various components of an electronic system [8]. The ever-growing need for denser integration or more input/outputs (I/Os) requires a scaling down of the dimensions and pitch of microbumps. Scaling the microbump pitch from hundreds to a few tens of microns is not straightforward because ultrafine microbumps require a stricter control of the electrical, thermal, and mechanical properties [9][10][11][12][13]. Ensuring the required reliability of the interconnect is thus critical [2]. Electromigration reliability requires improved interface control; dielectric reliability depends on the barrier-/low-k and chemical mechanical polishing (CMP)-interface control. Stress-induced voiding is driven by stress gradients, and chip-package interactions are increasingly difficult to minimize when porous dielectrics with a low Young's modulus are used [2].
An interconnect's lifetime is approximately proportional to the product of the line width and line height. When the dimensions are reduced, the reliability margin also decreases [2]; therefore, to ensure reliable 3D integrations, investigations are required into the microstructure of ultrafine interconnects. Unfortunately, there is a significant shortage of systematic studies on the microstructural formation and consequent microstructure-reliability relationship of ultrafine interconnects in the existing literature. The microstructure of microbumps, which can consist of Sn-based solder and Cu, Au, or In and Au, changes and for some cases intermetallic compounds (IMCs), e.g., Cu 6 Sn 5 , can essentially occupy the complete bump volume, e.g., through solid-liquid interdiffusion bonding (SLID) [14,15]. Intermetallic bumps have a higher melting temperature, thus making them a preferred candidate for 3D stacking without the danger of re-melting the previously made interconnections. It has been reported that SLID bonding provides a superior reliability performance [15]. In terms of growth kinetics, the activation energy for growing IMCs has been extensively reported in the literature (e.g., [16]). However, the material composition and processing can play a role in diffusion and IMC growth kinetics. Knowing only the activation energy and pre-exponential factor for the growth of IMCs is not sufficient for ensuring the reliability of microbumps in 3D integration [17]. For microbumps with dimensions of less than 10 μm, critical data, regarding the process and mechanism of the multiscale microstructural formation, the morphology and size of the phases, and the microstructural effects on reliability, need to be collected and understood to assist in a data-driven design for reliability.
This study aims to systematically investigate multiscale microstructural formation and subsequent microstructural effects on the reliability of microbumps for 3D integration using a modeling and simulation approach. Emphases are placed on the size and geometry effects on the microstructural formation and the relationship between microstructure and reliability of the microbumps. In particular, hourglass-shaped microbumps, which can be formed through a special fabrication process [18], are considered since it has been reported that an hourglass shape could improve joint fatigue lifetime [19]. The structure of the paper is as follows. Section 2 introduces the modeling methodology, including the phase field and phase field crystal models for solidification, the coupling between microstructural evolution and mechanical modeling, and the microstructure-based electromigration modeling. Section 3 is the main body of the paper; Section 3.1 presents the mesoscale microstructural formation during solidification of the microbumps using a phase field method. Atomic-scale microstructural formation during solidification of ultrafine microbumps using a phase field crystal method is presented in Section 3.2. Sections 3.3 and 3.4 present the microstructural effects on the mechanical and electrical reliability of the microbumps, with the stress effects on the phase separation and coarsening process in Section 3.3 and the coupling of microstructure, elastic mechanics, and electromigration in Section 3.4. Section 4 draws conclusions based on the modeling and simulation results.

Phase Field Model of Solidification
A phase field model introduced by Kim et al. [20] is used to simulate the microstructural evolution during solidification within one half of the cross-section of a microbump. The governing equations of the order parameter, ϕ, and the temperature field, T, are: where D is the thermal diffusivity; ∆H is the latent heat per unit volume; C p is the heat capacity per unit volume; L is the phase field mobility; ε is the gradient energy coefficient; θ is the angle between the normal vector of the solid-liquid interface and the x-axis; w is the height of the double-well potential; and T m is the melting temperature. The gradient energy coefficient, ε is assumed to be six-fold anisotropic [20]: where δ is the anisotropy and θ 0 is the orientation corresponding to the maximum anisotropy. In this study, θ 0 is set to 60°, and the six directions, which have angles of 30°, 90°, 150°, 210°, 270°, and 330° with respect to the x-axis, are called the anisotropic directions. The thermal diffusivity of the system is set as: where D S = 1.55 × 10 −5 m 2 /s and D L = 5 × 10 −6 m 2 /s are the thermal diffusivities of the solid and liquid phases, respectively. Note that this system is different from the model of Kim et al. [20], in which the thermal diffusivities are the same in both phases. The other parameters are based on the previous work by Kim et al. [20]. The boundary conditions (BCs) for the phase field model are set as follows. Zero-flux BCs are set at all of the boundaries for the phase field. For the temperature field, the curved boundaries are convectional with a coefficient of 2000 W·m −2 ·K −1 . The other boundaries are adiabatic. It is noted that this phase field model was originally used to simulate the solidification process of an undercooled pure material, and as such, it is a comparatively crude model for simulating the microstructure of microbumps. Despite the simplified nature of the simulated microstructure, the general trends of the size and geometry effects on the microstructural formation are expected to provide insight for the study of more realistic systems in which more complex microstructures can form.

Phase Field Crystal Model of Solidification
A phase field crystal (PFC) model introduced by Elder et al. [21,22] is used to simulate the solidification process in ultrafine interconnects. The governing equation of the PFC model is: where ψ is the time-averaged atomistic density, r = −0.25 is the undercooling, and ξ is the stochastic noise, which is not considered in this study. The solidification occurs from the nuclei, whose densities are set as: where A = −0.47034, q = 3 / 2 , and 0.285 ψ = . The nuclei are placed at different locations in the interconnects so that the effect of the nucleus location on the atomic arrangements can be studied. In the two-dimensional PFC model, the solid phase is assumed to have a hexagonal lattice, and the density of the liquid phase is set as ψ . The variables in the PFC model used in this study are dimensionless.

Coupling Microstructural Evolution with Elastic Mechanics
The model of Müller et al. [23] is used to simulate the phase separation and coarsening process in binary Sn-Pb alloys. This model incorporates strain energy into the phase separation process as a driving force for microstructural evolution. Therefore, the model can provide more accurate predictions regarding the microstructure of microbumps subjected to mechanical loads. The mass fraction of Sn, denoted as c, is used as an order parameter to distinguish different phases in the microstructure. The microstructural evolution is simulated by numerically solving the following equations (note that Einstein's summation convention is used): where ρ 0 is the mass density of the Sn-Pb alloy; J i is the diffusion flux; M ij is the mobility matrix; φ is the minimum between the Gibbs free energy of the α (Pb-rich) phase and the β (Sn-rich) phase; a kl is a matrix containing quantities related to surface tension; and W s is the elastic strain energy defined by: where ε kl is the total strain tensor; ε * kl is the thermal strain tensor; and C klrs is the stiffness tensor. It should be noted that all of the material properties mentioned above, i.e., ρ 0 , M ij , a kl , and C klrs , are obtained through experimental data except for C klrs of the β phase. The C klrs value of the β phase is set to 1/10 of its true value in this study to clearly demonstrate the influence of microstructure on the mechanical behavior of the microbumps. The stress reaches equilibrium much faster than the microstructure evolves; thus, the following static equation is used to solve for stress at each time step in the simulation: The above model is solved using the finite element method (FEM) in this study. The FEM enables a simulation of the phase separation process in microbumps of arbitrary geometries, thus incorporating the geometry effects on microstructural evolution. The geometries of the solder joints reported by Liu et al. [19,24], i.e., the barrel and hourglass shapes, are investigated. For the mechanical equations, i.e., Equations (9) and (10), the bottom boundaries are fixed, and the top boundaries are sheared by 0.01 μm to simulate the thermal stress caused by a mismatch in the coefficients of thermal expansion between the substrate and chip. For the microstructural evolution equations, i.e., Equations (7) and (8), a zero-flux condition for c is set to all boundaries. Initially, the microbumps are free of stress, and a fluctuation in c is imposed at the center of the domain. Microbumps with pad sizes of 1.5, 3, 4.5 and 8 μm in diameter are systematically investigated. The standoff height of each microbump is set to be equal to the pad size. The phase separation and coarsening process is simulated for 80 s, when the phase separation is complete for most of the microbumps studied. During the simulation, the shear modulus of the microbump is calculated according to the following equation: where τ is the average shear stress across the top boundary of the microbump; ∆x is shear displacement and set to be 0.01 μm; and h is the standoff height.

Microstructure-Based Electromigration Modeling
The electromigration model developed by Sukharev et al. [25] is extended for binary alloys in this study. The model incorporates the electric current and stress gradient as driving forces for vacancy diffusion and assumes that vacancy generation and annihilation occur at interfaces between different phases. The governing equation of vacancy diffusion is: where C v is the mole fraction of vacancy; J ch , J str , and J em are the vacancy flux driven by a gradient of the chemical potential, stress, and electric potential, respectively; and G is a term related to the generation/annihilation kinetics of the vacancy. Vacancy diffusion in the microbumps is a complicated process because it involves diffusion in multi-component systems. However, for simplicity in this study, the following effective vacancy fluxes based on the derivation of Darken's equation for interdiffusion [26] are used: where k is Boltzmann's constant; T is the absolute temperature; E f is the vacancy formation energy; Ω is the atomic volume; f is the ratio between the volume of the vacancy and the atom; σ h is the hydrostatic stress; e is the electron charge; V is the electric potential; C A , C B , D A , D B , Z A and Z B are the mole fractions, intrinsic diffusion coefficients, and effective charge numbers of components A and B in an A-B binary alloy, respectively; and D v is the interdiffusion coefficient expressed as: The diffusion coefficient at the interface is set as 100 times that of the bulk phase to simulate short-circuit diffusion.
The generation/annihilation term in Equation (12) is [25]: where: is the mole fraction of vacancy in the thermodynamic equilibrium state in a microbump subjected to stress. This mole fraction is dependent on the vacancy formation energy, E f , and the hydrostatic stress, σ h ; τ is the time related to the rate of vacancy generation and annihilation.
To predict a more accurate stress distribution at the interface, the atom exchange between the bulk phase and the interface is considered by introducing the concept of "a plated atom" [25]. It is assumed that plated atoms are generated or annihilated simultaneously with the vacancies at the interfaces and that the plated atoms cannot move due to their relatively small mobility compared to the vacancies. The following equation calculates the mole fraction of plated atom, C p [25]: Under a current stressing, the mole fractions of vacancy and plated atom change and thus result in a dilatation effect. Such a dilatation effect introduces an inelastic strain [25] defined as: where δ ij is the Kronecker delta and C v0 and C p0 are the mole fractions of vacancy and plated atom before a current is applied, respectively. When an inelastic strain, ε em kl , is generated, stress will develop in the microbumps. The stress reaches equilibrium much faster than the vacancy diffusion; thus, the following static equation is used [23]: where ε kl is the total strain. In this way, the vacancy diffusion and stress evolution, i.e., Equations (12)- (21), are coupled with the following continuity equation for a steady state current [25]: where ρ is the electrical resistivity. The aforementioned model is used to simulate the electromigration process in SnPb and SnCu microbumps at 150 °C for 500 h. In the simulations, the top and bottom boundaries of the joints are fixed, whereas the left and right boundaries can deform freely. The net flux of the vacancies, i.e., J ch + J str + J em , is assumed to be zero at all boundaries of the microbumps. An electric potential of +0.01 V is applied to the bottom boundaries of the microbumps, and the top boundaries are grounded. The initial conditions are set as C v0 = C eq v and C p0 = 0. The material properties used in the simulations are listed in Table 1, where C 11 , C 12 , and C 44 are the components of the stiffness tensor in the Voigt notation, and D is the intrinsic diffusion coefficient. Note that Sn exhibits a tetragonal crystal structure at 150 °C, and the material properties of the (110) crystal plane are used. In addition, the effective charge numbers of Sn and Pb in the Sn-rich and Pb-rich phases are assumed to be the same as those in pure Sn and Pb, respectively. Table 1. Material properties used in the electromigration model. Alloy Phase

Mesoscale Microstructural Formation during Solidification
This section discusses the mesoscale microstructural formation during the solidification of microbumps with an emphasis on the size and geometry effects on the microstructural formation. Note that a phase field model with a multi-component and multi-phase nature is beyond the scope of this study. Therefore, the solidification model used by Kim et al. [20] for pure materials is adopted. Despite this simplification, the features of microstructural formation during solidification and the size and geometry effects on microstructural formation in microbumps can be studied. Figure 1 demonstrates the microstructures of three microbumps with different sizes at t = 5 × 10 −8 s after solidification. Note that the ratios between the standoff height and the pad size (only half is shown in Figure 1) of the three microbumps are kept the same, i.e., 1:4. A nucleus is assumed to form at the bottom left corner of the microbump as an initial condition for solidification. Differences in the morphology of the dendrites in the three microbumps can be observed. In the 3.75-μm-standoff height microbump, the 90° primary dendritic arm (hereinafter referred to as the 90° arm) has already reached the top boundary, and the 30° primary dendritic arm (hereinafter referred to as the 30° arm) splits into two branches. In the 5-μm-standoff height microbump, the 90° arm splits into two branches. Note that only one branch of the 90° arm is shown in Figure 1b; the other branch is symmetric to that shown along the y-axis. The 30° arm also splits into two branches, but the branches exhibit obvious differences in both size and morphology, i.e., one branch dominates the other, showing a relatively fully developed dendritic structure. In the 6.25-μm-standoff height microbump, both the 90° and 30° arms split. This phenomenon of the dendritic arm splitting has been reported by other researchers (e.g., [34]) and is believed to be caused by numerical errors. In addition, the morphologies of the secondary dendritic arms between the 30° and 90° arms are slightly different among the three microbumps. Overall, there is no major difference in microstructure for the three microbumps studied here. However, the results indicate that a minute numerical error can propagate and thus influence the microstructural formation.

Size Effects on Microstructural Formation
To further investigate the size effect on the microstructural formation during solidification, Figure 2 shows the microstructure at t = 5 × 10 −7 s after solidification in three hourglass-shaped microbumps with substantial differences in standoff heights, i.e., 8, 16 and 24 μm. The effect of size on the microstructure of the microbumps can be clearly observed in Figure 2. It is obvious that a larger microbump contains longer primary dendritic arms. In addition, it is noted that the degree to which the 30° arms deviate from the 30° direction increases as the standoff height increases. This phenomenon can be understood as follows. As a larger microbump contains more dendritic arms, the competition among the dendritic arms is more intense than in a smaller microbump, resulting in the deviation of the 30° arm from its ideal direction. Furthermore, it can be observed from Figure 2 that, in general, the width of the dendritic arm decreases as the standoff height increases. In other words, a larger microbump contains finer dendrites. Moreover, the microstructure of a larger microbump, due to its larger dimensions, appears to be more homogeneous. The size effects on the microstructure of microbumps with different standoff heights, as discussed above, raise the question of why such size effects are not observed in Figure 1. Considering the microstructural formation process of microbumps of different sizes can shed some light onto this problem. The formation of dendrites goes through three stages: nucleation, growth, and coarsening. In this study, a solid nucleus is placed at the bottom left corner of the microbump as an initial condition for solidification; therefore, the nucleation process is ignored. The time span required for the growing dendrites to touch the boundaries of the bump is less in a smaller microbump. When dendrites fully occupy a microbump, the coarsening stage begins. As a result, the dendrites in a smaller microbump enter the coarsening stage earlier and, thus, exhibit a coarser microstructure. Therefore, the reasons that the size effects are not obvious in Figure 1 is two-fold. First, the difference in size of the three microbumps is small. Second, not enough time is allowed for solidification, and the formation of dendrites has not reached the coarsening stage. Figure 3 demonstrates the microstructure of four barrel-shaped microbumps with different geometries at t = 2 × 10 −7 s after solidification. Note that that the pad sizes are all set as 20 μm, but the standoff heights are 7.5, 10, 12.5 and 15 μm. In this case, when the geometry of the microbump changes, the standoff height varies. Therefore, the influence on microstructural formation comes from a combined effect of both the size and geometry. Differences in the morphology of the dendrites can be observed in the four microbumps shown in Figure 3. Apart from the 15-μm-standoff height microbump, the dendrites fully occupy the microbump. As discussed in the previous section on size effects, a larger microbump contains longer primary and secondary dendritic arms. The dendritic structures in the region between the 30° and 90° arms are different in the four microbumps. The fact that a larger microbump contains a finer dendritic structure can also be observed in Figure 3. The geometry influences the microstructural formation through the curved boundary. In the solidification model, heat from the microbump is exchanged with the environment through the curved boundary by convection, with a convection coefficient of 2000 W·m −2 ·K −1 . However, the cooling effect from the curved boundary is not obvious in Figure 3 because the dendritic structure does not obviously change near the curved boundary.

Geometry Effects on Microstructural Formation
To investigate this point further, Figure 4 demonstrates the resulting solidification microstructures when different convection coefficients are set for the curved boundary. Note that the cooling effect from the curved boundary occurs, as shown in Figures 4c,f, when a convection coefficient of 1 × 10 7 W·m −2 ·K −1 is used. In reality, the natural convection coefficient in air is much less than 2000 W·m −2 ·K −1 , explaining why the geometry effect is not obvious in Figure 4a,b,d,e. It is worth noting that the simulated dendritic microstructure is a simplified version of the real microstructure of microbumps. First, mass diffusion is not considered in the solidification model. It is expected that the curved boundary will influence the diffusion kinetics and, hence, the solidified microstructure of real microbumps. Second, an undercooling of 200 K is used in the model to accelerate the solidification process; therefore, the heat generated during solidification cannot leave through the curved boundary unless a high heat exchange coefficient is set. As such, the geometry effects on the microstructural formation of real microbumps may be more complicated.  Figure 5 demonstrates the atomic-scale microstructure inside the interconnects of different geometries at the same non-dimensional time, t = 150. In the simulations, three nuclei are placed at the same positions inside the interconnects, i.e., at the coordinates (−80,80), (80,80) and (0,320). The nuclei solidify to form grains 1-3, as labeled in Figure 5. The three grains are of different orientations, i.e., 0°, 15° and 30°, respectively. The yellow region in Figure 5 represents the undercooled liquid phase. The distance marked "d1" in each interconnect measures the distance between the solidification fronts of grains 2 and 3. The distance marked "d2" in the inserts of Figure 5c,d measures the distance between the bottom right corner, i.e., (0,150), and the solidification front, which is indicated by the solid line, of grain 3. Two inserts are included in Figure 5c,d to more clearly show the atomic arrangement at the bottom right corner. The results indicate that the solidification rates of the grains are influenced by the geometry of the interconnects. The grains in the barrel-shaped and rectangular interconnects grow faster than those in the two hourglass-shaped interconnects, as suggested by the "d1" values. However, there is only a slight difference between the "d1" values in the barrel-shaped and rectangular interconnects because the amounts of undercooled liquid in the interconnects are different due to their geometries. In the barrel-shaped and rectangular interconnects, the amount of undercooled liquid can provide a sufficient supply for the solidification process. The two hourglass-shaped interconnects, on the other hand, contain less undercooled liquid, and the solidification process slows down due to the limited supply. The geometry effect on the solidification rate is more obvious when the grains solidify along the curved boundary. For example, grain 3 in the thinner hourglass-shaped interconnect solidifies more quickly than that in the wider hourglass-shaped interconnect, as indicated by "d2" values shown in the inserts of Figure 5c,d. Note that the curvature of the curved boundary in the thinner hourglass-shaped interconnect is larger, suggesting that the solidification rate of grain 3 increases with the curvature of the boundary. (c) (d) Figure 6 demonstrates the atomic-scale microstructure of the interconnects of different geometries at a non-dimensional time, t = 1000. Three grains with different orientations are found in each interconnect. The orientations of the grains are shown by solid arrows in Figure 6a. Three grain boundaries (GBs), where dislocations aggregate, are also found in each interconnect. Differences can be found with respect to the GBs in the interconnects. The length of the GBs between grains 1 and 3, as well as the dislocation density, differs among the interconnects. In addition, the GB between grains 1 and 2 in the hourglass-shaped interconnect is shifted compared to that in the thinner hourglass-shaped interconnect, as indicated by the dashed lines in Figure 6c  (c) (d) Figure 7 shows the atomic arrangements in an hourglass interconnect at a non-dimensional time, t = 1200. Figure 7a shows the perfect lattice, which is generated using Equation (6) and then cropped by the hourglass geometry. Figure 7b displays the simulated atomic arrangements obtained by solving the PFC governing equation. In the simulations, only the 15° nucleus is placed at position (−80,80) in the interconnects in order to avoid interactions with other grains. The result from the PFC simulation shows that the atomic arrangement is influenced by the geometry. For example, the atomic arrangements in the black rectangles of Figure 7 differ. The simulated atomic arrangement in the black rectangle of Figure 7b forms a GB-like structure, which is not observed in the perfect lattice, as shown in Figure 7a. To more clearly demonstrate the differences in the atomic arrangements, stream-line plots of the displacements of the atoms in a PFC simulation relative to the corresponding atoms in a perfect crystal lattice are constructed according to the following procedure. First, the centers of the atoms in the two lattices are identified by an image processing method. Second, for every atom A in the perfect lattice, the corresponding atom B in the simulated lattice is found by a criterion stating that B is the atom closest to A in the simulated lattice. Finally, a vector pointing from A to B is drawn to show the displacement. Figure 8 shows stream-line plots of the atom displacement in the hourglass-shaped and barrel-shaped interconnects. The nucleus was placed at the coordinate (−80,80). Note that the vectors are scaled by a factor of two to more clearly show the displacements. In Figure 8a, it is obvious that the atoms above y = 250 in the hourglass-shaped interconnect have a relatively larger displacement. In addition, the direction of the displacement rotates clockwise from approximately 15° with respect to the horizontal line at the curve boundary to approximately −15° at the center line of the microbump. For the atoms below y = 250, the displacements are virtually zero. In contrast, the displacements of the atoms in the barrel-shaped interconnect exhibit a different behavior, as shown in Figure 8b. The displacements are smaller, and the stream-line plot does not form a global pattern. However, some local patterns do form, e.g., the pattern near the coordinate (−175,25).  Figure 9b. Such changes indicate that the position of the nucleus is a key factor in determining the pattern of the stream-line plot for atom displacement. Furthermore, the displacements in the barrel-shaped interconnect are smaller than those in the hourglass-shaped interconnect, similar to Figure 8. Therefore, the magnitude of the atom displacement depends on the geometry of the interconnect, while the pattern of the displacement distribution is influenced by the location of the nucleus.

Effects of Stress on Phase Separation
The microstructures formed in a barrel-shaped Sn37Pb microbump when under a shear stress and when free of stress are studied to investigate the sensitivity of the microstructure to an external stress. The initial microstructure and the microstructures formed at t = 2 s in the phase separation process are demonstrated in Figure 10. Note that only one quarter of the microbump is shown. Figure 10b shows that the phase separation proceeds faster when a shear stress is applied, in particular, at the corner of the microbump, where the stress concentrates. In addition, ordered arrangements of the phases are found when external stress is applied. That is, the phases tend to grow along the boundaries, as shown in Figure 10b. This ordered phase arrangement can be observed even if the phase separation is complete, as shown in Figure 11a. In contrast, Figure 11b shows the relatively disordered phase distribution that occurs when no external stress is applied. It should be noted that in Müller's work, an unrealistic stress level, i.e., 10 GPa, was required to observe a dramatic impact on the microstructural development [23]. However, the maximum stress applied for the result presented in Figure 11 is approximately 300 MPa. The stronger effects of stress on the microstructure in this work are most likely due to the incorporation of geometrical constraints on the diffusion process, which was not considered in Müller's work. The effects of stress on the microstructure of microbumps should not be neglected. All of the microstructures discussed in the following sections are simulated under a shear stress.

Size Effects on Microstructure
The size effect on the phase morphology in ultrafine microbumps is demonstrated in Figure 12. The results show the microstructures after the phase separation process has proceeded for 80 s in barrel-shaped Sn37Pb microbumps of different pad sizes, i.e., 1.5, 3, 4.5 and 8 μm. It is clear that larger phases with respect to the size of the microbump are formed in smaller microbumps. In particular, in the 1.5-μm microbump, the size of the β phase is nearly 1/5 of the standoff height, as shown in Figure 12. Thus, the microstructures of smaller microbumps are highly heterogeneous. In addition, the ordered phase arrangement induced by stress becomes less observable for microbumps with a pad size below 3 μm due to relatively larger phase sizes. A quantitative characterization method must be developed to further study the effects of stress on the microstructure of microbumps at such small scales.

Figure 12.
Microstructures of barrel-shaped Sn37Pb microbumps with different pad sizes after the phase separation process has proceeded for 80 s.

Geometry Effects on Microstructure
The effects of geometry on the microstructure of microbumps can be observed by comparing Figures 12 and 13. Figure 13 demonstrates the microstructures after 80 s of the phase separation process for hourglass-shaped Sn37Pb microbumps for the same pad sizes used in Figure 12. A similar relationship between the microstructural heterogeneity and the size of the microbump can be observed for the two geometries. However, the hourglass-shaped microbump contains a larger fraction of the ordered phase region than the barrel-shaped microbump for the same pad size and standoff height. Therefore, the microstructure of an hourglass-shaped microbump exhibits a higher degree of heterogeneity.  Figure 14a,d, respectively. To demonstrate the results more clearly, Figure 14b,e display the stress distribution in the α phase, and Figure 14c,f show the stress distributions in the β phase for the barrel-shaped and hourglass-shaped microbumps. It is obvious that stress tends to concentrate in the β phase because the Young's modulus of the β phase is larger. The maximum stress in the hourglass-shaped microbump is 56 MPa, which is lower than that in the barrel-shaped microbump, i.e., 87.2 MPa, implying that the probability for damage to form and propagate may be lower in an hourglass-shaped microbump. Stress concentrates at different locations in the two geometries. In Figure 14a, higher stress levels are observed at the four corners of the barrel-shaped microbump. In contrast, higher stresses occur near the left and right boundaries and in the center of the hourglass-shaped microbump, as shown in Figure 14d.  Figure 14a. That is, the microstructure primarily influences the distribution of stress at low levels, i.e., below 20 MPa, in the barrel-shaped microbump. However, the region of stress greater than 30 MPa is larger in the hourglass-shaped microbump, as shown in Figure 14d, where the spatial distribution of stress is primarily determined by the microstructure. The above results suggest that the relationship between the microstructure and von Mises stress in the microbump is also related to the geometry of the microbump.

Shear Moduli of Microbumps
During the phase separation and coarsening process, the microstructural evolution can influence the mechanical properties of the microbumps. In particular, the shear moduli of the microbumps change with the microstructural evolution. In this study, the microbumps are under a shear load and the microstructures of the microbumps evolve in order to minimize the strain energy. Figure 15a clearly shows that the volume fraction of the β phase is almost the same for the microbumps of different sizes and geometries. Therefore, the strain energy must be minimized through a change in morphology, i.e., forming the ordered phase structures as shown in Figures 12 and 13. The presence of the ordered phase structures makes the microbumps less stiff under the shear load as indicated by a decrease of the shear modulus as a function of time as shown in Figure 15b. The size and geometry effects on the shear moduli of the microbumps can also be observed from Figure 15b. First, the shear moduli of the microbumps of different sizes and geometries at t = 0 s, when microstructure has not yet formed, suggest that the shear moduli of the microbumps are independent of their sizes and determined by their geometries only. The hourglass-shaped microbumps exhibit smaller shear moduli than the barrel-shaped ones throughout the microstructural evolution process from t = 0 s to 80 s. Second, the size effect on the shear moduli of microbumps can be observed when microstructure forms and evolves. The shear moduli of the smaller microbumps decrease faster with time. The difference in the shear modulus between the microbumps of the same geometry but with a pad size of 1.5 µm and 8 µm can reach 0.5 GPa at t = 80 s.

Composition Effects on Microstructure
The phase separation process in microbumps of three different compositions, i.e., Sn47Pb, Sn37Pb, and Sn27Pb, is studied using an hourglass geometry in this section. According to the above-mentioned results, the interaction between the microstructure and the von Mises stress is stronger in the hourglass geometry. Therefore, the composition effects should be more clearly demonstrated in the hourglass-shaped microbumps. The microstructures of the 8-μm-pad size microbump of the three compositions after 80 s of the phase separation are shown in Figure 16. It can be observed that the volume fraction of the β phase dramatically increases with a decrease in Pb content, which is consistent with the Sn-Pb binary phase diagram; thus, the effect of composition on the microstructure of the microbumps is significant. In the Sn47Pb microbump, the α phase and β phase are separated into small regions, while a completely connected region of β phase is formed in the Sn27Pb microbump. The morphology of the phases also changes with the composition. The β phase in the Sn47Pb microbump is the finest because the mobility of Sn in the Pb-rich regions is lower than that in the Sn-rich regions by three orders of magnitude. Thus, the kinetics of phase separation is relatively slow in the Sn47Pb microbump.

Electromigration in Sn27Pb and Sn15Cu Microbumps
The vacancy diffusion and stress evolution induced by electromigration in microbumps with compositions of Sn27Pb and Sn15Cu are compared in this section. The effects of microstructure on the distributions of vacancy and stress are also studied. In particular, the differences between the Pb-rich phase and the Cu 6 Sn 5 IMC phase in determining the electromigration behavior of the microbump are emphasized. The microstructure shown in Figure 18 is used for the electromigration simulations for both the Sn27Pb and Sn15Cu microbumps to eliminate the effects of phase morphology. Further investigations on the effects of morphology on electromigration will be discussed in the following sections. In Figure 18, the light phase represents the Sn-rich phase, whereas the dark phase represents the Pb-rich phase for the Sn27Pb microbump and the Cu 6 Sn 5 phase for the Sn15Cu microbump. First, the microstructure influences the vacancy diffusion by affecting the distribution of the electric field intensity, as shown in Figure 19. Figure 19 demonstrates that the electric field intensity tends to concentrate at the four corners of the barrel-shaped microbump. However, compared to the Sn-rich phase, a stronger electric field is found in the Pb-rich phase and the Cu 6 Sn 5 IMC phase due to their higher electrical resistivities. In particular, the intensity of the electric field in the Pb-rich phase is approximately twice as large as that of the Cu 6 Sn 5 IMC phase, which indicates that the Pb-rich phase provides a larger driving force for the current-induced vacancy flux than the Cu 6 Sn 5 IMC, according to Equation (15). Second, the direction of the vacancy flux is slightly influenced by the phase morphology, which is indicated by the arrows shown in Figure 19. Such an effect is more obvious in the Sn27Pb microbump, as shown in Figure 19a, because the difference in electrical resistivity between the Sn-rich phase and the Pb-rich phase is approximately six times greater than that between the Sn-rich phase and the Cu 6 Sn 5 IMC phase.  Other microstructural parameters related to the vacancy diffusion include the vacancy formation energy and parameter CDZ, i.e., the product of the mole fraction, diffusivity, and effective charge number, as introduced in Equation (15). First, the vacancies are predicted to accumulate in the phase with the lower vacancy formation energy, i.e., the Pb-rich phase shown in Figure 20a and the Sn-rich phase shown in Figure 20b. Second, the vacancies in the Sn27Pb microbump move more easily under the electric current than those in the Sn15Cu microbump because the CDZ in the Pb-rich phase is larger. Figure 20a shows that the vacancies in the Sn27Pb micro joint tend to migrate towards the cathode. In contrast, the mole fraction of vacancy is almost uniformly distributed throughout the phases of the Sn15Cu microbump, as shown in Figure 20b. To clarify this point further, the distributions of the vacancies along the center line shown in Figure 18a before and after current stressing are plotted in Figure 20c. Note that the interfaces in the microstructure studied here are not interconnected. The vacancies must diffuse across the bulk phase to reach the cathode. As a result, fast diffusion of the vacancies along the interfaces does not dominate the electromigration process for the microstructure studied here.
Accompanying the vacancy diffusion process, the internal stress in the microbump evolves due to the dilatation effect. Figure 21 demonstrates that stress in both the Sn27Pb and Sn15Cu microbumps concentrates at the interfaces. As a consequence, the risk of void formation is higher at the interfaces, possibly explaining the experimental observation of voiding at the interfaces in SnPb solder strips [35]. Figure 22 shows the distributions of plated atoms and von Mises stress along the center line shown in Figure 18. It is noted that higher stresses can be found in regions where more plated atoms are generated or annihilated, as shown in Figure 22. Due to this fact, stress concentration is most likely caused by the strain resulting from the generation or annihilation of vacancy-plated atom pairs. In the Sn27Pb microbump, the vacancies are predicted to accumulate at the cathode and deplete at the anode, as shown in Figure 20c. Therefore, more vacancies are annihilated at the cathode and generated at the anode, resulting in higher stress concentrations at both electrodes, as shown in Figure 21a. In contrast, the vacancy distribution in the Sn15Cu microbump remains almost unchanged before and after a current is applied, as shown in Figure 20c. Therefore, the electromigration-induced stress is lower in the Sn15Cu microbump.

Effects of Phase Coarsening on Electromigration
The results discussed above suggest an intimate relationship between the microstructure and electromigration behavior in microbumps. This section further investigates the morphological effects on vacancy diffusion and stress evolution during electromigration. The phase coarsening process in a Sn37Pb microbump is studied first.
The microstructures formed in a Sn37Pb microbump after aging at 125 °C for 300, 900 and 3000 s are demonstrated in Figure 23. The volume fractions of the Pb-rich phase are nearly the same for the three microstructures, i.e., 40.09%, 39.93% and 39.81%. The major changes in phase morphology during a coarsening process lie in the size of the phase and the amount of interface. Figure 24 shows the temporal evolution of the maximum von Mises stress in the Sn37Pb microbump during an electromigration simulation based on a microstructure aged at 125 °C for 300 s to 3000 s. The amount of interface is high at the initial stage of the coarsening process, as shown in Figure 23a.
At this time, the vacancies can reach vacancy sources or sinks after diffusing a relatively short distance. As a result, vacancy generation and annihilation occur more quickly, mitigating vacancy accumulation and releasing stress concentration. Therefore, the maximum stress in the microbump aged for 300 s is the lowest. After aging for a longer time, the amount of interface decreases, as shown in Figures 22c and 23b. In these cases, the vacancies move longer distances before encountering vacancy sinks, resulting in vacancy accumulation and stress concentration. Voiding likely occurs in the microbump when the stress caused by electromigration exceeds the yield stress of the microbump [36]. In this study, the void nucleation time is defined as the time when the von Mises stress reaches the yield stress of the Sn37Pb solder, i.e., 18.1 MPa [37]. The void nucleation time is found to decrease as the aging time increases, as shown in Figure 25; however, the microstructure in the microbump does not remain constant during the electromigration process. The electric current applied to the microbump can also accelerate the phase coarsening process by increasing the system temperature through Joule heating. Therefore, the void nucleation time may be overestimated in this study. Further investigations into the coupling between electromigration and microstructural evolution are currently ongoing to provide more accurate predictions on the lifetime of electromigration.

Effects of Composition on Electromigration
The effects of composition on microstructure are more complicated than the phase coarsening case because the amount of each phase and the amount of interface are both related to the composition. The microstructures of microbumps with compositions of Sn47Pb, Sn37Pb, and Sn27Pb are demonstrated in Figure 26. The corresponding volume fractions of the Pb-rich phase are 52.38%, 40.56% and 29.04%, respectively. As shown in Figure 26, the change in the volume fraction of the Pb-rich phase is accompanied by a change in the connectivity of the Pb-rich phase, which can affect the diffusion path for the vacancies. Figure 26a shows that the Pb-rich phase in the Sn47Pb microbump is completely interconnected. An interconnected Pb-rich phase can provide a fast diffusion path for vacancies from the anode to the cathode under an electric current because the vacancy can move more easily in the Pb-rich phase than in the Sn-rich phase. In contrast, the Pb-rich phase in the Sn27Pb microbump is scattered in the Sn-rich matrix, as shown in Figure 26c. In this case, vacancies driven by the electric current from the anode to the cathode must travel through the Sn-rich phase. As a result, vacancy accumulation at the cathode and depletion at the anode are slower in the Sn27Pb microbump. This result also indicates that the maximum stress in the three microbumps should be in the order of Sn47Pb > Sn37Pb > Sn27Pb.
The opposite conclusion may be drawn if we consider the effect of vacancy generation and annihilation at the interface. That is, a relatively larger amount of interface in the Sn-47Pb microbump, as shown in Figure 26a, can reduce the maximum stress by rapid generation or annihilation of vacancies, as discussed above. Figure 27 shows that the combined effects of the connectivity of the Pb-rich phase and the amount of interface result in a non-monotonic relationship between the maximum von Mises stress and the Pb content. The calculated maximum stresses are in the order of Sn37Pb > Sn47Pb > Sn27Pb, but it should be noted that the maximum stress in the Sn37Pb microbump is only slightly larger than that in the Sn47Pb microbump. This result may explain the experimental observation of rapid hillock formation in eutectic SnPb solder strips [35].

Conclusions
This paper simulates mesoscale and atomic-scale microstructural formation during solidification and discusses the size and geometry effects on the microstructures of ultrafine microbumps. In addition, simulations on coupled microstructural evolution, mechanical stress, and electromigration have been performed, and the microstructural effects on the reliability of microbumps have been highlighted. The following conclusions can be drawn: 1. The size and geometry of the microbumps can influence the mesoscale microstructural formation during solidification; a smaller microbump solidifies more quickly than a larger counterpart. The length of the dendritic arms increases with the size of the microbump, and a larger microbump contains finer dendrites and a relatively more homogeneous microstructure. 2. The simulated atomic arrangements obtained during solidification using the phase crystal method exhibit a displacement when compared to the ideal crystal lattice due to the combined effects of the location of the nucleus and the geometry of the microbump. The magnitude of the displacement depends on the geometry of the microbump, while the distribution of the displacement is controlled by the location of the nucleus. 3. Ordered phase growth along the boundaries of the microbump can be observed due to influences from the imposed external shear stress and the geometrical constraints on diffusion. The formation of the ordered phase structures reduces the shear moduli of the microbumps to minimize the strain energy. The shear moduli of the hourglass-shaped microbumps are smaller than the barrel-shaped ones and decrease faster in the smaller microbumps as the microstructures evolve with time. 4. More vacancy accumulation and stress concentration occur in the Sn27Pb microbump compared to the Sn15Cu microbump during electromigration. Hence, a lower electromigration reliability of the Sn27Pb microbump can be expected. It is also found that phase coarsening in the Sn37Pb microbump reduces the amount of interface and accelerates vacancy accumulation at the cathode and depletion at the anode. Therefore, stress can quickly increase during electromigration in Sn37Pb microbumps aged for longer times, thus increasing the risk of voiding. 5. A combined effect caused by the connectivity of the Pb-rich phase and the amount of interface should be considered when comparing the electromigration behaviors of the Sn47Pb, Sn37Pb, and Sn27Pb microbumps. The maximum von Mises stress present in the Sn37Pb microbump increases with the electromigration time and depends on the composition of the microbump, which could result in nonuniformity in bump heights if the process variables are not tightly controlled.