Localized heating element distribution in composite metal foam‐phase change material: Fourier's law and creeping flow effects

A numerical parametric study is presented of a domestic thermal storage heat exchanger to explore the effect of highly localized positive temperature coefficient cylindrical heating elements in a phase change material (PCM) with conductive enhancement by open‐pore metal foam. By using 90 L of commercially available Rubitherm RT70HC wax, 5.7 kWh of thermal energy is captured by the unit. The discharge is via a central convective air channel. The constant low‐temperature heating elements are inherently safe for combustible PCM. The heat distribution by Fourier's law and the creeping flow is investigated using the local thermal equilibrium assumption between the PCM and metal foam. Heating element position, diameter, and temperature are varied to optimize charge time and exit air temperature. Two heating elements of 1 cm diameter and constant temperature of 90°C produce a suitable performance for overnight store charging of 7.23 hours. Discharge via the air channel provides an average temperature of the output air over 30°C. The results indicated that the PCM inside metal foam almost follows Fourier's law. The creeping flow of molten PCM inside the pores of the porous medium (free convection heat effect) has an inconsiderable influence on heat transfer in the domain.


| INTRODUCTION
Energy demand has been improved due to the use of more energy by increasing the world population. Thermal energy storage (TES) can save energy by connecting the disparity between the provided energy and demand by better utilizing renewable power resources. 1 Latent heat thermal energy storage systems (LHTES) operate at low temperatures and have a high energy mass. 2 Latent heat storage in PCMs has found application in several areas. 3 The application of solar energy for solar water thermal units, solar air thermal units, solar cookers, and high performance buildings are summarized in the review by Sherma et al. 4 For building applications, PCMs have been classified as paraffin wax, hydrated salts, fatty acids, and eutectics. 5 The drawback is that all these materials have a weak thermal conductivity, which limits the thermal efficiency. Enhancing thermal conductivity could enhance the performance of the unit by shorting the charging time. 6 Different methods have been developed to improve PCM characteristics, including fins deployment, 7,8 incorporating design features in air conditioning ducts, 9 heat pipes, 10 multi-tube usages in shell and tube heat exchangers, 11 encapsulated PCMs in microencapsulation 12 in larger-scale capsules, 13 molten salt capsules, 14 direct contact heat transfer enhancement by convective transport. 15 Some other improvement methods are enlarging the temperature variation between the heating source and the PCMs, 16 localized heat generation inside the PCM, 17 using nanomaterials in the PCM, 18 and deployment of the extended surface heat transfer by inserting fibrous materials such as carbon fiber brushes and carbon nanofibre, 19,20 use of expanded graphite foam, 21,22 and the use of metal foams. 23,24 The high thermal conductive porous medium approach significantly improves PCM-based heat exchangers' thermal performance as the porous medium has great thermal conductivity, high porosity, and high thermal surface area with negligible effect on the storage capacity 25 using metal or graphite foams. Different investigations display the enhancement capability of composite porous medium/PCM systems. Yang et al. 26 experimentally investigated the influence of porous material in an LHTES system and found that the phase change duration can be reduced to one-third of that for the case without the foam. Abujas et al. 27 studied the influence of graphite foam PCM-based application numerically. They found that the foam with approximately 60% accessible porosity reduced the charging time by a factor of up to 100. The charging process rate of paraffin with high PPI (pore per inch) aluminum foam was observed by Zhu et al. 28 using a numerical method, which showed that the porous medium improved the thermal efficiency by up to 83%. Because of the limited features of experimental studies such as cost and time considerations, the range of the input parameters, and the limitation of the post-processing, the development of numerical models for PCM-based systems has resulted in different commercial software such as ESP-r, 29 TRNSYS, 30 and ANSYS Fluent. 31 Mahdi et al. 32 examined the charging procedure of nanoparticle-enhanced PCM in a porous foam triplex-tube using Fluent, demonstrating the significant effect of porous metal foam. PCM infiltrated microcellular metal porous medium was evaluated numerically by Sundarram and Li 33 with extreme heat production and minimal cooling terms. They found that a larger pore dimension affects a higher phase change rate, and the effective thermal conductivity became twice by decreasing the pore dimension from 0.1 to 0.025 mm. Hossain et al. 34 simulated the thermal efficiency of a rectangular TES unit filled with a porous medium drenched by a nano-PCM (C 6 H 12 + CuO nanoparticle). They used the outcomes to produce a dimensionless consideration of the main physical aspects affecting the situation.
Regarding the PCM-based TES applications in buildings, the review by Zhou et al. 35 inspect the state of PCM impregnation of building materials to save waste energy. Souayfane et al. 36 explored the application of building cooling utilizing PCM. Relating to the building's context, PCMs could be placed in compensation for the absence of thermal material in the building shield, which would be recharged after usage. In the cooling application, night ventilation is an ordinary technique for recharging the PCM, absorbing the heat stored during the daytime as defined in the practical study of Zalba et al. 37 Incorporating PCM in buildings for controlling the thermal benefits has been commonly investigated previously in both active and passive storage systems. 38,39 There are some limitations in the PCM-based applications in the buildings, including the small area, construction, architecture, location, and weather. [40][41][42] As an example of the active system, ice storage needs a loop of fluid to charge/discharge the storage tank, filling with the PCM, and passing the cooling load to the water loop. 43 A directly attached PCM-water system was studied by Martin et al. 44 They reported the storage power of 30-80 kW/m 3 . Different experimental and numerical studies have been demonstrating the heating application used for solar-PCM storage systems. Passive thermal energy storage uses the energy density of the building to save heat, 45 and the advanced format of this system is to use the PCM in the shell of the building. 35 PCM can be used in various building components such as window glass, 46 concrete, 47 gypsum boards, 48 wallboards, 49 roof, 50 floors, 51 and walls, 52,53 making the use of the surface area to save space. 54 The benefits of PCM layers on buildings materials have also been studied to enhance thermal efficiency. 55 This study investigates the design factors of a wallmountable heat energy storage unit for domestic space heating applications incorporating PCM with metal foam under the Fourier Law and creeping flow effects. The model assumes the melting process using positive thermal coefficient (PTC) type electrical heating components for charging the PCM during off-peak hours using the local thermal equilibrium assumption between the PCM and metal foam. PTC heating elements provide a constant temperature heat source preventing the dangerous presence of high-temperature fluid as a space heater. Moreover, the use of rod heating elements helps overcome the PCM's shrinkage problem through the phase change procedure. After the heat is stored, the thermal energy is subsequently delivered to the air for space heating purposes aiming to provide uniform temperature due to metal foam. The heating elements' characteristics are studied comprehensively during the charging process to optimize the melting process during the required melting hours (8 hours), including the diameter, number, and surface temperature. The subsequent discharging time to the air and the stability and uniformity of the output air temperature are also investigated during the discharge process. The study of metal foam-PCM-based energy storage using localized cylindrical heating elements has not been investigated thus far, which is included comprehensively in this paper toward an innovative application for energy consumption reduction and heat decarbonization.

| PROBLEM DESCRIPTION
Two geometrical arrangements are considered. Figure 1 shows schematic diagrams of the numerical volumes modeled for the proposed units using two and four constant temperature heating elements with 25 cm long. The computational domain is considered symmetric in two planes in x and z directions, as indicated in Figure 1 (dimensions 15 × 30 × 50 cm), which represents onefourth of the physical thermal-storage unit. A metal foam-PCM is inserted in same size rectangular compounds adjacent to adiabatic walls airflow path within the current unit. The metal foam is considered copper metal with a porosity of 0.95 and a pore density of 30 PPI. RT70HC (RUBITHERM) is applied as the PCM material, an organic PCM with a high energy density and its thermophysical specificities are shown in Table 1 along with the properties of copper as the metal foam material. 56 For the system using two heating elements, only half of the heating element is modeled, as shown in Figure 1A. The placement of the heating elements in the center of the PCM is intended to solve the physical problem of shrinkage of PCM when it solidifies (although the numerical model does not represent the volume change of the PCM) and to cause a uniform distribution of thermal energy in the domain. The whole unit is primarily at a temperature of 22 C. For the molten PCM drive in the system, the boundary layer develops as the PCM melts. Thus, a creeping flow of the molten PCM inside the pores of metal foam can be applied, and the convection process influences the fluid dynamic model. For the ideal transmission of heat from the PCM during discharge, a 1 mm copper layer is placed between the HTF path and PCM enclosure. The air-flow is set at 40 g/s in discharge. During the solidification process, the cold air moves through the tube to gain thermal energy from the stored thermal energy inside the PCM. In other words, the air enters the channel from the bottom and then, after gaining heat, exits from the top. The width of the air-flow path is 2 cm. The heating element's diameter and temperature are varied throughout the melting process to approach the ideal outcome of charging between 7 and 8 hours and discharging in 16 hours. In the full storage unit of four repeated volumes shown in Figure 1, the anticipated thermal storage is approximately 5.7 kWh. The mass is approximately 79 kg, and the total airflow rate is 0.04 kg/s. Note that as shown in Figure 1B, for a system with two heating elements, only half of the heating element is modeled due to the symmetrical condition of the domain.

| MATHEMATICAL MODELING
The enthalpy porosity method is applied for the phase change model, and the heat transfer model in the metal foam is based on the thermal equilibrium method. The governing equations in a 3D Cartesian coordinate system, assuming the incompressible Newtonian fluid flow of liquid PCM, ignoring viscous dissipation and a homogeneous and isotropic open-cell porous medium, are given as follows 31,57 : where A m is considered 10 5 . 58 The Boussinesq approximation is added to the momentum equation for the buoyancy effect. The reference temperature is considered as the melting point temperature, which is the average of solidus and liquidus temperatures. k e is the effective thermal conductivity of the combination of porous medium and PCM together, and given as 59 : λ (in Equation (2)) is the liquid fraction, and it can be found through: where ΔH is the PCM fractional latent heat (0 (solid) < ΔH < 1 (liquid)).
For the airflow simulation, the effects of metal foam and PCM are excluded from the governing formula. The conventional Navier-Stokes equations are solved using the continuity and energy equations for laminar fluid flow in a duct. 60 The conductive process is the main way of heat transfer in this unit. The heat flux following from thermal conduction is directly proportional to the temperature gradient's quantity and opposite to it in sign (Fourier law).
where the vector _ q x is the heat flux (W/m 2 ) in the positive x-direction, dT/dx is the (negative) temperature gradient (K/m) in the direction of heat flow. 61 Due to the presence of a porous medium, the molten PCM flow is presented by the creeping flow, which has a very low Reynolds number Re < 1. 62 The momentum is transported by viscous diffusion (not by convection). This means there is no wake on the downstream side of an object; pressure scales on viscosity, not on kinetic energy; when they occur, eddies are as likely upstream as downstream of a blunt body. 62 Here, the low Reynolds number means the low velocity of the molten PCM.

| NUMERICAL PROCESS AND VALIDATION
ANSYS-FLUENT solves the main equations using SIM-PLE pressure-velocity coupling together with PRESTO for momentum and QUICK for energy equation. The fluid motion in the porous medium is expected to be a creeping flow. The convergence values are set to 10 −4 for continuity and 10 −6 for the momentum and energy equations. Figure 2 shows the generated meshes for two and four heating elements. The number of cells for the computational domain using four and two heating elements shown in Figure 2 are 538 154 and 469 132, respectively. Denser meshes are also generated to determine mesh independence based on the PCM charging time, and the results showed no considerable effect. Therefore, these grids are selected for further investigations. For the time step size, two different sizes of 0.1 and 0.2 second are examined to find the independent outcomes from the size of the time-step, and the results are almost identical. Hence, the time step size of 0.2 second is chosen to have a lower computational time. Table 2 presented the charging time for the case using four heating elements for the different grid and time-step sizes for the independence analysis test. In this case, the heating elements' diameter is 2 cm with a constant temperature of 95 C. As shown, the results of different cell numbers, as well as time-step sizes, are almost similar.
The experimental works of Zhao et al. 63,64 and Liu et al. 65 were used to validate the numerical analysis. They studied a rectangular LHS system with a heated lower side utilizing a metal foam (Cu) embedded in the PCM (RT-58) employing the thermal equilibrium model. For the studied geometry, the heat loss is assumed for all the sides of the enclosure except for the base side. Figure 3 illustrates the variation of temperature at the distance of 8 mm from the base side applying both thermal equilibrium and nonequilibrium methods. The experimental work of Zhao et al. 62 showed agreement with the maximum percentage error of 3% with practical results. It worth to mention that, the thermal equilibrium form is used in the current work because of the limited ability and costly arithmetic process of the thermal nonequilibrium form in 3D simulations, also, the tiny difference between applied the mentioned two models presented in the authors' previous work. 59 Figures 4 and 5 show the contours of liquid fraction and temperature profile, respectively, in the melting processes in 2-hour time periods for the system using four heating elements (2 cm diameter) and fixed temperature of 95 C. When the PCM temperature achieves the solidification point temperature (69 C), melting begins, and thermal power is collected in the PCM as latent heat. Laminar flow with a very low Reynolds number is created during the charging of the PCM showing the generation of creeping movement of molten PCM in the pores of the porous medium. The length of the element is 25 cm (5 cm less than the PCM height), and the lower part of the PCM melts later than the upper part, as shown in 900 seconds. The effect of density variation and gravity produces a cone shape around the heating element. The conductive metal foam causes various heat transfer rate at the top and bottom of the heating element to be low since the angle of the temperature cone boundary is not sharp. After 3600 seconds, more than 85% of the PCM is converted to liquid, and by 5400 seconds, most of the domain reaches the equilibrium temperature. Note that in Figures 4 and 5, one-quarter of the domain is displayed due to the symmetrical condition of the problem. Figure 6A shows the PCM's liquid development for the diameter of the heating element 0.5, 1, and 2 cm at a constant element temperature of 95 C using four heating elements. Melting time reduces from 3.18 to 1.7 hours (reduces 45%) as the elements' diameter increases from 0.5 to 2 cm. PCM temperature increases and reaches the melting temperature presented in Figure 6B. These observations are directly related to increasing the surface area. As shown, the difference between the liquid fraction of cases with the element diameter of 2 and 1 cm almost equals the liquid part difference for the cases with an element diameter of 1 and 0.5 cm. This means that the thermal behavior inside the composite follows the Fourier's law, showing that the behavior of the composite block is more like a solid block for heat storage rather than a liquid block and the effect of creeping flow of molten PCM inside the pores of metal foam, which is the reason for the difference of the rate of liquid part increment between the various sizes of heating elements is small. Note that the velocity of molten PCM is in the order of 0.01 mm/s results in the Reynolds number of less than 1 in the PCM pores, which shows the creeping flow in the system. Figure 7 illustrates the liquid phase's contours and the temperature distribution for two different diameters (0.5 and 1 cm) after 1800 seconds. Liquid PCM volume is larger for the 1 cm element than the 0.5 cm case, again because of the heater's larger surface area. Liquid gathers  Figure 6. It is surprising at 1800 seconds (shown in Figure 6) the doubling of diameter (4 times surface area and volume) results in half the rate increase of quadrupling the diameter (16 times surface and volume). This indicates that it is effective to use the 1 cm diameter from a mass (and presumably cost) of heating element consideration.

| Element temperature effects
The higher temperature of the four heating elements, 2 cm diameter case, accelerates phase change, as presented in Figure 8. The total PCM converts to liquid in 1.7 hours for the element temperature of 95 C, which increases to 3.23 hours for the element temperature of 85 C. As shown, similar to the effect of heating element diameter, the difference between the liquid fraction of cases with the element temperatures of 95 and 90 C almost equals the liquid fraction difference for 90 and 85 C. This means that the heat transfer inside the composite follows the Fourier's law, showing that the behavior of the composite block is more like a solid block for heat storage rather than a liquid block and the effect of the creeping flow of molten PCM inside the pores of metal foam is small. Figure 9 illustrates the molten development and the local temperature in the system at 1800 seconds for heating elements at 85 C and 95 C. Inspecting the detail at 1800 seconds in Figure 8 shows that changing the temperature by 5 C produces a 0.08 change of volume fraction in a linear relationship in these conditions. The charging time reduces by 47.5% using the heating elements with a 95 C compared with that using 85 C.

| The number of elements effect
Based on the conduction process of Fourier's Law, there are a couple of options to improve thermal performance in conduction regimes. Presuming a fixed temperature difference, the thermal conduction enhances due to either the mean surface contributing in a thermal behavior or the considerable value of the effective thermal conductivity. Figure 10 shows the effect of using four 1 cm diameter heating elements at 95 C compared to two elements, as expected since the surface area of the four elements is double, and the heat distributes in four different points rather than two. Four elements cause the PCM to melt in 3.18 hours, whereas two elements melt in 8.1 hours, a factor of 2.55 acceleration. After 3.18 hours, when the PCM melts completely using four heating elements with a mean temperature of 90.5 C, while it is 73.9 C in using two elements and the melt fraction is 0.5. Based on the conduction process of Fourier's Law, heat moves from hot to cold. Figure 11 shows the liquid fraction and temperature distribution of half of the system using two 0.5 and 1 cm heating elements comparable with Figure 7 using four elements with a diameter of 1 cm. Using the four elements case, higher heat transferred from the heat source results in a higher amount of molten PCM, which is also illustrated in the temperature distribution pictures, showing that the red color (higher temperature) rings cover a wider area. Note that in Figure 11, half of the domain is illustrated for each figure due to placing the heating element in the middle of the PCM domain; however, the computational domain is one-quarter of the unit, described in the system description. Figure 12 shows the difference in heat flux rate of the heating element for different numbers and the diameter of the heating element. When the heating element is turned on in the simulation, the heat flux changes from zero to the maximum power output to generate the element's constant defined temperature. After that, thermal diffusion in the PCM reduces the rate of heat flux. With two heating elements, the maximum heat flux for diameters of 1 and 0.5 cm is 1107 and 720 W/m 2 . Four heating elements with a diameter of 0.5 cm produce a maximum power output of 1358 W/m 2 , which is 1.89 times higher than using two 0.5 cm heating elements. The distribution of the heating elements in the system helps by inducing a uniform heat distribution all around the PCM, increasing the heat transfer in the system. The energy output rate is similar for different diameters of the elements at a constant element number after 1800 seconds. Increasing the number of elements caused the rate of energy output to always be higher till the end of the charging procedure, which is built on the estimation of thin layer flow (film movement) and creeping (large-viscous) flow.

| Summary of the charging process
The aim of this work is to optimize the arrangement of the rod heating elements for charging in less than 8 hours. Table 3 summarizes the charging times and the average rate of charging during the melting process. As presented, applying double heating bars (1 cm diameter) with a temperature of 95 C, the charging process takes almost 7.23 hours, which is desirable in this study. Therefore, this system is chosen for further analysis of the discharging process. By using a higher temperature of the heating element, more energy can be stored by sensible heat, and lower charging time can be achieved. Since the maximum PCM operating temperature is 100 C, lower constant temperatures are selected for the heating elements.

| Discharging process
Once the power cut off on the system after 8 hours, discharging (solidification process) begins by warming the F I G U R E 1 0 Development of the PCM (A) liquid fraction and (B) average temperature for different heating element diameters at elements temperature of (95 C) for the heating element diameter of 1 cm air inside the central tube. The cold HTF enters the tube from the bottom of the heat exchanger and then, after gaining heat from the stored thermal energy inside the PCM, which results in a higher temperature of the HTF, exits to the channel from the top of the heat exchanger. Figures 13 and 14 illustrate how the molten part shrinks and the temperature distributes in the system using two heating elements with 1 cm diameter and 95 C temperature at the unit's mid-section. The air channel is shown in the middle of each contour while PCM is placed at the sides. The average temperature of the PCM after 8 hours of charging is almost 90 C before starting the discharging process. The contours of the molten PCM and temperature are almost unchanged in the unit's length direction due to the discharge mechanism. 66 Since the air enters the system from the lower side, solidification starts from the bottom since air has a lower temperature. As shown in the temperature distribution, an almost uniform temperature occurs in the PCM domain by the attendance of high conductivity metal foam. In the solidification process, after 16 hours, the liquid fraction drops to approximately 0.05 (as shown in Figure 13); therefore, 5% of the stored energy remains unused. Figure 15 illustrates the variation of PCM liquid fraction and average PCM and exit air temperatures F I G U R E 1 2 The variation of heat flux rate generated from the heating elements for different diameters and numbers of the elements in the discharging case (from 8 to 24). Because of the existence of the porous medium, the average temperature of PCM does not vary significantly during the solidification after a sharp decrease in the sensible heat section from the temperature of 90 C (initial temperature in the discharging) to the liquidus temperature (71 C). Furthermore, a uniform output temperature is achieved with an average temperature of almost 30.9 C. The mean output temperature reduces from almost 34.8 C to 32.4 C (start of solidification) and then reduces to 29.9 C after 24 hours.

| Impact of the metal foam
The unit without using the copper foam is also studied in this study to highlight the advantages of the presence of metal foams. The outcomes show that, during the melting process, all the PCM melts in 49 hours. However, this stage is achieved within 1.7 hours with metal foam with the same condition as shown in Figure 16A (taking into account the presence of the foam's solid content, the size of the system increasing to reimburse the lost volume of PCM to foam). The difference in the charging time illustrates the influence of large thermal conductivity starkly and enhanced thermal behavior surface area. Figure 16B illustrates the development of the mean PCM temperature with and without foam. The system reaches 90 C within 1.7 hours, while for the non-porous unit, it reaches 92 C without foam in 49 hours (once the molten of the whole PCM). The mean temperature in the case of PCM only is a little higher because, with foam, the temperature distributes more uniformly, and the PCM reaches the melting point uniformly. In contrast, the case without foam has a local temperature rise around the heating element only. Note that as mentioned for the melting, the composite unit follows Fourier's law, while the PCM individual case follows the free convection impact.

| CONCLUSION
A porous medium-PCM TES unit with energy exchanging to air was developed numerically for domestic space heating using localized rod constant temperature heating elements for charging. The local thermal equilibrium assumption is employed between the PCM and porous medium. During the discharge, thermal energy delivers to air via a tube in the mid of the system. There is a significant improvement of the performance of heat transfer in the presence of high conductivity porous medium so that for the case of using four heating rods (diameter of 2 cm) with the temperature of 95 C, the charging time reduces to 1.7 hours for the porous-PCM case from almost 49 hours for the PCM only case. By utilizing two heating rods (diameter of 1 cm) with a temperature of 90 C, the melting time is almost 7.2 hours, which is sufficient for dwelling application. During the solidification procedure, the uniform output temperature of 31 C is approximately achieved, showing the significant advantage of the proposed system for domestic space heating while delivering heat to the room. The results indicated that the PCM inside metal foam acts similar to a composite solid block and almost follows Fourier's law. The creeping flow of molten PCM inside the porous medium's pores (free convection effect) has a tiny impact on increasing the heat transfer in the domain. The outcomes of this work deliver rules for the application of PCM energy storage heat exchangers in domestic space heating.