Best Practices for Thermal Modeling in Microelectronics with Natural Convection Cooling: Sensitivity Analysis

Mamadou Kabirou Touré^*, Papa Momar Souaré, Julien Sylvestre

Institut Interdisciplinaire d’Innovation Technologique (3IT), Department of Mechanical Engineering, University of Sher-brooke, Sherbrooke, Canada.

DOI: 10.4236/jectc.2021.102002                

Abstract

A detailed sensitivity study was carried out on various key parameters from a high precision numerical model of a microelectronic package cooled by natural convection, to provide rules for the thermal modeling of microelectronic packages subjected to natural convection heat transfer. An accurate estimate of the junction temperature, with an error of less than 1˚C, was obtained compared to the experimental data for the vertical and horizontal orientations of the test vehicle in the JEDEC Still Air configuration. The sensitivity study showed that to have an accurate estimate of the temperature, the following elements should be present in the thermal model: radiation heat transfer in natural convection cooling; a computational fluid dynamics analysis to find realistic convection coefficients; detailed models of the high conductivity elements in the direction of the heat flow towards the environment; and finally precise values for the thicknesses of layers and the thermal properties of the substrate and the printed circuit board.

Keywords

Computational Fluid Dynamics, Computational Heat Transfer, Microelectronic Packaging, Natural Convection, Radiation, Thermal Analysis, Thermal Management

1. Introduction

The performance of microelectronic systems deteriorates rapidly when their temperature exceeds a certain limit. For some JEDEC standards used in qualification and reliability tests, operating temperature limits are set, for example, to 125˚C for processors and 85˚C for memory chips. The purpose of thermal management is to maintain throughout the equipment a temperature distribution with limited variations above the recommended limits. The advancement in IC package development driven by the increase of transistor density and miniaturization makes the task of the packaging designer more challenging.

The key to successful thermal management is the ability to obtain complete and accurate temperature data under realistic operating conditions. Collecting temperature data by direct measurements (e.g. thermocouples) is the most commonly used method, but it is limited by the number of measurement points and the small size of the components. Infrared Thermal Imaging (IR) alleviates some of these problems by providing complete external temperature mappings. However, this approach can be time-consuming and can only be applied after the design is complete and the parts have been fabricated. It is therefore important to have thermal models and numerical simulation tools that can provide the temperature profile of the chips to achieve high-quality designs before prototyping or manufacturing.

Several studies have been conducted to thermally model microelectronic packages [1]. Analytical methods can be used in simple cases such as isothermal or adiabatic flat plates in horizontal or vertical orientations, but these analytical methods have obvious limitations for complex geometries or thermal boundary conditions [2]. Numerical methods are frequently used for thermal simulations of electronic packages [3]. Most of these numerical methods use the finite element method [4], the finite volume method [5], or the finite difference method [6]. Many of these studies address only the conduction in the solid and assume or estimate the convection by empirical methods [7] [8], which are generally obtained by only considering the average temperature difference between the surfaces and the ambient environment. The variations of the cooling airflow velocity due to geometric changes are often neglected because it is difficult to take them into consideration. In addition, the radiation is most of the time neglected for low power applications [9], although there is little published information about the level of power, or about the temperature difference between surface and ambient environment, such the radiation has to be taken into account.

The objective of this study was to perform a sensitivity analysis on our coupled (conjugate) methodology between a highly accurate conduction model and a computational fluid dynamics (CFD) model, which was reported in [10] [11] to be highly accurate when compared to temperature measurements on actual hardware. The sensitivity analysis was performed using a one-factor-at-a-time approach to evaluate the influence of model features or parameter values on the accuracy of the simulation results. Our aim is to identify the key factors to consider for accurate thermal modeling in microelectronics packages with natural convection cooling.

We present in this paper a method that enables an accurate prediction of the temperature field in a microelectronic package mounted on a printed circuit board (PCB) when it is cooled by natural convection to a thermal steady state. This configuration allows the estimate of the standard junction-to-ambient (θJA) thermal resistance in natural convection, which is useful to evaluate the relative thermal performance of microelectronic packages [12]. Our method is based on the coupling of a conduction model and of a fluid flow model. It performs well because of the high level of geometrical details in the conduction model, and because of the high resolution of the airflow and heat transfer models around the package. The accuracy of the conduction model within the solid domain is improved by a fine spatial discretization, while the precision of the convection model in the fluid domain results from a good spatial discretization of the thermal boundary layer. Using this high precision model, we evaluate the impact of the variation of several parameters on the temperature field in the package assembly, in comparison to experimental measurements. These parameters include the lid conductivity; the thermal interface material (TIM) thickness and conductivity; the underfill conductivity and the constriction resistance in the calculation of the equivalent resistance of the underfill layer with the interconnects (flux tube model); the laminated substrate copper distribution and thickness variability; the substrate core and dielectric conductivity; the presence of vias in the substrate dielectric and in the substrate core; the ball grid array (BGA) conductivity; the PCB thickness variability and emissivity; the PCB dielectric material in-plane and out-of-plane conductivity; the presence of thermal vias in the PCB; the use of empirical formulas for convection modelling; the consideration of radiation heat exchange and the still air chamber inclination. This study highlights the most important features to take into account for the construction of precise thermal models.

The paper is divided into eight sections. A review of the literature on the sensitivity studies around thermal issues in microelectronics packaging is presented in Section 2. Section 3 describes the detailed methodology of the two-way coupling between the conduction model and the fluid flow model. The acquisition of experimental data for the validation of the methodology is presented in Section 4. The results from the detailed numerical models are discussed in Section 5. Section 6 presents the validation of the numerical models, by comparing the results of the simulations from Section 4 with the experimental data from Section 5. The sensitivity analysis is discussed in Section 7 (key parameters) and in Section 8 (results and discussion).

2. Related Work

Sensitivity studies have been carried out by some authors for thermal issues in microelectronic packaging, but the list of parameters studied is not as exhaustive as what is targeted in the present study. Reference [13] presents an efficient lid design for a 3D stacked flip-chip package that increases the contact surface of the TIM between the top chip and the lid for a better heat conduction transfer. A sensitivity analysis is presented for the thermal resistance of five key structures (TIM, underfill/C4 composite between the two chips and the lower chip to substrate, substrate and seal band), to identify the most critical structure for the thermal performance of the package. This study concluded that the chip-to-chip thermal resistance in the 3D package was most critical and could be reduced by improving the thermal conductivity of the underfill material, increasing the C4 density and by decreasing their height. Therefore, an accurate assessment of the chip-to-chip thermal resistance and of the TIM thermal resistance is required for accurate temperature predictions in similar 3D packages. However, in this study only the conduction was treated with a compact model, which has the advantage of producing a fast approximate solution with very few details but leaves out many important details, such as the effect of radiation, for instance. The validation of the model was performed by comparing the experimental thermal resistivity value from sensors at different positions in the chip (center, near edge and edge) and the data from numerical model. However, the numerical model was solved by setting a uniform convection coefficient on the top of the lid (1000 W/m²K) and the bottom of the substrate (10 W/m²K), neglecting edge effects and fluid velocity fluctuations on surfaces.

The heat dissipation of a microelectronic device mounted on a four-layer PCB cooled by natural convection was studied in references [14] and [8]. In the first study [14], a finite element conduction model was used, with empirical formulas for estimating the fluid heat transfer coefficients. In the second study [8], a CFD model was used to solve the conjugated solid/fluid heat transfer problem, with sensitivity analysis on the use of empirical model for the convective heat transfer, and on the PCB features (such as the number of copper layers, thermal vias and thermal BGA). These results showed that using a CFD methodology leads to better predictions compared to using empirical model to estimate fluid heat transfer coefficients, and that the PCB structures such as thermal vias, thermal BGAs and copper layers have a considerable influence on the equilibrium temperature field.

Reference [15] has presented a sensitivity analysis of conjugate heat transfer for cooling electronic devices for heat sink characterization and system-level thermal design optimization. The parameters of this study comprised: the distance between the top of the heat sink and the enclosure (the gap), the angle of the flow (0˚ to 60˚), the conductivity of the heat sink (25 to 430 W/mK) and the fluid flow velocity. The goal of this study was to minimize the thermal resistance of a heat dissipating component. These results showed that the minimal thermal resistance of the heat sink (8.1 K/W) was obtained with a minimum gap and a forced convection flow angle of 30˚ at a power of 9 W. The authors also showed that the resistance could be further reduced using a heat sink material with higher thermal conductivity.

Reference [16] has studied the conjugate heat transfer of a computer cooled by liquid. The study showed the variation of the maximum temperature with variations in the space between the chip and the heat sink, as well as the operating power between 15 W and 40 W. Hence, the maximum temperature was decreased by 38% using water in the enclosure, relative to the air cooling. Increasing the space between chip and heat sink reduced the temperature. However, beyond 50 mm, this increase in space had an impact of less than 1% on the temperature variation.

Table 1 summarizes the impact of the package thermal resistance (junction to cooling fluid) by varying by +10% the nominal value of continuous parameters and by taking into account or not in the simulations the discrete parameters, for the references cited above.

These studies show that some parameters (such as the number of copper layers, and the presence of thermal vias) have a much larger impact on the estimated thermal resistance than other parameters. While leaving important details out of the simulations can reduce accuracy considerably, adding unnecessary details consume more computing resources, sometimes with minimal benefits on accuracy. Therefore, there is a need for identifying key features for the construction of precise thermal models, which is the main goal of our study.

Table 1. Summary of the impact of the variation of some parameters in the literature.

3. Methodology

Our study uses the finite element method (FEM) for the heat conduction inside the package and the finite volume method (FVM) for the air convection in the volume surrounding the package. The FEM enables more precise solutions by increasing the order of the elements and by refining the mesh locally, while the FVM is robust for solving the non-linear conservation laws which appear in fluid transport problems.

The full heat transfer problem was simulated using a two-way (conjugate) coupling methodology between a FEM conduction model and an FVM fluid model. The detailed numerical models are presented in Section 5. The interaction between the two models was implemented through exchanges of temperature fields and convection coefficients at the boundaries between the solid domain and the fluid domain.

The two-way coupling methodology involved iterating through the following steps: first, the heat transfer in the conduction model was solved with Ansys APDL [17] by applying uniformly an initial estimate of the heat transfer coefficients (obtained for instance from empirical formula [18] ) at each solid-fluid interface. Second, the temperature field at the solid-fluid interfaces calculated in the conduction model was applied as a boundary condition in the FVM simulation (see Equation (2)), which was solved with Ansys Fluent [19] to determine the velocity and pressure fields around the package. For each solid-fluid interface FVM element, the total heat transfer coefficient h_tot was calculated using

$h_{\text{tot}}=\frac{q_{\text{wall}}}{T_{\text{element}}-T_{\text{ambient}}}$ (1)

where q_wall is the total heat flux including convection and radiation, T_element is the temperature of each element at the solid-fluid interface and T_ambient is the fluid temperature at infinity. The heat transfer coefficients were transferred to the FEM model (see below). Third, the heat transfer coefficients were applied to all elements in the conduction model, which was solved again. Several iterations of the second and third steps were executed, until the calculation was terminated after the third step when the difference of the maximum temperature between two successive iterations was less than 0.2 K. Figure 1 shows typical convergence data for two orientations of a package (horizontal and vertical).

The FEM and FVM meshes were non-conformal. We found that an effective way of transferring the temperature field from FEM to FVM was to fit the temperature distribution on each surface with a Gaussian function. The fitting equation of the FEM temperature field for each surface was

$T\left(x,y\right)=A_i{\text{e}}^{-\left(\frac{{\left(x-x_{0i}\right)}^2}{2{\sigma }_{xi}^2}+\frac{{\left(y-y_{0i}\right)}^2}{2{\sigma }_{yi}^2}\right)}+B_i$ (2)

where $A_i$, $B_i$, $x_{0i}$, $y_{0i}$ and ${\sigma }_{xi}$ and ${\sigma }_{yi}$ are the fit constants for each surface i (the top of the lid and the top and the bottom of the PCB).

The FVM convection coefficients were transferred to the FEM by projecting

Figure 1. Convergence of maximum temperature versus iteration number, for 10 W power consumption, for two configurations: horizontal (solid line) and vertical (dotted line).

the convection coefficients from the FVM mesh to the FEM mesh at the solid-fluid interfaces using linear interpolation, by using a Delaunay triangulation of the FVM data, and performing linear barycentric interpolation at each FEM nodes.

4. Temperature Measurements

The validation of our calculation methodology was performed using experimental data obtained under the Integrated Circuits Thermal Test Method Environmental Conditions-Natural Convection (Still Air) (JEDEC JESD51-2) [12] standard test in natural convection heat transfer in a still-air ambient environment. The still air chamber had a volume of 0.028 m³. The 6.35 mm thick plastic walls were fabricated from transparent acrylic with low thermal conductivity and high emissivity. The chamber included a fixed thermocouple probe for measuring the internal air temperature. The standard JESD51-2 recommends an increase of the box volume if during the tests the temperature in the chamber increases by 10% or if the power exceeds 3 W, which was the case in our tests. This recommendation is justified by the fact that there could be convection on the exterior surfaces of the box impacting the thermal characterization. To verify this, type-T thermocouples were bonded internally to the corners of the box and at the center of each wall of the box to measure an eventual heat flux outside the box for heating powers above 3 W. Based on these results, it was decided to increase the domain of the CFD model to take into account the air volume outside of the box.

The test vehicle consisted of a test module mounted on a JEDEC standard PCB of dimensions 127 × 139.7 mm² by an array of BGA solder balls. The test module was composed of a silicon chip of dimensions 12.57 × 12.57 mm², flip-chip mounted on an organic substrate of size 55 × 55 mm². The interconnections between the chip and the substrate were a matrix of solder bumps (C4) filled with an underfill polymeric material to increase the mechanical integrity of the solder bumps. The module was stiffened by a 1 mm-thick copper lid that also served as a first-level heat sink. The lid was thermally coupled to the chip by a thin thermal interface material (TIM, approximately 15 µm thick) and attached to the substrate with a silicone adhesive (seal band). There was an excess of the TIM during the attachment process that filled the gap between the lid and the underfill (Figure 2).

The silicon chip had nine heating elements and six in-situ resistance temperature detectors (RTD). These RTD were located at different positions on the chip: one at each corner and two at the center of the chip (Table 2). To better capture convection on external surfaces and thus validate our numerical model, eighteen type-T thermocouples were bonded to some exterior surfaces (top and bottom of the PCB and top of the lid). Their positions are listed in Table 2.

5. Numerical Models

In this paper, numerical simulations were performed using a specialized cloud software infrastructure (PACK [20] ). Pack is high performance numerical simulation software used as an interface to other commercial numerical simulation software (such as Ansys APDL and Ansys Fluent). PACK allows the efficient

Figure 2. Optical micrograph of a cross section of the vehicle test, showing lid, die, substrate, BGA and PCB with a focus on C4, underfill and TIM bond line.

Table 2. Thermocouples locations.

generation of parametric numerical models that use advanced pre- and post-processing capabilities using an object-oriented programming model.

5.1. FEM Model

The PACK infrastructure was used to build the conduction model with the Ansys APDL finite element software [17]. In this model, we chose to use the finite element method especially for the accuracy and complexity of the geometry of our test module (Figure 3). The conduction model represents the test module with the PCB, as explained in Section 3.

The properties of the materials were assigned individually to each element according to the position of the element in a given volume. The C4 and BGA were modeled with a homogenized layer with equivalent thermal properties. The equivalent out-of-plane resistance of these layers was the sum of two resistances in parallel, for the interconnection and the surrounding material. In-plane, it was assumed to equal the surrounding material thermal resistance (underfill for the C4 and air for the BGA). The thermal resistances of the C4 and BGA interconnections was calculated with the flux tube formula [21], given by

$R_{\text{ball}}=\frac{{\left(1-\epsilon \right)}^{1.5}}{2d_{\text{Ball}}}\cdot \left(\frac{1}{{\lambda }_1}+\frac{1}{{\lambda }_2}\right)+\frac{4h_{\text{Ball}}}{{\lambda }_{\text{Ball}}\pi d_{\text{Ball}}^2}$, (3)

with,

$\epsilon =\sqrt{\frac{\pi d_{\text{Ball}}^{2}N}{4A}}$ (4)

where ${\lambda }_1$, ${\lambda }_2$ and ${\lambda }_{\text{Ball}}$ are the thermal conductivities of the top layer, bottom layer, and the C4 or BGA, h_Ball is the C4 or BGA height, d_Ball is the C4 or BGA diameter, N is the number of C4 or BGA and A is the total area of the chip or substrate. Equation (3) is valid for 0 < ε ≤ 0.3. Table 3 summarizes the dimensions and thermal properties of the materials used in the different layers.

Figure 3. Description of the FEM Model with a focus on the specific mesh of 3,324,148 elements and 3,402,667 nodes.

Table 3. Material properties.

One of the key elements in the realization of the high precision conduction model was the detailed modeling of the organic substrate. Substrate layer design files were used to model accurately the copper distribution and the vias of each layer in the organic substrate by superimposing the substrate design file and the mesh. In this approach, a rectangular background grid (pixels) was constructed based on black and white images of each layer of the organic substrate. The pixel count within each finite element is used to calculate the effective isotropic material properties based on the local concentration of copper and dielectric. This, effectively, forms a map of the material properties across each layer of the organic substrate [21]. The local thermal conductivity properties ${\lambda }_{eff,n}$ of each finite element in the substrate were calculated using

${\lambda }_{eff,n}={\beta }_n\cdot {\lambda }_{\text{Copper}}+\left(1-{\beta }_n\right)\cdot {\lambda }_{\text{Dielectric}}$ (5)

where ${\beta }_n$ is the area fraction of element n covered by copper, ${\lambda }_{\text{Copper}}$ and ${\lambda }_{\text{Dielectric}}$ are the thermal conductivity of copper and dielectric material. This method allows for a more accurate representation of the organic substrate material distribution than the standard homogenization approach. Each layer of the PCB was modeled with homogeneous layers and the thermal vias were explicitly represented. The heat flux was applied uniformly at the bottom of the die as a source term. The simulations were performed in the steady state.

5.2. FVM Model

Ansys Fluent was used to build the FVM model [19]. The FVM model was composed of the air inside the chamber, the chamber itself, and the air outside the still air chamber (Figure 4). The dimensions of the outer domain of the still air chamber were taken large enough to put the temperature in the laboratory as boundary conditions (0.61 m away of each side of the still air chamber except the face that sits on the lab table, which had a gap of 25.4 mm).

A surface-to-surface radiation model [19] for the resolution of radiation heat transfer as well as the k-ε turbulent model for the transport problem [19], were chosen. The resolution of the energy equation was done in our model in conjunction with the k-εturbulent model derived from the Navier-Stokes equations [19]. For more details on numerical models and hypotheses, see reference [10].

6. Comparison of Numerical Results and Experimental Data

The validation of our numerical model was done by comparing the results of the simulations with the experimental data (Figure 5). The comparison was performed on the diagonal line of three surfaces (the top of the lid, the die junction and the bottom of the PCB) and for two orientations (horizontal and vertical).

The results presented in Figure 5 show that our two-way coupled methodology provides accurate estimates of the temperature fields for both the vertical and the horizontal orientations of the test vehicle. For more details on the validation of our numerical model, see reference [11].

7. Sensitivity Study

In order to contribute to the definition of best practices in the modeling of packages cooled by natural convection, a sensitivity study was performed on key parameters of our numerical model in order to evaluate their impact on the accuracy of the simulation results. These parameters can be classified into two types: continuous parameters where the impact of a variation of +10% of their nominal value was evaluated (e.g. variations in thickness and thermal conductivity) and discrete parameters where the impact of their taking into account or not in the simulations was evaluated (e.g. the presence of vias and the use of copper distributions). See Table 4 for the full list.

Figure 4. Description of the FVM model, excluding the volume outside of the chamber (box).

Figure 5. Comparison of the simulated temperature (round markers) and measurements (square makers) in the horizontal ((a), (c) and (e)) and vertical ((b), (d) and (f)) orientations at the die junction (a and b), the top of the lid (c and d) and the bottom of the PCB ((e) and (f)).

The “Underfill flux tube” parameter refers to taking into consideration or not the heat flux contraction in the calculation of the equivalent resistance of the underfill layer (Equation (3)). “Laminated substrate copper distribution” refers to the modeling of the copper distribution on the thermal properties of the substrate (Equation (5)). The parameters “Substrate core PTH”, “Substrate dielectric vias”, and “PCB thermal vias” evaluate the taking into account the presence of these structures in numerical models. When not included, a homogenized layer with uniform thermal properties was modeled. The use of empirical formulas for convection modeling and the consideration of radiation heat exchange are evaluated in the “CFD” and “Radiation” parameters, respectively. Finally, the impact of still air chamber inclination is also evaluated in the “Still air chamber inclination” parameter. The impact of all these parameters on the temperature was evaluated at different locations (junction, top lid, bottom PCB and corner PCB) and in the two configurations (horizontal and vertical orientations).

Table 4. List of parameters.

8. Results and Discussions

Table 5 summarizes the sensitivity study on the parameters for 10 W of applied power. This table presents the variation of the temperature at different locations relative to the temperature at the same locations when all parameters are at their nominal value. For discrete parameters, the variation is the absolute value of the temperature difference when each discrete parameter is deactivated. For continuous parameters, the variation is the average of the absolute value of the temperature difference for each parameter set at −10% and at +10% of its nominal value.

The parameters were classified according to four categories: those that we consider negligible (maximum variation of less than 1% on the temperature at either the junction of the chip, the top of the lid, the bottom of the PCB, or the corner of the PCB); influential (variations between 1% and 5%) and critical (variations of more than 5%). Three parameters have been identified as critical: “Radiation”, “CFD” and “substrate dielectric vias”. In many thermal studies in microelectronics, radiation heat transfer is neglected because of the underestimation of its importance in the cooling of electronics. Some authors neglect its

Table 5. Summary of sensitivity study.

importance in conditions where the difference in temperature is small enough to simplify their analysis. In our setting, this would lead to errors of more than 40% at the junction of the chip and up to 80% at the corners of the test card because of its high emissivity (0.85). For a more detailed investigation on the importance of radiation heat transfer in simulations, the contribution of radiation compared to that of convection in natural convection heat transfer as a function of the power dissipated in the horizontal and vertical configurations is shown in Figure 6 and Figure 7.

It can be seen in both the horizontal and vertical configurations that the radiation heat transfer is responsible for more than half of the total heat flux dissipated, especially at the lower power levels, where up to 60% of the heat can be dissipated by radiation. Radiation heat transfer is therefore very important in natural convection cooling and must not be neglected.

Figure 6. Contribution of radiation (dotted line) and convection (solid line) in natural convection heat transfer according to the power dissipated, in the horizontal configuration.

Figure 7. Contribution of radiation (dotted line) and convection (solid line) in natural convection heat transfer according to the power dissipated, in the vertical configuration.

Also, using empirical formulas to calculate convection and radiation heat transfer coefficients leads to errors of more than 8% at the junction of the chip and up to 25% at the corners of the test card in the vertical configuration. The approximations introduced by the empirical relations (for example, neglecting the non-uniformity of the velocity on the interfaces) can explain these errors. In this regard and in order to better understand the differences in temperatures agreement between the CFD and the empirical relations method, a comparison between their corresponding heat transfer coefficients has been performed in references [10] [11]. This showed that using CFD to find realistic convection coefficients is necessary to provide an accurate estimate of the temperature in the microelectronic packages.

In our previous studies [10] [11], we have shown that in natural convection most of the heat flow passes through the PCB. The vias in the substrate dielectric layers act as thermal bridges between the substrate and the PCB. Not including them in the thermal simulations leads to significant errors, about 5% at the junction of the chip. The presence of the vias in the substrate dielectric layers results in an increase by 2.25 times of the effective conductivity of the dielectric.

Four parameters were found to be influential (variation between 1% and 5%). Like the vias in the dielectric layers, the PTH in the core layer of the substrate also plays the role of thermal bridges. Their presence leads to an increase by 3 times on the effective conductivity of the substrate core. The thermal resistances of the substrate and of the PCB also had an impact on the validity of the simulations results. As shown above, radiation is important in natural convection heat transfer, as confirmed by the impact of the PCB Emissivity. It can be concluded that a good characterization to find their precise values is necessary for these four parameters.

Negligible parameters (variation less than 1%) include the copper distribution in the substrate, the inclination of the still air chamber, the thermal vias in PCB (~0.7% of the PCB thermal resistance), the thermal conductivity of the dielectric layers, core layer, lid, underfill, BGA layer, layers of FR4 in the PCB, and the thickness and thermal conductivity of the TIM, as well as the emissivity of the still air chamber (which has a temperature that is close to ambient). Not taking them into account or varying their nominal value up to 10% did not affect much the accuracy of our simulations.

9. Conclusions

An accurate estimate of the temperature field, with an error of less than 1˚C compared to experimental measurements for the vertical and horizontal orientations of the test vehicle, was obtained through a conjugate methodology that combines a high-resolution conduction model and a CFD-based model to more precisely simulate the natural convection and radiation heat transfers in microelectronic packages. A sensitivity study was performed using this high precision model for various parameters. The sensitivity study has revealed that in order to have an accurate solution for the temperature field, the following features must be incorporated into the model:

Using CFD to find realistic convection coefficients: the use of empirical formulas to calculate heat transfer coefficients can lead to large errors (8% at the junction of the chip and 25% at the corners of the test card in the vertical configuration). These error results are mainly from the assumptions of a uniform velocity field, which is not observed when edge effects and fluid flow in a confined box are present.

Including radiation heat transfer: radiation can be responsible for more than half of the total heat flux dissipated in the studied natural convection cooling configuration, even at low power. In our case, almost all of the radiation heat transfer passes through the PCB (98%) because of its high emissivity (0.85), which is 17 times greater than that of the lid (0.05), and its surface which is 8 times larger than the surface of the lid. Given the importance of radiation in natural convection heat transfer, proper knowledge of the precise values of the emissivity of the PCB and the still air chamber walls is also important.

Model the vertical interconnections along the heat dissipation path: the heat being dissipated mainly through the PCB in the two orientations (86% for horizontal, 84% for vertical), the vias in the dielectric layers as well as the PTH in the substrate core play the role of thermal bridges and can increase significantly the substrate out-of-plane conductivity.

Obtain the correct thicknesses and thermal properties for the substrate and the PCB: the thermal resistance of the substrate and of the PCB can have a significant impact on the simulation results. Variations that may come from manufacturing tolerances should be properly understood.

Other parameters could vary by up to ±10% or not be taken into account and not affect the calculated temperatures by more than 1%. These parameters include the copper distribution mapping in the substrate, the inclination of the still air chamber, the thermal vias in PCB, the thermal conductivity of the dielectric layers, the core layer, lid, underfill, BGA layer, layers of FR4 in the PCB, and the thickness and thermal conductivity of the TIM, as well as the emissivity of the still air chamber.

Acknowledgements

The authors thank Éric Duchesne, Benoît Foisy, and Michel Levesque from the IBM Corporation for useful discussions, technical support, as well as the provision and maintenance of the laboratory equipment.

Funding

This project was financially supported by the IBM Corporation, the Natural Sciences and Engineering Research Council of Canada and Prompt.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Garimella, S.V. (2006) Advances in Mesoscale Thermal Management Technologies for Microelectronics. Microelectronics Journal, 37, 1165-1185. https://doi.org/10.1016/j.mejo.2005.07.017
[2]	Chiang, T.-Y., Banerjee, K. and Saraswat, K.C. (2002) Analytical Thermal Model for Multilevel VLSI Interconnects Incorporating via Effect. IEEE Electron Device Letters, 23, 31-33. https://doi.org/10.1109/55.974803
[3]	Minkowycz, W.J. and Sparrow, E.M. (2000) Advances in Numerical Heat Transfer. Vol. 2, CRC Press, Taylor & Francis, Boca Raton.
[4]	Rzepka, S., Banerjee, K., Meusel, E. and Hu, C.M. (1998) Characterization of Self-Heating in Advanced VLSI Interconnect Lines Based on Thermal Finite Element Simulation. IEEE Transactions on Components, Packaging, and Manufacturing Technology: Part A, 21, 406-411. https://doi.org/10.1109/95.725203
[5]	Mizunuma, H., Lu, Y. and Yang, C. (2011) Thermal Modeling and Analysis for 3-D ICs with Integrated Microchannel Cooling. IEEE Transactions on Computer-Aided Design of Integrated Circuits and Systems, 30, 1293-1306. https://doi.org/10.1109/TCAD.2014.2334321
[6]	Liu, Z., Swarup, S., Tan, S.X., Chen, H. and Wang, H. (2014) Compact Lateral Thermal Resistance Model of TSVs for Fast Finite-Difference Based Thermal Analysis of 3-D Stacked ICs. IEEE Transactions on Computer-Aided Design of Integrated Circuits and Systems, 33, 1490-1502. https://doi.org/10.1109/TCAD.2014.2334321
[7]	Chen, W.-H., Cheng, H.-C. and Shen, H.-A. (2003) An Effective Methodology for Thermal Characterization of Electronic Packaging. IEEE Transactions on Components and Packaging Technologies, 26, 222-232. https://doi.org/10.1109/TCAPT.2002.806180
[8]	Hwang, P.-W., Cheng, H.-C., Fang, J. and Li, J.-H. (2007) CFD-Based Thermal Characterization of Board-Level Microelectronic Devices under Natural Convection Cooling. 2007 International Microsystems, Packaging, Assembly and Circuits Technology, Taipei, 1-3 October 2007, 79-82. https://doi.org/10.1109/IMPACT.2007.4433572
[9]	Hwu, F.-S., Sheu, G.-J. and Chen, J.-C. (2006) Thermal Modeling and Performance of LED Packaging for Illuminating Device. Proceedings of the 6th International Conference on Solid State Lighting, 6337, 327-333. https://www.spiedigitallibrary.org/conference-proceedings-of-spie/6337/1/Thermal-modeling-and-performance-of-LED-packaging-for-illuminating-device/10.1117/12.678949.short?SSO=1
[10]	Toure, M.K., Momar Souare, P., Allard, S., Foisy, B., Duchesne, E. and Sylvestre, J. (2018) CFD-Based Iterative Methodology for Modeling Natural Convection in Microelectronic Packages. 2018 19th International Conference on Thermal, Mechanical and Multi-Physics Simulation and Experiments in Microelectronics and Microsystems (EuroSimE), Toulouse, 15-18 April 2018, 1-8. https://doi.org/10.1109/EuroSimE.2018.8369862
[11]	Souare, P.M., Toure, M.K., Allard, S., et al. (2018) High Precision Numerical and Experimental Thermal Studies of Microelectronic Packages in Still Air Chamber Tests. 2018 7th Electronic System-Integration Technology Conference (ESTC), 18-21 September 2018, 1-8. https://doi.org/10.1109/ESTC.2018.8546346
[12]	JEDEC (2008) JEDEC Standard JESD51-2A, Integrated Circuits Thermal Test Method Environmental Conditions—Natural Convection (Still Air). JEDEC Solid State Technology Association, Arlington.
[13]	Sikka, K., Wakil, J., Toy, H. and Liu, H. (2012) An Efficient Lid Design for Cooling Stacked Flip-Chip 3D Packages. 13th InterSociety Conference on Thermal and Thermomechanical Phenomena in Electronic Systems, San Diego, 30 May-1 June 2012, 606-611. https://doi.org/10.1109/ITHERM.2012.6231484
[14]	Cheng, H.-C., Chen, W.-H. and Cheng, H.-F. (2008) Theoretical and Experimental Characterization of Heat Dissipation in a Board-Level Microelectronic Component. Applied Thermal Engineering, 28, 575-588. https://doi.org/10.1016/j.applthermaleng.2007.04.013
[15]	Dogruoz, M.B., Sathyamurthy, P. and Mathur, S. (2012) Sensitivity Analysis in Conjugate Heat Transfer for Electronics Cooling. 13th InterSociety Conference on Thermal and Thermomechanical Phenomena in Electronic Systems, San Diego, 30 May-1 June 2012, 1219-1224. https://doi.org/10.1109/ITHERM.2012.6231561
[16]	Gdhaidh, F.A., Hussain, K. and Qi, H.S. (2014) Numerical Study of Conjugate Natural Convection Heat Transfer Using One Phase Liquid Cooling. IOP Conference Series: Materials Science and Engineering, 65, 012012. https://doi.org/10.1088/1757-899X/65/1/012012
[17]	ANSYS (2013) ANSYS Mechanical APDL Basic Analysis Guide, Release 15.0. ANSYS Inc., Canonsburg.
[18]	Huetz, J. and Petit, J.-P. (1990) Notions de transfert thermique par convection. Techniques de l'ingénieur, A1540, 47. https://www.techniques-ingenieur.fr/base-documentaire/archives-th12/archives-plastiques-et-composites-tiaam/archive-1/notions-de-transfert-thermique-par-convection-a1540a/
[19]	ANSYS (2017) Ansys Fluent Theory Guide, Release 18.0. ANSYS Inc., Canonsburg.
[20]	Sylvestre, J. (2007) Integrated Modeling of C4 Interconnects. 2007 Proceedings 57th Electronic Components and Technology Conference, Sparks, 29 May-1 June 2007, 1084-1090. https://doi.org/10.1109/ECTC.2007.373932
[21]	Lau, J.H., Wong, C.P., Prince, J. and Nakayama, W. (1999) Electronic Packaging: Design, Materials, Process & Reliability. McGraw-Hill, New York.