Understanding Inconsistencies in Thermohydraulic Characteristics Between Experimental and Numerical Data for DI Water Flow Through a Rectangular Microchannel

The behavior of fluid flow in microchannels has been a longstanding topic of interest in the field of fluid dynamics and heat transfer [1]. Specifically, researchers investigated whether the fluid flow in microchannels behaves similarly to a conventional laminar flow or whether there are deviations from this theoretical description and what the reasons for these possible deviations might be. This question is of great importance because of the huge potential of microchannel applications in various engineering systems, such as cooling systems for electronic devices or in chemical process engineering to provide high-quality steam. To better understand the behavior of fluid flow and heat transfer in microchannels, researchers have conducted numerous experimental studies [2,3]. The results of these studies vary, with some showing deviations from the laminar flow theory [4,5] and others supporting its validity [6,7].

Some studies have shown that flow and heat transfer in the microscale often deviate from conventional laminar flow theory due to scaling effects such as surface roughness and entrance effects. For example, Peng et al. [8] investigated the flow characteristics of water flowing through rectangular microchannels and found that the friction factor (⁠$f$⁠) deviates significantly from the conventional laminar flow theory. Later, Peng and Peterson [4] studied the effect of channel size on single-phase flow in heated microchannels using water and methanol as working fluids. They claimed that laminar-to-turbulent flow transition occurs at a Reynolds number (⁠$Re$⁠) of $Re≥$ 300 and a fully developed turbulent flow regime was first obtained at $Re>$ 1000. In another study, Pfund et al. [9] conducted an experiment to measure pressure drop across a microchannel and identified different flow regimes for water as a working fluid. The onset of laminar-to-turbulent flow transition was found at a $Re$ range of 1500 − 2200, and the Poiseuille number (⁠$Po=fRe$⁠) was found to be significantly higher than the theoretical value for fully developed laminar flow, but the authors remained uncertain about which parameter, channel geometry or surface roughness, had a stronger effect on $Po$ due to experimental uncertainty.

Contrary to the observations discussed above, some studies have found no significant deviations from conventional laminar flow theory. For instance, Judy et al. [6] found no significant deviations from conventional laminar flow theory when investigating single-phase pressure drop in circular and square microchannels with diameters ranging from 0.015 to 0.15 mm and lengths ranging from 36 to 300 mm. Mokrani et al. [10] also conclude that conventional laws and correlations are applicable to low aspect ratio rectangular microchannels with hydraulic diameters less than 0.1 mm. They also report no effect of hydraulic diameter on the Nusselt number (⁠$Nu$⁠) in their study. Rosa et al. [7] investigate scaling effects on single-phase flow in microchannels and conclude that macroscale theory and correlations are valid at the microscale if measurement uncertainty and scaling effects are carefully considered. These scaling effects include entrance effects, viscous heating, conjugate heat transfer, electric double-layer effects, surface roughness, and temperature-dependent properties [3].

As described in the aforementioned studies, the uncertainty associated with experimental measurements is a crucial factor that needs to be considered in microchannel research, since even small errors in measurements might lead to significant deviations in results, making it difficult to draw meaningful conclusions. This highlights the importance of utilizing numerical simulations as a tool to better understand fluid flow and heat transfer in microchannels [2]. Numerical simulations allow for the configuration of different operational conditions and provide insight into the behavior of fluids at microscales. By incorporating the relevant physics and geometry of the microchannel, numerical simulations can help researchers interpret their experimental observations and provide a deeper understanding of the underlying mechanisms that govern fluid flow and heat transfer in microchannels. For example, Lee et al. [11] found no significant difference between the predictions made by the fully conjugated model and the thin-wall model for microchannels with a rectangular cross-sectional area. On the other hand, Sahar et al. [12] found that the uniform heat flux assumption was not valid, and the deviation between the experimental results and the 3D fully conjugated model was attributed to the nonuniform distribution of heat flux along the channel. Dharaiya and Kandlikar [13] report that the effects on $Nu$ in the entry region of the developing flow are insignificant. Additionally, Gunnasegaran et al. [14] report that the rectangular channel with the smallest hydraulic diameter had the highest heat transfer coefficient, while Sahar et al. [15] showed the thermal performance index should be taken into consideration in the analysis of the thermohydraulic performance of microchannels heat exchangers. Pan et al. [16] report the existence of the optimal aspect ratio in which the heat transfer performance of the microchannel heat sink reaches its peak. The optimal aspect ratio is found to be different for various working fluids and solid materials. However, it is worth noting that accurately modeling the behavior of fluid flow and heat transfer in microchannels can be challenging. This is mostly linked to the very small length scales and the fact that the flow can be influenced by various factors such as entrance effects, viscous heating, conjugate heat transfer, heat losses, surface roughness, inlet/outlet restrictions and pressure sensors placed outside the microchannel [17,18]. Additionally, concerting experimental and numerical studies with the same conditions inside a microchannel might be difficult due to the uncertainty in the geometry of the channel.

In the present study, we investigate the pressure drop, temperature profile and heat transfer for De-ionized (DI) water flow in a rectangular microchannel with a hydraulic diameter of 850 μm in the Reynolds number range up to $Re≈$ 1100 by means of experiments and fully resolved 3D simulations of the microchannel system including inlet and outlet manifolds and conjugate heat transfer. Our main objective is to contribute to the understanding of the behavior of fluid flow and heat transfer in microchannels through a concerted experimental-numerical study, by identifying the reasons for the observed discrepancies and the ways to resolve them.

In the experimental part of this study, isothermal pressure drop at 22.7 °C and temperature profile at a heating power of 10 W in a 492 μm high, 1.5 mm wide, and 65 mm long rectangular stainless steel (SS) microchannel with a hydraulic diameter of 848.9 μm is measured. The mass flowrate ranges from 1 to 68 gmin⁻¹ for pressure drop measurements in isothermal cases and from 1 to 49 gmin⁻¹ for temperature profile measurements in diabatic cases. This corresponds to Reynolds numbers between 16 and 1144. De-ionized (DI) water is used as the working fluid.

Figure 1 shows the considered microchannel heat sink with integrated resistance temperature detectors (RTDs) to measure the temperature profile along the microchannel and cartridge heaters [19] to heat up the microchannel. The same experimental apparatus has been used in a set of studies reported in Refs. [20–23]. The heat sink housing is specifically designed for mechanical fixation of the microchannel, RTDs, and cartridge heaters, for hermetic sealing of the microchannel, and for thermal insulation from the environment. A fluid reservoir storing DI water at room temperature is connected to a micro-annular gear pump mzr-4622 from HNP Mikrosysteme GmbH [24]. This pump supplies the DI water with low pulsation and at well-defined flow rates to the heat sink with the SS microchannel. The low pulsation of the micropump intends to reduce pressure oscillations during pressure drop measurements. A 10 μm particle filter is placed in front of the micropump to prevent contaminants from entering the microchannel heat sink. The heat sink outlet is connected to a wastewater reservoir with a high-precision balance to accurately determine the mass flowrate. The wastewater reservoir is primarily important during temperature profile measurements, since in closed-loop configurations, where the heated fluid at the outlet of the heat sink flows back into the fluid reservoir, the fluid temperature at the heat sink inlet would increase over time.

The rectangular microchannel (492 ± 5 μm × 1.5 mm × 65 mm) was fabricated in a 1.4404/316 L SS block (5 mm × 8 mm × 68.6 mm) by machining and mounted in the heat sink housing as shown in Figs. 1(a) and 1(b). The DI water is supplied and discharged via two holes at the ends of the machined microchannel, each leading to an inlet and outlet on the heat sink housing. The heating power required for the temperature profile and heat transfer analysis is applied via eleven resistive cartridge heaters located directly below the microchannel. The cartridge heaters are mounted in 4 mm deep holes along the bottom of the microchannel plate and fixed in the bottom plate of the heat sink housing with thermally insulating PEEK screws. They carry integrated type K thermocouples (T/C) in their top tips measuring the temperatures at the points of contact between the T/C and the microchannel block. The RTDs shown in Fig. 1(c) are fabricated in a clean room process on a Pyrex glass lid (2 mm × 25 mm × 75 mm) [20], which is pressed onto the microchannel plate from above by means of the heat sink housing lid. An O-ring in a groove surrounding the microchannel provides a hermetic seal between the microchannel and the glass lid. The heat sink housing surrounding the microchannel is fabricated by laser stereolithography from a high-temperature stable, inert resin with a low thermal conductivity of 0.62 Wm⁻¹K⁻¹ (at 23 °C) [25]. On both sides of the resin housing, parallel to the longitudinal axis of the microchannel, there are 17 holes each, in which a total of 32 low-resistance spring contacts are anchored [26]. These spring contacts establish electrical contact with up to 17 RTDs when the glass lid is pressed on. This allows a high spatial resolution of the inside glass lid temperature along the microchannel, which could be compared with the temperature measurements of the T/C at the tips of the cartridge heaters and the numerical simulation data (see Sec. 3). The temperature at the outer surface of the heat sink was measured using the benchtop thermal camera FLIR A65SC from Teledyne FLIR with a spatial resolution of 1.31 mrad and an accuracy ± 5% of reading.

Figure 2 shows the details of utilized microchannel. Due to the corner radius of 0.2 mm of the milling tool used for channel milling, the geometry of microchannel is not exactly rectangular. This is clearly visible in Fig. 2(c) in the rounded corners between the channel wall and the channel bottom surface. Based on the measurements using atomic force microscopy, the surface roughness of the channel bottom is 55 ± 23 nm. The fit between the glass lid and the heat sink cavity for the glass lid is estimated to be 100 ± 50 μm: the depth d of the rectangular cavity for the glass lid in the 3D printed reactor case is 2.1 mm (see Fig. 1) with an uncertainty of ±50 μm caused by the 3D printing process with the glass lid being 2.0 mm thick. As shown in Fig. 2(a), this clearance ensured that the RTDs did not directly rest on the top of the channel block. This is because direct contact between glass lid and the SS channel block could cause an electrical short circuit of the RTDs and mechanically damage these structures as well as the glass lid, leading to a distortion of the measured resistance and leakage. It has to be noted that the introduced clearance height is unknown and may vary among different experimental campaigns. We estimate that this additional clearance increases the height of the microchannel from 492 ± 5 μm to 592 ± 50 μm. More information on the properties of RTDs can be found in Ref. [20].

The isothermal pressure drop from the inlet to the outlet of the microchannel is determined using a differential pressure transducer with compensated line pressure and temperature dependency [27]. The compensation for line pressure is important as the line pressure changes over time with the level of the DI water and wastewater reservoir altering the pressure drop measurement. The differential pressure transducer is connected to fluidic ports located at the bottom of the heat sink housing as shown in Figs. 1(a) and 5. The heat sink inlet and outlet meet through-holes inside the heat sink which act as fluidic tees and lead the DI water to the pressure transducer ports and the microchannel inlet/outlet. The indicated mass flows are determined with the precision balance at the wastewater reservoir. For this purpose, the micro-annular gear pump was run for 60 s and the weight difference was determined. The mass flowrate corresponds to the weight difference divided by the time. The effect of evaporation from the free liquid surface of DI water in the wastewater reservoir was investigated by measuring the change in weight during the test period of 60 s with the pump turned off. There was no change in the measured weight, so the influence of evaporation on the measured weight difference during the test period is confirmed to be negligible. This is also consistent with the result reported in a similar study [6].

The pressure losses are computed based on the model by Lee and Garimella [28]. The detailed description for estimation of pressure losses and material properties can be found in Appendix  A.

The temperature profile along the microchannel was determined in steady-state conditions at a total heating power of 10 W in three different ways. Steady-state conditions are assumed after the measured temperature profile remains unchanged for a certain prescribed flowrate. For most flow rates, this occurred no later than two hours after the start of the experimental measurement. The first temperature measurement method utilized the T/C at the tip of the cartridge heaters which are recessed at the bottom of the channel. The second temperature measurement is conducted with the RTDs along the glass lid on the top wall of the microchannel and the third measurement method uses a thermal imaging camera FLIR A65SC above the glass lid. In addition, the fluid temperature is measured directly at the heat sink outlet (approx. 18 mm downstream of the microchannel outlet) with a type T thermocouple.

The total heating power of 10 W is supplied by the cartridge heaters at the bottom of the microchannel for all measurements with cartridge heaters connected in series. Due to the production-related slight variations in the internal resistance of the cartridges (standard deviation of ±4.4%), the exact heat output of each cartridge is calculated by its weighted internal resistance relative to the applied 10 W. The cartridges are positioned according to increasing heating output (860–984 mW) from the inlet to the outlet of the system. This helps to reduce the fluctuations in the temperature profile along the channel linked to variations of cartridge heating powers. The exact heat outputs of the cartridge heaters are used during the numerical simulation as data reduction method to reconstruct the average microchannel wall temperatures given in Fig. 18. The first and last heating cartridges were positioned 9.25 mm from the microchannel inlet and outlet, respectively. Each heating cartridge has a diameter of 3.1 mm and a distance of 4.5 mm from its center to the center of an adjacent cartridge.

where $TRT$ is the room temperature, $RRT$ is the resistance at room temperature, $ΔR=R–RRT$ with the measured resistance $R$ and $α=$ 2.98$×$10⁻³ °C⁻¹ [20].

The thermal imaging camera has been calibrated through following procedure. Firstly, the temperature of 26 °C reflected from the laboratory environment has been measured using a planar crumpled aluminum foil placed directly over the glass cover. Next, the emissivity $ε2$ of the glass lid has been determined as shown in Fig. 4 through local application of a black coating spray with an emissivity $ε1$ of 0.98 on a glass lid. First, the emissivity of the thermal imaging camera was set to 0.98 the emissivity of the black coating. Then, the glass lid was heated to approximately 100 °C on a hotplate and the temperature $T1$ of the black coating was measured with the thermal imaging camera, so a subsequent adjustment of the camera's emissivity until the full temperature correspondence $T2=T1$ could be performed with $T2$ being the temperature of a glass lid spot without black coating immediately adjacent to the coating. In this process, the glass lid emissivity $ε2$ was determined to be 0.89. This value depends on the exact composition of the borosilicate glass and is comparable to emissivities found in the literature, for example, 0.82 for Pyrex [29].

where $uy$ is the estimated uncertainty for $y$⁠, $σi$ represents the measurement values used to calculate $y$ and $uσi$ stands for the uncertainties of measurement values. Since the mathematical formulations for the measurement uncertainties derived from Eq. (3) are quite extensive, they are listed in Appendix  B, which is intended to serve to the reader as a reference for identical or similar experimental trials.

In this section, we provide a description of the models employed as a numerical counterpart for the investigation of the heat transfer and fluid dynamics within the microchannel heat sink.

To simplify the computational process and reduce computational costs, the geometry considered in the numerical simulation consists solely of the core of the heat sink, excluding other components present in the experimental setup described earlier. For isothermal simulations, the model considers only the liquid region in the heat sink, i.e., microchannel, including cylindrical inlet and outlet manifolds and a semirectangular microchannel as shown in Fig. 5(a). The fluid enters and leaves the channel vertically (y-axis), in a direction normal to the channel axis (x-axis). The inlet/outlet manifolds and channel dimensions correspond to the ones described in the experiment section. For diabatic simulation, the model considers the whole metal SS block region with heaters placed in the cylindrical blind holes at the bottom of the metal block as shown in Fig. 5(b). Here, the fluid region includes the microchannel and inlet/outlet plenums without pressure ports. The glass lid is modeled through the application of a special boundary condition placed at the top wall of the microchannel as described in the following section.

It has to be noted that the true shape of the cross-sectional area within the microchannel is unknown. As depicted in Fig. 6, we expect a deviation from the nominal shape to be present. While the nominal cross section of the microchannel is typically assumed to be rectangular, in reality, the bottom of the microchannel has curved corners as a result of the milling process, and the height of the channel is slightly greater due to the clearance gap between the glass lid and the heat sink block. Therefore, different levels of abstraction for the cross-sectional geometries were examined to select the optimal combination for the numerical setup by comparing the numerical result with the experiment (see Sec. 4.1.1). The final configuration incorporates both the curved corners and the additional gap size $Hfit$ (corresponding to the corner radius with gap cross section in Fig. 6) to better represent the real-world conditions and improve the accuracy of our numerical simulations. By incorporating these considerations, we aim to accurately capture the heat transfer and fluid dynamics phenomena within the microchannel heat sink while managing computational complexity and ensuring computational efficiency.

Both isothermal and diabatic simulations were carried out using the open-source computational fluid dynamics (CFD) package openFOAM-7 [31]. The geometric models for the simulations were generated using openFOAM's mesh generators, namely, blockMesh and snappyHexMesh. The number of cells for isothermal and diabatic simulations is around 1.6 × 10⁶ and 2.2 × 10⁶ cells, respectively, resulting in a 12 × 24 cell configuration for $H$ × $Wmc$ for the microchannel. To ensure the accuracy of the simulation results, a grid independence study is conducted. This configuration guarantees that the dimensionless wall distance $y+$ for the first cell of the microchannel never exceeds the value of 2.5. Figure 7 depicts the grid used in the diabatic simulations for the solid and the fluid region.

where the coefficients $aθ,i$ are summarized for the quantity $θ$ in the Table 2.

The flowRateInletVelocity boundary condition is applied to impose a constant velocity value (⁠$|v|=M/ρAmc$⁠) at inlet that matches the specified mass flowrate $M$ according to the experiment. At the walls, a noSlip condition (⁠$v=(0,0,0)$⁠) is employed. The outlet uses the pressureInletOutletVelocity condition, which applies a zero-gradient condition for outflow (⁠$∂v/∂n=0$⁠) or assigns a velocity based on the flux in the patch-normal direction in the case of inflow, the pressure is specified with totalPressure (⁠$p=p0−|v|2/2$⁠) with $p0=$ 0. For the inlet and walls, a zero-gradient condition is applied for the pressure $(∂p/∂n=$ 0). Here, $n$ represents the normal direction of the patch.

The thickness, $li$⁠, and thermal conductivity, $κi$⁠, of surrounding materials are considered. The housing has a thermal conductivity of 0.621 Wm⁻¹K⁻¹ (at 23 °C) [25] and a thickness of 15.6 mm except at the top which is 3.5 mm. The same boundary condition is applied at the microchannel top to account for the presence of the glass lid which has a thickness of 2 mm and a thermal conductivity of 1.2 Wm⁻¹K⁻¹. An assumed ambient air heat transfer coefficient of 25 Wm⁻²K⁻¹ is used. The impact of this assumption is examined in Sec. 4.2.2.

In this section, the numerical and experimental results related to the effects of mass flow variation on the pressure drop under isothermal conditions and on the temperature profile and local heat transfer along the microchannel under diabatic conditions are compared and discussed. The focus of the discussion is on identifying and resolving possible causes of inadequacies in experimental measurements and numerical simulations.

Figure 8 shows the pressure drop from numerical simulations evaluated for the different cross section abstraction levels (Fig. 6) of the experimental microchannel. By comparing the microchannel pressure drop curve for the nominal cross section with the pressure drop curve where the corner radius of 0.2 mm at the microchannel bottom is considered, it can be seen that the corner radius slightly increases the numerically simulated microchannel pressure drop. If the microchannel clearance of 100 μm between the microchannel and the glass lid is added (corner radius with gap) the simulated pressure drop reduces significantly. However, considering the small fluid volume above the microchannel walls caused by the 100 μm clearance between the microchannel and the glass lid (realistic) has no significant effect on the microchannel pressure drop and can be neglected. The consideration of geometric details of the microchannel such as the corner radius and the gap is thus essential to obtain a match between numerically simulated and experimentally measured pressure drop. However, the experimentally measured pressure drop $Δptot$ still deviates from the numerical pressure drop of the corner radius with the gap domain and the realistic domain with increasing mass flowrate, although geometric details are considered in the simulation. By subtracting the pressure losses at the 90 deg bends of the microchannel inlet and outlet (Fig. 3) from the experimentally measured pressure drop $Δptot$ using Eq. (1), this deviation can be eliminated. This is discussed in detail in the next section.

In Fig. 9(a) the experimentally measured pressure drop $Δptot$ and the pressure drop $Δpmc$ experienced in the microchannel and calculated according to Eq. (1) are compared with the pressure drop extracted from the corresponding numerical simulations. All numerical simulations are plotted for the corner radius with gap domain shown in Fig. 6.

When comparing the profile of the measured pressure drop $Δptot$ with the simulated pressure drop $Δpmc$⁠, an increasing deviation can be observed with increasing flowrate. This can be attributed to the growing pressure losses at the 90 deg-bends of the inlet and outlet manifolds, which are also included in the measurement due to the position of the pressure transducer connection ports on the bottom side of the heat sink (Fig. 1). A look at Fig. 9(a) shows, however, that these can be reliably estimated using the data reduction method by Lee and Garimella [28] to calculate $Δpin$ and $Δpout$ (see Appendix  A). It is important to note that such pressure losses lead to significant discrepancies between numerical and experimental results and should therefore be carefully considered. The influence of the fluidic tees at the inlet and outlet of the heat sink, which divide the DI water flow to the microchannel and to the pressure ports, the fluidic path to the pressure ports (Fig. 5), and the piping from the pressure ports to the pressure transducer were negligible in both experiment and simulation and did not lead to any significant deviation.

The magnitude of the relative deviation of the numerical simulation from the experimental pressure drop $Δpmc$ is shown in Fig. 9(b). At the lowest flowrate of 1 gmin⁻¹, where the measured differential pressure is only 38 Pa, the measurement approaches the uncertainty of ±35 Pa of the pressure transducer used (Table 1). This results in an uncertainty of the relative deviation magnitude that is at least 92% (if only the uncertainty of the pressure transducer is considered). However, the relative influence of this measurement uncertainty decreases rapidly with increasing mass flowrate. For example, at 68 gmin⁻¹ it constitutes only 0.01%. At the same time, we observe an increase in the absolute measurement uncertainty for $Δpmc$⁠, since the largest influencing factors are now the cross-sectional microchannel area $Amc$⁠, the microchannel height $Hmc$⁠, the microchannel width $Wmc$ and it scales with the squared mass flux $G$ (see Eqs. (B7)–(B9)). Nevertheless, at higher flow rates, we observe an overall reduction in the deviation between simulation results and the experimental data. This can be attributed to the fact that the relative uncertainty in the experimental data becomes smaller compared to the absolute value of pressure. As a result, the experimental data becomes more reliable and less prone to errors. However, it is important to note that even at high flow rates, the uncertainty related to the microchannel cross-area shape, especially due to the uncertainty of the microchannel clearance $uHfit$ (Table 1), might cause discrepancies between the simulation and experimental results. Compared to that, the influence of all other uncertainties in Eqs. (B7)–(B9), for example, the uncertainty of the water density $uρ$ and the uncertainty of the mass flowrate $uM$ are negligibly small. The equations for all measurement uncertainties are given in Appendix  B.

Fig. 10(a) illustrates the contribution of the inlet and outlet pressure loss $(pin+pout$⁠, referenced as inlet/outlet) estimated from the correlations and the microchannel pressure drop (⁠$pmc$⁠, referenced as microchannel) to the total measured pressure $Δptot$ with increasing mass flowrate. As described above, any change in the flow cross section, such as contraction at the microchannel inlet or expansion at the microchannel outlet, creates a local pressure loss that is also measured in the experiment and has to be considered when compared to the simulation data. In the heat sink system at all considered flow rates, the main pressure drop (minimum 85% of $Δptot$⁠) occurs in the microchannel.

The flow in the microchannel can be divided into a flow-developing zone, where the velocity profile keeps evolving, and a fully-developed zone where the velocity profile becomes invariant to the streamwise location. The flow-developing zone is often referred to as the hydrodynamic entrance zone and is characterized by its length. The channel entrance zone contributes to the additional pressure loss as the growth in the boundary layer accelerates the flow inside the inviscid central region. Unlike macroscale channels, in microscale devices, the entrance losses can become significant since the developing zone might even occupy the entire length of the channel. The hydrodynamic entrance length has traditionally been defined as the distance from the channel inlet to the location where the velocity profile reaches 99% of the fully developed velocity profile [33]. In the present study, this was approximated as the location where the centerline velocity $Uc$ of a developing flow reaches 99% of the centerline velocity expected in the fully developed profile $Ufd$⁠. In our present study, Fig. 10(b) shows the fraction of entrance region length relative to the overall length of the microchannel. It covers the entire microchannel length for the two highest mass flow rates considered. The length of the fully-developed zone is a little shorter than what you see in the picture. This is because we didn't exclude the part where the flow is about to leave the channel.

In the present study, the hydraulic diameter is $Dh=$ 846 μm. It has to be noted that the influence of channel aspect ratio on the dimensionless entrance length is negligible for Reynolds numbers $Re>$ 50 [33]. The coefficients $c1$⁠, $c2$⁠, and $c3$ are provided in Table 3 for various aspect ratios based on the findings of Galvis et al. [33]. In numerical studies, the inlet velocity profile is often assumed to be flat, and the choice of inlet conditions significantly affects the entrance length [38]. In practical devices, micro/mini channels are typically preceded by a fluidic element such as a tank or a plenum. Lobo and Chatterjee [37] recently investigated the entrance length of a microchannel connected to a plenum with varying aspect ratios. The authors presented a correlation that accounts for the effect of aspect ratio on the entrance length; however, their correlation is not applicable for low Reynolds numbers (especially for $Re<$ 20).

The hydrodynamic development of flows has often been related to the growth of the boundary layer along the channel's walls. However, in the present study, the flow is streaming in/out of the microchannel with a vertical pipe or manifold on each end, which causes a sudden change in the cross section following a 90 deg-bend (Fig. 3) from the pipe to the microchannel. This results in flow redirection and separation and generates a secondary flow along the microchannel, especially for high mass flow rates (⁠$M>$ 20 gmin⁻¹). Figure 11 shows the velocity profile in the entrance region of the microchannel, where both the streamwise (⁠$Ux$⁠) and the in-plane (⁠$Us=Uy2+Uz2$⁠) velocities are normalized with the average velocity, which is computed based on the microchannel mass flux (⁠$Uave=Gmc/ρ$⁠). At low mass flow rates, the velocity profile is rather uniform, while at higher mass flow rates, a stronger secondary flow emerges leading to flow asymmetry for the streamwise velocity component. The main streamwise velocity is, however, still much more influential than the secondary flow. The inlet condition of a microchannel has a significant impact on the length of the entrance region, especially at low mass flow rates. The flow undergoes a bending process as it enters the microchannel, and this bending contributes significantly to the developing zone in the entrance region. The size of the separated flow area varies depending on the mass flowrate, and it becomes more pronounced at higher flow rates. In Fig. 12, we present a comparison between the measured entrance lengths obtained from our numerical simulations and those reported in previous studies. To capture the relationship between mass flowrate and entrance length, we fitted a curve based on Eq. (8) using the selected mass flow rates, and the corresponding fitting parameters are provided in Table 3. Our findings reveal that the entrance region in our experimental setup is slightly longer than what has been reported in the existing literature. However, as the mass flowrate increases (⁠$Re>$ 50), the influence of the entrance condition on the length of the entrance region becomes insignificant compared to the impact of the boundary layer. Consequently, for high mass flow rates, the estimated entrance length in our study aligns well with the values documented in the literature.

where for square ducts (⁠$AR=$ 1) $Po$ is 56.9. In the present study, the nominal aspect ratio of the rectangular cross section is 0.393, which gives $Po=$ 65.76 (respectively $fd=$ 65.76$/Re$⁠) according to Eq. (11). In Fig. 13(a) the Poiseuille number $Po$ and in Fig. 13(b) the Darcy friction factor $fd$ is computed for both the numerical simulation and the experiment. They agree well with each other when the measurement uncertainty is considered. For both, $Po$ increases and $fd$ decreases with increasing flowrate, thus appearing to deviate from Stokes flow macroscale theory.

The experimental deviation from $Po=fdRe=const$ for various mass flow rates is a phenomenon that is also reported in previous studies [8,41] for microchannels. It is often linked to viscous heating and fluid polarity [6]. In a study by Gian Luca Morini [18], it is reported, that viscous dissipation produces a non-negligible effect for liquid flows in microchannels with a hydraulic diameter smaller than 100 μm and leads to a decreasing Poiseuille number with increasing Reynolds number. However, in the isothermal experimental pressure measurements of this study in a microchannel with a hydraulic diameter of 848.9 μm, no viscous heating of the DI water is observed. The fluid temperature remains constant from the inlet to the outlet of the heat sink for all considered flow rates. Our numerical study, however, shows that in spite of the absence of viscous dissipation effects the Poiseuille number can increase with increasing mass flowrate (Reynolds number), which is related to the length of the developing flow region. If the entrance length (⁠$Le$⁠) is excluded from the total microchannel length $Lmc$ and only the fully-developed region in the microchannel is considered, then $Po$ remains constant at about 63.3. This is clearly visible in Fig. 13 from the numerical simulations that only consider the developed region. Hence, the consideration of flow development at the microchannel entrance might be crucial for an explanation of result deviations from the Stokes flow theory. The minor deviation from the theoretical value of $Po=$ 65.76 calculated with Eq. (11) could be due either to the measurement uncertainty $uPo$ (see Eq. (B12)) mainly dominated by the uncertainty of the fit between heat sink and microchannel $uHfit$ (Table 1) or the shape of the microchannel cross section, which is not perfectly rectangular due to the corner radii at the bottom edges. The uncertainty of the friction factor $uf$ as well as the uncertainty of the Poiseuille number $uPo$ (calculated with Eqs. (B10)–(B12)) at the lowest flowrate of 1 gmin⁻¹ are estimated to be very high (Fig. 13). This is linked to the fact that they are in the first place dominated by the uncertainty of the pressure transducer $uptot$ (Table 1). However, at a flowrate of 5 gmin⁻¹ or higher, the influence of $uΔptot$ on the uncertainties $uf$ and $uPo$ becomes negligibly small so they are rather dominated by $uHfit$ and slightly increase with the squared mass flux $G$ as explained in Sec. 4.1.2.

Fig. 14(a) shows the temperature profile of the numerically simulated glass lid temperature profile together with the experimentally measured RTD temperatures from the inlet center (0 mm) to the outlet center (63.5 mm) of the microchannel. Fig. 14(b) shows the experimentally determined T/C temperatures at the tips of the cartridge heaters and the numerically simulated temperatures along the heater domain matching the positions of the T/C tips. The plotted standard deviation of the measured RTD and T/C temperatures is in the range of 3$u$⁠, so the measured values are expected to be within the plotted measurement uncertainty range with a probability of 99.7%. It can be observed that for mass flow rates lower than 49 gmin⁻¹ the RTD and T/C temperatures mostly agree with the numerically simulated temperature along the glass lid and heater domain within an uncertainty of 3$u$⁠. However, in the range above 48 mm downstream of the microchannel inlet, the measured RTD temperature seems to drop, whereas the simulated glass lid temperature continues to increase with a slightly reduced slope. This is particularly evident when comparing the RTD-measured temperature values with the simulated temperatures along the glass lid domain at 2 gmin⁻¹. The experimentally measured T/C temperatures increased already from the beginning at the microchannel inlet with a significantly lower slope than the simulated heater domain temperatures for a flowrate of 2 gmin⁻¹. It is hypothesized that the axial heat losses to the heat sink housing and the environment increase close to the microchannel outlet, which is underestimated in the numerical simulation due to the uncertainty of the introduced air heat transfer coefficient (more in Sec. 4.2.2). A similar heat loss is observed in a related heat sink design with a single stainless steel microchannel by Talebi et al. [42]. In the work of Talebi et al. a significant drop in channel temperature of up to 10 °C occurs from the center to the outlet of the microchannel, which is linked to an increased heat loss at the outlet of the microchannel. To support this argument and to rule out a possible defect in the RTD structure at the microchannel outlet, in Fig. 15 the RTD measurement at the outlet of the microchannel was compared with the T/C type T fluid temperature measurement at the outlet of the heat sink housing 18 mm downstream of the microchannel outlet. The relative deviation between the two temperatures was on average only 0.9%. It can therefore be assumed that the RTD temperature measured at the microchannel outlet approximates the fluid temperature at the heat sink outlet and correctly captures the data.

In the present study, an air heat transfer coefficient (HTC) of $h=$ 25 Wm⁻²K⁻¹ is initially assumed as suggested in Ref. [43] for the surrounding air in contact with the heat sink. However, it is important to note that the HTC for air is not fixed nor uniform and is primarily influenced by the surface temperature. In general, the HTC increases with a larger temperature difference between the air and the surface. This variability in HTC presents a significant source of uncertainty in heat transfer simulations.

As depicted in Figs. 14(a) and 14(b), assuming a fixed HTC can lead to an underestimation of heat loss, especially in scenarios with low mass flow rates inside the heat sink where the surface temperature of the heat sink is higher. To further emphasize the impact of HTC variation, Figs. 16(b) and Fig. 16(a) illustrate the numerical temperature profile of the heater and glass lid domain for three different mass flow rates. The figure confirms that at lower mass flow rates the HTC has a significant impact on the result of the numerically simulated heater and glass lid domain temperatures. Conversely, at higher mass flow rates, where the surface temperature is lower, the uncertainty in the HTC has less influence on the overall result. Generally, the HTC uncertainty has a greater influence on the simulated heater domain temperature than on the glass lid domain temperature. Above a mass flowrate of 10 gmin⁻¹, the HTC of the simulated glass lid domain has a negligible influence on the temperature curve. This is confirmed by the surface temperature line profiles measured with the infrared camera along the external glass lid surface from the microchannel inlet to the microchannel outlet (Fig. 17(c)) and along the external heat sink surface from the first to the last heater cartridge (Fig. 17(d)). Figs. 17(a) and 17(b) show the corresponding infrared camera images of the heat sink top and heat sink front where the respective positions of the temperature line profiles are marked. A closer look at Figs. 17(c) and 17(d) confirms that at a low flowrate of 2 gmin⁻¹, the plotted external glass lid and heat sink surface temperatures are the highest and are therefore most likely to be affected by uncertainties of the estimated HTC. The measured surface temperatures of the glass lid and the heat sink front side drop considerably between 2 and 5 gmin⁻¹. From 5 gmin⁻¹, the surface temperature of the glass lid remains constant up to a mass flowrate of 50 gmin⁻¹, whereas the surface temperature of the heat sink front continues to decrease, but at a rate that decreases with the mass flowrate. Based on a set of observations, it has to be concluded that proper estimation of HTC is crucial for obtaining accurate and reliable results in heat transfer simulations.

To check the effect on variation in the provided heat source at particular cartridges, a nonuniform heat distribution is used based on the exactly measured heat output of each heater cartridge given in Sec. 2. Figure 18 shows the numerically simulated average microchannel wall temperature from inlet to outlet with an additional uniform heat source temperature profile at 2 gmin⁻¹. Due to the thermal conductivity of the SS microchannel (15 Wm⁻¹K⁻¹), the difference in the temperature profile for a uniform and a nonuniform heat source was only marginal across all flow rates. This finding is interesting as it suggests that the standard deviation of ±4.4% in heat power of the cartridge heaters does not significantly impact the temperature distribution in the microchannel and therefore does not leading to a disagreement between experiment and numerical simulation.

where $cth$ is the empirical coefficient and is suggested to be chosen in the range from 0.028 to 0.116 for different microchannel aspect ratios [44]. Based on correlations from Ref. [44], for the present configuration at $AR =$ 0.393, the coefficient is estimated to be $cth=$ 0.0719. This corresponds to thermal development length $Lth=$ 13.8, 35.1, 70.1, 140.6, 350.1 mm for the mass flow rates of $M=$ 2, 5, 10, 20, 50 gmin⁻¹, respectively. However, due to the fact, that we never approach the thermally developed state in the considered system, we cannot confirm these estimations through numerical simulations or experiments.

and the fully developed Nusselt number $Nufd=$ 4.5 for constant heat flux in a rectangular channel at $AR=$ 0.39 is assumed [45]. The observed deviation from the correlation at the outlet of the microchannel is attributed to the nonuniform heat flux conditions in our setup.

where $Pe*=RePrDh/Lmc$⁠. The deviation observed between the presented data and previous studies can be attributed to differences in experimental/numerical setup. In our present investigation, it is noteworthy that the heat flux distribution is nonuniform, which deviates from the assumptions of perfect uniform heat flux made in prior numerical studies. Furthermore, our system exhibits a nonrectangular cross section (small rounded corners at the bottom of microchannel wall), in contrast to idealized geometries assumed in previous research. Additionally, we consider the minor heat losses through the glass lid, a factor that was rarely considered or considered as adiabatic boundary condition in previous studies.

The pressure drop, temperature profile, and heat transfer of a heat sink with a single stainless steel microchannel with a hydraulic diameter of 848.9 μm has been studied experimentally and numerically for DI water flow at various flow rates. The inconsistencies between experiments and numerical results have been examined in order to understand their roots and develop a set of rules for a successful handshake between simulation and experiment.

The numerically simulated pressure drop is found to be well within the range of experimental measurement uncertainties after the nominal geometry has been adjusted to consider the realistic details of the setup—the corner radius of 0.2 mm at the edges of the microchannel bottom and the 100 μm larger height representing the unavoidable gap between the glass lid and the metal block. Furthermore, the total measured pressure loss has been corrected considering the additional pressure losses at the inlet and outlet manifold using correlation-based loss factors. A closer look at the laminar flow development from the inlet to the outlet of the microchannel demonstrates that already from a flowrate of 20 gmin⁻¹ half of the channel length is needed until the laminar flow reaches a fully developed state. The development of the flow is mainly governed by the presence of the channel bends at the inlet of the microchannel system, which induces three-dimensional pattern in the flow field and promotes a formation of a growing boundary layer downstream of the bend. The presence of developing regions leads to an increasing Poiseuille number at higher mass flow rates and hence introduces a disagreement compared to the laminar flow theory, where $Po$ is supposed to remain constant. Considering only the fully developed microchannel region, the Poiseuille number is demonstrated to be constant at 65.8 and shows no deviation from the laminar flow theory. Overall, it can be concluded, that the experimental measurement uncertainties of the pressure drop and the Poiseuille number are largely determined by the appropriate prescription of the microchannel cross section and its geometrical details.

The agreement between simulated and measured temperature profiles at the heater and glass lid improves as the flowrate increases. The largest deviation between numerically simulated and experimentally measured temperatures is found in the area of the microchannel outlet. Presumably, this is linked to the estimation of the prescribed heat transfer coefficient of the ambient air in the simulation, as the coefficient increases with increasing temperature gradients. Since the heat sink surface temperature is highest at low flow rates and increases from the microchannel inlet to the outlet, the largest temperature gradient occurs at the lowest flowrate of 2 gmin⁻¹ in the area of the microchannel outlet. The chosen heat transfer coefficient of 25 Wm⁻²K⁻¹ thus underestimates the heat transfer between ambient air and heat sink of the experiment at low flow rates in the outlet area of the microchannel, while it reliably estimates the heat transfer in the remaining area of the microchannel at higher flow rates.

The determined local Nusselt number compares well to the correlations for the laminar flow of DI water, except for the values near the channel outlet, which is again related to the choice of heat transfer coefficient at the boundaries of the simulation domain. This underscores the importance of prescribing accurate nonhomogeneous thermal boundary conditions in numerical simulations as a crucial factor for faithfully capturing all the experimental aspects of the process.

Concluding, we identify the geometrical description of the microchannel geometry including plenum design of the inlet/outlet of the system together with the knowledge about the flow state (developing or developed) as the most significant factors for the experimental-numerical handshake in terms of pressure drop estimation. Additionally, the proper choice of the thermal boundary condition at the outer boundary remains the most important factor for a proper prediction of the heat transfer behavior in microchannel systems.

• Deutsche Forschungsgemeinschaft (Grant Nos. STR 1585/2-1 and WO 883/24-2; Funder ID: 10.13039/501100001659).

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

The pressure loss due to flow contraction $Δpin$ at the microchannel inlet 90 deg bend (with $Amc<Ain$⁠) and flow expansion $Δpout$ at the microchannel outlet 90 deg bend (with $Aout>Amc$⁠) as depicted in Fig. 3, are calculated using the data reduction method by Lee and Garimella [28]

• pressure loss due to contraction at microchannel inlet

  $Δpin=[1−(AmcAin)2+Kcont]Gmc22ρ$

  (A1)

  with the loss coefficient

  $Kcont=0.0088(AR)2−0.1785(AR)+1.6027$

  (A2)
• pressure recovery due to expansion at microchannel outlet

  $Δpout=12KexpG22ρ$

  (A3)
  with the loss coefficient

  $Kexp=−2⋅1.33(AmcAout)[1−(AmcAout)]$

  (A4)

The following dimensions and properties are considered for evaluation:

• height of the microchannel $H=Hmc+Hfit$ including the gap to the glass lid: 592 ± 50 μm;

• microchannel width $Wmc$⁠: 1500 ± 1 μm;

• aspect ratio $AR=H/Wmc$⁠: 0.395 ± 0.033;

• cross-sectional area of the microchannel $Amc=(Hmc+Hfit)Wmc$⁠: 0.888 ± 0.075 mm²;

• area of the microchannel inlet and outlet $Ain=Aout=π(din/2)2=π(dout/2)2$⁠: 1.767 mm²;

• mass flux $Gmc=M/[(Hmc+Hfit)Wmc]$ of the DI water inside the microchannel in kgm⁻²s⁻¹;

• density of water at 22.7 ± 1.1 °C $ρ$⁠: 997.6 ± 0.3 kgm⁻³.

with $a1=$ 999.8531, $a2=$ 6.3269, $a3=$ 8.5238, $a4=$ 6.9432, $a5=$ 3.8212.

with the corresponding coefficients for water $A=$ 0.02939 Pas, $B=$ 507.88 K, $C=$ 149.3 K.

with $a=$ 4.2174356, $b=$ −0.0056181625, $c=$ 0.0012992528, $d=$ −0.00011535353, and $e=$ 4.14964 $× 10−6$⁠.

The following estimation for measurement uncertainties is utilized for the metrological characterization:

• microchannel cross-sectional area

  $uAmc=Amc⋅(uHH)2+(uWW)2$

  (B1)
• microchannel hydraulic diameter

  $uDh=2⋅AmcH+W⋅(uAmcAmc)2+uH2+uW2(H+W)2$

  (B2)
• microchannel inlet/outlet area

  $uAin=π2⋅din2⋅udin or uAout=π2⋅dout2⋅udout$

  (B3)
• microchannel mass flux

  $uG=G×(uMM)2+(uHmcHmc)2+(uWmcWmc)2$

  (B4)
• density of water

  $uρ=|(a2 × 10−2−2a3T × 10−3+3a4T2 × 10−5−4a5T3 × 10−7)⋅uT|$

  (B5)
• dynamic viscosity of water

  $uμ=(−A⋅B⋅eB/(T−C)1000⋅(T−C)2⋅uT)2$

  (B6)
• pressure loss due to contraction at inlet

  $uΔpin=((AmcAin)2+Kcont)·G22ρ… ⋅[4⋅(AmcAin)2⋅((uAmcAmc)2+(uAinAin)2)+(0.0176(HW)2−0.1785(HW))2⋅((uHH)2+(uWW)2)((AmcAin)2+Kcont)2… +4⋅((uMM)2+(uHH )2+(uWW)2)+(uρρ(T))2]1/2$

  (B7)
• pressure recovery due to expansion at outlet

  $uΔpout=12⋅Kexp⋅G22ρ⋅[(−2.66+5.32(AmcAout))2⋅(AmcAout)2⋅((uAmcAmc)2+(uAoutAout)2)Kexp2…+4⋅((uMM)2+(uHH )2+(uWW)2)+(uρρ(T))2]1/2$

  (B8)
• pressure drop of microchannel

  $uΔpmc=uΔptot2+uΔpin2+uΔpout2$

  (B9)
• friction factor of microchannel

  $uf=2⋅Δp⋅DhL⋅ρ⋅(AmcM)2⋅[(uΔpmcΔp)2+(uDhDh)2+(uLL)2+(uρρ)2+4⋅((uAmcAmc)2+(uMM)2)]0.5$

  (B10)
• Reynolds number of microchannel

  $uRe=Re⋅(uMM)2+(uAmcAmc)2+(uDhDh)2+(uμμ)2$

  (B11)
• Poiseuille number of microchannel

  $uPo=Po⋅(uff)2+(uReRe)2$

  (B12)

It is widely recognized that first-order numerical schemes provide stability, albeit at the cost of increased numerical diffusion. In order to obtain higher accuracy, second-order schemes such as linearUpwind for the momentum equation and limitedLinear for the energy equation have been employed for two distinct flow regimes: low-speed flow at 2 gmin⁻¹ and high-speed flow at 49 gmin⁻¹. Despite this, the application of second-order schemes did not yield a significant improvement in the temperature profiles at the locations of the T/Cs and RTDs as depicted in Fig. 20.