Numerical Simulation of the Heat Transfer Behavior of a Zigzag Plate Containing a Phase Change Material for Combustion Heat Recovery and Power Generation

1State Key Laboratory of Multiphase Complex Systems, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100190, China 2Beijing Shenwu Environment & Energy Technology Co., Ltd., Beijing 102299, China 3Global Energy Interconnection Research Institute, State Grid Corporation of China, Beijing 102211, China 4School of Engineering and Technology, China University of Geosciences, Beijing 100083, China 5Birmingham Centre for Thermal Energy Storage, University of Birmingham, Edgbaston, Birmingham B15 2TT, UK


Introduction
With the global development of heavy industrial sectors, an increasing consumption of fossil fuel which resulted in the risk of environmental degradation is becoming increasingly serious every day.The global energy consumption has been growing more rapidly than the economy over the past decades and in 2010 its growth reached 5.6% leading to a highest record since 1973 [1].In 2012, China's energy consumption grew by 11.2% overtaking the US as the world's largest energy consumer [2].But world primary energy consumption grew by 1.8% in 2012, well below the 10-year average of 2.6%, which benefited from the utilization of renewable energy especially the solar energy and the energy conversation.However, solar energy is susceptible thus irregular in nature and unavailable during cloudy time and nighttime.As a result, to collect and store the solar energy in the daytime for continual supply or evening use is critical.In China, the electricity generated by coal-fired power plant accounts for 77% of the total quantity of power industry.Therefore, in perspective of the energy utilization, it is important to make full use of the low-grade heat from combustion process and improve the total thermal efficiency of the whole plant.Because of its simple structure and high thermal efficiency to run in terms of the recovery of low temperature waste heat, phase change materials (PCMs) thermal heat storage technology is an ideal way to build a bridge between the energy consumption and utilization.Thermal energy can be stored in the form of latent or sensible heat or both.Compared with sensible heat thermal energy storage (SHTES) material, PCMs is particularly attractive since it can provide high storage capacity within a relatively small temperature range.Thus, the latent heat thermal energy storage (LHTES) systems have better prospect to fulfill the needs of the energy hungry society and improve the performance and reliability of the energy system [3], which is more attractive in the field of solar energy collecting [4], industrial waste heat recycling [5], and peak electricity shifting [6].However, because of the low thermal conductivity of PCMs (usually inorganic salts or organic materials), the performance of the LHTES systems is drastically affected during energy retrieval or withdrawal, which leads to low charging and discharging rates of the LHTES systems [7].
In order to overcome the limitations of the LHTES units, a large number of heat transfer enhancement techniques have been approached from different aspects.Fukai et al. [8], Mettawee and Assassa [9], and Sarı and Karaipekli [10] applied, respectively, carbon-fiber, aluminum powder and expanded graphite into the PCMs to improve the thermal conductivity of the PCMs themselves.Another key aspect is to optimize the performance of the heat exchanger, such as introducing fins into the shell-and-tube exchanger [11][12][13] and heat sink exchanger [14] and utilizing coil pipe [4].In addition, other designs for enhancing the heat transfer, such as spinning tube [15], microchannel [16], and multiple tube [17], have also been studied.
Because of the compact structure and large heat transfer area, the plate heat exchanger filled with PCMs is an ideal design.However, most studies on the LHTES unit employed the heat exchangers of shell-and-tube type or tube-tank type, thus the plate type LHTES unit was rarely studied.Besides, to the authors' best knowledge, no studies have been found in the literature about the impact of the zigzag surface topography on the heat transfer of plate type LHTES unit.Therefore, this work concentrates on the impact of zigzag configuration on the heat transfer of a plate type LHTES unit.A new case of a fixed input power with various combinations of inlet velocity and inlet temperature during the charging process was also studied.In order to mitigate the heavy demand on computing resources, the heat exchange unit employing zigzag configuration was assumed to be a two-dimensional model and the calculation was focused on the melting process.Moreover, an experimental system was designed accordingly to validate the model.flows between the two PCMs plates, leading to the melting or solidification of the PCMs upon absorbing/releasing heat from/to the HTF.For the smooth plate (a), the length (L) and thickness (2ℎ  ) of the computation domain are 360 and 3.75 mm; the geometries of the smooth plate are given as ℎ 1 = ℎ 2 = 4 mm; and the thickness of the plate wall is 1 mm, as shown in Figure 1(a).For the zigzag plate, the geometries of the zigzag plate are given as  1 =  2 = 2 mm;  3 =  4 = 4 mm;  5 = 2 mm and the thickness of the plate wall is 1 mm, as shown in Figure 1(b).Since the configuration of the heat exchanger is periodic, only one segment was employed for simulation (see Figure 1).Each LHTES unit, either smooth or zigzag plates type, contains the same amount of PCMs which is salt mixture composed of sodium chloride and magnesium chloride (NaCl/MgCl 2 ) with mass ratio 4 : 6.The physical properties of the PCMs and the HTF are tabulated in Table 1.

Mathematical Model.
The continuity equation is as follows: where  is density and ⃗  is the fluid velocity.The momentum equation is as follows: where  and V are the velocity in x and y directions;  is viscosity;  is pressure; t is time;   and   are the momentum sources in the enthalpy-porosity formulation in the form of   =  mush   (1 − ) 2 /( 3 + ) in which  mush is mushy zone constant that measures the amplitude of the damping and  is a small number (e.g.0.0001) to prevent division by zero; and   = (ℎ − ℎ ref )/  is the buoyancy source term.
The conservation equation of energy is as follows: where   = (Δ)/ + ∇( ⃗ Δ) is a phase related source term and the enthalpy of the material is the sum of the sensible enthalpy, ℎ, and the latent heat, Δ, in (3): where ℎ ref is the enthalpy value at the reference temperature  ref and  is the liquid fraction which is defined as follows: In (2),  is viscosity;   = (ℎ − ℎ ref )/  is the buoyancy source term;   and   are the momentum sources in the enthalpy-porosity formulation in the form of   =  mush   (1 − ) 2 /( 3 + ) in which  mush is mushy zone constant that measures the amplitude of the damping and  is a small number (e.g.0.0001) to prevent division by zero.In (3), ⃗  is the fluid velocity,  is the thermal conductivity,  and  are viscosity and density, respectively, and   = (Δ)/ + ∇( ⃗ Δ) is a phase related source term.In order to reduce the computational cost of the simulations, the geometry of the plate is considered to be sufficiently wide and the 2-dimensional model is used to simulate the thermal behavior of the PCMs-plate.Besides, the geometry of the PCMs-plate is periodic and only one channel surrounded by two halves of the PCMs plates is simulated.The computational domain is shown in Figure 1.

Boundary and Initial Conditions.
The interface between the PCMs and the inner tube was coupled-wall boundary condition for heat transfer.The inlet was assumed as velocity inlet with uniform velocity and temperature (which was higher than the melting point of the PCMs).The outlet was set as outflow boundary.The upper boundary and lower boundary of the simulation domain are assumed to be periodic.The front and back walls are considered as adiabatic boundary conditions.
The entire calculation domain was set initially at a zero velocity and a temperature of   (which was lower than the melting point of the PCMs).

Experimental Validation
3.1.Experiment System and Procedure.The numerical model was validated by experimental work.Figure 2 shows the schematics of the high temperature heat exchange system, consisting of an electric furnace, a high temperature blower, a LHTES unit, and a number of data acquisition equipment.During the charging process, the air was driven to circulate in the loop by the high temperature blower where it is heated by the electric furnace.The heat was transferred to and stored within the PCMs when the heated air flowed through the latent heat exchanger.The temperature of PCMs and outlet is measured by K-type thermocouples (the measurement range is 273.15K∼1273.15K and error range is ± 0.2%) and the flow rate of the air is measured by a orifice meter (error range ± 2%).Upon fully charged condition (judged from the temperature of the outlet air), the charging step would be ended by switching off the electric furnace and the high temperature blower.During the discharging process, the ambient air was switched on to deliver an air flow at a preset rate to extract the stored heat from the charged LHTES unit.Temperature of the phase change materials was measured by thermocouples (see Figure 3 for the thermocouple  arrangement).As shown in Figure 3, five K-type thermocouples connecting to a PC-based temperature data acquisition system inserted into the PCMs to measure the temperature of the PCMs.The exterior of the heat exchanger is shown in Figure 4.

Validation of the Simulation Models.
The computational domain is represented by a 2D base grid as shown in Figure 5.The finer grid is created by increasing the mesh density near both sides of the plate walls.In order to check the grid independence, grids of three sizes were checked, that is, 0.3 mm, 0.4 mm, and 0.5 mm with the cell amount of 84932, 63410, and 49163, respectively.No obvious differences in the liquid fraction were found among the three different grids.Therefore, the base grid with grid size of 0.5 mm was chosen in the following simulation, composed of 49163 cells including 19560 cells in the PCMs domain, 13770 cells in HTF domain, and the rest in the zigzag plate.
An initial run was performed and compared with the experimental data to validate the computational modeling of the melting process.In the experiment, the inlet temperature and the inlet velocity were 863.15K and 3.1 m/s and the NaCl-MgCl 2 mixture with a 4 : 6 mass ratio was used, presenting the phase change temperature of 713.15 K. Other thermophysical properties of the PCMs are shown in Table 1.The specific heat and thermal conductivity of HTF (air) within range of 200∼1000 K were given as follows.
As shown in Figure 6, the simulated value agrees reasonably well with the experimental data.The measured outlet temperature as a function of time is compared with the calculated result and is found to be in good agreement.However, the deviation of thermocouple 1 between the simulation and experiment in Figure 6 could be caused by several issues.One is that the PCMs with pores in the form of powder were used in the experiment, which reduces the conductivity of PCMs and more heat is delivered to the downstream region; other reasons may include that the properties of the PCMs are correlated with the varying temperature and the heat transfer fluid flow is not identical in all different plates, which were not considered in this simulation.7), respectively.The charging process, during which the heat is absorbed and stored, can be divided into three steps.At the beginning, the heat is stored in the PCMs as sensible heat because the temperature of the PCM is below its melting temperature and the temperature difference between the HTF and the PCMs is large.Next, part of the PCMs starts to melt as its temperature increased and the heat is stored as latent heat during the phase change process, while the readily melted part and the solid part of the PCMs are absorbing and storing heat as sensible heat without phase change.At the end of the charging process, PCMs in the exchangers are melted completely, so the heat is stored in the sensible way.

Results and Discussion
From Figure 7 it can be seen that the zigzag surface topography (with the length of 360 mm) can enhance the heat transfer of the LHTES system remarkably, since it costs about 3750 s to complete the whole melting process compared to that of about 4200 s by the system employing smooth plate of the same length with the zigzag plate.This effect can be attributed to two reasons.On one hand, the zigzag shape surface has the ability to change the flow form from laminar to turbulent by disturbing the flow direction, which can decrease the boundary layer thickness resulting in the decrease of thermal resistance at the HTF plate interface; on the other hand, the zigzag configuration presents larger heat transfer area than the smooth configuration, enhancing the heat transfer rate during the charging and discharging process.However, the smooth plate, which has the same heat transfer area with the zigzag plate, costs about 3850 s to complete the whole melting process.The zigzag surface topography does not have much effect on the heat transfer enhancement if the zigzag plate has the same heat transfer area with the smooth plate.But the zigzag topography has the ability of shortening the plate length.

The Effect of Tooth Shape on the Charging Process.
In this section, zigzag plates of four distinct zigzag lengths, that is,  = 12, 18, 24, and 30 mm, are studied and the simulation results of the liquid fraction as a function of time is shown in   8 that the four liquid fraction curves are almost coincident with each other over the entire charging process, implying that the zigzag length almost has no influence on the heat transfer capability of the LHTES system.It is well known that the heat transfer area and the flow distribution have great impact on the heat transfer capability of plate heat exchanger.With the increase of the tooth length, the heat transfer area decreases, but the adherence character of the fluid strengthens.As shown in Figure 9, the mainstream region is near the middle part of the passageway in the plates of small tooth length, but as the tooth length gets bigger, the mainstream region gets closer to the plate wall, where in the latter case the heat transfer is enhanced.Thus, the joint effect from the above two aspects leads to the tooth length almost having no impact on the melting process.
Furthermore, by changing the individual plate lengths, the impact of tooth rates on the heat transfer capability of the zigzag plate is explored.The parameter, tooth rate, is defined as  =  1 / 2 along the flow direction as shown in Figure 10. Figure 11 shows that the temporal variation of liquid fraction for different tooth rates are almost coincided with each other under the condition of V in = 6 m/s and  in = 923.15K. Thus, the tooth rate has little influence on the heat transfer capability, where two aspects should be considered.On one hand, the total tooth length ( 1 +  2 ) of the five plate heat exchangers with different tooth rates are equivalent but their individual heat transfer areas are different, which are 760.3mm, 747.984 mm, 745.032 mm, 747.984 mm, and 760.3 mm for the tooth rates of 1/5, 1/2, 1/1, 2/1, and 5/1, respectively.However, this difference is too small to produce noticeable impact on the heat transfer process.On the other hand, although the flow fields within the different plates are different with each other, they caused undistinguishable effect on the heat transfer capability of the heat exchangers.It is well known that the near-wall flow field formed near the upwind wall can enhance heat transfer more effectively than that formed near the lee side wall.The near-wall regions of the five cases showed distinct characteristic.When  < 1/1, all the near-wall regions are formed near the upwind wall at the two sides of the HTF passageway.However, when  ⩾ 1/1, the near-wall regions are formed, respectively, on the upwind wall and the lee side wall at the same side of the HTF passageway, as shown in Figure 12, which is the schematic plot of flow field near the wall.Figure 13 shows the velocity vector for different tooth ratio and, as shown in Table 2 and Figure 13, there is nearly no difference in   (  =   - +   -) between the first two cases with  < 1/1, but   of the later three cases with  ⩾ 1/1 increases with the increase of U. Besides, for the later three cases,   - (near lee side wall) increases with increasing U but   - (near upwind wall) decreases with increasing U. Therefore, the length of near-wall regions (  ) of the first two cases is smaller than that of the later three cases, but the near-wall regions of the first two cases locate close to the upwind wall while the two near-wall regions of the later three cases locate close to, respectively, the upwind wall and the lee side wall.Therefore, by considering the compromise between the length and the location of the near-wall regions, the change of flow filed caused by tooth rate almost has no impact on the heat exchange capability of the plates.
The heat stored by PCMs can be separated into sensible heat and latent heat due to the temperature rising and the reversible phase change between liquid and solid.In order to evaluate the ratio of latent heat (  ) in total stored heat (Q), the value of   / is studied in this section and the sensible heat (  ), latent heat, and total heat can be expressed as follows: where  is the total heat stored in the PCMs and  is the liquid fraction of the PCMs.
The pressure drop can be calculated as follows: where   and  in are the area of the outlet and inlet of computational domain, respectively.
Figure 13: Velocity vector for different tooth ratio.
Figure 14 shows the ratio of the latent heat (  ) in the total heat (, composed of latent heat and sensible heat) for various zigzag lengths.As shown in Figure 11, the value of   / rises quickly and gets to the first peak value at c.a. 500 s followed by a small decrease before c.a. 1000 s.This first   / peak can be attributed to the large temperature difference between the HTF and the PCMs at the beginning, resulting in the melting of the PCMs at the inlet region which rapidly absorbs and stores the heat of HTF as latent heat.After certain time, once the PCMs at the inlet region are completely melted, the PCMs continue absorbing heat and store it in sensible form leading to the decrease of   /.As the above process proceeds, the temperature difference between the PCMs and the HTF at the inlet region decreases and more heat transforms to the downstream region.As a result, the value of   / increases with the melting of the PCMs in the downstream region, causing the second peak in   /.Finally, since most of the PCMs are melted, more heat is continuously absorbed and stored in sensible form and the value of   / goes down gradually.
In order to evaluate the heat exchange unit comprehensively, the pressure drop is also investigated for different tooth length and tooth rate plates.Figure 15( tooth length increases from 12 mm to 39 mm.Nevertheless, from Figure 15(b) it can be seen that the pressure drop decreases from 51 Pa to 32 Pa while the tooth rate increases from 1/5 to 5/1.Therefore, the increase of tooth length or tooth rate can reduce the pressure drop effectively, enhancing the comprehensive performance of the heat exchange unit.

Effect of Different Working Conditions.
In order to investigate the effect of inlet velocity on the melting process, the temporal variation of the liquid fraction at different Re numbers with the same Ste number was calculated.The Reynolds number is generally defined as where   is the hydraulic diameter defined as   = 2  /( +   ) ≈ 2  , while  is the plate width inside the gasket and   is the mean spacing between plates., , and ] are the density, viscosity, and velocity of the HTF, respectively.Figure 16 presents the liquid fraction of the PCMs as a function of time obtained by employing different Re numbers.Obviously the slope of the curves is increasing while the Re numbers are increasing from 98 to 392.This results in a shortened melting time for the relatively larger Re numbers.This variation trend is attributed to two reasons.Firstly, higher Re corresponds to faster flow velocity which prompts the heat transfer.Besides, increasing Re number is accompanied with increasingly strong turbulence, which results in the HTF temperature being more uniform in each zigzag tooth.Secondly, the inlet heat energy and the temperature of HTF in the passageway also increase with the increase of Re number, enhancing the heat transfer.
As definition of Ste in (11), higher Ste brings about a greater temperature difference between the phase-changing region and the HTF.Thus, it is expected that the melting rate of PCMs increases with the increase of Ste number.This speculation is verified effectively through the results shown in Figure 17.A low Ste value means a low inlet temperature.This also means a small temperature difference between the HTF and the PCMs, and thus the heat transfer is lower at the smaller Ste values.Besides, Figure 17    implies that the heat transfer enhancement decreases with the Ste number.Figure 18 represents temporal variation of liquid fraction and outlet temperature under different conditions of the given fixed input power.
The input power can be expressed as is the reference temperature and the conditions of V in = 6 m/s,  in = 873.15K, and  ini = 663.15K are chosen as the base conditions.Five combinations of inlet condition at an identical input power are chosen, as shown in Table 3, and the impact of Ste and Re numbers on the energy storage is studied.
As shown in Figure 18, the variation trend is quite similar among the five outlet temperature curves, but overall the outlet temperature rises more quickly with small Ste number (i.e., larger Re number).The outlet temperature (  ) remains in low constant value in the early stage which is equal to the  ini and then rises slowly at moderate rates.When the liquid fraction approaches around 0.4, the   reaches a constant value which is close to the phase change temperature of the PCMs, but it can be seen that the onset of the constant temperature period becomes more delayed with decreasing Re numbers.At the end of the charging process, the outlet temperature rises quickly.The tendency of the outlet temperature can be explained as follows.At the beginning, the temperature of the PCMs is low and the heat is absorbed and stored in sensible form, leading to the constant   equaling  ini .At the phase change temperature of the PCMs, the material begins to melt while the heat is absorbed and stored in sensible and latent form, and due to the lower temperature of PCMs in the downstream region and the buffer effect of the melting process,  out rises slowly but is still lower than the melting point of PCMs at this stage.In the next stage, the temperature of PCMs in the downstream region rises and   is controlled by the phase change temperature of PCMs, which leads to   equaling the melting temperature of PCMs.Because most of PCMs melt and the heat is stored mainly in sensible form in the PCMs,   rises quickly at the end of the charging process.
Furthermore, the reaching time of constant   is delayed for smaller Re numbers.This is because the temperature difference between the HTF and the PCMs increases with decreasing Re numbers, which made more thermal energy get stored in the sensible way in the inlet region, the melting speed is lower in the downstream area and the reaching time of constant value is delayed.But the time duration of constant value decreases with increasing Ste number (decreasing Re number).Also, Figure 18 shows the liquid fraction as a function of time under different conditions for fixed input power.One can see that the overall melting time decreases with increasing Re number (decreasing Ste number), but the increase tendency slows down.This is because, for a given input power, an increase in the Ste number results in an increase in the inlet temperature and a decrease in the inlet velocity.The difference between inlet temperature and initial temperature increases with decreasing Re number according to (12).This implied that more heat is stored in the inlet region as sensible heat and less heat flows to the downstream region to melt the PCMs there, which leads to a lower melting rate in the condition of smaller Re for a given inlet power.

Conclusions
A numerical modeling of the melting process of PCMs within zigzag plate is performed using the solidification/melting model validated by experimental data.The validity of the model is investigated numerically in this study.According to this work, the following conclusions have been drawn: (i) A finite volume method with LES model is used in this study.The modeling results showed reasonable well agreement with the experimental data.
(ii) Zigzag shape can improve the heat transfer performance of the heat exchange plate compared to the smooth shape.
(iii) The tooth length and tooth rate almost have no impact on the heat transfer capability of the unit, but they both influence the pressure drop.Pressure drop decreases with the increase of tooth rate or tooth ratio.

Figure 1 :
Figure 1: The schematic of the LHTES unit.

2. 1 .
Physical Model.The configuration of the smooth plate (a) and zigzag plate (b) is shown in Figure 1.The HTF

Figure 6 :
Figure 6: Comparison between the measured and the calculated temperature as a function of time.

Figure 7 :
Figure 7: Liquid fraction as a function of time for the LHTES systems employing smooth and zigzag plates.

Figure 8 :Figure 9 :
Figure 8: Time variation of the liquid fraction at various zigzag lengths.

Figure 8 .
Figure 8.All cases are simulated under the condition of V in = 6 m/s and  in = 873.15K with the same amount of PCMs in the zigzag plates.It can be seen from Figure8that the four liquid fraction curves are almost coincident with each other over the entire charging process, implying that the zigzag length almost has no influence on the heat transfer capability of the LHTES system.It is well known that the heat transfer area and the flow distribution have great impact on the heat transfer capability of plate heat exchanger.With the increase of the tooth length, the heat transfer area decreases, but the adherence character of the fluid strengthens.As shown in Figure9, the mainstream region is near the middle part of the passageway in the plates of small tooth length, but as the tooth length gets bigger, the mainstream region gets closer to the plate wall, where in the latter case the heat transfer is enhanced.Thus, the joint effect from the above two aspects

Figure 10 :Figure 11 :
Figure 10: Schematic illustration of the shape of a tooth.

Figure 12 :
Figure 12: Schematic plot of flow field near the wall.

Figure 14 :
Figure14shows the ratio of the latent heat (  ) in the total heat (, composed of latent heat and sensible heat) for various zigzag lengths.As shown in Figure11, the value of   / rises quickly and gets to the first peak value at c.a. 500 s followed by a small decrease before c.a. 1000 s.This first   / peak can be attributed to the large temperature difference between the HTF and the PCMs at the beginning, resulting in the melting of the PCMs at the inlet region which rapidly absorbs and stores the heat of HTF as latent heat.After certain time, once the PCMs at the inlet region are completely melted, the PCMs continue absorbing heat and store it in sensible form leading to the decrease of   /.As the above process proceeds, the temperature difference between the PCMs and the HTF at the inlet region decreases and more heat transforms to the downstream region.As a result, the value of   / increases with the melting of the PCMs in the downstream region, causing the second peak in   /.Finally, since most of the PCMs are melted, more heat is continuously absorbed and stored in sensible form and the value of   / goes down gradually.In order to evaluate the heat exchange unit comprehensively, the pressure drop is also investigated for different tooth length and tooth rate plates.Figure15(a) shows the effects of zigzag tooth length on the pressure drop.The pressure drop decreases dramatically from 98 Pa to 49 Pa while the

Figure 15 :Figure 16 :
Figure 15: Pressure drop calculated for the heat exchanger unit.

Figure 17 Figure 17 :
Figure 17: The liquid fraction of the PCMs as a function of time during the charging process at different Ste number.
also indicates the change trend of the liquid fraction with increasing Ste number, which

Figure 18 :
Figure 18: The liquid fraction and outlet temperature under different Ste or Re numbers conditions for fixed input power.

Table 1 :
Thermophysical properties of the NaCl-MgCl 2 salt in the validation experiment.

Table 3 :
The working conditions at the same inlet power.
(iv) The melting time decreases with increasing Re and Ste numbers but decreases with increasing Re number (decreasing Ste number) when employing a fixed input power.Width of the inlet or outlet [m].