Conjugate Heat Transfer Analysis of PCM Suspensions in a Circular Tube under External Cooling Convection: Wall Conduction E ﬀ ects

Featured Application: heat exchanger design; thermosiphon design. Abstract: In this paper, a numerical method is used to investigate the conjugate heat transfer of a phase change material (PCM) suspension in a circular tube under external cooling convection. The following parameters and ranges were considered: dimensionless tube wall thickness, t w (0–0.5); wall-to-ﬂuid thermal conductivity ratio, k ∗ wf (0.1–10); volumetric fraction of PCM particles, c v (0.1); Biot Bi (1); Stefan number, Ste (0.1); and Peclet number, Pe (1000). The results show that the wall thermal conductivity considerably a ﬀ ects the outer / inner wall temperature of the tube, the average temperature of the working ﬂuid, and the volumetric liquid fraction of PCM particles. Thus, wall conduction e ﬀ ects must be properly accounted for to model heat transfer in a PCM suspension in tube ﬂow.


Introduction
Fine solid-liquid phase-change material (PCM) particles can be dispersed in a carrier or suspended in a working fluid to form PCM slurries or suspensions. The latent heat of absorption and release of the PCM can be used to reduce the maximum temperature and increase heat transfer in a loop. PCM particles can be prepared by microencapsulation technology [1][2][3][4][5] or emulsion technology [6][7][8].
In microencapsulation, polymer materials cover the PCM to create a barrier between the dispersed phase (PCM) and the dispersion medium (such as water and glycerol) that results in the suspension of PCM particles in the dispersion medium. The prepared PCM particles have a relatively fixed size and shape and relatively stable physical properties compared to particles prepared by other methods. However, these materials have the disadvantages of high cost, a difficult manufacturing process, reduced heat transfer efficiency from the polymer surface coating, and expansion and cracking after repeated heating. In emulsion technology, the addition of appropriate emulsifiers reduces the surface tension between the dispersed phase and the dispersion medium to produce a stable emulsion of PCM particles uniformly suspended in the dispersion medium.
Heat transfer enhancement using PCM suspensions as the working fluid has been demonstrated in forced convection. Inaba et al. [9] proposed a numerical simulation and an experiment on an isothermal and electronic cooling systems. The physical modeling of this study was inspired by the cooling section of the thermosiphon. The physical configuration of the conjugate heat transfer problem considered, as illustrated schematically in Figure 1a, is a circular pipe of inner radius r + i and wall thickness t + w (= r + o − r + i ). A tube is divided into three sections: (1) an upstream adiabatic section of length L + u ; (2) a heat exchange area with a convection boundary (T ∞ , h ∞ ) and a length of L + h , and (3) a downstream adiabatic section of length L + d . The PCM suspension enters the upstream adiabatic section at a fully developed velocity u + m , and the fluid has a uniform temperature T in . The following basic assumptions are made in the model.

1.
The working fluid is a suspension of PCM particles. The PCM particles can be completely solid phase, completely liquid phase, or a mixture of solid and liquid phases.

2.
A uniform diameter is assumed for the PCM particles suspended in the working fluid. 3.
All of the physical properties of the PCM suspension are assumed to be constant.

4.
The PCM suspension is regarded as a thin mixture.

5.
The PCM suspension is regarded as an incompressible Newtonian fluid. 6.
The PCM particles are neutral and evenly distributed in the working fluid. 8.
The PCM suspension enters the tube in fully developed laminar flow. 9.
The flow and temperature distribution of the PCM suspension in the tube are radially axisymmetric. 10. The change in the density of the PCM particles during solid-liquid phase changes is neglected.
Appl. Sci. 2020, 10, x FOR PEER REVIEW  3 of 15 and wall thickness (= -). A tube is divided into three sections: (1) an upstream adiabatic section of length ; (2) a heat exchange area with a convection boundary ( , ℎ ) and a length of , and (3) a downstream adiabatic section of length . The PCM suspension enters the upstream adiabatic section at a fully developed velocity , and the fluid has a uniform temperature . The following basic assumptions are made in the model.
1. The working fluid is a suspension of PCM particles. The PCM particles can be completely solid phase, completely liquid phase, or a mixture of solid and liquid phases. 2. A uniform diameter is assumed for the PCM particles suspended in the working fluid. 3. All of the physical properties of the PCM suspension are assumed to be constant. 4. The PCM suspension is regarded as a thin mixture. 5. The PCM suspension is regarded as an incompressible Newtonian fluid. 6. Viscous dissipation is neglected. 7. The PCM particles are neutral and evenly distributed in the working fluid. 8. The PCM suspension enters the tube in fully developed laminar flow. 9. The flow and temperature distribution of the PCM suspension in the tube are radially axisymmetric. 10. The change in the density of the PCM particles during solid-liquid phase changes is neglected.  Based on the foregoing assumptions, the following dimensionless governing equations and boundary conditions are used in this study.
(1) Governing equation The fluid energy equation and the tube wall heat conduction equation are as follows: where for the fluid region (0 < r < 1), u = u + u + m = 2 1 − r 2 , a = 1, and b = ρ * b f /ρ * p f ; and for the wall region The governing equation for the volumetric liquid fraction of a PCM particle ξ b is as follows: A detailed treatment of the solid-liquid phase change of particles can be found in Ho et al. [10].
(2) Boundary conditions The temperature and heat flux are assumed to be continuous at the inner wall of the tube: Thus, the mathematical formulation of the conjugate forced-convection problem consists of Equation (3), the dimensionless differential Equations (1) and (2), and the boundary conditions. The dimensionless parameters for the problem include the bulk fluid Pe; the Stefan number, Ste; the inlet dimensionless temperature, θ in ; the dimensionless particle diameter, d p ; the volumetric fraction of PCM particles, c v ; the property ratios for the PCM particles with respect to the suspending fluid, p,w f , and c * p,p f ; and the dimensionless tube wall thickness, t w . The particle Biot number, Bi p and the property ratios of the bulk with respect to the fluid, k * b f , ρ * b f , and c * p,b f , which appear in the foregoing formulation, are evaluated using the relationships given in Charunkyakorn et al. [16].

Dimensionless Thermal Quantities of Interest
Several dimensionless thermal quantities of interest can be derived from the solution to the governing differential equations.
(1) Local mean temperature of the suspension fluid at a given cross section of the tube (2) Mean liquid-phase or melted fraction of the PCM particles at a given cross section of the tube (3) Heat transfer on the inner wall surface The energy dissipated outward from the tube is shown in Figure 1b. The dimensionless local heat flux at the inner tube wall can be obtained as follows: Integrating the heat flux along the heat exchange section yields the heat transfer through the inner wall surface: Neglecting the axial heat conduction, the following equation can be deduced from the equation for the fluid energy: is the dimensionless total heat transfer rate to the tube wall from fluid convection heat transfer that has two components: The axial dimensionless heat transfer rate in the fluid is as follows: (4) Heat transfer for the outer wall The total heat transfer dissipated through the outer tube wall of the heat exchange section is as follows: where (5) Axial heat conduction in the tube wall (6) Nusselt number (Nu) of the inner wall surface (7) Nusselt number of the outer wall surface where the overall heat transfer coefficient is

Numerical Method
Spatial discretization of the differential equations was performed on a uniform mesh for the fluid and solid domains, incorporating the second-order central difference scheme for the diffusion terms and the Quadratic Upstream Interpolation for Convective Kinematics (QUICK) scheme [17] for the convective terms. Temporal derivative terms in the differential equations are treated using an implicit approximation with second-order accuracy. The algebraic discretization equations were solved through the line-by-line application of the tri-diagonal algorithm. Steady-state solutions to the governing equations were obtained numerically by a pseudo-transient approach. The temperature field was considered convergent when the maximum variation between the value computed at the present and previous iterations was less than 10 −5 . All computations were performed with double-precision arithmetic. Moreover, the global energy balance over the computational domain was checked for each simulation and found to always be within 0.1%.
In this paper, the grid convergence index (GCI) supplemented by the extrapolated relative error is used to search for the grid point and determine the data error. The grid point test is performed for pure water and a fluid with 20% PCM particles by volume. The various errors in the local Nu of the inner and outer wall surfaces, the local heat fluxes, the average cross-sectional solid-liquid volumetric fraction, and the average cross-sectional temperature are calculated using different numbers of grid points to determine the number of grid points. Figure 2 shows one of the obtained results. For a grid point test, obtained with c v = 0, Bi o = 1, Pe = 1000, Ste = 0.1, t w = 0, k * w f = 10, and θ in = 0, three cases have been tested, and the recommended numbers of points are 501 in the axial direction and 50 in the radial direction. For c v = 0.2, Bi o = 1, Pe = 1000, Ste = 0.1, t w = 0, and θ in = 0, four cases have been tested, and the recommended numbers of grid points are 2001 in the axial direction and 250 in the radial direction, which are used in this paper.
To validate the present formulation and numerical method, simulations were performed for selected cases of the forced convection of PCM suspension in a circular tube for comparison with the results that were presented in previous studies [10]. The results showed that there was no significant difference; thus, the reliability of the simulation results was confirmed.
points to determine the number of grid points. Figure 2 shows one of the obtained results. For a grid point test, obtained with 0, 1, 1000, 0.1, 0, * =10, and 0, three cases have been tested, and the recommended numbers of points are 501 in the axial direction and 50 in the radial direction. For 0.2, 1, 1000, 0.1, 0, and 0, four cases have been tested, and the recommended numbers of grid points are 2001 in the axial direction and 250 in the radial direction, which are used in this paper.

Results and Discussion
Numerical results have been obtained for the parameters pertinent to the present problem in the following ranges: the dimensionless tube wall thickness, t w = 0-0.5, and the wall-to-fluid thermal conductivity ratio, k * w f = 0.1-10. Other parameters were kept constant (c v (0.1), Bi o (1), Ste (0.1), and Pe (1000)).

Tube Wall and Fluid Temperature
The wall conduction effect is illustrated in Figure 3 Figure 3a shows that there is little change in the respective temperatures and ξ b between t w = 0.2 and 0.5. Thus, t w does not directly affect the heat transfer performance of the heat exchange section. By contrast, the temperature distributions are quite sensitive to changes in k * w f . The heat transfer enhancement by the suspension fluid at fixed Ste or c v is markedly impeded by a lower k * w f , resulting in considerably higher inner and outer wall temperatures.
In Figure 3b, the post-cooling effect is increasingly discernible with increasing k * w f . For k * w f = 10, the PCM completes heat exchange when the dimensionless parameter x + + L + u /(D i Pe) = 0.16. The length required for the entire PCM to dissipate heat is 91% of the length of the heat exchange section. The axial propagation of conduction heat transfer along the tube wall results in an over twofold decrease in the inner and outer wall temperatures at the exit of the cooled section. The PCM suspension microparticles still undergo an exothermic reaction in the downstream adiabatic section, as shown by ξ b (the yellow line) in Figure 3. Thus, the temperature increases the fluid temperature and inner and outer walls, as reflected by the change in θ w,i , as illustrated by the yellow region.

Heat Transfer Characteristics of Inner and Outer Tube Wall Surfaces
Figures 4 and 5 show the local Nu and local heat flux of the inner and outer wall surfaces along the dimensionless distance / , respectively. The wall surface is the medium for heat exchange with the external environment. Within the considered parameter range, has a nonsignificant effect on the heat transfer compared to * . In Figures 4a and 5a, for * of 10, the heat transfer and Nu are affected only by at the entrance of the heat exchange region and slightly affected at the inner wall. For k * w f = 1, the reaction in the PCM suspension microparticles remains exothermic in the downstream adiabatic section, as shown by the red line for the change in ξ b in Figure 3. Even at the terminal exit of the downstream adiabatic section, the PCM still contains unsolidified liquid, which increases the fluid temperature and thus the temperature of the inner and outer walls, as shown by the red region in Figure 3 for the change in θ w,i . Additionally, when the working fluid exits the adiabatic section, the temperatures are not yet steady (the steadiness corresponds to no further increase in the temperature). Thus, the temperature at the downstream outlet and its steadiness depend on whether the PCM latent heat is completely released in the heat exchange section. Figures 4 and 5 show the local Nu and local heat flux of the inner and outer wall surfaces along the dimensionless distance x + + L + u /(D i Pe), respectively. The wall surface is the medium for heat exchange with the external environment. Within the considered parameter range, t w has a nonsignificant effect on the heat transfer compared to k * w f . In Figures 4a and 5a, for k * w f of 10, the heat transfer and Nu are affected only by t w at the entrance of the heat exchange region and slightly affected at the inner wall.   However, increasing k * w f markedly affects the extent of the postcooling regions. Substantial heat transfer can occur along the indirectly cooled regions of the inner and outer walls, in accordance with the foregoing postcooling effect from combined axial conduction along the tube wall and fluid stream. In Figure 4b, the post-cooling region significantly increases in size with k * w f up to k * w f = 10. For k * w f = 1, the working fluid still dissipates heat after flowing out of the heat exchange area, and the heat exchange is completely terminated only when

Heat Transfer Characteristics of Inner and Outer Tube Wall Surfaces
x + +L + u D i Pe = 0.24. For k * w f = 10, the working fluid continues to dissipate heat even in the adiabatic section, confirming the significant effect of wall conduction on the heat dissipation capacity. Moreover, an iso-flux heat transfer limit along the cooled region can be increasingly achieved by decreasing k * w f , as demonstrated by the curves for k * w f = 0.1.

Analysis of Various Heat Transfer Contributions
Further insight into the heat transport mechanism of the suspension fluid flowing in the tube can be gained by analyzing the following quantities: the dimensionless outer wall heat transfer ( , * ),

Analysis of Various Heat Transfer Contributions
Further insight into the heat transport mechanism of the suspension fluid flowing in the tube can be gained by analyzing the following quantities: the dimensionless outer wall heat transfer (q * w,o ), the dimensionless inner wall heat transfer over the cooled section ( q * w,i L h ), the dimensionless sensible heat transfer ( q * w,i L h , sen ), the dimensionless latent heat transfer ( q * w,i L h , lat ), and the dimensionless wall conduction heat transfer ((q * w ) cond ). The respective plots for these heat transfer contributions are analyzed to identify the general features of the competing mechanisms of sensible and latent heat transfer.
First, for the considered parameters, the axial fluid conduction over the cooled section accounts for approximately 0.2% of the total heat transport in tube flow (data not shown). In Figure 6a, the individual heat transfer contributions do not change significantly with t w but increase with increasing k * w f , albeit at different increments. At low k * w f values of 0.1 and 1, the latent heat transfer is the dominant mechanism for convective transport in the suspension flow, and this contribution decreases considerably as k * w f increases. For k * w f = 0.1, each of the contributions to the heat transfer rate is quite low. When k * w f = 1, the latent heat transfer ( q *

Conclusions
In this paper, numerical results have been obtained to investigate the conjugate heat transfer of a PCM suspension in a circular tube under an external cooling convection boundary. Of particular emphasis, in this study are the wall conduction effects on the heat transfer behavior. The parameters pertinent to the present problem are considered in the following ranges: the dimensionless tube wall thickness, = 0-0.5, and the wall-to-fluid thermal conductivity ratio, * = 0.1-10. Other parameters are kept constant ( (0.1), (1), Ste (0.1), and Pe (1000)). The following conclusions are drawn from the results of the study.
(1) The most significant of the considered parameters appears to be * , which affects the

Conclusions
In this paper, numerical results have been obtained to investigate the conjugate heat transfer of a PCM suspension in a circular tube under an external cooling convection boundary. Of particular emphasis, in this study are the wall conduction effects on the heat transfer behavior. The parameters pertinent to the present problem are considered in the following ranges: the dimensionless tube wall thickness, t w = 0-0.5, and the wall-to-fluid thermal conductivity ratio, k * w f = 0.1-10. Other parameters are kept constant (c v (0.1), Bi o (1), Ste (0.1), and Pe (1000)). The following conclusions are drawn from the results of the study.
(1) The most significant of the considered parameters appears to be k * w f , which affects the temperature of the outer tube wall, the temperature of the inner tube wall, the average temperature of the fluid, and the volumetric liquid fraction in a PCM particle. The various contributions to the heat transfer and t w have a comparatively minor effect on the aforementioned variables.
(2) For high k * w f , the tube wall and fluid are significantly affected by postcooling in the downstream adiabatic section. Therefore, a substantial increase in k * w f creates an indirectly cooled region downstream of the heat exchange section where the PCM particles have completed solidification. Thus, there is a significant increase in the latent heat released in the solidification process and thus the efficacy of wall temperature control over the cooled section or even the heated section. Biot number at the outer surface of the tube (=h ∞ r + o /k w ) Bi p Biot number of a PCM particle (=Ur + p /k p ) C p specific heat (kJ/kg K) c * p,b f specific heat ratio (=c p,b /c p, f ) c * p,p f specific heat ratio (=c p,p /c p, f ) c * p,w f specific heat ratio (=c p,w /c p, f ) c v volumetric fraction of PCM particles d p dimensionless diameter of a PCM particle (=d + p /r + i ) d + p diameter of a PCM particle (m) Peclet number (=2r Nusselt number (=h 2r + i /k b ) q heat transfer rate (W) q local heat flux (W/m 2 ) q * local dimensionless heat flux (=q /q o ) q * dimensionless heat transfer rate (=q/q o ) r + radial coordinate r + p radius of a PCM particle (m) r dimensionless coordinate (=r + /r + i ) Ste modified Stefan number (=c p,b (T m − T ∞ )/h ls ). t time (sec) t w dimensionless tube wall thickness (=t + w /r + i ) t + w tube wall thickness (=r + o -r + i ) (m) T temperature (K) u dimensionless axial velocity (=u + /u + m ) u + axial velocity (m/s) U overall heat transfer coefficient (W/m 2 K) x dimensionless axial coordinate (=x + / r + i Pe ) outer surface of the tube p particle pf particle-to-fluid ratio sen sensible heat u upstream w tube wall wf wall-to-fluid ratio ∞ ambient Superscripts surface-averaged quantity * ratio of quantities + dimensional quantity