PERFORMANCE OPTIMIZATION OF PHASE CHANGE HEAT ACCUMULATOR WITH TURN-OVER TRIANGULAR TUBE

The performance optimization of the triangular-tube phase change heat accumulator is carried out from three aspects. The first is to compare with three different types tubes: the traditional round tube, the square tube, and the elliptical tube. Results show that equilateral triangle tube has the best heat storage and release characteristics. The second is to study the sub-cavity heat accumulator model based on multiphase change materials on the basis of the triangular tube. The third is to study the turn-over triangular tube heat accumulator. Both the flip function through discrete point fitting and the best reversal time are obtained.


INTRODUCTION
Phase change heat storage technology uses the characteristics of phase change materials (PCMs) to store/release heat. It has the advantages of simple device design, high heat storage density, and relatively constant heat storage temperature. However, the thermal conductivity of the phase change materials is low, so the structure optimization of the phase change heat accumulator must be carried out. As a result, many relevant research work has been developed, and they may be classified in two main categories: (i) the use of composite material with enhanced thermal conductivity, and (ii) the rationalized geometrical structure of the PCMs heat exchanger involving integration of finned tube and the employment of multiple melting temperatures.
Because the use of fossil fuels is limited, the development of the PCMs heat exchanger has attracted wide attentions in the past decades, due to a desire to improve the heat storage and release efficiency (Groll, 2020). A great deal of literature relevant to the enhancement of phase change heat storage has concentrated on the structural optimization of phase change heat accumulators. The simplest design of phase change heat accumulators is with rectangular or tubular tank containing the PCMs. Other configurations of the PCMs heat exchanger such as the tube heat exchanger with a bundle of tubes, shell-and-tube heat exchanger, and tube with external and internal finned surfaces, are also explored experimentally and numerically by various researchers (Kumara et al., 2021). Sciacovelli et al. (2015) studied the effect of Y-type ribs and Y-type bifurcation ribs on the thermal storage performance by numerical method. The results showed that Y-shaped ribs can effectively improve the heat storage efficiency. Ho et al. (2016) used the finite difference method to numerically study the fluid flow and heat transfer of a rectangular natural circulation loop that contains phase change material. Suspensions is presented to investigate how geometric parameters affect the thermal performance. The results show that in several of the geometric configurations, the heat transfer performance of the loop is significantly affected by the latent heat contribution associated with the melting/freezing of the PCM particles. Yuan et al. (2016) investigated the melting characteristics of PCM in an annular with different installation angle of fins ( θ ) via a numerical simulation method, which was based on an enthalpy-porous medium model to study the effect of fins installed at different locations on the heat transfer, and that results are verified via experimentation. Cao et al. (2018) investigated the effect of natural convection during melting process in eccentric horizontal shell and tube storage unit in order to improve melting rate of latent heat storage system. Air-side thermal and hydraulic performance of the finned flat-tube heat exchanger with single row is experimentally investigated by Zhang and Tang (2020).
Most recently, Mousavi et al. (2019) numerically studied melting thermal performance of a vertical cylindrical thermal energy storage unit in presence of a phase change material (PCM) including nanoparticles (nano-PCM) and also horizontal radial copper fins .The results indicate that the best result corresponds to the case with simultaneous use of nano-PCM (5% of Al2O3) and fins (three fins) that could improve the melting time up to 28.3%. Vogel et al. (2020) proposed a simplified and fast model for large-scale industrial storage systems. The model was implemented with a self-tailored MATLAB code and consists of two coupled parts: the heat transfer fluid region and the storage region, which includes the phase change material and the heat exchanger. Moreno et al. (2020) presented a transient numerical and experimental study to focus on describing the heat transfer in a cavity filled with air, which has one vertical wall with a phase change material (PCM). The wall with PCM was exposed to constant heat flux, whereas the opposite vertical wall was kept at a constant temperature. Elsanusi et al. (2020) investigated the performance of multiple PCMs for different arrangements in a horizontally positioned heat exchanger during melting of the PCMs. The results identified the effects of conduction and natural convection heat transfer for different PCM arrangements in the heat exchanger. Huang et al. (2020) conduct an experimental study on the thermal performance of the finned metal foam (FMF) heat sinks with phase change material. The dynamic temperature response of an FMF heat sink is analyzed and compared with the corresponding finned heat sink. Buonomo et al. (2020) numerically studied a latent heat thermal energy storage system with a highly conductive metal foam phase change material called Nano-PCM to enhance the heat transfer inside the system. Kukreja et al. (2020) performed an intensive numerical study inside the shell and tube type heat exchanger to find out the melting performance of a Phase Change Material. An axis symmetric virtual model of the thermal storage unit is created and various types of numerical schemes are implemented over it, to validate the suggested experimental results.
Several studies consider the effect of shell shape on heat transfer mechanism. Wu et al. (2020) found that all special-shaped tubes have better heating performance than ordinary round tubes for aircraft deicing fluid. Wang et al. (2020) used numerical simulation to study the falling film evaporator with 5 kinds of heat exchange tubes (round tube, flat tube, elliptical tube, dropper tube and egg tube). The egg-shaped tube is found to be the best film-forming tube among the five types of heat exchanger tubes. Yan et al. (2010) studied the heat transfer performance of the special-shaped tube in the heat pipe and found that the average heat transfer coefficient of the triangular heat pipe is the highest, and the overall synergy angle in the tube is large. Shahsavar et al. (2020) presented an experimental assessment on the effects of employing serpentine tubes with three different cross-sections of circular, triangular, and rectangular, on the characteristics of a photovoltaic/thermal (PV/T) unit in terms of energy and exergy efficiencies compared with a conventional PV system. Chen et al. (2019) developed one innovative heat-transfer enhancement tubevortex-shaped tube with the aid of roughing surface, and the geometrical model of the vortex-shaped tube, for the first time, is built by modeling the fabrication process based on the structural mechanics.
Among the studies that have considered different shell shapes few studies have concentrated on PCMs heat exchanger. Rabienataj Darzi et al. (2016) investigated numerically the melting and solidification of a PCM within three various horizontal annulus configurations including two circular cylinders, one elliptical cylinder in a circular cylinder and one finned cylinder in a circular cylinder in terms of the aspect ratio and the orientation of the ellipse and the number of fins. Faghani et al. (2018) simulated and analyzed the melting process between two elliptical cylinders. The shapes of the shell and the central tube are assumed to be circular or elliptical. Jourabian et al. (2020) studied constrained melting of ice as a PCM in inclined elliptical annulus. Zhao et al. (2019Zhao et al. ( , 2020 recently studied geometrical effects on energy storage tank with multiple tubes, and the fins are added by topology optimization method. In this paper, the special-shaped tube design is introduced into the structural optimization of the phase change heat accumulator, and the structure optimization of the PCM heat accumulator is studied by numerical methods. The research is mainly based on three points: First, the triangular tube heat accumulator is compared with three different types tubes. The second is the optimized model of the sub-cavity triangular-tube heat accumulator that can be filled with multiphase change materials. Three different melting temperature combinations of RT50, RT60, and RT35 are selected, and a combination of two filling forms is proposed. The third is to study the heat storage characteristics of the reversible triangular tube heat accumulator.

Physical model
The physical model diagram is shown in Fig. 1. The heat accumulator unit is a sleeve type, and model 1-3 is three comparison models for the triangular tube. Model 1 is a traditional round tube with an inner radius R in of 20 mm and an outer radius R out of 40 mm. The internal tube of model 2 is a square design, and the internal tube of model 3 is an elliptical tube design. Model 4 is an equilateral triangle designed with special-shaped tubes, with the sharp corners facing down. Model 5 is the same triangle with the sharp corners facing upwards (model 5 can be regarded as a flip model of model 4). Model 6 is based on model 4. The method of sub-cavity is optimized, and the same or different PCMs are filled in the sub-cavities for the model 6. The mass flow rates of the heat exchange fluids of all models are equal, and the cross-sectional area of the channel tube in the two-dimensional cross-sectional view is equal.

Mathematical model
Considering the influence of the natural convection on the melting of the phase change material, the continuity equations, momentum equations and energy equations should be established together in the cross section. The numerical simulation is performed by ANSYS Fluent which adopts the enthalpy-porosity method. The word "enthalpy" in this method is embodied in the treatment of the energy equation, that is, the enthalpy is directly used as the variable to be solved. Instead of tracking the liquid-solid front explicitly, phase change moving interface is determined after obtaining enthalpy and temperature fields for the enthalpy formulation. So the value of enthalpy of PCMs must be input into Fluent. The word "porosity" in the method essentially deals with the unity of the solid-liquid phase and solid-liquid interface in the momentum equation, and the source term in the momentum equation embodies the connotation of "porosity". In this way, the "two-phase fuzzy zone" of the phase change interface is considered as the porous medium zone whose porosity is equal to the liquid fraction. For a solid phase with zero porosity and zero velocity, the source term is infinite, which is equivalent to treating a solid as a liquid with high flow resistance, mathematically speaking, by the action of the source term. For the liquid phase, the porosity is 1, the source term is 0, and for the fuzzy zone, the porosity varies from 0 to 1, this leads to a decrease in momentum as the porosity decreases as the fuzzy zone solidifies. Therefore, through the action of the source term, the momentum equation which can be applied to the liquid phase region is applied to the solid phase and the "two-phase" region, thus the equation is unified and convenient for numerical calculation. Because the temperature of pipe along the axial direction is almost constant and the effect of straight pipe at the exit and bend pipe can be ignored, the physical model can be simplified as twodimensional heat transfer problem and the following assumptions are made: (1) The subcooling effect of PCMs when melting is ignored; (2) The laminar flow model is used in the phase change process of PCMs, which is a low-speed non-steady incompressible flow process; (3) PCMs and tube wall material are homogeneous and isotropic; (4) The temperature of pipe wall is constant in heat release and storage process; (5) Thermal conductivity on the thin outer wall of storage tank is ignored.
Based on the above analysis, the governing equations can be written as follows: Continuity equation: Where v  is the velocity vector; t is the time; ρ is the density of PCM; p is the pressure; μ is the dynamic viscosity of PCM; α is the expansion coefficient of PCM; g  is the gravity acceleration; T is the temperature of PCM; T ref is the reference temperature; λ is the thermal conductivity of PCM; C p is the specific heat of PCM; and H is the total enthalpy of PCM, which can be expressed as the sum of the sensible heat enthalpy (h) and the latent heat enthalpy (Δh) : h ref is the reference enthalpy; Γ p is the latent heat; β is liquid fraction that can be expressed as follows: Where T s is the solid temperature; T l is the liquid temperature. Some PCMs have same value; g S  is added into the momentum conservation equation as a source term, which is defined as: A mush is the fuzzy zone constant, generally, the value is between 10 4 -10 8 . The higher the value is, the faster the speed of solidification of the phase change material will fall, and the wave and divergence of the solution will be easily caused. The value is chosen 10 5 . ε is the constant 0.0001, the point of this value is to provide a decimal number to avoid the denominator of 0. Boundary conditions are set as follows: the inner boundary T W is 353.15K, the initial temperature T init is 293.15K, and the outer boundary is an adiabatic boundary condition.

NUMERICAL PROCEDURE AND CODE VALIDATION
The Solidification & Melting model of ANSYS Fluent based on the enthalpy-porosity method is introduced to deal with the material phase change process for numerical simulation, and the change physical parameters such as the density are realized by UDF. The primary physical properties of PCM are listed in Tab. 1. The velocity and pressure solution method is the unsteady implicit SIMPLEC algorithm, the pressure term is discrete using the PRESTO method, the convection term is discrete using the QUICK format, and the energy is the first-order upwind style. The pressure, density, velocity, liquid fraction, and energy sub-relaxation factor are set to 0.4, 1, 0.7, 0.4, and 0.7, respectively. The grid-independent solution has been performed for the cases of thermal performance of model 1(round tube design), the non-uniform structured grid is used to divide the calculation area. Three grid systems with the 1.1 inward growth rate in the radial direction, encrypting the boundary wall of the inner near wall (including ribs) are used, and the number of grids is 22,990, 44,551 and 59,501(see Fig. 2).
(a) Average temperature at different grid numbers (b) liquid phase fractions at different grid numbers Fig. 2 Grid size assessment As seen from Fig.2, the maximum calculation error of liquid fraction and PCMs average temperature is 0.11% and 0.74%, respectively. Therefore, the grid number 44,551 is retained after tradeoff between accuracy and computational cost.
The time-step-independent solution is obtained by trial calculation. Three different time steps, 0.2, 0.5, and 0.7 s, are considered for the case of thermal performance of model 1 (round tube design) at the grid number is 44,551. As seen from Fig. 3, the maximum calculation error of liquid fraction and PCM average temperature is 0.43% and 0.52%, respectively. Therefore, the time step is determined to be 0.5 s.
(a) Average temperature at different time steps (b) liquid fraction at different time steps

Fig. 3 Time step assessment
The validation of the computer code is performed by comparing with the experimental results provided by Yuan et al. (2016) and Cao et al. (2018). The number of Fo is introduced, which is defined as follows: The ∂ is the thermal diffusion coefficient, which is equal to λ ρ /( ) p c . τ is the time used for calculation. D is the equivalent diameter, which is equal to 2R out -2R in. The liquid fraction calculated by model 1 is compared with the experimental data in reference Yuan et al. (2016). The comparison result is shown in Fig.4(a), and the maximum data error is 3.7%. The liquid phase diagram after post-treatment is compared with the experimental data from reference Cao et al. (2018), and the result is shown in Fig.4(b). The results both show that the simulated data and the experimental data are in good agreement.

Comparison of heat storage and release processes of Model 1, 2, 3, 4 and 5
Numerical simulation of the heat storage process for comparative analysis is performed. The simulation results of Model 1, 2, 3, 4 and 5 are shown in the figure 5 and figure 6 below. As can be seen from Fig. 5 and Fig. 6, the melting process of the phase change material in the heat storage process for all models is mostly similar. The heat storage process can be divided into the following three stages. Stage 1: A thin layer of liquid phase PCMs fluid is generated on the heating wall(see Fig. 5), and the liquid phase area increases in the radial direction. Stage 2: Due to buoyancy, the liquid PCMs gathers on the upper part of the heat accumulator. At this time, the natural convection is very strong(see Fig. 6). The overall temperature in figure 6 decreases from top to bottom. After the upper part of the PCMs is completely melted, the liquid phase area diffuses downward from both sides. Stage 3: Due to the low thermal conductivity of the phase change material, when the natural convection is further reduced, the heat is difficult to diffuse, and the liquefaction speed of the lower part of the phase change area slows down, and the bottom hard-to-melt area is formed in model 1. Since the bottom of model 4 has the shortest relative distance, the lower area melts fastest(see Fig. 5(b)), and the other special-shaped tubes also melt faster than round tubes.
By analyzing the numerical simulation results of the exothermic process of model 1, 2, 3, 4 and 5, the liquid fraction and temperature field are shown in Fig. 7 and Fig. 8.   figure 7(a), the liquid PCMs around the inner wall surface heats up to form a solid thin layer. At this time, the PCMs around the inner wall mainly transfers heat in the form of heat conduction, and the solid thin layer grows in the radial direction. Affected by gravity, solid materials gather in the lower part. Due to the low thermal conductivity of PCMs, it is difficult for heat to pass through the solid PCMs to the cooling wall. It can be obtained from Fig.  7(b) that the melting times of models 1, 2, 3, 4, and 5 are 8244s, 6139s, 5848s, 3051s, 7371s, respectively. Compared with the round tube, the heat storage time of models 2, 3, 4, and 5 are shortened by 22.8%, 29.1%, 62.9% and 9.44%, respectively. Comparing Model 4 and Model 5, the efficiency of the two is almost the same in the exothermic process, indicating that the direction of gravity has little effect on the heat release efficiency when the triangular heat accumulator releases heat. The time for model 1 to completely solidify after 27000s, the time for model 2, 3, 4, 5 to completely solidify are 409.7min, 436min, 75min and 74.6min, respectively. For model 1, 2 and 3, a hard-to-solidify area is formed on the upper part.
(a)Liquid fractions (b) Average liquid fraction Fig. 7 Exothermic processes Figure 8 shows that the temperature gradient is mainly reflected in the radial direction during the heat release process, which is consistent with the above analysis. It can be concluded that the triangular tube significantly enhances the utilization efficiency of energy storage.

Cavity model based on triangular tube
Based on the above research on the smooth-tube heat accumulator, the model 4 is evenly divided into six equal-area regions with six fins, the fin thickness is 1mm, and they are named S1, S2, S3 from top to bottom. RT50, RT60 and RT35 are chose to propose a combination of different filling forms (listed in Tab. 2), the model 6 (sub-cavity model)with two different configurations is obtained.  Figure 9 shows the liquid phase diagram during heat storage process for the two configurations. It can be seen from Fig. 9 that the two configurations which adopt the radial fin sub-cavity design both significantly improve the heat storage efficiency. Due to the effect of buoyancy, the liquid PCMs flows upward along the walls and fins. In the S2 and S3 regions, compared with Fig.5 (a) of model 4，the downward diffusion speed of the liquid region is significantly improved.

Fig. 9
Liquid phase diagram during heat storage process of model 6 Figure 10 shows the temperature diagram during heat storage for the two configurations. According to Fig.10, the temperature field initially distributes in the radial direction and after a short time of heat exchange, the temperature field finally distributes in the vertical direction, which is consistent with the conclusion of the finless pipe. The velocity field covering the S1, S2 and S3 regions are shown in Fig.11. It can be observed that the fins have a very obvious drainage effect on the liquid PCMs. In each cavity area, the liquid PCMs first flows upward along the heating surface, and then returns to the lower part of the heating wall around the cavity. Therefore, the upper part of each cavity has a higher temperature than the lower part, and before the PCMs in the heat accumulator is completely melted, the disturbance of natural convection has always affected the distribution of the entire temperature field. It can be concluded that in the sub-cavity during the heat storage process of the special-shaped tube, natural convection and heat conduction take effect at the same time.  Fig. 12, it can be observed that the heat storage efficiency of the two configurations are both significantly improved. In the S2 and S3 regions, the diffusion speed of the liquid region is significantly higher than that in figure 5(b). The average temperature of each regions of model 6 are shown in Fig. 13. It can be seen from Fig. 13 that the heating rate of configuration 1 and configuration 2 in each area is almost the same, but in S3, the melting speed of configuration 2 is greatly improved, which is close to the melting speed of S1 and S2. In summary, the thermal storage time of configuration 1 is shortened by 27.9% compared with model 4, and configuration 2 is shortened by 46.47% compared with model 4, and 25.7% compared with configuration 1. Configuration 2 basically meets the requirement of simultaneous melting of the overall PCMs, which improves the heat storage efficiency. Figure 14 shows the average liquid fraction during the exothermic process of configuration 1 and 2. The time for complete solidification 7 of program 1 and 2 is 2362s and 2076s respectively. Compared with model 4, the time has been shortened by 47% and 53.9% respectively. Configuration 2 is 12.1% shorter than configuration 1. The sub-cavity model also has an optimized effect on the heat release process of the heat accumulator, and configuration 2 is better. Fig. 14 The average liquid fraction of model 6 program 1 and 2 during the exothermic process

Fitting of flip-type phase change heat storage
The design of the flip model: when the upper part of the PCMs is basically melted, and the lower part is solid, the original model is flipped 180° (that is, from sharp corners downwards to sharp corners upwards) for using natural convection to melt the originally slower melting area, thus it can form a velocity field similar to model 5 at the initial stage of melting, and shorten the melting time.
After a large number of numerical calculations at different turningover time, it is found that within 660s~780s, the PCMs complete melting time has a minimum value. The following Tab. 3 shows the final melting time at different turning-over times. Table 3 Final melting time at different turning points The function fitting curve is similar to a parabola with an opening upward, and with the lowest point between 700s and 720s. Take the data in Tab. 3 to perform polynomial interpolation fitting, and the fitting curve is shown in the following Fig. 15.

Fig. 15 Polynomial fitting curve diagram
The fourth degree polynomial fitting can best show the influence of the flip time on the final melting time. From the calculation of the fitted curve, we can obtained that the minimum value of y is 1088s, and x is 708s (after take an integer), that should be the best flip time.
The fitting curve equation obtained according to the data is:  (10): x is the turn-over time, and y is the total time for the final melting.
The correlation coefficient R of the equation is 0.99806. Set the turn-over time to 708s, perform numerical simulation to verify the equation, and the result is shown in the following Fig. 16 and Fig. 17.  with an error of about 0.18%, which proved the reliability of the fitting polynomial From Fig. 17, after the 708s reversal, the flow direction of the liquid phase region has been reversed, the liquefaction rate has a slight increase, and the Nu number has increased. In the overturned heat accumulator, the high-temperature liquid PCM is mainly concentrated in the lower half, and the liquid PCM is rapidly gathered upwards driven by the heat. After the flip, the average temperature drops first and then quickly rises to the level before the flip and accelerates the temperature rise. Finally, the PCMs melt completely in 1090s, shortening the melting time by 64.27%.

CONCLUSIONS
The performance optimization study of the triangular-tube phase change heat accumulator is carried out and the conclusions are drawn as follows: (1) Compared with the traditional round tube heat accumulator, models 2, 3, 4, and 5 shorten the heat storage time by 22.8%, 29.1%, 62.9% and 9.44% respectively. The triangular tube has the highest heat storage efficiency.
(2) Both the configuration 1 and 2 of the sub-cavity model (model 6) can significantly improve the heat storage and release efficiency, and the configuration 2 is better. The thermal storage time of configuration 1 is shortened by 27.9% compared with model 4, and program 2 is shortened by 46.47% compared with model 4, and 25.7% compared with configuration 1. During the exothermic process, the time taken for the complete solidification of program 1 and 2 is reduced by 47% and 53.9% compared with model 4. Configuration 2 is 12.1% shorter than configuration 1.
( x is the turn-over time, and y is the total time for the final melting.