Applied Thermal Engineering 24 (2004) 2501–2509 www.elsevier.com/locate/apthermeng
Numerical model for two-phase solidification problem in a pipe R. Conde, M.T. Parra *, F. Castro, J.M. Villafruela, M.A. Rodrıguez, C. Mendez Dpto. Ingenierıa Energetica y Fluidomecanica, E. T. S. Ingenieros Industriales, Universidad de Valladolid, Paseo del Cauce s/n, Valladolid 47011, Spain Received 21 January 2004; accepted 22 April 2004 Available online 19 June 2004
Abstract A two-dimensional model of heat transfer and solidification of a laminar flow inside a tube has been completed using the commercial fluid dynamics code, Fluent . The non-equilibrium balance of heat flux on the interface fluid/solid provides the blockage of the pipe. The model is validated with experimental results of blocking lengths. In addition a correlation is obtained to predict the blocking lengths of different fluids and operating conditions such as pipe diameter, mean velocity, wall and liquid temperatures. 2004 Elsevier Ltd. All rights reserved. Keywords: Cooling process; Solidification; Numerical model
1. Introduction Solidification process is specially important in many fields such as the manufacturing industry because of the alloy casting in moulds [1] or the material phase-change in storage system that uses paraffins, salts or fatty acids because of their latent heat [2]. It is commonly known that the classic Stefan problem (19th century) looks for the temperature distribution and freezing-front movement of a solidifying fluid. During the phase change in a pipe, the solid–liquid interface moves from the wall toward the interior of the pipe. While the solid phase increases its volume, the surface heat flux decreases due to the increasing thermal resistance of the no-molten medium. There is a blockage when the pipe length is large enough to avoid that
*
Corresponding author. Tel.: +34-983-42-33-13x4414; fax: +34-983-42-33-63. E-mail address:
[email protected] (M.T. Parra).
1359-4311/$ - see front matter 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.applthermaleng.2004.04.010
2502
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509
the interface reaches an equilibrium balance between the heat release from the fluid and the conduction from the solid phase. The aim of this study is to identify the blocking length as well as to analyse the influence of different parameters. These are geometric characteristics of the duct such as the diameter, operating conditions as the mean velocity and initial temperature of the flow as well as the temperature of the duct wall. Also interesting is the behaviour of different fluids. There is a wide variety of fluids that suffer solidifications on industrial processes ranging from liquids metals (Prandtl number around 102 ) to oils (Prandtl number around 103 ). In this paper, the validity of a correlation to predict the blocked location in a tube as a function of the previous parameters is analysed.
2. Model description, assumptions and boundary conditions A numerical model has been built up to simulate the two-dimensional, axisymmetric solidification of a liquid laminar flow in a cylindrical pipe. The simulations have been carried out with the Fluent v. 5.4.8 code. The geometrical model, presented in Fig. 1, reproduces a straight tube of 3 m length with a constant transversal area of diameter D. The surface of the pipe has a temperature Tw lower than the melting point of the fluid to provide its change of phase. The velocity profile is the characteristic one of a perfect developed laminar flow with a mean velocity V . The initial fluid temperature is T0 higher than the melting point. The independence of results to the grid refinement provided a computational domain of 6000 · 20 rectangular cells which means a spatial resolution Dz ¼ 0:5 mm/cell and Dr ¼ 0:12 mm/cell. Thermal profiles in the melt and throughout the solid were calculated as well as the liquid fraction and the velocity profiles of the fluid flow. Properties of the melt Eq. (1) and solid phase (s) materials are given in Table 1. The main numerical difficulty is that this version of the code does not have the facility to adapt the mesh to the interface movement and, therefore, it is unable to analyse the freezing of materials with strong difference of solid/liquid density. An assumption made on the numerical model is that the change of density between liquid and solid phase was not considered, [3]. Every mass flow or pipe diameter requires two kinds of simulations. The first one is a steady simulation with uniform temperature to provide the velocity profile of laminar flow. The second one is a transient simulation with heat transfer that uses the former velocity field as initial condition and a pressure inlet as boundary condition. When the solid phase is forming, the open Tw < Tf V
Solid r D/2
T0 > T f z
Fig. 1. Scheme of the physical model.
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509
2503
Table 1 Thermal properties of materials used in this study Material
Density q (kg m3 )
Water
998
Olive oil
914
Propulsion fuel Aluminium
8700 2690
Specific heat c (J kg1 K1 )
Thermal conductivity K (W m1 K1 )
Latent heat L Melting (kJ kg1 ) point Tf (K)
Viscosity l (kg m1 s1 )
4182 (l) 2037 (s) 2257 (l) 2048 (s) 502 (l) 544 (s) 1200 (l) 1060 (s)
0.6 (l) 2.2 (s) 0.2 (l) 0.2 (s) 2.9 (l) 2.5 (s) 100.0 (l) 200.0 (s)
333
273
0.001
209
266
0.084
280
3043
0.005
391
933
0.010
transversal area reduces and, therefore, the mass flow decreases. If the interface does not reach a heat transfer balance, the flow is blocked at a certain location. In that case, the solid phase will attain the axis of the pipe. The blocking length, Z, is determined as the z-position where the liquid fraction on the axis of the pipe is lower than 0.4. The conservation equations of mass and momentum are decoupled from the one of thermal energy. These equations are solved using a segregated solver with a second order upwind scheme. The transient model uses an implicit algorithm. 2.1. Modelling of the phase change The phase change model makes use of the melt/freeze enthalpy-porosity technique described in Refs. [4–7]. This facility is available in Fluent v. 5.4.8, [8,9]. The melt/solid interface is not explicitly tracked; instead a liquid fraction value is used for each cell in the computational domain. Hence, the mushy zone is a region where liquid-fraction values lie between zero and unity and its temperature ranges between the ones of the liquid and the solid. The energy equation used in heat transfer calculations, is expressed in terms of the enthalpy h at the temperature of the melt: Z T cp dT ð1Þ h ¼ href þ Tref
where href is a reference enthalpy, T the temperature, cp the specific heat capacity at constant pressure and Tref is the reference temperature. The liquid fraction can be defined as: 8 T > Tl <1 ð2Þ fl ¼ ðT Ts Þ=ðTl Ts Þ Ts 6 T 6 Tl : 0 T < Ts where T is the mushy zone temperature. The phase change takes place on the range of temperatures between Tl and Ts .
2504
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509
For a phase change problem, the energy conservation equation can be written as: o ðqh þ qfl LÞ þ divðqvhÞ ¼ divðkrT Þ ð3Þ ot Two unknowns are encountered in Eq. (3): enthalpy and liquid fraction. The last one is calculated as a function of the temperature according to the relation given by Eq. (2). The solution for T is obtained by iteration between Eqs. (2) and (3). In the solid region, the temporal evolution of the energy can be expressed as a function of conduction heat flux. Assuming that the liquid fraction and the velocity are zero, Eq. (3) is valid.
3. Validation of the numerical model In order to validate the performance of the numerical model, the conditions of the experimental work provided by Cheung and Baker [10] to predict the blocking length under different operating conditions have been reproduced. Fig. 2 shows the comparison between the blocking length of numerical results from our model and experimental results from Ref. [10]. These are obtained using water on different operating conditions such as mean velocity of 0.6, 0.95 and 1.4 m/s; pipe diameters of 4.7 and 7.3 mm and wall temperature of 77 and 195 K. Also an olive oil has been simulated. Both, experimental and numerical results match reasonably well providing numerical errors lower than 11% in the case of water. The biggest error occurs when the initial temperature of the liquid is high. With reference to the olive oil, the errors are important, up to 40%. The authors consider that there is an important uncertainty with reference to the properties of the olive oil, which are not specified on the Cheung and Baker’s work [10] and can suffer important variations from some sources to others.
Num;Wat;V=0.60m/s;D=4.7mm;Tw=77K 4.5 Num;Wat;V=0.95m/s;D=4.7mm;Tw=77K Num;Wat;V=1.40m/s;D=4.7mm;Tw=77K 4.0 Num;Wat;V=0.60m/s;D=7.3mm;Tw=77K Num;Wat;V=0.60m/s;D=4.7mm;Tw=195K Num;Oil;V=0.17m/s;D=4.7mm;Tw=195K 3.5
Exp;W at;V=0.60m/s;D=4.7mm;Tw=77K Exp;W at;V=0.95m/s;D=4.7mm;Tw=77K Exp;W at;V=1.40m/s;D=4.7mm;Tw=77K Exp;W at;V=0.60m/s;D=7.3mm;Tw=77K Exp;W at;V=0.60m/s;D=4.7mm;Tw=195K Exp;Oil;V=0.17m/s;D=4.7mm;Tw=195K
3.0 2.5
Z (m)
2.0 1.5 1.0 0.5 0.0 270
290
310
330
350
370
To (K)
Fig. 2. Comparison of numerical results and experimental results from Ref. [10].
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509
2505
As it was expected, a larger blocking length, Z, is provided by a higher value of: the initial temperature of the liquid phase, T0 ; the mean velocity flow, V ; the diameter of the pipe, D and the wall temperature, TW .
4. Numerical correlation for the blocking length On this section, the results of the previous section are treated according the procedure proposed by Cheung and Baker [10] to obtain a correlation that was able to predict the blocking length under different operating conditions. The blocking length on its non-dimensional form, Z , is a function of the following parameters: Z ¼
Z ¼ f ðWl ; Re; Ws ; Pr; as =al Þ D
ð4Þ
where Wl is the ratio between the effective heat, C, and the latent heat, L. The effective heat is a linear combination of the latent heat and the liquid sensible heat, being the proportionality factor c a characteristic of every fluid. This parameter depends on the initial temperature of the fluid, T0 . It is known that a hypothetic liquid with zero latent and specific heats will have a zero blocking length. Wl ¼
C L þ ccl ðT0 Tf Þ ¼ L L
Ws ¼ L=cs ðTf TW Þ
Ws is the Jakob number, which considers the relative importance of the latent heat with reference to the heat transfer on the solid phase. This parameter is a function of the wall temperature. Re ¼
qVD l
Pr ¼
lcl Kl
a¼
K qc
The remaining parameters are the Reynolds number involved with the mean velocity of the flow, the Prandtl number that characterises the fluid and finally, the ratio of solid and liquid thermal diffusivities as =al . It is experimentally verified that the blocking length is lineally proportional to the effective heat and have a potential relationship with the velocity, the wall temperature, the Prandtl number and the ratio of melt and solid thermal diffusivities. Therefore, the correlation to build up will have the form of Eq. (5). Z ¼ eWl Ren Wsm Pra ðas =al Þb
ð5Þ
First step is to obtain the value of the factor c for every fluid in order to calculate the effective heat C. Therefore, the linear correlation of the numerical blocking length as a function of the initial temperature of the fluid when the rest of parameters (V , D, TW and fluid) remain constants, let us identify the factor values: cwater ¼ 0:696; colive oil ¼ 0:698 and cpropulsion fuel ¼ 1:205. Fig. 3 shows the variation of the numerical blocking length versus the effective heat under different operating conditions. Table 2 has a summary of the linear correlation obtained from Fig. 3. The lower regression coefficient corresponds with the olive oil. The correlation analysis of the slopes corresponding with three velocities of water provides the exponent of the Reynolds number:
2506
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509 4
Wat;V=0.60m/s;D=4.7mm;Tw=77K Wat;V=0.95m/s;D=4.7mm;Tw=77K Wat;V=1.40m/s;D=4.7mm;Tw=77K Wat;V=0.60m/s;D=7.3mm;Tw=77K Wat;V=0.60m/s;D=4.7mm;Tw=195K Olive oil;V=0.17m/s;D=4.7mm;Tw=195K Fuel;V=1m/s;D=4.7mm;Tw=1660K
3.5
Z Num (m)
3 2.5 2 1.5 1 0.5 0 40
60
80 100 120 Effective Heat Γ (cal/g)
140
160
Fig. 3. Effects of the velocity, diameter and wall and liquid temperatures on the blocking length, Z.
Table 2 Linear correlation of the numerical blocking length versus the effective heat for different velocities, diameters, wall temperatures and fluids Material
V (m/s)
D (mm)
TW (K)
Znum (m) versus C (cal/g)
Water
0.6 0.95 1.4 0.6 0.6 0.17 1
4.7 4.7 4.7 7.3 4.7 4.7 4.7
77 77 77 77 195 195 1660
Z Z Z Z Z Z Z
Olive oil Fuel
¼ 0:0099 C ¼ 0:0149 C ¼ 0:02 C ¼ 0:0231 C ¼ 0:0189 C ¼ 0:025 C ¼ 0:007 C
R2 0.91 0.98 0.60 0.98 0.76 0.53 0.85
n ¼ 0:846. The modification of the diameter of the pipe is a way to verify the exponent of the Reynolds number. The variation of the wall temperature for the water cases allows to obtain the exponent of the parameter Ws , resulting, m ¼ 0:717. Finally, the influence of the remaining parameters requires the use of at least three fluids: water, olive oil and a propulsion fuel. The resolution provides the exponents a ¼ 0:896 and b ¼ 0:07 for the Prandtl number and the ratio of thermal diffusivities, respectively, and the proportionality constant e ¼ 0:049. The complete correlation obtained from the numerical blocking lengths is: Z ¼ 0:049Wl Re0:846 Ws0:717 Pr0:896 ðas =al Þ0:07
ð6Þ
5. Results and discussion of the parametric analysis The capability of the former correlation to predict the blocking length under different operating conditions and fluids from those used on previous sections is analysed.
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509 Temperature of aluminium (K) 1800 2000
1600
3
2200
2507
0.4
2.5
1.5
0.3
1 0.25
0.5 Wat;V=0.60m/s;D=4.8mm;Tw=50K Wat;V=0.60m/s;D=4.8mm;Tw=100K Wat;V=0.60m/s;D=4.8mm;Tw=150K Wat;V=0.8 m/s;D=4.8mm;Tw=77K Wat;V=1.2m/s;D=4.8mm;Tw=77K Oil;visc=0.03kg/ms;V=0.21m/s;D=4.7mm;Tw=195K Oil;visc=0.6kg/ms;V=0.17m/s;D=4.7mm;Tw=195K Aluminium;V=1.2m/s;D=4.7mm;Tw=830K
0 -0.5 -1 278
293 313 333 Temperature of water and oils (K)
Z num of aluminium (m)
Z num of water and oils (m)
0.35 2
0.2
0.15
363
Fig. 4. Numerical results from the parametric analysis.
The parametric analysis considers the solidification of water with velocities of 0.8 and 1.2 m/s; pipe temperatures of 50, 100 and 150 K. Also other fluids have been used as melt aluminium (caluminium ¼ 2:17) and two olive oils of viscosities 0.03 kgm1 s1 (coil1 ¼ 1:54) and 0.6 kgm1 s1 (coil2 ¼ 2:01). In both cases, the velocity value assures laminar flow regime. Fig. 4 indicates the results of the numerical model. Exp;Wat;V=0.6m/s;D=4.8mm;Tw=77K 4.5
Exp;Wat;V=1.4m/s;D=4.8mm;Tw=77K Exp;Wat;V=0.6m/s;D=4.8mm;Tw=195K 4 Num;Wat;V=0.8m/s;D=4.8mm;Tw=77K Num;Wat;V=0.6m/s;D=4.8mm;Tw=50K Num;Wat;V=0.6m/s;D=4.8mm;Tw=150K 3.5 Num;Oil;V=0.17m/s;D=4.8mm;Tw=195K
Exp;Wat;V=0.95m/s;D=4.8mm;Tw=77K Exp;Wat;V=0.6m/s;D=7.3mm;Tw=77K Exp;Oil;V=0.17m/s;D=4.8mm;Tw=195K Num;Wat;V=1.2m/s;D=4.8mm;Tw=77K Num;Wat;V=0.6m/s;D=4.8mm;Tw=100K Num;Alum;V=1.2m/s;D=4.8mm;Tw=830K
y=x
Z Corr (m)
3 2.5 2 1.5 1 0.5 0 0
0.5
1
1.5 Z (m)
2
2.5
3
Fig. 5. Verification of the correlation. Experimental values from Ref. [10].
2508
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509 1 Tw = 700 K
0.9 0.8
Tw = 600 K Tw = 500 K Tw = 400 K
Z Corr (m)
0.7
Tw = 300 K
0.6 0.5 0.4 0.3 0.2 0.1 0 900
950
1000
1050 To (K)
1100
1150
1200
Fig. 6. Prediction of the correlation for aluminium, D ¼ 0:07 m; V ¼ 0:1 m/s.
Fig. 5 shows the performance of the correlation with regard to the experimental results from Cheung and Baker’s work [10] and the numerical results from the parametric analysis. The solid dots correspond with the cases of the parametric analysis, whereas the rest of the dots represent the experimental values. The arithmetic averaged error of the correlation to the numerical result is 5.7% for the aluminium, 10.2 % for the water and 39% for the olive oil.
6. Exploitation of the correlation The obtained correlation is used to predict the blocking length of a laminar flow of melt aluminium under different wall temperatures and initial temperatures. Results are presented in Fig. 6. It was chosen a pipe diameter of D ¼ 0:07 m and a mean velocity of V ¼ 0:1 m/s because these are common values used down stream of a magnetic pump working to cast in moulds.
7. Summary A numerical model has been developed to analyse the freezing of water, olive oil or aluminium in cylindrical shape ducts under laminar flow conditions. A wide range of velocities, temperatures and diameters are imposed to analyse their influence. The model is validated with experimental results from other sources. The numerical model is able to provide good accuracy according to the dispersion of numerical and experimental blocking lengths.
R. Conde et al. / Applied Thermal Engineering 24 (2004) 2501–2509
2509
The exploitation of the numerical results provides useful information to build-up a correlation, which is able to predict the blocking length as a function of Prandtl, Reynolds and Jakob numbers as well as the thermal properties of the fluid and the liquid temperature.
References [1] Z. Huan, G.D. Jordaan, Galerkin finite element analysis of spin casting cooling process, Applied Thermal Engineering 24 (2004) 95–110. [2] P. Lamberg, Approximate analytical model for two-phase solidification problem in a finned phase-change material storage, Applied Energy 77 (2004) 131–152. [3] R. Conde, Simulaci on numerica de solidificaci on de fluidos en conductos circulares, Ph.D. Dissertation, University of Valladolid, 2004. [4] M. Cross, W. Voller, Assessment of weak solution techniques for solving Stefan problems, Numerical Methods in Thermal Problems (1979) 172. [5] V.R. Voller, M. Cross, Accurate solutions of moving boundary problems using the enthalpy method, International Journal of Heat and Mass Transfer 24 (1981) 545–556. [6] N. Shamsunder, E.M. Sparrow, Effect of density change on multidimensional conduction phase change, Transactions of ASME 98 (1976) 550–557. [7] A.B. Crowley, Numerical solution of Stefan problems, International Journal of Heat and Mass Transfer 21 (1978) 215–219. [8] Fluent 5.0 User’s Guide, Graphics Version 4.10-44, Fluent Inc., 1993. [9] V.R. Voller, C. Pankrash, International Journal of Heat Mass Transfer 30 (1987). [10] F.B. Cheung, L. Baker Jr., Transient freezing of liquids in tube flow, Nuclear Science and Engineering 60 (1976) 1–9.