Analysis of Heat Transfers inside Counterflow Plate Heat Exchanger Augmented by an Auxiliary Fluid Flow

Enhancement of heat transfers in counterflow plate heat exchanger due to presence of an intermediate auxiliary fluid flow is investigated. The intermediate auxiliary channel is supported by transverse conducting pins. The momentum and energy equations for the primary fluids are solved numerically and validated against a derived approximate analytical solution. A parametric study including the effect of the various plate heat exchanger, and auxiliary channel dimensionless parameters is conducted. Different enhancement performance indicators are computed. The various trends of parameters that can better enhance heat transfer rates above those for the conventional plate heat exchanger are identified. Large enhancement factors are obtained under fully developed flow conditions. The maximum enhancement factors can be increased by above 8.0- and 5.0-fold for the step and exponential distributions of the pins, respectively. Finally, counterflow plate heat exchangers with auxiliary fluid flows are recommended over the typical ones if these flows can be provided with the least cost.


Introduction
Counterflow plate heat exchangers are widely used in various engineering applications especially preheat, chemical, pharmaceutical, and food processing applications [1]. This is because both hot and cold fluids within the plate heat exchanger are exposed to a much larger surface area per unit volume than that in the conventional (double pipe) heat exchanger [2]. Also, plate heat exchangers can have hydraulic diameters smaller than 2 mm. This can lead to having larger heat transfer coefficients. Thus, plate heat exchangers have larger effectiveness compared to conventional counterflow heat exchangers. Additionally, many of the passive heat transfer enhancement tools like fins and rough surfaces [3][4][5] can easily be installed in the plate heat exchanger as compared to the conventional heat exchanger. This is why finned plate heat exchangers [6] and gasketed plate heat exchangers [7] are widely spread in many industrial applications.
The most recent literature reviews on passive heat transfer enhancements in heat exchangers [8,9] show that the major analyzed enhancement methods are the following: (1) twisted tape, (2) wire coil, (3) swirl flow, (4) conical ring, and (5) ribs. All of these devices augment heat transfer because they tend to disturb the fluid flows [3,10]. Therefore, it can be concluded that enhancing heat transfer in plate heat exchangers under laminar flow conditions did not receive much attention by researchers. Perhaps the most recent proposal for heat transfer enhancement in heat exchangers under laminar flow conditions is the use of nanofluids [11][12][13]. However, not all nanofluids can be adequate for processing special products like pharmaceutical and food products. This is because the commonly used nanoparticles can be harmful to human body [14,15]. Consequently, the present work aims to propose and analyze a new method for enhancing heat transfer in plate heat exchanger without altering either the velocity profiles or compositions of both hot and cold fluids.
The proposed plate heat exchanger is composed of hot and cold fluid channels separated by an auxiliary fluid channel. This auxiliary channel may contain as many passive enhancement tools as possible. Accordingly, both the velocity profile and the composition of the hot and cold fluids are preserved. The heat transfer enhancement in the proposed system is due to the following combined effects: (1) convection of the auxiliary fluid and (2) passive enhancement mechanisms in the auxiliary channel. In the present work, transverse pins connecting the facing boundaries of both 2 The Scientific World Journal hot and cold fluid channels are considered as one of passive enhancement mechanisms [16,17]. Moreover, the auxiliary fluid is considered to flow in the direction cross to both hot and cold fluid flow directions. Accordingly, the auxiliary channel length (hot/cold channels width) can be selected to be small enough to have boundary layer flows [18,19]. Hence, convection thermal resistances between the auxiliary fluid and both hot and cold fluids are minimized. The heat transfer rates within the present system are expected to be higher than those in conventional system for specific auxiliary flow conditions. Accordingly, the present work additionally aims to identify some trends of parameters that cause enhancement ratios to be above unity. Modeling laminar flow and heat transfer inside two dimensional channels including auxiliary channels is well established in the literature [18][19][20][21][22].
In the present work, heat transfer inside plate heat exchanger with auxiliary fluid channel separating the hot and cold fluid channels is modeled and analyzed. Both hot and cold fluid flows are considered to be laminar under hydrodynamically fully developed condition. The energy equations of the hot and cold fluids are coupled with the energy equations of the auxiliary fluid boundary layers. The solution of the momentum and energy equations within the boundary layers is well established [18][19][20]. Accordingly, both coupled hot and cold fluid energy equations are solved numerically using finite difference methods. Approximate analytical solutions for the heat transfer rates under the fully developed flow and very long pins conditions are derived. A number of heat transfer performance ratios including the heat exchanger effectiveness ratios are computed. A parametric study for heat transfer enhancement is made to recognize the conditions of controlling parameters that produce favorable enhancement factors.

Modeling of Flow and Heat Transfer inside the Hot and Cold Fluid Channels.
Consider two parallel channels of length and width . The first channel confines the hot fluid flow while the second one contains the cold fluid flow. The solid boundaries of these channels facing each other are perfectly connected together via cylindrical pins of diameter and length as shown in Figure 1. These pins are surrounded by an auxiliary fluid stream of free stream temperature ∞ and convection heat transfer coefficient ℎ . The convection heat transfer coefficient between the auxiliary fluid stream and the channel boundaries facing that stream is ℎ . Accordingly, heat transfers between the channels and the auxiliary fluid are due to convection over the unfinned surfaces and conduction via the connecting pins. In contrast, the outermost solid boundaries of the heat exchanger are considered to be adiabatic so that heat transfer rates in the heat exchanger are maximized.
The dimensionless momentum and energy equations of the hot and cold fluids are [21] 2 ℎ, 2 ℎ, where ℎ and are the dimensionless axial velocity fields for the hot and cold fluids, respectively. ℎ and are the hot and cold fluid dimensionless temperatures, respectively. ℎ and are the dimensionless axial positions of the hot and cold fluids, respectively. ℎ and are the dimensionless transverse positions of the hot and cold fluids, respectively. The channels aspect ratios ℎ, as well as the hot and cold flow Péclet numbers Pe ℎ and Pe are given by where ℎ and are the heights of the hot and cold fluid channels, respectively. ( ℎ , ), ( ℎ , ), and ( ℎ , ) are the density, specific heat, and thermal conductivity pairs of the hot and cold fluids, respectively. ℎ and are the mean axial velocities of the hot and cold fluids, respectively.
where ℎ and ℎ are the axial and transverse positions of the hot fluid, respectively. and are the corresponding positions of the cold fluid, respectively. ℎ and start from zero at the fluid inlets while ℎ and start from zero at the adiabatic boundaries. ℎ and are the axial velocity fields for the hot and cold fluids, respectively. ℎ and are the hot and cold fluid temperatures, respectively. ℎ1 and 1 are the inlet temperatures of the hot and cold fluids, respectively.
The boundary conditions of (2) are given by where ℎ and are the conduction heat fluxes through the pin bases at the hot and cold surfaces, respectively. However, ℎ and are the convection heat fluxes at the unfinned portions of the hot and cold surfaces, respectively. The local pins base area concentration denoted by can be calculated from the following expression: where / is the local axial gradient of the number of pins ( ).

Modeling of the Pins Conduction and Auxiliary Fluid
Convection Heat Fluxes. The one-dimensional fin equation [18,19] can be used to model the conduction heat transfer through the pins. This fin equation has the following dimensionless form: where the dimensionless pin distance, , the dimensionless pin local temperature, , and the dimensionless pin thermal length, * , are given by where Bi = ℎ / is the pin Biot number and = / is the pin aspect ratio. The boundary conditions of (9) are given by where ℎ, = ℎ ( ℎ , ℎ = 1) and , = ( = 1 − ℎ , = 1) are the dimensionless temperatures of the hot and cold boundaries, respectively. is the dimensionless cold excess temperature. It is equal to where 0 ≤ ≤ 1 as auxiliary fluids are usually hotter than the cold reservoir. The solution of (9) is Therefore, the conduction heat flux at the pin bases is equal to The Scientific World Journal Note that ℎ = cond | =0 and = − cond | = . Recall , where ℎ and are the temperatures at the hot and cold boundaries facing the auxiliary fluid, respectively. As such, (7c) and (7d) change to where Bi ℎ = ℎ ℎ / ℎ is the hot fluid Biot number.

The Heat Transfer Rates through the Heat Exchanger.
The heat transfer rate per unit width from the hot fluid ( ℎ ) and that to the cold fluid ( )can be calculated from where ℎ2 and 2 are the mean bulk temperatures at the hot and cold fluid exit ports, respectively. ℎ2 and 2 are the dimensionless values of ℎ2 and 2 , respectively. In terms of the dimensionless parameters, ℎ and are given by The integral form of (2) can be expressed as 2.4. Modeling of the Pin Base Area Distribution over the Channel Boundary. Different distributions for pins will be analyzed in this work. These distributions allow having concentrated pin distribution either far from the middle section of the channels or about this section. Two families of distributions are considered. They are the step function and exponential distributions. The step function distribution that has concentrated pins far from the channels midsection (see Figure 2) has the following mathematical form: where 0 ≤ ( / ) ≤ ( / ) = 1 − ( / ). is the maximum value of that produce 99% of the upper limit of pins base area concentration ( ).
The exponential distribution of the pins has the following functional form: The Scientific World Journal 5 where | | < . is the upper limit value that makes = . It can be accurately correlated to through the following correlation: The pins are more concentrated near the channels inlet/exit sections when > 0. However, they are more concentrated around the channels midsection when < 0. These trends are seen in Figure 2.

Hot and Cold Fluid Flow Nusselt Numbers.
The convection heat transfer coefficient for hot and cold fluid flows ℎ ℎ and ℎ , respectively, is defined as Thus, the local Nusselt numbers Nu ℎ and Nu are equal to

Heat Exchanger Effectiveness
Ratios. The maximum heat transfer rate per unit width from the hot fluid ℎ Max and that to the cold fluid Max are obtainable when ℎ2 = 1 and 2 = ℎ1 . Using (16), ℎ Max and Max are equal to Define the heat exchanger effectiveness factor ℎ as the ratio of heat transfer rate from the hot fluid to ℎ Max . Also, is defined as the ratio of heat transfer rate to the cold fluid to Max . As such, ℎ and are mathematically equal to

Heat Exchanger Second Set of Performance Ratios.
Let the reference case for the second performance ratios be the counterflow heat exchanger with perfect indirect contact between the hot and cold fluids. For this case, the boundary conditions given by (7c) and (7d) change to The heat transfer rate between the two fluids for this case, , is equal to where [ ℎ2 ] and [ 2 ] are the dimensionless exit mean bulk temperatures of hot and cold fluids, respectively, for the reference case. Define the heat exchanger performance indicators ℎ and as the ratio of heat transfer rate from the hot fluid and that to the cold fluid, respectively, to the reference heat transfer rate. Mathematically, they are equal to . (29)

Heat Exchanger Set of Performance Ratios due to Stratified Pin Distribution.
The last set of performance indicators for the present heat exchanger denoted by ( ℎ , ) are defined as the ratios of the heat transfer rate from the hot fluid and that to the cold fluid to the corresponding quantities when = , respectively. Mathematically, they are equal to 2.6. Analytical Model. Utilizing (15a), (15b), and (24), it can be shown that where ℎ ℎ, ≪ ℎ . This condition is necessary as the aim of introducing the auxiliary fluid flow is to enhance heat transfer 6 The Scientific World Journal inside the hot and cold fluid channels. The coefficients 1 and 2 are given by Using (31), the heat transfer rates at the hot and cold differential boundary elements are given by where is equal to (34) Both sides of (33) can be arranged by separation of variables to the following differential equation: Integrating (35) over the heat exchanger length results in the following solution: Using (26), ℎ and can be shown to be equal to where (Bi ℎ ) eff and (Bi ) eff are equal to The solutions of (2) for this case under fully developed flow condition can be numerically obtained with high accuracy by following the methodology described in Section 3 of this work. For this case, it can be shown that the following correlations of Nusselt numbers where ℎ = ( / ℎ )( ℎ / ) and = ( / )( / ).
The Scientific World Journal 7

Upper Limits of Heat Fluxes in Absence of Pins at
Fully Developed Condition. The dimensionless maximum convection heat fluxes in absence of pins can be obtained when ( ℎ , ) = ( ℎ1 , 1 ) and → 0. They are equal to where Bi = Bi ℎ ( ℎ / )( / ℎ ).

Modeling of Auxiliary Fluid Flow Convection Coefficients and Minimum .
Let the auxiliary stream be laminar flow along the direction of the channels width axis. Therefore, the convection heat transfer coefficients for the fin and unfinned surfaces can be computed using the following correlations [19]: where , Pr , and ] are the auxiliary fluid thermal conductivity, Prandtl number, and kinematic viscosity, respectively.
∞ is the auxiliary free velocity. Re = ∞ /] and Re = ∞ /] are the Reynolds numbers for the streams across the pins and along the unfinned surface, respectively. Thus, Bi ℎ , * and the relationships between the latter Reynolds numbers are equal to The minimum requirements for average pin base area concentration can be obtained when the heat transfer rates between (hot, cold) fluids and auxiliary stream are equal to those between the hot and cold fluids under perfect indirect contact between the fluids. These quantities can be obtained analytically for the ideal cases: (1) ℎ Pe ℎ and are very large, while Pe is very small and (2) Pe and are very large, while ℎ Pe ℎ is very small. For these conditions, the minimum average pin base area concentration denoted by ℎ Min, Min , respectively, can be shown to be equal to

Numerical Methodology and Results
Equation (2) is coupled via the boundary conditions given by (15a) and (15b). These equations can be solved by iterations using the implicit finite difference method discussed by Khaled and Vafai [23]. Equation (2) was discretized by employing three-point central differencing quotients for the first and second derivatives with respect to ℎ, directions. Furthermore, two-point backward differencing quotients were used in the discretization of the first derivatives with respect to ℎ, directions. For (15a) and (15b), three-point central differencing quotients were used to discretize the first derivatives with respect to ℎ, directions. The finite difference equations of (2) are given by The pairs ( , ) and ( = − + 1, ) represent the location of the discretized points in the numerical grids of the hot and cold fluid domains, respectively.
is the total number of either or sections. ℎ and are the total number of the discretized points per and sections, respectively. , ℎ , and were taken to be equal to = 1001 and ℎ = = 201. The applications of (46) for all discretized points at given and sections result in ℎ and tridiagonal systems of algebraic equations. These equations can easily be solved using the Thomas algorithm [24], if internal boundary temperatures are known. The numerical solution procedure is summarized in the following steps.

8
The Scientific World Journal (4) The corrected cold boundary temperatures ( ) corrected were found using the finite difference equation of (15b).
(5) Steps (2)-(4) were repeated by replacing ( ) assumed with the corrected ( ) corrected until the following condition is satisfied: Using doubled mesh sizes results in less than 0.3% error in the calculated parameters for moderate ℎ Pe ℎ and Pe values. This ascertains the grid-size independent results. The numerical results shown in Figures 3-12 were generated for the hot, cold, and the auxiliary stream fluids shown in Table 1. These selections correspond to an important application which is oil cooling cold water stream augmented by an air stream. The numerical results were compared with the analytical solution given by (35), (36), and (40) as shown in Table 2.
The numerical results are seen to have a good agreement with the analytical solution as ℎ, Pe ℎ, ≪ 1. The latter constraint is the major assumption used to generate (40).

The Role of Internal Flow Reynolds Numbers in the Performance Ratios.
As Re ℎ increases, both hot fluid mean bulk and heated boundary temperatures decrease due to the increase in advection and the widening effect of the thermal entry region. As a result, the convection heat transfer rate to the auxiliary stream and conduction through the pins increase due to the increase in the heated boundary and pins excess temperatures ( ℎ − ∞ , ℎ − ), respectively. These excess temperatures increase as decreases. Accordingly, the heat transfer rate from the hot fluid increases which causes ℎ to increase as both Re ℎ and (1− ) increase as seen in Figure 3. This figure shows that ℎ decreases as Re increases. This indicates that the increase in the reference case heat transfer rate ( ) is larger than that for the present system ( ℎ ). Also, it is noticed from Figure 3 that ℎ can be larger than one at both smaller and Re values and larger Re ℎ values. Similarly, can be larger than one for smaller (1 − ) and Re ℎ values and larger Re values.
The increase in Re ℎ decreases the hot fluid effectiveness ℎ as shown in Figure 4, since it majorly results in reduction of the hot fluid mean bulk temperature. However, a decrease in is noticed to cause an increase in ℎ due to the associated Table 2: A comparison between the numerical solution and the analytical one given by (37)  increase in ℎ − ∞ . Figures 7 and 9 show that an increase in Re causes an increase in ℎ . This is because both ∞ − , ℎ − increase as Re increases; thus, the conduction through the pins is enhanced. The latter enhancement cannot be clearly identified from Figure 4 as the pin aspect ratio is very small for this figure which is = 0.05. As indicated earlier, when → 0 the heat transfer rate in both channels will be uncoupled. In a similar manner, it can be concluded that decreases as Re increases and as (1 − ) decreases and it increases as Re ℎ increases as noticed in Figure 4.

The Role of Pins Aspect Ratio in the Performance Ratios.
As increases ( decreases), the pin surface area increases causing the fully developed maximum fluxes Θ ℎ and Θ to increase until they reach their asymptotic values as shown in Figure 5. This causes ℎ to increase as decreases when = 0.25 as shown in Figure 6. When = 0.5, Θ ≫ Θ ℎ ; thus, sharply increases causing sharp reduction in ℎ − . Therefore, ℎ decreases as decreases when = 0.5. The cases considered in Figures 6 and 7 have ℎ Pe ℎ ≫ Pe > 1. Thus, the hot fluid flow is dominated by the thermal entry region. For this condition, ℎ is less sensitive to , while increases apparently as increases. As a result, ∞ − , and ℎ − decrease as increases. Since pins  conduction is linearly dependent on ∞ − , and ℎ − while it is less sensitive to , decreases as decreases as seen in Figure 6. Due to the previous analysis, ℎ and increase as increases except when = 0.25 where ℎ is noticed in Figure 7 to decrease as increases. Also, it is shown from this figure that ℎ > 1 when = 0.25 and

The Role of Pins Base Area Concentration in the Performance Indicators.
Two limiting cases can be encountered in the present heat exchanger. They are as follows: (1) pure convection between the channels and the auxiliary stream when → 0 and (2) pure conduction between the channels when → . As seen in Figure 5, (Θ ℎ , Θ ) > (Θ ℎ , Θ ) when = 0.2 which represents the condition for Figures 8  and 9. Accordingly, ( ℎ , ) pair is expected to decrease as increases. This is shown in Figure 8 except for the ℎ plot when = 0.25. For this case, ℎ − ∞ , and ℎ − are very close to their upper limit ( ℎ1 − 1 ) while this limit is reduced to ( ℎ1 − ∞ ) for the pure convection condition as → 0. And since the heat flux is linearly proportional ℎ − ∞ , and ℎ − as can be seen in (15a) and (15b), ℎ is increased when is increased for = 0.25. Because of the previous facts, ( ℎ , ) pair decreases as increases except when = 0.25 where ℎ is noticed in Figure 9 to increase as increases. Furthermore, it is shown from this figure that ℎ > 1 when = 0.25 and Re = 10 while > 1 when = 0.5. This demonstrates the superiority of the present heat exchanger over the conventional counterflow plate heat exchanger.

The Role of Pins Distribution in the Performance
Indicators. Since = 0.25 and ℎ Pe ℎ ≫ Pe > 1 in Figures 10, 11, and 12, ℎ − ∞ , and ℎ − are very close to the maximum value ( ℎ1 − 1 ). These excess temperatures become closer to that limit near the hot fluid inlet. Far from this region, they tend to apparently decrease as decreases since Θ ℎ > Θ ℎ . Consequently, the uniform distribution of pins reveals the maximum ℎ in which ℎ = 1 as seen from these figures. On the other hand, ∞ − turns out to be much smaller than ( ∞ − 1 ) particularly near the cold fluid exit Increases in ℎ Pe ℎ , and Pe reduce ℎ , and , respectively. As a result, ℎ − ∞ , and ∞ − increase causing {( ℎ ) max , ( ) max } to increase as ℎ Pe ℎ , and Pe increase, respectively. These trends are seen clearly in Figures  13 and 14. These maximum ratios are obtained based on {Nu * ℎ , Nu * } approximations. Some critical parameters that produce the maximum quantities {( ℎ ) max , ( ) max } are listed in Table 3. Finally, {( ℎ ) max , ( ) max } as seen from Figures  13 and 14 can be much larger than the one indicating that properly distributing the pins is an efficient mechanism for heat transfer enhancement under fully developed laminar flow condition.

Conclusions
Heat transfer inside counter-flow plate heat exchanger subject to internal convections with an auxiliary fluid was investigated in this work. The auxiliary fluid passage is surrounded by the hot and cold fluid channels and it is supported by highly conductive pins connected to both channels. Good   agreement was noticed between the numerical solution and an approximate analytical solution based on fully developed flow and very long pin conditions. The results of the current study can be summarized by the following concluding remarks.
(1) The heat transfer rate from/to the hot/cold fluid of the present system can be higher than that for the conventional counterflow plate heat exchanger under the following conditions: (1) large hot/cold flow Reynolds number, (2) small cold/hot flow Reynolds number, and (3) large hot/cold fluid excess temperature.  (2) Increasing either the number pins or pins length may increase the hot/cold heat transfer rate above that for the conventional counterflow plate heat exchanger under the following conditions: (1) effective hot/cold fluid Biot number due to only pins conduction being larger than that due to only convection with the unfinned surface and (2) large hot/cold fluid excess temperature as when having large hot/cold flow Reynolds number.