A review of phase change heat transfer in shape-stabilized phase change materials (ss-PCMs) based on porous supports for thermal energy storage

Latent heat thermal energy storage (LHTES) uses phase change materials (PCMs) to store and release heat, and can effectively address the mismatch between energy supply and demand. However, it suffers from low thermal conductivity and the leakage problem. One of the solutions is integrating porous supports and PCMs to fabricate shape-stabilized phase change materials (ss-PCMs). The phase change heat transfer in porous ss-PCMs is of fundamental importance for determining thermalfluidic behaviours and evaluating LHTES system performance. This paper reviews the


Introduction
Renewable energy has been attracting increasing interest in the past decades due to its significant potentials for addressing growing energy demand and greenhouse gas emissions. According to "World Energy Outlook 2018" by International Energy Agency (IEA), the growing world economy and newly added 1.7 billion people are projected to push up global energy demand by 1/4 to 2040. And more energyrelated CO2 will be emitted, which is estimated to exceed efforts taken to tackle climate change [1]. Renewable energy provides an attractive solution to alleviate growingly global energy demand and achieve a low-carbon future.
Among various renewable energy sources, solar energy continues to be a promising candidate to produce thermal energy for domestic applications, industrial manufactures and buildings. However, there is a temporal and geographic mismatch between thermal energy supply and demand. As a result, thermal energy storage (TES) is proposed to address this problem [2].
Currently, there are mainly three thermal energy storage technologies: sensible heat storage, latent heat storage and thermochemical reaction storage [3]. Table 1 summarizes the principle, advantages, typical materials and application areas of these TES technologies. As shown in this table, latent heat materials, or phase change materials (PCMs) has advantages of high energy storage density, high latent heat and the ability to maintain an almost constant temperature, and thus it is most widely used [4]. These advantages contribute to not only reduce equipment ability required and cost but also improve thermal storage performance. As a consequence, PCMs are most investigated among these three TES technologies [5]. Table 1 Principle, advantages, typical material and application areas of three TES technologies [6,7]  Nevertheless, pure PCMs suffer from two problems: low thermal conductivity (≈ 0.1W m -1 K -1 ) and leakage, which limits its utilization in many sectors [8]. For example, in some electronic devices, the chip generates heat transiently or periodically, which requires an efficient cooling system to dissipate heat. However, the poor conductivity of pure PCMs decreases the heat transfer rate and increases the possibility of chip exposure to an extremely high-temperature environment [9]. One of the solutions is integrating porous supports with PCMs to fabricate shape-stabilized phase change materials (ss-PCMs) [10]. According to the pore size, porous supports are classified as macroporous (>50nm), mesoporous (2~50nm), microporous (<2nm) and hierarchical porous (ranging from marco to micro) materials [11]. The widely-used support materials are shown in Fig.1. Different types of support materials have varying properties and functions: metal foams, a typical macroporous material, possess high thermal conductivity and thus are implemented as heat delivery promoter; due to the small pore and large surface area, meso-and microporous materials show a strong guest-host interaction and therefore avoid leakage; in hierarchical porous materials, macropores act as the cavity to store PCMs, mesopores provide transport pathway and micropores give capillary force to immobilize PCMs [12]. Porous support materials of ss-PCMs: graphite foam [13], mesoporous silica [12], carbon nanotube (CNT) [14], metal-organic framework (MOF) [12], porous coordination polymer (PCP) [15], hierarchical porous polymer (HPP) [16], hierarchical porous carbon (HPC) [17].
Many researchers have reviewed the progress on porous ss-PCMs. Rehman et al. [18] reported a review on heat transfer augmentation of PCMs using porous metal foams and carbon materials. Kenisarin et al. [19] introduced the fabrication of porous ss-PCMs with expanded materials such as perlite and vermiculite in their review paper.
Umair et al. [13] reviewed researches on fabrication of porous shape-stabilized organic PCMs. Huang et al. [8] summarized the key studies on porous materials used as the supports of porous ss-CPMs. Feng et al. [11] presented a review on fabrication, characterization, enhancement and molecular simulation of nanoporous ss-PCMs.
Gao et al. [12] summarized recent studies on nanoconfinement effects on thermal properties of nanoporous ss-PCMs. Zhang et al. [5] reported a review on porous ss-PCMs fabricated by using metal foams and carbon materials, from the aspects of production, characterization, application as well as mathematical models describing phase change heat transfer in the composite PCMs. Phase change heat transfer is of remarkable importance for determining the thermal-fluidic behaviours and evaluating the performance of LHTES systems [20]. A good understanding of flow and heat transfer characteristics during the phase change process helps to realise the interactions between porous supports and PCM cores, and thus to design LHTES systems with higher loading, longer durability and higher thermal performance. A great number of researchers have carried out studies and significant advancements have been achieved on phase change heat transfer in porous ss-PCMs. However, to the best of our knowledge, there are very limited reviews on the solid-liquid phase-change heat transfer in porous ss-PCMs [5]. Since the phase change heat transfer in porous ss-PCMs has attracted increasing interest, this paper reviews the recent experimental

Experimental study on solid-liquid phase-change heat transfer in porous ss-PCM
The major advantage of experimental studies is that they are able to provide "directly interpretable" and reliable results [21]. In addition, data collected from experiments provide the validation source for numerical studies. Hence, a great number of experiments are performed by researchers to study the solid-liquid phasechange heat transfer in ss-PCMs based on porous supports.

Materials, methods and experimental apparatus in
phase-change heat-transfer investigations 2

.1.1 Materials and methods
To best of our knowledge, the first experimental investigation on the solid-liquid phase-change heat transfer of porous ss-PCM was conducted by Weaver and Viskanta in 1986 [22]. In their research, water and glass beads were employed as PCM and porous support, respectively. Since then, it has been over 30 years for investigations on phase-change heat transfer of porous ss-PCMs, and materials have changed a lot.  The metal foam is the most used porous support in the experimental studies.
There are mainly two methods to fabricate metal foam based ss-PCMs. The first one is the directly pouring method, which is adopted by Righetti et al. [24], Mallow et al. [25], etc. The procedure of this method is that the hot liquid PCM is directly poured into a container and mixed with the porous support. This approach is simple and does not need extra equipment, however, due to the air existing in the porous support, liquid PCM cannot infiltrate the support completely. Thus more researchers employed the vacuum impregnation (or vacuum infiltration) method, such as Zhang et al. [20], Jin et al. [26], etc. A typical flow chart is presented in Fig. 3. The procedure consists of six steps: firstly, the solid PCM, metal foam and a mesh which is used to support the foam are placed into a container. A vacuum pump is connected to the container and switched on to evacuate the air. Secondly, the container is heated and the metal foam sinks into the molten PCM. After the porous support is fully impregnated by liquid PCM, the heating process is ended and the vacuum pump is switched off at the same time.
Then the container is cooled. When the PCM is completely solidified, the container is reheated slightly to make it easy to withdraw the specimen. Finally, the composite PCM is taken out and the surplus PCM on the surface is removed.

Experimental apparatus
The experimental apparatus is relatively simple in early researches. For example, Weaver et al. [22] just employed a conventional camera to study the freezing process of water. Later, more and more high-tech measurement devices are utilized to investigate the solid-liquid phase-change heat transfer in porous ss-PCMs. Table 3 gives a summary of the measurement device and type of measured data in literature from 1986 to 2019. An example of experimental apparatus is also shown in Fig. 4.  [42] Generally, temperature and solid/liquid interface are two commonly monitored parameters. To record the temperature variation, the thermocouple is used by most researchers, as Table 3 shows. Thermocouples are inserted into the ss-PCM and linked with data acquisition. Hence, the temperature at selected points can be monitored. The number of thermocouples varies in different researches. To our best knowledge, Martinelli et al. [59] employed the maximum number of thermocouples, 40, to investigate the phase-change heat-transfer process in a 400mm tube-and-shell test rig. Most researchers utilized thermocouples to monitor the temperature variation of PCM. However, since there is an apparent difference in thermal-physical properties between porous supports and PCM, the temperature of the two components may be not the same [60]. To explore this thermal phenomenon, Zhang et al. [20] built a set of experimental apparatus where two T-type thermocouples were set into the ss-PCM.
One was inserted into PCM in the support material pore to record the PCM temperature, while the other was soldered on the metal skeleton surface to measure the skeleton temperature. This experimental setup is first built by Zhang et al. [20] and it can be employed to experimentally investigate the local thermal non-equilibrium and provide data for the comparison with simulation results. However, the thermocouple has an intrinsic disadvantage: it only records the temperature variation at certain points, rather than the whole temperature field. To solve this problem, Jackson et al. [32], Zhang et al. [20] and Yao et al. [55] applied the infrared camera to capture the evolution of the whole temperature field.
With regard to the phase interface, the conventional camera was widely utilized to capture the solid/liquid interface. Some researchers also employed a high definition (HD) camera to obtain higher quality images [36,41,48,55]. After images of phase interface being processed, some other parameters, such as melting/solidification fraction and rate, can be obtained. It is noteworthy that Chen et al. [31] employed an optical microscope to study the melting evolution of PCM at the pore scale, which provided new insights into the measurement of phase interface.  (1) conduction-dominated stage: the melting interface was parallel to the heating wall and the thermal energy was transferred to solid PCM in the form of sensible heat ( Fig. 5(b)).
(2) conduction-convection-mixed stage: liquid PCM flowed upwards due to the buoyancy force and a wide circulatory region was formed in the upper part, curving the interface; in the lower part, the interface was still vertical due to the conduction  melting fraction) [42] Diani et al. [45] used three kinds of paraffin with different melting temperatures as PCM and aluminium foam as porous support to study the phase change heat transfer in a 20mm(length) x 100mm(width) x 100mm(height) cavity. Venkateshwar et al. [46] investigated the melting process of aluminium foam/ coconut oil with CuO nanoparticles. The solid/liquid interface evolution obtained by these researchers is similar to that by Al-Jethelah et al. [42].
It should be noted that the above researches are under the condition of the left heating wall. Different heating position may lead to different solid/liquid interface. To explore the mechanism, Zheng et al. [36] conducted an experiment where the paraffin/copper foam ss-PCM was heated by the left, bottom, and top wall, respectively.
The interface propagation in the case of left heating wall was similar to Al-Jethelah et al. [42]'s results, while the other two cases were different. For the condition of the bottom heating wall, the melting front at the early stage was parallel to the heating boundary. However, at 4.5h, the PCM in the middle melted faster than on both sides because of the natural convection, and two symmetrical annular flows occurred. For the case of top heating, the phase interface was always parallel to the heating wall because the influence of natural convection was relatively insignificant.
The cylindrical container is another widely used enclosure for porous ss-PCMs [61]. In 2008, Siahpush et al. [47] investigated the melting/solidification process of eicosane/copper foam composite PCM in a 155.5mm(inner diameter) x 304.8mm(height) copper tube. During the test process, the outer wall temperature of the container was kept constant. 81 thermocouples were employed to monitor temperature variation. It was found that the curvature of solid/liquid interface in the case of metal foam was not as pronounced as the case of without metal foam because the metal foam enhanced the effective conductivity of ss-PCMs. Yang et al. [48] studied the phase change process in a tube-and-shell unit. Copper foam with/without a bottom fin was compounded into the paraffin to enhance heat transfer. For the PCM with copper foam, the interface was cone-shaped and developed from the inside to outside, while for the case with copper foam and bottom fin, the melting interface was inversed funnel-shaped. Recently, Yang et al. [41] performed a visual experiment to investigate the melting process of PCM/metal foam in a tube-and-shell unit. The solid/liquid interface propagation is presented in Fig. S1. It was found that at the initial stage (60min), the inner interface was vertical, indicating that the conduction dominated the heat transfer. With the elapse of time, more PCM melted and the natural convection remarkably contributed to the heat transfer. The hot molten PCM was pushed upwards by buoyancy force and accelerated the phase change in the upper region, thereby forming a funnel-shape interface after 60min. For the outer solid/liquid interface, it was horizontal during the whole phase change process.

Temperature distribution and variation
Thermocouples are widely utilized to record temperature distribution and variation during the phase change process. For instance, Zhou et al. [49] arranged four thermocouples at different positions in a rectangular cavity. They found that by adding copper foam into paraffin, the heater temperature was reduced dramatically; during the melting process, the temperature of PCM with metal foam was higher than that without metal foam, implying that the heat was conducted to PCM faster with the assistance of metal foam. Mancin et al. [52] investigated the solid-liquid phase-change process of paraffin/copper foam composite ss-PCM. A lower surface temperature was observed compared with no-foam case. Wang et al. [53] studied the paraffin/aluminum foam composite ss-PCM in a Li-on battery. Their results showed that the existence of aluminum foam improved the temperature uniformity of PCM. Zhang et al. [23] experimentally investigated the phase change process of eutectic salt in a metal foam.
Thermocouples were fixed along the axis and wall of a cylindrical container (see Fig.   6(a)). During the test procedure, the container was heated in the thermostatic oil bath.
The temperature variation is shown in Fig. 6(b). For the pure PCM, the temperature at point D was higher than that at point B and C, while for copper foam/salt composite, the temperature at point B was highest. This is because the natural convection dominated the melting process of pure PCM while it was suppressed in the case of copper foam/salt composite. process [23] Recently, the infrared camera was used by some researchers to capture the entire temperature field. In Zhang et al. [20]'s experiment, the cavity was heated by the left wall and the results was presented in Fig.7. It can be seen that at the initial stage, the temperature contours were almost vertical, indicating that the heat conduction dominated the heat transfer. With the elapse of time, more PCM melted and the local natural convection became gradually notable, which accelerated the interface propagation in the top domain. The trend of temperature fields was similar to that of the solid/liquid interface in Fig. 5. In Fan and Jin's research, an infrared camera was equipped to record the temperature field during the melting process [51]. It was found that the temperature difference between the porous support and PCM became more pronounced as the phase change proceeded. Yao et al. [55] investigated the melting process of paraffin/copper foam with the assistance of an infrared camera. Their results showed that the copper foam improved the temperature distribution uniformity and increased the melting rate by 2 times.

Comparison with pure PCM
The temperature field of PCM with/without metal foam is presented in Fig. 8. It can be seen that after adding the metal foam, the heat can be transfer effectively to the region away from the heating source and the temperature distribution is more uniform. It is generally believed that the addition of porous media promotes the heat conduction while it hampers the natural convection [18,53]. The effects of the two aspects on the phase-change process are the opposite. Table 4 lists outcomes of adding porous media into PCM. All the researches illustrate that the phase change heat transfer is enhanced by porous support, indicating that the improvement of heat conduction is more prominent than the suppression of natural convection.  The effective thermal conductivity was increased from 0.423 W/ (m·K) to 3.06 W/(m·K); Consuming time was reduced from 375min to 85min for solidification process, and from 500min to 250min for melting process. [24] Paraffin/ aluminium foam od: 62mm, h:800mm Tube The time of melting/solidification was shortened by more than 2h.

Effect of porous support configuration
The porous support configuration, i.e. porosity and pore density, has a remarkable effect on phase change heat transfer characteristics of porous ss-PCMs. Many researchers conducted experiments to explore it. A summary of relevant studies is given in Table 5. The heat transfer enhancement is characterized using different parameters, such as melting/solidification rate, temperature distribution uniformity, decrease in base temperature, etc. All the researches illustrate that the decrease in porosity enhances the phase change heat transfer because a small porosity indicates a large high-thermal-conductivity volume in porous ss-PCMs, thereby improving the whole heat conductivity capacity of ss-PCM.
In contrast, the effect of pore size varies in different researches. As Table 5 shows, the small pore size improved the thermal performance of porous ss-PCM in Zhao et al. [30] and Mallow et al. [25]'s researches, while it had the negative effect in Li et al.
[9] and Allen et al. [56]'s studies. The influence of pore size is attributed to two aspects: on the one hand, the decrease of pore size enhances the thermal conduction because it leads to a larger interactive surface area between PCM and support material; on the other hand, it hampers the natural convection of molten PCM because smaller pores limit the motion of liquid PCM. The two aspects compete and the final result is determined by the more significant one. The mechanism can be explained by Jin et al.
[26]'s experiment. At 20℃wall superheat, the influences of the two aspects are equal, and thus the melting rates are the same. When the wall superheat is increased to 30℃, the natural convection is more intensive and the confine on convection exerted by small pores is more obvious. Hence, the melting rate in the case of 50PPI becomes smaller than that of 30PPI.

Numerical investigation on solid-liquid phasechange heat transfer in porous ss-PCM
Although experimental investigations can provide "directly interpretable" results of the phase change process, they are time-consuming and cannot depict some detailed flow and heat transfer characteristics, such as the flow field. Besides, it would be very expensive to experimentally investigate all the parameters influencing the phase change process of porous ss-PCMs [21,64]. The numerical simulation provides an effective solution to address the problems and thus, a great number of researchers conducted numerical investigations [61].
Generally, the simulation methods for solid/liquid phase-change heat transfer in porous ss-PCMs can be classified into two categories: representative elementary volume (REV)-scale method and pore-scale method [65]. The schematic of the two methods is presented in Fig. 9 and the comparison is listed in Table 6. The REV-scale simulation treats the porous ss-PCM as a uniformly mixed medium (see the typical computational domain in Fig. 10(a)) and does not require an accurate description of the support structure. In other words, it ignores the complex geometric information of the medium. Instead, this method utilizes some statistical parameters, such as porosity, permeability, and effective thermal conductivity, to characterize the porous structure.   Comparison between REV-scale and pore-scale simulation

REV-scale simulation
According to the local thermal/non-thermal equilibrium, models used in the REVscale method are classified into one-temperature model and two-temperature model [5]. The one-temperature model assumes there is a thermal equilibrium between PCM core and porous support. This model was first applied to investigate the melting of ice in the porous medium by Weaver et al. [28]. Due to its simple computation, several researchers continued to employ the one-temperature model to perform simulations on the phase-change heat transfer in porous ss-PCMs. The relating studies are listed in Table 7.
The one-temperature model is applicable when the thermal equilibrium between the PCM and porous supports can be achieved. However, since there is an apparent difference between the thermo-physical properties of the PCM core and support material, the macroscopic temperature of the two components is not always the same [70]. Hence, the one-temperature energy model is not valid in many cases. In contrast, the two-temperature model takes into account the temperature difference between the PCM and porous support, and thus it is more accuracy [9]. As a result, it is widely employed by researchers. This section will introduce the twotemperature model and its applications in simulating phase-change heat-transfer process of porous ss-PCMs. Forchheimer-extended Darcy model [76] and the Brinkman-extended Darcy model [77]. Taking into account the above models, the generalized non-Darcy momentum equations for liquid PCM are given by [5] where , and is velocity in , and direction, respectively; the second term on the right side of above three equations accounts for the Darcy effect; the third term explains Forchheimer-extended Darcy effect; the fourth term represents the Brinkmanextended Darcy effect; the fifth term in Eq.(5), ( − 1 ), denotes the natural convection driven by temperature difference; the last term is the superficial velocity source term to discriminate the solid-liquid region [78]. is defined as in which ℎ is the mushy constant and commonly valued between 10 3 and 10 9 [79,80]. ∅ is set to 0.01 in order to prevent the division by zero.
It is noteworthy that Eqs.

Energy equation
In the two-temperature model, the energy equations for PCM and porous support are formulated separately.
For PCM [9,78]: For porous skeleton: where is thermal conductivity induced by dispersion; and are the effective thermal conductivity of PCM and the solid support, respectively; ℎ is interfacial heat transfer coefficient; is interfacial surface area. The definitions of ℎ and can be found in [78]. The last term, − , in Eq. (7) acts as the source term to calculate the liquid fraction . It denotes that the phase change is taken into account at each iteration. Readers could refer to Ref. [84] to find the detailed calculation method of .
When the local thermal equilibrium is assumed, the two-temperature energy equations reduce to the one-temperature energy equation [83]: Where ̅̅̅ and are volume-averaged parameters and expressed as ̅̅̅ = (1 − ) + (10) There are also some different forms of the two-temperature model. For example, in order to simulate the melting process in the paraffin/metal foam composite ss-PCM, Zhao et al. [30] proposed the following two-temperature energy equations: Compared with Eqs. (7) and (8) where ℎ is the volumetric heat transfer coefficient; is the thermal conductivity of PCM. In the above formulation, the effect of PCM thermal conductivity was characterized as , rather than + in Eq. (7). Hu et al. [71] employed the twotemperature model to simulate the phase change process in a micro-foam. In their model, the term + in Eq. (7) was replaced by .

Correlations of interfacial heat transfer coefficient
Effective thermal conductivity and interfacial heat transfer coefficient are two critical parameters for modelling phase change heat transfer in porous ss-PCMs.
Some researchers have summarized correlations of effective thermal conductivity in their reviewer papers [5,18]. Thus, in this section, we reviewed correlations of interfacial heat transfer coefficient used in numerical investigations on phase change heat transfer in porous ss-PCMs.
Interfacial heat transfer coefficient, ℎ , describes the thermal interaction between ligaments of the porous support and PCM. However, it is very difficult to measure the coefficient by experiments due to the complex support structure and complicated flow scenario in porous ss-PCMs. Therefore, most researchers used empirical correlations as a substitution. The commonly used correlations are listed in Table 8. Žukauskas In this correlation, the effect of porosity is taken into account. The formula is adopted by some researchers to estimate ℎ , such as Chen et al. [86] and Zhu et al. [82]. Di Giorgio et al. [93] proposed a new method to directly calculate ℎ by using the Kelvin 3D model of the foam structure. They divided the melting process into three phases: solid-solid, melting and convective and thus three coefficients were determined. They found ℎ did not vary remarkably in the solid-solid and melting conditions and therefore they assumed ℎ was constant in the two phases. Later Caliano et al. [38] employed this approach in their numerical study on the cold thermal energy storage unit. simulations using different numerical methods is given in Table 9. In addition, a statistics of numerical tools used in REV-scale simulation is presented in Fig. 11. The earliest numerical simulation was performed by Weaver and Viskanta in 1986. They utilized FDM to solve the model. Nevertheless, it is obvious from Fig. 11 that FVM is the most commonly used method. This may be attributed to the rapid development of commercial CFD software, such as Fluent because Fluent is based on FVM. As Fig.   11(b) shows, Fluent is widely used for simulation, promoting the application of FVM.

Applications
A summary of REV-scale simulations based on the two-temperature model is given in Table 9. It covers literature from 2000 to 2019. Fig.12 gives the statistics of types of porous supports and pore size. It is found that most support materials are metal foam, and the pore size is over 0.1mm. In other words, all ss-PCMs studied are macroporous.   considering natural convection [68] To compare effects of heat conduction and natural convection on the phasechange process, Yang et al. [68] conducted a simulation in a 76.2mm(length) × 152.4mm(height) × 25.4mm(width) domain. Their results are presented in Fig. 15.
Without considering the natural convection, the melting front is parallel to the heated wall during the whole phase-change process because the liquid PCM is assumed to not flow. In contrast, taking into account the natural convection, a sloped-shape interface is formed. Yang et al. also compared contributions of the two heat transfer mechanisms to PCM melting and the results in Fig. S3 showed that, although the influence of natural convection was unignorable, heat conduction played a more important role. From this point of view, it also explains the reason why the addition of high thermal-conductivity supports enhances the global heat transfer rate although it suppresses the natural convection.
Pore size and porosity of the porous support also exert effects on the phasechange heat-transfer performance of ss-PCM. In Tian et al. [98]'s study, metal foams with smaller porosity and pore size were recommended to improve the heat transfer rate because smaller pore size increased the contact area between PCM and the porous support and smaller porosity indicated more support materials embedded. Zhu et al. [96] investigated the thermal performance of paraffin/aluminum foam composite in a rectangular container and they found that using metal foam with smaller pore size could improve the thermal response of PCM. Zhao et al. [90] studied the solid-liquid phase-change heat transfer of succinonitrile/copper foam in a square container and analyzed the effect of porosity and pore size. The predicted results are shown in Fig.   S4. As porosity decreased, the melting process was accelerated, while the difference in melting time tended to smaller. The phase change process could also be promoted by decreasing pore size because the larger interfacial area can be generated. The similar laws can be found in Sundarram et al. [73]'s study. Moreover, the laws are consistent with the experimental finding by Zhao et al. [30]. Yang et al. [99] numerically analysed the solid-liquid phase change of ss-PCM in a 10cm(width) x 30cm(height) domain. It was found that the heat transfer was enhanced and the total melting time was shorted by increasing porosity linearly from bottom to top, compared with that for a constant porosity, owing to the enhanced natural convection.

Cylindrical container
An example of physical model of ss-PCM in cylindrical container is shown in Fig.   S5. In 2008, Lafdi et al. [100] numerically investigated the phase change process of PCM/graphite foam in a cylindrical domain. They found that the heat transfer performance of the PCM system was enhanced significantly owing to the addition of high thermal-conductivity graphite foam. The average output power of the system enhanced by graphite foam was increased more than 8 times for the space application, while approximately 5 times for the terrestrial application. Liu et al. [83] performed a 3D simulation in a 40mm (inner diameter) x 82mm (inner diameter) x 250mm (length) computational domain. In their study, the solid/liquid interface propagation and temperature distribution etc. were predicted. As Fig. 16 shows, the temperature along the central tube decreased because the heat was absorbed by the surrounding PCM as the heat transfer fluid flowed. A 335K iso-surface moved from the inlet to the outlet during the melting process. The right side of the iso-surface was completely molten PCM while the left represented where the melting was proceeding. Later, Zhang et al.
[23] studied the melting/solidification of ss-PCM in a vertical cylindrical storage unit.
They used NaNO3 and KNO3 as PCM and found that the solid/liquid interface curved gradually during the melting process, indicating increasing natural convection; for the solidification process, the phase interface was almost parallel to the cooling wall, illustrating that the heat conduction was dominant. Fig. 16. Temperature field evolution during melting process in a cylindrical container [83] Some researchers optimized the configuration of enhancers to improve thermal performance. For instance, Nithyanandam et al. [85] conducted a numerical study on the metal foam and heat pipe enhanced TES unit. They simulated the melting/solidification process of four different configurations: (1) no heat pipe; (2) no heat pipe-metal foam; (3) 2 horizontal heat pipes-metal foam; (4) 2 vertical heat pipesmetal foam. Among, the 2 vertical heat pipes with metal foam was recommended to improve melting/solidification efficiency. Later, Xu et al. [101] evaluated and optimized the melting performance of a cylindrical TES unit partially filled with a porous medium.
The predicted temperature and flow field are presented in Fig. 17. Their research showed that, to make the most of natural convection and save cost, the porous medium should be placed in the lower part of the tube and the optimal filling height ratio was 0.7. 3.2 Pore-scale LB simulation The above simulations are continuum-based. Lattice Boltzmann method (LBM) is a relatively new approach which is particle-based and employs a simple kinetic model [102]. The principle of this numerical method is that, fluid is discretized into small particles and macroscopic heat and mass transfer characteristics are predicted by evolving thermal motion of fluid particles. Compared with the above numerical method, LB method has three advantages: (1) simple calculation procedure; (2) parallel computation; (3) robust ability to handle complex geometries [102,103]. (1) PCM and porous supports are homogeneous and isotropic; (2) Liquid PCM is incompressible; (3) Thermo-physical properties of PCM and porous supports are constant; (4) Porous supports are rigid.
Based on the above assumptions, the continuity equation of PCM is simplified as [31]: Due to considering the real pore structure, the semi-empirical models are not needed, and the momentum equation can be expressed as: The energy equation for the PCM is given by: For the porous support: The total enthalpy H in Eq. (17) includes both sensible and latent enthalpy and is given by: where and are the density distribution function and equilibrium density distribution function, respectively. The equilibrium density distribution function is given by: is the forcing term and expressed as in which is the total body force: It should be noted that the total body force includes the Darcy (the first term on the right hand) and Forchheimer drag force (the second term on the right hand). The evolution equation in Guo and Zhao's model can recover the macroscopic generalized non-Darcy equation through the Chapman-Enskog expansion [66].
However, for the pore-scale simulation, the semi-empirical correlations are not needed. Hence, many researchers simplified Guo and Zhao's BKG-LB model [109,110]. Huo et al. [109] only considered the buoyancy force and defined as: This expression is also employed in Li et al. [110]'s model.
In addition, as the porosity is a statistic parameter for the REV-scale simulation, it is not needed in the pore-scale simulation. The equilibrium density distribution function and the forcing term are rewritten as [109] = [1 + Li et al. [69] utilized a more simplified form to define : For the two dimensional (2D) cases, D2Q9 lattice scheme is commonly used for discretizing the velocity space [103] and the nine velocities in the D2Q9 lattice ( Where is the streaming direction, and is the streaming speed defined as = ∆ /∆ , in which ∆ and ∆ are the lattice cell size and the lattice time step, respectively.

MRT-LB model
The MRT-LB model was proposed by d'Humières et al. [113]. It has a different collision term from the BGK-LB model. The standard evolution equation of the MRT-LB model (without a forcing term) is expressed as [114,115]: where and are the transition matrix and relaxation matrix, respectively. and are the velocity moments of and their equilibria: linearly transforms the distribution function into the velocity moment : Although the MRT-LB model is widely applied in the REV-scale simulation [105,115,116], few studies employ the MRT-LB method to model the flow field in the porescale simulation.
where and is the enthalpy distribution function and the equilibrium enthalpy distribution function: in which ̃ is the weight coefficients, is the sound speed of lattice, ∈ (0, 1).
For the porous support, due to the absence of flow, the equilibrium temperature distribution function is simplified as [109,117] =̃ is the relaxation time related to the thermal diffusivity coefficient : It is noteworthy that Jiaung et al. [121]'s and Chatterjee et al. [122]'s models are suitable for the conduction-dominated heat transfer. In 2008, Huber et al. [123] developed an LB model which incorporated the natural convection. They analyzed the transition from the conduction-dominated heat transfer to the fully developed convection. Their predicted results were in good consistency with scaling laws calculated by Jany et al. [124].
Huber et al. [123]'s thermal BGK-LB model was utilized by Song et al. [125] to perform the pore-scale simulation of the freezing process in the soil. Chen et al. [31] employed Eq. (40) to simulate the melting process of paraffin/aluminum foam ss-PCM .
The numerical results were in good agreement with experimental observations. Li et al. [69] used the thermal BGK-LB model to investigate the pore-scale heat-transfer performance in a 3D geometry. An interesting finding is that the secondary convection was captured in the transverse direction and neglecting the secondary convection in the 2D model led to a significant error in predicting heat transfer performance.
The transformation matrix relates the moment space and the velocity space In the pore-scale simulation using thermal MRT-LB model, the D2Q9 lattice scheme is commonly employed [107,128]. The transformation matrix can be found in Ref. [129][130][131] while the relaxation matrix is expressed as The equilibrium distribution function is given by Eq. (36) or Eq. (38). The equilibrium moment is given by = For the enthalpy-based MRT-LB model, the thermal conductivity is expressed as = 2 ( − 0.5)∆ (46) while for the temperature-based MRT-LB model, the thermal diffusivity is given by The application of the thermal MRT-LB model in the pore-scale simulation is relatively less compared with that of the thermal BGK-LB model. Ren et al. [128] employed the thermal MRT-LB model to perform the pore-scale simulation of the melting process in a metal foam. One year later, they utilized this model for another time to numerically investigate the phase-change process in a heat pipe-assisted metal foam.

Reconstruction of porous structure
To perform the pore-structure simulation, the first step is describing the pore structure of the ss-PCM. Currently, there are four methods to reconstruct the pore structure: (1) X-ray computed tomography (CT) [132,133] (2) quartet structure generation set (QSGS) method [134,135] (3) hard-sphere Monte-Carlo method [136] (4) random obstacle location method [137]. The comparison of these methods is presented in Table S1.
Among the four approaches, X-ray CT and QSGS methods are commonly utilized to construct the computational domain in simulating the phase-change heat-transfer process of porous ss-PCMs. For instance, Ren et al. [128] employed the QSGS method to reconstruct the microstructure of 20mm × 20mm metal foam (Fig. 18). The pore size ranged from 0.5mm to 1.0mm. One year later, they utilized this approach for another time to construct the porous structure of a heat pipe-assisted TES unit [107].
In contrast, Liu et al. [111] employed the X-ray CT approach and reconstructed a 5mm × 5mm × 5mm geometry. Li et al. [110] used X-ray CT to obtain the detailed geometric information of the ss-PCM structure. As a great number of computational resources would be consumed if the entire reconstructed domain was considered, they extracted a 20mm × 20mm computational domain to perform the simulation , as shown in Fig.   18(b1) and (b2). Other applications of methods of reconstructing ss-PCM structure are list in Table 10.  [20]'s investigation whose physical model was also heated by the left wall. As Fig. 19 shows, the overall distribution of the temperature field predicted by pore-and REV-scale method is similar: at the initial stage, the temperature contours were approximately parallel to the heating wall. With the elapse of time, more PCM melted and the effect of the natural convection was increasingly significant, leading to slope-shape isotherms.
Compared to the REV-scale simulation, the pore-scale simulation can reflect the influence of the porous structure. As seen from Fig. 19, the isotherms predicted by the pore-scale method were not as smooth as those by the REV-scale method. In the pore-scale simulation, some parts of the isotherms were parallel to the skeleton interface (see the marked area in Fig. 19(a)), indicating that the porous structure exerted an effect on the temperature distribution. However, this phenomenon cannot be simulated by the REV-scale method. In fact, the pore-scale simulation result is more reasonable. The pore distribution is ununiform, and some support protrudes at the interface; due to the high thermal conductivity, the support temperature is higher, leading to the melting of PCMs surrounding the support and forming an irregular interface. Fig. 19. Comparison of temperature fields predicted by different scale methods [20,110] The temperature evolution leads to the occurrence of solid/liquid phase change.
The effect of pore structure is presented more clearly through the comparison of solid/liquid interfaces in Fig. 20. As this figure shows, in the pore-scale simulation, the phase interface is zig-zag-shaped and many parts of the interface are parallel to the skeleton, while in REV-scale simulation, this detailed information cannot be captured. Ren et al. [128] simulated the melting process in a 20mm × 20mm square tank.
They found that PCM close to the skeleton melted faster than that far away. Huo et al. [109] investigated the melting process in a paraffin/aluminum porous ss-PCM with various porosities. They generated the porous structure via QSGS method. As more porous medium was located in the right region, PCM at this area melted more quickly.
The above pore-scale investigations illustrate that the phase change heat transfer in porous ss-PCMs is affected by the porous structure, while the REV-scale simulation neglect the influence of the porous structure.

Flow field
Ren et al. [128] simulated the flow field of the phase change process as shown in Fig. 21(a). It can be clearly seen that the flow of liquid PCM was driven by the buoyancy force and passed through the gap between the support. At = 0.06, the natural convection was further developed and a large vortex was formed in the middle field. The typical flow field simulated by the REV-scale method [68] is also presented in Fig. 21(b). It can be observed that the detailed flow through pores was ignored and the flow field was simplified greatly. Later, Ren et al. [107] studied the melting process of nanoparticle-PCM in a heat pipe-assisted TES unit. They found that the liquid PCM flowed through the gap between the metal foam and formed a small vortex at the early stage while the vortex did not occur in the case of high porosity and nanoparticle volume fraction. In Li et al. [110]'s pore-scale simulation, a large vortex was formed with a clockwise direction while some small vortices were generated due to the shear stress.
At the top area, the circulation drove the high-temperature PCM to the low-temperature region while at the bottom area, it dragged the low-temperature PCM to the hightemperature region. The detailed description of the flow field can also be found in Huo et al. [109]'s investigation. The flow patterns are shown in Fig. 22. The flow was weakened in the nearly closed pores while it was accelerated in the pore throats. Fig. 22. Streamlines in pores simulated by the pore-scale method [109] In summary, the pore-scale simulation can provide the overall phase-change heat-transfer characteristics like the REV-scale simulation. More importantly, it is able to display the flow and heat transfer in pores. This advantage is of great importance especially for mesoporous, microporous and hierarchical porous materials because the pore effect is significant in these supports [11,12]. In other words, when investigating phase change heat transfer in these support materials, the support structure must be considered. Hence, the pore-scale simulation enables researchers to investigate the influence of pores on the phase change heat transfer in porous ss-PCMs and this method should be paid much more attention.
Both experimental and numerical studies on phase change heat transfer in porous ss-PCMs are counted and the results are shown in Fig. 23(a). We also count investigations on support material preparations for LHTES in the recent 20 years and the outcomes are presented in Fig. 23 To make the best of new materials, more efforts should be devoted to investigating phase change heat transfer in mesoporous, microporous and hierarchical porous ss-PCMs. (1) The coverage of PCMs used in current studies is very narrow. As Fig. 2 shows, the phase change temperature of overwhelming PCMs is from 0℃ to 70℃, within the range of low-temperature PCMs. For the middle/high-temperature LHTES applications, such as industrial waste heat recovery, researchers should broaden the selection of PCMs with higher phase change temperature.
(2) There exists a research gap between phase change heat transfer and material preparation. Both experimental and numerical investigations focus on macroporous ss-PCMs. However, the preparation of nanoporous materials, especially mesoporous and hierarchical porous materials, for LHTES has developed rapidly in recent years.
In order to make the best of new materials, it is essential to investigate the phase change heat transfer in mesoporous, microporous and hierarchical porous ss-PCMs in the future.
(3) The pore-scale method exhibits great potential for the simulation of mesoporous, microporous and hierarchical porous materials and LBM is a suitable simulation tool due to its transient inherence. However, according to previous researches, the thermodynamic properties of PCMs in nanoporous supports changes a lot due to the size effect and confinement effect. Thus, to explore the phase change heat transfer in nanoporous ss-PCMs, researchers may adopt LBM with other advanced methods, such as molecular dynamics simulation.

Conflicts of interest
There are no conflicts to declare.