1. Introduction
Phase change material (PCM) can absorb and release a large amount of latent heat during the melting or solidification process. It is therefore widely used in energy storage devices, such as the photovoltaic (PV) devices [
1,
2], domestic buildings [
3], waste heat recovery [
4], and the thermal management of battery [
5], as well as the garment [
6]. In these applications, it is crucial to fulfilling the energy storage and management in a short time because PCM is widely applied in time-dependent environments, such as the harvest of intermittent solar and wind energy, the charging and discharging processes of lithium-ion batteries of electric vehicles, and so on. Therefore, there is a persistent demand on the understanding of heat transfer and melting rate for the PCM to improve its thermal management ability and satisfy its diversity application.
The rectangular enclosure is the most common container of thermal energy storage devices. For instance, the rectangular heatsink was used for PV [
7,
8], the rectangular collector device was used for solar power [
9], and PCM filled in the rectangular containers integrated into the building constructs was used to store or release the thermal energy [
10,
11]. Due to the low heat conductivity of PCMs especially for nonmetallic PCM, many researchers presented useful technologies to improve the melting rate when the PCMs were melting in the rectangular cavities. For instance, an inverted technique to change the direction of gravitational acceleration [
12] and lots of pin fins to expand the heat transfer area of the cavity [
13] are proposed to enhance the melting rate of PCM in a rectangular enclosure. Indeed, among these methods of enhancing the melting rate of PCM, increasing the heat transfer area is an effective method [
14,
15,
16]. Beyond that, adding the nanoparticle into phase change material in the rectangular cavities can also improve the melting rate of PCM effectively [
17,
18].
Most researchers focus on improving the low thermal conductivity of PCM to enhance melting rates. On the other hand, the geometrical characteristic of containers also influences the melting rate of pure PCM [
19]. In Kamkari et al. [
20] and Zhao et al. [
21], a series of experiments and simulations were performed to study the influence of tilt angle of PCM enclosures on melting by heating different sides. It can be found that the inclined angle of containers could affect the melting rate of PCM greatly. Moreover, the aspect ratio of the rectangular enclosure also has a contribution to the melting rate of PCM [
22,
23,
24]. Besides the influence of orientation and aspect ratio of containers, heating the different sides of the rectangular enclosure also plays an important role in the melting process of PCM. For instance, by using a heat flow meter apparatus, Fantucci et al. [
25] have measured a significant change of equivalent thermal conductivity of the PCM between the upper and the below-heated samples. The difference was mainly attributed to the convection occurring in the below-heated sample. Therefore, it is necessary to study the influence of heated side for PCM melting in the rectangular containers. Because of the influence of natural convection, different heated sides could acquire different melting rate of PCM [
26,
27,
28,
29]. Konig-Haagen et al. numerically studied the melting of PCM with the influence of natural convection in a rectangular cavity with five different solvers by heating the lateral side [
30]. They had shown that the approach formulated in terms of enthalpy appear to be more robust than the other using temperature formulations. Hence, the enthalpy-porosity method is adopted in this study.
Although the influence of the aspect ratio of the rectangular enclosure and the heated sides on melting rate are well identified in the previous studies, these influence factors are more emphasized on the qualitative behavior. How to guide the engineer to design a suitable enclosure is still indistinct. Therefore, the first part of this paper is devoted to answering: How to design the rectangular cavity to achieve the highest melting rate of pure PCM with same thermal storage capacity through heating one side of the cavity. Therefore, heating from the top, right, and bottom side and the setting of different aspect ratios are discussed in detail. The second part of this paper is devoted to achieving the dimensionless scaling law or empirical relation for such a melting process in a rectangular enclosure with different aspect ratios. The proposed fitting correlations could predict the liquid fraction, Nusselt number, and the dimensionless specific enthalpy versus dimensionless time with the influence of aspect ratio and heated side. This study provides a fundamental understanding of the melting process of PCM with heating the different sides, and significant references for designing the rectangular container of PCM.
This paper is organized as follows. In
Section 2, the mathematical formulations and numerical methods, as well as numerical validations, are presented. In
Section 3, the numerical results of three different aspect ratios heated by different sides are reported to illustrate the influence of aspect ratio for the different heated side. Then the results with more different aspect ratios are presented in the original unit to exhibit the melting rate with different settings. The analysis and correlation fitting in the non-dimensional form are presented in
Section 4, followed by a concluding remark in
Section 5.
3. Results and Discussion
In this section, a parametric study on the aspect ratio of enclosure heated by different walls for a given cross-sectional area, i.e., the same thermal storage capacity, is presented to compare their melting rates. Considering the dimensionless melting time is related to the height of the enclosure, then for the cases with different aspect ratios, the unit of dimensionless melting time corresponds to different real-time. In order to make the comparison clear enough, the results are presented in a real physical unit in this section. Firstly, the three typical enclosures of different aspect ratios of 0.5, 1.0, and 2.0 are selected to illustrate the influence of aspect ratio and a different heated side.
Figure 7 shows the temperature contours and streamlines of the PCM melting in the three different horizontal enclosures heated from the top, right, and bottom sides to illustrate the influence of aspect ratio and different heated sides. It is noted that the heat transfer for the top heated case is obviously only dependent on conduction. The heat flux is doubtless the minimum among the three different heated walls and could even be calculated without simulation. However, this study still made the simulation in order to extract the critical time on the transition from conductive to the convective regime for the bottom heated case as it includes the both effective after a critical transitional time. It is clear that even though the mass of the PCM in these three different enclosures is the same, the melting front, as well as the melting of PCM at given instants, exhibit non-trivial differences.
At the time of 30 s, the interfaces of PCM in these enclosures heated from the top or bottom side is approximately planar and has equivalent distance to the heated wall, which suggests that the heat transfer is dominative in conductive regime, although there are even number of rolls at the bottom of enclosures when heated from below, as shown in
Figure 7. The rolls are approximately square which suggests that the influence of buoyancy is working but weak so that the width of the roll is as large as the height of rolls. However, when the right side is heated, the buoyancy of PCM plays a more important role since the effective height is much larger than the other two cases. Consequently, the influence of natural convection is dominative, and the fluid with the higher temperature near the top of the enclosure facilitates the melting of PCM and vice versa. Ultimately, an approximately triangle liquid cavity is produced. As the melting goes on, the height of the fluid zone increases for the bottom heated cases. Consequently, the corresponding effective Rayleigh number increases, which is defined based on the real averaged height of liquid instead of the height of the enclosure (
H). Then the bottom heated case switch to convective regimes. Meanwhile, the rolls deform to irregular shapes to facilitate heat transfer, as shown in
Figure 7.
Figure 8 displays the evolutional liquid fraction in the three different enclosures by heating the different sides. In
Figure 8a,b, the melting of PCM in the enclosure heated from below is faster than the other two situations. However, this rule can be not applied to the enclosure with
AR = 2.0, where the melting rate of PCM in the enclosure heated from the right side is the fastest (see
Figure 8c). It indicates that increasing the aspect ratio can be converted from bottom heating dominant to right heating dominant. Because not only the size of the right side becomes larger but also the effective Rayleigh number becomes larger, to enhance the heat transfer of the unit area. Then, it can be known that the right heated case has a higher melting rate as
AR = 2.0 and the higher melting rate of the bottom heated case occurs in the range of 0.5 ≤
AR ≤ 1.
When the melting process is dominated by the conductive regime, the corresponding Nusselt number has the scaling law as follow [
38]:
For the melting process in the quasi-steady convection regime, there is also a corresponding boundary layer convective scaling law [
38]. It well agrees with the 1/4 law of Rayleigh number as below:
Figure 9 shows the evolution of Nusselt numbers on the heated sides as well as a comparison to theoretical prediction [
38]. The best coefficients for the above-mentioned scaling law are fitted and listed in
Table 2 for these three enclosures. It is clear that the invariant coefficient 1.155 is satisfied for the conductive regime and the coefficients for the convective regime are a litter bit affected by the increased Rayleigh number range from 0.397 to 0.411. The reason could be explained by the comparable trivial contribution of conduction in these convection regimes. To simplify the analysis, an invariant coefficient of 0.4 as an acceptable value is proposed for these cases.
In
Figure 9, it can be observed that the Nusselt numbers are in a good agreement with the scaling law as stated in Equations (12) and (13) when the enclosures are heated from the top and right sides. It is also demonstrated that convection is the main factor during the late stage of the melting of PCM in enclosures heated from the right side. The melting of PCM in rectangular enclosures heated from below is affected by conduction first and then by convection through a clear transition.
Since the aspect ratio plays an important role in the melting of PCM in the rectangular enclosures, its influence on melting is checked in detail.
Figure 10 shows the evolutional liquid fraction for the PCM melting process in rectangular enclosures with many different aspect ratios by heating the different sides. The evolutional liquid fraction monotonously varies with the influence of aspect ratio. It is indicated that the higher melting rate is supported by a larger heated surface for all cases. From
Figure 10b, it is observed that there is a transition of the curves at a larger aspect ratio for the right heated case. It could be explained by the fact that the top side of the cavity is melted into the liquid so that the area of the melting front reduces. Furthermore, the transition comes early as the increase of aspect ratio.
Nusselt number of the PCM in the cavity is plotted in
Figure 11. For the right heated side (
Figure 11b), it is observed that the Nusselt number slightly decreases at the end of the melting process, and the drop in the Nusselt number increases as the aspect ratio increases. In fact, the melting front does not advance at the same pace so that there is also a residual solid zone as part of the top is liquid. In
Figure 11c, an interesting phenomenon is that there are two obvious peaks of Nusselt number during the melting process, especially for the small aspect ratio cases. One of the peaks could be explained by the transition from conduction to convection. The other peak of Nusselt number could be the reason due to the influence of the decreasing number of rolls in flow fields.
Figure 12 plots the variation of enthalpy along the melting periods. The curves of enthalpy have similar behavior as that of the liquid fraction of PCM. It is noted that the initial value of the enthalpy is negative because the PCM is undercooled initially. Although the trend of these curves is not affected by different aspect ratios, it is obvious that the melting rate is highly dependent on it. How to evaluate the influence of aspect ratio on the melting rate of PCM is essential to the design of some energy storage devices. The foremost works are to describe the problem in a non-dimensional way and identify the possible critical time from the conductive regime to the existing convective regime under the influence of the aspect ratio.
4. Analysis and Correlation
Considering that the melting of PCM is widely used in many significant industrial fields, its physical mechanism attracts many pioneer researchers. However, it is widely admitted that the melting phenomenon is quite complicated so that Viskanta thought that no unified theoretical treatment was achieved [
33]. The instantaneous shape of the solid-liquid interface of PCM is the results of the strong coupling between the complicated flows of liquid and the melting rate of solid. Therefore, this work does not aim to study the theoretical solution on a present problem but to present the semi-empirical mathematical expression based on the basic fundamental understanding of this problem as well as the correlation of present numerical simulations. The purpose is to guide the design of the rectangular enclosure with the right aspect ratio and heated wall.
One of the fundamental understandings of the melting of PCM is that it contains at least conductive and convective regimes for the bottom or lateral heated cases. In fact, more heat transfer regimes on the melting of PCM are proposed in the literature but their influence is comparatively slight and ignores here to simplify the analysis. Based on the comparison between top and bottom heated cases, one could extract a critical time on the transition of the heat transfer regime and it is supposed to be related to aspect ratio. For the bottom side heated cases, it is essential to achieve the critical time when natural convection takes over conduction in the melting process. Then the influence of aspect ratio on these critical time are plotted, as shown in
Figure 13. It is clear that the dimensionless critical time
θc is inversely proportional to the aspect ratio of rectangular enclosures, approximately.
By introducing the critical dimensionless time, it is reasonable to divide the melting process into two parts for the bottom heated cases. Another issue is to unify the influence of aspect ratio as shown by the scattering curves in
Figure 10,
Figure 11 and
Figure 12 in a dimensionless way. In fact, for the top heated case, the unified non-dimensional relation is immediately available when the physical parameter is non-dimensional. For the lateral and bottom heated cases, the non-dimensional melting time should be implemented to separate the influence between conduction and convection. By fitting the results above, the liquid fraction and Nusselt number are plotted, as well as dimensionless specific enthalpy in a consistent way, as shown in
Figure 14,
Figure 15 and
Figure 16.
Figure 14 shows the fitted correlations of the liquid fraction compared with the numerical results. And these correlations are proposed in Equations (14)–(18). When the top side is heated, the liquid fraction with different
AR collapse to a single curve as expected when these data are plotted versus dimensionless time
θ. The problem is essentially independent on the
AR, as shown in
Figure 14a. The liquid fraction (
LF) in this case can be calculated as follows:
When the right side is heated, the aspect ratio plays an important role. By using nonlinear regression analysis, the correlation of liquid fraction (
LF) could be approximated as follows:
where,
In the above equation, the influence of conduction and convection is distinguished with different exponential relation to non-dimensional melting time. Then the liquid fraction with different
AR is approximately convergent to the approximate correlation, as shown in
Figure 14b. When the bottom side is heated, the fitted correlation should be divided into two parts by critical dimensionless time
θc, which could be calculated as follows:
where
c is a coefficient dependent on the variation of the cross-sectional area ((
A/
Aref)
1.5), in order to consider the geometrical contribution on Rayleigh number, in this study, it is 0.05436. By separating the different regimes, the correlation of liquid fraction (
LF) can be approached by following expressions according to the dimensionless melting time:
It is found that the liquid fraction has an exponential relation to AR in the convective regime. The above mathematical expressions on the melting process provide the ability to compare and quantify the melting rate by heating different sides of the rectangular enclosure with different aspect ratios. Especially, in the range of aspect ratio that could provide a higher melting rate by heating the bottom side than that by heating the right side. The answer is given in the last part of this paper by plotting the neutral curve, which provides the critical aspect ratios in the Rayleigh number versus aspect ratio diagram when the bottom and lateral heated wall need the same melting time.
In the same way, the correlation of Nusselt number for the melting process is also proposed based on different heated side.
Figure 15 shows the fitted correlations of Nusselt number compared with the numerical results. And these correlations are shown in Equations (19)–(21). In
Figure 15a, when the top side is heated, the fitted correlation is as follows:
When the right side is heated, Jany and Bejan [
38] proposed a unified mathematical expression for the Nusselt number based on conduction, convection, and mixing regimes of heat transfer. They also presented the necessary empirical constants based on their numerical simulation. In this work, the influence of the aspect ratio is also considered. When the aspect ratio is implemented in their expression and the empirical constant is slightly adjusted, it can be found that the Nusselt number can be calculated as follows:
For the bottom heated cases, it is noticed that the evolutional Nusselt number is irregular by the influence of irregular solid-liquid interface. It is, therefore, unnecessary to fit this irregular curve because it is sensitive to any random process. However, the whole process corresponds to a useful averaged Nusselt number. Hence, here is the above-mentioned correlation as follows:
Figure 16 shows the proposed correlations of dimensionless specific enthalpy (normal to latent heat) and its comparison with numerical results. And these proposed correlations, are shown in Equations (22)–(24). In
Figure 16a, when the top side is heated, the fitted correlation is as follows:
In
Figure 16b, when the right side is heated, the fitted correlation is as follows:
In
Figure 16(c1,c2), when the bottom side is heated, the fitted correlation is as follows:
It is obvious that the aspect ratio plays an important role in the convective regime. With above-mentioned correlations, these results could collapse to almost one line.
In order to validate the above-mentioned correlations, further evaluation is necessary. Here several independent experimental data are selected by heating a different side [
20,
28,
36,
39,
40,
41] and the available correlation in reference [
38] as a comparison is shown in
Figure 17. It seems that our proposed correlations are in good agreement with several different independent experimental results. Compare to the previous prediction, the present correlation exhibits a fidelity accuracy. Furthermore, the Rayleigh number in these experiments are in a wide range from 10
5 to 10
8 as shown in
Table 3, it suggests that our proposed correlation could extend as high Rayleigh number as 10
8. Hence, these correlations could be used to predict the melting of PCM in different rectangular enclosures by heating the different sides.
At last, it comes to the theme of this work about how to design a rectangular enclosure of PCM with the same heat storage capacity (the same cross-sectional area) by heating one sole side to achieve higher melting rate. It is clear that the top heated scheme is never possible to be the optimum option. The melting rate of PCM could be enhanced by a small aspect ratio for the bottom heated case or a large one for the lateral wall heated situation. However, the melting rate of PCM could be at the same level by choosing a suitable aspect ratio. By setting the same liquid fraction and non-dimensional melting time in Equations (18) and (21), the neutral curves of critical aspect ratio and its relationship to Rayleigh number can be plotted, as shown in
Figure 18.
From the Rayleigh number adopted in the present simulation, the critical aspect ratio is 1.25 when all of the PCM is changed to liquid, i.e., LF = 1.0. In fact, PCMs are not always expected to be melted out in their designing condition of many applications. Therefore, the neutral curves of the critical aspect ratio for LF = 0.5 and LF = 0.75 are also plotted in the same way. For the present numerical work, the corresponding critical aspect ratios are 0.82 and 0.93. It is noted that as the aiming liquid fraction reduces, the critical aspect ratio decreases. Based on the previous numerical results, it is clear that the melting rate is low at the beginning and increases in the evolutional process for the bottom heated cases. On the contrary, the melting rate reduces at the last stage of phase-changing for the lateral heated case. Consequently, the critical aspect ratio reduces as the aiming liquid fraction decreases.
It should be pointed out that the mushy zone constant in the numerical simulations is set to 105, which is the most selected constant, may not be suitable for all kinds of PCM. The reference density is set to the density of liquid PCM. Consequently, the empirical constants of correlation are supposed to be slightly changed based on the variation of mushy zone constant, different PCM, or variable density. The influences of the mushy zone constant and the variable density on the critical aspect ratio are not the scope of this study but would be the future work.
5. Conclusions
Considering the rectangular enclosure is frequently adopted in phase change energy storage, the melting process in the rectangular enclosures with different aspect ratios heated from the different sides are numerically studied in this work. The Rayleigh numbers in these cases are in the order of 105. This study compared the melting rate for the enclosure with the same sectional area but a different aspect ratio with different heated sides. It is clear that the top heated enclosure only affected by the conductive regime, while the other heated sides could also benefit from the convective regime. For the bottom heated side case, the difference to the top heated side is compared to the extract time transition from conductive to the convective regime at a critical time. This transitional time is observed as inversely proportional to the aspect ratio of the rectangular enclosure.
As the aspect ratio increases, the lateral heated side provides a more obvious advantage. On the contrary, the bottom heated side may result in the best melting rate for the small aspect ratio cases. However, the choice of aspect ratio is sometimes limited by many complex factors, such as the cost, the area of the heated surface, available space, and so on. It is obscured to know which one is better between the bottom and laterally heated schemes by their competitive melting rate of PCM as the rectangular enclosure close to be a square. Therefore, it is desirable to know a critical aspect ratio for one scheme better than the other for the design of the phase change energy storage.
In order to recommend the best aspect ratio of rectangular enclosure quantitatively for energy storage, a series of correlations in dimensionless style for the liquid fraction and Nusselt number, as well as the specific enthalpy, are proposed. These correlations are fitted from numerous numerical results and validated for several independent experimental results in a broad range of Rayleigh numbers. It is found that the proposed correlation could also predict the melting process much better than the other available correlation in reference to the rectangular enclosure with different heated sides. The present correlation could also be used for a higher Rayleigh number as large as 108. Base on the comparison of different schemes on the heated side by using the proposed correlation expression, the neutral curves of critical aspect ratio versus Rayleigh number is obtained in a diagram. The diagram shows that when AR < 1.25, the total melting time of PCM in the cavity heated from below is the smallest (LF = 1.0). The critical aspect ratios for LF = 0.5 and LF = 0.75 are respectively, 0.82 and 0.93. This study is helpful for the selection of the appropriate aspect ratio and heating method to achieve excellent comprehensive energy storage performance of PCM.