NUMERICAL STUDY OF SOME INFLUENCING PARAMETERS ON MELTING PROCESS OF PHASE CHANGE STORAGE UNITS INTEGRATED WITH A SOLAR WATER HEATER

This work theoretically analyzes the melting process of a phase change material (PCM) packed in thermal energy storage units. The units are filled with paraffin wax and inserted in the water storage tank of a solar water heater to improve the storage capacity. A computer program was built to simulate the performance of the whole system during a typical day between 8 a.m. and 3:45 p.m.. The simulation period was divided into small time intervals of 10 seconds. The Finite Difference Method (FDM) was used to analyze the melting process of the PCM units. Each unit was discretized into a large number of nodes forming a two dimensional grid. The effect of three PCM parameters were studied which are the shape, volume and surface area. The cross-sectional shape included square, annular square, circular and circular annulus bars. It was found that the circular sections provide higher melt ratio than the square sections and that the annular sections are more effective than the solid sections. Increasing the PCM surface area at constant PCM volume significantly increases the melt ratio for all PCM shapes. The system equipped with PCM having the least volume operates more effectively than a system with larger PCM volumes.


Introduction
The storage tank in a solar hot water system incorporates water as a storage medium.Heat is stored sensibly in the water while its temperature increases.When a Phase Change Material (PCM) is inserted in the storage tank along with the water, it stores energy as latent heat of fusion.The PCM undergoes melting while absorbing heat at approximately constant temperature.
The first application of Phase change material (PCM) described in the literature was studied by Telkes [1,2] in the 1940s for heating and cooling in buildings, a survey of various applications of PCM for thermal storage was given by Lane [3].Telkes [4] published the idea of using PCMs in walls, better known as Trombe walls.
Yanadori and Masuda [5] used calcium chloride hexahydrate as a latent heat storage material.The material was filled in a vertical cylindrical heat storage container.Water flowed in a vertical single copper pipe placed at the center of the container and transfers heat to the PCM around it.An experimental equation for thermal design of PCM containers was developed.
Ghoneim [6] developed a numerical model to simulate the performance of cylindrical PCM units packed vertically in a vertical vessel.The working fluid flows inside the vessel parallel to the PCM tubes.The Finite Difference Method was used to solve the thermal model.Three types of phase change materials were used to obtain the annual solar fraction as a function of PCM storage mass.
The work of modeling the melting phenomenon was developed further with the work of Darling et al. [7].They included in their model other parameters such as buoyancy, surface tension, viscosity, density difference and void growth.These parameters determine the location of the phase change interface.The energy, momentum and continuity equations were solved using a finite volume formulation.Velraj et al. [8] performed experimental and numerical analyses on a cylindrical PCM unit enhanced by internal longitudinal fins.They studied the inward solidification process of paraffin RT60 inside a vertical cylinder equipped with vertical fins.The cylinder was cooled from outside by water.
Lamberg and Siren [9,10] presented simplified analytical models which predict the solid-liquid interface location of the PCM and temperature distribution of the adjacent straight fin.The first paper tackled the melting process continued by the second paper that tackled the solidification process.The results were compared to numerical results of the two dimensional case with good agreement.A new factor called the fraction of solidified PCM was also introduced.
Liu et al. [11] experimentally studied the melting processes of stearic acid in an annulus cylindrical container.
Vyshak and Jilani [12] presented a comparative study of the total melting time of a phase change material packed in three containers of different geometric configurations, viz.rectangular, cylindrical and cylindrical shell, having the same volume and surface area of heat transfer.Employing a slightly modified enthalpy method, it was concluded that cylindrical shell containers take the least time for melting.
Bony and Citherlet [13] developed a numerical model to simulate heat transfer in phase change materials plunged in water storage tank.Comparison of the model with experimental results showed good agreement.
Gunther et al. [14] presented a new algorithm to simulate the PCM subcooling problem.Enthalpy approach and finite difference implicit scheme were used and said to be particularly useful in simulating applications with low cooling rates such as in building applications.
Shukla and Singh [15] conducted theoretical analysis of a thermal storage unit containing paraffin wax.The unit was intended for absorbing night coolness and provides cooled air at comfort temperatures during daytime in summer season.A MATLAB simulation tool was used to compute the variation of air temperature with various parameters affecting the system which included; location, time in year, air flow rate, PCM mass as well as charging and discharging time.
Lopez Navarro et al. [16] employed experimental results for characterization of a versatile latent heat storage tank capable of handling an organic PCM having melting temperature in the range between -10 o C and 100 o C. The enthalpy-temperature curve, the specific heat and density were measured for the paraffin wax being tested.
Antony Aroul Raj et al. [17] developed a transient numerical model to study the melting and solidification processes of paraffin wax in an annular cylindrical container.Enthalpy method was adopted to model the system and solved by Finite Difference Method in its implicit scheme.The temperature variation of the PCM was calculated along the two axes of the polar coordinates (r, z) along with the time required for melting and solidification.
Another numerical study was conducted by Yosr Allouche et al. [18] who simulated the heat transfer processes of a PCM inside a 100 l cylindrical tank, horizontally placed.The PCM used was microencapsulated slurry in 45% w/w concentration.ANSYS/FLUENT commercial software was used to solve the governing equations and the results were also validated experimentally.
Abdel Illah Nabil Korti and Fatima ZohraTlemsani [19] experimentally studied three different types of PCMs for latent heat thermal storage.Water was used as the heat transfer fluid.The temperatures of the PCM and water, solid fraction and thermal effectiveness were analyzed.The effects of the inlet temperature and the flow rate of water and the type of PCM used on the time required for charging and discharging were also discussed.It was concluded that the flow rate of water is more effective in charging than in discharging of heat while the water inlet temperature had more pronounced effect in this regard.
A similar experimental study was done by Meng and Zhang [20], however, the PCM used was thermally enhanced by combining it with copper foam.The study was backed by a three dimensional numerical model employing enthalpy-porosity method.The addition of copper foam caused a significant improvement in PCM thermal conductivity which yielded a decrease in charging and discharging times.
In all previous investigations, the effects of the shape, volume and surface area of PCM units were not accounted for.The present work numerically simulates the behavior of PCM storage units of various shapes, volumes and surface areas, surrounded by the storage water of a solar water heater.The varying water temperature throughout the day constitutes the boundary condition influencing the phase change process of the PCM units.A computer simulation procedure was developed to estimate the storage water temperature surrounding the PCM storage units.The varying tank temperature constitutes the driving heat that controls the phase change.The details of the simulation procedure are given by Al-Tabbakh [21].The thermal models of the PCM units in the simulation program are solved numerically using a modified form of the transient Finite Difference Method.The Melt Ratio of the PCM at each time interval is estimated using this technique.The PCM is studied in terms of cross sectional shape, volume and surface area.The effect of each parameter is studied by fixing the remaining two.

PCM Units Configuration
The phase change storage units considered in the present study take the shape of rectangular and cylindrical rods immersed in the water storage tank (Fig. 1).The PCM rod height is 0.707 m which equals the tank height.Both heights are kept constant in the simulation.The cross sectional area of each rod ranges between (0.0025 -0.01) m 2 depending on the PCM volume and surface area of heat transfer.Heat is transferred from the surrounding water to the PCM rod equally from all its lateral sides.Therefore, the analysis adopted is two dimensional.The two-dimensional conduction problem with a change of phase does not have a known analytical solution [9].It will be solved numerically using the transient Finite Difference Method (FDM) in implicit scheme.The method will be adapted to take into account the change of phase.In the case of square and annular square sections, Cartesian coordinates are used and the section is discretized into 441 nodes, 21 in each direction.The nodes or more precisely, the rectangular PCM lumps surrounding the nodes will encounter sensible heating until reaching the melting point.At the melting temperature, a heat balance is initiated for each node that reached the melting temperature.The heat balance will compare the required amount of heat to melt all the mass of the PCM lump with the actual amount of heat transferred to the lump from the neighboring lumps.If this amount of heat equals or exceeds that required to complete melting in the considered interval, then that lump is considered in a liquid state and the next calculation will be switched back to the sensible heat scheme, whereby the temperature of the liquid PCM will keep increasing.If the heat received by the lump at the considered interval is not enough to complete melting, which is the most probable case, the program will calculate the part of that lump that turned to liquid.Successive additions of heat to the lump will increase the molten PCM until the whole lump becomes liquid.The overall melt ratio is the summation of the molten parts of all the 441 lumps divided by the total mass of the rod.In the case of circular and circular annulus sections, cylindrical coordinates will be used to analyze the heat transfer problem.The section will be discretized into 21 annular concentric lumps and the same procedure followed in the Cartesian coordinates is applied in the cylindrical coordinates.

Thermal Analysis
The PCM is analyzed as a two-dimensional transient conduction medium governed by the following equation given by Chapman [22]:- When the PCM sections are discretized, several locations of nodes will appear.Some at the outer boundaries, others at the inner boundaries, and at interior and exterior corners.Each type of node requires a different finite difference formulation.The general Finite Difference equation applicable to all nodes is [22] Eq. ( 2) takes a slightly different form with different notation at each node.R ij is the thermal resistance between the node in question and its neighboring ones.It takes the following form:- Where l ij is the distance between adjacent nodes in the x or y directions.A ij is the area normal to the heat transfer direction which is normal to the length l ij .C i is the thermal capacitance defined as follows:- Where ΔV i is the lump volume that surrounds node i.When the calculations for each 10 seconds interval are finished, the result will be a set of nodal temperatures in a specific section.These temperatures are used to estimate the heat transferred to the PCM unit from the surrounding water.Thus, a heat balance is initiated for each node.The node, or the PCM lump surrounding the node, will receive heat from the higher temperature nodes, and deliver it to the nodes at the lower temperature.This process is carried out for each of the 441 nodes and at each time interval in the calculation procedure.When the heat received is greater than that delivered, then, the PCM is absorbing heat from the water and vice versa.The net sensible heat transferred to or from any node at each time interval is calculated by the following equation:- When melting commences, the conventional heat diffusion equation (eq.( 1)) is no longer valid.The governing energy equation in the presence of phase change is [9]:- This equation is difficult to solve numerically.However, the set of nodal temperatures given by the solution of eq. ( 1) along with the amount of heat transferred to each node are used to calculate the molten PCM at each node.The sum of all molten PCM of all nodes is used to estimate the total melt ratio at the specific time interval which is an implicit solution of eq. ( 6).When a node reaches the melting temperature, the heat transferred will cause a change of phase.The heat necessary to melt the PCM mass in the specific lump is found by the following equation:- Where h fu is the latent heat of fusion of the PCM.So, if q i,s >q i,m then all the lump mass is considered to be converted to liquid.If q i,s <q i,m then the amount of PCM mass converted to liquid within the PCM lump is found by the following equation:- The molten PCM for the whole 441 nodes is found by summing up the molten parts of each PCM node as follows:- The total melt ratio of all the nodes for that time interval, with reference to the total PCM volume V p , is calculated by the following equation:- Similarly, the total amount of heat absorbed by the 441 PCM nodes is:- The solid-liquid interface moves continuously during the phase change and takes a varying contour during the process.The present computer program is capable of giving all the nodal temperatures of the cross section at each 10 seconds.The location of the interface is considered as the locus of the nodal points that are at the melting temperature.If there are no nodes at exactly the melting temperature, an interpolation can be used to locate the contour between the nodes that are above and those that are below the melting point.As a result, contour plots can be drawn at successive time intervals which give a clear picture of the shape of the interface and the speed of its progression.

Results and Discussion
The typical results presented here were for a typical clear day on the 21 st of December for Baghdad, Iraq (33 o N).The simulation period extended from 8 a.m. to 3:45 p.m. Four PCM rod cross sections were studied (Fig. 1).To ensure constant PCM volume and surface area for all cross sections, the number of PCM rods must be changed.The ratio between the volume and surface area of each rod remained constant and the small and large rods remained dimensionally similar.
Fig. (2) shows the discretization of a square and an annular square sections.The combined cross sectional area and the perimeter of the two solid squares equal the exposed cross sectional areas of the inner and outer squares and their perimeters.However, the number of nodes in the annular section is approximately equal to the number of nodes in only one of the squares.The near equality of the number of nodes ensures dimensional similarity, where the ratio of each PCM cell length to the section side length is nearly equal for two sections.The equality of the cross sectional areas and perimeters ensures the equality of PCM volumes and surface areas because the length of the rods is the same.Fig. 3 shows the variation of the melt ratio for the four sections.Melting starts after the required period to bring the PCM from its initial temperature to the melting point.This period depends on the initial water and PCM temperatures and the PCM geometry and volume.The four sections start melting at the same time but afterwards differ in the rate of melting.The circular annulus section gives the highest melt ratio.The annular square section exhibits a melting curve very close to the circular section, while the square section has the lowest melt ratio of the four.
The results for increased surface area are shown in Fig. 4. Increasing the surface area of heat transfer results in an increase in the melt ratio for all the sections and complete melting was achieved for all sections.
The circular annulus section takes the least melting time (about 1:45 hours) followed by the annular square, circular and square sections respectively.For the circular annulus section, Fig. 5 shows that the surface area greatly influences the melt ratio and melting time.Increasing the surface area in the present work was achieved by increasing the number of PCM units and decreasing the size of each unit while keeping the total PCM volume constant.Fig. 6 shows the variation of PCM melt ratio for three PCM volumes, which are; 0.028 m 3 , 0.042 m 3 , and 0.056 m 3 .Increasing PCM quantity leads to earlier melting.The tank with more PCM contains less water and reaches the melting temperature faster.However, early melting does not mean an early melting completion, because the melting rate becomes slower with increased PCM quantity.So, larger PCM quantity takes more melting time despite the early start of melting.Figs.7 and 8 show contours for square and annular square sections at different hours in the simulation period.The low thermal conductivity of the PCM results in a slow penetration of heat through the material.The contours are at four consecutive times two hours apart.The corners of the sections tend to receive heat and melt faster than the core regions.A temperature difference of about 10 o C is observed between the boundaries and the core at 10:00 a.m. and 12 p.m. for a solid square section of 7.07 cm side length.The temperature difference increases towards the end of the simulation period.The temperature contours for the annular square section, having a side length of 10.6 cm and an inner side square length of 3.53 cm. in Fig. 8 show a smaller difference between the highest and lowest temperatures in the annular square section than in the solid square section.Higher melt ratios were reached at 2:00 p.m. and 3:45 p.m. in the annular sections indicating more effective heat transfer.
The shape of the PCM section has a small effect on the tank temperature for all the studied PCM volumes and surface areas as shown in Figs. 9 to 12.In other words the same share of heat is received by all sections regardless of their shapes.The effect of shape is confined in changing the melt ratio and each section takes a different time for melting.Each section has a different effective penetration depth for the flow of heat.The circular sections (solid or annular) utilize the same amount of heat more effectively than the square sections (solid or annular).Also, the annular sections (square or circular) are more effective than the solid sections (square or circular).This fact may be further investigated in future works.Fig. 11 shows that increasing PCM surface area of heat transfer makes the tank temperature closer to the melting temperature.The higher heat transfer associated with the increased area utilizes the collected heat in melting to a larger extent than to increase the water temperature.When the PCM volume increases, its effect on the storage water temperature becomes more pronounced and the water temperature follows the PCM melting temperature until complete melting is attained as shown in Fig. 12.
Table 1 shows the values of the final tank temperature, the PCM stored heat, the total stored heat (water and PCM) and the final melt ratio for all the cases investigated for the PCM having a melting temperature of 46 o C in the absence of tank loading.Clearly, complete melting is always accompanied by a decrease in the final tank temperature and an increase in the stored heat.The PCM Heat Ratio is the ratio of the heat stored in the PCM to the total heat stored.

Conclusions
a) Circular annulus PCM sections provide the best melt ratio.The maximum increase in PCM melt ratio was 28% over other geometries.b) Increasing the PCM surface area of heat transfer at constant PCM volume significantly increases the melt ratio for all PCM shapes.c) A significant temperature differences exists between the core and the boundaries of a square PCM section.The other sections showed lower temperature differences between the core and the boundaries.d) PCM melting caused a lower water tank temperature because the received heat is utilized for phase change rather than increasing the tank temperature sensibly.e) Increasing the PCM volume 100% for a constant total tank volume, decreases the melt ratio about 25% and decreases the rate of melting by 55%.

Fig. ( 2 )
Fig. (2) The discretization of two sections:-(a) annular square and (b) solid square.The combined two square sections in fig.(b) have surface area and volume that are equal to their counterparts for one annular square in fig.(a).(All dimensions in meters).

Table (
* per square meter of collector area.