Melting characteristics of organic phase change material in a wavy trapezoidal cavity

Journal Pre-proof Melting characteristics of organic phase change material in a wavy trapezoidal cavity

Zoubida Haddad, Farida Iachachene PII:

S0167-7322(19)34881-0

DOI:

https://doi.org/10.1016/j.molliq.2019.112132

Reference:

MOLLIQ 112132

To appear in:

Journal of Molecular Liquids

Received date:

31 August 2019

Revised date:

9 November 2019

Accepted date:

12 November 2019

Please cite this article as: Z. Haddad and F. Iachachene, Melting characteristics of organic phase change material in a wavy trapezoidal cavity, Journal of Molecular Liquids(2019), https://doi.org/10.1016/j.molliq.2019.112132

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Published by Elsevier.

Journal Pre-proof

Melting characteristics of organic phase change material in a wavy trapezoidal cavity Zoubida Haddad1, Farida Iachachene2 Institute of Electrical and Electronic Engineering, University M’Hamed BOUGARA of Boumerdes, 35000 Boumerdes, Algeria 2 Department of physics, faculty of sciences, University M’Hamed Bougara Boumerdes, Algeria 1

Abstract Melting of n-eicosane in a wavy trapezoidal cavity heated from below is investigated

of

numerically based on the enthalpy-porosity technique. The physical properties of n-eicosane are examined in details, and the impact of the amplitude of the wavy wall on the heat transfer

ro

and flow is discussed. The wave amplitude ranging is from 0 mm to 1.5 mm and the temperature difference range of 20-50°C. The obtained results show that during the melting

-p

process, at low temperature difference, the rolls merge together, leading to a pair of larger and

re

counter-rotating rolls. However, increasing the temperature difference for a = 1 mm and a = 1.5 mm causes the flow to become chaotic in time and space, where the rolls become

lP

irregular. The study also shows that the number of convection cells is time and wave dependent. Furthermore, it is observed that the melting time in the flat surface is about 2.5

na

times faster than that in the wavy surface at ΔT = 50°C.

Jo ur

Keywords: Melting, Rayleigh-Bénard convection, wavy trapezoidal cavity, heat transfer, phase change material, numerical simulation.

1

Journal Pre-proof

Nomenclatures wave amplitude (mm)

Amush

Mushy zone constant (kg/m3s)

Cp

Specific heat (J/kg·K)

f

Liquid fraction

h

Height of the cavity (m)

H

Sensible enthalpy (J/kg)

h

Latent heat ([J/kg)

k

Thermal conductivity (W/m·K)

lf

Latent heat of fusion ( J/kg)

P

Pressure (Pa)

Q

Heat flux (W/m2)

Ra

Rayleith number (-)

S

Source term

t

Time (s)

T

Temperature (K)

Tf

Melting temperature (K)

x,y

Cartesian coordinates (m)

U,V

Velocity components (m·s-1)

ro -p re

lP

na

Greek symbols

Dynamic viscosity (Ns/m2)



Fluid density (kg·m-3)

 Subscript PCM f

Jo ur

 ΔT

of

a

Temperature difference (°C) thermal expansion coefficient (k-1)

Phase change material Fusion

s

Solid

l

Liquid

2

Journal Pre-proof 1. Introduction Nowadays, several technologies are emerged to reduce energy usage, hence, create an effective cooling system such a development of modern electronic devices, thermal insulation materials in the building, and energy apparatus [1-7]. One type of material that meets the aforementioned purpose is Phase Change Materials (PCMs). The interest in PCMs derives from their interesting properties against conventional materials [8]. PCMs have many sciences and engineering applications such as electronics cooling, heating, thermal storage systems, and textiles [9-13]. Various numerical and experimental studies have been devoted to convection

of

(conduction)-dominated melting (solidification) of PCMs in different shapes like rectangular, spherical, cylindrical, triangular, trapezoidal and annular containers [14]. In a pioneer work

ro

done by Gau and Viskanta [15], the effect of buoyancy-driven gallium flow in a rectangular

-p

cavity is experimentally investigated. They found that natural convection, which is present during phase change, accelerates the melting rate. Later this work was validated numerically

re

by Brent et al. [16], using the Enthalpy-porosity technique. Both studies have been taken as a basis by many researchers in their experimental and theoretical studies. El Omari et al. [17]

lP

analyzed numerically a passive cooling system using enclosures with different geometries filled with thermal conductivity-enhanced PCM. Five geometries containing the same volume

na

of PCM are compared while cooling the same surface. Their results showed that the best efficiency is obtained for an enclosure shifted vertically relative to the cooled surface.

Jo ur

Hosseinizadeh et al. [18] explored numerically and experimentally the unconstrained melting of n-octadecane paraffin wax as PCM encapsulated in a hot spherical container with different operating conditions, and various diameters. It was found that although the Stefan number is an important parameter in the melting process, the results indicate that the geometrical parameter such as the sphere diameter has more influence on the melting rate and heat transfer compared to the operating conditions. Shmueli et al. [19] performed a numerical investigation of PCM melting process in vertical circular tubes, insulated from bottom and heated from sides, while the upper part is exposed to air. They demonstrated that at the beginning of the process, heat transfer is governed by conduction. However, as melting progresses, the natural convection effect becomes more dominant, Darzi et al. [20] presented a numerical study of PCM melting inside concentric and eccentric horizontal cylindrical annulus. The inner cylindrical tube is heated isothermally, while the outer tube is insulated. They found that melt rate in the top half becomes faster than the bottom half of the horizontal annulus. Furthermore, their results revealed that when the inner cylinder tube is moved downward, the 3

Journal Pre-proof melting rate increases significantly, which is attributed to the effect of eccentricity. Sharma et al. [21] simulated the solidification process of copper–water nanofluid in an isosceles trapezoidal cavity. The results revealed that the heat transfer performance of Nano-Enhanced PCM is significantly enhanced with trapezoidal cavity when compared to a square cavity having the same internal area. Arici et al. [22] undertook a numerical study on the melting of paraffin wax with Al2O3 nanoparticles in a partially heated and cooled square cavity. They illustrated that the heat energy stored by PCM can be enhanced by changing orientation of the cavity, dispersing nanoparticles or applying both simultaneously. Zeng et al. [23] investigated the melting process of lauric acid within a rectangular geometry at five different orientations.

of

It was found that the melting time could be significantly affected by changing orientation of the heat exchange surface. Dhaidan [24] studied numerically the melting of n-eicosane filled

ro

in a triangular container with two different orientations, where one storage cavity is a lower

-p

base container and the other is the upper-base container. The results unfold that melt rate is higher in the upper base than that in the lower base. In addition, it was found that melt rate of

re

the upper part of the container was higher than that recorded in the lower part for both orientations. Iachachene et al. [25] investigated the orientation and nanoparticles effects on

lP

the melting of paraffin wax filled in a trapezoidal cavity. It was observed that both effects are beneficial for heat transfer enhancements in PCMs. However, Nano-Enhanced PCM lead to a

na

lower heat transfer performance when the nanofluid thermal conductivity enhancement was less than 80%. Hoseinzadeh and Chamkha [26] studied numerically the melting process of

Jo ur

various PCMs embedded in a rectangular cavity. Air is used as the heat transfer fluid and CaCl2. 6H2O and RT25 as the PCMs. Their results showed that when the flow rate is increased, the melting rate and outlet air temperature are increased. Gong and Mujumdar [27] studied numerically the melting of PCM in a rectangular cavity heated from below. It was observed that complex and time-dependent flow patterns are obtained at high Rayleigh numbers. Fteiti and Nasrallah [28] investigated the melting process in a rectangular enclosure heated from below and cooled from above. They found that the melting process is more rapid as the aspect ratio of the enclosure is small, but the volume fraction melted at the steady state is less important. Kousksou et al. [29] conducted a numerical study to analyze the melting in a rectangular closed enclosure by subjecting the bottom wavy surface to a uniform temperature. The effect of the amplitude of the wavy surface on the flow structure, and heat transfer characteristics is investigated at Rayleigh number of 6×105. It was found that melting rate increases with increasing the amplitude of wavy surface.

4

Journal Pre-proof To the best of the authors’ knowledge, the impact of wave amplitude for different temperature differences (i.e, Rayleigh numbers) on natural convection with PCM in cavities cavity has not been addressed yet. The motivation of the present research work is first, to better understand the Rayleigh Bénard with phase change for the design and optimisation of phase change systems for energy storage or cooling devices, and second, to find a relationship between the wave amplitude and temperature difference. In order to perform accurate numerical simulations of the melting of n-eicosane in a trapezoidal cavity, a comprehensive examination of its physical properties has been performed. Correlations are presented along with comments about their validity. The effects of different parameters, like the temperature

of

difference and wave amplitude, on the melting rate, heat transfer characteristics, and flow structure evolution during the melting process are examined in details.

ro

2. Physical Model

-p

Schematic view of the physical model is shown in Fig. 1. It consists of a trapezoidal enclosure with bottom wavy wall. The length of the top wall and the height of the enclosure

re

are both equal to h =1 cm. The left and right walls are each inclined at an angle of 45° from the y-axis. The wavy wall has length 3h, and the wave amplitude is varied in the range of a =

lP

0 mm to a = 1.5 mm. The sinusoidal shape of the wall is described by the following equation: 0  x  3h

(1)

na

   h  x    f ( x)  a cos  2        h  2 

The cavity is filled with n-eicosane, which is assumed to be initially at a temperature below

Jo ur

the fusion point. The bottom wall is kept at constant temperature TH higher than melting temperature, while the other walls are maintained adiabatic. y

A

c ati ab i d

Adiabatic

h

Ad

iab

ati c

Isotherme surface

h x

Figure 1: Physical model

3. Mathematical formulation The heat transfer in the cavity takes place with the effect of conduction as well convection in solid and liquid phase, respectively. The melting of n-eicosane is assumed to be Newtonian and incompressible. The flow caused due to the melting is laminar and the viscous

5

Journal Pre-proof dissipations,

thermal

radiation,

and

three-dimensional

convection

are

negligible.

Thermophysical properties of PCM are both constant and temperature dependent. Considering the above mentioned assumptions, the governing equations for conservation of mass, momentum and energy can be written as [25]:  ( U )

 ( U ) t  ( V ) t  ( H ) t

t   



 ( V ) t

 (  UU ) x  ( UV ) x

(2)

 ( UV ) y

 ( VV )



 ( UH ) x



0

y







 ( UH ) y

P y



P x 



U    U    x  x  y  y  

 

V 

 

   Sx 

V 

     (  ) g  T  Tm   S y x  x  y  y 

  T    T  k   k  x  x  y  y 

of

t



(3) (1)

) (4) (5)

ro



Sx 

(1  f )2 f 3 

Amush U , S y 

(1  f )2 f 3 

Amush V

-p

In equations (3) and (4), S represents source term which is defined as [25]: (6)

re

where ε=0.001 to prevent division by zero, Amush is the mushy zone constant and f is the liquid volume fraction.

T

c

Tref

dT

(7) (8) (9)

Jo ur

H  fL f

p

na

h  href 

lP

Enthalpy, H is the sum of sensible enthalpy, h , and the latent heat of fusion, H : H  h  H

The enthalpy formulation requires a single domain in which the same set of governing equations are used to model both solid and liquid phases of a PCM. The transition from solid to liquid, and vice versa, occurs over a finite temperature range (ΔTf = 1°C) generating an artificial mushy region at the solid-liquid interface. The fluid velocity within the mushy region varies from zero (at the solid boundary) to the natural convection velocity (at the liquid boundary) as the melt fraction varies from 0 to 1. In both cases, phase change is quantified through the following equation for the melt fraction [25]:  0   T  Ts f   Tl  Ts  1

, T  Ts , Ts  T  Tl

(10)

, T  Tl

6

Journal Pre-proof Different parameters are used in Eqs. (2)-(5) to represent the thermophysical properties of liquid and solid phases of PCM. Detailed descriptions and expressions of these parameters are described in the following subsection. 3.1. Physical properties Accurate knowledge of thermo-physical properties of phase change materials in the desired operation temperature is a prerequisite for any reliable utilization of the materials. However, some difficulties have been encountered in the search of thermophysical properties of neicosane due to dispersion and lack of data in both liquid and solid phases. The thermophysical properties were taken from different sources. The significant difference between

of

experimental and numerical results leads to the suspicion of being using wrong properties. Therefore, a thorough search of physical-properties has been carried out. It was found that

ro

very little information has been found for some of them, such as density and thermal

presented in the literature is reviewed.

re

3.1.1. Dynamic viscosity

-p

conductivity of solid n-eicosane. In what follows, the physical properties of n-eicosane

To the best of our knowledge, there exists only one dynamic viscosity correlation reported

9.2095  1.82  103 for 310 k < T < 767 K T  0.0168 T  1.29  105 T 2

(11)

na

log10  

lP

in the literature, which is expressed as follows [30]

The dynamic viscosity is in units of centipoise, while the temperature in K.

Jo ur

3.1.2. Thermal conductivity

Nabil and Khodadadi [31] observed that the thermal conductivity of n-eicosane is independent of temperature in the range from 283.5 K to 306 K with an average value of 0.4235 W/m·K. However, an increase by 30% was observed at 308.15 K. This result matches the experimental data of Stryker and Sparrow [32] very closely, where they found that the thermal conductivity is independent of temperature for the 283.15 K–305.65 K range. The average value for these points is 0.4234 W/m·k. However, it was found that the conductivity drops off at T > 305.65 K. Velez et al. [33] measured the thermal conductivity of n-eicosane in the range 258–348 K. They found that the thermal conductivity is maintained nearly constant for temperatures much lower than the melting temperature and decreased gradually with temperatures near the melting temperature. Based on the experimental data of Velez et al. [33], a thermal conductivity correlation of neicosane in its solid and liquid states was derived. A two-dimensional regression is performed to develop such correlations. 7

Journal Pre-proof 

Solid thermal conductivity

The R2 of the regression is 99% and the maximum relative error is 2.9%. The correlation for solid phase is given as: k s  0.45964  7.83234  1012 et / 12.86634 for 258


(12)

Liquid thermal conductivity

The R2 of the regression is 99.4% and the maximum relative error is 0.25%. The correlation for liquid phase is given as: kl  0.45964  7.83234  104 T for 308.84
(13)

of

Yaws [30] also proposed a thermal conductivity correlation of n-eicosane at temperatures between 308.8 K and 620 K, and it is expressed as

ro

kl  0.20681  1.9990  104 T  2.6667  108 T 2

(14)

-p

Where k is expressed in W/m·K, and the temperature in K. 3.1.3. Density Solid density

re



A very few data were reported on the density of solid eicosane. Stryker and Sparrow [32]

lP

observed that the density is independent of the solidification temperature. The mean value is 840 kg m-3 with maximum variation of 1.1%. This result is in the same range as those

na

reported in the literature, which includes 856 kg·m at 307.15 K [34]. 910 kg·m-3 at 283.15 K and 856 kg·m-3 at 307.15 K [35], and 848.95 kg·m-3 at 309.54 K [36]. Liquid density

Jo ur



Velez et al. [33] measured the density of liquid n-eicosane for different temperatures. A gradual decrease of the density is observed as the temperature increases. They proposed the following expression for the density up to 348 K:

l  0.549T  945.7220

(15)

Yaws [30] also proposed a correlation for the density at temperatures between 309.58 K and 788.59 K, and it is expressed as 0.28571

T   1   788.59 

l  243.2  0.26098

(16)

where  is expressed in kg·m-3, and the temperature in K. 3.1.4. Heat capacity 

Liquid heat capacity

8

Journal Pre-proof Miltenburg et al. [37] measured the heat capacity of n-eicosane from 310 K to 390 K with an adiabatic calorimeter. The heat capacity data of the liquid were fitted to the following polynomial function form: Cpl  2057.0872  1.3782T  0.006554T 2

(17)

Yaws [30] presented a heat capacity correlation at temperatures between 311 K and 390 K, and it is given by: Cpl  187.2093  12.0994T  25.916  103 T 2  22.7294  106 T 3

(18)

It is worthy to note that the correlation of Yaws [30] agrees very well with the 

of

experimental data from [34]. The difference is around 5%. Solid heat capacity

ro

A heat capacity correlation in solid state was developed from the experimental data of Miltenburg et al. [37]. The R2 of the regression is 99.99% and the maximum relative error is

-p

1.28%. The present developed correlation for 10 < T <309.65 K can be expressed as:

re

Cps  246.64867  15.86901 T  0.07488 T 2  1.78646  104 T 3  9.02412  10 8 T 4

(19)

The heat capacity correlation developed by Yaws [30] is described by Cps  480.623  1.798T  11.007  103 T 2 for 94
lP

(20)

The heat capacity correlations are expressed in J/Kg·K , and the temperature in K.

na

3.1.5. Melting temperature and latent heat of fusion Table 1 summarizes some melting temperatures Tf and latent heat of fusion L f reported in

Jo ur

other studies, which are found to be very close to each other, except the data reported by Jiang et al. [40], which showed greater deviations of 4.2% and 1.2% for L f , and Tf, respectively. This may be attributed to the higher heating rate used (i.e., 10 oC/min). Table 1 Reported melting temperature (Tf) and heat of melting ( L f ) of n-eicosane. References Tf (K) L f (kJ/kg)

[33] 308.84 247.05

[37] 309.65 247.47

[38] 311.60 244.31

[39] 312.32 237.1

[40] 309.75 248.0

[41] 309.50 246.8

[42] 310.0 247.32

[43] 310.0 247.6

4. Numerical method The numerical solution of the melting process is obtained using the software ANSYS 16.0Fluent, which uses finite volume method along with the enthalpy porosity formulation [44]. Since all of the numerical simulations conducted in this study involved a phase change process, and the fluids in these simulations were assumed to be incompressible, the pressurebased model was utilized. The energy and momentum equations were discretized, using a 9

Journal Pre-proof second-order upwind scheme. The under-relaxation factors for the liquid fraction, thermal energy, velocity components and pressure correction are 0.9, 0.95, 0.7, and 0.3, respectively. The convergence criterion is fixed at 10-10 for the energy equation and at 10-6 for the momentum and continuity equations, with number of iterations for every time step is considered 100. The time step used for calculations is fixed at 0.1 s. The typical CPU time for each modeling is about 35 h. 4.1. Mesh independency test The enclosure is meshed with a non-uniform rectangular grid with a very fine spacing near the walls, as shown in Fig.2. An extensive mesh testing procedure was conducted to guarantee

of

a grid-independent solution. Various mesh combinations of 2501, 4941, 8181 and 12221 cells were explored for a=0.5 mm and ΔT = 30°C. The present code was tested for grid

ro

independence by calculating the melting rate in the cavity. The complete melting time for all

-p

grids was numerically computed and shown in Table 2. It can be observed that the relative error between the grid size of 8181 and 12221 is almost the same, therefore the grid size of

Jo ur

na

lP

re

8181ensures a grid-independent solution, and will be adopted in the present study.

Fig. 2. Grid structure.

Table 2 Mesh dependency test. Test 1 2 3 4

Mesh size 2501 4914 8181 12221

Total melting time (s) 534 540 542 543

errors 1.12 0.37 0.18

4.2. Code validation Due to the lack of data for natural convection with PCM in trapezoidal cavities, the experimental data of Gau and Viskanta [15] and the numerical results of Brent et al. [16] for a melting of PCM in a rectangular cavity are compared with the present data. The present results show a reasonable agreement with the published data, which indicates that the present

10

Journal Pre-proof model is adequate for convective flow with phase change materials. The details of the validation are described in [25]. 2 min

6 min

17 min

10 min

0,06

0,04

Y(m) 0,02 Gau and Viskanta [15] Brent et al.[16] Present study 0,04 X(m)

0,06

0,08

ro

0,02

of

0,00 0,00

Fig.3. Comparison of melting front at various times obtained by the present study and previously

-p

reported works.

5. Results and discussion

re

Using the above-described numerical model, simulation runs were carried out to study the melting of n-eicosane in a wavy trapezoidal cavity heated from below, the other walls are

lP

assumed to be adiabatic. The properties for the computed problem are both constant and temperature dependent, and are listed in Table 3. The results are governed by heat transfer in

na

solid and liquid regions, and fluid flow in liquid region. Particular efforts have been focused on the effects of the temperature difference, and wave amplitude which is in the range of 0 ≤ a

Jo ur

≤ 1.5 mm on the transport characteristics during the melting process. The results are presented and discussed in terms of isotherms, streamlines, heat flux distribution, and liquid fraction at various times.

Table 3 Thermophysical properties of n-eicosane. Solid -3

Density( kg·m ) Thermal conductivity (W/m·K) Specific heat (kJ/kg·K)

Liquid   0.549T  945.7220

856 0.551

4

k  0.20681  1.9990  10 T 8

2.6667  10 T Cp  246.64867  15.86901T  0.07488T 4

8

1.78646  10 T  9.02412  10 T

Melting temperature (K) Latent heat (kJ/kg)

3

310 247.6

4

2

2

Cp  2057.0872  1.3782T  0.006554T

2

-

Fig. 4(a) depicts the effect of the temperature difference, and wave amplitude on constant temperature profiles in a trapezoidal cavity for selected six time instants. As expected, during 11

Journal Pre-proof the early stage of melting, the buoyancy forces are unable to conquer the viscous resistance. Hence, the prevailing heat transfer mechanism is conduction as the temperature contours are parallel to the heated wall for a=0 mm, and fit with the wavy surface of the enclosure for all values of wave amplitude (i.e. isothermal stratifications occur). As time elapses to t = 200 s, the strength of the buoyancy-driven flow in the melted PCM is enhanced appreciably in the cavity. As shown from these isotherms, the basic feature of Rayleigh-Bénard problem is depicted accurately, such as the appearance of the thermal plumes. The formation of seven thermal plumes whose appearance are associated with the presence of the isothermal bottom wall are observed. The melted PCM is heated to the highest temperature at different positions

of

on the bottom and then floats up, reaches the phase change interface and splits into different thermal plumes. The melt is cooled as it flows through the phase change interface. This

ro

explains why the phase change interface has a peak at the center of each thermal plume. As

-p

time is more progressed, the natural convection becomes more significant. One can notice that the number of thermal plumes is time dependent due to convection circulation flow; the

re

thermal plumes that start appearing near the cavity sides collide and convolute with the adjacent thermal plumes and move toward the core. This process is repeated along the heated

lP

wall resulting into one main large thermal plume at t =500 s. and t = 650 s. Moreover, it is noted that the amplitude of the wavy surface impacts the thermal plumes along the wavy wall;

na

the number of thermal plumes is decreased with an increase of the wave amplitude up to a = 1 mm (the number of thermal plumes decreases from 7 to 3 at t = 200s). The same number of

Jo ur

thermal plumes is observed for the case of a = 1 mm and a = 1.5 mm. However, stronger natural convection effect is observed for a = 1.5 mm. It is further noticed that when the temperature difference increases from ΔT = 20°C to ΔT = 50°C, as shown in Fig.4(b), The isotherms show a similar trend to those presented in Fig.4a for ΔT = 20°C, but natural convection develops much earlier. As the melts heats up and the wave amplitude increases, the natural convection is intensified and hence the melted PCM increases rapidly. In addition, the thermal plumes are no longer identified and completely destroyed.

12

Journal Pre-proof T(K) 50s

200s

350s

500s

650s

a=0mm

a=0.5mm

a=1mm

a=1.5mm

T(k) 100s

150s

a=0mm

250s

-p

a=0.5mm

b) ΔT=50°C

na

lP

re

a=1mm

a=1.5mm

200s

ro

50s

of

a) ΔT=20°C

Jo ur

Fig.4. Isotherms for different wave amplitudes and time instants for a) ΔT = 20°C, b) ΔT = 50°C.

The effect of the wave amplitude on the streamlines for several times at ΔT = 20°C is illustrated in Fig.5. The displayed streamlines show the development of Rayleigh–Bénard convection cells along the bottom wall. Buoyancy forces induce the melted liquid on the bottom wall to rise vertically toward the upper cool surface. The temperature of the melted liquid cools down in the liquid-solid interface and the flow is driven downward which eventually results in a multicellular flow in the melted region. It is of interest to note that the number of convection cells is time dependent. As time proceeds, the melt height increases, the size of the convection cells increases, and the number of cells decreases. This dependence is in agreement with the result from some previous investigations [27-29]. Due to the low buoyancy during the early stages of the melting process (t <50 s), the convection cells are not fully established especially in the cavity core for a = 0 mm. After the conduction dominating stage, the formation of fourteen convective cells is observed at t = 200 s. As melting progresses, the cells are gradually getting larger and the convection in the enclosure is 13

Journal Pre-proof strengthened; the cells near the enclosure sides eventually collide with an adjacent cell and slowly drift horizontally toward the center, leading to a pair of larger and counter-rotating cells, asymmetric about the horizontal midplane of the enclosure at t = 500 s and t = 650 s. To further reveal the structure of the induced convection cells, the local heat flux at the bottom wall for a = 0 mm at different time instants are presented in Fig. 6. The heat flux distribution is wave-like corresponding to the multiple convection cells. The distribution and number of peaks are related to the velocity distribution in the fluid generated by the ascent of hot fluid plumes under the action of buoyancy forces. At t = 200 s, the evolution has six crests and seven troughs. These latter correspond to the thirteen junctions of the fourteen convection

of

cells in the streamlines shown in Fig. 5. The first crest from left corresponds to the junction of the first and second convection cells. The flow direction of the first circulation is clockwise

ro

and the second circulation is anticlockwise. The liquid layers from the two circulation cells

-p

are cooled after passing the phase change interface and then reach the junction of the two circulation zones at the bottom. This causes a low temperature zone to develop near the

re

junction at the bottom surface of the container. The low temperature zone can be seen in the isotherms in Fig. 4(a). Since the bottom surface of the container is isothermal, a low

lP

temperature near the bottom isothermal surface means a large temperature difference for heat transfer. This results in higher heat flux. At t = 650 s, the heat flux curve presents one peak

na

which corresponds to a pair of asymmetric rolls. The peak has low intensity due to the fact that the melting temperature increases, which reduces the thermal gradients and noticeably

1

a =0mm :

1

Jo ur

reduces the heat transfer at the bottom wall. 50s

200s

350s

500s

a=0.5mm :

a=1mm

a=1.5mm

Fig.5. Streamlines for different wave amplitudes and time instants for ΔT = 20°C.

14

650s

Journal Pre-proof a=0 mm

5

200s

350s

650s

2

Q (kw/m )

4 3 2 1 0 0,0

0,5

1,0

1,5

2,0

2,5

3,0

x/h

of

Fig.6. Heat flux distribution along the bottom surface at different time instants for ΔT = 20°C.

ro

A close inspection of the convective cells in Fig. 5, however, reveals that as the wave amplitude increases, the convective rolls strengthen due the enhancement of the buoyancy

-p

convective effect, and less convective cells are induced. More specifically, there are 14, 10, 6

re

rolls in the cavity, respectively, for a = 0 mm, a =0.5 mm, and a = 1 mm at t = 200 s, nevertheless, the number of convection cells is the same for a = 1 mm and 1.5 mm. In

lP

addition, the same number of rolls is observed for all the wave amplitudes at t > 200s. Note that this result contradicts the earlier result of melting in a rectangular closed enclosure by

na

subjecting the bottom wavy surface to a uniform temperature [29], where they found that the number of convective cells is independent of the wave amplitudes for high temperature

Jo ur

difference. It is further observed that under the continuing action of buoyancy at increasing time, the convective cells reduce for all the four cases, and later, the cells merge together into two main large rolls for a = 0 mm and a = 0.5 mm, but, in addition to the two main cells, a third relatively weak roll at the bottom wall of the cavity starts showing up for a = 1 mm and a = 1.5 mm.

Fig. 7 shows the heat flux at the bottom wall with respect to the wave amplitudes at t = 200 s and t = 650 s. It can be seen that with the increase of melting time, the heat flux is decreased due to the merging of the rolls. Moreover, a lower heat flux is observed for higher wave amplitude. Similar explanation applies to the crests and troughs in the heat flux distributions.

15

Journal Pre-proof 200s 6 a=0,5 mm

a=0 mm

5

a=1,5 mm

a=1 mm

3,0

3

4

2,5

3

2

2

Q (kw/m )

4

Q (kw/m )

2

Q (kw/m )

5

2 2

1 1 0 0,0

1,0

1,5 1,0 0,5

0 0,0 0,5

2,0

0,5

1,5

1,0

2,0

2,5

1,5

2,0

0,0 0,0

x/H 3,0

3,0 1,0

of

x/h a) t=200s

2,5 0,5

1,5

2,0

2,5

3,0

x/h b) t=650s

-p

ro

Fig.7. Heat flux distribution along the bottom surface for ΔT = 20°C at a) t = 200 s, b) t = 650 s.

It is of interest to unveil how the convective cells are formed in the cavity when the

re

temperature difference is raised from ΔT = 20°C to ΔT = 50°C. This is illustrated in Fig. 8 by showing the evolution of flow patterns at different time instants. The result indicates that at

lP

even higher buoyancy, the cavity is dominated by time and wave amplitude dependent rolls. At the early stage of melting, more regular convection cells are induced. As time proceeds,

na

convection cells strengthen, and become curved and are somewhat irregular to a certain degree. As the melting process continue, the convection cells on the sides of the cavity merge

Jo ur

gradually with their neighboring cells and two asymmetric cells are finally formed in the core of the cavity for a = 0 mm and a = 0.5 mm. Moreover, the convection cells for a= 1 mm and a = 1.5 mm continue their evolution without merging with each other, where the RayleighBénard convection becomes more complex, and the flow with multiple cells appears in the entire enclosure. Fig. 9 presents the heat flux distribution for different wave amplitudes at t = 100 s and t = 200 s. It can be observed that raising the temperature difference for a = 1 mm and a = 1.5 mm causes the flow to become chaotic in time, where the rolls become irregular in time and space.

16

Journal Pre-proof 50s

100s

150s

200s

250s

a =0mm a=0.5mm

a=1mm

a=1.5mm

200s

Fig.8. Streamlines for different wave amplitudes and time instants for ΔT = 50°C. 6 5

a=1,5 mm

a=1 mm

14

ro -p

3 2

8

8

2 6

1

4

6 4

re

2

Q (kw/m )

10

4

Q (kw/m )

2

Q (kw/m )

10

12

of

a=0,5 mm

a=0 mm

0,5

1,0

0,5 1,5

x/h

2,0

1,0 2,5

1,5

x/H 3,0

na

0 0,0

0 0,0

lP

2

2

a) t=100s

02,0 0,0

2,5 0,5

3,0 1,0

1,5

2,0

2,5

3,0

x/h b) t=250s

Jo ur

Fig.9. Heat flux distribution along the bottom surface for ΔT = 50°C at a) t = 100 s, b) t = 250 s.

It is important to elucidate the melting rate of n-eicosane in the cavity at increasing buoyancy and wave amplitude. This is illustrated in Fig.10, by showing the liquid fraction evolution for wave amplitude varying between a = 0 mm and a = 1.5 mm for ΔT = 20°C and ΔT = 50°C. As expected, whatever the wave amplitude, with increasing the temperature difference, the convective heat transfer regime advances faster. Consequently, the melting progresses faster. For example, for a = 0 mm, the complete melting time of the PCM for ΔT = 20°C and ΔT = 50°C is 853 s and 307 s, respectively. It can be seen that when the temperature difference is increased about 150%, the complete melting time is 64% increased. Moreover, increasing the wave amplitude ameliorates the heat transfer and leads to an increase in the melting rate of the PCM. However, it is further observed that the wave amplitude has high effect when the temperature difference increases. When the wave amplitude is increased from

17

Journal Pre-proof a = 0 mm to a = 1.5 mm, the complete melting time is increased by 25% for ΔT = 50°C and 18% for ΔT = 20°C. 0

0

T=50 C T= 40 C

0

T= 30 C

0

T= 20 C

1,0

0,6 a=0 mm

0,4

of

a=0,5 mm a=1 mm

0,2

a=1,5 mm

ro

Liquid Fraction

0,8

0,0 150

300

450

600

-p

0

750

900

t(s)

re

Fig.9. Liquid fraction evolution for different wave amplitudes and temperature differences.

6. Conclusion

lP

In the present work, a numerical study of Rayleigh-Bénard convection with phase change in a wavy trapezoidal Cavity is presented. Finite volume method for solving the

na

governing equations and enthalpy-porosity formulation to study the melting process are used. The effects of the temperature difference and wave amplitude on the performance of PCM

Jo ur

melting are investigated. The wave amplitude ranging is from 0 mm to 1.5 mm and the temperature difference up to ΔT = 50°C. With the intention of accurately simulating the melting of n-eicosane inside a wavy trapezoidal cavity, a thorough search for thermo-physical properties is performed. The correlations used have been presented, indicating their origin and accuracy. It was found that very little information has been found for some of them, such as density and thermal conductivity of the solid phase. The effect of the amplitude of the wavy surface on the flow structure and heat transfer characteristics is investigated in details, and the following conclusions can be drawn: 

The number of convection cells is time and wave dependent within the selected temperature difference and wave amplitude ranges.



Increasing the temperature difference for a = 1 mm and a = 1.5 mm causes the flow to become chaotic in time, where the rolls become irregular in time and space.

18

Journal Pre-proof 

Whatever the wave amplitude, with increasing the temperature difference, the convective heat transfer regime advances faster. Consequently, the melting progresses faster; the complete melting time is 64% increased at ΔT = 50°C. Increasing the wave amplitude increases the melting rate but not as higher as increasing the temperature difference; the complete melting time for a = 1.5 m is increased by 25% at ΔT = 50°C. Therefore, it is suggested to change the orientation of

na

lP

re

-p

ro

of

the cavity to significantly increase the melting rate.

Jo ur



19

Journal Pre-proof

References [1] Zietlow DC. Optimization of cooling systems. New York: Momentum Press; 2016. [2] S. Hoseinzadeh, S.A.R. Sahebi, R. Ghasemiasl et al, Experimental analysis to improving thermosyphon (TPCT) thermal efficiency using nanoparticles/based fluids (water), Eur. Phys. J. Plus (2017) 132: 197. [3] A. Nematpour Keshteli, M. Sheikholeslami, Nanoparticle enhanced PCM applications for intensification of thermal performance in building: A review, Journal of Molecular Liquids

of

274 (2019) 516–533.

[4] S. Hoseinzadeh, M.H. Zakeri, A. Shirkhani, and A.J. Chamkha, Analysis of energy

ro

consumption improvements of a zero-energy building in a humid mountainous area, J.

-p

Renewable Sustainable Energy 11, 015103 (2019); doi: 10.1063/1.5046512. [5] Tomasz Kisilewicz , Ma lgorzata Fedorczak-Cisak , Tamas Barkanyi , Active thermal

re

insulation as an element limiting heat loss through external walls, Energy & Buildings (2019), doi: https://doi.org/10.1016/j.enbuild.2019.109541

lP

[6] S. Hoseinzadeh, P. S. Heyns, A. J. Chamkha, A. Shirkhani, Thermal analysis of porous fins enclosure with the comparison of analytical and numerical methods, J Therm Anal

na

Calorim (2019) 138: 727. https://doi.org/10.1007/s10973-019-08203-x. [7] S. Hoseinzadeh , P.S. Heyns , H. Kariman , Numerical investigation of heat transfer of

Jo ur

laminar and turbulent pulsating Al2O3/water nanofluid flow", International Journal of Numerical Methods for Heat & Fluid Flow, https://doi.org/10.1108/HFF-06-2019-0485. [8] A.K Mishra, B.B. Lahiri, J. Philip, Thermal conductivity enhancement in organic phase change material (phenol-water system) upon addition of Al2O3, SiO2 and TiO2 nanoinclusions, doi:10.1016/j.molliq.2018.08.001. [9] Raoux S, Wuttig M. Phase change materials: Science and applications. Berlin: Springer; 2009. [10] Xu B, Li P, Chan C. Application of phase change materials for thermal energy storage in concentrated solar thermal power plants: A review to recent developments. Appl Energy 2015;160:286-307. [11] M.J. Muthoka, Z. Xuelai, Particle size effect on thermophysical properties of nanofluid and nanofluid based phase change materials: A review, doi:10.1016/j.molliq.2018.05.129

20

Journal Pre-proof [12] A. Nematpour Keshteli, M. Sheikholeslam, Nanoparticle enhanced PCM applications for intensification of thermal performance in building: A review, J. Mol. Liq. 274 (2019) 516– 533. [13] M. Sheikholeslami Numerical modeling of nano enhanced PCM solidification in an enclosure with metallic fin, J. Mol. Liq. 259 (2019) 424–438. [14] N.S. Dhaidan, J.M. Khodadadi, Melting and convection of phase change materials in different shape containers: a review, Renew. Sustain. Energy Rev. 43 (2015) 449–477. [15] C. Gau, R. Viskanta, Melting and Solidification of a Pure Metal on a Vertical Wall, J. Heat Transfer 1986;108:174-181.

of

[16] Brent AD, Voller VR, Reid KJ . Enthalpy-porosity technique form odeling convectiondiffusion phase change: application to the melting of a pure metal. NumerHeat transfer 1988;

ro

13: 297–318.

-p

[17] K. El Omari1 T. Kousksou Y. Le Guer, Impact of shape of container on natural convection and melting inside enclosures used for passive cooling of electronic devices,

re

Applied Thermal Engineering, Vol.31, Issue 14-15, October 2011, Pages pp.3022-3035, 2011. [18] Hosseinizadeh SF, Darzi AR, Tan FL, Khodadadi JM. Unconstrained melting inside a

lP

sphere.Int J ThermSci 2013; 63:55–64.

[19] Shmueli H, Ziskind G, Letan R. Melting in a vertical cylindrical tube: numerical

na

investigation and comparison with experiments. Int J Heat Mass Transfer 2010; 53: 4082–91. [20] Darzi AR, Farhadi M, Sedighi K. Numerical study of melting inside concentric and

Jo ur

eccentric horizontal annulus. Appl Math Modell 2012; 36: 4080–6. [21] R.K. Sharma, P. Ganesan, J.N. Sahu, H.S.C. Metselaar, T.M.I. Mahlia, Numerical study for enhancement of solidification of phase change materials using trapezoidal cavity, Powder Technol, 268 ,2014, 38–47.

[22] M. Arıcı , E. Tütüncü, M. Kan, H. Karabay,Melting of nanoparticle-enhanced paraffin wax in a rectangularenclosure with partially active walls, Int.J. Heat and Mass Transfer, 104, 2017, 7–17. [23] L. Zeng, J. Lu, Y. Li, W. Li, S. Liu, J. Zhu, Numerical study of the influences of geometry orientation on phase change material’s melting process, Advances in Mechanical Engineering 2017; 9(10): 1–11. [24] Nabeel S. Dhaidan Melting phase change of n-eicosane inside triangular cavity of two orientations, Journal of renewable and substainable energy, 054101 (2017). [25] F. Iachachen, Z. Haddad, HF. Oztop, E. Abu-Nada, Melting of phase change materials in a trapezoidal cavity: Orientation and nanoparticles effects. J Mol Liq 318:441–450 (2019). 21

Journal Pre-proof [26] S. Hoseinzadeh, R. Ghasemias, D. Havaei, A. J. Chamkha, Numerical investigation of rectangular thermal energy storage units with multiple phase change materials, Journal of Molecular Liquids, 271, 655-660 (2018). [27] Gong ZX, Mujumdar AS. Flow and heat transfer convection-dominated melting in a rectangular cavity heat from below. Int J Heat Mass Transfer 1998;41:2573-80. [28] Fteiti M, Ben Nasrallah S. Numerical study of interaction between the fluid structure and the moving interface during the melting from below in a rectangular closed enclosure. ComputMech 2005;35:161-9. [29] T. Kousksou, M. Mahdaoui, A. Ahmed, A. Ait Msaad, Melting over a wavy surface in a

of

rectangular cavity heated from below, Energy 2014;64:212-219.

[30] Yaws CL. Handbook of thermodynamic and physical properties of chemical compounds.

ro

Knovel; 2003.

-p

[31] Nabil M, Khodadadi JM. Experimental determination of temperature dependent thermal

Heat Mass Transfer 2013;67:301–10.

re

conductivity of solid eicosane-based nanostructure enhanced phase change materials. Int J

[32] P.C. Stryker, E.M. Sparrow, Application of a spherical thermal conductivity cell to solid

lP

n-eicosane paraffin, Int. J. Heat Mass Transfer 1990;33:1781–1793. [33] C. Velez, M. Khayet, J.M. Ortiz de Zarate, Temperature-dependent thermal properties of

na

solid/liquid phase change even-numbered n-alkanes: n-Hexadecane, n-octadecane and neicosane, Applied Energy 2015;143:383–394.

Jo ur

[34] D. V. Hale, M. J. Hoover and M. J. O'Neill, Phase change materials handbook, NASA Contractor Report, NASA'CR-61363 (197l). [35] Research and development study on thermal control by use of fusible materials, Northrop Corporation Interim Report NSL 65-16 to NASA-MSFC on Contract No NAS 8-11163 (February 1965).

[36] J. A. Broadbent, Melting and freezing in a vertical cylinder, M.S. Thesis, Department of Mechanical Engineering, University of Minnesota. Minneapolis. Minnesota (1982). [37] van Miltenburg JC, Oonk HAJ, Metivaud V. Heat capacities and derived thermodynamic functions of n-nonadecane and n-eicosane between 10 K and 390 K. J Chem Eng Data 1999;44:715–20. [38] K. Khimechea, Y. Boumraha, M. Benziane, A. Dahmanib, Solid–liquid equilibria and purity determination for binary n-alkane+naphthalene systems, ThermochimicaActa 2006;444:166–172.

22

Journal Pre-proof [39] Jiang F, Wang X, Wu D. Design and synthesis of magnetic microcapsules based on neicosane core and Fe3O4/SiO2 hybrid shell for dual-functional phase change materials. Appl Energy 2014;134:456–68. [40] Himran S, Suwono A, Mansoori GA. Characterization of alkanes and paraffin waxes for application as phase change energy storage medium. Energy Source 1994;16:117–28. [41] V. Chevalier, D. Petitjean, M. Bouroukba, M. Dirand, Mixtures of numerous different nalkanes: 2. Studies by X-ray diffraction and differential thermal analyses with increasing temperature, Polymer 1999;40:2129, 2137.

at high pressuresFluid Phase Equilib 2004;218:57–68.

of

[42] U. Domanska, P. Morawski, Solid+liquid equilibria of (n-alkane+cyclohexane) mixtures

[43] Gonzalez JA, Zawadzki M, Domanska U. Thermodynamics of mixtures containing

ro

polycyclic aromatic hydrocarbons. J Mol Liq 2008;143:134–40.

Jo ur

na

lP

re

-p

[44] ANSYS FLUENT Theory Guide Release 15.0, ANSYS Inc., 2013.

23

Journal Pre-proof Highlights 

The number of convection cells is time and wave dependent within the selected temperature difference and wave amplitude ranges.



As the temperature difference increases, the melting advances faster whatever the wave amplitude.



The wave amplitude has high effect when the temperature difference increases.



The melting time in the flat surface is about 2.5 times faster than that in the wavy

Jo ur

na

lP

re

-p

ro

of

surface.

24

Journal Pre-proof Conflict of interest

Jo ur

na

lP

re

-p

ro

of

Authors declare that there is no conflict of interest.

25