International Communications in Heat and Mass Transfer 48 (2013) 108–115
Contents lists available at ScienceDirect
International Communications in Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ichmt
A 3D numerical study of heat transfer in a single-phase micro-channel heat sink using graphene, aluminum and silicon as substrates☆ Ahmed Jassim Shkarah a,b,⁎, Mohd Yusoff Bin Sulaiman a, Md Razali Bin Hj Ayob a, Hussein Togun b,c a b c
Faculty of Mechanical Engineering, Universiti Teknikal Malaysia Melaka (UTeM), Melaka, Malaysia Department of Mechanical Engineering, Thi-Qar University, 64001 Nassiriya, Iraq Department of Mechanical Engineering, University of Malaya, 50603 Kuala Lumpur, Malaysia
a r t i c l e
i n f o
Available online 7 September 2013 Keywords: Micro-channel Single phase Graphene Thermal resistance
a b s t r a c t The study describes numerical simulations conducted on micro-channel heat sinks. Three different shapes related to the micro-channel depth and width is chosen for examination. Silicon, aluminum, and graphene are used as substrate materials for this study. The overall heat sink consisted of an array of rectangular micro-channels. Three different surface heat fluxes and three different volumetric flow rates are used for three cases. Water with nontemperature-dependent thermal properties is used as a coolant for steady-state, fully developed laminar flow in the micro-channels. From a heat transfer (thermal performance) perspective, it is found that graphene most effectively reduce the thermal resistance. Based on these results, graphene was further studied as a substrate material for a micro-channel heat sink. © 2013 Elsevier Ltd. All rights reserved.
1. Introduction Over the last decade, micromachining technology has been increasingly used to develop highly efficient heat sink cooling devices due to advantages such as lower coolant demands and smaller machinable dimensions. One of the most important micromachining technologies is the ability to fabricate micro-channels. Hence, the studies of fluid flow and heat transfer in micro-channels, which are two essential parts of such devices, have attracted attention due to their broad potential for solving both engineering and medical problems. Heat sinks are classified as either single-phase or two-phase according to whether liquid boiling occurs inside of the micro-channels [1]. The primary parameters that determine the single-phase and two-phase operating regimes are the heat flux through the channel wall and the coolant flow rate. For a fixed heat flux (heat load), the coolant may maintain its liquid state throughout the micro-channels. For a lower flow rate, the liquid coolant flowing inside the micro-channel may reach its boiling point, causing flow boiling to occurs, resulting in a two-phase heat sink.
1.1. Micro-channel and its use Tuckerman and Pease [2] first made use of miniaturization for the purposes of heat removal within the scope of a Ph.D. study in 1981. Their publication, entitled “High Performance Heat Sinking for VLSI,” was the ☆ Communicated by W.J. Minkowycz. ⁎ Corresponding author at: Faculty of Mechanical Engineering, Universiti Teknikal Malaysia Melaka (UTeM), Melaka, Malaysia. E-mail address:
[email protected] (A.J. Shkarah). 0735-1933/$ – see front matter © 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.icheatmasstransfer.2013.08.006
first study on micro-channel heat transfer. This pioneering work motivated many researchers to focus on micro-channel flow and led to its recognition as a high-performance heat removal tool. Before further discussion micro-channel flow and heat transfer, it is appropriate to introduce a definition for the term “micro-channel.” The scope of the term is a topic of debate among researchers in the field. Mehendale et al. [3] used the following classification based on the manufacturing techniques required to obtain various ranges of the smallest channel dimension, D: 1 μm b D b 100 μm micro-channels 100 μm b D b 1 mm minichannels 1 mm b D b 6 mm compact passages 6 mm b D conventional passages Kandlikar and Grande [4] adopted a different classification based on the rarefaction effect of gases in various ranges of the smallest channel dimensions, D: 1 μm b D b 10 μm transitional micro-channels 10 μm b D b 200 μm micro-channels 200 μm b D b 3 mm mini-channels 3 mm b D conventional passages A simpler classification was proposed by Obot [5] and was based on the channel hydraulic diameter rather than the smallest channel dimension. Obot classified channels with a hydraulic diameter less than 1 mm (Dh b 1 mm) as micro-channels. This standard was also adopted by other researchers, such as Bahrami and Jovanovich [6], and is considered to be more appropriate for the purposes of this work. The higher
A.J. Shkarah et al. / International Communications in Heat and Mass Transfer 48 (2013) 108–115
Nomenclature A Ac Cp D Dh H h K L P q” Q Rth T t u v wc ww w x y z ΔP P ρ μ
Total surface area of micro-channel (m2) Area of cross‐section of micro-channel (m2) Specific heat capacity (J/kg K) Diameter (m) Hydraulic diameter (m) Height of channel (m) Convective heat transfer coefficient Thermal conductivity (W/m K) Heat sink length (m) Total pressure (Pa) Heat flux (W/cm2) Volumetric flow rate (cm3/s) Thermal resistance (K/W) Temperature (K) Wall thickness at bottom (m) Fluid x-component velocity (m/s) Fluid y-component velocity (m/s) Channel width (m) Wall thickness at bottom (m) Fluid z-component velocity (m/s) Axial coordinate Vertical coordinate Horizontal coordinate Pressure drop (Pa) Pressure Density (kg/m3) Dynamic viscosity (m2/s)
Subscripts: o output in input s substrate w water m mean value max maximum
volumetric heat transfer densities required for advanced manufacturing techniques have led to more complex manifold designs [7]. Many of the same manufacturing techniques developed for the fabrication of electronic circuits are being used in the fabrication of compact heat exchangers. Micro-channel heat sinks constitute an innovative cooling technology for the removal of large amounts of heat through a small area. This technology is one of the potential alternatives for replacing conventional finned tube heat exchangers, which are primarily used in industries where automobiles, air conditioning and refrigeration are prevalent. A heat sink is usually made from a high thermal conductivity material, such as silicon or copper [8].The micro-channels are fabricated into its surface by either precision machining or micro-fabrication technology. A micro-channel heat sinks typically contains a large number of parallel micro-channels. Coolant is forced to pass through these channels to carry heat from a hot surface. In micro-channel heat exchangers, the flow is typically laminar and the heat transfer coefficients are proportional to the velocity. Micro-channel heat sinks provide very high surface area to volume ratios and have large convective heat transfer coefficients, a small mass and volume, and low coolant requirements. These attributes allow these heat sinks to be suitable for cooling devices such as high-performance microprocessors, laser diode arrays, radars, and high-energy-laser mirrors. Additionally micro-channel heat exchangers may be easier to repair than their conventional counterparts. Micro-channel heat exchangers offer other benefits, including an increased latent capacity for micro-channel evaporators. Micro-channel
109
heat exchangers improve heat transfer in two ways. First, the smaller dimensions in the refrigerant flow passages increase the refrigerant-side heat transfer rate. Second, the flat tube orientation reduces the airside flow resistance, leading to either an increase in airflow or the option to reduce the fan power. A micro-reactor is a device, typically with lateral dimensions of less than 1.0 mm, in which confined chemical reactions take place. If the reactor needs to be maintained at a particular temperature, a micro-reactor may also act as a micro-channel heat exchanger. Micro-reactors used for these applications are typically continuous flow reactors rather than batch reactors. Continuous flow models provide better performance in material synthesis than is possible with batch reactors. These reactors offer many advantages over conventionally sized reactors, including vast improvements in energy efficiency, reaction speed and yield, safety, reliability, scalability, on-site/on-demand production, and a much finer degree of process control. These improvements enable miniaturization of the fuel processor to minimize heat and mass transfer resistance. A micro-reactor may also be used to perform small scale reactions to determine if the potential for dangerous situations exists such as runaway reactions or the generation of excessive levels of heat. 1.2. Thermal properties of graphene High-performance and low-cost composites are highly desirable materials for mechanical, civil, and aerospace applications. Carbon fibers, first created in 1950s, are still major components in high-performance composites due to their remarkable mechanical properties, relatively easy and cheap fabrication process and low density. Recently, due to the development of nano-scale synthesis, engineering technologies, and incorporation of the hierarchical structures in biological materials, nano-composites featuring superior stiffness, strength, and energy dissipation capacities are being touted as the next-generation of multifunctional super materials [9–11]. For a long time, the performance of integrated circuits has been improved by geometrically scaling the silicon. In recent years, several severe scaling limitations have led to the introduction of novel materials via CMOS technologies. Recently, gate oxide scaling limits were overcome through the introduction of novel gate stack materials [REFERENZ]. Now, the MOSFET channel seems to be one of the limiting factors to continued scaling. One major impediment to further shrinking the channel length is the occurrence of short channel effects that result in higher off-state leakage currents [12]. The authors observed a thermal conductivity of K = 3410 W/m*K [9], which is clearly higher than the bulk graphite limit of 2000 W/m*K and is in agreement with the experiments performed by Ghosh et al. [13–15]. 2. Mathematical formulation A 3D conjugate heat transfer approach was used to simulate the micro-channel flow, incorporating both the conduction from the substrate and the convection from the channel fluid. The physical model and coordinate system is shown in Fig. 1. The simulations were performed using the FLUENT 14 [16] code. The micro-channel was simulated for several different dimensions, as shown in Table 1. The following assumptions were made when modeling the heat transfer characteristics in the rectangular cross-section of the microchannels: (i) The fluid properties are constant. (ii) The fluid is laminar, steady, and incompressible. (iii) The heat loss to the environment is negligible. In this study, the length L and the width LZ of the micro-heat sink were chosen to be 10 mm × 10 m (1 cm2) (see Fig. 1). A full 3-dimensional (3D) conjugate heat transfer model of the micro-channel heat sinks was created. As shown in Fig. 1, this numerical
110
A.J. Shkarah et al. / International Communications in Heat and Mass Transfer 48 (2013) 108–115
Fig. 1. A schematic model of micro-channel heat sink.
model utilizes a cell with one rectangular channel as the computational domain, due to the symmetrical structure of the heat sink. If the gravitational forces and heat dissipation from fluid viscosity are neglected for steady, incompressible, laminar flow, the governing equations for the liquid flow are as follows: The continuity equation:
Rth ¼
∂u ∂v ∂w þ þ ¼0 ∂x ∂y ∂z
ð1Þ
z ∂u ∂u ∂u ∂P ∂ u ∂z u ∂z u þ þ ¼− þμ ρ u þv þw ∂x ∂y ∂z ∂x ∂xz ∂yz ∂zz
ð2Þ
ð3Þ
z z z ∂w ∂w ∂w ∂P ∂w ∂w ∂w þ þ ρ u þv þw ¼− þμ ∂x ∂y ∂z ∂z ∂xz ∂yz ∂zz
ð4Þ
The hydrodynamic boundary conditions: at the channel wall surface ðno‐slipÞ; at the inlet;
z ¼ L; P ¼ P out at the outlet;
ð5Þ ð6Þ ð7Þ
the energy equation for the liquid in the channel: z ∂T ∂T ∂T ∂ T 1 ∂z T 1 ∂z T 1 þ þ ρl C P u 1 þ v 1 þ w 1 ¼ λ l ∂x ∂y ∂z ∂xz ∂yz ∂zz
ð8Þ
and the heat conduction in the solid: ∂ ∂T ∂ ∂T ∂ ∂T λs s þ λs s þ λs s ¼ 0 ∂x ∂y ∂z ∂x ∂y ∂z
ð9Þ
For the thermal boundary conditions, adiabatic boundary conditions are applied to all the region boundaries except for heat sink bottom wall. At this boundary, a constant heat flux is assumed: ks
ð11Þ
∂T s ¼ q″ ; for 0 ≤ z ≤ L; 0 ≤ x ≤ W; and y ¼ H: ∂z
ΔT
max
¼ T s;0 −T f;i
ð12Þ
The pumping power, P, defined as
z ∂v ∂v ∂v ∂P ∂ v ∂z v ∂z v ρ u þv þw þ þ ¼− þμ ∂x ∂y ∂z ∂y ∂xz ∂yz ∂zz
z ¼ 0; P ¼ P m ; u ¼ constant ; T ¼ T m
ΔT max q}As
where As is the area of the substrate subjected to the heat flux, and ΔTmax is the maximum temperature rise in the heat sink defined as:
The momentum equations:
u¼v¼w¼0
It should be noted that in actual applications, the heat loss from the solid substrate to the environment may be substantial and would need to be considered. The overall thermal resistance can be defined by the following equation:
ð10Þ
P ¼V˙ ΔP ¼ N um Ac ΔP
ð13Þ
where V˙ is the total volume flow rate, N is the number of channels, um is the mean channel velocity, Ac is the cross‐sectional area of the channel, and ΔP is the pressure drop across the channel. This expression and approach has been widely used by previous researchers [17–19]. Uniform heat fluxes, representing heat sources, were applied to the base of the main domain in the simulations. The side walls and the top wall were set to be adiabatic. The interior walls of the substrate (the channel surface) were set to be “coupled walls.” A uniform velocity profile was specified at the micro-channel inlet for defining and simultaneously developing the flow conditions. A “pressure outlet” boundary condition was applied at the outlet. The substrate material used in the simulations was a combination of aluminum, silicon, and grapheme, with water as the coolant liquid. A finite volume method (FVM) was used to convert the governing equations to algebraic relations via a first order upwind scheme. The SIMPLE algorithm was used in FLUENT 14 to enforce mass conservation and obtain the pressure field. The segregated solver was used to solve for the conservation of mass, momentum, and energy governing integral equations. From this analysis, the distribution of velocity, pressure, and temperature in the micro-channel heat sink was determined. 3. Results and discussion To characterize any numerical and modeling errors, three different meshes were generated, with the ANSYS ICEM-CFD mesh generator to solve for the property distributions. A grid-independent study was performed using the three different hexahedral element
A.J. Shkarah et al. / International Communications in Heat and Mass Transfer 48 (2013) 108–115
111
Table 1 Dimensions of all the simulated micro-channel cases. Case
Wc(μm)
Ww(μm)
H(μm)
1 2 3
56 55 50
44 45 50
320 287 302
Table 2 Results from the grid-independent study. Mesh size
Temperature at Z = 0.01 mm (K)
10,20,400 30,40,200 60,80,400
343.176 343.669 343.794
Table 3 Numerical validation of the thermal resistances using both experimental data [20] and numerical data [21]. Cases
q" (W/cm2)
Q(cm3/s)
Rth(oC/W), experimental data [20]
Rth(oC/W), numerical [21]
Rth(oC/W), present work
1 2 3
181 277 790
4.7 6.5 8.6
0.110 0.113 0.090
0.157 0.128 0.105
0.145 0.129 0.109
Fig. 3. Temperature distribution for case 3 in the bottom of the heat sink (using graphene) for different heat fluxes at Q = 6.5 cm3/s.
mesh sizes. The results from the different meshes are listed in Table 2. Additionally, Table 2 shows that the solution is independent of the mesh size, and a further increase in the mesh size will not have a significant effect on the accuracy of the solution. Therefore, obtain good accuracy, (10, 20, 400), (30, 40, 200), and (60, 80, 400) sized meshes are used. The convergence criterion for the momentum and energy equations was set to be less than 10− 6. The PC used for these analysis possessed Intel ® Core™ i7-3770, a 3.4-GHz processor, and a 32-GB RAM.
To ensure the accuracy of the numerical model, results were compared with experimental results from reference [20] and other numerical results from reference [21]. Temperature-dependent fluid properties, as shown in Table 3, and silicon were used. The model used in this work shows good agreement with the experimental results, even for higher heat fluxes. Table 3 shows the agreement between the model results and the results from references [20,21]. From these results, it can be concluded that the present model can be used to simulate a micro-channel heat sink with water as a cooling fluid.
Fig. 2. Temperature distribution for case 1 in the bottom of the heat sink (using graphene) for different fluid flow rates at q" = 181 W/cm2.
Fig. 4. Temperature distribution for case 3 in the bottom of the heat sink (using graphene) for different heat fluxes at Q = 8.6 cm3/s.
112
A.J. Shkarah et al. / International Communications in Heat and Mass Transfer 48 (2013) 108–115
Fig. 5. Temperature distribution for case 2 in the fluid (using graphene) for different heat fluxes at Q = 8.6 cm3/s.
Fig. 7. Thermal resistance for case 2 at different volume flow rates and q" = 277 W/cm2.
Fig. 2 shows the temperature distribution at the bottom of the heat sink using graphene as the substrate material. The temperature decreased as the volume flow rate increased (for same identical heat flux values) due to the increase in flow velocity. The temperature distributions at different flow rates have similar trends, as shown in Fig. 2. Figs. 3 and 4 show the temperature distributions for the heat fluxes of q" = (181, 277, and 790), W/cm2 and volume flow rate of Q = (6.5 and 8.6) cm3/s. These figures illustrate that the temperature increased from the entrance of the micro-channel to the exit
of micro-channel similarly for the different cases. The difference between the input and the output temperatures increased with the amount of surface heat flux applied at the bottom heat sink was changed. Fig. 5 shows the fluid temperature distribution (using graphene) for different heat fluxes at a flow rate Q = 8.6 cm3/s and for simulation case 2. Fig. 5 also illustrates that the fluid temperature increased from the inlet of the micro-channel to the output of the micro-channel. The input temperature is same for the all applied heat fluxes. The trend is same for all the heat fluxes.
Fig. 6. Thermal resistance for case 1 at different volume flow rates and q" = 790 W/cm2.
Fig. 8. Thermal resistance for case 3 at different volume flow rates and q" = 181 W/cm2.
A.J. Shkarah et al. / International Communications in Heat and Mass Transfer 48 (2013) 108–115
Fig. 9. Thermal resistance for Silicon at different volume flow rates and q" = 181 W/cm2.
Figs. 6–8 show that the thermal resistance
113
Fig. 11. Thermal resistance for aluminum at different volume flow rates and q" = 181 W/cm2.
max deRth ¼ ΔTq}A s
creased as the volume flow rate increased when using identical substrate materials. This is expected because when the volume flow rate increased the heat transmitted by fluid also increased so that the maximum temperature at the solid (Ts,0). decreased. Additionally, the thermal resistance for a graphene substrate material was the lowest examined. Figs. 6 and 7 show the same decreasing thermal resistance trend for silicon, aluminum, and grapheme. Additionally, the deference in the thermal resistance between silicon, aluminum, and graphene did not change. From case 1 to case 2 and from q" = 277 W/cm2 to
q" = 790 W/cm2, no difference were observed in the trends, except for the low heat flux case (q" = 181 W/cm2). Fig. 8 shows that the difference in the thermal resistance at the input of the micro-channel was greater than at the output. Figs. 6–8 show that the graphene is the best material to use as a substrate. Figs. 9–13 show the effects of the micro-channel shape on the thermal resistance for a controlled variation in the volume flow rate. Fig. 9 shows that the thermal resistance in case 2 is greater than in case 1 and 3.
Fig. 10. Thermal resistance for graphene at different volume flow rates and q" = 181 W/cm2.
Fig. 12. Thermal resistance for silicon at different volume flow rates and q" = 277 W/cm2.
114
A.J. Shkarah et al. / International Communications in Heat and Mass Transfer 48 (2013) 108–115
Fig. 15. Temperature contours for case 1 using aluminum at q" = (a: 181 b: 277 c: 790) W/cm2 and Q = 4.7 cm3/s. Fig. 13. Thermal resistance for graphene at different volume flow rates and q" = 790 W/cm2.
Fig. 9 also shows that case 2 has the greatest thermal resistance. Additionally, Fig. 9 shows that there is not a significant difference in the thermal resistances between cases 1, 2, and 4 for silicon with an applied heat flux of q" = 181 W/cm2. Figs. 10 and 13 show the thermal resistance for all three cases with a graphene substrate; the thermal resistance for case 2 is the greatest one. In Fig. 10, the value of the thermal resistance is almost equal to that of an aluminum substrate with volume flow rates greater than 6.5 cm3/s and a heat flux of q" = 181 W/cm2. Additionally, in Fig. 12, the thermal resistance is nearly equal for small volume flow rates, but the difference in thermal resistance increased as the volume flow rate increased. Figs. 14–16 show the temperature contours for silicon, aluminum, and graphene substrate for case 1 with different volumetric
Fig. 16. Temperature contours for case 1 using graphene at q" = (a: 181 b: 277 c: 500 d: 790) W/cm2 and Q = 4.7 cm3/s.
flow rates and surface heat fluxes. These figures show that the maximum temperature is at the end of the micro-channels heat sink bottom. References
Fig. 14. Temperature contours for case 1 using (a) silicon, (b) aluminum, and (c) graphene at q" = 790 W/cm2 and Q = 4.7 cm3/s.
[1] A. Jassim Shkarah, M.Y. Bin Sulaiman, R.H. Ayob, Numerical analysis of two phase flow heat transfer in rectangular micro channel heat sink, Arch. Sci. 65 (2012) 602–613. [2] D.B. Tuckerman, R.F.W. Pease, High-performance heat sinking for VLSI, IEEE Electron Device Lett. 2 (1981) 126–129. [3] S.S. Mehendale, A.M. Jacobi, R.K. Shah, Fluid flow and heat transfer at micro- and meso-scales with application to heat exchanger design, Appl. Mech. Rev. 53 (2000) 175–193. [4] S.G. Kandlikar, W.J. Grande, Evolution of microchannel flow passages-thermohydraulic performance and fabrication technology, Heat Transfer Eng. 24 (2003) 3–17. [5] N.T. Obot, Toward a better understanding of friction and heat/mass transfer in microchannels—a literature review, Microscale Thermophys. Eng. 6 (2002) 155–173.
A.J. Shkarah et al. / International Communications in Heat and Mass Transfer 48 (2013) 108–115 [6] M. Bahrami, M.M. Yovanovich, J.R. Culham, Pressure drop of fully developed, laminar flow in rough microtubes, J. Fluids Eng. 128 (2006) 632–637. [7] S.G. Kandlikar, Chapter 5—flow boiling in minichannels and microchannels, Heat Transfer and Fluid Flow in Minichannels and Microchannels, Elsevier Science Ltd, Oxford, 2006. 175–226. [8] A. Jassim Shkarah, M.Y. Bin Sulaiman, R.H. Ayob, Two-phase flow in micro-channel heat sink review paper, Int. Rev. Mech. Eng. 7 (2013) 231–237. [9] Y. Liu, B. Xie, Z. Zhang, Q. Zheng, Z. Xu, Mechanical properties of graphene papers, J. Mech. Phys. Solids 60 (2012) 591–605. [10] F. Barthelat, H. Tang, P.D. Zavattieri, C.M. Li, H.D. Espinosa, On the mechanics of mother-of-pearl: a key feature in the material hierarchical structure, J. Mech. Phys. Solids 55 (2007) 306–337. [11] P. Fratzl, S. Weiner, Bio-inspired materials—mining the old literature for new ideas, Adv. Mater. 22 (2010) 4547–4550. [12] M. Baus, T.J. Echtermeyer, B.N. Szafranek, M.C. Lemme, H. Kurz, Device architectures based on graphene channels, Integrated Circuit Design and Technology and Tutorial, ICICDT 2008 IEEE International Conference on2008, 2008, pp. 269–272. [13] S. Ghosh, S. Subrina, V. Goyal, D.L. Nika, E.P. Pokatilov, A.A. Balandin, Extraordinary thermal conductivity of graphene: possibility of thermal management applications,
[14]
[15]
[16] [17] [18] [19] [20] [21]
115
Thermal and Thermomechanical Phenomena in Electronic Systems (ITherm), 12th IEEE Intersociety Conference on 2010, 2010, pp. 1–5. A.A. Balandin, S. Ghosh, W. Bao, I. Calizo, D. Teweldebrhan, F. Miao, et al., Superior thermal conductivity of single-layer graphene, Nano Lett. 8 (2008) 902–907. S. Ghosh, I. Calizo, D. Teweldebrhan, E.P. Pokatilov, D.L. Nika, A.A. Balandin, et al., Extremely high thermal conductivity of graphene: prospects for thermal management applications in nanoelectronic circuits, Appl. Phys. Lett. 92 (2008) 151911–151913. C. ANSYS Inc., ANSYS FLUENT Theory Guide, Release 14.0, November 2011. A. Weisberg, H.H. Bau, J.N. Zemel, Analysis of microchannels for integrated cooling, Int. J. Heat Mass Transf. 35 (1992) 2465–2474. J.H. Ryu, D.H. Choi, S.J. Kim, Numerical optimization of the thermal performance of a microchannel heat sink, Int. J. Heat Mass Transf. 45 (2002) 2823–2827. L. Ji, G.P. Peterson, Geometric optimization of a micro heat sink with liquid flow, IEEE Trans. Compon. Packag. Technol. 29 (2006) 145–154. D.B. Tuckerman, Heat transfer microstructures for integrated circuits, Stanford University, 1984. K.C. Toh, X.Y. Chen, J.C. Chai, Numerical computation of fluid flow and heat transfer in microchannels, Int. J. Heat Mass Transf. 45 (2002) 5133–5141.