Abstract

We numerically investigate the melting and solidification behavior of phase-change materials (PCMs) encapsulated in a small-radii cylinder subjected to a cyclic convective boundary condition (square-wave). First, we explore the effects of the Stefan and Biot numbers on the nondimensionalized time required for a PCM initially held at Tcold to melt and reach the crossflow temperature Thot (i.e., reference Fourier number T̃ref). The increase in either Stefan or Biot number decreases T̃ref which can be predicted accurately using the correlation developed in this work. The variations of the PCM melt fraction, surface temperature, and heat transfer rate as a function of Fourier number are reported and analyzed. We further study the effect of the cyclic Fourier number T̃ on the periodic melting and freezing process. The melting or freezing front initiates at the outer periphery of the PCM and propagates toward the center. At higher frequencies, multiple two-phase interfaces are generated (propagating inward), and the surface temperature oscillates in the vicinity of the melting temperature. This increases the effective temperature difference with the crossflow and leads to a higher overall heat transfer.

Introduction

Phase-change materials (PCMs) are widely used in energy storage applications because of their high latent heat [13]. The temperature and density of PCMs remain nearly constant during the phase change [4], which makes them suitable for a variety of thermal management applications including heat exchangers [1,5], electronic cooling [6,7], building applications [8], solar energy storage [9], and spacecraft thermal control [10,11]. In these applications, PCMs are often encapsulated in a variety of enclosures (e.g., cylinders, spheres, annular cavities, etc.) to reduce their reactivity with the environment, control the volume change, and provide a supporting structure [12].

Phase-change processes in the encapsulated geometries are inherently complex and transient, since the propagation of the solid–liquid interfaces depend on the rate and direction of heat transfer. In larger geometries, the flow induced by the buoyancy can be important. Apart from experimental and analytical approaches, different numerical methods based on enthalpy formulation models, effective heat capacity models, and temperature-transforming models have been used in the literature to simulate the complex phase-change process [1319]. Analytical solutions are obtained by tracking the two-phase interface and solving each phase independently while coupling the heat transfer at the interface. However, the enthalpy formulation [20,21] or the effective heat capacity [22] models utilize a single differential equation that governs both phases, and hence there is no need to track the interface. In the enthalpy formulation model, the total enthalpy of the PCM is decomposed into latent and sensible heat to obtain a differential equation that includes the melt fraction. On the other hand, in the effective capacity model, a large heat capacity is used in the phase-change temperature range to account for the latent heat effects [17,22,23].

Many experimental, analytical, and numerical studies on the melting and solidification of the PCM in encapsulated geometries have been reported in the literature by imposing constant temperature or heat flux at the encapsulation boundary and for a variety of encapsulation geometries such as spheres [2426], rectangular containers [2730], and horizontal [3133] and vertical cylinders [3443]. These studies provide a good understanding of the heat transfer and thermal storage characteristics of the encapsulated PCMs (EPCMs) under fixed temperature or heat flux boundary conditions.

In many practical applications such as PCM heat exchanger design [1], buildings exposed to solar radiation with embedded PCMs for energy management, and waste heat thermal storage, the encapsulated PCM is exposed to cyclic thermal conditions. Due to the transient nature of the melting front and the low thermal conductivity associated with the PCMs, the heat transfer rate could vary significantly under various cycling frequencies. Moreover, multiple phase bands inside the encapsulated domain may form [44], which are not studied extensively in the literature. Among the limited number of studies that considered a cyclic boundary condition, Ho and Chu [45] imposed a sinusoidal surface temperature to study a natural convection-dominated melting process. Casano and Piva [46] performed numerical and experimental studies of the melting and freezing process in a plane PCM slab under cyclic temperature boundary conditions and studied the effects of Fourier and Stefan numbers on the energy storage capacity. They observed that the mean value of the energy storage in the system decreases at higher values of Fourier and Stefan numbers. Mazzeo et al. [47] solved a Stefan problem exposed to a periodic boundary condition analytically by considering a single melting interface. More recently, Kant et al. [48] studied the melting behavior of PCMs under sinusoidally varying heat flux conditions and reported that the time of melting reduces under a variable heat flux when compared to an equivalent constant heat flux.

This paper focuses on the thermal analysis of PCM encapsulated in a long cylindrical polymer tube of a small diameter (larger aspect ratio) under a cyclic convective boundary condition. The effects of natural convection are neglected because of the small tube diameters [49]. The motivation for the research stems from the need to accurately capture the transient movement of the phase-change interface and heat transfer rate to and from an encapsulated PCM exposed to cyclic heating and cooling [1] over a broad range of heating and cooling times that feature multiple phase-change interfaces. The paper employs a numerical approach and characterizes the effect of PCM thermophysical properties, tube radius, ambient air temperature, convective heat transfer coefficient, and the duration of heating and cooling on the phase-change behavior and thermal performance of the encapsulated PCM. Results are presented using dimensionless parameters (i.e., Fourier, Stefan, and Biot numbers) to provide a general framework for utilizing PCM in applications where the cyclic heating and cooling is required or beneficial. For example, Karmakar et al. [1] proposed an EPCM heat exchanger designed to cool the water flow exiting a steam condenser in a steam power plant. In the EPCM heat exchanger, warm water from the condenser flows over an array of EPCM tubes and melts the PCM. This is followed by a period where ambient air flows over the tubes to refreeze the PCM and make it available for the next period of melting under the warm water flow. Thus, the cylindrical EPCM tubes were periodically subjected to heating and cooling. In addition to cooling in power plants [1], there is a potential application of similar EPCM heat exchangers in automotive and building heating and cooling applications, data center cooling, and a variety of other industrial applications.

Computational Details

The conservation of total enthalpy of the above-described system with negligible convection can be written as [50]
ρht=·(kT)
(1)

where h is the specific enthalpy, and T is the temperature.

To isolate the nonlinearity involved in the enthalpy–temperature relation, Voller et al. [50,51] proposed a new method which separates the sensible and latent heat components, and introduced the latent heat component as a nonlinear source term [50,52] in the above equation. This methodology has been widely used in the numerical simulations of the phase-change process and, in general, is referred to as the enthalpy method [5358]. The change in the enthalpy H of a PCM element with mass m over a small change in temperature dT can be written as the change in the sum of the enthalpy of liquid and solid portions
dH=d(Hs+Hl)
(2)
d(mh)=d(mshs+mlhl)
(3)
where m=ms+ml, h denotes the specific enthalpy, and subscripts s and l indicate solid and liquid components, respectively. Using the conservation of mass, Eq. (3) can be expanded and simplified as
(ms+ml)dh=msdhs+hsdms+mldhl+hldml
(4)
Further, assuming a constant density of solid and liquid ρ=ρs=ρl, we have
ρ(Vl+Vm)dh=ρ(Vsdhs+hsdVs+Vldhl+hldVl)
(5)
Rearranging the terms, Eq. (5) can be written as
dh=(1γ)dhshsdγ+γdhl+hldγ
(6)
where γ=Vl/(Vm+Vl) is the liquid volume fraction. By definition, the change in the specific enthalpy of each phase in a reversible constant pressure process can be expressed as
dhs=csdT
(7)
dhl=cldT
(8)
Substituting Eqs. (7) and (8) in Eq. (6), we have
dh=[(1γ)cs+γcl]dT+Ldγ
(9)
where L=hlhs is the specific enthalpy of fusion. Using Eq. (9), one can write Eq. (1) as
ρcTt=·(kT)ρLγt
(10)
or
ρcTt=·(kT)ρhlatt
(11)
where
c=(1γ)cs+γcl
(12)
In this paper, Eq. (11) is solved using a fully implicit method. The time derivatives are discretized using a first-order scheme, and the spatial derivatives are discretized using the second-order central difference scheme. At each time-step, Eq. (13) estimates the heat transfer to each discretized volume as
Q=ρV·(kTnt+Δt)Δt
(13)
where V is the cell volume, Δt is the time-step, and the subscript n represents the subiteration number at each time-step. Using the total enthalpy–temperature curve, we estimate the rate of change in the latent portion of the enthalpy (last term in the right-hand side of Eq. (11)). We then solve Eq. (11) implicitly to calculate Tn+1t+Δt and repeat the process until we achieve convergence such that
max|Tn+1t+ΔtTnt+Δt|<ϵ
(14)

with the convergence threshold of ϵ=106. Accuracy of the numerical simulations depends on the convergence threshold, especially for steep enthalpy curves in the melting temperature range as the error in temperature propagates into the error in enthalpy with a factor equal to the rate of change of enthalpy with temperature (i.e., eh=eT×dh/dT). Hence, it is imperative to use lower threshold levels for the steep melting curves. To start the subiteration process, the temperature at the previous time-step estimates the initial temperature of T0t+Δt, which results in the explicit estimation of heat transfer to each cell. It is important to note that the above procedure requires the knowledge of a functional relationship between the enthalpy and the temperature of the PCM material.

The above-described numerical method is validated with the analytical solutions as well as with the experimental measurements of the phase-change process in a vertical cylinder. Paterson [59] derived an analytical solution of the isothermal melting process inside a cylinder with a line heat source at r = 0. For the isothermal melting process, there is no mushy region, and both the solid Rs(t) and liquid Rl(t) interfaces coincide. According to the analytical solution [59], the location of the liquid–solid interface is given by
Rs(t)=Rl(t)=λt
(15)
where λ is a coefficient that needs to be determined using a transcendental equation (Eq. (29) in Ref. [59]) and is a function of the thermophysical properties of the substance, strength of the line heat source, and the initial temperature. For an ice–water system with the thermophysical properties reported in Table 1, a uniform temperature of 2°C (i.e., ice) at t = 0 and a line heat source of 995.8 J m−1 at r = 0, we obtain λ=6.637×104 m s1/2 [59]. For these conditions, the following latent enthalpy curve is considered for the numerical simulations:
hlat(T)=Lγ(T)={0T<TsLTTsTlTsTsTTlLT>Tl
(16)

where Ts=Tm0.01 and Tl=Tm+0.01. Numerical simulations involving Eqs. (11) and (16) are conducted for a one-dimensional axisymmetric cylindrical domain of radius 1×103 m, where a uniformly distributed grid of 500 cells and a time-step of Δt=1×105 s are employed. Figure 1 shows the analytical solution and numerical predictions of the location of the solid–liquid interface as a function of time. On average, the predictions deviate by 0.7% from the analytical solution.

Fig. 1
Analytical solution [59] and numerical predictions of the melt front as a function of time for an ice–water system with a line heat source
Fig. 1
Analytical solution [59] and numerical predictions of the melt front as a function of time for an ice–water system with a line heat source
Close modal
Table 1

Thermophysical properties of water and ice [59]

ParameterValueParameterValue
ρ (kg m−3)1000cs (J kg−1 K−1)1430.97
ks (W m−1 K−1)2.218cl (J kg−1 K−1)4180.56
kl (W m−1 K−1)0.602Tm (°C)0
L (kJ kg−1)307.94
ParameterValueParameterValue
ρ (kg m−3)1000cs (J kg−1 K−1)1430.97
ks (W m−1 K−1)2.218cl (J kg−1 K−1)4180.56
kl (W m−1 K−1)0.602Tm (°C)0
L (kJ kg−1)307.94
The exact solution of a nonisothermal phase change in cylindrical coordinates with a line heat sink at r = 0 is obtained by Özişik and Uzzell [60]. Due to the presence of a large phase-change temperature range, a mushy region exists between the solid and liquid interfaces which are determined analytically by
Rs(t)=λst
(17)
Rl(t)=λlt
(18)
The current validation involves the freezing process of aluminum–copper alloy containing 5% copper in a one-dimensional axisymmetric cylindrical domain of radius 5 m with a line heat sink at the center [60]. The thermophysical properties of both phases are shown in Table 2 and are an eutectic with the solid fraction at the solidus front of γsu=0.8952. Here, the following enthalpy function is considered for the numerical simulations [50]:
hlat(T)=Lγ(T)={0T<TsLγsuT+(1γsu)TlTsTlTsTsTTlLT>Tl
(19)
Table 2

Thermophysical properties of aluminum–copper alloy (5% copper) [60]

ParameterValueParameterValue
ρ (kg m−3)2723.2Ts (°C)547.8
ks (W m−1 K−1)197.3Tl (°C)642.2
kl (W m−1 K−1)181.7γsu (°C)0.8952
cs (J kg−1 K−1)1046.7L (J/kg)395403
cl (J kg−1 K−1)1256
ParameterValueParameterValue
ρ (kg m−3)2723.2Ts (°C)547.8
ks (W m−1 K−1)197.3Tl (°C)642.2
kl (W m−1 K−1)181.7γsu (°C)0.8952
cs (J kg−1 K−1)1046.7L (J/kg)395403
cl (J kg−1 K−1)1256

Numerical simulations are performed for a cylindrical domain of radius of 5 m for which Eqs. (11) and (19) are solved. The numerical domain is discretized uniformly in the radial direction with 4000 cells, and a time-step of Δt=1s is utilized. Figure 2 shows the analytical and predicted liquid and solid fronts as a function of time. Current numerical predictions agree well with the analytical solution with an average deviation of 2.5% which validates the accuracy of the numerical approach.

Fig. 2
Analytical [60] and predicted liquid and solid fronts as a function of time for a nonisothermal phase change with a line heat source
Fig. 2
Analytical [60] and predicted liquid and solid fronts as a function of time for a nonisothermal phase change with a line heat source
Close modal

Encapsulated Tube With Convective Boundary Condition.

This section presents the ability of our solver to predict the cyclic melting and freezing process of a phase-change material encapsulated in a small diameter vertical cylinder. The experimental test setup by Kapilow et al. [49] was designed to evaluate the melting and freezing characteristics of a commercial PCM (RT35HC) [61] encapsulated inside a cylindrical high-density polyethylene (HDPE) tube subjected to convective boundary condition. They reported the measured PCM volume phase fraction and the overall heat transfer coefficient of the entire system over time. Figure 3 shows the schematic of the encapsulated PCM cross section where the air flows from left to right. The inner and outer tube diameters are 2.78×103 m and 3.18×103 m, respectively. Based on the outer tube diameter as the length scale, the Reynolds number of the air flowing over the tube is 1797. The reported Rayleigh number of the experimental tests is very low (50–1200) due to the small tube diameter and thus suggests insignificant effects of natural convection [49]. Hence, the natural convection effects are also neglected in the current study.

Fig. 3
Encapsulated PCM inside a polymer cylinder
Fig. 3
Encapsulated PCM inside a polymer cylinder
Close modal
The thermophysical properties of the air, PCM, and HDPE tube are reported in Table 3. The reported total specific enthalpy curve of RT35HC by the manufacturer is shown in Fig. 4 (open symbols) [61] assuming Tref=0. The cs and cl are calculated as the slope of htot in the solid and liquid states, respectively (Table 3). The latent heat hlat is then calculated by subtracting the sensible heat from the total specific enthalpy curve. Finally, the latent heat distribution with temperature is fitted using the following equation:
hlat=aerf(Tbc)+d
(20)

where a = 0.45775, b = 34.86900, c = 1.213776, and d = 0.454482. The fitted hlat,hsen, and htot are shown with lines in Fig. 4. It should be noted that the reported cooling and heating enthalpy curve in the manufacturer's data sheet slightly differs [61], and hence the averaged enthalpy curve is considered here. Utilizing the reported enthalpy curve of the RT35HC eliminates any assumptions regarding the enthalpy function, melting temperature range, and potentially leads to more accurate numerical simulations.

Fig. 4
Specific enthalpy of RT35HC as a function of temperature: Symbols [61], lines Eq. (20) and Table 3.
Fig. 4
Specific enthalpy of RT35HC as a function of temperature: Symbols [61], lines Eq. (20) and Table 3.
Close modal
Table 3

Thermophysical properties of air, PCM (RT35HC) [61], and HDPE

ParameterAirRT35HCHDPE
ρ (kg m−3)1.1845825970
ks (W m−1 K−1)0.360.49
kl (W m−1 K−1)0.14
kg (W m−1 K−1)2.597 ×102
cs (J kg−1 K−1)3357.12250
cl (J kg−1 K−1)3228.6
cg (J kg−1 K−1)1006
ν (m2 s−1)1.5571 ×105
ParameterAirRT35HCHDPE
ρ (kg m−3)1.1845825970
ks (W m−1 K−1)0.360.49
kl (W m−1 K−1)0.14
kg (W m−1 K−1)2.597 ×102
cs (J kg−1 K−1)3357.12250
cl (J kg−1 K−1)3228.6
cg (J kg−1 K−1)1006
ν (m2 s−1)1.5571 ×105
In the experiments [49], an encapsulated PCM tube, installed vertically inside a test section, was exposed to an air flow with different temperatures to alternatively melt and freeze the PCM. Initially, cold ambient air of Tcold=25.3°C which is lower than the melting temperature range of the RT35HC flows over a fully solid encapsulated PCM. Air temperature is then sharply increased to Thot=42.9°C which initiates the melting process. The air temperature T(t) can be represented as
T(t)=12(TcoldThot)[sign(sinπtT)+1]+Thot
(21)

where T is the half cycle. In Kapilow's experiment [49], T is large enough that allows PCM to fully melt and freeze, and the experiment is only conducted over one cycle of melting and freezing. Equation (21), however, provides a general cyclic square-wave function which alternates the air temperature between Tcold and Thot.

To model the effects of the external convective heat transfer and the thermal resistance of HDPE tube, a convective boundary condition with an overall heat transfer coefficient for a unit length of the tube is applied to the cylindrical surface of the PCM as
Uoverall(TiT(t))=kPCMriTir
(22)
Uoverall=1Rair+Re
(23)
Rair=1roUair
(24)
Re=ln(ro/ri)ke
(25)

where Rair is the thermal resistance of the air, Re is the thermal resistance of the HDPE tube, ri and ro are the inner and outer radius of the HDPE tube, respectively, and hair is the heat transfer coefficient around the outer periphery of the tube.

To validate the numerical solver with the experimental data, we perform simulations at a low frequency to allow full melting and solidification. The average liquid phase fraction and the predicted overall heat transfer are compared with the measurements of Kapilow et al. [49]. Figure 5 shows the average liquid melt fraction γ of the PCM material with a cyclic convective boundary condition with T=400s. To ensure a grid-independent solution, we perform numerical simulations with three different grids of 10 × 15, 20 × 30, and 40 × 60 grids in radial and azimuthal directions. The dependency of the solution on the time-step is also evaluated by conducting the simulations with the highest grid resolution with Δt= 5, 1, and 0.1 s, which corresponds to Δt/tmelting of 3.37×102,6.76×103, and 6.76×104, respectively (Fig. 6). The computed average liquid fraction over time is grid-independent for all evaluated grids (Fig. 5) and does not change for Δt less than 1 s. Numerical predictions agree well with the measurements during the melting period (200s<t<400s), but the rate of freezing (400s<t<600s) is slightly overpredicted. As the numerical procedures are validated with the analytical solutions, this discrepancy may be associated with:

Fig. 5
Grid independency results for the validation case [49] in terms of average liquid melt fraction over time
Fig. 5
Grid independency results for the validation case [49] in terms of average liquid melt fraction over time
Close modal
Fig. 6
Time-step independency results for the validation case [49] in terms of average liquid melt fraction over time
Fig. 6
Time-step independency results for the validation case [49] in terms of average liquid melt fraction over time
Close modal
  1. The difference between the RT35HC enthalpy curve during melting and freezing [61].

  2. Uncertainty in the thermophysical characteristics of RT35HC (see different values reported in Refs. [6163]).

  3. Variation in the air velocity along the axis of the vertical cylinder in the experiments leads to a lower overall heat transfer coefficient.

This emphasizes that attention needs to be made in characterizing the thermophysical properties of the PCM materials.

Figure 7 shows the heat transfer rate predictions as a function of time during the melting and freezing process. The integration of the rate of change of the latent enthalpy over the numerical domain determines the latent portion of the total heat transfer rate. Predictions of the latent portion of the heat transfer rate agree well with the measurements [49]. The PCM temperature rises from the initial temperature after the air temperature increases at t=200s. The phase change only begins after the PCM temperature in the outer portion of the PCM reaches the melting temperature range, and hence the qlat increases sharply for t>200s and reaches a maximum value which corresponds to the lowest possible resistance of the PCM to the heat transfer. During the melting process, the melt front progresses inward, and hence PCM exerts an extra resistance to the heat transfer, which is associated with the decrease in the latent heat transfer rate as observed in Ref. [49]. The freezing process behaves similarly with the only difference that the freezing rate is faster than the melting due to the higher thermal conductivity of the PCM in its solid phase and more substantial temperature difference with the ambient (see Table 3).

Fig. 7
Grid independency results for the validation case [49] in terms of heat transfer rate as a function of time
Fig. 7
Grid independency results for the validation case [49] in terms of heat transfer rate as a function of time
Close modal

Nondimensionalized Governing Equations.

To broaden the applicability of the current study, Eq. (11) can be nondimensionalized for encapsulated PCMs exposed to convective boundary conditions. Without loss of generality, we consider a PCM exposed to an external convective boundary condition with ambient temperature varying according to Eq. (21) and an overall heat transfer coefficient U that incorporates any tube-wall conduction resistance (see Eq. (23)). A nondimensionalized temperature can be defined as
θ=TTmThotTcold
(26)
where Thot and Tcold are the upper and lower limits of the ambient temperature. Assuming constant thermophysical properties, Eq. (11) can be rewritten as
θt~=~2θ1Steγt~
(27)
where thermal diffusivity α=k/ρc, Fourier number t̃=tα/l2 with l being an appropriate length scale, ~=l, and Ste is the Stefan number, defined as
Ste=c(ThotTcold)L
(28)
and represents the ratio of the sensible to the latent heat transferred to a PCM initially held at Tcold<Tm and heated to reach Thot>Tm (i.e., complete heating), or vice versa. The total heat transferred per unit mass during the complete heating or cooling process is Qmax=L+c(ThotTcold) and can be written in a nondimensional format as
Q̃max=QmaxL=1+Ste
(29)
The specific enthalpy of the PCM, assuming c=cs=cl, reads
h̃=hL=θrefθStedθ+γ
(30)
The convective boundary condition can be written as
Bi(θθ)=n~θ
(31)
where n is the unit vector normal to the surface, and Bi=Ul/k is Biot number. Similar to the numerical approach described earlier, Eq. (27) can be solved using Eqs. (31) and (32) and a known γ(θ) function. We limit the current study to an isothermal melting and freezing process encountered for pure materials (Eq. (16)). A convective boundary condition with a uniform heat transfer coefficient (Bi) and an ambient temperature governed by Eq. (21) is imposed on an encapsulated PCM cylinder of radius R. Hence, θ is only a function of radius and Fourier number θ=θ(r̃,t̃) for the specified Bi and Ste. Assuming an equal ambient temperature difference with respect to the melting temperature (i.e., ThotTm=TmTcold), Eq. (21) can be simplified to
θ(t̃)=0.5[sign(sinπt̃T̃)+1]+0.5
(32)
The goal is to evaluate the time-averaged heat transfer rate of the described system as a function Ste, Bi, and T̃
q̃=f(Bi,Ste,T̃)
(33)
This is not readily available in the literature due to the complex behavior of the PCM under cyclic melting and freezing, especially if the T̃ is not large enough to allow complete heating or cooling. Under this condition, one can postulate that the PCM is expected to reach a cyclic equilibrium behavior with the same cyclic characteristics of the ambient temperature variations after a long time, and this equilibrium state is not known a priori. Furthermore, as the boundary condition varies at the outer periphery of the PCM tube cyclically, it is expected that the melting and freezing interfaces initiate at the outer periphery location and progress toward the core of the tube. However, if the melting or freezing of the entire PCM is not completed before the boundary condition changes from heating to cooling or vice versa, then a single unique interface between the solid and liquid phases may not exist inside the tube. To better understand the above scenario, it is important to first identify the Fourier number in which the heating or cooling process of the entire PCM is completed (i.e., the reference Fourier number) as a function of Stefan and Biot numbers, i.e., T̃ref=g(Ste,Bi). Therefore, we first perform numerical simulations for different values of Ste and Bi by initializing the PCM at θ(r̃,0)=θcold=0.5 with θ=θhot=0.5 and predict the Fourier number for which the temperature of the PCM becomes uniform and reaches θhot. In this process, PCM undergoes both sensible and latent heat transfer. Initially, the PCM is heated to reach the melting temperature (i.e., θ = 0), followed by the melting process and the second phase of sensible heat transfer until the PCM temperature approaches the ambient temperature. The total heat transfer during this process is Qmax=1+Ste, which is the sum of the total sensible and latent heat transfer. The simulations are conducted for 0.01Ste1 and 0.1Bi10 to characterize the effect of Biot and Stefan numbers on T̃ref. For the second part of the study, we define a parameter ζ=T̃/T̃ref and perform simulations for a cyclic ambient temperature governed by Eq. (32) and different values of Ste, Bi, and ζ to evaluate the heat transfer rate as
q̃=f(Bi,Ste,ζ)
(34)

Simulations are performed for the same range of Ste and Bi number and for 0.1ζ1 and continue until the PCM shows a cyclic behavior. Due to the uniform heat transfer coefficient assumption around the circular tube, all simulations are conducted using a one-dimensional grid with 40 grid points uniformly distributed from the center to the outer radius of the PCM using Δt̃=0.01. The selected Δt̃ provides Δt̃/T̃ref<5.1×103 for all simulations.

Results and Discussion

Figure 8 shows the variations of the predicted T̃ref as a function of Ste and Bi. The nature of the plot shows that T̃ref rapidly decreases with increasing values of Ste and Bi numbers. The predicted T̃ref values from the simulations are used to obtain a correlation as a function of Bi and Ste values as follows:
T̃ref(Ste,Bi)=1+1+3.5Ste+0.5BiSteBi
(35)

with a mean and maximum deviation of 2.3% and 3.9% from the predictions, respectively. Figure 9 compares the predicted T̃ref using numerical simulations with the ones calculated using the proposed correlation Eq. (35) at selected Biot numbers. The applicable range of the above correlation is 0.01Ste1 and 0.1Bi10. Although the effect of the Ste and Bi on the T̃ref is clear, it might be also necessary to understand the effect of each individual physical quantity (including the thermophysical properties of the PCM, and tube radius) on the actual time required to complete the heating or cooling process Tref. The Appendix provides more insights on this matter. However, to simplify the interpretation of the Ste and Bi on the PCM behavior, in the context of the following analysis and without loss of generality, we consider constant thermophysical properties and a fixed tube radius. Hence, TrefT̃ref and the variation of the Ste and Bi numbers are limited to the changes in the values of temperature difference (Thot − Tcold) and the overall heat transfer coefficient U, respectively.

Fig. 8
T̃ref as a function of Ste and Bi
Fig. 8
T̃ref as a function of Ste and Bi
Close modal
Fig. 9
T̃ref as a function of Ste at selected Bi values. Symbols represent numerical simulations, and solid lines depict the proposed correlation Eq. (35).
Fig. 9
T̃ref as a function of Ste at selected Bi values. Symbols represent numerical simulations, and solid lines depict the proposed correlation Eq. (35).
Close modal

To understand the effect of Ste, Fig. 10 shows the temporal variation of the melt fraction and the dimensionless heat transfer rate of the system as a function of Fourier number t̃ for Ste values of 0.01, 0.1, and 1 with a fixed Biot number of one. Under a fixed Bi number, the corresponding reference Fourier number T̃ref for these three cases are 155, 20, and 6, respectively. This process consists of three phases: an initial sensible heat transfer where the value of γ remains unchanged, followed by a combined latent and sensible heat transfer where the value of γ changes from 0 to 1, and then again a pure sensible heat transfer toward the end with a constant value of γ. Varying Ste only affects the second phase of heat transfer as when the value of γ changes over time, and thus the initial and final sensible heat transfer phases (i.e., γ/t̃=0) are unaffected by Ste (see Eq. (27)). By definition, the Ste number represents the ratio of sensible heat to the latent heat. Consequently, at very low values of Ste that corresponds to a small temperature difference, the PCM melts at a very low heat transfer rate (see q̃ for Ste = 0.01 in Fig. 10), and thus the duration of PCM melting is the highest. Figure 11 shows this behavior where the surface temperature takes the longest (i.e., larger Fourier number) to reach the ambient condition. On the other hand, when the Ste number increases to 1, the rate of heat transfer in the PCM system is the highest (see q̃ for Ste = 1 in Fig. 10) due to the highest temperature difference, and thus the duration of melting turns out to be the least. The plot of the surface temperature at Ste = 1 in Fig. 11 also confirms this behavior as it reaches the ambient temperature at the fastest rate.

Fig. 10
Effect of Stefan number on melt fraction and heat transfer rate as a function of Fourier number
Fig. 10
Effect of Stefan number on melt fraction and heat transfer rate as a function of Fourier number
Close modal
Fig. 11
Effect of Stefan number on temperature at r̃ = 1 as a function of Fourier number
Fig. 11
Effect of Stefan number on temperature at r̃ = 1 as a function of Fourier number
Close modal

Figure 12 shows the effect of the Bi number on the temporal variation of melt fraction and the heat transfer rate of the system. Under a fixed Ste number, the corresponding reference Fourier number for these three cases are T̃ref = 141, 20, and 7, respectively. Similar to the previous scenario, the entire process of heat transfer in the system consists of two purely sensible heat transfer portions and a combined latent and sensible heat transfer. The rise in the Bi number corresponds to an increase in the values of the convective heat transfer coefficient, and as such, the reference Fourier number decreases due to the increased heat transfer from the ambient to the PCM system. At a small value of the Bi number, the convective resistance dominates the resistance due to conduction, and hence the heat transfer rate is not significantly affected by the variation of the conductive resistance inside the tube during the phase-change process and remains almost constant (see Fig. 12)). The heat transfer rate is also the lowest at this Bi number due to the small heat transfer coefficient. The surface temperature for Bi = 0.1 (shown in Fig. 13) remains close to the PCM melt temperature over most of the period. When the Bi number rises to 1, the heat transfer rate decreases almost linearly with time as the convective resistance is no longer dominant, and hence the heat transfer rate is affected by the variations in the conductive thermal resistance inside the tube. Due to the higher Bi number, the melting rate also increases. Finally, at Bi = 10, the conductive resistance becomes more dominant such that the heat transfer rate in the system shows the strongest dependence on the conductive resistance of the PCM during the phase-change process, and hence, decreases very rapidly with the progression of the phase-change interface toward the center of the tube. The rate of melting is the highest under this condition as the value of the external heat transfer coefficient is the highest. The surface temperature at Bi = 10 approaches the ambient temperature very quickly compared to the previous conditions.

Fig. 12
Effect of Biot number on melt fraction and heat transfer rate as a function of Fourier number
Fig. 12
Effect of Biot number on melt fraction and heat transfer rate as a function of Fourier number
Close modal
Fig. 13
Effect of Biot number on temperature at r̃ = 1 as a function of Fourier number
Fig. 13
Effect of Biot number on temperature at r̃ = 1 as a function of Fourier number
Close modal

The discussion below presents the effect of the cyclic exposure of the PCM cylinder to hot and cold ambient temperature under various fractions of T̃ref. For a fixed Ste = 0.1 and Bi = 1, Fig. 14 shows the variation of melt fraction as a function of Fourier number for ζ = 0.1, 0.5, and 1. This plot corresponds to the Fourier number interval of 2T̃ref after each of the cases has reached a cyclic state; ζ = 0.1 has 10 cycles, ζ = 0.5 has 2 cycles, and ζ = 1 has 1 cycle of alternating hot and cold ambient conditions. At ζ = 1, the heating and cooling process completes in each half period. However, reducing the value of ζ to 0.5 prohibits the completion of the heating and cooling process, and the melt fraction approximately varies between 0.1 and 0.9. For the lowest value of ζ, the melt fraction varies between 0.4 and 0.6 as the heating and cooling periods shorten. It is important to note that at the cyclic equilibrium condition, the melt fraction swings around 0.5 as both the heat transfer coefficient, and the difference between the melting temperature and the ambient temperature are identical for these heating and cooling processes. It is clear that having different values for these parameters during heating and melting could alter the cyclic equilibrium behavior of the PCM.

Fig. 14
Effect of ζ on melt fraction as a function of Fourier number at Ste = 0.1 and Bi = 1
Fig. 14
Effect of ζ on melt fraction as a function of Fourier number at Ste = 0.1 and Bi = 1
Close modal

Figure 15 shows the normalized heat transfer rate distribution for the same Fourier number intervals as in Fig. 14. The time-averaged heat transfer rate for ζ = 0.1 is around 0.95, whereas it decreases to values of 0.81 and 0.57 when ζ is increased to 0.5 and 1, respectively. Thus, under the current conditions, exposure of the PCM to alternate heating or cooling for a shorter duration than T̃ref is more effective in enhancing heat transfer. This interesting behavior can be explained in terms of the surface temperature at r̃ = 1 for various ζ values as shown in Fig. 16. At ζ = 0.1, the surface temperature oscillates closely around the PCM phase-change temperature (θ = 0), and therefore the rate of heat transfer is greater due to the higher temperature difference between the PCM outer surface and the ambient conditions. When the value of ζ increases to 0.5, the temporal variation of the PCM surface temperature moves away from the PCM melt temperature toward the ambient temperature, which decreases the effective temperature difference and the time-averaged heat transfer rate. At the highest value of ζ = 1, the PCM surface temperature reaches the ambient temperature, which further reduces the heat transfer rate, and thus the average heat transfer rate is lower.

Fig. 15
Effect of ζ on heat transfer rate as a function of Fourier number at Ste = 0.1 and Bi = 1
Fig. 15
Effect of ζ on heat transfer rate as a function of Fourier number at Ste = 0.1 and Bi = 1
Close modal
Fig. 16
Effect of ζ on temperature at r̃ = 1 as a function of Fourier number at Ste = 0.1 and Bi = 1
Fig. 16
Effect of ζ on temperature at r̃ = 1 as a function of Fourier number at Ste = 0.1 and Bi = 1
Close modal

As discussed earlier, the cyclic variation of the ambient temperature with a period lower than T̃ref could potentially generate multiple phase-change interfaces. The outer portion of PCM is always exposed to this variation in the boundary condition, and hence the melting or freezing fronts are generated at this outer location. Here, we define the phase-change interface as any location in which the sign of θ changes (since θ = 0 represents the melting temperature). For significantly longer periods, a phase-change front initiates at r̃=1 and progresses toward the core of the tube (r̃=0) such that eventually, the whole PCM material is either liquid or solid. Figure 17 shows this behavior. When ζ is reduced to 0.5, a new phase-change front will be initiated at r̃=1 before the former front reaches the core (see Fig. 18). Hence, the case of ζ=0.5 shows a maximum of two phase-change fronts. The predictions at the ζ value of 0.1 (Fig. 19) shows a maximum of three simultaneous phase-change fronts. It is interesting to note that the speed at which these phase fronts move toward the core is also affected by the value of ζ and increases for a lower values of ζ.

Fig. 17
Interface location as a function of Fourier number for ζ = 1, Ste = 0.1, and Bi = 1
Fig. 17
Interface location as a function of Fourier number for ζ = 1, Ste = 0.1, and Bi = 1
Close modal
Fig. 18
Interface locations as a function of Fourier number for ζ = 0.5, Ste = 0.1, and Bi = 1
Fig. 18
Interface locations as a function of Fourier number for ζ = 0.5, Ste = 0.1, and Bi = 1
Close modal
Fig. 19
Interface locations as a function of Fourier number for ζ = 0.1, Ste = 0.1, and Bi = 1
Fig. 19
Interface locations as a function of Fourier number for ζ = 0.1, Ste = 0.1, and Bi = 1
Close modal

Conclusions

In this study, we numerically investigate the melting or solidification behavior of phase-change materials encapsulated in a small-radii cylinder. We further study the case of the encapsulated PCM undergoing a cyclic heating and cooling process through a convective boundary condition with variable ambient temperature (square-wave) at the cylinder surface. A new correlation is proposed to estimate the reference Fourier number that corresponds to complete melting or freezing as a function of Ste and Bi. Higher Bi and Ste accelerates the freezing/melting process which is completed at smaller Fourier numbers. Studying the effects of Bi and Ste numbers on the rate of heat transfer and surface temperature helps us understand the impact of these parameters on the reference Fourier number.

The effects of the time period of cyclic heating and cooling on the time-averaged heat transfer rate are studied. Shorter periods of cyclic heating enhance the time-averaged heat transfer rate due to lower surface temperature variations. The tube consists of multiple phase-change fronts simultaneously; for lower ζ values (ζ=T̃/T̃ref), the corresponding phase-change fronts move faster toward the center of the tube.

The proposed correlation and the effect of Bi and Ste numbers on overall heat transfer and melting time can be used for applications that utilize small-radii cylindrical encapsulated PCM exposed to cyclic convective boundary conditions [1]. The significant changes observed in the thermal behavior of the encapsulated PCM under different cyclic conditions indicate the importance of conducting such studies for PCMs encapsulated in other geometries such as spheres or rectangular enclosures.

Acknowledgment

This work was supported by a grant from the Advanced Research Project Agency-Energy (DE-AR0000572) under the Advanced Research in Dry cooling (ARID) program with the Electric Power Research Institute (EPRI) as the prime contractor. This financial support is gratefully acknowledged.

Funding Data

  • Advanced Research Project Agency-Energy (DE-AR0000572; Funder ID: 10.13039/100000015).

Nomenclature

     
  • c =

    specific heat (J kg−1 K−1)

  •  
  • D =

    cylinder diameter (m)

  •  
  • e =

    error

  •  
  • h =

    specific enthalpy (J kg−1)

  •  
  • H =

    enthalpy (J)

  •  
  • k =

    thermal conductivity (W m−1 K−1)

  •  
  • L =

    specific enthalpy of fusion (J kg−1)

  •  
  • m =

    mass (kg)

  •  
  • q =

    heat transfer rate (W)

  •  
  • Q =

    heat transfer per unit mass (J kg−1)

  •  
  • Q =

    heat (J)

  •  
  • r =

    radial distance from cylinder center (m)

  •  
  • R =

    thermal resistance (W−1 m2 K)

  •  
  • R =

    melt front location (m)

  •  
  • t =

    time (s)

  •  
  • T =

    temperature (K)

  •  
  • T =

    time period of the half cycle of melting/freezing (s)

  •  
  • U =

    heat transfer coefficient (W−1 m2 K)

  •  
  • V =

    volume (m3)

  •  
  • V =

    cell volume (m3)

  •  
  • Δt =

    time-step (s)

Greek Symbols

    Greek Symbols
     
  • α =

    thermal diffusivity (m2 s−1)

  •  
  • γ =

    liquid volume fraction

  •  
  • ϵ =

    convergence criteria

  •  
  • ζ =

    fraction of the dimensionless half cycle time period of melting/freezing

  •  
  • θ =

    dimensionless temperature

  •  
  • λ =

    empirical constant

  •  
  • ν =

    kinematic viscosity (m2 s−1)

  •  
  • ρ =

    density (kg m−3)

Superscripts and Subscripts

    Superscripts and Subscripts
     
  • cold =

    lower limit of the ambient temperature

  •  
  • e =

    HDPE tube

  •  
  • hot =

    upper limit of the ambient temperature

  •  
  • i =

    inner

  •  
  • l =

    liquid

  •  
  • lat =

    latent heat

  •  
  • m =

    melting

  •  
  • n =

    subiteration number

  •  
  • o =

    outer

  •  
  • ref =

    reference

  •  
  • s =

    solid

  •  
  • su =

    eutectic

  •  
  • =

    ambient

Dimensionless Groups

    Dimensionless Groups
     
  • Bi =

    Biot number

  •  
  • q̃ =

    dimensionless heat transfer rate

  •  
  • r̃ =

    dimensionless radial distance from cylinder center

  •  
  • Ste =

    Stefan number

  •  
  • t̃ =

    Fourier number, /L2

  •  
  • T̃ =

    dimensionless time period of the half cycle of melting/freezing (s)

Appendix: Effect of Physical Quantities on the Complete Heating and Cooling Time

To study the effect of PCM thermal conductivity k, convective heat transfer coefficient h, and the tube radius R on the reference time Tref, we first rewrite Eq. (35) as
T̃ref(Ste,Bi)=ABi+B
(A1)
where A and B are only a function of the Stefan number and are defined as
A=3.5Ste+1Ste,B=Ste+0.5Ste
then
Tref=(ABi+B)R2ρck=(kAhR+B)R2ρck=(RAh+R2Bk)ρc
(A2)

Hence, for a constant Stefan number, Eq. (35) indicates that the increase in Biot number monotonically decreases the reference Fourier number T̃ref. The increase in the Biot number can be due to higher values of h or R or a lower value of k. Equation (A2), however, indicates that the reference time Tref increases for larger length scales or lower values of the thermal conductivity, but decreases for higher heat transfer coefficients.

To study the effect of ΔT, PCM heat capacity c, and the latent heat of the fusion L on the reference time Tref, Eq. (35) can be rewritten as
T̃ref(Ste,Bi)=CSte+D
(A3)
where C and D are only a function of the Biot number and are defined as
C=0.5Bi+1Bi,D=Bi+0.5Bi
then
Tref=(CSte+D)R2ρck=(LCcΔT+D)R2ρck=(LCΔT+cD)R2ρk
(A4)

At a constant Biot number, T̃ref monotonically decreases with Stefan number (see Eq. (35)) which can be due to higher ΔT or thermal capacity, or lower latent heat of fusion. However, although lowering L or increasing ΔT still decreases Tref, an increase in the thermal capacity (i.e., higher Stefan number) actually increases the reference time. Furthermore, Eqs. (A2) and (A4) show that for a fixed Bi and Ste number, increasing the density of the PCM increases the time required to complete the heating or cooling process.

References

1.
Karmakar
,
A.
,
Kanani
,
Y.
,
Bhattacharya
,
A.
,
Acharya
,
S.
,
Taghizadeh
,
S.
, and
Ling
,
K.
,
2019
, “
Optimization and Analysis of a Heat Exchanger With Encapsulated Phase Change Material
,”
J. Thermophys Heat Trans.,
33(4), pp.
1161
1175
.10.2514/1.t5720
2.
Song
,
S.
,
Qiu
,
F.
,
Zhu
,
W.
,
Guo
,
Y.
,
Zhang
,
Y.
,
Ju
,
Y.
,
Feng
,
R.
,
Liu
,
Y.
,
Chen
,
Z.
,
Zhou
,
J.
,
Xiong
,
C.
, and
Dong
,
L.
,
2019
, “
Polyethylene Glycol/Halloysite@Ag Nanocomposite PCM for Thermal Energy Storage: Simultaneously High Latent Heat and Enhanced Thermal Conductivity
,”
Sol. Energy Mater. Sol. Cells
,
193
, pp.
237
245
.10.1016/j.solmat.2019.01.023
3.
Wu
,
X.
,
Zhu
,
Z.
,
Zhang
,
H.
,
Xu
,
S.
,
Fang
,
Y.
, and
Yan
,
Z.
,
2020
, “
Structural Optimization of Light-Weight Battery Module Based on Hybrid Liquid Cooling With High Latent Heat PCM
,”
Int. J. Heat Mass Transfer
,
163
, p.
120495
.10.1016/j.ijheatmasstransfer.2020.120495
4.
Castell
,
A.
, and
Solé
,
C.
,
2015
, “
An Overview on Design Methodologies for Liquid-Solid PCM Storage Systems
,”
Renewable Sustainable Energy Rev.
,
52
, pp.
289
307
.10.1016/j.rser.2015.07.119
5.
Miers
,
C. S.
, and
Marconnet
,
A.
,
2021
, “
Experimental Investigation of Composite Phase Change Material Heat Sinks for Enhanced Passive Thermal Management
,”
ASME J. Heat Transfer-Trans. ASME
,
143
(
1
), p.
013001
.10.1115/1.4048620
6.
Ali
,
H. M.
, and
Arshad
,
A.
,
2017
, “
Experimental Investigation of n-Eicosane Based Circular Pin-Fin Heat Sinks for Passive Cooling of Electronic Devices
,”
Int. J. Heat Mass Transfer
,
112
, pp.
649
661
.10.1016/j.ijheatmasstransfer.2017.05.004
7.
Colla
,
L.
,
Ercole
,
D.
,
Fedele
,
L.
,
Mancin
,
S.
,
Manca
,
O.
, and
Bobbo
,
S.
,
2017
, “
Nano-Phase Change Materials for Electronics Cooling Applications
,”
ASME J. Heat Transfer-Trans. ASME
,
139
(
5
), p.
052406
.10.1115/1.4036017
8.
Konuklu
,
Y.
,
Ostry
,
M.
,
Paksoy
,
H. O.
, and
Charvat
,
P.
,
2015
, “
Review on Using Microencapsulated Phase Change Materials (PCM) in Building Applications
,”
Energy Build.
,
106
, pp.
134
155
.10.1016/j.enbuild.2015.07.019
9.
Mofijur
,
M.
,
Mahlia
,
T.
,
Silitonga
,
A.
,
Ong
,
H.
,
Silakhori
,
M.
,
Hasan
,
M.
,
Putra
,
N.
, and
Rahman
,
S.
,
2019
, “
Phase Change Materials (PCM) for Solar Energy Usages and Storage: An Overview
,”
Energies
,
12
(
16
), p.
3167
.10.3390/en12163167
10.
Kim
,
T. Y.
,
Hyun
,
B.-S.
,
Lee
,
J.-J.
, and
Rhee
,
J.
,
2013
, “
Numerical Study of the Spacecraft Thermal Control Hardware Combining Solid–Liquid Phase Change Material and a Heat Pipe
,”
Aerosp. Sci. Technol.
,
27
(
1
), pp.
10
16
.10.1016/j.ast.2012.05.007
11.
Desai
,
A.
,
Singh
,
V.
, and
Bhavsar
,
R.
,
2017
, “
Numerical Investigation of PCM Based Thermal Control Module for Space Applications
,” Proceedings of the 24th National and Second International ISHMT-ASTFE Heat and Mass Transfer Conference (
IHMTC-2017
), Hyderabad, India, 2017, Begel House, pp.
621
628
.10.1615/IHMT C-2017.870
12.
Sarbu
,
I.
, and
Dorca
,
A.
,
2019
, “
Review on Heat Transfer Analysis in Thermal Energy Storage Using Latent Heat Storage Systems and Phase Change Materials
,”
Int. J. Energy Res.
,
43
(
1
), pp.
29
64
.10.1002/er.4196
13.
Sharma
,
A.
,
Tyagi
,
V.
,
Chen
,
C.
, and
Buddhi
,
D.
,
2009
, “
Review on Thermal Energy Storage With Phase Change Materials and Applications
,”
Renewable Sustainable Energy Rev.
,
13
(
2
), pp.
318
345
.10.1016/j.rser.2007.10.005
14.
Lohrasbi
,
S.
,
Sheikholeslami
,
M.
, and
Ganji
,
D. D.
,
2017
, “
Multi-Objective RSM Optimization of Fin Assisted Latent Heat Thermal Energy Storage System Based on Solidification Process of Phase Change Material in Presence of Copper Nanoparticles
,”
Appl. Therm. Eng.
,
118
, pp.
430
447
.10.1016/j.applthermaleng.2017.03.005
15.
Dutil
,
Y.
,
Rousse
,
D. R.
,
Salah
,
N. B.
,
Lassue
,
S.
, and
Zalewski
,
L.
,
2011
, “
A Review on Phase-Change Materials: Mathematical Modeling and Simulations
,”
Renewable Sustainable Energy Rev.
,
15
(
1
), pp.
112
130
.10.1016/j.rser.2010.06.011
16.
Ghani
,
F.
,
Waser
,
R.
,
O'Donovan
,
T. S.
,
Schuetz
,
P.
,
Zaglio
,
M.
, and
Wortischek
,
J.
,
2018
, “
Non-Linear System Identification of a Latent Heat Thermal Energy Storage System
,”
Appl. Therm. Eng.
,
134
, pp.
585
593
.10.1016/j.applthermaleng.2018.02.035
17.
Liu
,
S.
,
Li
,
Y.
, and
Zhang
,
Y.
,
2015
, “
Review on Heat Transfer Mechanisms and Characteristics in Encapsulated PCMs
,”
Heat Transfer Eng.
,
36
(
10
), pp.
880
901
.10.1080/01457632.2015.965093
18.
Zeng
,
X.
, and
Faghri
,
A.
,
1994
, “
Temperature-Transforming Model for Binary Solid-Liquid Phase-Change Problems Part I: Mathematical Modeling and Numerical Methodology
,”
Numer. Heat Transfer, Part B
,
25
(
4
), pp.
467
480
.10.1080/10407799408955931
19.
Huo
,
Y.
, and
Rao
,
Z.
,
2018
, “
The Enthalpy-Transforming-Based Lattice Boltzmann Model for Solid–Liquid Phase Change
,”
ASME J. Heat Transfer-Trans. ASME
,
140
(
10
), p.
102301
.10.1115/1.4040345
20.
Meyer
,
G. H.
,
1973
, “
Multidimensional Stefan Problems
,”
SIAM J. Numer. Anal.
,
10
(
3
), pp.
522
538
.10.1137/0710047
21.
Shamsundar
,
N.
, and
Sparrow
,
E. M.
,
1975
, “
Analysis of Multidimensional Conduction Phase Change Via the Enthalpy Model
,”
ASME J. Heat Transfer-Trans. ASME
,
97
(
3
), pp.
333
340
.10.1115/1.3450375
22.
Bonacina
,
C.
,
Comini
,
G.
,
Fasano
,
A.
, and
Primicerio
,
M.
,
1973
, “
Numerical Solution of Phase-Change Problems
,”
Int. J. Heat Mass Transfer
,
16
(
10
), pp.
1825
1832
.10.1016/0017-9310(73)90202-0
23.
Yao
,
M.
, and
Chait
,
A.
,
1993
, “
An Alternative Formulation of the Apparent Heat Capacity Method for Phase-Change Problems
,”
Numer. Heat Transfer, Part B
,
24
(
3
), pp.
279
300
.10.1080/10407799308955894
24.
Fan
,
L.-W.
,
Zhu
,
Z.-Q.
,
Liu
,
M.-J.
,
Xu
,
C.-L.
,
Zeng
,
Y.
,
Lu
,
H.
, and
Yu
,
Z.-T.
,
2016
, “
Heat Transfer During Constrained Melting of Nano-Enhanced Phase Change Materials in a Spherical Capsule: An Experimental Study
,”
ASME J. Heat Transfer-Trans. ASME
,
138
(
12
), p.
122402
.10.1115/1.4034163
25.
Zhu
,
Z.-Q.
,
Liu
,
M.-J.
,
Hu
,
N.
,
Huang
,
Y.-K.
,
Fan
,
L.-W.
,
Yu
,
Z.-T.
, and
Ge
,
J.
,
2018
, “
Inward Solidification Heat Transfer of Nano-Enhanced Phase Change Materials in a Spherical Capsule: An Experimental Study
,”
ASME J. Heat Transfer-Trans. ASME
,
140
(
2
), p.
022301
.10.1115/1.4037776
26.
Hu
,
N.
,
Zhu
,
Z.-Q.
,
Li
,
Z.-R.
,
Tu
,
J.
, and
Fan
,
L.-W.
,
2019
, “
Unconstrained Melting Heat Transfer of Nano-Enhanced Phase-Change Materials in a Spherical Capsule for Latent Heat Storage: Effects of the Capsule Size
,”
ASME J. Heat Transfer-Trans. ASME
,
141
(
7
), p.
072301
.10.1115/1.4043621
27.
Safdari
,
M.
,
Ahmadi
,
R.
, and
Sadeghzadeh
,
S.
,
2020
, “
Numerical Investigation on PCM Encapsulation Shape Used in the Passive-Active Battery Thermal Management
,”
Energy
,
193
, p.
116840
.10.1016/j.energy.2019.116840
28.
Raj
,
A.
,
Srinivas
,
M.
, and
Jayaraj
,
S.
,
2019
, “
A Cost-Effective Method to Improve the Performance of Solar Air Heaters Using Discrete Macro-Encapsulated PCM Capsules for Drying Applications
,”
Appl. Therm. Eng.
,
146
, pp.
910
920
.10.1016/j.applthermaleng.2018.10.055
29.
Orozco
,
A.
,
Hinojosa
,
J.
, and
Amaya
,
K.
,
2021
, “
The Effect of a Segmented Wall Filled With PCM on Heat Transfer and Airflow in a Closed Cavity
,”
ASME J. Heat Transfer-Trans. ASME
,
143
(
9
), p.
093001
.10.1115/1.4051600
30.
Rao
,
Z.
,
Huo
,
Y.
, and
Li
,
Y.
,
2018
, “
The Lattice Boltzmann Investigation for the Melting Process of Phase Change Material in an Inclined Cavity
,”
ASME J. Heat Transfer-Trans. ASME
,
140
(
1
), p.
012301
.10.1115/1.4037908
31.
Hu
,
N.
,
Li
,
Z.-R.
,
Zhang
,
R.-H.
,
Liu
,
J.
, and
Fan
,
L.-W.
,
2021
, “
An Experimental Investigation of Constrained Melting Heat Transfer of Nano-Enhanced Phase Change Materials in a Horizontal Cylindrical Capsule Using Thermochromic Liquid Crystal Thermography
,”
ASME J. Heat Transfer-Trans. ASME
,
143
(
1
), p.
012401
.10.1115/1.4048471
32.
Ghalambaz
,
M.
,
Mehryan
,
S.
,
Mozaffari
,
M.
,
Zadeh
,
S. M. H.
, and
Pour
,
M. S.
,
2020
, “
Study of Thermal and Hydrodynamic Characteristics of Water-Nano-Encapsulated Phase Change Particles Suspension in an Annulus of a Porous Eccentric Horizontal Cylinder
,”
Int. J. Heat Mass Transfer
,
156
, p.
119792
.10.1016/j.ijheatmasstransfer.2020.119792
33.
Balikowski
,
J. R.
, and
Mollendorf
,
J. C.
,
2007
, “
Performance of Phase Change Materials in a Horizontal Annulus of a Double-Pipe Heat Exchanger in a Water-Circulating Loop
,”
ASME J. Heat Transfer-Trans. ASME
,
129
(
3
), pp.
265
272
.10.1115/1.2426359
34.
Mahamudur Rahman
,
M.
,
Hu
,
H.
,
Shabgard
,
H.
,
Boettcher
,
P.
,
Sun
,
Y.
, and
McCarthy
,
M.
,
2016
, “
Experimental Characterization of Inward Freezing and Melting of Additive-Enhanced Phase-Change Materials Within Millimeter-Scale Cylindrical Enclosures
,”
ASME J. Heat Transfer-Trans. ASME
,
138
(
7
), p.
072301
.10.1115/1.4033007
35.
Allen
,
M. J.
,
Bergman
,
T. L.
,
Faghri
,
A.
, and
Sharifi
,
N.
,
2015
, “
Robust Heat Transfer Enhancement During Melting and Solidification of a Phase Change Material Using a Combined Heat Pipe-Metal Foam or Foil Configuration
,”
ASME J. Heat Transfer-Trans. ASME
,
137
(
10
), p.
102301
.10.1115/1.4029970
36.
Wickramaratne
,
C.
,
Moloney
,
F.
,
Pirasaci
,
T.
,
Kamal
,
R.
,
Goswami
,
D. Y.
,
Stefanakos
,
E.
, and
Dhau
,
J.
,
2016
, “
Experimental Study on Thermal Storage Performance of Cylindrically Encapsulated PCM in a Cylindrical Storage Tank With Axial Flow
,”
ASME Paper No. POWER2016-59427
.10.1115/POWER2016-59427
37.
Archibold
,
A. R.
,
Bhardwaj
,
A.
,
Rahman
,
M. M.
,
Goswami
,
D. Y.
, and
Stefanakos
,
E. L.
,
2016
, “
Comparison of Numerical and Experimental Assessment of a Latent Heat Energy Storage Module for a High-Temperature Phase-Change Material
,”
ASME J. Energy Resour. Technol.
,
138
(
5
), p.
052007
.10.1115/1.4033585
38.
Pirasaci
,
T.
,
Wickramaratne
,
C.
,
Moloney
,
F.
,
Goswami
,
D. Y.
, and
Stefanakos
,
E.
,
2017
, “
Dynamics of Phase Change in a Vertical PCM Capsule in the Presence of Radiation at High Temperatures
,”
Appl. Energy
,
206
, pp.
498
506
.10.1016/j.apenergy.2017.08.187
39.
Raj
,
A.
,
Srinivas
,
M.
, and
Jayaraj
,
S.
,
2019
, “
CFD Modeling of Macro-Encapsulated Latent Heat Storage System Used for Solar Heating Applications
,”
Int. J. Therm. Sci.
,
139
, pp.
88
104
.10.1016/j.ijthermalsci.2019.02.010
40.
Xu
,
T.
,
Chiu
,
J. N.
,
Palm
,
B.
, and
Sawalha
,
S.
,
2019
, “
Experimental Investigation on Cylindrically Macro-Encapsulated Latent Heat Storage for Space Heating Applications
,”
Energy Convers. Manage.
,
182
, pp.
166
177
.10.1016/j.enconman.2018.12.056
41.
Izgi
,
B.
, and
Arslan
,
M.
,
2020
, “
Numerical Analysis of Solidification of PCM in a Closed Vertical Cylinder for Thermal Energy Storage Applications
,”
Heat Mass Transfer
,
56
(
10
), pp.
2909
2922
.10.1007/s00231-020-02911-z
42.
Selimefendigil
,
F.
, and
Öztop
,
H. F.
,
2021
, “
Analysis of Hybrid Nanofluid and Surface Corrugation in the Laminar Convective Flow Through an Encapsulated PCM Filled Vertical Cylinder and POD-Based Modeling
,”
Int. J. Heat Mass Transfer
,
178
, p.
121623
.10.1016/j.ijheatmasstransfer.2021.121623
43.
Kumar
,
A.
, and
Saha
,
S. K.
,
2021
, “
Performance Analysis of a Packed Bed Latent Heat Thermal Energy Storage With Cylindrical-Shaped Encapsulation
,”
Int. J. Energy Res.
,
45
(
9
), pp.
13130
13148
.10.1002/er.6639
44.
Mazzeo
,
D.
, and
Oliveti
,
G.
,
2018
, “
Thermal Field and Heat Storage in a Cyclic Phase Change Process Caused by Several Moving Melting and Solidification Interfaces in the Layer
,”
Int. J. Therm. Sci.
,
129
, pp.
462
488
.10.1016/j.ijthermalsci.2017.12.026
45.
Ho
,
C. J.
, and
Chu
,
C. H.
,
1993
, “
Periodic Melting Within a Square Enclosure With an Oscillatory Surface Temperature
,”
Int. J. Heat Mass Transfer
,
36
(
3
), pp.
725
733
.10.1016/0017-9310(93)80048-Y
46.
Casano
,
G.
, and
Piva
,
S.
,
2002
, “
Experimental and Numerical Investigation of the Steady Periodic Solid-Liquid Phase-Change Heat Transfer
,”
Int. J. Heat Mass Transfer
,
45
(
20
), pp.
4181
4190
.10.1016/S0017-9310(02)00122-9
47.
Mazzeo
,
D.
,
Oliveti
,
G.
,
De Simone
,
M.
, and
Arcuri
,
N.
,
2015
, “
Analytical Model for Solidification and Melting in a Finite PCM in Steady Periodic Regime
,”
Int. J. Heat Mass Transfer
,
88
, pp.
844
861
.10.1016/j.ijheatmasstransfer.2015.04.109
48.
Kant
,
K.
,
Shukla
,
A.
,
Sharma
,
A.
, and
Biwole
,
P. H.
,
2018
, “
Melting and Solidification Behaviour of Phase Change Materials With Cyclic Heating and Cooling
,”
J. Energy Storage
,
15
, pp.
274
282
.10.1016/j.est.2017.12.005
49.
Kapilow
,
D.
,
Hsuan
,
Y. G.
,
Sun
,
Y.
, and
McCarthy
,
M.
,
2018
, “
Convective Melting and Freezing of Phase Change Materials Encapsulated Within Small Diameter Polymer Tubes
,”
Exp. Therm. Fluid Sci.
,
92
, pp.
259
269
.10.1016/j.expthermflusci.2017.11.012
50.
Voller
,
V. R.
, and
Swaminathan
,
C. R.
,
1991
, “
General Source-Based Method for Solidification Phase Change
,”
Numer. Heat Transfer, Part B
,
19
(
2
), pp.
175
189
.10.1080/10407799108944962
51.
Voller
,
V. R.
,
Swaminathan
,
C. R.
, and
Thomas
,
B. G.
,
1990
, “
Fixed Grid Techniques for Phase Change Problems: A Review
,”
Int. J. Numer. Methods Eng.
,
30
(
4
), pp.
875
898
.10.1002/nme.1620300419
52.
Voller
,
V. R.
, and
Prakash
,
C.
,
1987
, “
A Fixed Grid Numerical Modelling Methodology for Convection-Diffusion Mushy Region Phase-Change Problems
,”
Int. J. Heat Mass Transfer
,
30
(
8
), pp.
1709
1719
.10.1016/0017-9310(87)90317-6
53.
Rösler
,
F.
, and
Brüggemann
,
D.
,
2011
, “
Shell-and-Tube Type Latent Heat Thermal Energy Storage: Numerical Analysis and Comparison With Experiments
,”
Heat Mass Transfer
,
47
(
8
), pp.
1027
1033
.10.1007/s00231-011-0866-9
54.
König-Haagen
,
A.
,
Franquet
,
E.
,
Pernot
,
E.
, and
Brüggemann
,
D.
,
2017
, “
A Comprehensive Benchmark of Fixed-Grid Methods for the Modeling of Melting
,”
Int. J. Therm. Sci.
,
118
, pp.
69
103
.10.1016/j.ijthermalsci.2017.04.008
55.
Rakotondrandisa
,
A.
,
Danaila
,
I.
, and
Danaila
,
L.
,
2019
, “
Numerical Modelling of a Melting-Solidification Cycle of a Phase-Change Material With Complete or Partial Melting
,”
Int. J. Heat Fluid Flow
,
76
, pp.
57
71
.10.1016/j.ijheatfluidflow.2018.11.004
56.
Zhao
,
W.
,
Elmozughi
,
A. F.
,
Oztekin
,
A.
, and
Neti
,
S.
,
2013
, “
Heat Transfer Analysis of Encapsulated Phase Change Material for Thermal Energy Storage
,”
Int. J. Heat Mass Transfer
,
63
, pp.
323
335
.10.1016/j.ijheatmasstransfer.2013.03.061
57.
Muhammad
,
M. D.
,
Badr
,
O.
, and
Yeung
,
H.
,
2015
, “
Validation of a CFD Melting and Solidification Model for Phase Change in Vertical Cylinders
,”
Numer. Heat Transfer, Part A
,
68
(
5
), pp.
501
511
.10.1080/10407782.2014.994432
58.
Longeon
,
M.
,
Soupart
,
A.
,
Fourmigué
,
J. F.
,
Bruch
,
A.
, and
Marty
,
P.
,
2013
, “
Experimental and Numerical Study of Annular PCM Storage in the Presence of Natural Convection
,”
Appl. Energy
,
112
, pp.
175
184
.10.1016/j.apenergy.2013.06.007
59.
Paterson
,
S.
,
1952
, “
Propagation of a Boundary of Fusion
,”
Proc. Glasgow Math. Assoc.
,
1
(
1
), pp.
42
47
.10.1017/S2040618500032937
60.
Özişik
,
M. N.
, and
Uzzell
,
J. C.
,
1979
, “
Exact Solution for Freezing in Cylindrical Symmetry With Extended Freezing Temperature Range
,”
ASME J. Heat Transfer-Trans. ASME
,
101
(
2
), pp.
331
334
.10.1115/1.3450969
61.
Rubitherm Technologies GmbH
,
2018
, “
Datasheet, RT35HC
,”
Rubitherm Technologies GmbH
,
Berlin, Germany
, www.rubitherm.eu/media/products/datasheets/Techdata_-RT35HC_EN_09102020.PDF
62.
Valentova
,
K.
,
Pechackova
,
K.
,
Prikryl
,
R.
,
Ostry
,
M.
, and
Zmeskal
,
O.
,
2017
, “
Study of the Thermal Properties of Selected PCMs for Latent Heat Storage in Buildings
,”
AIP Conf. Proc.
,
1866
, p.
040042
.10.1063/1.4994522
63.
Kraiem
,
M.
,
Karkri
,
M.
,
Nasrallah
,
S. B.
,
Sobolciak
,
P.
,
Fois
,
M.
, and
Alnuaimi
,
N. A.
,
2019
, “
Thermophysical Characterization and Numerical Investigation of Three Paraffin Waxes as Latent Heat Storage Materials
,”
Preprints 2019030034.
10.20944/preprints201903.0034.v1