Numerical modelling of phase change material melting process embedded in porous media: Effect of heat storage size

The aim of this paper is to study the influence of enclosure size in latent heat thermal energy storage systems embedded in a porous medium for domestic usage of latent heat thermal energy storage heat exchangers. A 2-D rectangular enclosure is considered as the computational domain to study the heat transfer improvement for a phase change material embedded in a copper foam considering a constant heat flux from the bottom surface. Different dimensions of the composite system are examined compared with a system without a porous medium. The thermal non-equilibrium model with enthalpy-porosity method is employed for the effects of porous medium and phase change in the governing equations, respectively. The phase change material liquid fraction, temperature, velocity, stream lines and the rate of heat transfer are studied. The presence of a porous medium increases the heat transfer significantly, but the improvement in melting performance is strongly related to the system's dimensions. For the dimensions of 200 × 100 mm (W × H), the melting time of porous-phase change material with the porosity of 95% is reduced by 17% compared with phase change material-only system. For the same storage volume and total amount of thermal energy added, the melting time is lower for the system with a lower height, especially for the phase change material-only system due to a higher area of the input heat. The non-dimensional analysis results in curve-fitting correlations between the liquid fraction and Fo.Ste.Ra −0.02 for rectangular latent heat thermal energy storage systems for both phase change material-only and composite-phase change material systems within the parameter range of 1.16 < Ste < 37.13, 0 < Fo < 1.5, 2.9 × 104 < Ra < 9.5 × 108, 0 < L f  < 1 and 0 < Fo.Ste.Ra −0.02  < 0.57. Over a range of system's volume, heat flux and surface area of the input heat flux, the benefit of composite phase change material is variable and, in some cases, is negligible compared with the phase change material-only system.


Introduction
Latent heat thermal energy storage (LHTES) systems are used due to high capacity of heat storage, typically 5-14 times higher than sensible heat storage systems, with the added advantage of almost constant temperature during the solidification/melting process. 1,2 LHTES systems have been employed in both domestic and industrial sectors in order to reduce the energy demand. 3 For domestic heat exchangers, LHS can be employed in order to reduce the consumed energy as well as help peak-shaving in the buildings as the main sector of energy consumption. 3,4 Phase change materials (PCMs) are used in LHTES systems to store heat during the melting process and then release heat during the solidification process. The main challenge of efficient PCM is the long melting time due to low thermal conductivity and low thermal diffusivity within the bulk PCM; these disadvantages lead to limited use of LHTES systems. 5 Methods in the literature that help to overcome the weak characteristics of PCM include the use of a more complex configuration of heat exchangers, [6][7][8][9] use of extended metal surfaces, 10,11 use of encapsulated PCM [12][13][14] and use of high conductive nanoparticles, [15][16][17][18] and the most recent way is the use of a conductive porous medium. [19][20][21][22] The use of porous media in LHTES systems enhances the heat transfer rate inside the PCM. Attention has focused on investigating the use of porous media on different heat exchangers to enhance the thermal performance of the system. Py et al., 23 as one of the pioneer researchers, studied a composite paraffin-graphite matrix impregnated by the capillary forces. They presented that the equivalent thermal conductivity of the composite is ranged from 4 to 70 W/m K based on the percentage of paraffin in the matrix instead of 0.24 W/m K for the paraffin. Reducing the porosity of the matrix enhances the thermal conductivity of the composite. They showed a reduction in the solidification time and higher stability of the TES system with a composite structure. Mesalhy et al. 24 performed a numerical analysis on the melting process of PCM embedded in a high conductivity porous matrix using a thermal nonequilibrium model. There is a significant effect of the porous matrix on the rate of heat transfer and melting time. They claimed that although the melting rate increases by using low porosity medium, due to reducing the convection effect, the PCM storage with high porosity and high thermal conductivity is the best technique. Zhao et al. 25,26 studied on heat transfer enhancement in PCM embedded in a copper metal foam experimentally and numerically. They worked on the effect of porosity and pore density of a metal foam-PCM compared with a PCM-only system. They showed significant increase of heat conduction rate by using a metal foam. However, in comparison with the PCM-only system, the natural convection effect reduces significantly. Lower porosity and higher pore densities can increase the rate of heat transfer in the domain. They showed that the heat transfer increases 5-20 times in the solid phase zone and 3-10 times in the charging process with metal foam compared with the PCM-only system, with similar effect during discharge. Tay et al. 27 performed a CFD simulation using ANSYS FLUENT software accompanied with experimental validation for a complex PCM-based heat storage unit. A complicated and unique tube arrangement for heat transfer fluid was designed to charge and discharge the PCM. They showed that the capability of CFD code on simulating phase changes in a complex heat storage system. They presented the details of the CFD simulation which can accurately predict the behaviour of LHTES systems. Using CFD, Liu et al. 28 studied the melting of phase change in a porous medium in a shell and tube LHTES system in 2-D and 3-D cases. They claimed that the melting time can be reduced by more than seven times in comparison with PCM-only system. In the study of the pore size, they showed that, in general, the pore size has a very small effect on the melting process in the porosity of 95%. For very high pore densities (60 PPI), the natural convection is almost negligible due to high flow resistance in the melting process. Decreasing the pore density (from small sizes such as 10 PPI to the moderate sizes such as 30 PPI), higher performance is achieved due to a higher rate of heat transfer. Sciacovelli et al. 29 performed a CFD simulation on the phase-change process in a vertical shell-and-tube LHTES unit. They verified the CFD code using a 2D axisymmetric model with experimental data from the literature. They added copper nanoparticles to the PCM to increase the thermal conductivity and showed that the nano-enhanced PCM with 4% concentration of nanoparticles reduced the melting time by 15% and increase the heat flux rate by 16%. Gupta et al. 30 performed an experimental study and CFD simulation of an active cylindrical LHTES tank for domestic applications compared with water-filled tank due to the high heat capacity of suitable PCMs. They used hydrated salt as the PCM and used OpenFOAM software for CFD simulation. They showed accurate results of the CFD code compared with the experimental data. They reported the advantages of PCM heat storage unit than the water tank in the initial moments; however, the discharging heat was reduced due to the high viscosity of liquid PCM. Therefore, they recommended an improvement in heat exchanger design and the use of heat transfer enhancement methods. Hossain et al. 31 studied a 2-D rectangular enclosure heat storage system using both nanoparticles and porous medium. CuO was considered as the nanoparticles. They considered the source of thermal energy at the upper surface of the enclosure in order to reduce the effect of natural convection in the presence of a porous medium. Nanoparticles increase the thermal conductivity of PCMs, and the porous medium increases the rate of heat transfer through the PCM. They showed that the melting process improved more significantly by the presence of porous medium rather than nanoparticles and a faster rate of melting happens at lower porosity and higher volume fraction of nanoparticles. For the porosity of 85%, the melting time reduced from 3100 s for the nanoparticle concentration of 10% to 3440 s without nanoparticles. Zhang et al. 32 studied numerically different metal foams of copper and nickel in a latent heat solar energy storage system using molten salt PCM. They showed the reduced effect of natural convection in the presence of composite metal foam due to the flow resistance. They showed that due to the large difference between the thermal conductivity of metal foam and PCM, a considerable temperature gradient exists between the porous medium and PCM, and therefore the thermal non-equilibrium model should be correctly considered in the simulation. They showed that the copper metal foam has a better performance than nickel, having the maximum temperature difference of 6.8 C between the foam and PCM and the reduction in total melting time of 28.3%. Mahdi and Nsofar 33 studied a triplex-tube LHTES system using both Al 2 O 3 nanoparticles and copper foam to increase the melting time of the heat exchanger. They showed that simultaneous use of porous medium and nanoparticles can improve the melting of PCM significantly. They showed that by a reduction of PCM volume due to the presence of metal foam and nanoparticles, high porosity metal foam with a low volume fraction of nanoparticles is recommended. By employing 95% porosity copper metal foam, the melting time reduces from 162 min for the non-porous case min to 18 min; however, employing 5% nanoparticle just reduced the melting time to 130 min.
Based on the presented literature review, researchers employed metal foams to increase the rate of heat transfer in PCM and reduce the melting/solidification time compared with the PCM-only systems. However, the benefit of the composite depends on the volume of the PCM which is not mentioned in the literature. Therefore, in this paper, different sizes of the storage system and amounts of heat flux are examined for the copper foam LHTES system compared with the PCM-only case to comprehensively study the effect of the porous medium in relation to natural convection in the domain. The thermal non-equilibrium model is employed to model a composite metal foam with PCM in a 2-D rectangular enclosure, representative of domestic thermal storage units, with heat flux from the lower surface. The effect of solidliquid phase change is considered using the enthalpyporosity method. The liquid fraction, temperature and velocity distributions, streamlines and melting time are all studied for the composite PCM in comparison with the PCM-only system in this paper. The results of this paper provide guidelines for the better design of domestic thermal energy storage systems.

Mathematical description
In a PCM-only LHTES system, heat is transferred to the solid PCM by conduction and then spread out in the liquid PCM by convection and conduction. To increase the heat transfer inside the PCM, a porous foam is prepared, and the PCM is injected into the porous media. 23 The heat is more transferred by the high conductivity solid porous medium through the PCM. In the modelling of PCM, to consider the effect of phase change, the enthalpy-porosity technique is used where the porosity is set equal to the liquid fraction in each cell. Then, by applying Darcy's law for fluid flow in the porous medium and employing the Kozeny-Carman equation, the flow is modelled in the mushy zone where the liquid fraction is varied between 0 and 1. 34 In the presence of a porous medium, in addition to the pressure drop caused by solidified materials, a pressure loss is considered due to viscous and inertial losses in the momentum equation employing the Brinkman-Forchheimer equation. 35 For simulating heat transfer process in the presence of porous media, two thermal models can be applied: i.e. the thermal equilibrium and non-equilibrium models. In the equilibrium model, the liquid PCM and the porous media have the same temperature, and the Brinkman-Forchheimer equation is just added to the momentum equation as a new source term (fourth term on the right-hand side of equation (2) below). 36 However, in the thermal non-equilibrium model, the porous medium and fluid are not in thermal equilibrium, and the temperatures of the solid and the fluid are solved separately using a semi-heuristic approach based on a local thermal non-equilibrium model. This model is more accurate and is employed in this study 24,32 derived first by Vafai and Tien 37 and then modified by Cheng and Hsu. 38 .The following assumptions are also assumed in the present numerical investigation: 28 Based on these assumptions, the set of governing equations for the Brinkman-Forchheimer-extended Darcy model can be formulated as follows: 32

Continuity equation
The continuity equation is defined as follows based on the mentioned assumptions 35

Momentum equation
The momentum equations in x and y-directions are given as follows: . Thermal non-equilibrium model 35 . Thermal equilibrium model 36 On the right side of the momentum equation, the second terms belong to the viscous resistance, and the third term is the Kozeny-Carman equation, where A m is assumed equal to 10 5 kg/m 3 s as recommended in several studies. [39][40][41][42][43][44] The fourth and fifth terms account for the extension of Darcy's law to explain the non-Darcy effects. Due to the direction of gravitational acceleration in the negative y-direction, the Boussinesq approximation is added as the sixth term by f g" T À T ref À Á . 33 Note that the difference between the momentum equations for the non-equilibrium model compared with the equilibrium one is the presence of ð 1 " 2 Þ in the first term on the left and right sides of the equations and also ð 1 " 2 Þ in the second term on the left-hand side of the momentum equation.
The first term related to the presence of the porous medium is the viscous loss term, where K is calculated based on the empirical equation of Calmidi and Mahajan 45 which was obtained from experiments as follows 46 The second term is an inertial loss term, where C is determined as 46 and d l is the ligament or cell diameter given as 33 The characteristics of the porous medium are defined by three parameters including the ", ! and d l , and by knowing two of these parameters, the other one can be calculated from equation (6), where d p is given as 33 where PPI means part per inch.

Energy equation
. Thermal non-equilibrium model: The energy equations for the solid and fluid are given as 35 For the PCM For the porous medium In the energy equation, l is defined based of the fluid, solidus and liquidus temperatures as 34 ð9Þ and the total enthalpy is the summation of sensible and latent heat given as where h is defined as and ÁH is the fractional latent heat of the PCM that may vary between zero for solid and L (latent heat of fusion) for liquid.
In the thermal non-equilibrium model, the effective thermal conductivity of the fluid and solid should be determined and employed independently. There are different theoretical models in the literature based on the porous characteristics, the thermal conductivity of solid and fluid. 45,[47][48][49][50] The models are classified based on the unit cell. The 'tetrakaidecahedron' cell model of Boomsma and Poulikakos is employed in this study which is first introduced in 2001 48 and then corrected in 2011. 51 In this model, the general form of the effective thermal conductivity is defined as follows 48,51 where where e ¼ 0:16 and Then the effective thermal conductivity of fluid and solid is calculated from equation (12) as Due to lack of an accurate model to describe the local heat transfer between the solid porous medium and the fluid, one simplification is to calculate it based on theoretical studies of other geometries or fitting methods from numerical or experimental results. [52][53][54] The pores are usually treated as another geometry such as cylinders, spheres or even vertical plate. 53,[55][56][57] All the equations can be found in Zhang et al. 32 Among different methods, the porous structure is usually considered as cylinders, and the laminar flow of liquid PCM in porous structure is considered similar to the flow around a cylinder. 16,28,32 Therefore, the laminar fluid flow in the pores is assumed as the flow around a cylinder where the following empirical correlation is used given as 53 Re and A sf is given as 16 . Thermal equilibrium model: The energy equation is given as follows which is the same for the fluid and porous medium 36 where k e is calculated as the volume average of the thermal conductivities of porous matrix material and PCM is given as 28 Note that for the PCM-only system, in the governing equations, the porosity is equal to 1, and the source terms in the momentum equations due to the presence of the porous medium are omitted.

Problem description
The important effect of the presence of a high conductivity porous medium in the LHTES systems is to reduce the charging time. Therefore, in this paper, to study the effect of presenting a porous medium, different rectangular geometries with various areas are considered for both porous-PCM and PCM-only with various heat fluxes to see the benefits of the system.

Geometry and boundary conditions
The schematic of the physical domain studied in this paper is illustrated in Figure 1. In the system, porous-PCM/PCM-only is embedded in a rectangular enclosure with adiabatic walls, and a heat flux is implemented at the bottom of the system. Different widths and heights are examined ranging from 12.5 to 200 mm to completely investigate the effect of a porous medium in the various geometries using PCM. The material of the container as well as the porous medium is considered to be copper with a thermal conductivity of 387.6 W/m K. The heat flux is varied from 800 to 6400 W/m 2 regarding the size of the storage system. Table 1 lists the studied dimensions of the system in order to study the effects of system width, height and area. Note that the heat flux magnitude is considered variously for different simulations which is assigned according to the area of the heat storage system and also the bottom surface.
RT-58 (RUBITHERM) is considered as the PCM, and the physical properties of the RT-58 are listed in Table 2. 28 The reason for using RT-58 is that the melting point is suitable for domestic usage of LHS units in air and water heat exchangers.
Note that PCM is initially at 23 C to study the solid zone as well as the melting and liquid zones during the heat storage process.

Numerical modelling
To solve the governing equations, ANSYS-FLUENT software is used in this study with the aid of user-defined functions for calculating the interfacial heat transfer coefficient between the fluid and porous medium at each location and each time. The governing equations are discretized using the SIMPLE algorithm with Presto scheme for pressure correction equation and QUICK scheme for the momentum and energy equations. The values for under-relaxation factors are set to 0.3, 0.6, 1 and 0.9 for the pressure, velocity, energy and liquid fraction, respectively. The convergence criteria for continuity and momentum equations are set to 10 À6 and 10 À9 for the energy equation. Note that the higher values are also checked for the convergence criteria, and no change is seen in the results.
Different density of mesh is studied at the beginning, and, finally, the grid length is considered 1 mm in the x-direction and 0.5 mm in the y-direction due to the existence of natural convection in the y-direction. Therefore, for example, for the dimensions of 100 Â 25 mm, the number of nodes in x and y-direction are 101 and 51, respectively. The mesh density is equivalent to others in the literature in the rectangular geometry. 26,28 Furthermore, the time step size is considered 0.5 s for the porous-PCM and 0.1 s for the PCM-only. It should be noted that a lower time step size is also studied, and the results are almost the same, since the velocity is low, and the thermal front movement is similarly low, reflecting a low Peclet number and Courant number situation.
In the Kozeny-Carman equation, A m is the mushy-zone constant usually within the range of 10 4 -10 7 kg/m 3 s which should be determined based on an experimental study and varied with the geometrical and operational parameters of the system. 32 However, in the presence of a porous medium, due to the small effect of natural convection, A m has a small effect on PCM behaviour. 43 Figure 2(a) and (b) illustrates the variation of liquid fraction and mean temperature in terms of time for three different values of Am. As explained and shown, Am has small effects especially during the phase change when the PCM places between the solidus and liquidus temperature. Therefore, the values of 10 5 are considered for A m in this study.  5  800  25  1600  50  3200  100  6400  Area  50  100  6400  100  50  3200  200 25 1600

Model verification
Two different geometries are studied to verify the model for both porous-PCM and PCM-only cases in this study. For a porous-PCM, Zhao et al. 25,26 studied a metal-foam PCM with the porosity of 95% and pore density of 10 PPI using RT-58 numerically and experimentally. In Figure 3(a), the geometry as well as the boundary conditions can be seen. Liu et al. 28 also used the data from Zhao et al. 25,26 to verify their simulations.
The results of the present data in comparison with Zhao et al. 25,26 and Liu et al. 28 are displayed in Figure 4. As shown, an excellent agreement can be found between the non-equilibrium thermal model and the experimental results of Zhao et al. 25,26 and numerical results of Liu et al. 28 Furthermore, as mentioned in Liu et al., 28 regarding the numerical results of Zhao et al., 25,26 since a constant melting temperature is considered in their paper, a small variation was observed with the experimental results. However, in this study and the study of Liu et al., 28 different liquidus and solidus temperatures are considered for the simulations. Note that the temperature at the height of 8 mm from the bottom is presented in Figure 4.
To verify the code for a PCM-only system, the experimental study of Mat et al. 39 is simulated, and then the results are compared with each other. In the experiment of Mat et al., 39 a triple tube heat storage system with internal-external fins filled with RT-82 PCM was investigated, and the temperature is measured at different locations to find the average    temperature of the PCM. In Figure 3(b), a schematic of the simulated domain is displayed. Figure 5 shows the results of the current simulation in comparison with the results of Mat et al. 39 Temperature and liquid fraction are shown in Figure 5(a) and (b), respectively. Excellent agreement is achieved in this study in comparison with the experimental and numerical results of Mat et al. 39 Results and discussion In the following, first, for the geometry with the dimension of 200 Â 50, the results are presented completely for both porous and non-porous simulations, and then the results of other simulations are discussed.

Effect of porous-PCM compared with PCM-only
The presence of the porous medium inside the LHTES system results in the enhancement of heat transfer rate in the domain. Therefore, the heat is moved from the external source faster in the entire domain which can cause a shorter melting time. However, due to making a flow resistance, natural convection also has a small effect on the melting process. The other main advantage of using porous media is that the PCM melts and solidifies more uniformly than the case of PCM-only system. In this section, a fundamental discussion is presented by contour plots to better describe the difference between the PCM-only and porous-PCM systems. Figure 6 shows the contour plot of liquid fraction for the cases of porous-PCM (on the left) and PCMonly (on the right) at different times for the dimensions of 200 Â 100 and heat flux of 1600 W/m 2 . As shown, at the beginning (before 20 min), both systems melt from the bottom and the liquid is formed. After that, in the porous-PCM, the liquid cannot move properly in the domain due to flow resistance by the solid surface, and the heat is transferred by conduction mechanism through the high conductivity metal porous media through the PCM. However, in the PCM-only system, natural convection gradually becomes dominant for transferring heat to the top layers and by circulating a low conductivity PCM, and as a result of natural convection enhancement, the rate of heat transfer increases. Therefore, at the time of 40 min, PCM is melted more in the PCM-only system rather than composite PCM system. Another issue is that, as can be seen easily after the time of 60 min, since the heat is transferred to the top layer by the porous medium very quickly, the top layers of the PCM melt at almost the same time similar to the bottom layers, however, with a lower melting rate. However, for the PCM-only case, it can be seen that the top layers are solid until the moving interface of liquid-solid region rises to the top layers. This is a big advantage of porous PCM system especially in discharge time when the entire PCM solidifies almost at the same time. Figure 7 shows the contour plot of temperature at the times of 40, 80 and 120 min for the systems of porous-PCM (on the left) and PCM-only (on the right). As shown at the time of 40, due to absorbing heat by the porous medium and transferring it to the top layers due to a higher conductivity, the heat is transferred faster to the top, and higher temperature can be found at the top layers. However, for the PCM-only case, the temperature is at the initial temperature at the top layers. Furthermore, since the natural convection dominates, the Benard formation can be seen in the domain. However, for the porous-PCM, just the conduction heat transfer can be seen. As time passes, natural convection starts affecting the porous-PCM to a very small extent which is the reason for a higher temperature in a circulating shape in the presence of porous-PCM. However, for the PCM-only system, the generated vortices due to natural convection start joining together and forming bigger vortices as can be seen at the time of 80 min. This is established in the literature. 58 Furthermore, in the entire system during the melting process, a lower temperature gradient happens in the presence of porous medium; the magnitude of maximum temperature is less than that for the PCM-only case; the magnitude of minimum temperature is much higher than PCM-only. For example, the maximum and minimum temperatures are 357 and 327 K for the porous-PCM case and 372 K and 311 K for the PCM-only case at 80 min. Note that the variation of the temperature of porous medium is very small in the entire domain and changes from 329.9 to 332.2 K. It means that the maximum temperature difference between the liquid PCM and porous medium is 40 K. It shows the necessity of using the non-equilibrium thermal model on tackling porous-PCM. At the time of 120 min, when the porous PCM melts completely, and the PCM-only system is still in melting process, again, a low temperature gradient is seen in the porous-PCM and a high temperature gradient in the PCMonly case.
To show the effect of natural convection in the domain, Figure 8 displays the velocity contours as well as the streamlines for both systems of porous-PCM and PCM-only at different times of 40, 80 and 120 min. In the presence of porous media, the effect of natural convection is much lower which is well established in the literature. 26 However, in a porous-PCM, when a noticeable amount of liquid is made in the domain, the temperature rises significantly, and natural convection affects the distribution of temperature. At the beginning, several vortices are generated in the small melted area, and the Benard formation can be seen in natural convection. Then, as time passes and the melting area increases, the vortices merge together and create a bigger recirculation area. 58 Furthermore, the magnitude of velocity is lower by three orders of magnitude for the porous-PCM than the PCM-only, which proves small effect of natural convection in the presence of porous media. However, for the PCM-only system, the main mechanism of heat transfer is natural convection, which causes higher values of velocity magnitude in the domain. The average velocity magnitude inside the PCM is 2.8 Â 10 À7 m/s for the porous-PCM and 2.9 Â 10 À4 m/s for the PCM-only at 80 min. Figure 9 displays the average temperature of the porous medium, the PCM embedded in the porous medium and the PCM-only system as a function of time. As shown, the porous medium temperature is higher than the PCM temperature in the composite PCM foam until the time of almost 40 min, and then the temperatures are almost equal. Furthermore, the temperature of PCM in the presence of porous medium is higher than PCM-only system until the time of 40 min due to a higher rate of heat transfer through the metal foam. After that, as shown in the contours of temperature, the temperature in the PCM-only system is higher than the porous-PCM. At last, after the entire PCM is melted, the temperature rises and is equal for both cases of the porous-PCM and PCM-only system. Figure 10 presents the average temperature at different locations for both systems of PCM-only and porous-PCM. For the PCM-only system, the difference between the temperatures of different lines is high since enough time is needed that the heat is transferred to the top layers. For the mean temperature of the first line (Y ¼ 12.5 mm), at almost time of 40 min, the PCM starts melting at the temperature of 48 C. However, for the middle line (Y ¼ 25 mm), almost after 70 min, the PCM starts melting. Furthermore, for the porous PCM, due to a high rate of heat transfer, the temperatures of the lines are too close to each other, and for all of them, the PCM starts melting almost from the time of 40 min until the time of around 120 min from the beginning. For PCM-only system, the wall temperature increases with a high rate and then decreases sharply at about 20 min. This can also be seen for the porous-PCM case, but with only a small reduction. The reason is that at the melting point, PCM absorbs the heat from the bottom wall. By natural convection and recirculation of the flow for each of the generated vortices, cooler fluid moves to the bottom layers which absorbs more heat from the bottom and therefore causes a temperature reduction on the surface. Note that the melting point of PCM is a lot lower than the temperature of the bottom plate especially in the case of PCM-only system. This issue can be seen in the literature in case of heat flux boundary condition. 25  porous-PCM over the PCM-only system. In the solid zone, especially, and also the liquid zone, significant heat transfer enhancement occurs in the presence of the porous medium as indicated by low ÁT. Note that a reduction of temperature difference occurs due to the enhancement of natural convection for the PCM-only system. 25 Regarding the average temperatures in different phase zones (Figure 11(b)), it can be seen that the overall thermal diffusivity in the case of the porous medium is almost 2.5 times in the solid zone, 1.8 times in the mushy zone and 1.5 times in the liquid zone which is varied for different dimensions of the system. Note that in Figure 11(b), the mean temperatures at each phase zones are shown by colour lines for porous-PCM and PCM-only cases.
A comparison of the liquid fraction for both studied cases is shown in Figure 12. As shown, at the beginning, the liquid fraction of PCM-only system is a little higher than the porous-PCM. The reason is that the porous medium transfers the induced heat flux to the top layers, and the PCM does not receive most of the heat for melting in comparison with the PCM-only system. In other words, there is not a build-up of heat in the lower levels as the heat is more effectively transferred to the top level by the porous medium. Hence, it takes longer for the lower levels to reach the melting point. Without the porous medium, the heat is trapped in the lower levels, and    hence the PCM in the lower levels reaches the melting point quicker.
After that, the rate of melting increases when the heat is transferred to the entire domain, and the entire PCM starts melting. The melting time is higher for the PCM-only system, as expected, but the difference between the systems of porous-PCM and PCM-only in the melting is not too significant due to the small dimensions of the system, and as a result, small meting times are observed.

Effect of heat storage height
By increasing the height of heat storage system with constant heat flux of 1.6 kW/m 2 from the bottom, the melting time increases. However, in the porous-PCM system, due to the presence of high conductivity porous medium, the heat is transferred faster to the top layers, and therefore, the presence of porous medium is more effective. Figure 13 displays the variation of the liquid fraction as a function of time for different heights of the storage system. As shown, by doubling the height and therefore volume of the system, double melting time almost occurs, i.e. the rate of heat transfer enhancement is almost constant in porous-PCM systems. Figure 14 shows the melting time for different heights for the constant heat flux for the porous-PCM in comparison with the PCM-only system. As shown, by increasing the height of the heat storage system, the melting time increases at almost constant rate. Furthermore, the rate of enhancement in the melting time is higher for the PCM-only system due to lack of porous medium with the produced-free convection. However, the melting time for the storage system with the height of 12.5 mm for the PCM-only system is almost equal to the porous-PCM. After that, the difference of melting time between the PCM-only system and porous-PCM increases for higher heights of the storage system. Note that for very small heights, the PCM melts from the bottom entirely, and the presence of porous medium has no effect on the melting time. In other words, before that the porous medium is able to affect the PCM, all the PCM melts, and therefore the presence of porous medium has no effect. Figure 15 illustrates the melting time of both porous PCM and PCM-only cases considering different heat flux from the bottom according to Table 1. The heat flux is increased according to the increased volume of the system. As shown, for the porous-PCM, due to increasing the rate of heat transfer by increasing the surface area of the system boundaries, the melting time is almost equal for different heights of the storage system. Just for the height of 100 mm, a little increase in the melting time can be seen due to the large size of the storage system. However, for the PCM-only system, the melting time increases with increasing the height with much higher rate than the porous-PCM.

Effect of heat storage width
By increasing the width of the storage system for a constant height, if a constant heat flux is induced to the bottom wall, it is obvious that the melting time is almost constant. Therefore, to better understand the effect of heat storage width, the storage systems with various widths and heat fluxes are examined as presented in Table 1. Figure 16 displays the variation of the liquid fraction for different widths and heat fluxes of the storage system. As shown, similar to the height study, due to the presence of the porous medium, the melting time almost increases in a constant rate. Figure 17(a) and (b) displays the variation of melting time for both cases of PCM-only and porous-PCM for constant and variable-induced heat fluxes according to Table 2 from the bottom, respectively. As mentioned, for a constant heat flux, the melting times for different widths of the storage system are almost constant for both cases of PCM-only and porous-PCM. Just for the PCM-only system, the melting time is a little higher than that for the porous-PCM. As shown in Figure 17(b), by applying different heat fluxes proportional to the system's width, the melting time increases with a constant rate for both systems of PCM-only and porous-PCM. In contrast to the height study, the same increment can be seen between the melting time of PCM-only system and porous-PCM system similar to the constant heat flux.

Effect of the heat storage dimensions
To study the effect of constant surface area of the system, three different geometries, with the dimensions of 200 Â 25, 100 Â 50 and 50 Â 100, are studied according to Table 1. Since the amount of input heat should be equal for better comparing the systems, the heat flux is changed for them related to the bottom surface which is equal to 1600, 3200 and 6400 W/m 2 , respectively, for more rate of heat input. Figure 18 shows the variation of the liquid fraction for three studied heat storage sizes with the same area through the time. As shown, due to the presence of porous medium and enhancement of heat transfer rate inside the domain, the difference on the melting time is small between the studied geometries. Just in the case of 50 Â 100 due to a large height of the system  Table 1 for both cases of PCM-only and porous-PCM.  which causes a larger time to transfer heat from the bottom to the entire area even at much higher heat fluxes. Therefore, it can be concluded that even in the presence of porous medium in LHTES systems, by considering the heat source from the bottom, a system with a lower height has a shorter melting time and therefore is recommended. Figure 19 displays the total melting time for the systems with and without considering the porous medium for different dimensions of the LHTES systems with the same area. As shown, the effect of heat storage dimensions is more pronounced when a PCMonly system is used in the system. It has a small effect when a porous PCM is employed.

Dimensionless analysis
To generalize the effect of different parameters, a nondimensional analysis is performed to normalize the liquid fraction. The first non-dimensional number utilized in unsteady heat transfer problems in Fourier number is given as 36 where is the thermal diffusivity (m 2 /s), t is time from the start point of the melting and H is unit height. Stefan number is employed in latent heat storage problems as the ratio of sensible to the latent heat given as where 60 Therefore, Stefan number is defined as 60 For natural convection problems, Grashof number or more comprehensively Rayleigh number is used indicating the ratio of natural convection to conduction heat transfer and for a uniform wall heating flux is defined as 61 The mentioned dimensionless parameters are employed in different studies in the literature for defining the variation of liquid fraction. 42,[60][61][62][63][64] To find the variation of liquid fraction with the mentioned dimensionless parameters, since the problem is transient, the effect of Fourier number is studied. Based on the literature and due to the fact that the liquid fraction is directly a function of time, it is proposed that the liquid fraction is varied with the Fourier number directly. Figure 20 illustrates the variation of liquid fraction in terms of Fourier number for different Stefan number for PCM-only system. Note that for different simulations, Stefan number is calculated based on the dimensions and heat flux and also properties of the PCM according to equation (24). For a constant Stefan number, the variation of    liquid fraction is almost similar for different cases (shown for the Ste ¼ 2.32 and Ste ¼ 9.28). However, as shown, liquid fraction is varied by changing the Stefan number, and therefore the Fourier number alone is not sufficient to generalize the liquid fraction due to the effect of latent heat and phase change. As shown, a higher liquid fraction is achieved by reducing the Stefan number at a constant Fourier number, or a constant liquid fraction is achieved at a higher Fourier number for a lower Stefan number which is also reported in the literature. 61 The reason can be explained based on equation (24). For constant properties of PCM and a constant height, Stefan number decreases by increasing q w , and therefore, at a same time or for equal furious numbers, a higher liquid fraction is achieved. Figure 20(b) illustrates the variation of liquid fraction for different porous cases showing a similar trend with non-porous cases. Therefore, it is proved that the liquid fraction is a function Stefan number. Based on the literature review and also due to the fact that the liquid fraction is changed directly with the inverse of latent heat of fusion based on the storage capacity of LHTES systems, the liquid fraction is directly a function of Stefan number. The independent parameter of Fo.Ste is almost enough for the generalization of the melting process; however, due to the effect of natural convection, especially for high Stefan number cases and especially for the case of PCM only, when the effect of natural convection is significant, Rayleigh number is better to be included to propose the effect of natural convection. After careful calculation to find the best combination of Fo.Ste with Ra, Figure 21 Since the simulation results correlate well with the proposed non-dimensional number, a curve-fitting analysis is performed between liquid fraction and x ¼ Fo.Ste.Ra À0.02 for both PCM-only and composite-PCM systems. The fitted polynomial curves are as follows PCM-only L f ¼ À8:06x 4 þ 6:8729x 3 À 0:8998x 2 þ 1:5144x À 0:0153 Composite-PCM L f ¼ À5:7659x 4 À 1:9444x 3 þ 6:1535x 2 À 0:0379x þ 0:0172 The R-square values for the PCM only and composite-PCM correlations are 0.9922 and 0.996, respectively. Note that the fitted equations are reasonable for the rectangular geometry within the parameter range: 1.16 < Ste < 37.13, 0 < Fo < 1.5, 2.9 Â 10 4 < Ra < 9.5 Â 108, 0 < L f < 1, 0 < Fo.Ste. Ra À0.02 < 0.57. The maximum deviation of liquid fraction from the curve fitted equations is 0.065 and 0.051, respectively, for the PCM-only and composite-PCM cases showing the good accuracy of the predicted correlations.

Conclusion
Previous research in the area of thermal energy storage with PCM has considered the application to a generic storage system, typically at a large scale for central delivery of thermal energy for industrial and large-scale architectural purposes. This work shows that if a PCM storage heater is required in domestic applications, with a unit on the wall of a domestic dwelling, then the surface/volume characteristics and the charging surface position and size are significantly important.
Different dimensions of the thermal storage system show the effect of a high conductivity porous medium in different heat storage system volumes on the melting time, temperature, velocity distributions and heat transfer rate. After verification of the code, the results showed increase of heat transfer rate by using the porous-PCM storage system. However, the effect depends on the size of the storage unit. Higher effect of porous medium on melting process is observed in the PCM's solid-phase zone when the thermal   diffusivity increases by almost 2.5 times greater than the PCM-only system. Furthermore, the presence of porous medium on the melting time is highly related to the aspect ratio (width and height) of the storage systems. For the PCM-only system, the height has a higher effect than the width due to the contribution of natural convection in the free flowing system. Considering equal input heat, by changing the dimensions of the storage system with the same surface area and volume, the melting time reduces by 18 The heat flux rate is therefore significant in determining the efficacy of the porous medium's influence on the PCM. The non-dimensional analysis results in curve-fitting correlations between the liquid fraction and Fo.Ste.Ra À0.02 for rectangular latent heat storage systems for both PCM-only and composite-PCM systems within the parameter range of 1.16 < Ste < 37.13, 0 < Fo < 1.5, 2.9 Â 10 4 < Ra < 9.5 Â 108, 0 < L f < 1 and 0 < Fo.Ste.Ra À0.02 < 0.57.

Declaration of Conflicting Interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was funded by the EPSRC (Engineering and Physical Sciences Research Council) via Supergen Energy Storage II, grant reference EP/P003435/1, titled 'nano-Structured PCM Composites for Compact Space Heating: n-CoSH.