- Email: [email protected]
Simulation research on PCM freezing process to recover and store the cold energy of cryogenic gas
International Journal of Thermal Sciences 50 (2011) 2220e2227
Contents lists available at ScienceDirect
International Journal of Thermal Sciences journal homepage: www.elsevier.com/locate/ijts
Simulation research on PCM freezing process to recover and store the cold energy of cryogenic gas Hongbo Tan a, Cui Li a, Yanzhong Li a, b, * a b
School of Energy and Power Engineering, Xi’an Jiaotong University, Xi’an 710049, People’s Republic of China State Key Laboratory of Multiphase Flow in Power Engineering, Xi’an Jiaotong University, Xi’an 710049, People’s Republic of China
a r t i c l e i n f o
a b s t r a c t
Article history: Received 18 August 2010 Received in revised form 7 March 2011 Accepted 20 April 2011 Available online 20 July 2011
This paper is concerned with the Phase Change Material (PCM) freezing process that is used for recovery and storage of the cryogenic gas cold energy. Water has been employed as the PCM to investigate the solidification phenomenon outside a tube in which cryogenic nitrogen gas flowed, and numerical study has been conducted on this two-dimensional transient freezing problem by using the Solidification and Melting model in the Computational Fluid Dynamics (CFD) code FLUENT. The calculation results, verified by experimentally measured data, showed that the dimensionless numbers, such as Biot number and Stefan number of PCM, and the Stanton number of the coolant flowing in the tube, have remarkable effects on the characteristics of the two-dimensional freezing problem. Larger Biot number is beneficial to the heat transfer between PCM and the coolant, and can promote higher freezing rate. Higher Stefan number appears to result in larger freezing rate in a fixed axial position. Moreover, the frozen layer slope along the tube is steeper at larger Stanton numbers. Ó 2011 Elsevier Masson SAS. All rights reserved.
Keywords: Cryogenic cooling recovery Two-dimensional freezing Numerical simulation LNG
1. Introduction Cryogenic fuels, such as Liquefied Natural Gas (LNG) and liquid hydrogen (LH2) have been deemed to be promising alternative fuels due to their high energy density. They are stored at extremely low temperatures to maintain the liquid state. In the fuel liquefaction process, a considerable amount of cold energy is transformed from electric energy and then stored within the cryogenic fuels. At consumers end, the cryogenic liquids need to be regasified to about the ambient temperature, and plenty of cryogenic cooling capacity would be released during this vaporizing process. Therefore, it is quite important and imperative to maximumly recover the cold energy of these cryogenic liquids. Up to present, theoretical researches on the cryogenic liquids cold energy and its highefficient utilization have been conducted rather frequently. Many power cycles, based on Rankine and Brayton cycles, were proposed to recover cooling capacity of cryogenic liquids [1e5]. In these studies, cryogenic liquids vaporizers are usually used as lowtemperature thermal sinks for the bottom Rankine cycle, with waste heat or low-grade thermal energy as thermal sources. For
* Corresponding author. State Key Laboratory of Multiphase Flow in Power Engineering, Xi’an Jiaotong University, Xi’an 710049, People’s Republic of China. Tel.: þ86 29 82668725; fax: þ86 29 82668789. E-mail addresses: [email protected] (H. Tan), [email protected] (Y. Li). 1290-0729/$ e see front matter Ó 2011 Elsevier Masson SAS. All rights reserved. doi:10.1016/j.ijthermalsci.2011.04.017
such cases, the performance of cold energy recovery is sensitive to the cold sink stability and its isothermal characteristic. Consequently, a Cold Storage System (CSS), an integrant part of the Cryogenic Fuel Cold Energy Recovery Systems (CFCERS), is an important facility to balance the discrepancy between supply and demand of the cold energy. On one hand, the superfluous cold energy can be conserved by CSS to avoid the unhelpful (even harmful) supercooling of the heat transfer medium. On the other hand, the cold supply fluctuation that is induced by the variation of the cryogenic liquid fuel consumption could be diminished. Therefore, there has been an increasing need for investigating the high-efficient cold storage to balance the discrepancy between supply and demand of cold energy. Due to its isothermal behavior during the charging/discharging process and large energy storage density, the Latent Heat Cold Storage System (LHCSS) is well acknowledged as an appropriate solution to energy conservation and peak cold load shifting for aforementioned CFCERS. The LHCSS techniques applied in energy management systems such as solar energy accumulators [6], energy-saving building envelopes [7], and cold-storage air conditioning systems [8], have been investigated intensively in the last decades. It is believed that the shell-and-tube type heat exchanger is the most promising device as a key apparatus of high efficiency LHCSS. In such an energy storage unit, Phase Change Material (PCM) is filled in the gap between the shell and tubes and heat transfer fluid (HTF) flows in the tubes. During the charging process,
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
Fig. 1. The schematic diagram of freezing process.
the PCM will undergo a liquidesolid phase change in the charging process to store the cold energy and a solideliquid phase change for the energy discharging. As an important characteristic of LHCSS, the PCM freezing process has attracted a great deal of research interest, and the ice formation phenomenon in the shell side was extensively investigated by using both numerical and experimental methods [9e13]. In these studies, the HTF (e.g. ethylene glycol based water solutions) flowing inside the tube was warmed by the latent heat of freezing of water and experienced inappreciable temperature increase because of its large specific heat. It is reported that the HTF temperature increases by about 1 C along a 12.3-m-long tube [13] and the axial variation of tube wall temperature is also negligible. Therefore, the tube wall could be simplified as a constant wall temperature thermal boundary. Moreover, due to relatively high convective heat transfer coefficient inside the tube, the main thermal resistance appears in PCM side for most ice storage systems. In this sense, the effective heat conductivity of the PCM should be enhanced, for example, by adopting finned tubes [11,12], or adding copper particles and graphite matrix impregnated to PCM [10], to promote the PCM freezing process. However, the LHCSS in the CFCERS is quite different from the conventional cold-storage units mentioned above. Firstly, cryogenic gas, characterized by weak convective heat transfer property, is used as the HTF and according the main thermal resistance of the conjugated heat transfer process turns to the HTF side. In this case adding fins to the PCM-side heat transfer surface is no longer an efficient way for the heat transfer enhancement in LHCCS. Secondly, because of the rather small specific heat, cryogenic gas is heated substantially along the tube and this brings about an obvious temperature gradient in the axial direction, which in turn leads to a strongly changed wall temperature distribution in the freezing process of water. Hence the PCM freezing on the surface of the tube may show a distinct pattern from conventional case. To the best of the authors’ knowledge, few studies were reported on the PCM freezing phenomenon in CFCERS in the open literatures, and further work is highly needed to obtain full understanding on this problem. In view of the particularity of the thermal resistance distribution and the varied wall temperature of the conjugated heat transfer with PCM freezing, the present study aims at illuminating the
2221
characteristics of the PCM freezing in the CFCERS and analyzing the effects of the thermal boundary and thermal resistance distribution on the freezing process. By using the Solidification and Melting model, the two-dimensional freezing process of PCM, with different volumetric flow rates and inlet temperatures of gaseous coolant taken into consideration, was simulated in the Computational Fluid Dynamics (CFD) code FLUENT. An experimental investigation of the cryogenic cooling capacity recovery process with LHCSS has been conducted to validate the numerical model. The frozen layer growth performance, the variation of the tube wall temperature, and the cold energy recovered in the CFCERS were obtained and their dependency on influencing factors such as the Biot and Stefan numbers of PCM, and the coolant Stanton number was demonstrated. Considering the physical meanings of the dimensionless numbers mentioned above, their effects on the phase-changing characteristics were also discussed in this paper. 2. Numerical modeling of the two-dimensional freezing problem 2.1. Physical model The physical model of the two-dimensional freezing problem is shown as Figs. 1 and 2. The test tube is made of copper and has a length of 1000 mm, an outer diameter of 20 mm and a thickness of 1 mm. It is placed inside a rectangular tank of 1000 mm 150 mm 200 mm (L W H), which is filled with water, that is, the PCM. Cryogenic gas flows in the copper test tube as the HTF and heat is transferred from PCM to HTF. With a thermal insulation layer the tank wall, the heat loss of which could be neglected comparing with the heat transferred from PCM to HTF, can be regarded as an adiabatic boundary condition. As shown in Fig. 1, the gaseous coolant is substantially heated inside the heat transfer tube due to its small specific heat capacity. Correspondingly, there is a remarkable rise in the coolant temperature along the tube and ice growth rate descends along the tube length obviously. The frozen layer not only increases in radial direction but also propagates in axial direction. If the fact that the initial temperature of water is near its freezing point is taken into account, the natural convection effect in liquid bulk region could be ignored. Therefore, the PCM freezing problem could be simplified to a two-dimensional (in the axial and radial directions) transient freezing problem. 2.2. Numerical modeling In this study, the Solidification and Melting model in FLUENT is adopted to simulate the two-dimensional freezing process. This model employs the enthalpy-porosity technique and treats the phase-changing region as a porous zone with porosity equal to the liquid fraction (b). Appropriate momentum sinks are added to
Fig. 2. The main solidification test unit and thermocouples locations.
2222
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
the momentum equation to account for the pressure drop caused by the presence of solid material. Sinks are also added to the turbulence equations to account for reduced porosity in the solid regions. 2.2.1. Governing equations for the PCM The energy equation of the enthalpy-porosity model is written as
The following hypotheses are made,
v y HÞ ¼ V$ðlVTÞ þ S ðrHÞ þ V$ðr! vs
(1)
H ¼ h þ DH
(2)
y is the fluid velocity, S is a source term, H is the enthalpy of where ! PCM. For the enthalpy-porosity technique, H is computed as the sum of the sensible heat, h, and the latent heat, DH,
where
ZT h ¼ href þ
kin ¼ 3(uinIin)2/2; where Iin represents the turbulence intensity of the velocity at the inlet; 3=4 3=2 3in ¼ ðCm kin =0:07DH Þ, where Cm is an empirical constant specified in the turbulence model (approximately 0.09), DH is the hydraulic diameter of the inlet section.
cp dT
(3)
Axial conduction of cryogenic gas and PCM is negligible; The internal gas flow has reached fully developed state in the test tube; The physical properties of PCM (both liquid phase and solid phase) are uniform and constant; As the initial temperature of PCM is considered to be close to its phase change temperature, the natural convection around the tube can be neglected; Thermal losses and conduction through the outer wall of PCM storage tank have been ignored, i.e., adiabatic condition is applied at the wall of the storage tank.
Tref
and href is the reference enthalpy, Tref is the reference temperature. The latent heat,
DH ¼ bL
(4)
where b is the liquid fraction which is defined as,
b¼
8 < :
T < Tsolid
0 ðT Tsolid Þ Tliquid Tsolid
Tsolid < T < Tliquid
1
T > Tliquid
(5)
Tsolid and Tliquid represent the solidification and melting temperature of PCM respectively. A momentum sink is added to the momentum equation to account for the velocity extinguishment due to the porosity reduction during solidification process. The momentum sink (Smom) resulting from the reduced porosity in the mushy zone takes the following form, 2
Smom
ð1 bÞ Amush ! y ! yp ¼ 3 b þ3
(6)
where 3 is a small number (0.001) to prevent division by zero, Amush y p is the solid velocity due to the is the mushy zone constant, and ! pulling of solid out the domain (also referred to as the pull velocity when the modeling is used to solve the continuous casting problem). In the present study, the solid phase of PCM is not being y p is set equal to zero, that is, ! y p ¼ 0. pulling from the domain, so ! Sinks are added to all of the turbulence equations in the mushy and solidified zones to account for the presence of solid matter. The sink term (Stur) is given as following,
ð1 bÞ Amush f Stur ¼ b3 þ 3 2
(7)
Since the geometrical model is symmetric, only half the domain was taken for simulation. In this study, structured quadrilateral grids were generated and local grid refinement was adopted to account for the complicated heat transfer and the intensive phase change process that exist in the bulk region of PCM near the tube surface. The computational domain of this study was divided into 15,750, 34,000 and 51,600 cells respectively, and the time step was set as 0.001 s. The simulation results showed that the maximum deviations of the internal fluid temperature and ice-layer thickness were just 1.9% and 1.1% for different grids. Evidently grid independent was satisfied for the grid size of 15,750 cells and therefore simulation results demonstrated in the following sections is conducted on this grid. 3. Verification and validation 3.1. Experiment apparatus and procedure An experimental investigation has been conducted on the cryogenic cold energy recovery process in a Latent Heat Cold Storage System (LHCSS) to validate the numerical model. The main solidification test unit is shown in Fig. 2, in which the locations of thermocouples are marked. The main geometrical parameters of the test unit are given in Section 2.1. The test procedure and the error analysis of the experiment system have been introduced in detail in the previous paper [14]. As shown in Fig. 2, the thermocouples were placed both on tube surface and inside the tube, separately, to measure the axial temperature distributions of the tube wall and the cryogenic fluid, and with an equal axial interval of 125 mm starting from the position x ¼ 250 mm. Temperature variations of PCM were also tested at three observed cross-sections, i.e., x ¼ 250 mm, 500 mm and 750 mm. At each cross-section, four thermocouples were symmetrically located with an interval of 5 mm in the PCM tank (see Fig. 2). The inlet and outlet temperatures of cryogenic nitrogen
where f represents the turbulence quantity being solved (k, 3, u). 2.2.2. Boundary conditions and hypothesis The standard ke3 model was adopted to calculate the internal forced convection of gaseous coolant. The numerical parameters for the internal forced convection of cryogenic gas are specified as Uniform distribution of inlet temperature (Tin) and velocity (uin) are assumed;
Table 1 The uncertainties of experimental parameters [14]. Property
Uncertainty
Range
Cryogenic gas flow rate PCM temperature Cryogenic gas temperature Tube wall temperature Ice radius
0.15 Nm3/h 0.02 C 0.02 C 0.02 C 0.35 mm
7.62e11.67 Nm3/h 30 to 10 C 90 to 5 C 55 to 10 C 10e48 mm
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227 Table 2 Experiment conditions summary. Test no.
Gas inlet temperature ( C)
Gas flow rate (Nm3/h)
Initial water temperature ( C)
1 2 3
46.0 61.1 77.9
7.82 11.18 11.67
6.8 6.6 7.1
gas inside tube, i.e., Tf,in and Tf,out, respectively, were measured by T-type thermocouples with a calibrated accuracy of 0.2 C. The flow rate of nitrogen gas was measured by a gas turbine flow meter with a measurement range of 2e20 Nm3/h and the uncertainty was 0.75%. A 40 channels HP type 34970A with the accuracy of 0.04% was used to measure the DCV output of thermocouples. The measured voltage was converted into temperature automatically and transmitted to PC. The 24-bit color images of ice formation at the size of 512 512 pixels were taken by a digital camera and transmitted to PC through a RS-232 interface. The images of ice layer were manipulated with special software so that the liquidesolid interface was clearly indicated and the ice thickness was able to be measured. The uncertainty of the ice thickness is estimated to be within 0.35 mm by employing the basic analysis of Holman and Gajda [15]. Table 1 shows the uncertainties of the measured parameters in this experiment. In the present studied Latent Heat Cold Storage System, water is used as PCM and the test tube is internally cooled by cryogenic nitrogen gas. Ice is formed on the external surface of the test tube to
2223
recover and store of cold energy of cryogenic nitrogen. As summarized in Table 2, several tests were performed under different nitrogen gas inlet conditions. 3.2. Verification and validation of simulation models As summarized in Table 2, three groups of experiments were conducted for CFD model validation. Comparison of the simulated ice-layer thickness and the tube wall temperature with the corresponding experiment data is shown in Fig. 3. Fig. 3(a)e(c) show the simulated time-wise variation of the frozen layer thickness in comparison with the experiment results of Test no. 1, 2 and 3, respectively. As shown in these figures, the freezing began in the position near the tube entrance and the ice growth rate decreased obviously along the tube length for all the cases. This is because the cryogenic gas is heated significantly due to its small specific heat capacity. This brings about an obvious axial temperature gradient along the copper tube, which in turn exerts a varied temperature condition at the tube wall surface (as shown in Fig. 3(d)) in the freezing process. This is probably the fundamental reason for the varied ice growth rates at different axial positions. The calculation results of the time-wise variation of the tube wall temperature are also compared to the experiment results. In view of the similarity of all these results, the case of Test no. 3 is shown in Fig. 3(d) as an example. As the experiment progresses, the tube wall temperature at different axial positions will decrease. That is because the thermal resistance between the coolant and the liquid PCM bulk region increases with the increase of ice layer
Fig. 3. Comparisons of numerical and experimental results.
2224
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
while the coolant temperature keeps almost constant for given test condition. Besides, the wall temperature along the test tube will increase at a certain time which may arise from the remarkable temperature gradient of the gaseous coolant that flows inside the heat transfer tube. As expected, it can be seen that there is a good agreement between experimental data and numerical calculation results. This means the model used in the present simulation is appropriate for the water freezing process and the predicted results are reliable. 4. Results and discussion A validated numerical model was employed to investigate the PCM solidification characteristics under different coolant inlet conditions. The inlet temperatures varied from 40 C to 120 C. Inlet velocities of coolant varied from 4 m/s to 20 m/s. The initial temperature of PCM was set to 7 C. The characteristic parameters, including the ice-layer thickness, the tube wall temperature, and the recovered cold energy, were obtained. The calculation results indicate that these characteristic parameters are influenced by three dimensionless numbers, that is, the coolant Stanton number, the Biot number and the PCM Stefan numbers, which reflect the effect of the internal convective heat transfer of gaseous coolant. The weak convective heat transfer leads to a large thermal resistance at the inner side of the test tube, which, in present study, is the dominate resistance of the conjugate phase change/convection heat transfer process and thus has significant effect on the characteristics of PCM freezing process. 4.1. Frozen layer thickness 4.1.1. The effect of Biot number on the frozen layer growth performance The frozen layer growth performance is acknowledged as the most important characteristics of the LHCSS systems. This section emphatically analyzes the effects of Biot number and Stefan number on the frozen layer growth performance. Biot number, Bi, is defined as Bi ¼ aro/l, where a is the convective heat transfer coefficient of the gaseous coolant (W/m2 K1); ro is the outer radius of the test tube (m); l is the conductivity of the solid phase of PCM (W/m K1). It can be regarded as a measure of the convective heat transfer coefficient of the coolant flow, a, for a fixed conductivity of the solid phase PCM. In other words, a larger Biot number means a relatively high internal convective heat transfer coefficient and enhanced heat transfer from the PCM to the coolant. As a result, a high freezing rate is observed at large Biot number, which is confirmed by the results of Fig. 4. Fig. 4 shows the growth of ice layer at different axial positions (X ¼ 250 mm, 500 mm, 700 mm) with different Biot numbers. As shown in the figure, the freezing begins in the position near the entrance and the ice thickness increasing rate decreases obviously along the tube for all the cases. Taken the case of Bi ¼ 0.217 for an example, at the end of the simulation the ice-layer thicknesses at X ¼ 250 mm, 500 mm, 700 mm are about 16 mm, 12 mm and 8 mm, respectively. This phenomenon is caused by the remarkable axial temperature increase of the gaseous coolant as mentioned in the previous sections. In Fig. 4, ice thickness is more sensitive to Bi at small times than at large times. Specifically, the difference of the ice thickness growth rates between the three different Biot numbers is greater at the beginning of the freezing process. That is because the ice layer is thin at the beginning of freezing process and the thermal resistance of ice is relatively small when compared with the internal thermal resistance of gaseous coolant. At large Biot numbers, the overall thermal resistance will be reduced more significantly in the early stage of the freezing process than that in
Fig. 4. The effect of Biot number on ice-layer growth performance at difference axial position; (a) X ¼ 250 mm; (b) X ¼ 500 mm; (c) X ¼ 750 mm.
the mid or later stage which means the ice layer is thicker and the ice thermal resistance plays an essential role in the overall thermal resistance. Due to the different ice-layer increasing rates at different axial locations, the liquidesolid phase interfaces show an approximately conical shape along the tube length. Fig. 5 shows the axial ice-layer distributions at s ¼ 16 min and s ¼ 120 min for different simulating conditions. It is indicated that the axial ice slope is sensitive to the Stanton number of the coolant. Stanton number measures the ratio
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
2225
Fig. 5. The axial variation of ice-layer thickness at different time with different Stanton number.
of heat transferred into a fluid to the thermal capacity of fluid and is defined as St ¼ a/(urcp), where a is the convective heat transfer coefficient of gaseous coolant (W/m2 K1); u, r and cp present the velocity, density and specific heat of coolant, respectively. A larger Stanton number means more heat is transferred into the coolant and a larger rise in coolant temperature occurs due to its small specific heat. So, larger St yields higher downstream coolant temperatures and smaller PCM-to-coolant temperature differences. Therefore, it can be expected that larger values of St is responsible for the greater decrease of frozen layer thickness along the tube. This expectation is verified by Fig. 5, where the axial slopes of the ice thickness curves are greater for St ¼ 8.58 than those for St ¼ 3.4. 4.1.2. The effect of Stefan number on the frozen layer growth performance The Stefan number, Ste, is defined as the ratio of the solid sensible heat to latent heat in the PCM freezing process. It is given by the formula Ste ¼ cp,ice(T* Tin)/L, where cp,ice is the specific heat of ice (J/kg K1); T* is the solidification temperature (K); Tin is the inlet temperature of coolant (K); L is the latent heat of PCM (J/kg). A high Ste denotes a comparatively low coolant inlet temperature, which means a greater temperature difference exists between the phase-changing region and the gaseous coolant inside the tube. Therefore, higher frozen layer growth rates are observed at larger Ste numbers and the results presented in Fig. 6 come to the same conclusion. As analyzed in the previous section, the enhancement of the freezing performance is more sensitive to Ste at the beginning of freezing process. In the early stage, relatively large temperature difference is found between the phase interface and the coolant when the Stefan number takes a high value, and this will results in a fast growth of ice layer. After the rapid growth period of ice layer, the heat transferred from PCM to the cryogenic gaseous coolant will reduce as the increase of the ice-layer thermal resistance, followed by a modest enhancement effect of Ste on the frozen layer growth rate. Fig. 7 shows the ice-layer thickness along the tube at different times and different Ste numbers for a fixed Biot number (Bi ¼ 0.423). The effect of Stefan number on the ice-layer growth performance is illustrated more intuitively. As shown, the ice layer for Ste ¼ 0.675 is thicker than that for Ste ¼ 0.225 all the time. In addition, the Stanton numbers varies slightly for these cases in Fig. 7 due to the change of physical properties with different Stefan numbers. Therefore in Fig. 7, the axial slopes of ice layers are not obvious and can not be identified.
Fig. 6. The effect of Stefan number on ice-layer growth performance at difference axial positions.
4.2. Tube wall temperature The tube wall temperature distribution reflects the axial variation of cryogenic gaseous coolant temperature and the characteristics of the conjugated heat transfer during the freezing process partially. Fig. 8 shows the effects of Biot number on the time-wise variation of tube wall temperature. As shown in this figure, at a fixed location and time, the tube wall temperature (e.g. X ¼ 250 mm) for Bi ¼ 0.476 is lower than that for Bi ¼ 0.217. That is because, the overall thermal resistance would be significantly decreased for large values of Biot number due to the
2226
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
Fig. 9. The effects of Stefan number on the tube wall temperature. Fig. 7. The axial variation of ice-layer thickness at different time with different Stefan number.
enhancement of the heat transfer inside the tube. Hence for large Bi, the tube wall temperature is much closer to the cryogenic gaseous coolant and the wall-to-coolant temperature difference is smaller when compared with small Bi number. Fig. 9 illustrates the effect of the Stefan number on the timewise variation of the tube wall temperature. As can be seen from the figure, the tube temperature curves for Ste ¼ 0.45 lie below those for Ste ¼ 0.225 and the tube temperature decrease in a fixed heat transfer tube section (e.g. from X ¼ 250 mm to X ¼ 500 mm) is greater for Ste ¼ 0.45 than that for Ste ¼ 0.225. Moreover, higher heat transfer rates appears for Ste ¼ 0.45 due to the lower gaseous coolant inlet temperature. 4.3. Cold energy recovered With respect to thermal storage applications, the amount of cold energy stored in the PCM freezing process is of great importance. By neglecting the heat losses from the surrounding into the test unit, the cold energy recovered by the LHCSS will be equal to the total enthalpy increase of the cryogenic gas, which is given by
_ p Tf ;out Tf ;in Q ¼ mc
(8)
_ is the mass flow rate of coolant (kg/s); cp is the specific where m heat of coolant (J/kg K1); Tf,in and Tf,out respectively represent the inlet and outlet temperatures of the coolant flow (K).
Fig. 8. The effects of Biot number on the tube wall temperature.
Fig. 10 shows the effects of Biot number on the cold per unit tube length for the same Stefan number. At a fixed Stefan number, the stored cold is seen to be quite sensitive to the Biot number and increases remarkably as the Biot number increases. Due to the growth of frozen layer, the ice thermal resistance increases as the freezing process progresses and the overall thermal resistance increases correspondingly. Hence, it might be expected that the heat transfer rate should decrease with the thickening of the ice layer. However, as shown in the figure, the decreasing trend of the heat transfer rate is moderate, and this can be explained by the fact that the overall thermal resistance is not so sensitive to the increasing frozen layer thermal resistance. Actually in this study, the ice-layer thermal resistance just accounts for less than ten percent of the overall resistance at the end of the freezing process [14]. Therefore, the cold energy storage rate can be simplified to a constant as marked in Figs. 10 and 11. The effects of Stefan number on the recovered cold energy per unit tube length for the same Biot number are illustrated in Fig. 11. The heat flux for Ste ¼ 0.675 is 1.25 times larger than that for Ste ¼ 0.217. It can be explained by the fact that for a large Stefan number, larger PCM-to-coolant temperature difference results from the lower cryogenic gaseous coolant inlet temperature. Correspondingly, the logarithmic mean temperature difference increases as the Stefan number increases at fixed Bi. Therefore, the heat transferred from PCM to coolant will increase accordingly.
Fig. 10. The effects of Biot number on the cold recovered for Ste ¼ 0.562.
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
2227
2) The Stanton number reflects the convective heat transfer capabilities of the gaseous coolant. Larger Stanton number means that more heat is transferred into the coolant and a larger rise occurs in the coolant temperature. This is responsible for the greater decrease of the frozen layer thickness along the tube length. 3) In LHCSS the stored cold energy per unit length can be simplified to a linear function of time. The cold-storage rate is seen to be quite sensitive to the Biot number and the Stefan number. Acknowledgments This work was supported by the High-Tech Research and Development Program under the contract No. 2007AA05Z216 in China. Fig. 11. The effects of Stefan number on the cold recovered for Bi ¼ 0.432.
5. Conclusions This paper presented a numerical analysis on the PCM freezing process used to recover and store cryogenic cold energy of gaseous coolant. Based on the Solidification and Melting model, the twodimensional phase change problem was simulated in the Computational Fluid Dynamics (CFD) code FLUENT. An experimental investigation was conducted to validate the CFD model and the simulation results showed a good agreement with experimental results. When a flow of cryogenic gas passes through a horizontal copper tube submerged in liquid PCM, solidification of PCM will occur on the outer surface of the tube. The major characteristics of PCM freezing in this investigation are the weak convective heat transfer inside the tube, the large thermal resistance of the gaseous coolant that dominates the overall thermal resistance, and the large temperature increase experienced by the cryogenic gas. Effects of Biot number, Stanton number, and Stefan number on the PCM freezing characteristics have been analyzed in detail. The main conclusions can be summarized as follows: 1) The frozen layer growth performance is sensitive to Biot number and Stefan number. A larger Biot number results in a relatively higher ice layer increasing rate at a fixed axial position due to an increase in the internal convective heat transfer coefficient. Higher values of Stefan number leads to lower coolant inlet temperature and greater temperature difference between the phase-changing region and the gaseous coolant inside the tube, and consequently fast growth rates of frozen layer.
References [1] Na Zhang, Noam Lior, A novel Brayton cycle with the integration of liquid hydrogen cryogenic exergy utilization, Int. J. Hydrogen Energ. 33 (1) (2008) 214e224. [2] Qiang Wang, Yanzhong Li, Jiang Wang, Analysis of power cycle based on cold energy of liquefied natural gas and low-grade heat source, Appl. Therm. Eng. 24 (4) (2004) 539e548. [3] Qiang Wang, Yanzhong Li, Chen Xi, Exergy analysis of liquefied natural gas cold energy recovering cycles, Int. J. Energ. Res. 29 (2005) 65e78. [4] Na Zhang, Noam Lior, A novel near-zero CO2 emission thermal cycle with LNG cryogenic exergy utilization, Energy 31 (10e11) (2006) 1666e1679. [5] Jan Szargut, Ireneusz Szczygiel, Utilization of the cryogenic exergy of liquid natural gas (LNG) for the production of electricity, Energy 34 (7) (2009) 827e837. [6] Belén Zalba, José M. Marín, Luisa F. Cabeza, et al., Review on thermal energy storage with phase change: materials, heat transfer analysis and applications, Appl. Therm. Eng. 23 (3) (2003) 251e283. [7] A.M. Khudhair, M.M. Farid, A review on energy conservation in building applications with thermal storage by latent heat using phase change materials, Energ. Convers. Manage. 45 (2) (2004) 263e275. [8] Akio Saito, Recent advances in research on cold thermal energy storage, Int. J. Refrigeration 25 (2) (2002) 177e189. [9] Babak K. Soltan, Morteza M. Ardehali, Numerical simulation of water solidification phenomenon for ice-on-coil thermal energy storage application, Energ. Convers. Manage. 44 (2003) 85e92. [10] L.F. Cabeza, H. Mehling, S. Hiebler, F. Ziegler, Heat transfer enhancement in water when used as PCM in thermal energy storage, Appl. Therm. Eng. 22 (2002) 1141e1151. _ [11] Aytunc Erek, Zafer Ilken, M. Ali Acar, Experimental and numerical investigation of thermal energy storage with a finned tube, Int. J. Energ. Res. 29 (2005) 283e301. [12] N. Kayansayan, M. Ali Acar, Ice formation around a finned-tube heat exchanger for cold thermal energy storage, Int. J. Therm. Sci. 45 (2006) 405e418. [13] Badr A. Habeebullah, An experimental study on ice formation around horizontal long tubes, Int. J. Refrigeration 20 (2007) 789e797. [14] Hongbo Tan, Yanzhong Li, Hanfei Tuo, Man Zhou, Baocong Tian, Experimental study on liquid/solid phase change for cold energy storage of LNG refrigerated vehicle, Energy 35 (5) (2010) 1927e1935. [15] J.P. Holman, W.J. Gajda, Experimental Methods for Engineers. McGraw-Hill, New York, 1994.
Contents lists available at ScienceDirect
International Journal of Thermal Sciences journal homepage: www.elsevier.com/locate/ijts
Simulation research on PCM freezing process to recover and store the cold energy of cryogenic gas Hongbo Tan a, Cui Li a, Yanzhong Li a, b, * a b
School of Energy and Power Engineering, Xi’an Jiaotong University, Xi’an 710049, People’s Republic of China State Key Laboratory of Multiphase Flow in Power Engineering, Xi’an Jiaotong University, Xi’an 710049, People’s Republic of China
a r t i c l e i n f o
a b s t r a c t
Article history: Received 18 August 2010 Received in revised form 7 March 2011 Accepted 20 April 2011 Available online 20 July 2011
This paper is concerned with the Phase Change Material (PCM) freezing process that is used for recovery and storage of the cryogenic gas cold energy. Water has been employed as the PCM to investigate the solidification phenomenon outside a tube in which cryogenic nitrogen gas flowed, and numerical study has been conducted on this two-dimensional transient freezing problem by using the Solidification and Melting model in the Computational Fluid Dynamics (CFD) code FLUENT. The calculation results, verified by experimentally measured data, showed that the dimensionless numbers, such as Biot number and Stefan number of PCM, and the Stanton number of the coolant flowing in the tube, have remarkable effects on the characteristics of the two-dimensional freezing problem. Larger Biot number is beneficial to the heat transfer between PCM and the coolant, and can promote higher freezing rate. Higher Stefan number appears to result in larger freezing rate in a fixed axial position. Moreover, the frozen layer slope along the tube is steeper at larger Stanton numbers. Ó 2011 Elsevier Masson SAS. All rights reserved.
Keywords: Cryogenic cooling recovery Two-dimensional freezing Numerical simulation LNG
1. Introduction Cryogenic fuels, such as Liquefied Natural Gas (LNG) and liquid hydrogen (LH2) have been deemed to be promising alternative fuels due to their high energy density. They are stored at extremely low temperatures to maintain the liquid state. In the fuel liquefaction process, a considerable amount of cold energy is transformed from electric energy and then stored within the cryogenic fuels. At consumers end, the cryogenic liquids need to be regasified to about the ambient temperature, and plenty of cryogenic cooling capacity would be released during this vaporizing process. Therefore, it is quite important and imperative to maximumly recover the cold energy of these cryogenic liquids. Up to present, theoretical researches on the cryogenic liquids cold energy and its highefficient utilization have been conducted rather frequently. Many power cycles, based on Rankine and Brayton cycles, were proposed to recover cooling capacity of cryogenic liquids [1e5]. In these studies, cryogenic liquids vaporizers are usually used as lowtemperature thermal sinks for the bottom Rankine cycle, with waste heat or low-grade thermal energy as thermal sources. For
* Corresponding author. State Key Laboratory of Multiphase Flow in Power Engineering, Xi’an Jiaotong University, Xi’an 710049, People’s Republic of China. Tel.: þ86 29 82668725; fax: þ86 29 82668789. E-mail addresses: [email protected] (H. Tan), [email protected] (Y. Li). 1290-0729/$ e see front matter Ó 2011 Elsevier Masson SAS. All rights reserved. doi:10.1016/j.ijthermalsci.2011.04.017
such cases, the performance of cold energy recovery is sensitive to the cold sink stability and its isothermal characteristic. Consequently, a Cold Storage System (CSS), an integrant part of the Cryogenic Fuel Cold Energy Recovery Systems (CFCERS), is an important facility to balance the discrepancy between supply and demand of the cold energy. On one hand, the superfluous cold energy can be conserved by CSS to avoid the unhelpful (even harmful) supercooling of the heat transfer medium. On the other hand, the cold supply fluctuation that is induced by the variation of the cryogenic liquid fuel consumption could be diminished. Therefore, there has been an increasing need for investigating the high-efficient cold storage to balance the discrepancy between supply and demand of cold energy. Due to its isothermal behavior during the charging/discharging process and large energy storage density, the Latent Heat Cold Storage System (LHCSS) is well acknowledged as an appropriate solution to energy conservation and peak cold load shifting for aforementioned CFCERS. The LHCSS techniques applied in energy management systems such as solar energy accumulators [6], energy-saving building envelopes [7], and cold-storage air conditioning systems [8], have been investigated intensively in the last decades. It is believed that the shell-and-tube type heat exchanger is the most promising device as a key apparatus of high efficiency LHCSS. In such an energy storage unit, Phase Change Material (PCM) is filled in the gap between the shell and tubes and heat transfer fluid (HTF) flows in the tubes. During the charging process,
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
Fig. 1. The schematic diagram of freezing process.
the PCM will undergo a liquidesolid phase change in the charging process to store the cold energy and a solideliquid phase change for the energy discharging. As an important characteristic of LHCSS, the PCM freezing process has attracted a great deal of research interest, and the ice formation phenomenon in the shell side was extensively investigated by using both numerical and experimental methods [9e13]. In these studies, the HTF (e.g. ethylene glycol based water solutions) flowing inside the tube was warmed by the latent heat of freezing of water and experienced inappreciable temperature increase because of its large specific heat. It is reported that the HTF temperature increases by about 1 C along a 12.3-m-long tube [13] and the axial variation of tube wall temperature is also negligible. Therefore, the tube wall could be simplified as a constant wall temperature thermal boundary. Moreover, due to relatively high convective heat transfer coefficient inside the tube, the main thermal resistance appears in PCM side for most ice storage systems. In this sense, the effective heat conductivity of the PCM should be enhanced, for example, by adopting finned tubes [11,12], or adding copper particles and graphite matrix impregnated to PCM [10], to promote the PCM freezing process. However, the LHCSS in the CFCERS is quite different from the conventional cold-storage units mentioned above. Firstly, cryogenic gas, characterized by weak convective heat transfer property, is used as the HTF and according the main thermal resistance of the conjugated heat transfer process turns to the HTF side. In this case adding fins to the PCM-side heat transfer surface is no longer an efficient way for the heat transfer enhancement in LHCCS. Secondly, because of the rather small specific heat, cryogenic gas is heated substantially along the tube and this brings about an obvious temperature gradient in the axial direction, which in turn leads to a strongly changed wall temperature distribution in the freezing process of water. Hence the PCM freezing on the surface of the tube may show a distinct pattern from conventional case. To the best of the authors’ knowledge, few studies were reported on the PCM freezing phenomenon in CFCERS in the open literatures, and further work is highly needed to obtain full understanding on this problem. In view of the particularity of the thermal resistance distribution and the varied wall temperature of the conjugated heat transfer with PCM freezing, the present study aims at illuminating the
2221
characteristics of the PCM freezing in the CFCERS and analyzing the effects of the thermal boundary and thermal resistance distribution on the freezing process. By using the Solidification and Melting model, the two-dimensional freezing process of PCM, with different volumetric flow rates and inlet temperatures of gaseous coolant taken into consideration, was simulated in the Computational Fluid Dynamics (CFD) code FLUENT. An experimental investigation of the cryogenic cooling capacity recovery process with LHCSS has been conducted to validate the numerical model. The frozen layer growth performance, the variation of the tube wall temperature, and the cold energy recovered in the CFCERS were obtained and their dependency on influencing factors such as the Biot and Stefan numbers of PCM, and the coolant Stanton number was demonstrated. Considering the physical meanings of the dimensionless numbers mentioned above, their effects on the phase-changing characteristics were also discussed in this paper. 2. Numerical modeling of the two-dimensional freezing problem 2.1. Physical model The physical model of the two-dimensional freezing problem is shown as Figs. 1 and 2. The test tube is made of copper and has a length of 1000 mm, an outer diameter of 20 mm and a thickness of 1 mm. It is placed inside a rectangular tank of 1000 mm 150 mm 200 mm (L W H), which is filled with water, that is, the PCM. Cryogenic gas flows in the copper test tube as the HTF and heat is transferred from PCM to HTF. With a thermal insulation layer the tank wall, the heat loss of which could be neglected comparing with the heat transferred from PCM to HTF, can be regarded as an adiabatic boundary condition. As shown in Fig. 1, the gaseous coolant is substantially heated inside the heat transfer tube due to its small specific heat capacity. Correspondingly, there is a remarkable rise in the coolant temperature along the tube and ice growth rate descends along the tube length obviously. The frozen layer not only increases in radial direction but also propagates in axial direction. If the fact that the initial temperature of water is near its freezing point is taken into account, the natural convection effect in liquid bulk region could be ignored. Therefore, the PCM freezing problem could be simplified to a two-dimensional (in the axial and radial directions) transient freezing problem. 2.2. Numerical modeling In this study, the Solidification and Melting model in FLUENT is adopted to simulate the two-dimensional freezing process. This model employs the enthalpy-porosity technique and treats the phase-changing region as a porous zone with porosity equal to the liquid fraction (b). Appropriate momentum sinks are added to
Fig. 2. The main solidification test unit and thermocouples locations.
2222
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
the momentum equation to account for the pressure drop caused by the presence of solid material. Sinks are also added to the turbulence equations to account for reduced porosity in the solid regions. 2.2.1. Governing equations for the PCM The energy equation of the enthalpy-porosity model is written as
The following hypotheses are made,
v y HÞ ¼ V$ðlVTÞ þ S ðrHÞ þ V$ðr! vs
(1)
H ¼ h þ DH
(2)
y is the fluid velocity, S is a source term, H is the enthalpy of where ! PCM. For the enthalpy-porosity technique, H is computed as the sum of the sensible heat, h, and the latent heat, DH,
where
ZT h ¼ href þ
kin ¼ 3(uinIin)2/2; where Iin represents the turbulence intensity of the velocity at the inlet; 3=4 3=2 3in ¼ ðCm kin =0:07DH Þ, where Cm is an empirical constant specified in the turbulence model (approximately 0.09), DH is the hydraulic diameter of the inlet section.
cp dT
(3)
Axial conduction of cryogenic gas and PCM is negligible; The internal gas flow has reached fully developed state in the test tube; The physical properties of PCM (both liquid phase and solid phase) are uniform and constant; As the initial temperature of PCM is considered to be close to its phase change temperature, the natural convection around the tube can be neglected; Thermal losses and conduction through the outer wall of PCM storage tank have been ignored, i.e., adiabatic condition is applied at the wall of the storage tank.
Tref
and href is the reference enthalpy, Tref is the reference temperature. The latent heat,
DH ¼ bL
(4)
where b is the liquid fraction which is defined as,
b¼
8 < :
T < Tsolid
0 ðT Tsolid Þ Tliquid Tsolid
Tsolid < T < Tliquid
1
T > Tliquid
(5)
Tsolid and Tliquid represent the solidification and melting temperature of PCM respectively. A momentum sink is added to the momentum equation to account for the velocity extinguishment due to the porosity reduction during solidification process. The momentum sink (Smom) resulting from the reduced porosity in the mushy zone takes the following form, 2
Smom
ð1 bÞ Amush ! y ! yp ¼ 3 b þ3
(6)
where 3 is a small number (0.001) to prevent division by zero, Amush y p is the solid velocity due to the is the mushy zone constant, and ! pulling of solid out the domain (also referred to as the pull velocity when the modeling is used to solve the continuous casting problem). In the present study, the solid phase of PCM is not being y p is set equal to zero, that is, ! y p ¼ 0. pulling from the domain, so ! Sinks are added to all of the turbulence equations in the mushy and solidified zones to account for the presence of solid matter. The sink term (Stur) is given as following,
ð1 bÞ Amush f Stur ¼ b3 þ 3 2
(7)
Since the geometrical model is symmetric, only half the domain was taken for simulation. In this study, structured quadrilateral grids were generated and local grid refinement was adopted to account for the complicated heat transfer and the intensive phase change process that exist in the bulk region of PCM near the tube surface. The computational domain of this study was divided into 15,750, 34,000 and 51,600 cells respectively, and the time step was set as 0.001 s. The simulation results showed that the maximum deviations of the internal fluid temperature and ice-layer thickness were just 1.9% and 1.1% for different grids. Evidently grid independent was satisfied for the grid size of 15,750 cells and therefore simulation results demonstrated in the following sections is conducted on this grid. 3. Verification and validation 3.1. Experiment apparatus and procedure An experimental investigation has been conducted on the cryogenic cold energy recovery process in a Latent Heat Cold Storage System (LHCSS) to validate the numerical model. The main solidification test unit is shown in Fig. 2, in which the locations of thermocouples are marked. The main geometrical parameters of the test unit are given in Section 2.1. The test procedure and the error analysis of the experiment system have been introduced in detail in the previous paper [14]. As shown in Fig. 2, the thermocouples were placed both on tube surface and inside the tube, separately, to measure the axial temperature distributions of the tube wall and the cryogenic fluid, and with an equal axial interval of 125 mm starting from the position x ¼ 250 mm. Temperature variations of PCM were also tested at three observed cross-sections, i.e., x ¼ 250 mm, 500 mm and 750 mm. At each cross-section, four thermocouples were symmetrically located with an interval of 5 mm in the PCM tank (see Fig. 2). The inlet and outlet temperatures of cryogenic nitrogen
where f represents the turbulence quantity being solved (k, 3, u). 2.2.2. Boundary conditions and hypothesis The standard ke3 model was adopted to calculate the internal forced convection of gaseous coolant. The numerical parameters for the internal forced convection of cryogenic gas are specified as Uniform distribution of inlet temperature (Tin) and velocity (uin) are assumed;
Table 1 The uncertainties of experimental parameters [14]. Property
Uncertainty
Range
Cryogenic gas flow rate PCM temperature Cryogenic gas temperature Tube wall temperature Ice radius
0.15 Nm3/h 0.02 C 0.02 C 0.02 C 0.35 mm
7.62e11.67 Nm3/h 30 to 10 C 90 to 5 C 55 to 10 C 10e48 mm
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227 Table 2 Experiment conditions summary. Test no.
Gas inlet temperature ( C)
Gas flow rate (Nm3/h)
Initial water temperature ( C)
1 2 3
46.0 61.1 77.9
7.82 11.18 11.67
6.8 6.6 7.1
gas inside tube, i.e., Tf,in and Tf,out, respectively, were measured by T-type thermocouples with a calibrated accuracy of 0.2 C. The flow rate of nitrogen gas was measured by a gas turbine flow meter with a measurement range of 2e20 Nm3/h and the uncertainty was 0.75%. A 40 channels HP type 34970A with the accuracy of 0.04% was used to measure the DCV output of thermocouples. The measured voltage was converted into temperature automatically and transmitted to PC. The 24-bit color images of ice formation at the size of 512 512 pixels were taken by a digital camera and transmitted to PC through a RS-232 interface. The images of ice layer were manipulated with special software so that the liquidesolid interface was clearly indicated and the ice thickness was able to be measured. The uncertainty of the ice thickness is estimated to be within 0.35 mm by employing the basic analysis of Holman and Gajda [15]. Table 1 shows the uncertainties of the measured parameters in this experiment. In the present studied Latent Heat Cold Storage System, water is used as PCM and the test tube is internally cooled by cryogenic nitrogen gas. Ice is formed on the external surface of the test tube to
2223
recover and store of cold energy of cryogenic nitrogen. As summarized in Table 2, several tests were performed under different nitrogen gas inlet conditions. 3.2. Verification and validation of simulation models As summarized in Table 2, three groups of experiments were conducted for CFD model validation. Comparison of the simulated ice-layer thickness and the tube wall temperature with the corresponding experiment data is shown in Fig. 3. Fig. 3(a)e(c) show the simulated time-wise variation of the frozen layer thickness in comparison with the experiment results of Test no. 1, 2 and 3, respectively. As shown in these figures, the freezing began in the position near the tube entrance and the ice growth rate decreased obviously along the tube length for all the cases. This is because the cryogenic gas is heated significantly due to its small specific heat capacity. This brings about an obvious axial temperature gradient along the copper tube, which in turn exerts a varied temperature condition at the tube wall surface (as shown in Fig. 3(d)) in the freezing process. This is probably the fundamental reason for the varied ice growth rates at different axial positions. The calculation results of the time-wise variation of the tube wall temperature are also compared to the experiment results. In view of the similarity of all these results, the case of Test no. 3 is shown in Fig. 3(d) as an example. As the experiment progresses, the tube wall temperature at different axial positions will decrease. That is because the thermal resistance between the coolant and the liquid PCM bulk region increases with the increase of ice layer
Fig. 3. Comparisons of numerical and experimental results.
2224
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
while the coolant temperature keeps almost constant for given test condition. Besides, the wall temperature along the test tube will increase at a certain time which may arise from the remarkable temperature gradient of the gaseous coolant that flows inside the heat transfer tube. As expected, it can be seen that there is a good agreement between experimental data and numerical calculation results. This means the model used in the present simulation is appropriate for the water freezing process and the predicted results are reliable. 4. Results and discussion A validated numerical model was employed to investigate the PCM solidification characteristics under different coolant inlet conditions. The inlet temperatures varied from 40 C to 120 C. Inlet velocities of coolant varied from 4 m/s to 20 m/s. The initial temperature of PCM was set to 7 C. The characteristic parameters, including the ice-layer thickness, the tube wall temperature, and the recovered cold energy, were obtained. The calculation results indicate that these characteristic parameters are influenced by three dimensionless numbers, that is, the coolant Stanton number, the Biot number and the PCM Stefan numbers, which reflect the effect of the internal convective heat transfer of gaseous coolant. The weak convective heat transfer leads to a large thermal resistance at the inner side of the test tube, which, in present study, is the dominate resistance of the conjugate phase change/convection heat transfer process and thus has significant effect on the characteristics of PCM freezing process. 4.1. Frozen layer thickness 4.1.1. The effect of Biot number on the frozen layer growth performance The frozen layer growth performance is acknowledged as the most important characteristics of the LHCSS systems. This section emphatically analyzes the effects of Biot number and Stefan number on the frozen layer growth performance. Biot number, Bi, is defined as Bi ¼ aro/l, where a is the convective heat transfer coefficient of the gaseous coolant (W/m2 K1); ro is the outer radius of the test tube (m); l is the conductivity of the solid phase of PCM (W/m K1). It can be regarded as a measure of the convective heat transfer coefficient of the coolant flow, a, for a fixed conductivity of the solid phase PCM. In other words, a larger Biot number means a relatively high internal convective heat transfer coefficient and enhanced heat transfer from the PCM to the coolant. As a result, a high freezing rate is observed at large Biot number, which is confirmed by the results of Fig. 4. Fig. 4 shows the growth of ice layer at different axial positions (X ¼ 250 mm, 500 mm, 700 mm) with different Biot numbers. As shown in the figure, the freezing begins in the position near the entrance and the ice thickness increasing rate decreases obviously along the tube for all the cases. Taken the case of Bi ¼ 0.217 for an example, at the end of the simulation the ice-layer thicknesses at X ¼ 250 mm, 500 mm, 700 mm are about 16 mm, 12 mm and 8 mm, respectively. This phenomenon is caused by the remarkable axial temperature increase of the gaseous coolant as mentioned in the previous sections. In Fig. 4, ice thickness is more sensitive to Bi at small times than at large times. Specifically, the difference of the ice thickness growth rates between the three different Biot numbers is greater at the beginning of the freezing process. That is because the ice layer is thin at the beginning of freezing process and the thermal resistance of ice is relatively small when compared with the internal thermal resistance of gaseous coolant. At large Biot numbers, the overall thermal resistance will be reduced more significantly in the early stage of the freezing process than that in
Fig. 4. The effect of Biot number on ice-layer growth performance at difference axial position; (a) X ¼ 250 mm; (b) X ¼ 500 mm; (c) X ¼ 750 mm.
the mid or later stage which means the ice layer is thicker and the ice thermal resistance plays an essential role in the overall thermal resistance. Due to the different ice-layer increasing rates at different axial locations, the liquidesolid phase interfaces show an approximately conical shape along the tube length. Fig. 5 shows the axial ice-layer distributions at s ¼ 16 min and s ¼ 120 min for different simulating conditions. It is indicated that the axial ice slope is sensitive to the Stanton number of the coolant. Stanton number measures the ratio
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
2225
Fig. 5. The axial variation of ice-layer thickness at different time with different Stanton number.
of heat transferred into a fluid to the thermal capacity of fluid and is defined as St ¼ a/(urcp), where a is the convective heat transfer coefficient of gaseous coolant (W/m2 K1); u, r and cp present the velocity, density and specific heat of coolant, respectively. A larger Stanton number means more heat is transferred into the coolant and a larger rise in coolant temperature occurs due to its small specific heat. So, larger St yields higher downstream coolant temperatures and smaller PCM-to-coolant temperature differences. Therefore, it can be expected that larger values of St is responsible for the greater decrease of frozen layer thickness along the tube. This expectation is verified by Fig. 5, where the axial slopes of the ice thickness curves are greater for St ¼ 8.58 than those for St ¼ 3.4. 4.1.2. The effect of Stefan number on the frozen layer growth performance The Stefan number, Ste, is defined as the ratio of the solid sensible heat to latent heat in the PCM freezing process. It is given by the formula Ste ¼ cp,ice(T* Tin)/L, where cp,ice is the specific heat of ice (J/kg K1); T* is the solidification temperature (K); Tin is the inlet temperature of coolant (K); L is the latent heat of PCM (J/kg). A high Ste denotes a comparatively low coolant inlet temperature, which means a greater temperature difference exists between the phase-changing region and the gaseous coolant inside the tube. Therefore, higher frozen layer growth rates are observed at larger Ste numbers and the results presented in Fig. 6 come to the same conclusion. As analyzed in the previous section, the enhancement of the freezing performance is more sensitive to Ste at the beginning of freezing process. In the early stage, relatively large temperature difference is found between the phase interface and the coolant when the Stefan number takes a high value, and this will results in a fast growth of ice layer. After the rapid growth period of ice layer, the heat transferred from PCM to the cryogenic gaseous coolant will reduce as the increase of the ice-layer thermal resistance, followed by a modest enhancement effect of Ste on the frozen layer growth rate. Fig. 7 shows the ice-layer thickness along the tube at different times and different Ste numbers for a fixed Biot number (Bi ¼ 0.423). The effect of Stefan number on the ice-layer growth performance is illustrated more intuitively. As shown, the ice layer for Ste ¼ 0.675 is thicker than that for Ste ¼ 0.225 all the time. In addition, the Stanton numbers varies slightly for these cases in Fig. 7 due to the change of physical properties with different Stefan numbers. Therefore in Fig. 7, the axial slopes of ice layers are not obvious and can not be identified.
Fig. 6. The effect of Stefan number on ice-layer growth performance at difference axial positions.
4.2. Tube wall temperature The tube wall temperature distribution reflects the axial variation of cryogenic gaseous coolant temperature and the characteristics of the conjugated heat transfer during the freezing process partially. Fig. 8 shows the effects of Biot number on the time-wise variation of tube wall temperature. As shown in this figure, at a fixed location and time, the tube wall temperature (e.g. X ¼ 250 mm) for Bi ¼ 0.476 is lower than that for Bi ¼ 0.217. That is because, the overall thermal resistance would be significantly decreased for large values of Biot number due to the
2226
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
Fig. 9. The effects of Stefan number on the tube wall temperature. Fig. 7. The axial variation of ice-layer thickness at different time with different Stefan number.
enhancement of the heat transfer inside the tube. Hence for large Bi, the tube wall temperature is much closer to the cryogenic gaseous coolant and the wall-to-coolant temperature difference is smaller when compared with small Bi number. Fig. 9 illustrates the effect of the Stefan number on the timewise variation of the tube wall temperature. As can be seen from the figure, the tube temperature curves for Ste ¼ 0.45 lie below those for Ste ¼ 0.225 and the tube temperature decrease in a fixed heat transfer tube section (e.g. from X ¼ 250 mm to X ¼ 500 mm) is greater for Ste ¼ 0.45 than that for Ste ¼ 0.225. Moreover, higher heat transfer rates appears for Ste ¼ 0.45 due to the lower gaseous coolant inlet temperature. 4.3. Cold energy recovered With respect to thermal storage applications, the amount of cold energy stored in the PCM freezing process is of great importance. By neglecting the heat losses from the surrounding into the test unit, the cold energy recovered by the LHCSS will be equal to the total enthalpy increase of the cryogenic gas, which is given by
_ p Tf ;out Tf ;in Q ¼ mc
(8)
_ is the mass flow rate of coolant (kg/s); cp is the specific where m heat of coolant (J/kg K1); Tf,in and Tf,out respectively represent the inlet and outlet temperatures of the coolant flow (K).
Fig. 8. The effects of Biot number on the tube wall temperature.
Fig. 10 shows the effects of Biot number on the cold per unit tube length for the same Stefan number. At a fixed Stefan number, the stored cold is seen to be quite sensitive to the Biot number and increases remarkably as the Biot number increases. Due to the growth of frozen layer, the ice thermal resistance increases as the freezing process progresses and the overall thermal resistance increases correspondingly. Hence, it might be expected that the heat transfer rate should decrease with the thickening of the ice layer. However, as shown in the figure, the decreasing trend of the heat transfer rate is moderate, and this can be explained by the fact that the overall thermal resistance is not so sensitive to the increasing frozen layer thermal resistance. Actually in this study, the ice-layer thermal resistance just accounts for less than ten percent of the overall resistance at the end of the freezing process [14]. Therefore, the cold energy storage rate can be simplified to a constant as marked in Figs. 10 and 11. The effects of Stefan number on the recovered cold energy per unit tube length for the same Biot number are illustrated in Fig. 11. The heat flux for Ste ¼ 0.675 is 1.25 times larger than that for Ste ¼ 0.217. It can be explained by the fact that for a large Stefan number, larger PCM-to-coolant temperature difference results from the lower cryogenic gaseous coolant inlet temperature. Correspondingly, the logarithmic mean temperature difference increases as the Stefan number increases at fixed Bi. Therefore, the heat transferred from PCM to coolant will increase accordingly.
Fig. 10. The effects of Biot number on the cold recovered for Ste ¼ 0.562.
H. Tan et al. / International Journal of Thermal Sciences 50 (2011) 2220e2227
2227
2) The Stanton number reflects the convective heat transfer capabilities of the gaseous coolant. Larger Stanton number means that more heat is transferred into the coolant and a larger rise occurs in the coolant temperature. This is responsible for the greater decrease of the frozen layer thickness along the tube length. 3) In LHCSS the stored cold energy per unit length can be simplified to a linear function of time. The cold-storage rate is seen to be quite sensitive to the Biot number and the Stefan number. Acknowledgments This work was supported by the High-Tech Research and Development Program under the contract No. 2007AA05Z216 in China. Fig. 11. The effects of Stefan number on the cold recovered for Bi ¼ 0.432.
5. Conclusions This paper presented a numerical analysis on the PCM freezing process used to recover and store cryogenic cold energy of gaseous coolant. Based on the Solidification and Melting model, the twodimensional phase change problem was simulated in the Computational Fluid Dynamics (CFD) code FLUENT. An experimental investigation was conducted to validate the CFD model and the simulation results showed a good agreement with experimental results. When a flow of cryogenic gas passes through a horizontal copper tube submerged in liquid PCM, solidification of PCM will occur on the outer surface of the tube. The major characteristics of PCM freezing in this investigation are the weak convective heat transfer inside the tube, the large thermal resistance of the gaseous coolant that dominates the overall thermal resistance, and the large temperature increase experienced by the cryogenic gas. Effects of Biot number, Stanton number, and Stefan number on the PCM freezing characteristics have been analyzed in detail. The main conclusions can be summarized as follows: 1) The frozen layer growth performance is sensitive to Biot number and Stefan number. A larger Biot number results in a relatively higher ice layer increasing rate at a fixed axial position due to an increase in the internal convective heat transfer coefficient. Higher values of Stefan number leads to lower coolant inlet temperature and greater temperature difference between the phase-changing region and the gaseous coolant inside the tube, and consequently fast growth rates of frozen layer.
References [1] Na Zhang, Noam Lior, A novel Brayton cycle with the integration of liquid hydrogen cryogenic exergy utilization, Int. J. Hydrogen Energ. 33 (1) (2008) 214e224. [2] Qiang Wang, Yanzhong Li, Jiang Wang, Analysis of power cycle based on cold energy of liquefied natural gas and low-grade heat source, Appl. Therm. Eng. 24 (4) (2004) 539e548. [3] Qiang Wang, Yanzhong Li, Chen Xi, Exergy analysis of liquefied natural gas cold energy recovering cycles, Int. J. Energ. Res. 29 (2005) 65e78. [4] Na Zhang, Noam Lior, A novel near-zero CO2 emission thermal cycle with LNG cryogenic exergy utilization, Energy 31 (10e11) (2006) 1666e1679. [5] Jan Szargut, Ireneusz Szczygiel, Utilization of the cryogenic exergy of liquid natural gas (LNG) for the production of electricity, Energy 34 (7) (2009) 827e837. [6] Belén Zalba, José M. Marín, Luisa F. Cabeza, et al., Review on thermal energy storage with phase change: materials, heat transfer analysis and applications, Appl. Therm. Eng. 23 (3) (2003) 251e283. [7] A.M. Khudhair, M.M. Farid, A review on energy conservation in building applications with thermal storage by latent heat using phase change materials, Energ. Convers. Manage. 45 (2) (2004) 263e275. [8] Akio Saito, Recent advances in research on cold thermal energy storage, Int. J. Refrigeration 25 (2) (2002) 177e189. [9] Babak K. Soltan, Morteza M. Ardehali, Numerical simulation of water solidification phenomenon for ice-on-coil thermal energy storage application, Energ. Convers. Manage. 44 (2003) 85e92. [10] L.F. Cabeza, H. Mehling, S. Hiebler, F. Ziegler, Heat transfer enhancement in water when used as PCM in thermal energy storage, Appl. Therm. Eng. 22 (2002) 1141e1151. _ [11] Aytunc Erek, Zafer Ilken, M. Ali Acar, Experimental and numerical investigation of thermal energy storage with a finned tube, Int. J. Energ. Res. 29 (2005) 283e301. [12] N. Kayansayan, M. Ali Acar, Ice formation around a finned-tube heat exchanger for cold thermal energy storage, Int. J. Therm. Sci. 45 (2006) 405e418. [13] Badr A. Habeebullah, An experimental study on ice formation around horizontal long tubes, Int. J. Refrigeration 20 (2007) 789e797. [14] Hongbo Tan, Yanzhong Li, Hanfei Tuo, Man Zhou, Baocong Tian, Experimental study on liquid/solid phase change for cold energy storage of LNG refrigerated vehicle, Energy 35 (5) (2010) 1927e1935. [15] J.P. Holman, W.J. Gajda, Experimental Methods for Engineers. McGraw-Hill, New York, 1994.