Abstract

The industry shift to multicore microprocessor architecture will likely cause higher temperature nonuniformity on chip surfaces, exacerbating the problem of chip reliability and lifespan. While advanced cooling technologies like two phase embedded cooling exist, the technological risks of such solutions make conventional cooling technologies more desirable. One such solution is remote cooling with heatsinks with sequential conduction resistance from chip to module. The objective of this work is to numerically demonstrate a novel concept to remotely cool chips with hotspots and maximize chip temperature uniformity using an optimized flow distribution under constrained geometric parameters for the heatsink. The optimally distributed flow conditions presented here are intended to maximize the heat transfer from a nonuniform chip power map by actively directing flow to a hotspot region. The hotspot-targeted parallel microchannel liquid cooling design is evaluated against a baseline uniform flow conventional liquid cooling design for the industry pressure drop limit of approximately 20 kPa. For an average steady-state heat flux of 145 W/cm2 on core areas (hotspots) and 18 W/cm2 on the remaining chip area (background), the chip temperature uniformity is improved by 10%. Moreover, the heatsink design has improved chip temperature uniformity without a need for any additional system level complexity, which also reduces reliability risks.

1 Introduction

On-chip hotspots, which are areas of elevated temperature that occur on the chip level, have intense adverse effects on electronic device performance and reliability. New techniques and remedies have been innovated by thermal researchers to alleviate these effects [1]. Much research has focused on microfabricated thermoelectric coolers (TEC), orthotropic thermal interface materials (TIMs), and anisotropic heat spreaders for hotspot remediation [26].

The hotspot remedy considered in this study is a modification of traditional microchannel liquid-cooled heatsinks, which adapts them to cooling nonuniform heat fluxes. This modification to heatsink architecture (channel network and/or manifolds) results in conducting more flow toward the hotspots and/or increasing convective heat transfer surface area on the top of the hotspots [711]. The former (conducting more flow toward the hotspots) is the core idea of this study. To the best of our knowledge, the work presented in this paper is a precursive study on the design and optimization of liquid-cooled heatsinks that prioritizes hotspots over background for the flow distribution. A few designs were found in literature that target areas of high heat flux by coolant flow rate [1214]. A microgap heatsink that operated with two separate unmixed fluids was studied by Green et al. [12]. In their design, one fluid is acted as a fluidic spreader dedicated to cooling the hotspot, while the second fluid served as both a coolant for the background heat flux and an on-chip regenerator for the hotspot fluid. Sharma et al. [13] designed a manifold microchannel liquid-cooled heatsink to target on-chip hotspots by conducting the inlet flow to the top of the hotspots. In their design, coolant impinged on the top of the hotspots after leaving the inlet manifold and was collected from the background to the outlet manifold. Kumar and Singh [14] managed chip temperature nonuniformity by taking advantage of flow maldistribution created by different types of heatsink manifolds. Using numerical modeling, they showed that different manifolds (I-Type, Z-Type, C-Type, and N-Type) could reduce the thermal resistance of heatsinks cooling chips with center, right, left, and dual streamlines heating.

A scientific optimization methodology based on design of experiments (DOE) was implemented to optimize this study's design. DOE is a statistical technique that examines the relationship between design and response parameters [15]. When combined with regression analysis and optimization algorithms, DOE is used extensively for optimizing the thermohydraulic performance of cooling devices such as gas and liquid-cooled heatsinks. Rao et al. [16] used DOE response surface method (RSM) and computational fluid dynamics (CFD) to calculate the regression functions for thermal resistance and pump power of a microchannel heatsink. Then, using an iterative JAYA algorithm, they calculated the optimal design parameters. In order to minimize the thermal resistance and pump power of a double-layered microchannel heatsink, Kulkarni et al. [17] performed a multi-objective optimization using DOE RSM and numerical simulation. Naqiuddin et al. [18] employed DOE Taguchi–Grey method to maximize the heat transfer performance and minimize the pressure drop and temperature gradient in a microchannel heatsink. In another study, RSM was utilized for the shape optimization of a microchannel heatsink employed to cool a multidie module (discrete nonuniform heat flux boundary condition), which maximized the heatsink thermal performance and minimized its pressure drop [19]. Using DOE, microchannel heatsinks with different wall surfaces including uncoated, hydrophobic, and super-hydrophobic, were optimized geometrically to minimize thermal resistance and pump power [20]. Nakhchi [21] optimized the heat transfer and pressure drop inside a liquid-cooled sinusoidal wavy channel. The author used DOE response surface method to calculate the optimal geometries that possessed minimum pressure drop and maximum Nusselt number. Chen et al. [22] used DOE Taguchi method to perform a screening experimental design that determined the crucial parameters in the design of a parallel-plain fin heatsink. Using DOE response surface method, they optimized the thermohydraulic performance of the heatsink.

Modern remedies for hotspots, such as TECs, anisotropic graphite heat spreaders, and orthotropic TIMs benefit high-tech designs in material science. However, there are serious limitations and drawbacks associated with them. TECs suffer from high cost and low thermal efficiency. Graphite heat spreaders offer low mechanical strength and flexibility. Orthotropic TIMs have operational temperature limits that make them inappropriate for high heat flux applications. The disadvantages of modern hotspot solutions were the motivation behind this study.

The objective of this work is to numerically demonstrate a novel concept to remotely cool chips with hotspots and to maximize chip temperature uniformity using an optimized flow distribution under constrained geometric parameters for the heatsink and constant total coolant flowrate. A system of manifolds was designed to concentrate the coolant flowrate to the channels over a hotspot with a heat flux eight times that of background. The hotspot, located on the center of the chip, yielded a heat flux of 145.5 W/cm2 compared to the background heat flux of 18.1 W/cm2. The effect of flow concentration toward the hotspot on the thermal and heat spreading resistivity of the heatsink, chip temperature nonuniformity, and pump power was investigated. Next, using DOE with a full factorial design (FFD), the distribution was numerically optimized to achieve a minimum chip temperature nonuniformity, heatsink thermal resistivity, and pump power. The obtained optimal designs were compared with a baseline uniform flow conventional liquid cooling design.

The studies available in the literatures on hotspot mitigation deal with various techniques to direct the coldest fluid to the hottest region at high flow rates or the cooling solution is embedded on to the die targeting the hotspot. Lee et al. [23,24] presented a modified fin structure that could be the closest study we could observe on modified microchannel design. We have sincerely made a literature search and to the best of our knowledge we would like to affirm that this is the foremost study on split flow multizone microchannels for hotspot mitigation.

2 Problem Description

An exploded schematic view of the heatsink design is shown in Fig. 1. As shown in the figure, the chip sits at the bottom of the heatsink with a hotspot located at its center. A thin layer of TIM is placed between the chip and the heatsink. There are two inlet flows of coolant at the inlet manifold of the heatsink (QI and QII). Each of them is supplied by separate pumps. That is because the channels in this heatsink design are classified as middle class (on the hotspot) (class I) and side class (class II) channels. This design enables different flows of the coolant to be channeled to each of these channel classes. In this case, QI is the flow supplied to the middle channels over the hotspot and QII is the flow directed towards the side channels. It is worth mentioning that both QI and QII flow rates and the number of channels in each class can be controlled based on requirements. Similarly, two ports in the outlet manifold direct the coolant out of the system. Lastly, a cover sits atop the entire heatsink.

Fig. 1
An exploded schematic of the module showing different parts of the heatsink as well as channel classes with their flow rate and number of channels. Channels in class I are over the hotspot.
Fig. 1
An exploded schematic of the module showing different parts of the heatsink as well as channel classes with their flow rate and number of channels. Channels in class I are over the hotspot.
Close modal

Limitations of the design are categorized as geometrical and operational. The server designers who determine the available space in a server and a heatsink location dictate the geometrical limitations of the heatsink design, which are considered fixed and unchangeable during the optimization process. The geometrical parameters of the module are shown in Fig. 2 and provided in Table 1. Additionally, server designers determine the operational conditions that must be considered in the heatsink design and optimization process. The operational parameters are provided in Table 2. As seen in the table, the hotspot area is less than 3% of the chip area, while the hotspot heat flux is almost eight times that of the background. It is worth mentioning that, in this study, the values of the geometric and operational fixed parameters (fin size, coolant flow rate, hotspot to background heat flux ratio, etc.) are selected based on the previous researches on liquid-cooled heatsinks used for high power applications [16,2527].

Fig. 2
Geometric parameters of the module (given in Table 1)
Fig. 2
Geometric parameters of the module (given in Table 1)
Close modal
Table 1

Geometric parameters of the module

SymbolDefinitionValue (mm)
hBBase height3.00
hFFin height4.00
hTTIM height0.10
lChipChip length43.40
lFFin length61.80
lHHotspot length6.50
tFFin thickness0.18
wCNChannel network width38.20
wCChannel width0.20
wChipChip width26.80
wHHotspot width4.76
SymbolDefinitionValue (mm)
hBBase height3.00
hFFin height4.00
hTTIM height0.10
lChipChip length43.40
lFFin length61.80
lHHotspot length6.50
tFFin thickness0.18
wCNChannel network width38.20
wCChannel width0.20
wChipChip width26.80
wHHotspot width4.76
Table 2

Operational parameters of the module

SymbolDefinitionValue
TInCoolant inlet temperature300 (K)
q˙ChipChip total power250 (W)
q˙HSHotspot power45 (W)
q˙BGBackground power205 (W)
AChipChip area1163.12 (mm2)
AHSHotspot area30.94 (mm2)
AHS/AChipArea ratio2.66%
qHSHotspot heat flux145.5 (W/cm2)
qBGBackground heat flux18.1 (W/cm2)
qHS/qBGHeat flux ratio8.03
QCoolant flow rate33.33 (cm2/s)
SymbolDefinitionValue
TInCoolant inlet temperature300 (K)
q˙ChipChip total power250 (W)
q˙HSHotspot power45 (W)
q˙BGBackground power205 (W)
AChipChip area1163.12 (mm2)
AHSHotspot area30.94 (mm2)
AHS/AChipArea ratio2.66%
qHSHotspot heat flux145.5 (W/cm2)
qBGBackground heat flux18.1 (W/cm2)
qHS/qBGHeat flux ratio8.03
QCoolant flow rate33.33 (cm2/s)

The hotspot creates a drastically nonuniform temperature distribution on the chip and heatsink base, which reduces their reliability and lifespan.

3 Optimization Design Parameters

There are two design parameters for the present optimization problem. The first one is the percentage of channels considered for the hotspot area in class I
λI=NIN
(1)
where N=101 is the total number of channels and NI is the number of channels concentrated on the hotspot (class I). The second design parameter is the percentage of flow that class I channels receive, given by
ωI=QIQ
(2)

in which, Q is the total heatsink coolant flow rate (a constant value of 2lit/min) and QI is the portion of this flow rate that travels through the NI channels. The number of channels in class II and the coolant flowrate through them depend on the dimensionless values of λI and ωI. Values of λI range from 0.2 to 0.5 (0.2λI0.5), and ωI can take values ranging between 0.5 and 0.8 (0.5ωI0.8). A uniform flow is achieved when λI is at its maximum (λI=0.5) and ωI is at its minimum (ωI=0.5), whereas a maximum flow concentration on the hotspot is achieved when λI is at its minimum (λI=0.2) and ωI is at its maximum (ωI=0.8).

4 Response Parameters

Heatsink thermal resistivity is a measure of the heatsink thermal performance and is considered a thermal response parameter of the system. It is commonly formulated as
RTh=TBMaxTInq˙ChipAChip
(3)
where TBMax and TIn are the maximum base temperature and coolant inlet temperature, respectively, and q˙Chip is the total chip power. Similarly, heatsink heat spreading resistivity is defined as follows:
RSp=TBMaxTBAveq˙ChipAChip
(4)

where TBAve denotes the heatsink average base temperature. With the presence of hotspots in a module, removing the excessive heat and temperature distribution nonuniformity are of higher priority.

The main issue with hotspots is that they generate elevated temperature regions that lead to a high on-chip temperature gradient. This significantly reduces the reliability and lifespan of the chip. Therefore, it is necessary to address the temperature distribution nonuniformity. This is performed during the optimization process by defining the parameter as
ψ=TChipMaxTInTChipMinTIn
(5)

in which, TChipMax and TChipMin represent the maximum and minimum chip temperatures, respectively. To achieve a close-to-uniform temperature distribution profile, it is desirable that ψ be close to unity. A perfectly uniform temperature distribution is achieved when ψ=1, which means a completely uniform on-chip temperature distribution is obtained (due to TChipMax=TChipMin=TChipAve).

The next response parameter is the total pump power, given by
W˙P=W˙P,I+W˙P,II
(6)
W˙P,I=ΔPI·QI
(7)
W˙P,II=ΔPII·QII
(8)

where W˙P,I and W˙P,II are the powers for pump I and II, respectively, and ΔPI and ΔPII are the average pressure drops across channels in class I and class II, respectively.

The ratio of the channel Reynolds number for class I to that of the class II is estimated for each design point using Eq. (9)
ReDhIReDhII=(QIQII)(NIINI)
(9)

5 Numerical Analysis

5.1 Basic Assumptions.

The following assumptions are considered in the numerical modeling:

  1. steady, incompressible, three-dimensional flow;

  2. fully laminar flow in all channels (ReC<1000);

  3. gravity is negligible;

  4. constant thermophysical properties of both solid and fluid phases due to minute temperature variations;

  5. negligible viscous heat generation;

  6. negligible radiation heat transfer.

5.2 Governing Equations.

The following are the governing equations for a three-dimensional conduction–convection heat transfer problem in an incompressible, steady laminar flow regime [28,29].

Conservation of mass for the coolant fluid is given by
·V=0
(10)

where V is velocity vector field.

Equations of motion are expressed as
ρf(V·)V=∇P+μf2V
(11)

In Eq. (11), ρf is the density of the coolant fluid, P is the pressure field, and μf represents the absolute viscosity of the fluid.

Energy equation of the coolant fluid is
ρfcp,f(V·)Tf=kf2Tf
(12)

In which kf is the conductivity, and Tf is the temperature of the fluid.cp,f denotes the fluid specific heat capacity at constant pressure.

Finally, the energy equation for the solid phases (metal and TIM) is given by
2Ts=0
(13)

where Ts is the temperature in the solid phase.

5.3 Numerical Domain and Boundary Conditions.

Previous studies were able to reduce the computational domain to one-half of a unit cell of the mini/microchannel array, given the simplicity of a uniform heat flux boundary condition [16,17,30,31]. In this study, the presence of a hotspot introduces a nonuniform heat flux boundary condition. With the geometric symmetry of both the heatsink and the heat flux boundary condition in mind, the only way to conserve computational resources is to use a half symmetry model (Fig. 1).

A schematic of the computational domain and the corresponding boundary conditions are illustrated in Fig. 3. It is important to mention that this figure is not drawn to scale but is drawn with exaggerated dimensions to provide a clear view. The thickness of the fin and width of the channel are significantly larger in the figure than their actual dimensions in the model. In Fig. 3, the coolant enters the channels from the “velocity inlet” and leaves the heatsink from the “pressure constant” surfaces. In order to employ symmetry in the model, zero-gradient boundary conditions are imposed on the symmetry plane, shown as the violet dashed area in Fig. 3. Specific heat fluxes are applied to the chip background and the hotspot. The heat flux value at the hotspot is eight times greater than that of the background to enforce the nonuniform temperature distribution that is typical in such cases. No-slip hydraulic and thermally coupled boundary conditions are applied to the surfaces of the channel walls that are in contact with the flow. An adiabatic boundary condition is applied to all the surrounding surfaces.

Fig. 3
A schematic view of the computational domain with boundary conditions. Symmetric geometry and thermal boundary conditions make it possible to use a half-symmetry analysis.
Fig. 3
A schematic view of the computational domain with boundary conditions. Symmetric geometry and thermal boundary conditions make it possible to use a half-symmetry analysis.
Close modal

5.4 Numerical Method.

The system of equations for motion and energy in the computational domain is discretized on a staggered grid. It is performed using the finite volume method (FVM). Second-order upwind and central difference schemes are applied to the spatial discretization of the convective and diffusion terms, respectively. The governing equations are solved by ansysfluent and the field functions (pressure, velocity, and temperature) calculated using the iterative SIMPLEC algorithm [32].

Compared to other discretization methods (finite difference and finite element), the finite volume method is more compatible with the transport equations governing fluid flow and heat transfer. FVM discretizes the governing equations in an integral form by dividing the physical space into several arbitrary polyhedral control volumes (cells). Since FVM is based on the direct discretization of the conservation equations, mass, momentum, and energy are conserved by the numerical scheme that results in high robustness and calculation efficiency [33]. FVM is straightforward to implement on nonuniform/unstructured grids and flexible enough for complex geometries.

A mapped grid, made of hexahedral cells with fine divisions at the TIM layer and near the channel walls, is applied to the computational domain. Regions with finer cell dimensions correspond to regions of higher gradients.

A grid sensitivity analysis is performed based on both temperature nonuniformity of the chip and pump power given by Eqs. (5) and (6), respectively. Figure 4 shows the results of this analysis. As seen in the figure, when the number of cells is greater than 9 × 106, the pump power and chip temperature nonuniformity variations are less than 1%. As a result, the model is discretized to almost 9 × 106 cells in this study.

Fig. 4
Grid sensitivity analysis based on pump power and temperature nonuniformity (from Eqs. (5) and (6))
Fig. 4
Grid sensitivity analysis based on pump power and temperature nonuniformity (from Eqs. (5) and (6))
Close modal

5.5 Coolant, Metal, and Thermal Interface Material Properties.

There are two distinct solid phases in the numerical simulation, which possess different thermophysical properties. These two phases are TIM and copper. The thermophysical properties of copper at room temperature (300 K) are used. As for the TIM, the material properties of ShinEtsuMicroSi, X23-7783D, which is a commonly used TIM in electronic packaging, is incorporated into the model. Water is chosen as the coolant, with an inlet temperature set to 300 K. It is assumed that the thermophysical properties of the coolant are constant. The thermophysical properties of the coolant and the solid phases are provided in Table 3.

Table 3

Thermophysical properties of coolant and solid phases

cp(Jkg.K)kf(Wm.K)ρ(kgm3)μ(kgm.s)
Coolant41790.6139970.000855
Copper381387.68978
TIM100062500
cp(Jkg.K)kf(Wm.K)ρ(kgm3)μ(kgm.s)
Coolant41790.6139970.000855
Copper381387.68978
TIM100062500

6 Hydraulic Verification of the Numerical Simulation

Using the analytical approach below, the hydraulic pump power results across the DOE design points from the numerical simulation are verified. It is important to mention that the analytical solution is obtained for an asymptotic (including both hydrodynamically developing and fully developed) laminar flow in rectangular channels [34].

The hydraulic power is found from the following equations:
WP˙=ΔP·Q
(14)
ΔP=12ρfflCDhVC2
(15)
In this equation
f=ϕReDh
(16)
and Reynolds number
ReDh=ρfVCDhμf=4ρfQCμfpC
(17)
in which, Dh is the channel hydraulic diameter, and
ϕ=4[(3.44l+)2+(ζReDh)2]12
(18)
l+=lCDhReDh
(19)
and
ζReDh=24(1+ϵ)2[1192ϵπ5tanh(π2ϵ)]
(20)

where ϵ is the aspect ratio of the channel, which is the ratio of the longer to shorter edges of the channel cross section.

As seen in Fig. 5, the simulation results match very well with the analytical solution for all the cases that are studied in the analysis. Design model numbers (Table 4) in this figure refer to the pump power for each case in ascending order. The maximum and average discrepancies are less than 1.6% and 2.2%, respectively.

Fig. 5
Comparison of analytically and numerically calculated pump power at different FFD design points (Table 4) ordered in an ascending trend
Fig. 5
Comparison of analytically and numerically calculated pump power at different FFD design points (Table 4) ordered in an ascending trend
Close modal
Table 4

Calculated response parameters at FFD design points

No.ωIλININIIQI(lit/min)QII(lit/min)ReDhIReDhIIW˙P(W)RTh(K.cm2/W)ψ
10.50.551501.01.01.00.220.6259.202
20.60.551501.20.81.50.230.6149.113
30.50.441601.01.01.50.230.6139.105
40.70.551501.40.62.30.260.6089.007
50.60.441601.20.82.20.260.6059.002
60.50.331701.01.02.30.260.6049.001
70.80.551501.60.43.90.300.6028.845
80.70.441601.40.63.40.300.6008.870
90.60.331701.20.83.40.320.6008.890
100.50.221801.01.03.80.340.6028.917
110.80.441601.60.45.90.370.5998.653
120.70.331701.40.65.30.390.5998.738
130.60.221801.20.85.70.440.6028.806
140.80.331701.60.49.00.490.6028.472
150.70.221801.40.68.90.560.6078.644
160.80.221801.60.415.20.720.6198.339
No.ωIλININIIQI(lit/min)QII(lit/min)ReDhIReDhIIW˙P(W)RTh(K.cm2/W)ψ
10.50.551501.01.01.00.220.6259.202
20.60.551501.20.81.50.230.6149.113
30.50.441601.01.01.50.230.6139.105
40.70.551501.40.62.30.260.6089.007
50.60.441601.20.82.20.260.6059.002
60.50.331701.01.02.30.260.6049.001
70.80.551501.60.43.90.300.6028.845
80.70.441601.40.63.40.300.6008.870
90.60.331701.20.83.40.320.6008.890
100.50.221801.01.03.80.340.6028.917
110.80.441601.60.45.90.370.5998.653
120.70.331701.40.65.30.390.5998.738
130.60.221801.20.85.70.440.6028.806
140.80.331701.60.49.00.490.6028.472
150.70.221801.40.68.90.560.6078.644
160.80.221801.60.415.20.720.6198.339

7 Interpretation of Computational Fluid Dynamics Results

The FFD with four levels is applied to the design parameters to generate a total of 16 points in the design space. The total pump power, thermal resistivity, and temperature nonuniformity of the chip are calculated using CFD at the design points.

Table 4 shows the response parameters (W˙P,ψ, and RTh) and ReDhI/ReDhII at the FFD points. For model No. 1 with a uniform flow distribution (ReDhI/ReDhII=1), the minimum pump power is achieved, but both thermal resistivity and temperature nonuniformity are at their peaks. The minimum value of ψ is obtained when pump power and the Reynolds number ratio (ReDhI/ReDhII=15) are at their maximum. This corresponds to the design point of model No. 16 in which 80% of the flowrate is pumped to 20% of the channels on the hotspot (class I). Increasing ReDhI/ReDhII causes the pump power to increase and the temperature nonuniformity to decrease.

Although the maximum thermal resistivity coincides with the maximum temperature nonuniformity and minimum pump power (model No. 1), it reaches a minimum value when pump power and temperature nonuniformity are at intermediate values (model No. 11). When 80% of the flowrate travels through 40% of the channels in the middle (class I), as is the case for model No. 11, the maximum chip/base temperature is at its lowest.

Figure 6 shows the variation of ψ, RTh, and RSp (Eq. (4)) with pump power at the FFD design points. As shown in the figure, ψ and RSp have similar trends and minimums at the same pump power. This implies that a more uniform chip temperature distribution is achieved when heat spreads more efficiently from the hotspot. However, RTh shows a different trend. The thermal resistivity of the heatsink is calculated by
RTh=RDir+RSp+RConv+RCal
(21)
where RDir is direct conduction resistivity across the heatsink base. RConv is convective resistivity and RCal, known as caloric thermal resistivity, is expressed as
RCal=AChipm˙cp,f
(22)
Fig. 6
Variation of (a) chip temperature nonuniformity, (b) thermal resistivity, and (c) heat spreading resistivity with pump power at fast Fourier transform design points (Table 4)
Fig. 6
Variation of (a) chip temperature nonuniformity, (b) thermal resistivity, and (c) heat spreading resistivity with pump power at fast Fourier transform design points (Table 4)
Close modal

By concentrating the flowrate to the channels above the hotspot (increasing W˙P) and deviating from a uniform flow distribution, the heat spreading resistivity (RSp), at the right-hand side of Eq. (21) decreases, which results in a reduction in RTh. This reduction in RTh initiates at uniform flow distribution of model No. 1 and continues to model No. 11 where 80% of the flow travels through 40% of the channels on the hotspot. At this point, RTh reaches its minimum. Concentrating more flow to the channels on the hotspot increases the RTh and changes its trend from descending to ascending (Fig. 6(b)). The reason is that concentrating most of the flow to the middle channels attenuates the thermal performance of the side channels. This results in an increase in the heatsink RConv. Figure 7 illustrates how the base of the heatsink spreads the heat generated by the chip hotspot and background to the channels with low flowrate at the sides. Using the temperature profile on the plane in the middle of the channels perpendicular to the flow direction, Fig. 8 compares the spread of the heat for two design models, No. 1 and No. 16, corresponding to maximum and minimum RSp (or ψ), respectively. Comparing Figs. 8(a) and 8(b), it is concluded that concentrating the flow on the hotspot improves the spread of the heat through the heatsink base and fins, which, in turn, results in more uniform chip temperature. Figure 9 shows the corresponding temperature profiles of the models No. 1 and No. 16 on a plane close to the root of the fins. These temperature profiles confirm there is an enhanced spread of heat in model No. 16 compared to that of model No. 1.

Fig. 7
A schematic of heatsink cross section when flow rate is concentrated on the top of the hotspot. A part of the heat conducts to side channels with low flow rate and heat transfer performance.
Fig. 7
A schematic of heatsink cross section when flow rate is concentrated on the top of the hotspot. A part of the heat conducts to side channels with low flow rate and heat transfer performance.
Close modal
Fig. 8
Temperature profile on the middle heatsink plane perpendicular to the flow direction. (a) Design model No. 1 (ωI=0.5 and λI=0.5) and (b) design model No. 16 (ωI=0.8 and λI=0.2).
Fig. 8
Temperature profile on the middle heatsink plane perpendicular to the flow direction. (a) Design model No. 1 (ωI=0.5 and λI=0.5) and (b) design model No. 16 (ωI=0.8 and λI=0.2).
Close modal
Fig. 9
Temperature profile on the plane close to the fin roots and parallel to the flow direction. (a) Design model No. 1 (ωI=0.5 and λI=0.5) and (b) design model No. 16 (ωI=0.8 and λI=0.2).
Fig. 9
Temperature profile on the plane close to the fin roots and parallel to the flow direction. (a) Design model No. 1 (ωI=0.5 and λI=0.5) and (b) design model No. 16 (ωI=0.8 and λI=0.2).
Close modal

8 Parametric Study

8.1 Pump Power.

Figure 10 shows the variation of pump power with the design parameters λI and ωI. As expected, design model No. 1 with a uniform flow distribution possesses the minimum pump power. Increasing ωI or decreasing λI concentrates the flow to the channels on the hotspot, which results in a higher pressure drop and total pump power.

Fig. 10
The effect of design parameters ωI and λI on the total pump power
Fig. 10
The effect of design parameters ωI and λI on the total pump power
Close modal

8.2 Chip Temperature Nonuniformity.

Figure 11 illustrates the effect of on-chip flow concentration on the chip temperature nonuniformity and heatsink heat spreading resistivity. The figure shows that both chip temperature nonuniformity and heat spreading resistivity follow the same trend with ωI at different values of λI. Focusing the flow to the channels on the hotspot enhances heat spreading from it and results in a more uniform chip temperature profile.

Fig. 11
The effect of design parameters ωI and λI on the (a) chip temperature nonuniformity and (b) heatsink heat spreading resistivity
Fig. 11
The effect of design parameters ωI and λI on the (a) chip temperature nonuniformity and (b) heatsink heat spreading resistivity
Close modal

8.3 Thermal Resistivity.

The effect of ωI and λI on the total heatsink thermal resistivity is demonstrated in Fig. 12. The graph for λI=0.5 and λI=0.4, which correspond to 50% and 40% of the channels on the hotspot, shows that RTh is reduced by increasing ωI. In these cases (λI=0.5 and λI=0.4), the low number of side channels (class II), makes it so that the concentration of the flow on the middle channels does not increase heatsink convective thermal resistivity remarkably. Consequently, by increasing ωI and decreasing RSp (Eq. (21)), the heatsink thermal resistivity is reduced. In other words, for λI=0.5 and λI=0.4, RSp is the only term on the right side of Eq. (21) that depends on ωI. As Fig. 12 shows for λI=0.3, the trend of thermal resistivity is almost flat. For λI=0.3, when 60% of the channels are placed on the sides, RSp andRConv in Eq. (21) depend on ωI. Increasing ωI decreases RSp, on the one hand, and increases RConv on the other hand, which results in nearly constant RTh. This is the reason that when λI=0.3, thermal resistvity (RTh) is a weak function of ωI. An ascending trend for thermal resistivity is recognizable at λI=0.2. As mentioned in Sec. 7, due to the large number of channels located on the sides (class II), increasing ωI weakens convective heat transfer and flow bulk motion in these channels which increases heatsink RConv and results in RTh growth.

Fig. 12
The effect of design parameters ωI and λI on the heatsink thermal resistivity
Fig. 12
The effect of design parameters ωI and λI on the heatsink thermal resistivity
Close modal

9 Optimization

Two multi-objective optimizations are performed based on the designed heatsink. First ωI and λI are optimized to minimize system pump power (W˙P) and chip temperature nonuniformity (ψ). This optimization problem can be expressed as
Minimize:W˙P=F1(ωI,λI)
(23)
Minimize:ψ=F2(ωI,λI)
(24)
Subjected to
0.5ωI0.8
(25a)
0.2λI0.5
(25b)
In the second optimization problem, pump power (W˙P) and thermal resistivity (RTh) are considered as the response parameters
Minimize:W˙P=F1(ωI,λI)
(26)
Minimize:RTh=F3(ωI,λI)
(27)
Subjected to
0.5ωI0.8
(28a)
0.2λI0.5
(28b)
Regression analysis is employed to estimate the functions of F1, F2, and F3.
W˙P=exp(1.377+1.524ωI2.738λI+1.489ωI2+5.184λI24.816ωIλI)
(29)
ψ=8.613+1.953ωI0.334λI3.360ωI2+2.403ωIλI
(30)
RTh=exp(0.4660.169ωI+0.079λI+0.297ωI2+0.584λI20.071ωIλI)
(31)

The calculated regression functions are verified both statistically, considering coefficients of determination (R2s), and numerically by modeling a few intermediate design points.

The obtained values for R2 and R2(Adjusted) for W˙P are 99.8% and 99.7%, respectively. Similar R2 values of 99.4% and 99.2% are calculated for ψ. Finally, the values of 96.8% and 94.6% are acquired for R2 and R2(Adjusted) of RTh.

In addition to the statistical validation, the accuracy of the regression functions in Eqs. (29)(31) is confirmed numerically. Three intermediate designs are chosen from the available range of design parameters. Table 5 compares the results of the numerical modeling with regression function predictions. As Table 5 shows the discrepancy is less than 1% at all intermediate points for all regression functions.

Table 5

Regression functions verification by CFD results in three intermediate points

Design ParametersW˙P(W)ψRTh (K.cm2/W)
ωIλISimulationRegressionDiscrepancySimulationRegressionDiscrepancySimulationRegressionDiscrepancy
0.550.440.23330.23350.09%9.09819.10530.08%0.61250.61150.16%
0.650.340.31600.31570.09%8.86768.88040.14%0.59960.59950.02%
0.750.240.53500.53720.40%8.55608.54010.19%0.60490.60670.35%
Design ParametersW˙P(W)ψRTh (K.cm2/W)
ωIλISimulationRegressionDiscrepancySimulationRegressionDiscrepancySimulationRegressionDiscrepancy
0.550.440.23330.23350.09%9.09819.10530.08%0.61250.61150.16%
0.650.340.31600.31570.09%8.86768.88040.14%0.59960.59950.02%
0.750.240.53500.53720.40%8.55608.54010.19%0.60490.60670.35%

This study also tried linear regression functions (quadratic forms) for all the response parameters. Still, in the case of pump power and thermal resistivity, they failed to satisfy linearity and homoscedasticity assumptions. To fix the problem, a logarithmic transformation is used for pump power and thermal resistivity (Eqs. (29) and (31)) [35]. A regression function is valid only if it satisfies the standard regression assumptions [36] expressed as:

  • Linearity assumption: The relation between response and design parameters can be estimated by standard multiregression when it is linear in nature.

  • Normality assumption: The regression function errors should be normally distributed.

  • Homoscedasticity (constant variance) assumption: The variance of the regression function errors should be constant and independent of their mean.

In addition to the above assumptions, the regression functions should be checked for outliers and influential points. A more detailed analysis of regression diagnostics can be found in Ref. [36].

Although logarithmic transformation increases the normality and decreases the effects of outliers, it is mostly used to enhance data linearity and remove heteroscedasticity. A regression diagnostics analysis is performed for the new nonlinear regression functions (Eqs. (29) and (31)) that shows they satisfy standard regression assumptions thoroughly.

A weight optimization approach is used to solve the above multi-objective optimization problems of W˙Pψ and W˙PRTh. The objective functions are defined as weighted combinations of the response parameters
ΞW˙Pψ=ξW˙P(W˙PW˙PMin)+ξψ(ψψMin)
(32)
ΞW˙PRTh=ξW˙P(W˙PW˙PMin)+ξRTh``(RThRThMin)
(33)

In Eqs. (32) and (33), W˙P, ψ, and RTh are substituted by Eqs. (29)(31), respectively. W˙PMin, ψMin, and RThMin are the minimum values of the regression functions of Eqs. (29)(31).ξW˙P, ξψ, and ξRTh are the weight factors assigned to W˙P, ψ, and RTh, respectively. Weight factors are chosen according to the designer's priority and the importance they place on the design parameters. The iterative algorithm of JAYA [37] is used for minimizing the objective functions expressed in Eqs. (32) and (33). Figure 13 shows the Pareto curve for the optimization problem of W˙Pψ associated with the values of thermal resistivity (RTh) at these points. Table 6 indicates the design parameters of these optimal designs. The Pareto curve for the optimization problem of W˙PRTh together with the values of ψ at these points are shown in Fig. 14. Table 7 displays ωI and λI at the obtained optimal points.

Fig. 13
The results of the multi-objective optimization problem of W˙P−ψ with corresponding thermal resistivity values
Fig. 13
The results of the multi-objective optimization problem of W˙P−ψ with corresponding thermal resistivity values
Close modal
Fig. 14
The results of the multi-objective optimization problem of W˙P−RTh″ with corresponding chip temperature nonuniformity values
Fig. 14
The results of the multi-objective optimization problem of W˙P−RTh″ with corresponding chip temperature nonuniformity values
Close modal
Table 6

Design and response parameters at the optimal points of the multi-objective optimization problem of W˙Pψ

ωIλIW˙P(W)ψRTh (K.cm2/W)
10.50.500.229.180.626
20.50.400.239.100.612
30.80.500.308.820.602
40.80.460.328.720.600
50.80.420.358.690.599
60.80.380.388.630.600
70.80.340.438.560.602
80.80.310.488.510.604
90.80.270.548.460.607
100.80.230.658.390.612
110.80.200.738.340.616
ωIλIW˙P(W)ψRTh (K.cm2/W)
10.50.500.229.180.626
20.50.400.239.100.612
30.80.500.308.820.602
40.80.460.328.720.600
50.80.420.358.690.599
60.80.380.388.630.600
70.80.340.438.560.602
80.80.310.488.510.604
90.80.270.548.460.607
100.80.230.658.390.612
110.80.200.738.340.616
Table 7

Design and response parameters at the optimal points of the multi-objective optimization problem of W˙PRTh

ωIλIW˙P(W)ψRTh (K.cm2/W)
10.50.500.2199.180.626
20.50.480.2199.160.623
30.50.440.2229.130.617
40.50.420.2269.110.613
50.50.370.2389.070.608
60.50.360.2439.060.607
70.540.360.2549.040.604
80.590.370.2708.990.602
90.650.370.2948.920.600
100.690.380.3128.860.599
110.720.200.3378.780.599
ωIλIW˙P(W)ψRTh (K.cm2/W)
10.50.500.2199.180.626
20.50.480.2199.160.623
30.50.440.2229.130.617
40.50.420.2269.110.613
50.50.370.2389.070.608
60.50.360.2439.060.607
70.540.360.2549.040.604
80.590.370.2708.990.602
90.650.370.2948.920.600
100.690.380.3128.860.599
110.720.200.3378.780.599

10 Sensitivity Analysis

Sensitivity analysis is a mathematical technique used to evaluate how the response parameters are influenced by uncertainties in the design parameters. The analysis can be performed both by taking partial derivatives of regression functions (Eqs. (29)(31)) with respect to independent (design) parameters [38] or the results of F-values obtained from the analysis of variance (ANOVA) [39,40]. The F-value is a non-negative real number defined as the ratio of the mean regression sum of squares to mean error sum of squares. The higher F-value that a term has in the analysis of variance, the more significant its effect on the response parameter. Using the results of F-values calculated through the analysis of variance for generating regression functions (Eqs. (29)(31)), a sensitivity analysis was done to estimate the influence of design parameters and their combinations (regression functions terms) on the response parameters. Figure 15 shows the results of this sensitivity analysis. As seen in the figure, quadratic and interaction terms have greater impacts on the response parameters than linear terms. The interaction term ω1λ1 is the most influential term for this problem. It is attributed to the fact that the flow concentration on the high heat flux area includes the contribution of both design parameters simultaneously. In other words, to create a concentrated flow distribution, one must specify both the number of channels and the flowrate in class I.

Fig. 15
A sensitivity analysis of the design parameters and their quadratic and interaction combinations (temperature nonuniformity is independent of λI2 (Eq. (30))
Fig. 15
A sensitivity analysis of the design parameters and their quadratic and interaction combinations (temperature nonuniformity is independent of λI2 (Eq. (30))
Close modal

11 Conclusion

In this study, the authors numerically demonstrated a novel concept for energy efficient, hotspot-targeted liquid cooling of multicore microprocessors. The multizone heatsink consisted of varyingly distributed channels and flow over the hotspots and the background as a potential design to maximize temperature uniformity over the chip surface. Uniform flow distribution over the entire chip surface was taken as the baseline for comparison. Detailed three dimensional CFD simulations were performed for a series of designs and then a surrogate model representing the actual physics was developed. The surrogate model was used for optimization. The optimization study results made it evident that concentrating the flow on the small areas of elevated heat flux improves the spread of heat from these regions and results in a more uniform chip temperature. An improvement of 10% was achieved in the chip temperature uniformity compared to the baseline case. However, the results showed that concentrating the flow on the hotspot decreased the system thermal resistivity up to some point (where RSp and RConv are minimum) after which, its effect followed an adversely reverse course. In other words, there was a relative minimum for thermal resistivity with respect to flow concentration on the hotspot.

Funding Data

  • NSF Industry/University Cooperative Research Center (Award No. IIP-1738793; Funder ID: 10.13039/100000001).

Nomenclature

A =

area (m2)

cp =

specific heat capacity (J/kgK)

Dh =

hydraulic diameter (m)

f =

friction factor (Eqs. (15) and (16))

F =

general function (Eqs. (23), (24), (26), and (27))

h =

height (m)

k =

thermal conductivity (J/mK)

l =

length (m)

m˙ =

mass flow rate (kg/s)

p =

perimeter (m)

P =

pressure (Pa)

Q =

coolant volumetric flow rate (m3/s)

q˙ =

power (W)

q =

heat Flux (W/m2)

R =

thermal resistivity (Km2/W)

R2 =

coefficient of determination

Re =

Reynolds number

t =

thickness (m)

T =

temperature (K)

V =

velocity (m/s)

w =

width (m)

W˙ =

pump power (W)

Greek Symbols

Greek Symbols
ε =

channel aspect ratio

λ =

normalized channel number

μ =

absolute viscosity (kg/ms)

Ξ =

objective function (Eqs. (32 and (33))

ξ =

weight factor (Eqs. (32) and (33))

ρ =

density (kg/m3)

ψ =

chip temperature nonuniformity

ω =

normalized flow rate

Subscripts

Subscripts
B =

base

BG =

background

C =

channel

Chip =

chip

Cal =

caloric

CN =

channel network

Conv =

convective

Dir =

direct

F =

fin

f =

fluid

HS =

hotspot

In =

inlet

P =

pump

S =

solid

Sp =

spreading

Th =

thermal

I =

class I

II =

class II

Subscripts

Subscripts
Ave =

average

Min =

minimum

Max =

maximum

I =

class I

II =

class II

References

1.
Bar-Cohen
,
A.
, and
Wang
,
P.
,
2010
, “
On-Chip Thermal Management and Hot-Spot Remediation
,”
Nano-Bio-Electronic, Photonic and MEMS Packaging
,
C. P.
Wong
,
K.-S.
Moon
, and
Y.
(Grace) Li
, eds.,
Springer US
,
Boston, MA
, pp.
349
429
.10.1007/978-1-4419-0040-1_12
2.
Bulman
,
G.
,
Barletta
,
P.
,
Lewis
,
J.
,
Baldasaro
,
N.
,
Manno
,
M.
,
Bar-Cohen
,
A.
, and
Yang
,
B.
,
2016
, “
Superlattice-Based Thin-Film Thermoelectric Modules With High Cooling Fluxes
,”
Nat. Commun.
,
7
(
1
), pp.
1
7
.10.1038/ncomms10302
3.
Hao
,
X.
,
Peng
,
B.
,
Xie
,
G.
, and
Chen
,
Y.
,
2016
, “
Efficient on-Chip Hotspot Removal Combined Solution of Thermoelectric Cooler and Mini-Channel Heat Sink
,”
Appl. Therm. Eng.
,
100
, pp.
170
178
.10.1016/j.applthermaleng.2016.01.131
4.
Marcinichen
,
J. B.
,
d'Entremont
,
B. P.
,
Thome
,
J. R.
,
Bulman
,
G.
,
Lewis
,
J.
, and
Venkatasubramanian
,
R.
,
2014
, “
Thermal Management of Ultra Intense Hot Spots With Two-Phase Multi-Microchannels and Embedded Thermoelectric Cooling
,”
ASME
Paper No. IPACK2013-73276.10.1115/IPACK2013-73276
5.
Wang
,
P.
,
Yang
,
B.
, and
Bar-Cohen
,
A.
,
2009
, “
Mini-Contact Enhanced Thermoelectric Coolers for on-Chip Hot Spot Cooling
,”
Heat Transfer Eng.
,
30
(
9
), pp.
736
743
.10.1080/01457630802678391
6.
Bachmann
,
C.
, and
Bar-Cohen
,
A.
,
2008
, “
Hotspot Remediation With Anisotropic Thermal Interface Materials
,”
11th Intersociety Conference on Thermal and Thermomechanical Phenomena in Electronic Systems
, Orlando, FL, May 28–31, pp.
238
247
.10.1109/ITHERM.2008.4544276
7.
Ansari
,
D.
, and
Kim
,
K.-Y.
,
2018
, “
Hotspot Thermal Management Using a Microchannel-Pinfin Hybrid Heat Sink
,”
Int. J. Therm. Sci.
,
134
, pp.
27
39
.10.1016/j.ijthermalsci.2018.07.043
8.
Ansari
,
D.
, and
Kim
,
K.-Y.
,
2019
, “
Hotspot Management Using a Hybrid Heat Sink With Stepped Pin-Fins
,”
Numer. Heat Transfer, Part A: Appl.
,
75
(
6
), pp.
359
380
.10.1080/10407782.2019.1599272
9.
Gonzalez-Hernandez
,
J.-L.
,
Kandlikar
,
S. G.
, and
Hernandez-Guerrero
,
A.
,
2016
, “
Performance Assessment Comparison of Variable Fin Density Microchannels With Offset Configurations
,”
Heat Transfer Eng.
,
37
(
16
), pp.
1369
1381
.10.1080/01457632.2015.1136146
10.
Rubio-Jimenez
,
C. A.
,
Kandlikar
,
S. G.
, and
Hernandez-Guerrero
,
A.
,
2012
, “
Numerical Analysis of Novel Micro Pin Fin Heat Sink With Variable Fin Density
,”
IEEE Trans. Compon., Packag. Manuf. Technol.
,
2
(
5
), pp.
825
833
.10.1109/TCPMT.2012.2189925
11.
Rubio-Jimenez
,
C. A.
,
Kandlikar
,
S. G.
, and
Hernandez-Guerrero
,
A.
,
2013
, “
Performance of Online and Offset Micro Pin-Fin Heat Sinks With Variable Fin Density
,”
IEEE Trans. Compon., Packaging Manuf. Technol.
,
3
(
1
), pp.
86
93
.10.1109/TCPMT.2012.2225143
12.
Green
,
C.
,
Fedorov
,
A. G.
, and
Joshi
,
Y. K.
,
2009
, “
Fluid-to-Fluid Spot-to-Spreader (F2/S2) Hybrid Heat Sink for Integrated Chip-Level and Hot Spot-Level Thermal Management
,”
ASME J. Electron. Packag
,
131
(
2
), p.
025002
.10.1115/1.3104029
13.
Sharma
,
C. S.
,
Tiwari
,
M. K.
,
Zimmermann
,
S.
,
Brunschwiler
,
T.
,
Schlottig
,
G.
,
Michel
,
B.
, and
Poulikakos
,
D.
,
2015
, “
Energy Efficient Hotspot-Targeted Embedded Liquid Cooling of Electronics
,”
Appl. Energy
,
138
, pp.
414
422
.10.1016/j.apenergy.2014.10.068
14.
Kumar
,
S.
, and
Singh
,
P. K.
,
2019
, “
A Novel Approach to Manage Temperature Non-Uniformity in Minichannel Heat Sink by Using Intentional Flow Maldistribution
,”
Appl. Therm. Eng.
,
163
, p.
114403
.10.1016/j.applthermaleng.2019.114403
15.
Montgomery
,
D. C.
,
2019
,
Design and Analysis of Experiments
,
Wiley
,
New York
.
16.
Rao
,
R. V.
,
More
,
K. C.
,
Taler
,
J.
, and
Ocłoń
,
P.
,
2016
, “
Dimensional Optimization of a Micro-Channel Heat Sink Using Jaya Algorithm
,”
Appl. Therm. Eng.
,
103
, pp.
572
582
.10.1016/j.applthermaleng.2016.04.135
17.
Kulkarni
,
K.
,
Afzal
,
A.
, and
Kim
,
K. Y.
,
2016
, “
Multi-Objective Optimization of a Double-Layered Microchannel Heat Sink With Temperature-Dependent Fluid Properties
,”
Appl. Therm. Eng.
,
99
(
25
), pp.
262
272
.10.1016/j.applthermaleng.2016.01.039
18.
Naqiuddin
,
N. H.
,
Saw
,
L. H.
,
Yew
,
M. C.
,
Yusof
,
F.
,
Poon
,
H. M.
,
Cai
,
Z.
, and
Thiam
,
H. S.
,
2018
, “
Numerical Investigation for Optimizing Segmented Micro-Channel Heat Sink by Taguchi-Grey Method
,”
Appl. Energy
,
222
, pp.
437
450
.10.1016/j.apenergy.2018.03.186
19.
Hadad
,
Y.
,
Ramakrishnan
,
B.
,
Pejman
,
R.
,
Rangarajan
,
S.
,
Chiarot
,
P. R.
,
Pattamatta
,
A.
, and
Sammakia
,
B.
,
2019
, “
Three-Objective Shape Optimization and Parametric Study of a Micro-Channel Heat Sink With Discrete Non-Uniform Heat Flux Boundary Conditions
,”
Appl. Therm. Eng.
,
150
, pp.
720
737
.10.1016/j.applthermaleng.2018.12.128
20.
Ermagan
,
H.
, and
Rafee
,
R.
,
2018
, “
Geometric Optimization of an Enhanced Microchannel Heat Sink With Superhydrophobic Walls
,”
Appl. Therm. Eng.
,
130
, pp.
384
394
.10.1016/j.applthermaleng.2017.11.039
21.
Nakhchi
,
M. E.
,
2019
, “
Experimental Optimization of Geometrical Parameters on Heat Transfer and Pressure Drop Inside Sinusoidal Wavy Channels
,”
Therm. Sci. Eng. Prog.
,
9
, pp.
121
131
.10.1016/j.tsep.2018.11.006
22.
Chen
,
H.-T.
,
Horng
,
J.-T.
,
Chen
,
P.-L.
, and
Hung
,
Y.-H.
,
2004
, “
Optimal Design for PPF Heat Sinks in Electronics Cooling Applications
,”
ASME J. Electron. Packag
,
126
(
4
), pp.
410
422
.10.1115/1.1826078
23.
Lee
,
P.-S.
, and
Garimella
,
S. V.
,
2008
, “
Hot-Spot Thermal Management With Flow Modulation in a Microchannel Heat Sink
,”
ASME
Paper No. IMECE2005-79562.10.1115/IMECE2005-79562
24.
Lee
,
Y. J.
,
Lee
,
P. S.
, and
Chou
,
S. K.
,
2013
, “
Hotspot Mitigating With Obliquely Finned Microchannel Heat Sink—an Experimental Study
,”
IEEE Trans. Compon., Packag. Manuf. Technol.
,
3
(
8
), pp.
1332
1341
.10.1109/TCPMT.2013.2244164
25.
Qu
,
W.
, and
Mudawar
,
I.
,
2002
, “
Experimental and Numerical Study of Pressure Drop and Heat Transfer in a Single-Phase Micro-Channel Heat Sink
,”
Int. J. Heat Mass Transfer
,
45
(
12
), pp.
2549
2565
.10.1016/S0017-9310(01)00337-4
26.
Yang
,
X.-H.
,
Tan
,
S.-C.
,
Ding
,
Y.-J.
, and
Liu
,
J.
,
2017
, “
Flow and Thermal Modeling and Optimization of Micro/Mini-Channel Heat Sink
,”
Appl. Therm. Eng.
,
117
, pp.
289
296
.10.1016/j.applthermaleng.2016.12.089
27.
Lin
,
L.
,
Zhao
,
J.
,
Lu
,
G.
,
Wang
,
X.-D.
, and
Yan
,
W.-M.
,
2017
, “
Heat Transfer Enhancement in Microchannel Heat Sink by Wavy Channel With Changing Wavelength/Amplitude
,”
Int. J. Therm. Sci.
,
118
, pp.
423
434
.10.1016/j.ijthermalsci.2017.05.013
28.
Bejan
,
A.
,
2013
,
Convection Heat Transfer
,
Wiley
,
New York
.
29.
Hahn
,
D. W.
, and
Qzisik
,
M. N.
,
2012
,
Heat Conduction
,
Wiley
,
Hoboken, NJ
.
30.
Lin
,
L.
,
Chen
,
Y. Y.
,
Zhang
,
X. X.
, and
Wang
,
X. D.
,
2014
, “
Optimization of Geometry and Flow Rate Distribution for Double-Layer Microchannel Heat Sink
,”
Int. J. Therm. Sci.
,
78
, pp.
158
168
.10.1016/j.ijthermalsci.2013.12.009
31.
Xie
,
X. L.
,
Liu
,
Z. J.
,
He
,
Y. L.
, and
Tao
,
W. Q.
,
2009
, “
Numerical Study of Laminar Heat Transfer and Pressure Drop Characteristics in a Water-Cooled Minichannel Heat Sink
,”
Appl. Therm. Eng.
,
29
(
1
), pp.
64
74
.10.1016/j.applthermaleng.2008.02.002
32.
Versteeg
,
H. K.
, and
Malalasekera
,
W.
,
2007
,
An Introduction to Computational Fluid Dynamics: The Finite Volume Method
,
Pearson Prentice Hall
,
London, UK
.
33.
Liu
,
R.
,
Wang
,
D.
,
Zhang
,
X.
,
Li
,
W.
, and
Yu
,
B.
,
2013
, “
Comparison Study on the Performances of Finite Volume Method and Finite Difference Method
,”
J. Appl. Math.
,
2013
, pp.
1
10
.10.1155/2013/596218
34.
Muzychka
,
Y. S.
, and
Yovanovich
,
M. M.
,
2009
, “
Pressure Drop in Laminar Developing Flow in Noncircular Ducts: A Scaling and Modeling Approach
,”
ASME J. Fluids Eng.
,
131
(
11
), p.
111105
.10.1115/1.4000377
35.
Graybill
,
F. A.
, and
Hariharan
,
K. I.
,
1994
,
Regression Analysis Concepts and Applications
,
Duxbury Press
,
Belmont, CA
.
36.
Chatterjee
,
S.
, and
Hadi
,
A. S.
,
2012
,
Regression Analysis by Example
,
Wiley
,
Hoboken, NJ
.
37.
Venkata Rao
,
R.
,
2016
, “
Jaya: A Simple and New Optimization Algorithm for Solving Constrained and Unconstrained Optimization Problems
,”
Int. J. Ind. Eng. Comput.
,
7
, pp.
19
34
.10.5267/j.ijiec.2015.8.004
38.
Zheng
,
N.
,
Liu
,
P.
,
Liu
,
Z.
, and
Liu
,
W.
,
2017
, “
Numerical Simulation and Sensitivity Analysis of Heat Transfer Enhancement in a Flat Heat Exchanger Tube With Discrete Inclined Ribs
,”
Int. J. Heat Mass Transfer
,
112
, pp.
509
520
.10.1016/j.ijheatmasstransfer.2017.05.019
39.
Chiang
,
K.-T.
,
Chou
,
C.-C.
, and
Liu
,
N.-M.
,
2009
, “
Application of Response Surface Methodology in Describing the Thermal Performances of a Pin-Fin Heat Sink
,”
Int. J. Therm. Sci.
,
48
(
6
), pp.
1196
1205
.10.1016/j.ijthermalsci.2008.10.009
40.
Chiang
,
K.-T.
,
2007
, “
Modeling and Optimization of Designing Parameters for a Parallel-Plain Fin Heat Sink With Confined Impinging Jet Using the Response Surface Methodology
,”
Appl. Therm. Eng.
,
27
(
14–15
), pp.
2473
2482
.10.1016/j.applthermaleng.2007.02.004