FACTA UNIVERSITATIS Series: Electronics and Energetics Vol. 31, N o 4, December 2018, pp. 519-528 https://doi.org/10.2298/FUEE1804519J CONSIDERATIONS ON THE IMPORTANCE OF PROPER HEAT TRANSFER COEFFICIENT MODELING IN AIR COOLED ELECTRONIC SYSTEMS  Marcin Janicki, Agnieszka Samson, Tomasz Raszkowski, Tomasz Torzewicz, Andrzej Napieralski Department of Microelectronics and Computer Science, Lodz University of Technology, Poland Abstract. This paper illustrates, based on a practical example of a hybrid circuit, the influence of proper heat transfer coefficient modelling in air cooled electronic systems on the accuracy of thermal simulations. This circuit contains a transistor heat source and a set of temperature sensors. The measurements of their temperature responses are taken in natural convection and forced air cooling conditions. The experimental data provide the information necessary to estimate the local heat transfer coefficient values in heat source and temperature sensor locations. Moreover, the experiments rendered possible the fitting of parameters of an empirical heat transfer coefficient model for different surface temperature rise values and cooling air velocities, and hence allowed significant improvement of thermal simulation accuracy. Key words: thermal modeling, air cooling, heat transfer coefficient. 1. INTRODUCTION Temperature is the most important factor influencing the operation of all electronic systems and, at the same time, it is the principal cause of their failures [1]. Thus, thermal analyses have become an indispensable stage of electronic system design. Unfortunately, usually in simulations it is assumed that thermal model parameter values are temperature independent. However, in most real cases this assumption may not be true since, as shown in [2], the heat transfer coefficient values depend strongly on surface temperature rise and cooling air velocity. Obviously, there exist thermal models, mostly empirical, which allow taking into account such dependencies, however typically for flat plates these relations were derived assuming the uniformity of surface temperature or heat flux [3]-[5]. Thus, considering that in electronic circuits heat sources usually occupy a small part of their surface, large temperature gradients occur and consequently local heat transfer coefficient Received September 9, 2018 Corresponding author: Marcin Janicki Lodz University of Technology, 116 Żeromskiego Street 90-924 Lodz, Poland (E-mail: janicki@dmcs.pl)  520 M. JANICKI, A. SAMSON, T. RASZKOWSKI, T. TORZEWICZ, A. NAPIERALSKI values differ considerably in various locations rendering those standard relations useless in practice [6]-[7]. This problem is illustrated in this paper based on the example of a real hybrid circuit containing a transistor heat source and several thermistor temperature sensors. During the experiments, the heat source and sensor temperature values were measured both in natural and forced convection cooling conditions with variable air velocity. Local heat transfer coefficient values were estimated owing to the efficient coupling of the forward thermal solver with different inverse algorithms implementing chosen methods described in [8]- [10] and a simple algorithm proposed by the authors. Moreover, the experiments rendered possible the fitting of empirical model parameters. The fitting was carried out employing various swarm intelligence algorithms [11]-[13]. The following section of this paper will introduce in detail the investigated test hybrid circuit and provide the results of its temperature measurements. Then, the estimation results of local heat transfer coefficient values are presented and discussed. Finally, these values are fitted to an empirical model allowing the computation of local heat transfer coefficient at a chosen location for any surface temperature rise and cooling air velocity. 2. EXPERIMENTAL RESULTS 2.1. Hybrid circuit description The investigations presented in this paper concern a test hybrid circuit, whose layout is presented in Fig. 1, containing a Bipolar Junction Transistor (BJT) and five platinum thermistors with a positive thermal coefficient (PT1-5). All these components can be used as temperature sensors and the BJT can serve at the same time as a heat source. The PCB is manufactured in the insulated metal substrate technology and is made of a 1 mm-thick aluminum alloy. The board is a long but narrow rectangle with the devices located in the middle of its width so that the heat diffusion along the substrate gradually becomes quasi one-dimensional. Each time the distances between individual sensors are doubled. Initially all devices were calibrated using a cold plate. The BJT was calibrated for the emitter current values, which were used later on for heating, i.e. the currents ranging from 250 mA to 1000 mA with the step of 250 mA, and the thermistor current value was equal to 750 A. As expected, the measured device temperature characteristics were linear. The BJT base-emitter voltage sensitivity decreased with the emitter current from -1.86 mV/K at 250 mA to -1.51 mV/K at 1000 mA and the thermistor sensitivity was 96 V/K. 33 5 5 1045 20 160 40 BJT P T 1 P T 2 P T 3 P T 4 P T 5 wind direction Fig. 1 Layout of the hybrid circuit (dimensions in millimeters). Considerations on The Importance of Proper Heat Tranfer Coefficient Modeling in Air Cooled Electronic ... 521 2.2. Temperature measurement results The experiments consisted in the measurements of transistor heat source and sensor dynamic thermal responses for various levels of dissipated power with natural convection cooling and with forced convection cooling in the wind tunnel at different air velocities. The transistor temperature values were recorded with the transient thermal tester T3Ster manufactured by Mentor Graphics and the thermistor temperature values were registered with a custom measurement board designed by the authors. The dependence of BJT and sensor steady state temperature rise values on dissipated power in natural convection cooling conditions for the four values of emitter current used during the calibration is plotted in Fig. 2. As can be seen, this dependence is not perfectly linear and for the power dissipation of 7 W the temperature rise values of all thermistors are at least 20% lower than the ones, which would be obtained by projecting the values measured for the lowest power dissipation. This clearly indicates the existence of a non- linearity related to the change of the heat transfer coefficient value. This effect apparently is not so strong for the BJT, however as demonstrated in [14] it results from the increase of the package thermal resistance what counterbalances the change of cooling conditions. The transistor and sensor steady state temperature rise values obtained at the heating current of 1000 mA for different air velocities are shown in Fig. 3. The wind direction, indicated in Fig. 1 by the arrow, was chosen expressly in order to avoid the sensor heating by the air passing over the board. Looking at the figure it could be concluded that owing to the forced air movement the temperature rise values are noticeably reduced, but for the higher air velocities this decrease is less important. Moreover, the transistor temperature rise was reduced by one fourth, whereas for the sensor locations this reduction is more significant and it amounts from nearly 50% for the sensor closest the source to 80% for the furthest one. The results obtained for the other heating current values are very similar, except for the fact that the temperature rise values are correspondingly lower. Finally, the recorded device dynamic thermal responses are presented in Fig. 4. For the sake of picture clarity, the results are shown only for the heat source and two selected sensors located 10 mm and 40 mm away from the BJT. As can be seen, there is a visible delay of a few seconds in the thermal response at the sensor locations due to the diffusion of generated heat through the board. Moreover, the change of cooling conditions affects the heating curves only after approximately 10 s. 3. HEAT TRANSFER COEFFICIENT MODELLING 3.1. Estimation of local heat transfer coefficient values For the purposes of the local heat transfer coefficient estimation, a three-dimensional thermal model of the circuit was created. Taking into account that the only nonlinearity is related to temperature dependent boundary conditions, the temperature solutions were obtained with the analytical Green‟s function solver because of short computation time. This forward solver was coupled with different inverse estimation algorithms. The current local heat transfer coefficient values were updated in each iteration step based on the comparison between measured and simulated data. The procedure was stopped when the difference was smaller than the assumed simulation accuracy. 522 M. JANICKI, A. SAMSON, T. RASZKOWSKI, T. TORZEWICZ, A. NAPIERALSKI 0 10 20 30 40 50 60 70 80 90 100 110 120 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0 6.5 7.0 T e m p e ra tu re r is e [ K ] Power [W] BJT PT1 PT2 PT3 PT4 PT5 Fig. 2 Device steady state temperature rise for different dissipated power values with natural convection cooling. 0 10 20 30 40 50 60 70 80 90 100 110 120 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 T e m p e ra tu re r is e [ K ] Air velocity [m/s] BJT PT1 PT2 PT3 PT4 PT5 Fig. 3 Device steady state temperature rise for different air velocities at the 1000 mA emitter current with forced convection cooling. 0 10 20 30 40 50 60 70 80 90 100 110 120 1E-5 1E-4 1E-3 1E-2 1E-1 1E+0 1E+1 1E+2 1E+3 T e m p e ra tu re r is e ( K ) Time (s) BJT 1.0A 0.0m/s BJT 1.0A 2.0m/s PT2 1.0A 0.0m/s PT2 1.0A 2.0m/s PT4 1.0A 0.0m/s PT4 1.0A 2.0m/s Fig. 4 Measured heating curves of selected devices for different air velocities. Considerations on The Importance of Proper Heat Tranfer Coefficient Modeling in Air Cooled Electronic ... 523 The effectiveness of tested inverse estimation algorithms was assessed using as the performance indicators the number of iteration steps and the simulation time required to attain the desired accuracy. These algorithms included the false position, the Newton- Raphson, the secant and the bisection methods. However, the best results were obtained employing a simple numerical algorithm proposed by the authors. This algorithm has the shortest computation time, mainly because of the smallest number of algebraic operations required and partly owing to the efficient storage of previously computed heat transfer coefficient values. Obviously, these results cannot be generalized directly since each time the computation time depends on the initial heat transfer coefficient value h0 or the step size of the coefficient value update h0. For more detailed description of the algorithm tests, please refer to [15]. The main idea of the proposed algorithm consists in the iterative improvement of the heat transfer coefficient estimates h. When the difference between the measured value and the estimated one is positive, the temperature is underestimated and the coefficient value has to be decreased, and when the difference is negative, this value is increased. If the difference changes its sign, the step size h is reduced and the search becomes more accurate. It is also worth mentioning that in each iteration step the current heat transfer coefficient values and their corresponding temperature values in the transistor and sensor locations are stored, what allows significant reduction of the computation time. The total value of the heat transfer coefficient is supposed to model both radiation and convection cooling mechanisms occurring at outer structure surfaces and its value could be split into two components, as shown in Eq. 1, which can be estimated independently. The first of them encompasses the phenomena which are dependent on the surface and the surrounding ambient temperature values, i.e. mainly radiation and natural convection, and the other one reflects the effects related to the cooling air velocity with forced convection. fcnc hhh  (1) Therefore, originally the first component hnc was estimated assuming that its initial value equalled 5 W/(m 2 K), what corresponds approximately to the theoretical value for the pure radiation cooling at room temperature. The initial step size for the update of the heat transfer coefficient value h was 1 W/(m 2 K) and the desired temperature simulation accuracy (stopping criterion) was set to 0.05 K. The estimation results obtained using the algorithm proposed by the authors for the four BJT heating current values are shown with different markers in Fig. 5. As can be seen, the estimated dependence of the local heat transfer coefficient value on circuit surface temperature rise is not a linear function and the heat transfer coefficient ranges from 9.7 W/(m 2 K) 12 K over the ambient temperature to 15.4 W/(m 2 K) when the temperature rise amounts to 91 K. Next, the second component of the heat coefficient value, due to the forced movement of the cooling air, was estimated. First, the total values of the coefficient were determined, similarly as in the case of the radiation and the natural convection component, and then the previously stored coefficient values dependent only on the surface temperature were subtracted from the total values. The estimation results obtained for the heat source are plotted in Fig. 7 again with different markers corresponding to the respective BJT heating currents. As can be seen, between the air velocities of 1 m/s and 3 m/s the values of this component are at least doubled. The results obtained for sensor locations are very similar. 524 M. JANICKI, A. SAMSON, T. RASZKOWSKI, T. TORZEWICZ, A. NAPIERALSKI Fig. 5 Estimated dependence of the radiation and natural convection component on the surface temperature rise. Fig. 6 Estimated dependence of the forced convection component on the air velocity. 3.2. Determination of model parameter values Once the values of the heat transfer coefficient were estimated, it was possible to fit these values to a simple parametric model allowing the computation of local heat transfer coefficient values for given temperature rise and air velocity values. The first component hnc, related to the surface temperature rise T, can be modeled by the relation presented in Eq. 2. When a surface is at ambient temperature, the constant „a‟ reflects the radiation cooling. The theoretical value of the exponent „c‟ for a surface of uniform temperature and with pure natural convection cooling is 0.25 [4]-[5], but in reality it is much higher because it depends also on the radiation and surface temperature gradients. Considerations on The Importance of Proper Heat Tranfer Coefficient Modeling in Air Cooled Electronic ... 525 c Δba Thnc  (2) The heat transfer coefficient component related to the forced convection cooling could be modeled by a power function of the air velocity v shown in Eq. 3. The theoretical value of the exponent „e‟ is in the range of 2/3 ÷ 3/4 [4]-[5]. e d vh fc  (3) In order to estimate the unknown model parameters contained in Eq. 2-3, a modified version of the bee colony optimization algorithm, described in [15], was used. The results of the fitting obtained in the case of the natural and forced convection cooling are plotted in Figs. 7-8 respectively. The thick black line in these figures represents the results of the global fitting carried out for all measurement data available. The final expression for the total value of heat transfer coefficient is expressed by the following formula: 778.00.395 67.885.183.4 vTh  (4) The fitted dependence of the first heat transfer coefficient component hnc shows a very good agreement with the measurements. The value of model parameter „a‟ is close to the theoretical one and the value of the exponent „c‟, as expected, is relatively high due to the small size of the heat source occupying only slightly more than 1% of the circuit surface area. On the other hand, the value of the exponent „v‟ in the forced convection component hfc is just slightly higher than the one reported in literature. The global fit results in Fig. 8 seem not to be accurate, but the black line represents, as already mentioned, the results produced for the entire measurement data set. More accurate results can be obtained when the fitting is performed individually for the particular temperature sensor locations. Then, the heat transfer coefficient values predicted for the transistor location are much higher due to its elevated temperature. Fig. 7 Fitting results of the heat transfer coefficient component modeling the radiation and natural convection cooling. 526 M. JANICKI, A. SAMSON, T. RASZKOWSKI, T. TORZEWICZ, A. NAPIERALSKI Fig. 8 Fitting results of the heat transfer coefficient component modeling the forced convection cooling. The final parametric model for the heat transfer coefficient value, expressed by Eq. 4, was implemented in the forward thermal solver. The simulations were repeated using the variable coefficient values and the constant ones equal to the lowest and the highest value used in the preceding simulation with the variable temperature and air velocity dependent value. The results of these simulation are compared in Fig. 9 with the measurement. The assumption of a low heat transfer coefficient value leads to the important overestimation of the heat source temperature (dotted line), by almost 30%. On the other hand, the high value of the coefficient assures the accurate result in the thermal steady state (solid line), but in the region when the heat diffuses through the board, between 1 s and a few minutes, the simulation errors exceed 10 K. 0 10 20 30 40 50 60 70 80 90 100 110 120 130 140 150 1E-5 1E-4 1E-3 1E-2 1E-1 1E+0 1E+1 1E+2 1E+3 T e m p e ra tu re r is e ( K ) Time (s) BJT 1.0A 0m/s MES BJT 1.0A 0m/s HIGH BJT 1.0A 0m/s VAR BJT 1.0A 0m/s LOW Fig. 9 Comparison of measured and simulated transistor heating curves. Considerations on The Importance of Proper Heat Tranfer Coefficient Modeling in Air Cooled Electronic ... 527 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 1.00 1.05 0 10 20 30 40 50 60 70 80 90 100 110 120 130 140 150 160 N o rm a li s e d te m p e ra tu re r is e (- ) Distance (mm) natural convection forced convection Fig. 10 Comparison of temperature profiles along the cross-section of the circuit board in different cooling conditions. Finally, it was instructive to compare the temperature profiles across the circuit board obtained in the case of the natural and forced convection cooling. These results, presented in Fig. 10, represent the temperature rise relative to the heat source temperature. As can be seen with the forced convection cooling the temperature profile is much steeper and consequently the differences in the heat transfer coefficient values are more important, but the thermal influence coefficients between the heat source and sensor locations decrease reducing the so-called „cooling circle. Similar conclusions were drawn also in [16]. 4. CONCLUSIONS This paper discussed the problem of proper modelling of the heat transfer coefficient value in electronic circuits. Currently, it is a common practice that for thermal simulations the same coefficient value, usually around 10 W/(m 2 K) for natural convection, is assumed for the entire surface of an electronic circuit. However, as it was demonstrated here based on a practical example, the use of constant coefficient values might lead to significant simulation errors, both in transient and thermal steady states. This is caused mainly by the fact that the heat transfer coefficient values can vary across the surface of a circuit even by an order of magnitude. The simple parametric model proposed by the authors allowed the computation of the local heat transfer coefficient values in function of surface temperature rise and cooling air velocity. This model, fitted to experimental data, can be used for the iterative updates of the local heat transfer coefficient values during thermal simulations. When included into standard FEM thermal analysis tools or lumped thermal models, such as the ones described in [17]-[18], it could improve significantly the thermal simulation accuracy. Additionally, it was shown that the parameter values in the proposed model differ substantially from the standard heat transfer textbook formulas, which are derived mostly for cases when surface temperature or heat flux is uniform. 528 M. JANICKI, A. SAMSON, T. RASZKOWSKI, T. TORZEWICZ, A. NAPIERALSKI REFERENCES [1] R. Ross, ed., Microelectronics Failure Analysis: Desk Reference. EDFAS, ASM International, 2011. [2] M. Janicki, Z. Sarkany and A. Napieralski, “Impact of Nonlinearities on Electronic Device Transient Thermal Responses”, Microelectron. J., vol. 45, pp. 1721-1725, December 2014. [3] F.P. Incropera and D.P. De Witt, Fundamentals of Heat and Mass Transfer. Wiley & Sons, 2002. [4] J.P. Holman, Heat Transfer. McGraw-Hill, 1986. [5] J.H. Lienhard IV and J.H. Lienhard V, A Heat Transfer Textbook. Phlogiston Press, 2012. [6] G. Ellison, Thermal Computations for Electronics: Conductive, Radiative, Convective Air Cooling. CRC Press, 2010. [7] Y. M. Li and A. Ortega, “Forced Convection from A Rectangular Heat Source in Uniform Shear Flow: The Conjugate Peclet Number in The Thin Plate Limit”, In Proceedings of the 6 th ITHERM. Seattle, WA: IEEE, 1998, pp. 284-294. [8] R. Aster, B. Borchers, and C. Thurber, Parameter Estimation and Inverse Problems. Elsevier, 2005. [9] M.N. Ozisik and H.R.B. Orlande, Inverse Heat Transfer. Taylor & Francis, 2000. [10] M. Janicki and S. Kindermann, “Recovering Temperature Dependence of Heat Transfer Coefficient in Electronic Circuits”, Inverse Probl. Sci. En., vol. 17, pp. 1129-1142, October 2009. [11] D. Karaboga, and B. Akay, “Survey: Algorithms Simulating Bee Swarm Intelligence”, Artif. Intell. Rev., vol. 31, pp. 68-85, June 2009. [12] D. Teodorovic, P. Lucic, G. Z. Markovic and M. Dell'Orco, “Bee Colony Optimization: Principles and Applications”, In Proceedings of the 8 th NEUREL. Belgrade: IEEE, 2006, pp. 151-156. [13] A. Colorni, M. Dorigo, and V. Maniezzo, “Distributed Optimization by Ant Colonies”, In Proceedings of the ECAL. Paris: Elsevier, 1991, pp. 134-142. [14] T. Torzewicz, A. Samson, T. Raszkowski, A. Sobczak, M. Janicki, M. Zubert, A. Napieralski, “Thermal Analysis of Hybrid Circuits with Variable Heat Transfer Coefficient”. In Proceedings of the 33 rd SEMI- THERM, San Jose, CA, 2017, pp. 19-22. [15] T. Raszkowski and A. Samson, “Application of Genetic Algorithm and Swarm Intelligence Algorithms to Heat Transfer Coefficient Estimation”, Bulletin de la Société des Sciences et des Lettres de Łódź. Série: Recherches sur les Déformations, vol. LXVII, pp. 103-125, October 2017. [16] G. A. (W.) Luiten, “Characteristic Length and Cooling Circle”, In Proceedings of the 26 th SEMI- THERM. Santa Clara, CA: IEEE, 2010, pp. 7-13. [17] R. Menozzi, P. Cova, N. Delmonte, F. Giuliani, and G. Sozzi, “Thermal and Electro-thermal Modeling of Components and Systems: Review of The Research at The University of Parma”, Facta Universitatis, Series: Electronics and Energetics, vol. 28, pp. 325-344, September 2015. [18] S. Stanisic, M. Jevtic, B. Das, and Z. Radakovic, “FEM CFD Analysis of Air Flow in Kiosk Substation with Oil Immersed Distribution Transformer”, Facta Universitatis, Series: Electronics and Energetics, vol. 31, pp. 411-423, September 2018.