Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 306 Effect of turbulent Prandtl number in the convective turbulent heat transfer modelling Ahmed Abed Al-Kadhem Majhool Mechanical Engineering Department, College of Engineering, Al-Qadissiya University ahmidkadhim@yahoo.com Received 1 September 2014 Accepted 19 October 2014 ABSTRACT The aim of this work is to investigate the capability of the numerical predications for both dynamic and thermal fields in a turbulent flow by using eddy-viscosity turbulence model. A special case was chosen to predict forced convection heat transfer in low turbulent Reynolds number and very low Prandtl number fluids and investigate the sensitivity of these predictions to the type of the in-house code and to several input parameters. Results obtained with a standard LRN (k-ε) turbulence model at relevant Reynolds number for the fully developed turbulent channel flow (Reτ = 180 and Pr= 0.025), are presented and discussed. The mean velocity and temperature profiles agree very well with that of a reference Direct Numerical simulation (DNS). Keywords: Turbulent Prandtl number, Turbulent heat flux, Channel flow, Turbulence modelling تاثير عدد براندال االضطرابي على نمذجة وسط انتقال حراري مضطرب الملخص ي في حيز مضطرب يعتمد على موديل لزوجة كي والحرارالتنبوء العددي لكال المجالين الحرالهدف من هذا البجث هو تحقيق امكانية النتقال حراري قسري باستخدام تيارات الحمل في وسط مضطرب عند الدوامة. الحالة الخاصة المختارة في هذا البحث هي التنبؤء ا بتغير في التنبؤات في برنامج صمم من قبل المؤلف متبوع هعدد رينولدز مضطرب قليل وعدد برندال قليل جدا وكذلك تحقيق هذ متغيرات االدخال. النتائج المستخلصة باستحدام موديل اضطراب يعتمد على نمذجة لزوجة الدوامة لعدد رينولدز منخفض لقناة في mailto:ahmidkadhim@yahoo.com Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 307 معدل توزيع السرعة ودرجة قد اظهرت ونوقشت. ,(Reτ = 180 and Pr= 0.025)حالة اضطراب كامل في ظروف جريان بالمقارنة مع المصدر من المحاكة العددية المباشرة. الحرارة يظهر نتائج ممتازة 1. Introduction The Prandtl number is defined as the ratio between the viscous diffusion rate to the thermal diffusion rate and low Prandtl numbers of fluids are fluids which have a Prandtl number much lower than unity. In this case, the thermal diffusivity of the fluid is dominant. This property affects the boundary layer thickness of the thermal and the momentum boundary layer of moving fluids. For low Prandtl number of fluid, the thermal boundary layer is much thicker than the momentum boundary layer. The transfer of heat in turbulent flow fields at low range of Prandtl numbers can be affected by the turbulence for two reasons. Firstly, the molecular diffusivity of fluid plays a more important direct role because it has the same order as the turbulent diffusivity. Secondly, molecular motion causes thermal tags to escape from eddies and therefore reduces the effectiveness of the turbulence in mixing. A modified low Reynolds number turbulence stress and heat flux equation model was developed by Nishimura [1]. The low Reynolds number turbulence model was implemented for numerical simulation of the turbulent-laminar transition. The old existing low Reynolds number turbulence models generally require very thin mesh width between a wall and the first computational node to assure the accuracy in calculation. The modified version was based on a criterion for the distance between the wall and the first computational node, in which non-dimensional distance (y + ) must be less than 0.5. A sequence of direct numerical simulations (DNSs) of turbulent heat transfer in a channel flow have been carried out by Kawamura et al. [2] with different boundary conditions for velocity and temperature fields. In the velocity field, the turbulent Poiseuille and Couette flows were used; while in the temperature field, the uniform heat flux heating and the constant temperature difference were employed. They observed that the turbulent Prandtl number is close to the unity in the near-wall region independent of the boundary conditions; whereas in the central region it depends mainly upon the conditions. Therefore the turbulent Prandtl number cannot be assumed to be a constant. Direct numerical simulations study for the velocity and temperature fields for turbulent flow in a channel have been carried out by Piller et al. [3] in order to examine the influence of Prandtl number on turbulent transport. They recorded the influence of Prandtl number on the eddy diffusivity and on statistical properties of the fluctuating temperature field. From their noticed that the ratio of the turbulent diffusivity at certain Prandtl number Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 308 to the turbulent diffusivity at a Prandtl number of unity is roughly constant over the whole cross- section of the channel. A numerical simulation of forced convection over a vertical backward-facing step was studied by Otic et al. [4] which used a hybrid between large eddy simulation and direct numerical simulation methods, where the study for analyzing turbulent liquid metal flow. The numerical simulation was conducted at a Reynolds number of 10,000 based on the step height and maximum mean inlet velocity, Prandtl number of 0.006, and an expansion ratio of 1.20. The modeling convective turbulent heat transfer of a liquid metal was presented by Bricteux et al. [5] with a Prandtl number down to 0.01, which is the order of magnitude of lead-bismuth eutectic in a liquid metal reactor. They observed that the heat transfer for a low Prandtl number fluid could be molecular even though the flow is turbulent, making the common Reynolds analogy unsuitable. Thermal fluctuations in a heated periodic channel have been modeled by Monod et al. [6] by using Large Eddy Simulation (LES). The fluid Prandtl number was chosen equal to 0.01 and the friction Reynolds number to 395. They showed that, in low Prandtl fluids, the characteristic scales of thermal fluctuations are much bigger than the dynamic ones. Both Isothermal and isoflux boundary conditions were considered in the work along with conjugate heat transfer boundary conditions. Bricteux et al. [7] used different approaches for the numerical prediction of heat transfer in a turbulent channel flow at very low Prandtl number and high Reynolds number to compare among them. Their investigations permitted to suggest that for the range of Reynolds numbers envisaged the Chien model can handle with to be the most proper one. They also showed that the correlation of Kays allows to obtain reliable prediction of the temperature profile. In this paper, standard LRN (k-ε) turbulence model at relevant Reynolds number for the fully developed turbulent channel flow for evaluating the turbulent Prandtl number and fixed values are compared and their performance are assessed against DNS data. 2. Mean Flow Equations The governing equations for fluid motion are the Navier-Stokes equations. The flow is assumed to be a steady, turbulent and incompressible flow. When the density of a viscous fluid is constant, the equations are sufficient to model the flow in general form can be described in terms of the conservation of mass, momentum and energy equations, which can be written, in Cartesian coordinates as: Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 309 where the total derivative is defined as: and where the over line symbol represents the mean quantity while the lowercase letters correspond to the fluctuating component of that quantity . where μ and σt are the fluid viscosity and turbulent Prandtl number respectively. 3. Constant Wall Heat Flux Model In fully developed heat transfer, the mean temperature is varying linearly which is increased as a function of, x, and the increased rate of heat transfer of increase can be determined by applying an energy balance to a differential element of the channel. In this case, the constant wall heat flux, qw, as shown in Figure 1 is used which the given condition yields Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 310 Where is the heat flux applied at the upper and bottom walls, k is the thermal conductivity, Cp is the specific heat and is the bulk mean velocity. The right hand side of equation (8) is not constant, due to the presence of bulk velocity and the varying distance along the computational domain. By using the definition of Prandtl number, Pr =ν/α, the averaged energy balance equation becomes The above equation contains two unknown terms. The first term is the functions of mean temperature and the second term is turbulent heat flux. These are the quantities of primary interest. 4. The proposed Turbulence Model In this section, the eddy-viscosity models used in this study will be presented. The standard (k-ε) model of Launder and Spalding [9] is a Two Equation model requires two transport equations, one for k (the turbulent kinetic energy) and the other for (its turbulent dissipation rate) ǫ, to describe the turbulence. The low-Reynolds-number (k- ε) model of Launder and Sharma [8] contains certain modifications which have to be made if one wants to apply this model to near wall regions, where the Reynolds number is low. Therefore the k and ε equations are thus: where Pk, is the production term created by mean shear, defined as : 1- The first modification is the presence of damping functions in order to account for the near-wall region and it is done through introducing a viscous damping function, fµ, into the turbulent viscosity equation. Therefore, Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 311 The function fµ is used to account for both the true viscous damping at low Reynolds number and then decreases across the viscous sub-layer and the preferential damping of the wall-normal fluctuations as the wall is approached. In this work, the damping function fµ is given by: where the turbulent Reynolds number is defined as: 2- The second modification is that the ε equation is solved for , where the latter represents the isotropic part of the turbulent eddy dissipation rate, therefore it can be defined as : 3- The third modification is to introduce and add a new source term, this term accounts for extra dissipation near the wall in order to modify the turbulent kinetic energy prediction in this region where velocity gradients are changing rapidly. The new source term is commonly known as E-term, which is defined as: 5. Case Studied The standard test case of a flow field within a two dimensional, steady - state, incompressible turbulent flow is used to validate the solver. The problem of the physical case of interest in this study is the developed turbulent flow of a low Prandtl number fluid subjected to a uniformly heat flux in a channel. The case studied here is that of fully-developed flow through a plane channel of half-width (h) with a Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 312 wall of non-zero thickness (d) uniformly heated at its outer surface (as shown schematically in Figure 2). The initial values which are used in the present work are listed in Table 1. 6. Numerical Treatment The one dimensional Finite-Volume Method (FVM) is used to discretize the governing equations on structured grids. The Navier - Stokes equations are solving using SIMP LE algorithm of Patankar and Spalding [10]. To discretize the governing equations for fluid flow, the cell-centered finite volume method is selected. In this algorithm the overall solution procedure is iterative and is based on a pressure-correction equation is derived from the discretized equations for continuity, momentum, energy, turbulent kinetic energy and turbulent dissipation rate. All simulated cases are carried out using non-uniform, with 210 nodes in the fluid region with an expansion ratio of 1.025. It is worth mentioning that in order to fix the Reynolds number the height of the half-channel is taken to be 0.5 m. In the context of turbulent heat transfer, the boundary condition of Isoflux where the wall is imposed to a constant heat flux. The iterative method which is implemented in the in-house code to solve the discretized equations is known as Tri-Diagonal Matrix Algorithm (TDMA). 7. Results and discussions The mean velocity profile is presented as a function of y+ in Figure 3 and compared with the numerical direct numerical simulation (DNS) data of a low Reynolds number flow by Kawamura et al. [2]for plane channel flow at Reτ = 180. The present result is in good agreement with the numerical data. By taking DNS data of [2] as a reference, the discrepancy of 6.5% on the bulk velocity is calculated which the present result is in good agreement with the numerical data. It is found a slight deviation from the reference data in the logarithmic region. This is due to a characteristic difference between the channel and the boundary layer flows. The logarithmic region can be found even in the case of the lowest Reynolds number. In this work the effect of turbulent Prandtl number cannot be affected on the mean velocity distribution. Therefore, it can be shown from the figure that one curve is present versus the DNS data. In the near-wall region, the model has the accurate asymptotic behavior and it is expected that the proposed model is applicable in turbulent channel flows. The dimensionless temperature profile Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 313 where Tw, is the wall temperature and normalized by the friction temperature is displayed in Figure 4a for various turbulent Prandtl numbers. The profile is obtained at low Prandtl number of (0.025) which matches very well to that of the DNS of Kawamura et al. [2], which validates the numerical methods. In order to figure out the obtainable results the figure has been magnified as reported Figure 4b, only in a certain place to investigate the effect of turbulent Prandtl number. There is no remarkable discrepancy of mean temperature exist for all fixed Turbulent Prandtl numbers [0.6,0.9,1.8,2and 4] except when the turbulent Prandtl number is variable in the linear thermal regime. This is due to the law (T + = Pr∗ y+). Consequently, the variable turbulent Prandtl number is set in the code as defined in the literatures in the form of, Prt = νt/αt. It can be observed that in the thermal Log- law that below of Prt = 1.8 the predicted profiles are underestimated the temperature over the channel height compared with the DNS data. In the case of variable Prt the overestimation of the turbulent heat diffusivity αt due to inaccurate modelling for this parameter because the large thermal diffusion which is caused by the much larger turbulent intensity in the central region. Figure 5 shows the local temperature distribution along the computational domain by choosing a fixed five cases and one case variable turbulent Prandtl number. From the results and since the momentum and thermal spectrum are not similar in the case of low Prandtl number fluids, it can be concluded that changes throughout the flow field, with the largest changes occurring close to the wall. Going farther, to investigate the effect of turbulent Prandtl number on the heated turbulent behavior, the turbulent heat flux results from the combination of temperature and velocity fields. Thus it is also important to recognize the convergence nature of velocity fields which in turns affects the convergence of turbulent heat fluxes. Figure 6 shows the distribution of turbulent heat flux normalized by the friction velocity and temperature for six different turbulent Prandtl numbers compared with DNS data. From this comparison, it is clear that the peak of the turbulent heat flux decreases with the increase of turbulent Prandtl number and moves towards the wall as the conductive sub-layer becomes thinner. In addition the most suitable value is when Prt= 0.9 which leads to the effect of turbulent has the larger effect than the thermal field. The efficiency of the heat transfer from the heated solid wall to the fluid flow is described by the Nusselt number. For constant wall heat flux boundary conditions the Nusselt number is computed as Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 314 For more assessment to the above predictions and investigation Figure 7a shows the variation of Nusselt number at different turbulent Prandtl numbers. All simulated cases show the same trend because the large thermal diffusion is caused by the much larger turbulent intensity in the central region of the channel flow. By using the same way the Figure 7a is enlarged to show the overestimated results for turbulent Prandtl numbers less than 1.8. Finally Figure 8 shows variation of the turbulent Prandtl number along the channel height. The investigation will be only to justify the variable turbulent Prandtl number. In the near wall region, is close to the molecular Prandtl number value which is almost independent of both the thermal and velocity boundary conditions. In the central region, however, depends appreciably on the thermal and turbulent velocity fields and offers very low values approximately equal to 0.345. 8. Conclusions Despite obtaining fairly promoted results during the course of the case studied on the topic of fully developed turbulent channel flow with heat transfer, further work on this topic can be done in order to obtain better results for turbulent heat flux which is in more agreement with the available DNS data. A fixed and variable values of turbulent Prandtl number is used in the present work, but studies have shown that these constant values are have under predicted results across the channel and the variable case has an over predicted results. This could introduce inaccuracy in the present results and an alternative method should be found to model this. From this work, it might be conclude that the turbulence model needs to handle with the thermal eddy viscosity term. Since it was taken to be constant in this work. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 315 9. References [1] M. Nishimura. Development of A Low Reynolds Number Turbulence Stress and Heat Flux Equation Model. 7th International Conference on Nuclear Engineering, Tokyo, Japan, April 19-23, 1999. [2] H. Kawamura, H. Abe and K. Shingai. DNS of turbulence and heat transport in a channel flow with different Reynolds and Prandtl numbers and boundary conditions. Proceedings of 3rd International Symposium Turbulence, Heat and Mass Transfer, editors Y. Nagano et al., 15–32, 2000. [3] M. Piller, E. Nobile and T. J. Hanratty. DNS study of turbulent transport at low Prandtl numbers in a channel flow. Journal of Fluid Mechanics. 458, 419–441, 2002. [4] I. Otic , A. G. Class and T. Schulenberg. Numerical Investigation of Turbulent Low-Prandtl- Number Forced Convection Over a Vertical Backward-Facing Step. 7th International Topical Meeting on Nuclear Reactor Thermal Hydraulics, Operation and Safety Seoul, Korea, October 5-9, 2008. [5] L. Bricteux , M. Duponcheel and Y. Bartosiewicz. Direct and Large Eddy Simulation of Turbulent Heat Transfer at Very Low Prandtl Number: Application to Lead-bismuth Flows. Proceedings of the International Conference Nuclear Energy for New Europe 2010, Portoro, Slovenia, Sept.6-9, 2010. [6] R. Monod , G. Brillant , A. Toutant and F. Bataille. Large Eddy Simulations of a turbulent periodic channel with conjugate heat transfer at low Prandtl number. Journal of Physics: Conference Series 395, 2012. [7] L. Bricteux, M. Duponcheel, M. Manconi and Y. Bartosiewicz. Numerical prediction of turbulent heat transfer at low Prandtl number. Journal of Physics: Conference Series 395, 012044, 2012. [8] B. E. Launder and B. I. Sharma. Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc. Letters in Heat and Mass Transfer, 1:131138, 1974. [9] B. E. Launder and D. B. Spalding. Mathematical models of turbulence, London academic press, 1972. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 316 [10] S.V. Patankar. Numerical heat transfer and fluid flow. Hemisphere publishing corporation, Taylor and Francis group, New York, 1980. Nomenclature Cp Specific heat capacity at constant pressure Cµ Constant in equation 13 Cε1 Constant in equation 12 Cε2 Constant in equation 12 D Total derivative d Thickness Eε Molecular effect in equation 12 fε Damping function in equation 12 fµ Damping function in equation 14 H Convection coefficient h Height of half channel k Turbulence kinetic energy Nu Nusselt number P Mean pressure Pk Production term of k qw Wall heat flux Pr Prandtl number Prt Turbulent Prandtl number Re Reynolds number Reτ Turbulent Reynolds number T Temperature T + Dimensionless Temperature uiθ Turbulent heat flux vector U + Dimensionless velocity UB Bulk velocity x Cartesian coordinates in stream-wise directions y Cartesian coordinates in wall-normal directions Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 317 Y + Dimensionless distance from the wall Greek Symbols ε Dissipation rate of the turbulence kinetic energy ε ~ Homogenous part of ε µ Dynamic viscosity µt Turbulent dynamic viscosity ν Kinematic viscosity νt Turbulent kinematic viscosity σt Turbulent Prandtl number Figure 1: Schematic showing the constant wall heat flux. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 318 Figure 2: Geometry of plane channel flow with heated wall. Figure 3: Mean velocity distribution. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 319 (a) (b) Figure 4: Profiles of mean dimensionless temperature. Figure 5: Local temperature profile at several turbulent Prandtl numbers. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 320 . Figure 7: Variation of Nusselt number along the channel height. Figure 6: Turbulent heat flux with different turbulent Prandtl numbers Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 4 ….2014 321 Figure 8: Variation of turbulent Prandtl number along the channel height.