# Approximate Analytical Solution for One-Dimensional Solidification Problem of a Finite Superheating Phase Change Material Including the Effects of Wall and Thermal Contact Resistances

Authors: Hamid El Qarnia; Fayssal El Adnani; El Khadir Lakhal
Journal: Journal of Applied Mathematics
 (2012)

---

## Abstract

This work reports an analytical solution for the solidification of a superheating phase change material (PCM) contained in a rectangular enclosure with a finite height. The analytical solution has been obtained by solving nondimensional energy equations by using the perturbation method for a small perturbation parameter: the Stefan number,ε. This analytical solution, which takes into account the effects of the superheating of PCM, finite height of the enclosure, thickness of the wall, and wall-solid shell interfacial thermal resistances, was expressed in terms of nondimensional temperature distributions of the bottom wall of the enclosure and both PCM phases, and the dimensionless solid-liquid interface position and its dimensionless speed. The developed solution was firstly compared with that existing in the literature for the case of nonsuperheating PCM. The predicted results agreed well with those published in the literature. Next, a parametric study was carried out in order to study the impacts of the dimensionless control parameters on the dimensionless temperature distributions of the wall, the solid shell, and liquid phase of the PCM, as well as the solid-liquid interface position and its dimensionless speed.


---

## Body

## 1. Introduction

Solid-liquid phase change process (Stefan problems) in enclosures is a field of interest because of its wide range applications in engineering and great importance in academic research. Applications of such process are encountered in design of thermal storage systems, thermal comfort, electronic cooling, metallurgy, crystal growths, ice formation, freezing food, and other technologies. In contrast to heat transfer in single phase, heat transfer with melting and solidification is complex due to the presence of solid and liquid phases which are separated by a moving solid-liquid interface. Further, the position and the shape of this interface are unknown. The solid-liquid phase change process is classified as moving boundary problems. Such problems are nonlinear due to the presence of the moving solid-liquid interface. The corresponding exact analytical solutions are not always easy to obtain; exception was made for special 1D case [1]. Although these exact analytical solutions were developed, they remain, however, a good approximation to predict certain 1D solid-liquid phase change problems occurring in engineering applications. The most well-known exact solution is the Neumann’s solution to the one-dimensional problem of solidification in a semi-infinite region. Indeed, these exact analytical solutions do not take into consideration the effects of wall conduction, wall PCM interfacial thermal resistance, the superheating of PCM, and the height of the enclosure. In fact, in some technologies, the above effects exist. Therefore, a deep understanding of the solidification-melting process including such effects requires the development of the solution of equations that govern such process. This solution can predict the temperature distributions in phases and solidification/melting rate in order to control the speed of fabrication. To obtain such solutions, approximate, semianalytical, and numerical methods must be used to solve the complex solidification and melting processes. One of the techniques that can be applied in heat transfer with solidification/melting processes is the perturbation methods. These methods were widely used for solving many transient heat transfer problems.Charach et al. [2] have also developed a new perturbation method for one-dimensional and conduction-controlled melting in a planar slab. The approach is suitable for boundary conditions with small amplitude but almost arbitrary time dependence. The results are uniformly valid in time. The method is applied to the case of time-dependent heat flux boundary conditions and illustrated by the example of a step-like heat pulse. Li and Barber [3] have also employed a linear perturbation method to solve the two-dimensional heat conduction problem in which a liquid, initially at the melting temperature, becomes solidified by heat transfer to a plane mold the temperature of which is approximately uniform, but contains a small sinusoidal perturbation in one-space dimension. Results are obtained for the consequent sinusoidal perturbation in the melt/solid boundary as a function of time and for the temperature distribution throughout the solid shell. These results are expressed in terms of a series of confluent hypergeometric functions which shows good convergence for practical values of material properties and temperatures. The inverse problem, in which the melt/solid boundary is prescribed and the mold temperature is to be determined, is also briefly discussed.Aziz and Benzies [4] consider the application of regular parameter perturbation technique to obtain approximate solutions of heat transfer problems with temperature-dependent thermal properties. For pure conduction, linearly varying thermal conductivity and heat capacity are considered and two examples are solved. Next, a conducting-convecting fin is treated. Finally, a two-parameter perturbation is used to solve for temperature distribution in a conducting-convecting-radiating fin with temperature-dependent thermal conductivity. On comparing the perturbation solutions with corresponding numerical solutions, the accuracy is found to be good. Caldwell and Kwan [5] have developed perturbation methods for Stefan problems with time-dependent temperature at the boundary. The methods are applied to melting of ice in the half-plane, outward spherical solidification, and outward cylindrical solidification of a saturated liquid. The results are shown to compare well with those obtained by other numerical studies. Recently, Yigit [6] used a linear perturbation method to solve two-dimensional heat conduction during the solidification of a liquid caused by heat transfer to a plane mold of finite thickness. Heat flux drawn from the lower surface of the mold is approximately uniform, but contains a small sinusoidal perturbation in one-space dimension. Numerical results are obtained for the consequent sinusoidal perturbation in the solid/melt boundary as a function of time. Analytical results are obtained for the limiting case in which diffusivities of the solidified shell and the mold materials are infinitely large. Hwang et al. [7] used a perturbation approach for the phase-change problem during solidification involving the wall conduction and wall-material interfacial thermal contact resistance. Approximate analytical solutions were expressed in terms of dimensionless temperature distribution of the wall and both PCM phases. Comparison is made between the obtained results and those obtained in the case with neglecting the effects of wall-material interfacial thermal resistance and/or wall conduction. It was found that the solidification rate increases with decreasing wall-material interfacial thermal resistance, or increasing Stefan number, or decreasing the ratio of solidified-material thermal conductivity to wall thermal conductivity. Stephan and Holzknecht [8] developed approximate solutions for determining temperatures and phase change in solidification problems, based on perturbation method for small values of the Stefan number, ε. These solutions are applied to geometrically simple forms of the solidifying body. A comparison with a numerical solution shows that the approximate solutions of solidification on a flat plate in the range of ε<2 are sufficiently accurate for technical applications. In solidification on the outside of a cylinder the error of the solidification speed remains less than 2.6% for ε≤1 and less than 1.2% for ε≤0.5, whereas in solidification on the outside of a sphere the error remains less than 3.8% for ε≤1 and less than 1.3% for ε≤0.5. The error decreases rapidly with decreasing the Stefan number. Kharche and Howarth [9] considered the inward solidification of liquids (initially at melting temperature) in cylinders whose boundary values vary slightly with position axially and periodically along the cylinder. The problems are solved analytically, in terms of line arisen perturbations of the slight axial variation, for large Stefan number. The same authors [10] considered the inward solidification of liquids in cylinders (initially at melting temperature), with constant boundary temperatures (or heat flux), but where the boundary profile of the cylinder geometry varies slightly and periodically with the axial coordinate. Analytical solutions were also obtained, in terms of linearised perturbations of the slight axial variation, for large Stefan number. The lecturer can also refer to works carried out by Pedroso and Domoto [11, 12], Ching-Lun and Yen-Ping [13], Yan and Huang [14], Aziz and Lunardini [15], Aziz and Na [16], and Farhad Najafi and Ahmadi [17].In the knowledge of the authors, the effects of the wall conduction and the wall-solid shell interface thermal resistances on 1D solidification problem were only considered in the study conducted by Hwang et al. [7]. However, in this latest study, the effect of the initial superheating of the liquid PCM and the height of the enclosure, containing the PCM, were not taken into account. In the current work, these two limitations are overcoming by examining the 1D solidification process of an initial superheating PCM contained in an enclosure with finite height, including boundary wall conduction and interfacial wall-solid shell thermal resistances. Besides the energy equations for wall and solid-shell and their corresponding initial and boundary conditions, another energy equation with its initial and boundary conditions is included to predict the time temperature distribution of the liquid PCM. The study of this problem consists in developing analytical forms for the temperature distributions of solid-liquid phases and wall, as well as the solid-liquid interface position and its speed, by using the perturbation method. The validation of the current investigation was made by analyzing the case of saturated liquid PCM and comparing the predicted results with those reported in the literature [7]. After having validated the obtained analytical solution, a parametric study was carried out in order to examine and discuss the effects of the control parameters, such as, the superheating parameter, the dimensionless interfacial wall- solid shell thermal resistance, the ratio of solidified material thermal conductivity to wall thermal conductivity, and the Stefan number on the dimensionless temperature distributions of the wall and two PCM phases, and the speed and position of the solid-liquid interface.


## 2. Mathematical Formulation

The physical model, as presented in Figure1, consists of a rectangular enclosure which contains a PCM slab, initially at a temperature, Ti, higher than the melting temperature of the PCM, Tm. Suddenly, the temperature of the bottom wall of the enclosure is held at temperature To which is lower than the melting temperature, Tm. The height of the enclosure is, w, and the bottom wall thickness is, l. In this study, the height, w, is assumed negligible as compared to the other dimensions of the enclosure, and the effects of the edge are not taken into account. With these assumptions, the transient heat transfer conduction-dominated solidification may be considered, with a good approximation, as one-dimensional in the x-direction. It should be noticed, that as time progresses, the solid layer forms and an interface solid-liquid appears. This interface is at melting temperature, Tm, which is lower than the temperature, Ti. Therefore, the liquid phase is stable and there is no fluid motion due to free convection. This explains why the conduction heat transfer prevails in the PCM liquid. The densities of two phases are also assumed equal and independent of the temperature.Figure 1
Geometric of the physical model.The 1D heat transfer conduction-dominated solidification is governed by the energy equations for solid and liquid phases, and wall as follows:(2.1)∂2Ts∂x2=1αs∂Ts∂t0≤x≤h(t),∂2Tl∂x2=1αl∂Tl∂th(t)≤x≤w,∂2Tw∂x2=1αw∂Tw∂t-l≤x≤0.
The initial and boundary conditions are given as follows:
(2.2)t=0,Ts=Tl=Tw=T1,t>0,x=-l,Tw=To,x=0,-kw∂Tw∂x=-ks∂Ts∂x=Ts-Twr,x=h(t),Ts=Tl=Tf,ρΔH∂h∂t=ks∂Ts∂x-kl∂Tl∂x,x=w,Tl=T1.
In order to develop the dimensionless governing equations, the following nondimensional variables, listed below, are used.(2.3)θs,w,l=Ts,w,l-ToTf-To,X=xl,H=h(t)l,W=wl,τ=εαstl2,ε=cp,s(Tf-To)ΔH,R=ksrl,β=αsαw,ξ=kskw,γ=αsαl,K=klks.
The nondimensional energy equations are given as follows:
(2.4)∂2θs∂X2=εg∂θs∂H,∂2θl∂X2=εγg∂θl∂H,∂2θw∂X2=εβg∂θw∂H,
with
(2.5)g=dHdτ.The nondimensional initial and boundary conditions of (2.4) are
(2.6)H=0,θs=θl=θw=λ=T1-ToTf-To,(2.7a)H>0,θW(X=-1,H)=0,(2.7b)∂θw∂X|X=0=ξ∂θs∂X|X=0=ξθs-θwR|X=0,(2.7c)θs(X=H,H)=θl(X=H,H)=1,(2.7d)g=[∂θs∂X-K∂θl∂X]X=H,(2.7e)θl(X=W,H)=λ.As it can be seen from the analysis of the above nondimensional equations, the dimensionless control parameters are: R, λ, ε, ξ, β, γ, K, and W.


## 3. Solution Procedure

The foregoing energy equations are solved analytically using the perturbation method. It should be noted that this method requires the presence of a small perturbation parameter. In this study, the Stefan number,ε, is selected to play such a role. This dimensionless group means the relative importance of stored sensible heat to the stored latent heat. Its magnitude depends on the material and the temperature difference involved. This quantity is less than unity, for water, while for metals it varies in the range of 1–10 [18]. For materials such as n-octadecane, for example, the Stefan number is very small (ε≪1) if the temperature difference involved is ≪125°C. This condition is satisfied for low temperature applications.The perturbation solutions to the above system of equations can be expressed as follows:(3.1)θs=θso+εθs1+ε2θs2+⋯,θl=θlo+εθl1+ε2θl2+⋯,θw=θwo+εθw1+ε2θw2+⋯,(3.2)g=go+εg1+ε2g2+⋯.
The above expressions (3.1) and (3.2) are substituted into (2.4)–(2.7e). After rearranging them order by order we obtain the following equations.ε
o:
(3.3)∂2θso∂2X2=0∂2θlo∂2X2=0∂2θwo∂2X2=0go=∂θso∂X-K∂θlo∂X|X=Hθso(X=H)=1θl0(X=H)=1θlo(X=W)=λ=T1-ToTf-Toθwo(X=-1)=0∂θso∂X|X=0=θso-θwoR|X=0∂θwo∂X|X=0=ξθso-θwoR|X=0.ε
1:
(3.4)∂2θs1∂2X2=go∂θso∂H∂2θl1∂2X2=goγ∂θlo∂H∂2θw1∂2X2=goβ∂θwo∂Hg1=∂θs1∂X-K∂θl1∂X|X=Hθs1(X=H)=0θl1(X=H)=0θl1(X=W)=0θw1(X=-1)=0∂θs1∂X|X=0=θs1-θw1R|X=0∂θw1∂X|X=0=ξθs1-θw1R|X=0.ε
2:
(3.5)∂2θs2∂2X2=go∂θs1∂H+g1∂θso∂H∂2θl2∂2X2=γ[go∂θl1∂H+g1∂θlo∂H]∂2θw2∂2X2=β[go∂θw1∂H+g1∂θwo∂H]g2=∂θs2∂X-K∂θl2∂X|X=Hθs2(X=H)=0θl2(X=H)=0θl2(X=W)=0θw2(X=-1)=0∂θs2∂X|X=0=θs2-θw2R|X=0∂θw2∂X|X=0=ξθs2-θw2R|X=0a=K(1-λ),η=R+ξ.The solution for the above set of equations is given in Table1.Table 1
Solution of (2.4)-(2.5).
θ
s
,
0
X
+
η
η
+
H
θ
s
,
1
A
(
X
3
6
+
η
X
2
2
+
B
X
+
η
B
-
β
ξ
3
)
θ
s
,
2
E
(
X
5
120
+
η
X
4
24
+
F
X
3
6
+
(
η
F
-
β
ξ
3
)
X
2
2
+
G
X
-
G
H
-
I
)
θ
l
,
0
(
R
1
-
1
)
X
+
W
-
R
1
H
W
-
H
θ
l
,
1
C
(
X
3
6
-
W
X
2
2
+
D
X
-
W
D
+
W
3
3
)
θ
l
,
2
J
(
X
5
120
-
W
X
4
24
+
L
X
3
6
+
(
W
3
3
-
L
W
)
X
2
2
+
S
(
X
-
W
)
-
2
W
5
15
+
L
W
3
3
)
θ
w
,
0
ξ
(
X
+
1
)
η
+
H
θ
w
,
1
ξ
A
(
β
X
3
6
+
β
X
2
2
+
B
X
+
B
-
β
3
)
θ
w
,
2
β
ξ
E
(
β
X
5
120
+
β
X
4
24
+
F
X
3
6
+
(
F
-
β
3
)
X
2
2
-
(
ξ
η
T
+
G
H
+
I
β
η
)
X
+
R
η
T
-
G
H
+
I
β
η
)
g
o
1
η
+
H
+
a
W
-
H
g
1
-
(
1
η
+
H
+
a
W
-
H
)
(
η
H
(
η
+
H
)
2
+
β
ξ
+
H
3
3
(
η
+
H
)
3
-
γ
a
3
)
g
2
(
(
E
-
K
J
)
H
4
24
+
(
η
E
+
W
K
J
)
H
3
6
+
(
F
E
-
L
K
J
)
H
2
2
+
(
E
{
η
F
-
β
ξ
3
}
-
K
J
(
W
3
3
-
W
L
)
)
H
+
(
E
G
-
K
J
S
)
)The different coefficients appearing in the solution are given in Table2.Table 2
Coefficients appearing in the solution.
A
-
(
1
(
η
+
H
)
3
+
a
(
W
-
H
)
(
η
+
H
)
2
)
B
2
β
ξ
-
3
η
H
2
-
H
3
6
(
η
+
H
)
C
γ
R
1
-
1
(
W
-
H
)
2
(
1
(
η
+
H
)
+
a
(
W
-
H
)
)
D
2
W
(
W
+
H
)
-
H
2
6
E
(
1
η
+
H
+
a
W
-
H
)
(
3
(
η
+
H
)
4
+
a
2
W
-
η
-
3
H
(
η
+
H
)
3
(
W
-
H
)
2
)
F
B
+
(
η
H
η
+
H
+
β
ξ
+
H
3
3
(
η
+
H
)
2
)
(
1
+
a
(
(
η
+
H
)
/
(
W
-
H
)
)
3
/
(
η
+
H
)
+
a
(
(
2
W
-
η
-
3
H
)
/
(
W
-
H
)
2
)
)
+
(
η
H
/
(
η
+
H
)
2
)
+
(
(
β
ξ
+
H
3
)
/
3
(
η
+
H
)
3
)
-
a
(
γ
/
3
)
3
/
(
η
+
H
)
2
+
a
(
(
2
W
-
η
-
3
H
)
/
(
η
+
H
)
(
W
-
H
)
2
)
G
-
I
+
β
ξ
T
η
+
H
I
H
5
120
+
η
H
4
24
+
F
H
3
6
+
(
η
F
-
β
ξ
3
)
H
2
2
J
γ
2
(
R
1
-
1
)
(
W
-
H
)
3
(
1
η
+
H
+
a
W
-
H
)
(
3
H
+
2
η
-
W
(
η
+
H
)
2
+
3
a
(
W
-
H
)
)
L
2
W
(
W
+
H
)
-
H
2
6
+
(
1
(
3
H
+
2
η
-
W
)
/
(
η
+
H
)
+
3
a
(
η
+
H
)
/
(
W
-
H
)
)
×
(
W
-
H
3
+
2
a
(
η
+
H
)
3
-
1
γ
(
η
+
H
)
(
η
H
+
H
3
+
β
ξ
3
(
η
+
H
)
)
)
S
1
W
-
H
(
H
5
-
W
5
120
+
W
5
-
W
H
4
24
+
H
3
-
W
3
6
L
+
(
W
3
3
-
W
L
)
H
2
-
W
2
2
)
T
2
β
15
-
F
3As it can be seen from the above analysis, the solution for the one-dimensional solidification problem of a finite superheating PCM including the effects of wall and thermal contact resistances is obtained by using the perturbation method. The solid-liquid interface position is estimated using the Runge-Kutta method from ((2.5) and (3.2)) with the dimensionless time step, Δτ*=0.01. The solid-liquid interface speed (or solidification rate) is also calculated using (3.2). The developed solution (Table 1) is easy and simpler to use for the valuation of all interested quantities (temperature distribution of wall, solidified material, liquid PCM, solid-liquid interface, and its speed).


## 4. Results and Discussion

### 4.1. Validation

The developed mathematical model was validated by comparing the instantaneous solid-liquid interface position calculated by the obtained analytical solution and that obtained by Hwang et al. [7], for the case of λ=1 (without superheating effect). The results of such validation are displayed in Figures 2(a), 2(b), and 2(c). As it can be seen from such figures, the instantaneous solid-liquid interface position was calculated for different values of the dimensionless parameters R, ε, and ξ, respectively. The values of the other parameters are β=1.2, γ=2.75, K=0.41, and W=20. The values of γ and K were evaluated for the case of n-octadecane as a PCM. The analysis of these figures shows clearly the good agreement between the present-predicted results and those predicted by Hwang et al. [7]. The values obtained by the two models match with each other quite well. Consequently, the analytical solution obtained in the present study, for small Stefan number ε, is available to study the 1D heat conduction during solidification of a superheating PCM in a rectangular enclosure, with finite height, including the effects of bottom wall and interfacial thermal resistance between bottom wall and solid phase of PCM.Instantaneous solid-liquid interface position: (a) (ε=0.02, ξ=1.2), (b) (R=1, ξ=1.2), (c) (ε=0.02, R=1) (full symbols: present study, empty symbols: study of Hwang et al. [7]).
(a)
(b)
(c)


### 4.2. Results and Discussion

In this section, only the effects of the parametersR, λ, ε, and ξ on the thermal behavior of the system depicted in Figure 1 will be studied and analyzed. The other parameters: β=αs/αw, γ=αs/αl, K=kl/ks, and W=w/l were fixed to their respective values: 1.2, 2.75, 0.41, and 20.


#### 4.2.1. Effect of the Superheating Parameter,λ

The effect of the superheating of the liquid phase of the PCM on the solidification process is showed in Figure3, which displayed the temperature distribution in the wall and both solid and liquid phases of PCM, at dimensionless time τ*=100 (τ*=τ/ε with τ=ε(αst/l2)), for various values of the superheating parameter λ: 1, 1.5 and 2. The other parameters R, ε, and ξ are set to 1, 0.02, and 1.2, respectively. The analysis of this figure shows that for the case of λ=1 (liquid PCM not superheated, Ti=Tf), the dimensionless temperature in the liquid phase is uniform and equal to 1. As the dimensionless superheating parameter λ increases, the temperature gradient at the solid-liquid interface, in the liquid side, increases and approaches that in the solid side one. In the wall and solid PCM, however, the superheating has no effect on the dimensionless temperature distribution. As a result, the extracted heat flux from the solid-liquid interface to the wall decreases, which will lead to a decrease of the solidification rate as the superheating parameter S increases. This result agrees well with (2.7d), which expressed the dimensionless solid-liquid speed, g (or the solidification rate).Figure 3
Effect of the superheating parameter,λ, on the dimensionless temperature distribution in the system, at τ*=100, for R=1, (ε=0.02 and ξ=1.2).


#### 4.2.2. Effect of the Dimensionless Interfacial Thermal Resistance,R

Figure5 illustrates the dimensionless temperature across the wall, liquid and solid phases, at τ*=100, for different values of the dimensionless interfacial thermal resistance, R. The values of the other control parameters are given as follows: λ=2, ε=0.02, and ξ=1.2. As it is expected, a temperature discontinuity at the wall-solid PCM interface occurs. This is caused by the presence of the interfacial thermal resistance, and it is much important for larger interfacial thermal resistance. The figure also shows that as the dimensionless interfacial thermal resistance increases, the dimensionless temperatures at the solid phase and wall increase and decrease, respectively. The reason behind this is that as R increases, the dimensionless heat flux extracted from the solid PCM-wall interface diminishes, while the heat flux extracted from the solid-liquid interface of PCM remains almost constant. As a result, the recovery rate of the sensible heat dissipated by solid PCM decreases with increasing the interfacial thermal resistance, R. This leads to an increase in the dimensionless temperature of the solid PCM. The increasing in the dimensionless interfacial thermal resistance, R, has almost no effect on the heat flux dissipated by the lower boundary of the wall. However, the dimensionless heat flux received by the wall, at its upper boundary (wall-solid PCM interface), decreases with increasing, R. Consequently, the wall dissipates more heat to its environment as the dimensionless interfacial thermal resistance, R, increases. This explains the dimensionless wall temperature drop with increasing the dimensionless interfacial thermal resistance, R.Another result that can be also concluded from the analysis of this figure is that the larger solid-liquid interface position is obtained forR=1, which also corresponds to the higher solidification rate (the solid-liquid interface position corresponds to θ=1); such a result was expected. Indeed, in (2.7d), the term (K(∂θl/∂X)|H), representing dimensionless heat flux density received by the solid-liquid interface from the liquid phase of PCM, is almost constant for all values of the parameter R used in this study, while the term (∂θs/∂X|H) (the dimensionless heat flux density transmitted by the solid-liquid interface to the solid phase) decreases with increasing the parameter R. This is the reason behind the decrease of the dimensionless rate of the solidification and so the solid-liquid interface speed. The effect of the dimensionless interfacial thermal resistance, R, on the dimensionless temperature is also depicted in Figure 4, for various superheating parameters, λ. As it can be observed, higher solidification rates and solid-liquid interface velocities were obtained for lower superheating parameters, λ.Figure 4
Effect of the superheating parameter,λ, on the instantaneous dimensionless solid-liquid interface position for R=1 and 2 (ε=0.02 and ξ=1.2).Figure 5
Effect of the dimensionless interfacial thermal resistance,R, on the dimensionless temperature of the system at τ*=100 for λ=2, (ε=0.02 and ξ=1.2).


#### 4.2.3. Effect of the Stefan Number,ε

Because the dimensionless time,τ, depends on the Stefan number, ε, the dimensionless distribution temperature in the system, for different Stefan numbers, should not be displayed at this dimensionless time. The dimensionless time, τ*=τ/ε=100, which is independent of the Stefan number was, however, chosen to display dimensionless temperature across the wall, the solid and liquid phases of PCM, for Stefan numbers ε=0.02, 0.04, and 0.06, corresponding to dimensionless times τ=2, 4, and 6. The other parameters are set to the following values: λ=2, R=1, and ξ=1.2.As it can be seen in Figure6, lower dimensionless thermal gradients, (∂θs/∂X|H) and (K(∂θl/∂X)|H), and temperature of the wall and solid phase of PCM were obtained for larger values of the Stefan number, ε. Indeed, although these dimensionless thermal gradients decrease with increasing ε, the net dimensionless heat flux extracted from the solid-liquid interface (also the dimensionless velocity of the solid-liquid interface position, εg), becomes larger for greater values of ε as illustrated in Figure 7. It should be also noticed that the decrease of the net dimensionless heat flux, extracted from the solid-liquid interface, with the Stefan number, ε, leads to an increase in the recovery rate of the sensible energy of the wall and the solid PCM. This is the reason behind the decrease of the dimensionless temperature of both wall and solid PCM as ε increases as displayed in  Figure 6. The dimensionless speed of the solid-liquid interface, εg=dH/dτ*=ε[(∂θs/∂X)X=H-(K-l(∂θl/∂X))X=H], also increases as ε increases. Such a result is clearly showed in Figures 7 and 8 illustrating the temporal variations of the dimensionless solid-liquid interface speed and position, respectively, for ε=0.02 and 0.06 at different values of the superheating parameter, λ (λ=1, 1.5, and 2). Because the dimensionless solidification rate is proportional to the dimensionless velocity of the solid-liquid interface, it varies in the same manner with ε.Figure 6
Dimensionless temperature distribution in the system, for different Stefan numbers,ε, with λ=2, R=1, ξ=1.2 at τ*=100.Figure 7
Time variation of the dimensionless speed of the solid-liquid interface,εg=dH/dτ*, for ε=0.02 and 0.06 (λ=2, R=1, ξ=1.2).Figure 8
Effect of the Stefan number,ε, on the time variation of the dimensionless solid-liquid interface position for different superheating parameters, λ, with R=1 and ξ=1.2.


#### 4.2.4. Effect of the Thermal Conductivity Ratio,ξ

The effect of the thermal conductivity ratio,ξ, on the dimensionless temperature distribution of the wall and the two PCM phases at τ*=100 is displayed in Figure 9, for different values of ξ. The other parameters are kept constant at the following values: λ=2,  R=1, and ε=0.02. Figure 9 shows that larger dimensionless net-extracted heat flux from the solid-liquid interface was obtained for small parameter ξ. Indeed, small values of the parameter ξ, correspond to larger values of the thermal conductivity of the wall that enhance heat transfer between the bottom wall and its environment. As a result, the recoveries rate of sensible heat from the wall and the PCM increase, which leads to a decrease in the temperature of the wall and PCM, and an increase in the solidification rate as illustrated in Figure 10.Figure 9
Dimensionless temperature distribution in the system, for different thermal conductivity ratio,ξ,  at τ*=100, with λ=2, R=1, ε=0.02.Figure 10
Effect of the thermal conductivity ratio,ξ, and superheating parameter, λ, on the time variation of the dimensionless thickness of the solid PCM layer (R=1, ε=0.02).


## 4.1. Validation

The developed mathematical model was validated by comparing the instantaneous solid-liquid interface position calculated by the obtained analytical solution and that obtained by Hwang et al. [7], for the case of λ=1 (without superheating effect). The results of such validation are displayed in Figures 2(a), 2(b), and 2(c). As it can be seen from such figures, the instantaneous solid-liquid interface position was calculated for different values of the dimensionless parameters R, ε, and ξ, respectively. The values of the other parameters are β=1.2, γ=2.75, K=0.41, and W=20. The values of γ and K were evaluated for the case of n-octadecane as a PCM. The analysis of these figures shows clearly the good agreement between the present-predicted results and those predicted by Hwang et al. [7]. The values obtained by the two models match with each other quite well. Consequently, the analytical solution obtained in the present study, for small Stefan number ε, is available to study the 1D heat conduction during solidification of a superheating PCM in a rectangular enclosure, with finite height, including the effects of bottom wall and interfacial thermal resistance between bottom wall and solid phase of PCM.Instantaneous solid-liquid interface position: (a) (ε=0.02, ξ=1.2), (b) (R=1, ξ=1.2), (c) (ε=0.02, R=1) (full symbols: present study, empty symbols: study of Hwang et al. [7]).
(a)
(b)
(c)


## 4.2. Results and Discussion

In this section, only the effects of the parametersR, λ, ε, and ξ on the thermal behavior of the system depicted in Figure 1 will be studied and analyzed. The other parameters: β=αs/αw, γ=αs/αl, K=kl/ks, and W=w/l were fixed to their respective values: 1.2, 2.75, 0.41, and 20.


### 4.2.1. Effect of the Superheating Parameter,λ

The effect of the superheating of the liquid phase of the PCM on the solidification process is showed in Figure3, which displayed the temperature distribution in the wall and both solid and liquid phases of PCM, at dimensionless time τ*=100 (τ*=τ/ε with τ=ε(αst/l2)), for various values of the superheating parameter λ: 1, 1.5 and 2. The other parameters R, ε, and ξ are set to 1, 0.02, and 1.2, respectively. The analysis of this figure shows that for the case of λ=1 (liquid PCM not superheated, Ti=Tf), the dimensionless temperature in the liquid phase is uniform and equal to 1. As the dimensionless superheating parameter λ increases, the temperature gradient at the solid-liquid interface, in the liquid side, increases and approaches that in the solid side one. In the wall and solid PCM, however, the superheating has no effect on the dimensionless temperature distribution. As a result, the extracted heat flux from the solid-liquid interface to the wall decreases, which will lead to a decrease of the solidification rate as the superheating parameter S increases. This result agrees well with (2.7d), which expressed the dimensionless solid-liquid speed, g (or the solidification rate).Figure 3
Effect of the superheating parameter,λ, on the dimensionless temperature distribution in the system, at τ*=100, for R=1, (ε=0.02 and ξ=1.2).


### 4.2.2. Effect of the Dimensionless Interfacial Thermal Resistance,R

Figure5 illustrates the dimensionless temperature across the wall, liquid and solid phases, at τ*=100, for different values of the dimensionless interfacial thermal resistance, R. The values of the other control parameters are given as follows: λ=2, ε=0.02, and ξ=1.2. As it is expected, a temperature discontinuity at the wall-solid PCM interface occurs. This is caused by the presence of the interfacial thermal resistance, and it is much important for larger interfacial thermal resistance. The figure also shows that as the dimensionless interfacial thermal resistance increases, the dimensionless temperatures at the solid phase and wall increase and decrease, respectively. The reason behind this is that as R increases, the dimensionless heat flux extracted from the solid PCM-wall interface diminishes, while the heat flux extracted from the solid-liquid interface of PCM remains almost constant. As a result, the recovery rate of the sensible heat dissipated by solid PCM decreases with increasing the interfacial thermal resistance, R. This leads to an increase in the dimensionless temperature of the solid PCM. The increasing in the dimensionless interfacial thermal resistance, R, has almost no effect on the heat flux dissipated by the lower boundary of the wall. However, the dimensionless heat flux received by the wall, at its upper boundary (wall-solid PCM interface), decreases with increasing, R. Consequently, the wall dissipates more heat to its environment as the dimensionless interfacial thermal resistance, R, increases. This explains the dimensionless wall temperature drop with increasing the dimensionless interfacial thermal resistance, R.Another result that can be also concluded from the analysis of this figure is that the larger solid-liquid interface position is obtained forR=1, which also corresponds to the higher solidification rate (the solid-liquid interface position corresponds to θ=1); such a result was expected. Indeed, in (2.7d), the term (K(∂θl/∂X)|H), representing dimensionless heat flux density received by the solid-liquid interface from the liquid phase of PCM, is almost constant for all values of the parameter R used in this study, while the term (∂θs/∂X|H) (the dimensionless heat flux density transmitted by the solid-liquid interface to the solid phase) decreases with increasing the parameter R. This is the reason behind the decrease of the dimensionless rate of the solidification and so the solid-liquid interface speed. The effect of the dimensionless interfacial thermal resistance, R, on the dimensionless temperature is also depicted in Figure 4, for various superheating parameters, λ. As it can be observed, higher solidification rates and solid-liquid interface velocities were obtained for lower superheating parameters, λ.Figure 4
Effect of the superheating parameter,λ, on the instantaneous dimensionless solid-liquid interface position for R=1 and 2 (ε=0.02 and ξ=1.2).Figure 5
Effect of the dimensionless interfacial thermal resistance,R, on the dimensionless temperature of the system at τ*=100 for λ=2, (ε=0.02 and ξ=1.2).


### 4.2.3. Effect of the Stefan Number,ε

Because the dimensionless time,τ, depends on the Stefan number, ε, the dimensionless distribution temperature in the system, for different Stefan numbers, should not be displayed at this dimensionless time. The dimensionless time, τ*=τ/ε=100, which is independent of the Stefan number was, however, chosen to display dimensionless temperature across the wall, the solid and liquid phases of PCM, for Stefan numbers ε=0.02, 0.04, and 0.06, corresponding to dimensionless times τ=2, 4, and 6. The other parameters are set to the following values: λ=2, R=1, and ξ=1.2.As it can be seen in Figure6, lower dimensionless thermal gradients, (∂θs/∂X|H) and (K(∂θl/∂X)|H), and temperature of the wall and solid phase of PCM were obtained for larger values of the Stefan number, ε. Indeed, although these dimensionless thermal gradients decrease with increasing ε, the net dimensionless heat flux extracted from the solid-liquid interface (also the dimensionless velocity of the solid-liquid interface position, εg), becomes larger for greater values of ε as illustrated in Figure 7. It should be also noticed that the decrease of the net dimensionless heat flux, extracted from the solid-liquid interface, with the Stefan number, ε, leads to an increase in the recovery rate of the sensible energy of the wall and the solid PCM. This is the reason behind the decrease of the dimensionless temperature of both wall and solid PCM as ε increases as displayed in  Figure 6. The dimensionless speed of the solid-liquid interface, εg=dH/dτ*=ε[(∂θs/∂X)X=H-(K-l(∂θl/∂X))X=H], also increases as ε increases. Such a result is clearly showed in Figures 7 and 8 illustrating the temporal variations of the dimensionless solid-liquid interface speed and position, respectively, for ε=0.02 and 0.06 at different values of the superheating parameter, λ (λ=1, 1.5, and 2). Because the dimensionless solidification rate is proportional to the dimensionless velocity of the solid-liquid interface, it varies in the same manner with ε.Figure 6
Dimensionless temperature distribution in the system, for different Stefan numbers,ε, with λ=2, R=1, ξ=1.2 at τ*=100.Figure 7
Time variation of the dimensionless speed of the solid-liquid interface,εg=dH/dτ*, for ε=0.02 and 0.06 (λ=2, R=1, ξ=1.2).Figure 8
Effect of the Stefan number,ε, on the time variation of the dimensionless solid-liquid interface position for different superheating parameters, λ, with R=1 and ξ=1.2.


### 4.2.4. Effect of the Thermal Conductivity Ratio,ξ

The effect of the thermal conductivity ratio,ξ, on the dimensionless temperature distribution of the wall and the two PCM phases at τ*=100 is displayed in Figure 9, for different values of ξ. The other parameters are kept constant at the following values: λ=2,  R=1, and ε=0.02. Figure 9 shows that larger dimensionless net-extracted heat flux from the solid-liquid interface was obtained for small parameter ξ. Indeed, small values of the parameter ξ, correspond to larger values of the thermal conductivity of the wall that enhance heat transfer between the bottom wall and its environment. As a result, the recoveries rate of sensible heat from the wall and the PCM increase, which leads to a decrease in the temperature of the wall and PCM, and an increase in the solidification rate as illustrated in Figure 10.Figure 9
Dimensionless temperature distribution in the system, for different thermal conductivity ratio,ξ,  at τ*=100, with λ=2, R=1, ε=0.02.Figure 10
Effect of the thermal conductivity ratio,ξ, and superheating parameter, λ, on the time variation of the dimensionless thickness of the solid PCM layer (R=1, ε=0.02).


## 4.2.1. Effect of the Superheating Parameter,λ

The effect of the superheating of the liquid phase of the PCM on the solidification process is showed in Figure3, which displayed the temperature distribution in the wall and both solid and liquid phases of PCM, at dimensionless time τ*=100 (τ*=τ/ε with τ=ε(αst/l2)), for various values of the superheating parameter λ: 1, 1.5 and 2. The other parameters R, ε, and ξ are set to 1, 0.02, and 1.2, respectively. The analysis of this figure shows that for the case of λ=1 (liquid PCM not superheated, Ti=Tf), the dimensionless temperature in the liquid phase is uniform and equal to 1. As the dimensionless superheating parameter λ increases, the temperature gradient at the solid-liquid interface, in the liquid side, increases and approaches that in the solid side one. In the wall and solid PCM, however, the superheating has no effect on the dimensionless temperature distribution. As a result, the extracted heat flux from the solid-liquid interface to the wall decreases, which will lead to a decrease of the solidification rate as the superheating parameter S increases. This result agrees well with (2.7d), which expressed the dimensionless solid-liquid speed, g (or the solidification rate).Figure 3
Effect of the superheating parameter,λ, on the dimensionless temperature distribution in the system, at τ*=100, for R=1, (ε=0.02 and ξ=1.2).


## 4.2.2. Effect of the Dimensionless Interfacial Thermal Resistance,R

Figure5 illustrates the dimensionless temperature across the wall, liquid and solid phases, at τ*=100, for different values of the dimensionless interfacial thermal resistance, R. The values of the other control parameters are given as follows: λ=2, ε=0.02, and ξ=1.2. As it is expected, a temperature discontinuity at the wall-solid PCM interface occurs. This is caused by the presence of the interfacial thermal resistance, and it is much important for larger interfacial thermal resistance. The figure also shows that as the dimensionless interfacial thermal resistance increases, the dimensionless temperatures at the solid phase and wall increase and decrease, respectively. The reason behind this is that as R increases, the dimensionless heat flux extracted from the solid PCM-wall interface diminishes, while the heat flux extracted from the solid-liquid interface of PCM remains almost constant. As a result, the recovery rate of the sensible heat dissipated by solid PCM decreases with increasing the interfacial thermal resistance, R. This leads to an increase in the dimensionless temperature of the solid PCM. The increasing in the dimensionless interfacial thermal resistance, R, has almost no effect on the heat flux dissipated by the lower boundary of the wall. However, the dimensionless heat flux received by the wall, at its upper boundary (wall-solid PCM interface), decreases with increasing, R. Consequently, the wall dissipates more heat to its environment as the dimensionless interfacial thermal resistance, R, increases. This explains the dimensionless wall temperature drop with increasing the dimensionless interfacial thermal resistance, R.Another result that can be also concluded from the analysis of this figure is that the larger solid-liquid interface position is obtained forR=1, which also corresponds to the higher solidification rate (the solid-liquid interface position corresponds to θ=1); such a result was expected. Indeed, in (2.7d), the term (K(∂θl/∂X)|H), representing dimensionless heat flux density received by the solid-liquid interface from the liquid phase of PCM, is almost constant for all values of the parameter R used in this study, while the term (∂θs/∂X|H) (the dimensionless heat flux density transmitted by the solid-liquid interface to the solid phase) decreases with increasing the parameter R. This is the reason behind the decrease of the dimensionless rate of the solidification and so the solid-liquid interface speed. The effect of the dimensionless interfacial thermal resistance, R, on the dimensionless temperature is also depicted in Figure 4, for various superheating parameters, λ. As it can be observed, higher solidification rates and solid-liquid interface velocities were obtained for lower superheating parameters, λ.Figure 4
Effect of the superheating parameter,λ, on the instantaneous dimensionless solid-liquid interface position for R=1 and 2 (ε=0.02 and ξ=1.2).Figure 5
Effect of the dimensionless interfacial thermal resistance,R, on the dimensionless temperature of the system at τ*=100 for λ=2, (ε=0.02 and ξ=1.2).


## 4.2.3. Effect of the Stefan Number,ε

Because the dimensionless time,τ, depends on the Stefan number, ε, the dimensionless distribution temperature in the system, for different Stefan numbers, should not be displayed at this dimensionless time. The dimensionless time, τ*=τ/ε=100, which is independent of the Stefan number was, however, chosen to display dimensionless temperature across the wall, the solid and liquid phases of PCM, for Stefan numbers ε=0.02, 0.04, and 0.06, corresponding to dimensionless times τ=2, 4, and 6. The other parameters are set to the following values: λ=2, R=1, and ξ=1.2.As it can be seen in Figure6, lower dimensionless thermal gradients, (∂θs/∂X|H) and (K(∂θl/∂X)|H), and temperature of the wall and solid phase of PCM were obtained for larger values of the Stefan number, ε. Indeed, although these dimensionless thermal gradients decrease with increasing ε, the net dimensionless heat flux extracted from the solid-liquid interface (also the dimensionless velocity of the solid-liquid interface position, εg), becomes larger for greater values of ε as illustrated in Figure 7. It should be also noticed that the decrease of the net dimensionless heat flux, extracted from the solid-liquid interface, with the Stefan number, ε, leads to an increase in the recovery rate of the sensible energy of the wall and the solid PCM. This is the reason behind the decrease of the dimensionless temperature of both wall and solid PCM as ε increases as displayed in  Figure 6. The dimensionless speed of the solid-liquid interface, εg=dH/dτ*=ε[(∂θs/∂X)X=H-(K-l(∂θl/∂X))X=H], also increases as ε increases. Such a result is clearly showed in Figures 7 and 8 illustrating the temporal variations of the dimensionless solid-liquid interface speed and position, respectively, for ε=0.02 and 0.06 at different values of the superheating parameter, λ (λ=1, 1.5, and 2). Because the dimensionless solidification rate is proportional to the dimensionless velocity of the solid-liquid interface, it varies in the same manner with ε.Figure 6
Dimensionless temperature distribution in the system, for different Stefan numbers,ε, with λ=2, R=1, ξ=1.2 at τ*=100.Figure 7
Time variation of the dimensionless speed of the solid-liquid interface,εg=dH/dτ*, for ε=0.02 and 0.06 (λ=2, R=1, ξ=1.2).Figure 8
Effect of the Stefan number,ε, on the time variation of the dimensionless solid-liquid interface position for different superheating parameters, λ, with R=1 and ξ=1.2.


## 4.2.4. Effect of the Thermal Conductivity Ratio,ξ

The effect of the thermal conductivity ratio,ξ, on the dimensionless temperature distribution of the wall and the two PCM phases at τ*=100 is displayed in Figure 9, for different values of ξ. The other parameters are kept constant at the following values: λ=2,  R=1, and ε=0.02. Figure 9 shows that larger dimensionless net-extracted heat flux from the solid-liquid interface was obtained for small parameter ξ. Indeed, small values of the parameter ξ, correspond to larger values of the thermal conductivity of the wall that enhance heat transfer between the bottom wall and its environment. As a result, the recoveries rate of sensible heat from the wall and the PCM increase, which leads to a decrease in the temperature of the wall and PCM, and an increase in the solidification rate as illustrated in Figure 10.Figure 9
Dimensionless temperature distribution in the system, for different thermal conductivity ratio,ξ,  at τ*=100, with λ=2, R=1, ε=0.02.Figure 10
Effect of the thermal conductivity ratio,ξ, and superheating parameter, λ, on the time variation of the dimensionless thickness of the solid PCM layer (R=1, ε=0.02).


## 5. Conclusions

The present study shows clearly the efficiency of the perturbation method to solve one dimensional solidification problem and to calculate quantities that characterize such process. The analysis carried out in this work allows the following concluding remarks.(i)
Higher liquid PCM dimensionless temperature gradient at the solid-liquid interface are obtained for larger superheating parameter,λ, and it approaches dimensionless solid PCM temperature gradient at the solid-liquid interface, as λ increases.(ii)
In the range of the superheating parameter,λ, explored in this study, the dimensionless temperature distribution of the solid PCM and wall are slightly affected by the increasing of the superheating parameter, λ.(iii)
Lower dimensionless solidification rate is obtained for larger superheating parameter,λ.(iv)
Larger dimensionless solidification rate is obtained for lower dimensionless interfacial thermal resistance,R.(v)
An increase in the dimensionless interfacial thermal resistance,R, leads to an increase of the dimensionless solid PCM temperature and a decrease of the dimensionless bottom wall temperature.(vi)
An increase in the Stefan number,ε, implies a decrease of the dimensionless temperatures of the wall and solid PCM.(vii)
Larger dimensionless solidification rates are obtained for higher Stefan number,ε.(viii)
A decrease in the dimensionless temperature of the wall and the solid PCM occurs as the thermal conductivity ratio,ξ, decreases.(ix)
Larger dimensionless solidification rates are obtained for lower thermal conductivity ratio,ξ.



---

*Source: 174604-2012-11-14.xml*