Interfacial thermal stresses in solder joints of leadless chip resistors

Microelectronics Reliability 46 (2006) 873–884 www.elsevier.com/locate/microrel

Interfacial thermal stresses in solder joints of leadless chip resistors H.R. Ghorbani, J.K. Spelt

*

Department of Mechanical and Industrial Engineering, University of Toronto, 5 King’s College Road, Toronto, Ontario, Canada M5S 3G8 Received 4 January 2005; received in revised form 15 April 2005 Available online 15 August 2005

Abstract A novel two-dimensional analytical model has been developed for the interfacial thermal stresses in the solder joints of a leadless chip resistor (LCR) assembly under both plane stress and plane strain conditions. Both global and local expansivity mismatches are incorporated into the model. Interfacial thermal stresses are approximated using elementary strength of materials theory. Governing differential equations are linearized through a finite difference discretization procedure. The conditions of zero shear stress at the free edges and self-equilibrated peel stresses are met. The model is more accurate than existing ones, is mathematically straightforward, and can be extended to include inelastic behavior. The results are compared with available data in the literature and finite element analysis, and a geometry optimization is performed for a sample LCR as an example of the use of the model.  2005 Elsevier Ltd. All rights reserved.

1. Introduction Leadless chip resistors (LCRs) are surface-mounted microelectronic assemblies in which the package-toboard interconnection at each end of the chip resistor consists of a solder joint between metallized areas on the resistor and printed circuit board (PCB) (Fig. 1). When an LCR assembly is subject to changes in temperature, complex stresses and strains are generated by the different thermal expansions of: (i) the resistor and the substrate (PCB) to which it is soldered (global mismatch), and (ii) the solder and the materials (resistor, PCB) to which it is bonded (local mismatch). The failure mechanism of primary concern for LCRs in particular,

* Corresponding author. Tel.: +1 416 978 5435; fax: +1 416 978 7753. E-mail address: [email protected] (J.K. Spelt).

and for surface mount assemblies in general, is thermal creep-fatigue due to cyclic differential thermal expansion and contraction. Most commonly, fatigue cracks are observed to originate at the interface between the solder joint and the resistor. Predicting the thermal fatigue life of a solder joint requires a model of the solder response to cyclic temperature and stress, and a thermomechanical model to calculate the required input data; e.g. stress, strain, or strain energy range per cycle. Finite element analysis (FEA) gives singular stress fields near the free ends of such solder joints, and consequently is strongly mesh sensitive [2]. As well, FEA can be time-consuming to use as a primary design tool since changes in the solder geometry require that a new model be built. Analytical modeling of LCR assemblies has therefore been investigated extensively [15]. In general, the available analytical models of leadless assemblies may be classified into two categories: global models and global–local models.

0026-2714/$ - see front matter  2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.microrel.2005.05.022

874

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

Fig. 1. Schematic of an LCR assembly. The large open and shaded arrows illustrate global and local thermal expansion, respectively.

In global modeling, each solder joint is replaced by an equivalent beam element with the properties of the solder. Only the global coefficient of thermal expansion (CTE) mismatch between the component (chip resistor) and the substrate (PCB) is considered. This means that the thermal expansion is measured only from the centre of the component (axis A–A in Fig. 1). In addition, most global models assume that bending stresses are insignificant and the bending moments serve to accommodate some of the global expansion mismatch, thus reducing the shear stress magnitude in the joint [3,4]. Vandevelde et al. [16] developed a more sophisticated global model that also includes the bending stresses. In reality, however, the interfacial shear and peel stresses have complex distributions that cannot be captured by global models. Clech [1] added local CTE mismatch in leadless assemblies to HallÕs global model. In this model, the local mismatches were assumed to be independent of the global mismatch, peel stresses were neglected and the interfacial shear stresses were taken to be equal and constant for the lower and upper interfaces of the solder joint. Building on the large literature in the field of thermal stress analysis in layered beams (see [2] for a review), Suhir [11,13] developed a model of bimaterial assemblies adhesively bonded at the ends; mimicking the geometry of LCR assemblies with the solder joints being represented by adhesive. The local and global stresses were derived separately and then superimposed, although the coefficient of thermal expansion of the adhesive layer did not enter the formulation. The assembly was treated as an adhesively bonded bimaterial structure subject to external (global) shearing load and bending moments due to the global rigidity of the structure. However, the effect of the global bending moment on the interfacial shear stress was not considered. In addition, these models assume a thin-elongated assembly so that the radius of curvature can be taken to be constant through the thickness and transverse shear deformations can be ignored. By assuming a constant radius of curvature through the thickness of the assembly, the interfacial shear and peel stresses were decoupled, a common simplification which leads to highly approximate solutions

for short assemblies. In reality, however, the solder joints in LCRs (Fig. 1) are usually relatively short trilayer structures with large aspect ratio (thickness/ length). As well, the shear stresses at upper and lower interfaces were assumed to be the same, an assumption that is not necessarily true as the interfacial shear stresses are most likely to have different distributions depending strongly on the geometry and the thermomechanical properties of LCR members. Although the condition of self-equilibrated peel stresses was imposed following earlier approaches in layered assemblies [9,10,12], the models predict large interfacial shear stresses at the free ends of the assembly violating the stress-free boundary condition. As well, the accuracy of the predicted stresses remains uncertain, as they were not compared with FEA results. The present analytical model for LCRs utilizes a recent model for thermal stresses in trilayer assemblies [2]. The new LCR model captures the complex interaction of the global and local thermal mismatches, rather than estimating the interfacial stresses due to global and local mismatches separately and then adding them. The present model builds on the concept of non-local deformation [2,8] and has the following characteristics which make it more accurate and versatile than existing models: i. interfacial shear stresses are zero at free edges; ii. interfacial peel stresses are self-equilibrating; iii. radius of curvature is not assumed to be constant through the thickness of the bonded region of the assembly; iv. it is equally applicable to LCRs with short or long solder joints; v. it accounts for transverse shear and throughthickness deformations; vi. compliance parameters are included for both plane stress and plane strain conditions; vii. the governing equations are solved using a relatively simple finite difference procedure. The last point makes the approach mathematically straightforward and facilitates its extension to inelastic behavior; a key factor in predicting thermal fatigue of solder joints.

2. Problem formulation 2.1. Free-body diagram Fig. 1 shows a plane view of a simplified LCR assembly with two solder joints of length L = 2l located at each end of the resistor. The distance from the centerline of the assembly (axis A–A) to the centerline of each solder joint (axis B–B) is L 0 . The Poisson ratio, elastic modulus,

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

875

Tj(x), due to the interlaminar shear stresses, act at each cross-section at the mid-plane of the component, substrate, and solder joint and must be zero at x = l. At the other interfacial end (i.e. x = l), however, only T2(x) is zero. The peel stresses, pm(x), due to differences in the flexural stiffness and vertical compliance of the adjacent layers, have a symmetric distribution with respect to axis A–A, and are, therefore, self-equilibrating over each joint. Transverse shear forces, Vj(x), due to the interlaminar peel stresses, act at each cross-section at the midplane of the component, substrate, and interlayer. These forces must be zero at both free ends of the interfacial regions. Bending moments, Mj(x), due to both shear and peel stresses, act at each cross-section and are zero at x = l. At x = l, however, only M2(x) is zero.

coefficient of thermal expansion, and thickness of the assembly members are, respectively, mj, Ej, aj, and hj, where the indices j = 1, 2, and 3 denote the PCB (1), solder (2), and resistor (3) materials. Under the influence of a uniform temperature change, Dh, the resistor and PCB expand or contract horizontally relative to the centerline of the assembly axis A–A, and the solder joints expand or contract horizontally relative to their own axis B–B. In reality, the ends of the solder joints would have a fillet whose geometry is governed by the solder-pad wetting condition during solidification. As reported by Nicewarner [5] and other investigators, larger fillets provide for longer fatigue lives because of the corresponding increase in load bearing and crack propagation areas and the reduction in the stress concentration associated with a straight edge. As a result, the simplification resulting from the neglect of the fillet in the present and the other analytical models mentioned above, will overestimate the severity of the loading resulting in conservative life predictions. Fig. 2 shows the free-body diagram of half of the LCR assembly indicating interfacial shear, sm(x), and peel stresses, pm(x), stresses, where m = 1, 2 refers to the interface number, as well as the cross-sectional loads and moments (i.e. Tj(x), Vj(x), and Mj(x)). The horizontal and vertical displacements of interfacial points are denoted by umj ðxÞ and wmj ðxÞ, respectively. The interlaminar shear stresses, sm(x), due to the expansivity mismatch amongst assembly members, will vanish at the ends of each bonded region and are antisymmetric with respect to axis A–A. Axial forces,

2.2. Equilibrium equations in bonded region For the region l < x 6 l, at the right-hand side of the LCR assembly, the equilibrium of horizontal forces at each cross-section implies that (Fig. 2): Z l s1 ðnÞdn ð1Þ T 1 ðxÞ ¼ Zx l T 2 ðxÞ ¼ ½s2 ðnÞ  s1 ðnÞdn ¼ ½T 1 ðxÞ þ T 3 ðxÞ ð2Þ x Z l T 3 ðxÞ ¼  s2 ðnÞdn ð3Þ x

Equilibrium of the transverse forces at each crosssection gives:

L l A

M3

#3

l

l

B

Resistor T3 Y

l

2

τ1(x) τ1(x) #1

A

M1 Printed Circuit Board (PCB)

T1

mid-plane

p (x)

#2 Solder Joint

X

M3 (x) T3 (x)

τ2(x) τ2(x)

Inner face

V3 (x)

y x

mid-plane Outer face

p1(x)

l l

V2 (x) T2 (x)

M2 (x)

V1(x) B

M1 (x) mid-plane T1 (x)

L´

Fig. 2. Internal forces and moments in right half of LCR assembly due to uniform temperature change. Z = 0 corresponds to the x–y plane of symmetry along the centreline of the resistor.

876

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

Z

V 1 ðxÞ ¼

l

p1 ðnÞdn

ð4Þ

x

Z

V 2 ðxÞ ¼

l

½p2 ðnÞ  p1 ðnÞ dn ¼ ½V 1 ðxÞ þ V 3 ðxÞ Z l V 3 ðxÞ ¼  p2 ðnÞdn

ð5Þ

x

ð6Þ

x

Finally, equilibrium of the bending moments at cross-section yields: Z l h1 M 1 ðxÞ ¼  T 1 ðxÞ þ V 1 ðnÞdn 2 x Z l h2 M 2 ðxÞ ¼ ½T 3 ðxÞ  T 1 ðxÞ  ½V 1 ðnÞ þ V 3 ðnÞdn 2 x Z l h3 M 3 ðxÞ ¼ T 3 ðxÞ þ V 3 ðnÞdn 2 x

each

ð7Þ ð8Þ ð9Þ

These bending moments may be related to the local radius of curvature as: M j ðxÞ ¼ Dj =qj ðxÞ

ð10Þ

Ej h3j

Ej

¼ Ej for plane stress and where Dj ¼ 12 , with E Ej ¼ 1mj 2 for plane strain conditions. j

2.3. Equilibrium equations in unbonded portions In the unbonded region, L 0 6 x < (L 0  l), from Eqs. (1) and (3) it can be seen that the axial forces are constant: Z l T 1 ðxÞ ¼ T 1 ðlÞ ¼ s1 ðxÞdx ð11Þ l Z l T 3 ðxÞ ¼ T 3 ðlÞ ¼  s2 ðxÞdx ð12Þ l

Knowing that T2(l) is zero, from Eq. (2) we find that: T 1 ðlÞ þ T 3 ðlÞ ¼ 0 or, Z l l

ð13Þ

Z l h1 T 1 ðlÞ þ V 1 ðxÞdx 2 l Z l h3 V 3 ðxÞdx M 3 ðxÞ ¼ M 3 ðlÞ ¼ T 3 ðlÞ þ 2 l M 1 ðxÞ ¼ M 1 ðlÞ ¼ 

s1 ðxÞdx ¼

ð18Þ

Lastly, M2(l) = 0 so that Eq. (8) yields:  Z l h2 ½s1 ðxÞ þ s2 ðxÞ þ V 1 ðxÞ þ V 3 ðxÞ dx ¼ 0 2 l

ð19Þ

This condition eliminates the possibility of solid body rotation of the solder joint. 2.4. Interlaminar shear stresses The longitudinal interfacial displacements of the resistor, solder joint, and PCB in the neighborhood of each interface may be written as: Z x u11 ðxÞ ¼ a1 DhðL0 þ xÞ þ k1 T 1 ðnÞdn þ g1 ½s1 ðxÞ l Z x h1 dn þ ug1 ð20Þ  d1 s001 ðxÞ  2 l q1 ðnÞ Z x T 2 ðnÞdn  g2 ½s1 ðxÞ u12 ðxÞ ¼ a2 Dhx þ k2 l Z h2 x dn ð21Þ  d2 s001 ðxÞ þ 2 l q2 ðnÞ Z x T 2 ðnÞdn þ g2 ½s2 ðxÞ u22 ðxÞ ¼ a2 Dhx þ k2 l Z h2 x dn ð22Þ  d2 s002 ðxÞ  2 l q2 ðnÞ Z x u23 ðxÞ ¼ a3 DhðL0 þ xÞ þ k3 T 3 ðnÞdn  g3 ½s2 ðxÞ l Z x h3 dn þ ug3 ð23Þ  d3 s002 ðxÞ þ 2 l q3 ðnÞ where kj are the axial compliances, and gj and gjdj are the two components of the shear compliance of the jth layer. In plane stress, aj ¼ aj , kj ¼

Z

ð17Þ

1m2j Ej hj

, and gj can be obtained

from Eq. (A.1). Under plane strain, aj ¼ ð1 þ mj Þaj ,

l

s2 ðxÞdx

ð14Þ

l

This condition ensures equilibrium of the solder joint in the horizontal direction. The transverse shear forces in the unbonded regions are zero because the peel stresses in the bonded areas are self-equilibrating; therefore, from Eqs. (4) and (6): Z l V 1 ðxÞ ¼ V 1 ðlÞ ¼ p1 ðxÞdx ¼ 0 ð15Þ l Z l V 3 ðxÞ ¼ V 3 ðlÞ ¼ p2 ðxÞdx ¼ 0 ð16Þ l

The bending moments in the unbonded regions are also constant and can be found from Eqs. (7) and (8) as:

kj ¼

ð1þmj Þð12mj Þ , ð1mj ÞEj hj

and gj is that of Eq. (A.2). Under both h2

j plane stress and plane strain conditions d j ¼ 1000p 2 . The derivation of these compliances can be found in Ghorbani and Spelt [2]. In Eqs. (20)–(23), the first terms are displacements due to the free thermal expansions, the second terms are longitudinal displacements due to the axial loads Tj(x), and the third and fourth terms are the interfacial shear displacements and the flexural deformations, respectively. As well, ug1 and ug3 are the displacements of the unbonded regions of the resistor and PCB, respectively, due to the interfacial shear and peel stresses in the bonded region. These deformations consist of both axial and bending components:

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

 h1 M 1 ðlÞ ðL0  lÞ 2D1   1 h3 M 3 ðlÞ ðL0  lÞ ug3 ¼  T 3 ðlÞ þ E3 h3 2D3 ug1 ¼



1

E1 h1

T 1 ðlÞ 

ð24Þ ð25Þ

The present model can be regarded as ‘‘non-local’’ since it includes the effect of s00m ðxÞ in the interfacial displacements (Eqs. (20)–(23)). This non-localization is essential if the model is to capture the large shear stress gradient near the assembly ends, enabling the interfacial shear stress to vanish at the free surfaces [2,8]. Substituting for the radii of curvatures from Eq. (10), using compatibility of the horizontal displacements at each interface (i.e. u11 ðxÞ ¼ u12 ðxÞ and u22 ðxÞ ¼ u23 ðxÞ), and differentiating leads to the following system of coupled integro-differential equations:  k 1 T 1 ðxÞ þ k 2 T 3 ðxÞ  k 3 s01 ðxÞ þ k 6 s000 1 ðxÞ Z l Z l  k8 V 1 ðnÞdn  k 9 V 3 ðnÞdn ¼ ða1  a2 ÞDh x

x

ð26Þ  k 2 T 1 ðxÞ þ k 4 T 3 ðxÞ  k 5 s02 ðxÞ þ k 7 s000 2 ðxÞ Z l Z l  k9 V 1 ðnÞdn  k 10 V 3 ðnÞdn ¼ ða2  a3 ÞDh x

x

ð27Þ where the stiffnesses k1–k10 are given in Appendix A. Substituting for axial forces and differentiating yields:   d4 d2 k 6 4  k 3 2 þ k 1 s1 ðxÞ þ k 2 s2 ðxÞ dx dx þ k 8 V 1 ðxÞ þ k 9 V 3 ðxÞ ¼ 0   d4 d2 k 2 s1 ðxÞ þ k 7 4  k 5 2 þ k 4 s2 ðxÞ dx dx þ k 9 V 1 ðxÞ þ k 10 V 3 ðxÞ ¼ 0

ð28Þ

ð29Þ

It is interesting to note that the differential equations (28) and (29) governing the interfacial shear stresses in LCR assemblies are the same as those obtained for the trilayer assemblies [2]. It can be seen that the interfacial shear and transverse shear forces (or peel stresses) are coupled. Two other differential equations must therefore be derived from compatibility of vertical deformations as discussed below. 2.5. Interlaminar peel stresses The vertical deflections of the resistor, solder joint, and PCB in the neighborhood of each interface may be written as: Z x Z xZ x dndn0 w11 ðxÞ ¼ l1 þ wg1 V 1 ðnÞdn þ d1 p1 ðxÞ þ l l l q1 ðnÞ ð30Þ

w12 ðxÞ ¼ l2 w22 ðxÞ ¼ l2 w23 ðxÞ ¼ l3

Z

x

Z

l x

V 2 ðnÞdn  d2 p1 ðxÞ þ V 2 ðnÞdn þ d2 p2 ðxÞ þ

l Z x

V 3 ðnÞdn  d3 p2 ðxÞ þ

l

877

Z

x

Z

l x

Z

l x

Z

Z Z

l

x

l x l x l

dndn0 ð31Þ q2 ðnÞ dndn0 ð32Þ q2 ðnÞ dndn0 þ wg3 q3 ðnÞ ð33Þ

where lj are the transverse shear compliances, and dj are the through-thickness tensile compliances of the jth layer. lj ¼ hj1Gj for both plane stress and plane strain conditions. dj can be obtained from Eqs. (A.3) and (A.4) for plane stress and plane strain conditions, respectively. The derivation of these compliances can be found in Ghorbani and Spelt [2]. In Eqs. (30)–(33), the first terms are due to the Timoshenko effect (i.e. transverse deflections caused by the shear loads Vj(x)), the second terms are local displacements due to the peel stresses, and the third terms are the flexural deformations. As well, wg1 and wg3 are the vertical displacements of the unbonded portion of the resistor and PCB, respectively. These deformations are due to the bending moments M1(l) and M3(l) as follows: ðL0  lÞ2 M 1 ðlÞ 2D1 ðL0  lÞ2 M 3 ðlÞ wg3 ¼ 2D3

wg1 ¼

ð34Þ ð35Þ

Substituting for the radii of curvature from Eq. (10), using compatibility of the vertical deflections at each interface (i.e. w11 ðxÞ ¼ w12 ðxÞ and w22 ðxÞ ¼ w23 ðxÞ), and differentiating twice leads to the following system of coupled integro-differential equations: k 8 T 1 ðxÞ  k 9 T 3 ðxÞ þ D1 V 01 ðxÞ þ D2 V 03 ðxÞ  D3 V 000 1 ðxÞ Z l Z l þ D6 V 1 ðnÞdn þ D7 V 3 ðnÞdn ¼ 0 ð36Þ x

x

k 9 T 1 ðxÞ  k 10 T 3 ðxÞ þ D2 V 01 ðxÞ þ D4 V 03 ðxÞ  D5 V 000 3 ðxÞ Z l Z l þ D7 V 1 ðnÞdn þ D8 V 3 ðnÞdn ¼ 0 ð37Þ x

x

where the stiffnesses D1–D8 are given in Appendix A, and p1(x) and p2(x) in Eqs. (30)–(33) were replaced by V 01 ðxÞ and V 03 ðxÞ, respectively. Finally, differentiating Eqs. (36) and (37) gives:   d4 d2 k 8 s1 ðxÞ þ k 9 s2 ðxÞ þ D3 4  D1 2 þ D6 V 1 ðxÞ dx dx   d2 ð38Þ þ D2 2 þ D7 V 3 ðxÞ ¼ 0 dx   d2 k 9 s1 ðxÞ þ k 10 s2 ðxÞ þ D2 2 þ D7 V 1 ðxÞ dx   d4 d2 ð39Þ þ D5 4  D4 2 þ D8 V 3 ðxÞ ¼ 0 dx dx

878

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

Again, it is noted that the differential equations (38) and (39) governing the transverse shear forces in LCR assemblies are the same as those obtained for the trilayer assemblies [2].

are determined, peel stresses p1(x) and p2(x) can be calculated from V 01 ðxÞ and V 03 ðxÞ, respectively.

2.6. Boundary conditions

Upon finding the shear and peel stresses, as discussed in Section 3 below, the axial forces Tj(x) and bending moments Mj(x) can be determined using Eqs. (1)–(3) and (7)–(9). The axial normal stresses rmxj ðxÞ in each LCR member (j = 1, 2, 3) at each interface (m = 1, 2) for both plane stress and plane strain conditions may then be approximated by [6,9]:

In order to solve the coupled ordinary differential equations (28), (29), (38) and (39), for s1(x), s2(x), V1(x), and V3(x), sixteen boundary conditions are needed. In addition to the two boundary conditions provided by Eqs. (14) and (19), the following eight boundary conditions are evident: s1 ðlÞ ¼ s2 ðlÞ ¼ V 1 ðlÞ ¼ V 3 ðlÞ ¼ 0

ð40Þ



k 5 s02 ðlÞ

þ

k 7 s000 2 ðlÞ

¼

ða2



a3 ÞDh

rmxj ðxÞ ¼

T j ðxÞ 6M j ðxÞ  hj h2j

ð49Þ

ð41Þ

where ( + ) and () signs refer to the lower (m = 1) and upper (m = 2) solder joint interfaces, respectively. Under plane stress condition, the interfacial von Mises equivalent stresses, rmeq j ðxÞ, at each solder interface for each material may be calculated as:

ð42Þ

rmeq j ðxÞ ¼ ½rmxj ðxÞ2 þ pm ðxÞ2  rmxj ðxÞ  pm ðxÞ þ 3sm ðxÞ2 1=2

Knowing that T1(l) = T3(l) = 0, six boundary conditions may be obtained from Eqs. (26), (27), (36) and (37) as follows:    k 3 s01 ðlÞ þ k 6 s000 1 ðlÞ ¼ ða1  a2 ÞDh

2.7. Other stresses

D1 V 01 ðlÞ þ D2 V 03 ðlÞ  D3 V 000 1 ðlÞ ¼ 0

ð43Þ

ð50Þ

D2 V 01 ðlÞ þ D4 V 03 ðlÞ  D5 V 000 3 ðlÞ ¼ 0

ð44Þ

Under the plane strain condition, however, out-ofplane interfacial stresses, rmzj ðxÞ, will also be present and serve to increase the interfacial equivalent stress as follows [14]:

The remaining two boundary conditions aim at connecting the global and local deformations in solder joints as: u22 ðlÞ  u12 ðlÞ ¼ u23 ðlÞ  u11 ðlÞ

ð45Þ

w22 ðlÞ  w12 ðlÞ ¼ w23 ðlÞ  w11 ðlÞ

ð46Þ

These two conditions imply that the solder joint does not stand stationary against global displacements imposed by the unbonded portions of the resistor and the PCB. Substituting from Eqs. (20)–(23) into (45) and from (30)–(33) into (46), and using Eqs. (17) and (18) yields the following two boundary conditions: Z l k 6 s001 ðlÞ þ k 7 s002 ðlÞ þ ½k 11 s1 ðxÞ  k 12 s2 ðxÞ l

þ k 13 V 1 ðxÞ þ k 14 V 3 ðxÞdx ¼ ða1  a3 ÞDhðL0  lÞ ð47Þ Z l h1 h3 D9 s1 ðxÞ  D10 s2 ðxÞ D3 V 01 ðlÞ  D5 V 03 ðlÞ þ 2 2 l  ð48Þ D9 V 1 ðxÞ þ D10 V 3 ðxÞ dx ¼ 0 where the global stiffnesses k11–k14 and D9 and D10 are given in Appendix A. Eqs. (28), (29), (38) and (39) with the boundary conditions (14), (19), (40), (42), (43), (44), (47) and (48) represent the complete mathematical formulation of the problem that governs the interfacial stresses in LCR assemblies. Note that once the transverse shear forces

rmzj ðxÞ ¼ mj ½rmxj ðxÞ þ pm ðxÞ  aj Ej DT rmeq j ðxÞ

¼

f½rmxj ðxÞ

2

 pm ðxÞ þ

þ ½pm ðxÞ 

½rmxj ðxÞ

rmzj ðxÞ2

ð51Þ 

rmzj ðxÞ2

pffiffiffi þ 6sm ðxÞ2 g1=2 = 2

ð52Þ

3. Solution method A closed-form solution for the two-point boundary value problem of Eqs. (28), (29), (38) and (39) would lead to nonlinear (eigenvalue) expressions which would have to be solved numerically. The usual numerical solution for such problems involves parameterization of the differential equations followed by an appropriate shooting or relaxation technique or a combination of both. However, integral-type coupled boundary conditions (14), (19), (47) and (48) and the divergence of the Newton–Raphson iteration for short trilayer assemblies (such as LCRs) make it difficult to apply on the present boundary value problem. In addition, such an approach would make it difficult to extend the model to include inelastic effects. Consequently, a finite difference procedure was adopted [2], with the bonded region discretized into n divisions (n odd) such that i = 1 for the innermost interval and i = n for the outermost interval. The derivatives of sm(x) are approximated as:

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

sðr1Þ ði þ 1Þ  sðr1Þ ðiÞ m m ‘r smðr1Þ ðiÞ  smðr1Þ ði  1Þ sðrÞ m ðiÞ ¼ ‘r

sðrÞ m ðiÞ ¼

ð53Þ ð54Þ

where r = 1–4, and ‘ = L/n is the pitch between the interfacial divisions. Similar expressions can be written for derivatives of Vj(x) as well as for the higher derivatives. Using the above definition, the governing equations (28), (29), (38) and (39) and the boundary conditions (14), (19), (40)–(44), (47) and (48) may be easily discretized (some examples are given in Appendix B). The system of linear equations in Appendix B can then be solved for sm(i), the interfacial shear stress at the middle of each discretization interval, and similarly for V1(i) and V3(i). Peel stresses p1(i) and p2(i) may then be determined knowing they are related to the first derivatives of V1(i) and V3(i), respectively (Eqs. (B.1) and (B.2)). Other stresses can also be defined from the formulas of Section 2.7.

4. Results and discussion For illustrative purposes, numerical results are presented for the LCR assembly given in Table 1. The resistor is connected to an FR-4 printed circuit board using two SnPb eutectic solder joints at either end. The length of the solder joint (L) is 0.76 mm and the distance between the centerline of solder joints and that of the component and substrate, L 0 is 2.87 mm. The governing equations in Appendix B were solved for the interfacial shear stresses and transverse shear forces under plane strain condition using the Compaq Visual FORTRAN 6.6 [17] International Mathematical and Statistical Library (IMSL) subroutine LSARG. It was found that the number of elements (n) should be greater than 301 in order to attain solutions that are independent of n. The peel, axial, and Von Mises stresses were then calculated as described in Appendix B. The predictions of the present model were compared with those of a 2D finite element model of right half of the example LCR assembly constructed using ANSYS 8.0 [18] as shown in Fig. 3. The assembly was meshed using PLANE183 elements having 8 nodes and 16 degrees of freedom and the symmetry boundary conditions were applied to the left end of the model. Progressive refinements in mesh density (e.g. 307, 1560, 6232,

Table 1 Dimensions and properties of LCR used in the example

Resistor Solder joint PCB

h (mm)

E (MPa)

m

CTE (ppm/C)

0.65 0.12 1.23

131 000 64 540 22 000

0.30 0.40 0.28

2.8 22 17

Fig. 3. 2D FEA model of right half of the LCR of Table 1 with 6232 elements.

7907, and 9601 elements with element size as small as 3.5 · 0.8 lm) were examined for both the lower and upper interfaces. Regardless of element size or shape, all of these models displayed a stress distribution in the interfacial region near the free ends of the solder (i.e. jxj > 0.75l) that was highly sensitive to the mesh density. Away from these regions, on the other hand, the stress distribution was insensitive to mesh density with 1506 elements and greater. The present model assumes that the shear and peel stresses are the same for adjacent materials across each interface (Fig. 2). The axial stress, however, is predicted to be discontinuous across each interface (Eq. 49), depending strongly on the materials bounding the interface. Consequently, rz and req will also be discontinuous across the solder interfaces. In order to compare the solder interface stress predictions of the present model with the FE results, it was necessary to extract the ANSYS stresses at each interface for the solder, rather than use the default interface node averaging for these interface stresses. Under a uniform temperature increase of 75 C, the shear, peel, axial, and equivalent stress

120 100

Present Model FEA Pao Model

80

Shear Stress (MPa)

nþ1 : 2 nþ1 iP : 2 i6

879

60

upper interface

40 20 0 lower interface

-20 -40 -60

lower interface upper interface

-80 -1.0 -1.0 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 Normalized Distance From Center (x/l)

0.8

1.0

Fig. 4. Interfacial shear stresses in the right solder joint of the example LCR of Table 1.

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

Peel Stress (MPa)

50 40 30 lower interface 20 10 0 -10 -20 -30 upper interface -40 upper interface -50 -60 -70 Present Model -80 FEA -90 -100 -1.0 -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 1.0 Normalized Distance From Center (x/l)

Axial Stress (MPa)

Fig. 5. Interfacial peel stresses in the right solder joint of the example LCR of Table 1.

20 0 -20 -40 Lower Interface -60 -80 -100 -120 Upper Interface -140 -160 Present Model -180 FEA -200 -220 -1.0 -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 1.0 Normalized Distance From Center (x/l)

Fig. 6. Interfacial axial stresses in right solder joint of the example LCR of Table 1.

-80

Out of Plane Stress (MPa)

-100 -120

Lower Interface

-140 Upper Interface

-160 -180 -200

FEA Present Model

-220 -240

-1

-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 Normalized Distance From Center (x/l)

0.8

1

Fig. 7. Interfacial out-of-plane stresses, rz, in the right solder joint of the example LCR of Table 1.

distributions in the sample LCR are shown in Figs. 4–8 under plane strain conditions along the XY plane of

von Mises Stress (MPa)

880

240 220 200 180 160 140 120 100 80 60 40 20 0

FEA Present Model

Upper Interface / Plane Strain

Lower Interface / Plane Strain Upper Interface / Plane Stress

Lower Interface / Plane Stress

-1

-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 Normalized Distance From Center (x/l)

0.8

1

Fig. 8. Interfacial von Mises stresses in example LCR of Table 1.

symmetry at Z = 0. The results for the present model are based on the 501 divisions (n), while the FEA results are those of the model with 6232 elements (element size of 31.7 · 7.5 lm close to the free ends). Fig. 4 shows the calculated interfacial shear stresses at the lower and upper interfaces. It can be seen that both the FE and present models predict shear stresses that are not antisymmetric with respect to the solder joint centerline and that the shear stresses differ significantly at the two interfaces. Both models predict that the magnitude of the shear stress at the lower interface is greatest at the inside end. The situation is reversed for the upper interface where the present model predicts a larger shear stress at the outer end while FE predicts this maximum to be at the inner corner. This difference is due mainly to the ability of the present model to capture the fact that the shear stresses vanish at the free surfaces-a condition not satisfied by the FE model. Overall, the shear stress predicted by the FE and present model are consistent both qualitatively and quantitatively, except near the free ends where the very large stress gradients make the FE stresses singular. Particularly of note is the characteristic that the shear stresses at the midpoint (x = 0) of both interfaces are equal but not zero, almost-8.5 MPa, unlike the trend observed for trilayer assemblies [2]. The magnitude of this non-zero shear stress at x = 0 mainly depends on the global expansivity mismatch between the component and substrate and on the length L 0 ; i.e. the right side of boundary condition (47). On the other hand, the overall distribution and magnitude of the shear stresses, is a function of the thermomechanical and geometrical properties of each layer. In other words, the contributions of both the global and local expansivity and rigidity mismatches in the assembly are significant. It is interesting to note that the global model of Jih and Pao [4] predicted a constant shear stress of 6.2 MPa for the entire length of both interfaces of the example LCR.

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

Fig. 8 shows the resultant equivalent (von Mises) interlaminar stresses predicted by the FE and present models for the right-hand solder joint under both plane strain and plane stress conditions. It is noted that the plane stress results for the shear, peel and axial stress had similar distributions to those shown in Figs. 4–6 and were consistent with the FE results, but their magnitude was 30–50% lower than those under the plane strain condition. As expected, the von Mises stresses at both interfaces were much larger under the plane strain condition than under the plane stress condition, due to the presence of rz under the plane strain condition. Away from the free edges, the present model results match the FE results both qualitatively and quantitatively. In the regions close to the inner and outer ends of the joint, the von Mises stresses predicted by the present model reach finite values while the magnitudes of the FE stresses were strongly mesh sensitive. As well, the FE equivalent stresses close to the free ends differ from those of the present model because both the FE shear and axial stresses violate the zero stress condition at the free ends and the FE peel stresses are not self-equilibrating. The present results can be compared with published experimental thermal fatigue observations made with the LCR example of Table 1. Both the FE and the present model suggest that the driving force for the failure of the solder joint is greater at the upper solder interface regardless of whether the stress state is plane stress or plain strain. This is consistent with the experimental observation of Qi et al. [7] who reported that cracks usually initiated at the inside free end of the LCR solder layer and propagated outward into the fillet. It might be anticipated that the crack propagation time would be relatively short [7] due to the fairly constant magnitude of the von Mises stress along this interface. Fig. 9 shows the effect of the resistor length on the shear stress distribution at the lower interface under

10 5

Trilayer Model [2]

0 -5

Present Model (L′- l) 0

-15 -20

L′

-10 increased

Shear Stress (MPa)

Comparing this value with the actual shear stress distribution indicated in Fig. 4 illustrates the large difference between the predictions of the global models and those of FEA and the present model. Similarly, the fact that the shear stresses have opposite sign and are not equal at each interface is contradictory to the recent models of Suhir [11,13] that assume the same shear stress distribution at the upper and lower interfaces. The FE and present models both predict peel stresses that are not symmetric with respect to the solder joint centerline (Fig. 5). As with the shear stresses, the peel stresses calculated by both models are non-zero at the centerline and are approximately constant over most of the interfacial regions. As well, the present model predicts that the largest peel stresses occur at the outer edge of the upper interface. The magnitudes and distributions of the FE peel stresses at both interfaces were highly sensitive to the mesh density in the region jx/lj > 0.75 and were not self-equilibrating, although they did change sign along each interface. Therefore, quantitative comparisons between the models were not possible. For all mesh densities that were examined it was found that FE peel stresses were positive at the free ends of the lower interface and negative at the free ends of upper interface, in agreement with the present analytical model. Fig. 6 shows that the interfacial axial stress distributions predicted by both the FE and the present models have a similar trend in the right-hand LCR solder joint. Near the inner and outer ends, however, the FE results differed from those of the present model. This was mainly because the FE axial stresses are connected to the peel stresses through equilibrium equations within the context of elasticity; hence, they also become singular and could not satisfy the zero stress condition at the free edges. This boundary condition was satisfied by the present model. Comparing Figs. 4–6, it can be seen that the magnitudes of the axial stresses about the midpoint of the solder are much greater than the corresponding peel and shear stresses, and that the axial stress at the upper interface is greater than that at the lower interface. As well, it can be seen that the role of the bending moment in the solder layer is relatively small, since the axial stress remains compressive over most of the joint length. The interfacial out-of-plane stresses (rz) in the righthand solder joint are depicted in Fig. 7. Due to their connectivity to the axial and peel stresses (Eq. (51)), they carry forward the same limitation as these stresses; i.e. singularity at the free surfaces at either end of the joint. It is interesting to note that the magnitudes of these stresses are relatively large compared to the other stresses and that they are not necessarily zero at the free ends of the joint at x = ±l. Therefore, these out-of-plane stresses can play a significant role in determining the von Mises stresses in LCR assemblies with sufficient width to generate plane strain conditions in the solder.

881

present example Present Model (L′- l)

10 L′

-25 -1.0 -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 Normalized Distance From Center (x/l)

0.8 1.0

Fig. 9. Effect of L 0 on interfacial shear stress at lower interface for the resistor assembly of Table 1.

882

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

the plane strain condition. It is seen that, as L 0  l, i.e. the unbounded length of component and substrate, tends to zero (Fig. 1), the interfacial shear stress distribution approaches that calculated by the trilayer model of Ghorbani and Spelt [2] with a length of L, where the shear stress has an antisymmetric distribution. The present model does not work in the limit for L 0  l = 0 because of the boundary conditions imposed in Section 2.6. Conversely, increasing L 0 , while keeping the solder joint length L constant, accentuates the effect of global CTE mismatch, increasing the magnitude of shear stresses in the middle of the solder joint and at the inner free end, and reducing the shear stress at the outer free end. Beyond a certain unbounded resistor length (approximately 10 times the value of L 0 used in this example) the shear stress distribution becomes relatively insensitive to increases in L 0 (Fig. 9). This limiting length is a function of both the thermomechanical and geometrical properties of all assembly members. Analogous trends were observed for the other interfacial stresses.

5. Conclusions A semi-analytical model has been developed for LCR assemblies under thermal loads that is more accurate than existing models and matches most FE results. In contrast to earlier analytical models, the present model includes both local and global mismatches interactively, satisfies the zero shear stress condition at the free ends, and predicts a peel stress distribution that is self-equilibrating and non-singular. The present model captures the very large shear stress gradients near the free ends in the vicinity of the maximum shear stress. Comparisons with finite element (FE) results show similar trends in the shear, peel, axial, out-of-plane, and equivalent stress distributions, both qualitatively and quantitatively; however, at the free ends of the solder joint, FE models fail to provide the zero shear and axial stress boundary conditions. Similarly, detailed comparisons are inconclusive for peel stresses because of the strong mesh sensitivity of the FE models near the free ends of the assembly and non-self equilibrating peel stresses predicted by FE models. The present model has the practical benefit of being flexible and straightforward to implement using a finite difference scheme. By changing the constitutive properties of the finite difference elements, it is possible to model plastic and creep behavior in LCR assemblies.

and Packaging (CMAP) and its member companies, particularly Celestica Inc. Appendix A. Assembly stiffness terms  1 gj ¼ uj coth  uj Þ coth  uj ð3  mj  ð1 þ mj Þ 4p  2ð1  mj Þ L þð1 þ mj Þ uj  ðplane stressÞ  uj Gj  1 ð3  4mj   gj ¼ uj coth  uj Þ coth  uj 4pð1  mj Þ  2ð1  2mj Þ L þ uj  ðplane strainÞ  uj Gj  1 þ mj  j Þ coth w j dj ¼ ð3  mj þ ð1 þ mj Þ wj coth w 4p  4 L ð1 þ mj Þ wj  ðplane stressÞ  j Ej w  1 þ mj  j coth w  j Þ coth w j dj ¼ ð3  4mj þ w 4ð1  mj Þp  4ð1  mj Þ L  wj  ðplane strainÞ j w Ej  j ¼ phj =l. where  uj ¼ phj =L and w  2  1 h1 h22 ; D1 ¼ l1 þ l2 þ k 1 ¼ k1 þ k2 þ 4 D1 D2 1 h22  k2 ; 4 D2 k 3 ¼ g1 þ g2 ; k2 ¼

k 4 ¼ k2 þ k3 þ

ðA:2Þ

ðA:3Þ

ðA:4Þ

ðA:5Þ

D2 ¼ l2

ðA:6Þ

D3 ¼ d1 þ d2

ðA:7Þ

1 ; D4 ¼ l2 þ l3 þ 4 D2 D3

ðA:8Þ

D5 ¼ d2 þ d3

ðA:9Þ



h22

k 5 ¼ g2 þ g3 ; k 6 ¼ g1 d 1 þ g2 d 2 ; k 7 ¼ g2 d 2 þ g3 d 3 ;   1 h2 h1 k8 ¼ ;  2 D2 D1 k9 ¼

ðA:1Þ

1 h2 ; 2 D2

 h23

1 1 þ D1 D2 1 D7 ¼ D2 1 1 D8 ¼ þ D2 D3 D6 ¼

D9 ¼

ðA:10Þ ðA:11Þ ðA:12Þ

1 ðL0  lÞ2 2 D1 ðA:13Þ



 1 h2 h3  ; k 10 ¼ 2 D2 D3  k 11 ¼

0

D10 ¼

1 ðL  lÞ 2 D3

 h21



2

 h23

ðA:14Þ

1 1 1 1 ðL0  lÞ; k 12 ¼  þ ðL0  lÞ þ E1 h1 4 D1 E3 h3 4 D3 ðA:15Þ

Acknowledgement k 13 ¼ The authors are grateful for the financial and technical support of the Center for Microelectronics Assembly

1 h1 0 ðL  lÞ; 2 D1

k 14 ¼

1 h3 0 ðL  lÞ 2 D3 ðA:16Þ

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

Appendix B. Sample discretized equations1 Governing equations2:

Forces: n X T 1 ðiÞ ¼ ‘ s1 ðkÞ



 k6 4k 6 6k 6 k 3 s1 ði  2Þ s ði  4Þ  s ði  3Þ þ  1 1 ‘4 ‘4 ‘4 ‘2     4k 6 2k 3 k6 k3 þ  4 þ 2 s1 ði  1Þ þ 4  2 þ k 1 s1 ðiÞ ‘ ‘ ‘ ‘ þ k 2 s2 ðiÞ þ k 8 V 1 ðiÞ þ k 9 V 3 ðiÞ ¼ 0

883

ð1Þ

k¼i

Bending moments: M 1 ðiÞ ¼  ð28Þ

D3 4D3 k 8 s1 ðiÞ þ k 9 s2 ðiÞ þ 4 V 1 ði  4Þ  4 V 1 ði  3Þ ‘ ‘     6D3 D1 4D3 2D1 V V 1 ði  1Þ  ði  2Þ þ  þ þ 1 ‘4 ‘2 ‘4 ‘2   D3 D1 D2 þ 4  2 þ D6 V 1 ðiÞ  2 V 3 ði  2Þ ‘ ‘ ‘   2D2 D2 ð38Þ þ 2 V 3 ði  1Þ þ  2 þ D7 V 3 ðiÞ ¼ 0 ‘ ‘

n X h1 T 1 ðiÞ þ ‘ V 1 ðkÞ 2 k¼i

ð7Þ

Axial stresses: rmxj ðiÞ ¼

T j ðiÞ 6M j ðiÞ  hj h2j

ð49Þ

von Mises equivalent stresses (plane stress): rmeq j ðiÞ ¼ ½rmxj ðiÞ2 þ pm ðiÞ2  rmyj ðiÞ  pm ðiÞ þ 3sm ðiÞ2 1=2 ð50Þ Out-of-plane stresses: rmzj ðiÞ

¼ mj ½rmxj ðiÞ þ pm ðiÞ  aj Ej DT

ð51Þ

Boundary conditions: s1 ð1Þ ¼ s2 ð1Þ ¼ V 1 ð1Þ ¼ V 3 ð1Þ ¼ s1 ðnÞ ¼ s2 ðnÞ ¼ V 1 ðnÞ ¼ V 3 ðnÞ ¼ 0 k6 3k 6  3 s1 ðn  3Þ þ 3 s1 ðn  2Þ ‘ ‘   3k 6 k 3 s1 ðn  1Þ ¼ ða1  a2 ÞDh   ‘ ‘3

References ð40Þ

ð41Þ

  D3 3D3 3D3 D1 V 1 ð2Þ V ð4Þ þ V ð3Þ þ  þ 1 1 ‘ ‘3 ‘3 ‘3 D2 V 3 ð2Þ ¼ 0 ð43Þ þ ‘   D3 3D3 3D3 D1 V 1 ðn  1Þ V ðn  3Þ  V ðn  2Þ þ  1 1 ‘ ‘3 ‘3 ‘3 D2 V 3 ðn  1Þ ¼ 0 ð43Þ  ‘ 

k 6 ½s1 ð3Þ  2s1 ð2Þ þ s1 ð1Þ=‘2 þ k 7 ½s2 ð3Þ  2s2 ð2Þ n X ½k 11 s1 ðiÞ  k 12 s2 ðiÞ þ k 13 V 1 ðiÞ þ s2 ð1Þ=‘2 þ i¼1

þ k 14 V 3 ðiÞ‘ ¼ ða1  a3 ÞDhðL0  lÞ

ð47Þ

Peel stresses: i 6 ðn þ 1Þ=2: p1 ðiÞ ¼ ðV 1 ði þ 1Þ  V 1 ðiÞÞ=‘

ðB:1Þ

i P ðn þ 1Þ=2: p1 ðiÞ ¼ ðV 1 ðiÞ  V 1 ði  1ÞÞ=‘

ðB:2Þ

1 The equation numbers match the original equation in the text before discretization. 2 Note that ‘‘+’’ and ‘‘’’ signs in the (±) refer to the divisions in the left and right side of i ¼ nþ1 2 , respectively.

[1] Clech JP. Solder reliability solutions: A PC-based designfor-reliability tool. In: Surface mount international proceedings, San Jose, CA, US, September 8–12, 1996. p. 136– 51. [2] Ghorbani HR, Spelt JK. Interfacial thermal stresses in trilayer assemblies. Trans ASME J Electron Pack, in press. [3] Hall PM. Forces, moments and displacements during thermal chamber cycling of leadless ceramic chip carriers soldered to printed boards on components, hybrids and manufacturing technology. IEEE Trans 1984;7(4):314–27. [4] Jih E, Pao YH. Evaluation of design parameters for leadless chip resistors solder joints. Trans ASME J Electron Pack 1995;117:94–9. [5] Nicewarner E. Historical failure distribution and significant factors affecting surface mount solder joint fatigue life. Solder Surf Mount Technol 1994;17:22–9. [6] Pao YH, Eisele E. Interfacial shear and peel stresses in multilayered thin stacks subjected to uniform thermal loading. Trans ASME J Electron Pack 1991;113:164–72. [7] Qi Y, Ghorbani HR, Spelt JK. Accelerated thermal cycling of tin–lead and lead-free solder joints. In: SMTA international annual conference, Chicago, IL, USA, 2003. p. 735– 40. [8] Ru CQ. Interfacial thermal stresses in bimaterial elastic beams: modified beam models revisited. Trans ASME J Electron Pack 2002;124:141–6. [9] Suhir E. Calculated thermally induced stresses in adhesively bonded and soldered assemblies. In: International symposium on microelectronics proceedings, Atlanta, GA, October 1986. p. 383–92. [10] Suhir E. An approximate analysis of stresses in multilayered elastic thin films. J Appl Mech 1988;55:143–8. [11] Suhir E. Predicted thermal stresses in a bimaterial assembly adhesively bonded at the ends. J Appl Phys 2001;89(1): 120–9.

884

H.R. Ghorbani, J.K. Spelt / Microelectronics Reliability 46 (2006) 873–884

[12] Suhir E. Analysis of interfacial thermal stresses in a trimaterial assembly. J Appl Phys 2001;89(7):3685–94. [13] Suhir E. Bimaterial assembly with a low modulus bonding later at the ends. J Appl Phys 2003;93(6):3657–61. [14] Timoshenko SP, Goodier JN. Theory of elasticity. New York: McGraw-Hill; 1970. [15] Tummala RR, Rymaszewski EJ. Microelectronics packaging handbook. New York: Van Nostrand Reinhold; 1997.

[16] Vandevelde B, Christiaens F, Beyne E, Roggen J, Peeters J, Allaert K, et al. Thermomechanical models for leadless solder interconnections in flip chip assemblies. IEEE Trans Comp Pack Manuf Tech Part A 1998;21(1):177–85. [17] Compaq Visual Fortran 6.6, Compaq Computer Corporation, Houston, Texas, USA, 2000. [18] ANSYS 8.0, ANSYS Inc., Canonsburg, Pennsylvania, USA, 2003.