Heat Transfer Coefficient Determination during FC-72 Flow in a Minichannel Heat Sink Using the Trefftz Functions and ADINA Software

Abstract

This work focuses on subcooled boiling heat transfer during flow in a minichannel heat sink with three or five minichannels of 1 mm depth. The heated element for FC-72 flowing along the minichannels was a thin foil of which temperature on the outer surface was measured due to the infrared thermography. The test section was oriented vertically or horizontally. A steady state heat transfer process and a laminar, incompressible flow of the fluid in a central minichannel were assumed. The heat transfer problem was described by the energy equations with an appropriate system of boundary conditions. Several mathematical methods were applied to solve the heat transfer problem with the Robin condition to determine the local heat transfer coefficients at the fluid/heated foil interface. Besides the 1D approach as a simple analytical method, a more sophisticated 2D approach was proposed with solutions by the Trefftz functions and ADINA software. Finite element method (FEM) calculations were conducted to find the temperature field in the flowing fluid and in the heated wall. The results were illustrated by graphs of local heated foil temperature and transfer coefficients as a function of the distance from the minichannel inlet. Temperature distributions in the heater and the fluid obtained from the FEM computations carried out by ADINA software were also shown. Similar values of the heat transfer coefficient were obtained in both the FEM calculations and the 1D approach. Example boiling curves indicating nucleation hysteresis are shown and discussed.

1. Introduction

The use of boiling phenomena is one of the ways of intensifying heat transfer. Boiling heat transfer during flow in small size channels has great potential for large heat flux transfer due to the fact that it enables the meeting of conflicting requirements such as a high heat flux at minor temperature differences among the heated wall and saturated liquid, for slight dimensions of the heat transfer device. There is increasing interest on heat exchangers with mini- or microchannels due to these advantages mentioned above, resulting from their high process efficiency and compactness of technological solutions. It is worth mentioning that miniaturization combined with the increasing amount of transported heat is currently needed in most technical applications. This topic is widely recognized to be the key to the information technology industry for which heat flux in integrated chips is limited. Increasingly, applications in advanced technology require heat fluxes as high as possible. Knowledge on a more efficient cooling technology based on mini- and microchannels is needed, and this is the reason for its rapid expansion.

On the other hand, the impact of reduction in channel size on heat transfer has not been recognized yet, therefore further experimental and theoretical studies are required. Most contemporary research has focused on heat transfer with phase changes in flow and pool boiling on a miniature scale. In a review of the literature, selected experimental works and mathematical modelling regarding boiling heat transfer during flow in channels of small dimensions are summarized.

Researchers use various cooling fluids in flow boiling experiments conducted on experimental stands with minichannels of different geometries. The outcomes of a boiling heat transfer study on the flow of ethanol in a rectangular minichannel asymmetrically heated were discussed in Ref. [1]. Due to infrared thermography, a two-dimensional temperature field of the heated wall during steady-state experiments was assessed. The proposed mathematical model involved assumed a laminar flow and that the physical parameters of the test section were time-independent. The 2D temperature distributions in its main elements were calculated due to application of the Trefftz method and the hybrid Picard-Trefftz method. Boiling during fluid flow in minichannels engraved in a thin plate under normal, hyper- and microgravity has attracted considerable attention from the authors of Ref. [2]. The working fluid was HFE-7100. According to the authors, the analysis of temperature as a function of gravity evidenced that gravity influences the flow. Testing the inverse methods, the authors concluded that the sensors induced lots of disturbances and accentuated the ill-posed character. Boiling heat transfer during the flow of refrigerant R245fa in a multi-microchannel heat sink and the comparison with the data collected for refrigerant R236fa are provided in Ref. [3]. Several heat transfer trends were underlined depending on heat flux: the heat transfer coefficient values increased with vapor quality, but they were independent of heat flux and mass velocity (at low heat flux), the heat transfer coefficient values increased with heat flux but they were independent of vapor quality (at medium heat flux), the heat transfer coefficient values decreased with increasing heat flux and vapor quality but increased with mass velocity (at high heat flux). The investigation on boiling heat transfer during HFE 7000 flow along a circular minichannel vertically positioned was discussed in Ref. [4]. The results indicated that the heat transfer coefficient increased with heat flux and saturation temperature, but the mass flux did not have impact on the coefficient values. The obtained results were compared with selected heat transfer correlations from the literature. In Ref. [5], the authors developed a five-equation model of dryout. The model, based on the liquid balance in the film and in the core and the vapor balance in the flow core, was dedicated for minichannels. Verification of the model was provided using the experimental data collected for minitubes. In Ref. [6], a cylindrical heat exchanger with minijets was used in the experiments, whereas a tubular reference heat exchanger was used for comparison. During the testing of both heat exchangers at the single-phase convection, hot air was used as a heating medium, while water was a heated medium. It was noticed higher values of the heat transfer coefficient for the minijets heat exchanger was observed in comparison with the reference one. The main aim of Ref. [7] was computation of the boiling heat transfer coefficient at the working fluid-heated wall contact surface, during HFE-649 flow in an annular minigap. The mathematical model was proposed assuming the laminar fluid flow and axisymmetric steady-state heat transfer process. The Trefftz method was used for solving the inverse heat problem, using the experimental data obtained for the saturated boiling region. According to the authors, the heat transfer coefficient decreased with the distance from the minigap inlet, with the increasing of the vapor phase share in the liquid–vapor mixture.

Many of the researchers used Fluorinert FC-72 due to its physical properties and boiling point. The results from experimental studies on FC-72 flowing in a minichannel or microchannel set were discussed in [8,9,10,11,12]. In Ref. [8], the results from experiments conducted in the test section with a circular minichannel were discussed. The main aim was to determine the local heat transfer coefficients for subcooled and saturated flow boiling. It was shown that a relationship between the heat transfer coefficient and heat flux occurred, especially at the subcooled boiling region. An increase in local heat transfer coefficients with the increase in heat flux, whereas a weak dependence on the vapor quality was observed. A comparison between our own results and those calculated according to selected heat transfer correlations from the literature was added. Studies on FC-72 flow boiling in a multi-minichannel copper heat sink comprising 14 rectangular minichannels were described in Ref. [9]. A clear relationship between the average heat transfer coefficient with mass flux and heat flux was noticed. Moreover, the results gained at sub-atmospheric pressure and atmospheric pressure showed the same trend qualitatively. The studies on heat transfer coefficient and pressure drop under different conditions of water were provided. The authors concluded that both coolants can be used to achieve a maximum heat dissipation rate of 100 W/cm², whilst maintaining the maximum temperature of the heating wall below 90 °C. An experimental study on FC-72 and water flow boiling heat transfer and pressure drop in heat sinks with rectangular fins (two types) was presented in [10]. Under the same mass flow conditions for both types of heat sinks, similar values of heat transfer coefficient and pressure drop were observed. The experimental data comparison showed that the FC-72 flow boiling cooling power was up to 330 percent higher than that of liquid water. The FC-72 pressure drop was substantially greater than water. The topics of Ref. [11] was FC-72 flow boiling heat transfer and pressure drop in microchannels. Multi-ported rectangular microchannels were studied. It was found that with increasing vapor quality and mass flux, the pressure drop increased. The authors also detected that under equal mass flux and vapor content conditions, heat flux had a marginal effect on the pressure drop. Furthermore, a comparison was made between the expected pressure drop and heat transfer coefficient by the correlations selected from the literature and the values measured. The results of research on FC-72 flow boiling in asymmetrically heated rectangular minichannels were discussed in [12]. Liquid crystal thermography was applied to measure the temperature distribution on the heating foil, while the optical observation technique was used to observe flow patterns. The analysis of flow structure images enabled the measurement of the void fraction. The heat transfer model during the flow of the boiling liquid, based on the Trefftz method, was proposed.

A broad experimental data analysis performed on heat transfer and pressure drop in compact heat exchangers were delivered in [13,14]. The authors of Ref. [14] collected the data from experiments using several refrigerants (R410A, R407C, R134a and R404A) for their comparative analysis in terms of test parameters. The study concerned the condensation in minichannels. It was noticed that when the refrigerant mass flux density increases, the average heat transfer coefficient also increases. Furthermore, the per-channel heat transfer intensity was lower in multiports than in a single minichannel of the same diameter. The calculation methods of heat transfer determination during two-phase flow in micro- and minichannel heat sinks are often based on the results from experiments.

The most popular numerical methods in computational fluid dynamics (CFD) analysis are the finite element method (FEM) and finite volume method (FVM). There are many commercial software packages to carry out such computations: ANSYS CFX/Fluent [15,16], COMSOL Multiphysics [17], STAR-CMM+ [18], ABAQUS [19], ADINA [20], and many more. Among the above-mentioned programs, the most popular is undoubtedly ANSYS Fluent. Comparative studies of some aspects of CFD calculations using ANSYS Fluent and ADINA software were investigated in [21].

A short literature review concerning numerical studies on flow boiling in minichannels is presented. The authors of Ref. [22] investigated fluid flow and heat transfer in microchannels. Their studies showed that conventional correlations give accurate predictions for the laminar flow in rectangular microchannels of hydraulic diameters in the range of 244–974 μm. Experimental and theoretical research on single-phase heat transfer in microchannels were discussed in [23]. The heat transfer data for circular, triangular, rectangular, and trapezoidal microchannels with hydraulic diameters in the range of 60–2000 μm were analyzed. The authors focused on the comparison of the experimental data obtained by a number of researchers with the conventional heat transfer theory.

The comparison of the results showed a major discrepancy between the theoretical predictions and the measurements. For conventional and the proposed interrupted microchannel heat sinks, 3D numerical simulations of conjugate heat transfer were provided in [24]. In Ref. [25], an experimental and numerical analysis of the characteristics of fluid flow and heat transfer in a complex structured microchannel heat sink were presented. The thermal behavior in the corrugation microchannel heat sink was described by numerical simulation. The main subject of Ref. [26,27,28], was that there was thermohydraulic maldistribution in mini heat exchangers, and the shape of collector optimization would improve this problem. A numerical investigation on flow maldistribution in 50 parallel 1 mm × 1 mm minichannels and a minigap and 1 mm depth minigap section comprising rectangular, trapezoidal, triangular, or concave manifolds in Z-type flow configuration, were described in [26]. Calculations were performed in ANSYS Fluent. Using the finite volume method, the conservation equations of mass, momentum and energy were solved. In Ref. [27], the work of heat exchangers with different numbers of minichannels were simulated. The mathematical model, proposed in the MATLAB, was numerically investigated using the ANSYS CFX code. In Ref. [28], numerical simulations in ANSYS Fluent were discussed for various mass flow rates and a constant heat flux. Modifications to the configuration of the inlet manifold channel and plenum were carried out to improve the fluid distribution. The first attempts at numerical calculations according to the mathematical model describing heat transfer during fluid laminar flow in minichannels using ANSYS CFX program and ADINA software, and the computations based on the finite element method and the Trefftz functions were performed by the authors of this work and shown in [29,30], respectively. The outcomes of the comparative analysis were discussed. In [31], CFD modeling of laminar forced convection in a nanofluid flowing in a minichannel heat sink was proposed according to four models: single phase, VOF (volume of fluid), mixture and Eulerian. Using the finite volume method, the 3D steady-state governing partial differential equations were discretized. The continuous flow formation in minichannels during two-phase flow contributes to the occurrence of instabilities of heat transfer and pressure drop, and oscillations in the flow pattern. This issue was investigated and modelled using the windowed recurrence quantification analysis [32]. Another interesting topic is image analysis during multiphase flow in minichannels. Such analysis, including the boiling synchronization of two-phase flow structure in two parallel minichannels, each of 1 mm diameter, was discussed in [33]. In Ref. [34], for an advanced thermal-FSI (fluid–solid interaction) minichannel heat exchanger model, selected numerical modeling problems were presented. The heat transfer between the separated mediums for different mass flows was given special attention. Thermal-FSI flow equations (for the fluid domain) were solved by a CFD (computational fluid dynamics) solver, and the response of the solid body was obtained by using a CSD (computational solid dynamics) solver. Heat and fluid laminar flow in rectangular channels equipped with longitudinal vortex generators were numerically investigated in [35,36]. As a working fluid, a shear-thinning non-Newtonian liquid was utilized [35]. It was performed using three-dimensional simulations on a channel with five pairs of vortex generators and a plain channel for comparison. The results demonstrated an enhancement of heat transfer in the range of 39–188% for a non-Newtonian liquid in rectangular channels with longitudinal vortex generators regarding a Newtonian liquid flow. The authors stated that equipping rectangular channels with longitudinal vortex generators resulted in an enhancement in heat transfer performance compared to a plain channel, although larger pressure losses occurred. In the numerical simulations presented in [36], deionized water was used as the working fluid, with temperature-dependent thermo-physical properties. Five microchannel configurations with various angles of attack were tested for the longitudinal vortex generators. With the Nusselt number, a 2–25 percent increase for microchannels with longitudinal vortex generators was reported. Experimental testing of the inverse method of heat transfer coefficient determination was the subject of Ref. [37]. The main aim was to propose a method that can be used to solve the problems of nonlinear inverse heat conduction problem and to determine the local heat transfer coefficients on the unknown boundary. The proposed inverse method has been experimentally confirmed using a horizontal, cylindrical thick-walled tank as a collector. A comparison of the determined and measured temperature transients at points inside the collector wall confirmed that the transient temperature distribution is defined by the proposed inverse method. The authors reported that good agreements between temperature transients were achieved.

To date, the authors of the present work have studied heat transfer enhancement during flow boiling in one [38,39], two [40,41], or three minichannels [42] of rectangular cross-section and in an annular gap [7]. Most of the previous research focused on experiments conducted in the stationary state, although some of the latest works concerned time-dependent studies [41,43]. Based on the experimental data, analytical and analytical-numerical methods using the Trefftz functions, the heat transfer coefficient values on the heating surface contacting the fluid flowing in a minichannel were determined. Only the first attempts have been made to predict the behavior of cooling fluids flowing in minichannels due to numerical calculations using ANSYS CFX [30] and ADINA software [29].

This paper reports the results of a study on subcooled boiling heat transfer during FC-72 flow in a group of parallel minichannels of rectangular cross-sections that form a minichannel heat sink with three or five minichannels of 1 mm depth and 42 mm length. The results of steady-state experiments were taken into account. The key objective of the study was to develop mathematical calculation methods for the identification of the boiling heat transfer coefficient based on the data from own experiments. A mathematical model for the considered phenomenon was proposed. ADINA software and the authors’ calculation procedure based on the Trefftz functions [44] were used to carry out the FEM simulations leading to the determination of the temperature distributions in the fluid and a heated wall of a minichannel. The most important property of the Trefftz functions [45,46,47,48] is that they strictly satisfy the governing differential equation. Two novelties should be emphasized:

(i)

The new geometric layout and dimensions of the test section: in place of the test section with minichannels or minigaps 180 mm long, a system of several short minichannels (three and five) constituting a minichannel heat sink was applied (a square test section with an area of 40 mm × 40 mm contains several short parallel minichannels oriented in two spatial arrangements: vertical or horizontal);

(ii)

Calculations according to the 2D approach were performed applying the Trefftz functions and in the commercial program ADINA, the 1D approach was also used in the computations; therefore, three different mathematical calculations were applied to validate the results.

2. Experimental Database

The main loops of the experimental stand and apparatus are illustrated in Figure 1. In the flow loop, a working fluid-Fluorinert FC-72 is recirculated. The flow loop consisted of the following elements: a test section with a group of three or five minichannels—a minichannel heat sink, a gear pump, a heat exchanger, a compensating tank, a Coriolis mass flow meter, a deaerator, and a filter. The system for acquiring data and images, designed to collect measurement data, involved data acquisition stations, a PC equipped with special software, an infrared camera and a high-speed camera. Using pressure meters, the pressure of the fluid at the inlet and outlet of the test section were monitored. The supply and control system contained a power supply unit. An infrared camera enabled the measurement of temperature distributions on the heated foil. A high-speed camera was used to identify flow structures during fluid flow in minichannels. The detailed information about the main elements of the equipment realized in the experimental setup are described in [39,42].

A test section with parallel minichannels was the vital part of the experimental stand, Figure 2. The base part of the test section was made of aluminum alloy. The test section comprised three or five minichannels of 1 mm in depth and 42 mm in length. In a group of three minichannels, each channel was 11 mm wide, whereas in a group of five minichannels, they were 6 mm wide. The test section with three minichannels was oriented vertically with fluid upflow. The test section with five minichannels was oriented horizontally with FC-72 flow above the heated wall. Minichannels were created using a Teflon exchangeable element. The heated element for cooling liquid flowing along the minichannels was a thin foil, made of Haynes-230 alloy. This alloy was selected for its electrical resistivity, assuming slight changes in temperature resistivity. In the test section with three minichannels, 0.1 mm thick foil was applied as the heated wall of the minichannels, whereas the foil of 0.45 mm thickness was mounted in the test section with five minichannels. The temperature of the outer side of the foil was monitored with infrared thermography. The inlet and outlet collectors were equipped with K-type thermocouples and pressure meters.

In each experimental series, after deaeration, there was a flow of the working fluid in the main loop. Fluorinert FC-72 at temperatures below its boiling point flowed into the minichannels through an inlet collector. Heat flux supplied to the heated foil (which constituted one wall of the minichannels) was increased gradually by adjustment of the supplied heating power. The power increased gradually until subcooled boiling occurred. The temperature of the outer side of the heated foil was measured by an infrared camera. Two data acquisition stations recorded the signals from thermocouples and pressure meters (at the inlet and outlet to/from the test section), mass flow rate, the current supplied to the heated wall, and the voltage drop across the foil.

3. Experimental Uncertainties

According to the specification data released by the manufacturers, the thermal accuracy of the FLIR A655SC infrared camera is ±2 K or ±2 percent, and the mass flow rate accuracy of Coriolis mass flow meter Promass 80A04 is ±0.15% of the full scale (when applied for liquids). The overpressure values were registered with PMP71 Cerabar S (Endress + Hauser) pressure meters with accuracy ±0.05% of the full scale. The nominal measurement accuracy for a K-type thermocouple was 1.5 K in the temperature range of the conducted experiment. In [49], it was proved that the uncertainty of temperature measurements can be assumed as 0.6 K or lower due to the fact that thermocouples were additionally calibrated using 9102 HDRC thermocouple calibration instrument, manufactured by Hart Scientific. The average uncertainties of the heat transfer coefficient and heat flux in the subcooled boiling region (when the specified infrared camera was used for surface temperature measurement in experiments) were estimated as lower than 20% and 10.5%, respectively [42].

4. Mathematical Methods of Calculations

4.1. Two-Dimensional Approach—Mathematical Model and the FEM Calculations

4.1.1. Governing Equations, Boundary Conditions and Basic Assumptions

The mathematical approach presented in this part of the paper was based on the following assumptions:

• Only the central line along the minichannel axis, consistent with the flow direction, was taken into consideration;
• The heat transfer in the minichannel wall and flowing fluid was stationary and bidirectional: the y direction coincided with the flow direction, the x direction perpendicular to the flow direction referred to the depth of the minichannel;
• The independence of material properties (of the fluid and the foil) on temperature;
• The volumetric heat flux generated uniformly by the heating element was modeled as an internal heat source in the Poisson equation;
• Temperature measurements of the outer surface of the heated wall and the fluid temperature at the inlet and outlet of the minichannel were taken into account in boundary conditions;
• The fluid flow in the minichannel was laminar and incompressible with a constant mass flow rate;
• Heat losses to the environment were not taken into account (in previous works [40,41]; it was checked that they are almost negligible);
• Only one component of the fluid velocity vector (parallel to the minichannel heated surface) was not zero;
• It was assumed that the stabilized fluid temperature distribution occurred in the section distant by 25% of the channel depth from the fluid–heated foil contact interface.

Based on the above assumptions, the heat transfer in the fluid flowing in the minichannel and the heated wall can be described by the system of the following differential equations:

$${\kappa }_f {\nabla }^2{T}_f-w_y\left(x\right)\frac{\partial T_f}{\partial y}=0 \mathrm{for} \left(x,y\right)\in {\mathsf{\Omega }}_M,$$

(1)

$${\nabla }^2{T}_F=-\frac{q^{″}}{{\delta }_F·{\lambda }_F} \mathrm{for} \left(x,y\right)\in {\mathsf{\Omega }}_F,$$

(2)

where ${ \mathsf{\Omega }}_M=\left\{\left(x,y\right)\in R^2:0\le x\le {\delta }_M, 0\le y\le L \right\},{ \mathsf{\Omega }}_F=\left\{\left(x,y\right)\in R^2:{\delta }_M\le x\le {\delta }_M+{\delta }_F, 0\le y\le L\right\}$, $T_f$—the fluid temperature, $T_F$—the foil temperature, ${\nabla }^2$—Laplacian, $w_y\left(x\right)$—component of fluid velocity vector, $q^{″}$—heat flux, $q^{″}=\frac{I\mathrm{∆}U}{A}$, $A$—cross-sectional area of the heated foil, I—the current, ΔU—the voltage drop across the foil, ${\lambda }_F$—the thermal conductivity of the foil, ${\kappa }_f$—thermal diffusivity coefficient, ${\kappa }_f=\frac{{\lambda }_f }{{\rho }_f c_{p,f }}$, ${\lambda }_f—$the thermal conductivity of the fluid, ${\rho }_f—$fluid density, $c_{p,f }$—specific heat of the fluid, ${\delta }_{F }$—the foil thickness, ${\delta }_{M }$—the minichannel depth, L—the minichannel length.

Adopted boundary conditions are shown in Figure 3.

To analyze the heat transfer intensity, the local values of heat transfer coefficient on the surface between the heated foil and the boiling liquid were used. They were determined from the Robin boundary condition:

$${\alpha }_{2D}\left(y\right)=q^{″}\left(y\right)/\left(T_F\left({\delta }_M,y\right)-T_f\left({\delta }_M-x_i,y\right)\right),$$

(3)

where $T_f\left({\delta }_M-x_i,y\right)—$the fluid temperature at the selected distance $x_i=\frac{1}{4} {\delta }_M$.

4.1.2. Analytical–Numerical Method—The FEM with the Trefftz-Type Basis Functions

To calculate the wall and fluid temperature distributions, like in [50], the FEM with the Trefftz-type basis functions was used. The domain ${\mathsf{\Omega }}_M\cup {\mathsf{\Omega }}_F$ was divided into rectangular four-node elements. In each element, like in [43], the temperature was approximated by a linear combination of the basis functions which strictly satisfy the given differential equations, Equations (1) or (2), respectively:

$$T_{f }^j\left(x,y\right)= {{\displaystyle \sum }}_{k=1}^N{\widehat{T}}_f^r\left(x_{jk},y_{jk}\right)f_{jk}\left(x,y\right){ \mathrm{in} \mathsf{\Omega }}_{j,M} \mathrm{for} \mathrm{j}=1,2,\dots ,{\mathrm{J}}_M,$$

(4)

$$T_{F }^j\left(x,y\right)=u\left(x,y\right)+{{\displaystyle \sum }}_{k=1}^N\left({\widehat{T}}_F^r\left(x_{jk},y_{jk}\right)-u\left(x_{jk},y_{jk}\right)\right)g_{jk}\left(x,y\right) {\mathrm{in} \mathsf{\Omega }}_{j,F} \mathrm{for} \mathrm{j}=1,2,\dots ,{\mathrm{J}}_F,$$

(5)

where: ${\mathrm{\cup }}_{j=1}^{{\mathrm{J}}_M}{\mathsf{\Omega }}_{j,M}={ \mathsf{\Omega }}_M$, ${\mathrm{\cup }}_{j=1}^{{\mathrm{J}}_F}{\mathsf{\Omega }}_{j,F}={\mathsf{\Omega }}_F$, ${ \mathrm{J}}_{\mathrm{M}}$—the number of elements in ${\mathsf{\Omega }}_M$, ${ \mathrm{J}}_F$—the number of elements in ${\mathsf{\Omega }}_F$, $u\left(x,y\right)$—the particular solution to Equation (2), $f_{jk}\left(x,y\right)$—the basis functions specific to the Fourier–Kirchhoff equation, $g_{jk}\left(x,y\right)$—the basis functions referred to the Laplace’s equation, ${\widehat{T}}_f^r$—the temperature value in the r-th node of domain ${\mathsf{\Omega }}_M\cup {\mathsf{\Omega }}_F$, ${\widehat{T}}_F^r$—the temperature value in the r-th node of domain ${\mathsf{\Omega }}_M\cup {\mathsf{\Omega }}_F$, r—node number in the entire domain ${\mathsf{\Omega }}_F\cup {\mathsf{\Omega }}_M$, j—element number, k—basis function number in j-th element, N—the number of nodes in the element.

For example, the basis functions $f_{jk}\left(x,y\right)$ were obtained using the Lagrange interpolation with the following Trefftz functions, described in [51]:

$$1,x,y-\frac{{\mathrm{w}}_{avg}{x}^4}{2{({\delta }_M)}^2{\kappa }_f}+\frac{{\mathrm{w}}_{avg}{x}^3}{{\delta }_M{\kappa }_f},xy-\frac{3{\mathrm{w}}_{avg}{x}^5}{10{({\delta }_M)}^2{\kappa }_f}+\frac{{\mathrm{w}}_{avg}{x}^4}{2{\delta }_M{\kappa }_f}.$$

(6)

The polynomials listed in Equation (6) satisfy the Fourier–Kirchhoff equation (Equation (1) with the given parabolic profile of velocity [52], $w_y\left(x\right)=\frac{6w_{avg}}{{({\delta }_M)}^2}\left({\delta }_{M}x-x^2\right)$, where $w_{avg}$—the average fluid velocity.

The coefficients ${\widehat{T}}_f^{r }$ and ${\widehat{T}}_F^r$ were calculated by minimizing the functional, like in [53], which represented the approximate solution’s mean square error on the domain boundaries and at the common edges of adjacent subdomains.

The calculated foil temperature function enables the determination of the heat flux from Fourier’s law:

$$q^{″}\left(y\right)=-{\lambda }_F\frac{\partial T_F\left({\delta }_M,y\right)}{\partial x}.$$

(7)

4.1.3. Numerical Method—ADINA Software

ADINA software version 9.2 was used to carry out computations. The energy equation (Equation (1)), the momentum equation, and the mass conservation equation were solved (for more details see [54]). The system of boundary conditions (see Figure 3) had to be completed with the following conditions:

$$\frac{\partial T_f}{\partial x}\left(0,y\right)=0,$$

(8)

$$w\left(0,\mathrm{y}\right)=\left[0, 0\right],$$

(9)

$$w\left({\delta }_M,y\right)=\left[0, 0\right],$$

(10)

$$w\left(x,0\right)=\left[0,w_{avg}\right],$$

(11)

$$p\left(x,L\right)=p_{out},$$

(12)

where w = [w_x,w_y]—the velocity vector with components in x and y directions (w_x and w_y respectively), p_out—the absolute fluid pressure at the outlet of the minichannel, and w_avg—the average velocity of the fluid, both measured during experiments.

It was assumed that no heat losses to the environment occurred, as shown by Equation (8). Equations (9) and (10) are standard no-slip conditions at boundaries perpendicular to the flow direction. The Equation (11) defines the velocity profile at the minichannel inlet, while Equation (12) describes the pressure at the outlet of the minichannel.

It should be emphasized that in the computations using ADINA, fluid velocity is one of the variables to be found, while in calculations with the aid of the Trefftz functions a fixed distribution of the fluid velocity was assumed (see Section 4.2).

Typical FE mesh is presented in Figure 4. In computations, four-node planar FCBI (flow condition-based interpolation) elements were used. The lines parallel to the y-axis were subdivided into segments with equal lengths. Lines in the fluid domain perpendicular to the fluid flow (parallel to the x-axis) were subdivided in such a way that the side lengths of the segments close to the fluid/foil interface were narrower. At the fluid/foil interface, double nodes were considered (for the fluid and the foil separately). Constraint equations had to be applied to ensure the equality of temperatures at the fluid/foil interface in both domains.

The calculations of the local heat transfer coefficients were the main aim of the computations. Mesh density was set on the basis of the values of this coefficient. To show the process of choosing the mesh density, the data from the experiment with the test section comprising five minichannels at fixed heat flux (q″ = 13.8 kW/m²) were selected.

The symbol n is denoted by the number of segments into which the minichannel along y-axis was subdivided (see Figure 4). Lines perpendicular to the flow direction were subdivided into numbers of segments dependent on n. These numbers are denoted as a(n) for the fluid domain, and b(n) for the foil domain. Functions a(n) and b(n) are defined as follows:

$$a\left(n\right)=floor\left(n ·\frac{{\delta }_M}{L}\right),$$

(13)

$$b\left(n\right)=max\left\{floor\left(n ·\frac{{\delta }_F}{L}\right),1\right\},$$

(14)

where floor function denotes the greatest integer less than or equal to the argument.

The values of the heat transfer coefficient were obtained for selected values of n (n = 100, 200, 300, 400, 500, 600) and are shown in Figure 5. This figure represents the heat transfer coefficient vs. the distance from the minichannel inlet, while the six listed numbers of elements were tested in the numerical procedure of calculations. Maximum values of the relative differences rd_max between the values of heat transfer coefficient for adjacent values of n (n = 100 vs. n = 200; n = 200 vs. n = 300, etc.) are listed in

Table 1. The values of the maximum relative differences (rd_max) between heat transfer coefficients obtained from ADINA computations for adjacent values of n, the data for the test section with five minichannels, q″=13.8 kW/m².

Maximum Relative Differences (rdmax) between Heat Transfer Coefficients
rdmax (n = 100, n = 200)	rdmax (n = 200, n = 300)	rdmax (n = 300, n = 400)	rdmax (n = 400, n = 500)	rdmax (n = 500, n = 600)
13.3%	3%	1.7%	1.5%	0.95%

.

It is noticed that rd_max values decrease with n increases. Moreover, the maximum difference rd_max between curves for adjacent values of n dropping below 1% was assumed as effective enough because further increasing of n has not improved the result greatly. It is worth mentioning that similar results of the analysis were gained where other values of q″ were taken into account in calculations. That is the reason that a number n equal to 600 was applied in the numerical procedure carried out in most computations. The mesh used in calculations for the test section with five minichannels for n = 600 counted 31,284 elements and 33,060 nodes.

The distributions of the y-component of the fluid velocity as a function of the distance from the adiabatic wall of the minichannel (at x = 0, Figure 3) are shown in Figure 6a,b, for the test section with three and five minichannels, respectively. These distributions have been obtained at six distances from the minichannel inlet along the flow: 1 mm, 2 mm, 3 mm, 4 mm, 5 mm, and 16 mm. Additionally, velocity profiles, assumed in the calculations based on the Trefftz functions, were shown (solid black line). According to the assumption of the laminar flow, by increasing the distance from the minichannel inlet, the parabolic velocity of the fluid was observed. Furthermore, it can be noticed that the velocity profile stabilizes as the distance from the minichannel inlet increases.

The relative difference between the maximum values of w_y for the test section with three minichannels (see Figure 6a) obtained from ADINA at y = 5 mm and from the computations using the Trefftz functions was equal to 7.3%, while at y = 16 mm it was less than 6.6%. Moreover, the relative difference values turned out to be smaller for the test section with five minichannels (see Figure 6b): at y = 5 mm the difference achieved 5.5%, whereas at y = 16 mm it was less than 2.5%.

The 2D distributions of the y-component of fluid velocity according to ADINA software are illustrated in Figure 7 for the test section with three minichannels (Figure 7a) and with five minichannels (Figure 7b).

4.2. One-Dimensional Approach—Analytical Method

A simple 1D approach of local heat transfer coefficient calculations according to Newton’s law was used to compare the results. The variation of heat flux produced inside the foil along the flow and the direction consistent with the width of the minichannel was ignored in this approach. The heat transfer coefficient at the fluid/foil interface was determined at subcooled boiling region from the following dependence:

$${\alpha }_{1D}\left(y\right)=\frac{q^{″}}{T_{1D}\left({\delta }_M,y\right)-T_f\left(y\right) },$$

(15)

where $T_{1D}\left({\delta }_M,y\right)=T_{F,IR}\left(y\right)- \frac{q^{″}{\delta }_F}{{}_F}$, $q^{″}$—heat flux defined as for Equation (2), the local fluid temperature T_f is calculated as linear dependence between the measurements realized at the inlet and the outlet of the minichannel and other designations, assumed as for the 2D approach.

5. Results and Discussion

During each experiment, the current supplied to the heated wall was increased. This caused the increasing of heat flux transferred to the boiling fluid flowing in the minichannel. The experimental results pertain to the steady state (stationary condition). The study involved the identification of this coefficient in the subcooled boiling region.

The analysis includes the data from experiments performed using a minichannel heat sink with:

• Three channels, the test section oriented vertically with fluid upflow—four values of increasing heat flux: q″ = 68 kW/m², q″ = 75.7 kW/m², q″ = 100.2 kW/m² and q″ = 123.5 kW/m²;
• Five channels, the test section oriented horizontally with fluid flow above the heated wall—six values of the increasing heat flux: q″ = 13.8 kW/m², q″ = 24.6 kW/m², q″ = 31.3 kW/m², q″ = 38.7 kW/m², q″ = 47.2 kW/m² and q″ = 55.8 kW/m².

Other main experimental parameters of the selected sets are presented in

Table 2. Main experimental parameters.

The Test Section, a Minichannel Heat Sink with:	Spatial Position of the Test Section	Heat Flux (Range) q″ (kW/m2)	Inlet Pressure (Average) (kPa)	Mass Flow Rate (Average, per a Channel) (kg/s)	Inlet Liquid Subcooling (Average) (K)
Three channels	vertical, fluid upflow	68.0–123.5	105.7	0.0020/3	33.7
Five channels	horizontal, fluid flow above the heated wall	13.8–55.8	104.0	0.0026/5	30.7

.

The results are graphically described as:

• The temperature vs. the distance from the minichannel inlet at the fluid/foil interface: determined from the 2D approach (computations using the Trefftz function and ADINA software and the 1D approach temperature measured by infrared thermography at the outer foil surface: Figure 8 (the test section with three minichannels) and Figure 9 (the test section with five minichannels);
• Example temperature distributions in the foil and in the fluid obtained from the FEM computations carried out in ADINA software: Figure 10 (the test section with three minichannels) and Figure 11 (the test section with five minichannels);
• Example distribution of heat flux obtained from the FEM computations carried out by ADINA software: Figure 12 (the test section with three minichannels) and Figure 13 (the test section with five minichannels);
• The heat transfer coefficient vs. the distance from the inlet of the minichannel: achieved from 2D (computations using the Trefftz functions and ADINA software) and 1D approaches: Figure 14 (the test section with three minichannels) and Figure 15 (the test section with five minichannels);
• Example boiling curves generated for two distances from the inlet of the minichannel (0.015 m and 0.025 m)—Figure 16, based on the data from the experimental set using the test section with five minichannels;
• Example distributions of the y-component of the fluid velocity: Figure 6 (vs. the distance from the adiabatic wall of the minichannel) and Figure 7, obtained from computations in ADINA software (part “a”—the test section with three minichannels and part “b”—the test section with five minichannels).

5.1. Temperature Distributions and Heat Transfer Coefficients

It should be emphasized that the two analyzed data sets differ in experimental thermal and flow parameters, spatial positions of the test section, the number of minichannels and their width, and the thickness of the heated foil. In previous works [38,39], it was pointed out that the angle of inclination of the channel to the horizontal plane and the direction of the flow, flow rate and type of flow, as well as the pressure and inlet liquid subcooling at the minichannel inlet have an impact on the initiation and development of the boiling process in the minichannels (the first two of them seem to be the most important).

Figure 8 and Figure 9 illustrate the comparison between the temperatures at the fluid/foil interface from computations based on the Trefftz functions, ADINA software, and according to the 1D approach. Additionally, the local values of foil temperature from the measurements by an infrared camera at the outer surface of the heated foil are presented. Dependences shown in Figure 8 refer to the test section comprising three minichannels, and in Figure 9, to the test section with five minichannels.

The dependencies concerning temperatures at the fluid/foil interface are similar for both methods: according to the 1D approach and to the FEM calculations (Figure 8 and Figure 9). Generally, the results differ very slightly, and a similar course of dependencies are noticed. Comparing the results presented in Figure 8 and Figure 9, higher discrepancies between the values are observed with increasing heat flux. The temperature distributions in the foil and in the fluid obtained from the FEM computations carried out by ADINA software (Figure 10 and Figure 11) confirm the assumption that the fluid temperature stabilizes in the section distant by 1/4 of the channel depth from the fluid/heated foil contact interface (this assumption was made in Section 4.1).

Example distributions of heat flux obtained from the FEM computations carried out by ADINA software are illustrated in Figure 12 and Figure 13, for the test section with three and five minichannels, respectively. The heat flux distributions confirmed the observations from the temperature distributions illustrated in Figure 10 and Figure 11. Changes in temperature and heat flux values occur near the foil/fluid contact surface while they are almost constant close to the adiabatic minichannel wall.

In

Table 3. The heat flux calculated on the basis of experimental measurements and determined from the 2D approach calculations (average) and relative differences between them, the data for the test section with three minichannels.

Heat Flux (kW/m2)
Experiment, q″	68.6	75.7	100.2	123.5
ADINA, q″avg	66.5	76.4	99.79	122.83
Trefftz, q″avg	68.66	75.71	100.23	123.55
Relative Differences (%)
Experiment vs. ADINA	2.2	0.9	0.4	0.5
Experiment vs. Trefftz	0.09	0.01	0.03	0.04

and

Table 4. The heat flux calculated on the basis of experimental measurements and determined from the 2D approach calculations (average) and relative differences between them, the data for the test section with five minichannels.

Heat Flux (kW/m2)
Experiment, q″	13.8	24.6	31.3	38.7	47.2	55.8
ADINA, q″avg	12.97	22.82	28.77	36.62	46.96	54.54
Trefftz, q″avg	13.79	24.62	31.31	38.72	47.25	55.82
Relative Differences (%)
Experiment vs. ADINA	6	7.2	8.1	5.4	0.5	2.3
Experiment vs. Trefftz	0.07	0.08	0.03	0.05	0.1	0.04

, heat flux values calculated directly on the basis of experimental measurements (named “Experiment”) and the average heat flux determined from the 2D approach computations applying the Trefftz functions (named “Trefftz”) and ADINA software (named “ADINA”) are presented. Moreover, the relative differences between heat flux values obtained from experiments and calculated due to Trefftz method are in the range 0.01–0.09%, whereas the experimental heat fluxes compared to the resulting heat fluxes gained from ADINA are higher, reaching 8.1%.

Local heat transfer coefficients are the main results from the FEM computations (using the Trefftz functions and ADINA software). For comparison, the results of calculations according to the 1D approach are also illustrated. The presented results correspond to the temperature data presented in Figure 8 and Figure 9. Values of the local heat transfer coefficients shown in Figure 14 (determined on the basis of experimental data gained for the three minichannels) are in the range 3–11 kW/(m²K). The coefficient increases with the increase in the supplied heat flux. The highest values of the heat transfer coefficient exceed 10 kW/(m²K) near the channel outlet. Analyzing the dependencies shown in Figure 15, it is noticed that the coefficient increases with the distance from the minichannel inlet. Furthermore, comparing the results according to both methods applied in the FEM showed that the heat transfer coefficient achieved similar values.

The coefficient calculated from the 2D approach also has similar values to those obtained from the computations according to the simple 1D approach. Moreover, higher discrepancies between the local values of the heat transfer coefficient are detected near the channel outlet. In turn, the relationship illustrated in Figure 15 (obtained for the test section with five minichannels) showed the coefficient diminishing with increasing the distance from the channel inlet. Maximum values of the heat transfer coefficient did not exceed 3.3 kW/(m²K) at the minichannel inlet, among all results from the computations. The highest discrepancies of the coefficient values obtained from the considered mathematical methods are clearly noticeable near the minichannel outlet. The values of the maximal relative differences (rd_max) and average differences (rd_avg) between heat transfer coefficients obtained using two approaches (1D and 2D) are listed in

Table 5. The values of the maximum (rd_max) and average (rd_avg) relative differences between heat transfer coefficients obtained using the 1D approach and 2D approach (applying the Trefftz functions and ADINA software), the data for the test section with three minichannels.

Relative Differences (rd) between Values of the Heat Transfer Coefficient (%)
Experiment, q″ (kW/m2)		68	75.7	100.2	123.5
1D approach vs. 2D approach, Trefftz functions	rdmax	24.6	18.7	25.4	22.6
	rdavg	12.8	5.4	7.24	7
1D approach vs. 2D approach, ADINA software	rdmax	22.7	31.6	43.5	40.8
	rdavg	6.4	8.7	12.8	11.4
2D approach, Trefftz functions vs. 2D approach, ADINA software	rdmax	16.5	12.1	17.4	14.5
	rdavg	9.4	3.3	5.8	6.5

and

Table 6. The values of the maximum (rd_max) and the average (rd_avg) relative differences between heat transfer coefficients obtained using the 1D approach and the 2D approach (applying the Trefftz functions and ADINA software), the data for the test section with five minichannels.

Relative Differences (rd) between Values of the Heat Transfer Coefficient (%)
Experiment, q″ (kW/m2)		13.8	24.6	31.3	38.7	47.2	55.8
1D approach vs. 2D approach, Trefftz functions	rdmax	13.4	9.9	10.2	9.6	12.6	16.3
	rdavg	4	4.1	4.2	4.6	5.2	6.5
1D approach vs. 2D approach, ADINA software	rdmax	8.1	10.8	11.7	8.1	7.8	7.4
	rdavg	2.6	2.2	3.3	2.3	2.8	3.3
2D approach, Trefftz functions vs. 2D approach, ADINA software	rdmax	8.3	11.4	8.3	13.9	13.7	16.2
	rdavg	2.3	3.9	2.8	6	8.1	9.9

, based on the data for the test section with three and five minichannels, respectively.

When analyzing the results listed in Table 5 and Table 6, it is seen that the average (rd_avg) relative differences between the heat transfer coefficients obtained using the 1D approach and the 2D approach differ from 5.4% to 12.8% (three minichannels) and from 2.2% to 6.5% (five minichannels). The reason for these computational discrepancies is probably due to assuming a linear dependence of the fluid temperature from the inlet to the outlet in the calculation procedure introduced in the 1D approach. According to the 2D approach, the fluid temperature results from the FEM calculations were in a strictly defined position of the reference line. Moreover, it is obvious that with increasing heat flux, the temperature gradient in the boundary layer increases, but the 1D approach does not take this into consideration. Therefore, the 2D approach is expected to give more accurate values of the heat transfer coefficient. Average relative differences between the FEM calculations using ADINA and the Trefftz functions did not exceed 9.4% (the test section with three minichannels) and 9.9% (the test section with five minichannels). In most of the analyzed cases, for the highest heat fluxes, the difference between the local heat transfer coefficients calculated from two mathematical approaches (1D and 2D) was higher.

It is worth mentioning that the differences in the values of the heat transfer coefficient determined from both the FEM calculations are usually greater in the channel outlet in comparison with the results from calculations according to the 1D approach, especially for higher values of heat flux. It may be explained when analyzing each calculation model applied in the FEM as follows: (i) in the model using the Trefftz functions only the energy equation is solved (Equation (1), the unknown function to be found is a temperature function while the velocity profile is assumed; (ii) in the FEM calculations using ADINA software, the full set of CFD equations is solved (mass, momentums and energy equations), and pressure, temperature and velocity are functions to be designated.

In conclusion, the resulting heat transfer coefficients from the 1D approach are usually lower than from the 2D approach; the maximum discrepancies are approximately 40% for higher values of the heat transfer coefficient. Moreover, it can be stated that similar dependencies and values of the heat transfer coefficient were obtained in both the FEM calculations (average relative differences between the values of the heat transfer coefficient ranged from a few to several percent).

5.2. Boiling Curves

Example boiling curves as heat flux dependence in function of the heated wall and fluid temperature difference are shown in Figure 16. The curves were constructed on the basis of the data from the experiment using the test section with five minichannels, for two distances from the minichannel inlet: 0.015 m and 0.025 m. The course of the boiling curves is typical for refrigerants and confirms the previous results of the authors, published, among others, in Ref. [39].

The boiling curves course given in Figure 16 indicates that while the subcooled liquid flows into a minichannel, an increase in heat flux causes boiling incipience (the bottom line of the graph, ONB—the onset of nucleate boiling). In the area adjacent to the heated wall of the minichannel, the fluid was superheated, but in the flow core it was subcooled. Impulsive nucleation causes a temperature drop of the heated wall surface, called “nucleation hysteresis”. The bubbles absorb a great amount of energy transferred to the liquid, acting as internal heat sinks. Analyzing the course of both boiling curves, it is observed that the highest temperature difference drop at ONB diminishes with increasing the distance from the minichannel inlet.

6. Conclusions

Boiling heat transfer during FC-72 flow in a minichannel heat sink oriented vertically or horizontally was the main aim of this work. The key objective was to develop mathematical calculation methods for the identification of the boiling heat transfer coefficient based on the data from own research. The experimental setup comprised a test section with three or five channels, oriented in a different angle from the horizontal plane, and were selected to provide a variety of databases. During experiments, temperature measurement of the outer foil surface constituting the common heated wall of minichannels was provided by an infrared camera. Selected sets differed in experimental thermal and flow parameters, number, width of minichannels, spatial position of the test section, and thickness of the heated foil.

It was proposed that the mathematical model assumed a steady state heat transfer process and a laminar, incompressible flow of the fluid in an asymmetrically heated central minichannel. The data necessary for the creation of this mathematical model (boundary conditions, etc.) were taken from the measurements. The data from experiments at the subcooled boiling region were used for calculations. The FEM calculations (carried out simultaneously by the Trefftz functions and ADINA software) were conducted to find the temperature field in the flowing fluid and in the heated wall. Numerical computations due to ADINA software included testing the effect of the mesh density on the values of the heat transfer coefficient. Distributions of the y-component of the fluid velocity were also examined. It was observed that with increasing the distance from the minichannel inlet, the parabolic velocity profile of the fluid stabilized, and its shape did not differ much from the distributions of the y-component of the fluid velocity assumed in the computations based on Trefftz functions.

In the two-dimensional approach, the heat transfer coefficient was determined using the Robin condition. Moreover, a simple analytical method (the one-dimensional approach) was used in the calculations for comparison. The results were illustrated by graphs of the foil and fluid temperatures versus the distance from the minichannel inlet. Local heat transfer coefficients were the main results from the FEM computations (using the Trefftz functions and ADINA software). For comparison, the results of calculations according to the 1D approach were presented. Temperature distributions in the heater and the fluid obtained from the FEM computations carried out by ADINA software were also shown. The values of the heat transfer coefficient obtained from the 2D approach were agreed satisfactorily (average relative differences did not exceed 9.9%). The resulting heat transfer coefficients from the 1D approach were usually lower than from the 2D approach. The relative differences between the 1D and 2D approaches were greater (average relative difference reached 12.8%). The differences in the values of the heat transfer coefficient determined from both the FEM calculations and according to the 1D approach resulted from its simplification (the full set of CFD equations was set in the FEM computations).

Example boiling curves indicating nucleation hysteresis were shown and discussed. The course of the boiling curves was typical for refrigerants indicating nucleation hysteresis.

In the future, research on flow boiling of other refrigerants during stationary and time-dependent experiments will be conducted. The analysis of the effects of geometry and spatial orientation of channels with selected experimental parameters and enhanced heater surface on boiling heat transfer are planned. Visualization and analysis of two-phase flow pattern during flow boiling will be studied. Computation will be carried out using the commercial programs (ANSYS CFX/Fluent and STAR-CCM+) and own calculation approaches based on hybrid methods combining the Trefftz functions with Picard’s iteration method, the homotopy perturbation method, and the Beck method. Verification of the results due to comparison with the data known from literature for the interpretation of own experimental results and suitability of calculation methods are planned.

The innovative construction of a compact heat exchanger with a group of rectangular mini- and microchannels, with enhanced surfaces, will improve the heat transfer process. Such devices, which enable effective cooling of components during operation, can be applied in a number of industries such as microelectronics. There is a large contemporary market of potential users interested in high-efficiency heat exchangers.