Numerical simulation of turbulent mixed convective ...

1 downloads 0 Views 926KB Size Report
axial velocity profile becomes asymmetric in the fully developed region for lower Re and higher Riq. However, with increasing Re and decreasing Riq, the flow is ...
Numerical simulation of turbulent mixed convective pipe flow Sayeed Mohammed, Anuruddha Bhattacharjee, and Sumon Saha

Citation: AIP Conference Proceedings 1980, 050029 (2018); doi: 10.1063/1.5044365 View online: https://doi.org/10.1063/1.5044365 View Table of Contents: http://aip.scitation.org/toc/apc/1980/1 Published by the American Institute of Physics

Numerical Simulation of Turbulent Mixed Convective Pipe Flow Sayeed Mohammed a), Anuruddha Bhattacharjee b) and Sumon Saha c) Department of Mechanical Engineering, Bangladesh University of Engineering and Technology, Dhaka, Bangladesh a) Corresponding author: [email protected] b) [email protected], c)[email protected] Abstract. In the present study, three-dimensional numerical simulation of turbulent mixed convection in a smooth horizontal pipe filled with air is carried out using finite volume method. The variation of flow and thermal characteristics in the fully developed turbulent region is observed with the variation of Reynolds number (Re = 5×74100, 33130, 23000, 10480, 7410) for modified Richardson number (0.1 ≤ Riq ≤ 10), while the modified Grashof number (Grq) is kept constant at 5.49×108. The mathematical model of the present problem consists of ReynoldsAveraged Navier-Stokes and energy equations along with the realizable k-ε turbulence model. At first, validation of the present results is performed in terms of the variation of isotherm contour, axial velocity contour, dimensionless mean axial velocity (u+) and mean temperature (T+) profiles in the inner and outer region of the pipe with the previous work. Secondly, the effect of the variation of Re and the corresponding Riq in the mixed convection regime is observed on variation of the axial velocity and isotherm contours, u+ and T+ profiles. The results show that the axial velocity profile becomes asymmetric in the fully developed region for lower Re and higher Riq. However, with increasing Re and decreasing Riq, the flow is dominated by forced convection and the velocity profile becomes symmetric again. On the other hand, u+ and T+ profiles show similar increasing trend with increasing Re and decreasing Riq from inner to outer region of the pipe.

INTRODUCTION A VHCS (ventilated hollow core slab) model is used to activate the thermal mass of concrete slab to decrease the energy demand for making comfortable building environment. The VHCS was first proposed by Swedish engineers in development of Thermo-Deck system [1]. During summer days, floor and ceiling slabs behave as heat sink, which are cooled down through the ventilation of cool outdoor air through hollow core at night [2]. However, in winter, the function of the hollow core is reversed [3]. This paper considers VHCS applications, which is idealized by a smooth horizontal pipe. Many researchers carried out experiments on different thermal performance analysis considering turbulent flow with Reynolds numbers varying in between 9000 to 45000 [4-6]. In this range of Re, turbulence is generated in the hollow core which enhances heat transfer rate between the hollow core wall and the working fluid through forced convection. With these conditions, the temperature gradient causes the buoyancy induced natural convection which distorts temperature field, flow field, and affects the heat transfer behavior of hollow core [7]. Therefore, the performance of hollow core is considered to be dominated by turbulent mixed convection. Very few investigations in literature considered buoyancy induced natural convection with heat transfer in smooth horizontal pipe [8-10]. Grassi and Testi [11] have investigated the ability of the Reynolds-averaged Navier–Stokes (RANS) turbulence models in predicting the heat-transfer rate in turbulent mixed convection with a high Prandtl number (Pr = 12.4) fluid in a circular pipe. In their study, they considered Re = 5.750 and Ri = 9, and thus, the mixed convection was predominantly dominated by buoyancy induced natural convection. They have found that, the realizable k-ε (REAL k-ε) model provided an adequate representation of the heat transfer in the entry region, while the renormalization group k-ε (RNG k-ε) model produced better results in the fully developed or quasi-developed flow. Moreover, Petukhov and Polyakov [8] have studied temperature and velocity profiles experimentally for fully developed mixed convection of air (with a Prandtl number close to unity) in a circular pipe at Ri = 1.08 (where, Re = 43,000), and the heat transfer rate (in terms of Nu) at Ri = 0.83 (where, Re = 49,000). The Re values adopted in their study were higher than those typically encountered in VHCS applications. Mundhe and Deshmukh [12] have studied heat transfer for turbulent flow in a circular tube International Conference on Mechanical Engineering AIP Conf. Proc. 1980, 050029-1–050029-9; https://doi.org/10.1063/1.5044365 Published by AIP Publishing. 978-0-7354-1696-3/$30.00

050029-1

with air as working fluid. They have found that the separation of boundary layers causes reduction in temperature and velocity gradient in the flow region, which makes significant enhancement of heat transfer in turbulent region. Petukhov [13] has investigated heat transfer and friction in turbulent pipe flow with variable physical properties. The results show that, with constant fluid properties analysis hold good, though, with considerable variation in the physical properties caused disagreement between theoretical and experimental data. McKeon et al. [14] repeated the measurements Zagarola and Smits [15], of mean velocity profile in fully developed turbulent pipe flow using smaller pitot probe. The study confirmed the power law region near the wall and, for Reynolds numbers greater than 230 × 103 (Reτ > 5 × 103), a logarithmic region further out, however, the limits of these regions and some of the constants differ from those reported by Zagarola and Smits [15]. Considering the above, this paper focuses on numerical investigation of air flow structures, temperature fields, and heat exchange that occur within the hollow cores under mixed convection. This paper is to numerically predict the flow characteristics in a smooth horizontal pipe, in the fully developed turbulent region under mixed convection idealizing hollow core slab for VHCS applications. The concrete component typically used for VHCS applications is not considered in this study, since prime focus of this paper is on the investigation and analysis of the behavior of the pipe flow under mixed convection with turbulence. However, this study demonstrates the effect of the variation of Reynolds and Richardson number in the mixed convection regime which will definitely provide new scopes for further research and development of models for VHCS applications.

MODELING AND SIMULATION Physical Model A 2D view of the 3D horizontal smooth pipe designed using ANSYS Design Modeler 14.5 is shown in Fig. 1, where the diameter (d = 2R) and the length (L) of the pipe are 0.18 m and 12 m respectively. Air with an inlet velocity of uz is entering the pipe which is varied from 0.638 to 6.38 m/s. The inlet temperature (Tin) of air is " kept constant at 293 K. A constant heat flux (𝑞𝑤 ) of 100 W/m2 is applied uniformly over the wall of the pipe. The outlet is set to the outflow condition.

FIGURE 1. Two-dimensional azimuthal and cross-sectional views of a smooth horizontal pipe along with physical dimensions and boundary conditions.

Mathematical Model The mathematical model considered for the present work includes Reynolds-Averaged Navier-Stokes (RANS) and energy equations along with realizable k-ε turbulence model equations. The governing equations can be expressed as follows: .u = 0, (1) 1 u u = − [p +  (  + t ) u  + F, (2) 

 C p  uT  =. ( eff T ) ,

eff =  f + t ,

050029-2

(3) (4)



 ( u   ) k =     + 



 ( u  )  =     + 

t 

t    k  + Gk + Gb −  ,  k  

(5)

  2  ,    + C1 k C3 Gb − C2  k +   

(6)

The eddy viscosity µt is defined using turbulent kinetic energy, k and turbulent dissipation rate,  as follows: k2  t =  C , (7)



The governing parameters of the present simulation i.e. Re, Grq, Riq, and Pr are defined as follows: Grq g  qw" d 4 ud  Re = z , Grq = , Riq = 2 , Pr = , (8) 2   Re  The value of the model constants are C1ε = 1.44, C2ε = 1.9, k = 1.0 and ε = 1.2. In the above equations, 𝐶𝜇 is eddy-viscosity whose value is computed using Reynolds [16] formula, C3ε is a constant calculated by the following relation provided by Henkes et al. [17], u is the mean velocity vector field, u is the averaged velocity field,  is the dynamic viscosity, t is the eddy viscosity, Gk, and Gb are the generation of turbulence kinetic energy due to the mean velocity gradients and buoyancy respectively, 𝜎𝑘 and 𝜎 are the turbulent Prandtl numbers for k, and  respectively, p is the mean pressure, ρ is fluid the density, λt is the turbulent thermal conductivity, λf is the thermal conductivity of fluid. Thermo-physical properties of air at 293 K are presented in Table 1, while boundary conditions and list of the governing parameters used in the present simulations are presented in Tables 2 and 3 respectively. TABLE 1. Thermo-physical properties of air at 293 K [16].

Symbol Cp 𝜌 λf 𝛽 𝜈 Pr

Properties Heat capacitance Density Thermal conductivity of fluid Thermal expansion coefficient Kinematic viscosity Prandtl Number

Unit J/kg.K kg/m3 W/m.K 1/K m2/s -

Value 1005 1.18 0.0261 0.00335 1.55×10-5 0.70

TABLE 2. Boundary conditions for the present simulation.

Simulation Run 1 Run 2 Run 3 Run 4 Run 5

Air Inlet Velocity, 𝒖𝒛 (m/s) 6.380 2.850 1.980 0.902 0.638

Air Inlet Temperature, Tin (K)

Applied Wall Heat Flux, 𝒒"𝒘 (W/m2)

293

100

TABLE 3. List of governing parameters used in the present simulation.

Simulation Run 1 Run 2 Run 3 Run 4 Run 5

Reynolds Number, Re 74100 33130 23000 10480 7410

Modified Grashof Number, Grq

Modified Richardson Number, Riq

5.49×108

0.1 0.5 1 5 10

Numerical Method and Mesh Generation The present problem is solved numerically using finite volume method by means of a commercial software ANSYS FLUENT 14.5. The entire domain of the 3D horizontal smooth pipe is discretized into computational

050029-3

grids consist of quadrilateral elements. Three different grid arrangements are tested as shown in Fig. 2 in order to find the optimum one. The grids are simultaneously composed of structured and unstructured mesh. The unstructured mesh is auto generated through FLUENT 14.5. The structured region is done using edge sizing and inflation layer. This method helps to have fine regions near the wall which is instrumental to this study. This simulation is performed in steady state. The upwind scheme to discretize the governing equation is of second order. SIMPLE method is carried to couple the velocity and pressure fields. A convergence criteria for all governing equations has been set to 10-6. TABLE 4. Characteristics of the different grid arrangements considered in the present study.

Inflated Region Grid G1 G2 G3

Growth Ratio 1.2 1.1 1.05

(a)

No. of layers 5 10 32

Face Size Min (mm)

Max (mm)

5.98×10 2.99×10-3 1.75×10-3

0.598 0.299 0.35

-3

Element Number 150239 152101 162542

(b)

Edge Sizing 30 50 100

Velocity 0.9 1.0 1.2

(c)

FIGURE 2. Cross-section view of grid arrangement of the computational domain for (a) G1, (b) G2, (c) G3.

MODEL VALIDATION Validation has been done in terms of dimensionless turbulent velocity profile (u+) and temperature profile (T+). Ahmed et al. [18] plotted numerical and empirical profiles of velocity and temperature in the fully developed region along the vertical diameter of the pipe. His results hold good agreement with the empirical relations proposed by Polyakov [19]. In this study, validation is done for the profiles with Real k-ε model and the results are precisely accurate with the data of Ahmed et al. [18]. In Figs. 3 (a) and (b) u+ and T+ profiles are plotted against the dimensionless distance from the wall (y+) along the vertical diameter of the horizontal pipe. The trends of the two graphs are observed same. ( R − y)u + uz + Tw − T  qw" (9) y+ = , u = ,T = , u = w , T = ,  u T   C p u Here, 𝑢𝜏 is the friction velocity, 𝑇𝜏 is the friction temperature, 𝑇𝑤 is the wall temperature, and 𝑇 is the fluid temperature. To further demonstrate the validity of the present study, various flow representation in the fully developed regions are showed in terms of isotherm and velocity contours. Figs. 4 (a), (b), (c), (d), show that the contours lines of Ahmed et al. [18] and present work are quite similar. The deviations may be attributed due to the major simplification of the computational model. Therefore, it can be demonstrated that, the present simulation model is reasonably accurate to analyze.

050029-4

25

25 20

Present (Lower Half) Ahmed et al. [18] (Lower Half) Present (Upper Half) Ahmed et al. [18] (Upper Half)

20 15

u+

T+

15 10

10

5

5

0 0 10

101

102

y+

103

0 0 10

101

(a)

y+

102

103

(b)

FIGURE 3. Comparison of the variation of (a) u+ profiles and (b) T+ profiles in the fully developed region at z/d = 61.1 with Ahmed et al. [18].

(a)

(b)

(c)

(d)

FIGURE 4. Variation of velocity contours (a) Ahmed et al. [18], (b) present simulation, and isothermal contours (c) Ahmed et al. [18], (d) present simulation, along with vertical diameter in the fully developed region at z/d = 61.1 for Re = 23000.

050029-5

RESULT AND DISCUSSION Mean Stream-wise Axial Velocity Distribution The variation of vertical axial velocity profile with the variation of Reynolds number is shown in the Fig. 5. From the figure, it can be observed that, the increase of Reynolds number for the pipe, causes decrease of the buoyancy effect. Since, the velocity is lower, the shear forces are greater than the momentum of the air flow itself causing the profile to be asymmetric. However, due to the addition of constant heat flux over the wall of the pipe, the buoyancy effect is induced inside the pipe. As, temperature rises, the density of air becomes smaller, so, the hot air moves upward, and the cold air moves downward. This also creates an asymmetric profile as shown in the Fig. 5. Here, when the velocity of air is low, the buoyant effect is dominant. However, with the increase of velocity, this nature of the flow diminishes. On the other hand, an interesting point can be observed for Re = 74100. The velocity profile comes up with a wedge in the middle for that higher Re. This behavior is due to fact that, the air near to the wall is gaining temperature from the wall heat flux and thus getting lighter and moving faster. Although, the air around the core region could not receive the heat due to lack of adequate time for heat transfer.

1.0

y/R

0.5

0.0

Re = 74100, Riq = 0.1 Re = 33130, Riq = 0.5 Re = 23000, Riq = 1 Re = 10480, Riq = 5 Re = 7410, Riq = 10 Ahmed et al. [18] (Re = 23000, Riq = 1)

-0.5

-1.0 0.0

0.5

uz/ub

1.0

1.5

FIGURE 5. Effect of Reynolds number on the vertical axial velocity profile of a horizontal smooth pipe in the fully developed turbulent region (z/d = 61.1).

In Figs. 6, the dimensionless turbulent velocity (u+) is plotted against the dimensionless distance from the wall (y+) for vertical diameter of the pipe. In the viscous sub layer for y+≤ 5 the profile are mostly linear. They follow the normalized law of the wall (10) u+ = y+ , + In the overlap layer i.e. 5 ≤ 𝑦 ≤ 30 the velocity profiles are proportional with the logarithmic law of the turbulent boundary layer which can be expressed by 1 uz+ = ln y + + B, (11)



( )

Here, 𝜅 and B are constant which can be derived experimentally and their values are 0.40 and 5.0 respectively. This equation is termed as the universal turbulent flow velocity profile.

050029-6

25 20

25 Re = 74100, Riq = 0.1 Re = 33130, Riq = 0.5 Re = 23000, Riq = 1 Re = 10480, Riq = 5 Re = 7410, Riq = 10

20

15

u+

u+

15

10

10

5

5

0 0 10

101

y+ (a)

102

0 0 10

103

101

102

y+

103

(b)

FIGURE 6. Comparison between the u+ profiles in the fully developed region at z/d = 61.1 on the radius at (a) lower half of vertical diameter, (b) upper half of vertical diameter for various Reynolds number.

Mean Temperature Distribution In Figs. 7, T+ is plotted against y+ for vertical diameter of the pipe. From the figure, we can observe that, in viscous sub layer for y+≤ 5, the mean turbulent temperature is higher for lower Reynolds number due to the lower velocity, which enhances the heat transfer between the layers of fluid. For, 𝑦 + ≥ 30, in the lower half diameter of the pipe, the mean temperature reaches maximum and then lowers. This signifies that, at z/d = 61.1 near the core region of the pipe flow, the heat transfer has not occurred completely and is still in mixing phase. This can be viewed from Fig. 7, where, there is discrepancy in temperature contours in the core region. Further, down the stream, the heat transfer will be complemented, and, the profile would be more of a linear nature. In viscous sub layer for y+≤ 5 the temperature profile follows a certain law which is (12) T + = Pr y + , + For 5 ≤ 𝑦 ≤ 30, which is the overlap region, the temperature profiles follows a logarithmic expressions which is 1 T+ = ln ( y + ) +  ( Pr ) . (13)

T

25 20

25

Re = 74100, Riq = 0.1 Re = 33130, Riq = 0.5 Re = 23000, Riq = 1 Re = 10480, Riq = 5 Re = 7410, Riq = 10

20

15

T+

T+

15

10

10

5

5

0 0 10

101

y+

102

0 0 10

103

(a)

101

y+

102

103

(b)

FIGURE 7. Comparison between the T+ profiles in the fully developed region at z/d = 61.1 on the radius at (a) lower half of vertical diameter, (b) upper half of vertical diameter for various Reynolds number.

050029-7

Visualization of Flow and Thermal Fields along Pipe Cross-section In Figs. 8, velocity and isotherm contours are presented for Re = 74100, 10480, 7410 and for Riq = 0.1, 5, 10 respectively, in the fully developed turbulent region of the pipe at z/d = 61.1. From the Figs. 8, it can be evidently said that, with the increase of Richardson number, the buoyant effect of the flow increases. Therefore, the core velocity region moves downwards. Similar pattern can also be seen in the temperature contour. Here, the reason for the core moving downwards is due to the effect of gravitational acceleration. As, the velocity of the air flow is decreasing near the wall, the air regions are getting more heated due to the constant heat flux. Thus, the density gets lower and they rises up inside the pipe, while, the cold air goes down.

(a)

(b)

(c)

(d)

(e)

(f)

FIGURE 8. Comparison of velocity and isotherm contours (a), (d) Re = 74100 and Riq = 0.1, (b), (e) Re = 10480 and Riq = 5 and (c), (f) Re = 7410 and Riq = 10, at the cross section of a horizontal smooth pipe in the fully developed region at z/d = 61.1 whereas, the top row presents velocity contours, and the bottom row presents isotherm contours.

CONCLUSIONS A numerical study of mixed convective turbulent flow inside a horizontal smooth pipe is carried out at different governing conditions. The characteristics of velocity and temperature contours have been studied in the fully developed turbulent region, as well as, the mean flow distribution in terms of mean velocity and temperature profiles for different Reynolds and Richardson number. Based on the results obtained, the following conclusions can be drawn. Some important observations are pointed below: (i) As Reynolds number increase the buoyant effect in the pipe flow becomes feeble hence attains a symmetric velocity profile. (ii) The core region of the velocity and temperature contours goes down along the diameter of the pipe as the inlet velocity lowers (iii) In the fully developed region for higher Reynolds number heat transfer is not completed yet at z/d = 61.1. (iv) In the viscous sub layer the turbulent velocity profile follow a linear law but in the overlap regions it follows a logarithmic law.

050029-8

ACKNOWLEDGEMENTS The authors gratefully acknowledge the logistic support provided by Department of Mechanical Engineering, Bangladesh University of Engineering and Technology (BUET) during this research work.

REFERENCES 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. 19.

R. Winwood, R. Benstead, M. Mashrae and R. Edwards, Build. Serv. Eng. Res. T. 18(1), 1–6 (1997). C. Cheung and R. Fuller, EcoLibrium 6(2), 24–32 (2007). M. Standeven, R. Cohen, B. Bordass and A. Leaman, Build. Serv. J. 20, 37–42 (1998). M. R. Shaw, K. W. Treadaway and S. T. P. Willis, Renew. Energ. 5(5–8), 1028–1038 (1994). N. Barnard, Build. Serv. Res. Inf. Assoc., 28 (1994). M. J. Ren and J. A. Wright, Build. Environ. 33(1), 43–52 (1998). R. Winwood, R. Benstead, M. Mashrae, R. Edwards and K. Letherman, Build. Serv. Eng. Res. T. 15(3), 171–178 (1994). B. S. Petukhov and A. F. Polyakov, Heat Transfer in Turbulent Mixed Convection (Hemisphere, New York, 1988). W. Grassi and D. Testi, Int. J. Comput. Fluid D. 21(7–8), 267–276 (2007). D. C. Wilcox, Turbulence Modelling for CFD (DCW Industries, Inc., California, 1998). W. Grassi and D. Testi, ASME J. Heat Transf. 128(10), 1103–1107 (2006). S. V. Mundhe and P. W. Deshmukh, Int. J. Innovative Res. Adv. Eng. 1(3), 70–75 (2014). B. S. Petukhov, Adv. Heat Transf. 6, 503–564 (1970). B. J. McKeon, J. Li, W. Jiang, J. F. Morrison and A. J. Smits, J. Fluid Mech. 501, 135–147 (2004). M. V. Zagarola, and A. J. Smits, J. Fluid Mech. 373, 33–79 (1998). W.C. Reynolds Fundamentals of turbulence for turbulence modeling and simulation. Lecture Notes for Von Karman Institute Agard Report No. 755 (1987). A. W. M. Henkes, F. F. van der Flugt, and C. J. Hoogendoorn, Int. J. Heat Mass Transfer, 34:1543-1557 (1991). F. Ahmed, R. Gianluca, F. Francesco and L. Chengwang, ASME J. Heat Transf. 138, 12501–12511 (2016). A. F. Polyakov, ASME J. Appl. Mech. Tech. Phys. 15(5), 632–637 (1974).

050029-9