. Theoretical modeling techniques on resolving turbulent heat flows in a non-dimensionalcircular tube mounted obstacle using the wall Y+ as guidance in selecting the appropriate grid configuration and corresponding turbulence models are investigated using CFD Code. The results obtained shows that the heat fluxes of 20, 23, 30 and 40 kW/m2, increases as the Y+wall profile moves away from the near – wall region, this is due to the effect of viscosity, buoyance, acceleration and the friction of the turbulence modification.The results also indicates that the low Y+wall treatment is suitable only for low Reynolds turbulence models in which it is assumed that the viscous sub-layer is properly resolved. The simulated results obtained in this research are in good agreements with the experimental results in the literature, even though they over predicted the observed heat transfer deterioration both quantitatively and qualitatively. Keywords: Wall Y+, Turbulence, Supercritical fluid, Heat flux Corresponding author: enarosehaomamoke@gmail.com; eomamoke@delsu.edu.ng 1.0 INTRODUCTION With its high density features like liquids, large diffusivity like gases and excellent dissolving power, supercritical fluids are playing a significant role in purification and extraction processes of various industries ( Kiran et al,2000; McHardy and Sawan, 1998; Wang et al, 2009). Over 50 years back in their seminal paper, M.E.Fisher and B.Widom discussed liquid and gaslike supercritical states by observing the decay behaviour of pair correlation function at large distances using linear continuum and lattice models and challenged the existing description of supercritical fluid as a single homogeneous phase like other states of matter (Fisher and Widom,1969). Since then, many experimental studies had been executed to validate the heterogeneous nature of the supercritical fluids (Nishikawa and Tanaka,1995; Morita et al, 1997; Gorelli et al, 2006; Simeoni et al, 2010) . Simeoni et al. (2010) carried out inelastic X-ray scattering and molecular dynamics simulation to find out a demarcation line between two dynamically different regime ("liquidlike" and "gaslike") in supercritical fluid around critical pointcalled Widom line. In a recent study, however, it hasbeen found that this Widom line doesn't obey the corresponding states principle and the transition lines differ with different fluids ( Banuti et al, 2017; Kanka and Krishnamurthy, 2018). Turbulence is the primary example of a highly nonlinear phenomenon. However, thereis evidence that some processes of shear turbulence are controlled by linear dynamics,in particular the mechanism by which energy is transferred from the mean velocitycomponent of the flow to the spatially and temporally evolving perturbations ( Farrell and Ioannou, P, 1998; Kim and Lim, 2000; Jimenez, 2013) . It is agreed that the stream – wise rolls and streaks are ubiquitous in wall – shear flow (Klebanoff, et al 1962; Kline et al 1967) and that they are involved in a quasi – periodic regeneration cycle (Panton, 2001; Adrian, 2007; Smits, et al, 2011; Jim ´enez, 2012; Jimenez, 2018). The space – time structure of rolls and streaks is believed to play an important role insustaining and carrying shear – driven turbulence (Smits, et al, 2011; Kim et al, 1971; Jim´enez and Moin, 1991; Hamilton et al, 1995). The ultimate cause maintaining this self – sustaining cycle, and hence turbulence, is the energyextraction from the flow mean shear. Within the fluid mechanics community, there havebeen several mechanisms proposed as plausible scenarios for how this energy extractionoccurs. Conceptually, we can divide these mechanisms into three categories: (i) modal inflectionalinstability of the mean cross-flow, (ii) non-modal transient growth, and (iii) non – modaltransient growth assisted by parametric instability of the time – varying mean cross – flow. The goal of this present work is to investigate the mechanism using theoretical modeling techniques on resolving turbulent heat flows in a non – dimensional circular tube that mounted obstacle in wall Y+. 2.0MODEL AND METHODOLOGY This section presents the methodology used in the study of heat transfer mechanism with supercritical carbon dioxide as the working fluid. The experimental data of Kim et al. (2005) is used as input data for both inlet and outlet boundary conditions for modeling in the STAR-CCM+ CFD code. 2.1STAR-CCM+ CFD code The phenomenon of heat transfer to supercritical fluid is not fully understood and the Simulation of Turbulent flow in Arbitrary Regions Computational Continuum Mechanics C ++ for Computational Fluid Dynamics (STAR-CCM+ CFD)based on code has not been fully explored in this area of science. For this reason an in-depth analysis is necessary to understand this phenomenon so as to perform accurate heat transfer predictions for use in designing this research. The STAR – CCM+CFD code is a “general multipurpose” commercial software for computational fluid dynamics. A software STAR – design as geometric modeling and grid generation tool to create geometry and surface/volume mesh. The code is characterized by the full interactivity between user and solution: the integrated analysis and visualization tool provide a live feedback on the progress of the simulation. This allows changing parameters without stopping the solution, observing directly the consequent effects. STAR – CCM+CFD code is an entire computational physics and engineering process for solvingproblems involving flow heat transfer and stress. It provides a suite of integrated components that combine to produce a powerful package that can address a wide variety of modeling needs. It is specifically designed to handle large models quickly and efficiently using a unique client–server architecture that seamlessly meshes and simultaneously solves and post – processes over multiple computing resources without requiring additional effort from the user. The object – oriented nature of the code can be seen in the user interface. An object tree is provided for each live simulation,(Enaroseha et al, 2021) containing object representations of all the data associated with the simulation. The objects presented on the simulation tree reside on the server, which can run as either a serial or a parallel process (Enaroseha and Omehe, 2020). 2.2Governing Equations in the Axi – symmetric Flow STAR-CCM+ This research work concerns axi – symmetric flows in circular tubes therefore a brief description of the governing equation in axi – symmetric form are presented below. The adopted notations in the equations are U represent the axial velocity and V represent the radial velocity; x represent the axial coordinate and r is the radial coordinate, ρ is thedensity, P the pressure, µ is the molecular viscosity, eis the effective dynamic viscosity and g the gravity. 2.2.1Continuity equation Fluid density is kept under the sign of differentiation to take into account its variation withpressure and temperature. In general, variation of fluid properties is accounted forbyemploying formulas in a separate module. All the flow variables appearing in the continuity equation and in all the other transport equations are considered as mean flow, time- (or ensemble-) averaged ones, according to the classical Reynolds-averaging approach. 10U U t x r r        (2.1) WhereU represents the axial velocity, ρ is density, r is the radial distance, x is the axial coordinates and t is the time.In the momentum equations, the Reynolds stresses that result from the process of averaging of the instantaneous Navier – Stokes equations are modelled using the Boussinesq approximation, i.e., assuming that they are proportional to mean rates of deformation, with the turbulent viscosity t as the proportionality factor. 2.2.2 Axial velocity (U) – The turbulent momentum equation. The Reynolds stresses that result from the process of the instantaneous Navier – Stokes equations are modeled using Bassinesq approximation, i.e., assuming that they are proportional to mean rates of deformation, with turbulent viscosity t as the proportionality factor, the equations 2.2 and 2.3 represent the turbulent momentum equations: 2 21 e e xU U prUV r t x x r r r rU g U                            (2.2) Where eis the effective dynamic, r is the radial distance, x is the axial coordinates, is the density, U is the axial velocity, xg is the gravity acceleration and V is the radial velocity. 2.2.3 V- Momentum equation  2 21 1 2 0e e e e eVV V V P U V UV r r r t x x r r r r x r r r rV r                                                 (2.3) Where eis the effective dynamic, r is the radial distance, x is the axial coordinates, is the density, U is the axial velocity, and V is the radial velocity. 2.2.4 Energy Balance equation Supercritical fluid properties vary dramatically as a function of temperature in the vicinity of the pseudo-critical temperature. These sharp variations can give rise to numerical instabilities that make convergence difficult. ∂ ∂x(pUT)+1 r(prVT)=∂ ∂x[(μ pr+μT σT)σT ∂x]+1 r∂ ∂r[r (μ pr+μT σT)∂T ∂r] (2.4) Where U is the axial velocity, T is the temperature, is the molecular viscosity, r is the radial velocity, V is the radial velocity, x is the axial coordinate, and  is the density. 2.2.5 The -equation: This is given by Durbin, (1986).  1t t bUk rVk rt x x r r r GP                                                  (2.5) Where kp is turbulent production term, х is the axial coordinate, μ is the molecular viscosity, bG is the production term, Ɛ is the energy dissipation rate and ρ is the density. The turbulent production term is calculated from: 2 2 2 2 2 2 2U V V U V x r r r xp                                     (2.6) In the equation above, U is the axial velocity, V is the radial velocity, r is the radial distance, kpis the turbulent production term, х is the axial coordinate, μ is the molecular viscosity, bG is the production term, Ɛ is the energy dissipation rate and ρ is the density. 2.2.6The-equation The - equation is the energy dissipation rate and it is the most widely used in the simulations for supercritical water due to its good performance in the free-shear layer flows. However damping functions has to be employed in the modelling of STAR – CCM+ CFD code in order to avoid inherent singular defects in the near wall region.    2 3221t t bU rV rt x x r r y I Ef p fc c G c                                                          (2.7) In the above equation,  is the kinetic energy,  is the energy dissipation rate, C is a constant, U is the velocity, μ is the molecular viscosity, r is the radial distance, V is the radial velocity, x represent the axial coordinate, y is the radial coordinate, ρ is the density, fis damping function, T is the temperature, bGis the production term, kpis the turbulent production term and g is the gravitational acceleration, the subscript e is the effective dynamic viscosity, and Ε is the enthalpy. 2.3 Turbulence Modelling Computational modelling play a vital role in improving the understanding of the physics of convective heat transfer to fluids at supercritical pressures and assist with the development of correlation for engineering application. Raynolds-Averaged Navier-Stokes (RANS) approach which uses time average quantities and has the advantage to shorten computing time was used. The two types ofturbulence models used in this work were the - and the -ω turbulence models. The Ɛ isthe energy dissipation rate of k, while the ω is the specific heat dissipation rate of  proposed by Wilcox in 1998 but put into use by Manceau and Hanjalic, (2002) . Six turbulence models were selected for this work to predict the heat transfer deterioration and enhancement, namely; -epsilon Abe-Kondoh-Naganoturbulence model (AKN), Elliptic blending turbulence model (EB), V2F -epsilon turbulence model, standard low-Reynolds turbulence model, -omega Shear-Stress Transport turbulence model (SST) and the Standard Wilcox -omega turbulence model: The basic transport equations of the six different turbulence models, their  terms and  terms as well as their constants are discussed below 2.3.1 V2F - turbulence model The V2F - turbulence model ( Abe et al, 1994; Davidson et al, 2003; Lien et al,1998) is known to capture near - wall turbulence effects, which are crucial for accurate prediction of heat transfer, skin friction and flow separation. The model solves two additional turbulence quantities which are the normal stress function and the elliptic function in addition to and. The model is designed to handle wall effects in turbulent boundary layers and to accommodate non-local effects. 3.3.1.1 Turbulence Equation for V2F turbulence model The  term for V2F model is given by  0k b k Mv tgAkkdVt k k v daG G S v                           (2.8) The  term for the V2F model is given by  1 220 1v tgAka dVt T v daGS v                                      (2.9) Where Gk is turbulent production term, M= dilation dissipation, t= turbulent viscosityG bis Production term,= kinetic energy ,= dissipation rateρ= density s And s are the user-specified source term, 0 is the ambient turbulence value in the source terms that counteract turbulence decay. 3.3.2 Standard low-Reynolds - turbulence model This model has identical coefficients to the Standard - model, but provides more damping functions. These functions let it be applied in the viscous-affected regions near walls. It involves transport equation for turbulent kinetic energy and its dissipation rate Ɛ. The transport equations are suggested by Jones and Launder (1972) with the coefficients suggested by Launder and Sharma, (1972). Some additional terms have been added to the model in STAR-CCM+ to account for effects such as buoyancy and compressibility. 3.3.2.1 Turbulence Equation for standard low-Reynolds turbulence model The  term for standard low-Reynolds is given by  0k b k Mv tgAkkdVt k k v daG G S v                          (2.10) Where ks= source term,  = turbulent kinetic energy, ρ= density, = velocity, g= velocity of gas, = rate of dissipation of , 0= ambient turbulence value, = viscosity, t= turbulent viscosity, k= turbulent kinetic energy, bG= production term, kG= kinetic energy production term, a= change in area, V = change in volume, 1f,2f and 3f are damping functions.The  term for standard low-Reynolds is given by    1 1 3 2 021 k n bv tgAkkG pY dV t T v daf C G G C G C S v                                   (2.11) WherekSand sare the user-specified source term, 0is the ambient turbulence value inthe source term that counteract turbulence decay and ' Gis an additional term. 22 22 222 exp Re =1-Cexp Re Re Re k y T y tk ky kG Df G E fy v v        2.4Physics Models The physics models define the primary variables of the simulation, including pressure, temperature and velocity and the mathematical formulation. In this research the flow is turbulent and compressive. The segregated flow models together with the default k – epsilon turbulent model were used. In STAR – CCM+ the physics models are defined on a physics continuum. The data in Table 2.1 were used for the simulation. Table 2.1: Fluid physics and model specification Models Model Specification Space model 2-D Time model Steady state Material medium Single phase carbon dioxide Equation of state Polynomial density Flow model Segregated flow Energy model Segregated fluid temperature Viscous regime Turbulent Turbulent model  – epsilon -epsilon model Realizable two layer Wall function Low y+ wall treatment 3.0Simulation The models described in section 2 above was run with a constant mass flow rate of 314 kg/m2s, inlet fluid temperature of 15 OC, pressure of 8 MPa and heat fluxes of 20, 23, 30 and 40 kW/m2. The simulation was run in a steady-state condition with segregated flow model and - turbulent model as default, the material medium was single phase carbon dioxide (gas) with segregated fluid temperature. The low y+ wall treatment was chosen because it is suitable only for low Reynolds turbulence models in which it is assumed that the viscous sub-layer is properly resolved. The simulation was run for 20,000 iteration as the maximum step or stopping criteria, it was observed in the residual as it converged to a reasonable level. 4.0 PRESENTATION AND DISCUSSION OF RESULTS 4.1 Fluid Velocity within the Circular Tube Figure 4.1 shows the velocity profile of the circular tube at heat flux of 23 kW/m2 and mass flux of 314 kg/m2s at a distance of 1.4 m along the tube. From the Figure 4.2, it can be observed that the various turbulence models predict the radial velocity profile. From Figure 4.1 we see that the velocity obtained from the simulation initially has a value around 0.38 m/s. As the simulation continues, the velocity increases steadily as it gets to the walls of the tube. The maximum of the velocity is achieved just before the prism layer of the tube was created. The profile from here shows a sharp downward trend. This downward trend shows that there truly is a laminar sub-layer just beneath the walls and mimicked by the prism layer created in the physical model. In this laminar sub-layer, it shows that fluid movement is very slow and can be described as the molecules of the CO 2 coming to a stop. Figure 4.1: Radial velocity at a distance of 1.4 m with heat flux of 23 kW/m2 It is worth noting that even though the various -ɛ models predicted the same shape of velocity, there seems to be some differences in the velocities obtained by the different models. Using the AKN model data at flux of 20 kW/m2, we obtained the Velocity profile details of some selected radial distances with different heat fluxes, these were presented in Table 1 of Appendix I. The results obtained from Table 1 in Appendix I shows that the value of the Position [0.0, 1.0, 0.0] (m) – 0.8m has a constant value of 0.0039 even with the continuous simulation but the Velocity Component [0 – 0.8 m (m/s)] decreases from 0.4337 to 0.0988. From Table 2 of Appendix I, the AKN model data at heat flux of 23 kW/m2 were used in the simulation. The result shows that the value of the Position [0.0, 1.0, 0.0] (m) – 0.8m was 0.0039 even with the continuous simulation but the Velocity Component [0 – 0.8 m (m/s)] increases from 0.0008 to 0.0594 with value of the position value remaining constant. 4.2Fluid Density within the Circular Tube Again only the -ɛ models used in this research would be considered since they gave appreciable results when compared with the wall temperatures obtained from the experimental data. Figure 4.2 shows the radial density profile of the flow with supercritical CO 2 at an axial distance of 1.4 m with a heat flux of 23 kW/m2. The -ɛ models of AKN, EB, Standard Low-Reynolds and V2F describes the density of supercritical CO 2 at 8 MPa at a radial distance of 1.4 m. from the graph we also observed that the supercritical fluid at the middle of the tube shows a relatively constant density. As the fluid come close to the walls of the tube, the density begins to drop very sharply. It can be said that its drop is almost instantaneous showing that the fluid has reached it pseudo-critical point. This is the region which has proven challenging for current correlations in the field of computational Physics engineering to deal with. Figure 4.2: The density at a distance of 1.4 m with heat flux of 23 kW/m2 Density profile details of some selected radial distances with different heat fluxes are also presented in Appendix II. The AKN model simulation result data for density at heat flux of 20 kW/m2 as shown in Table 1 of Appendix II shows that the value of the Position [0.0, 1.0, 0.0] (m) – 0.8m has a constant value of 0.0039 with Radial density data at different axial distances and heat fluxes that ranges from 868.4071 as the highest density and 760.3031 as lowest density. Table 2 of Appendix II shows the AKN model data for density at heat flux of 30 kW/m2, the results shows 868.4073 and 149.0942 as the highest and lowest densities respectively with a constant position of 0.0039 4.3 Y+ Wall Treatment The wall Y+is a non-dimensional distance and is often used in CFD to describe how coarse or fine a mesh is for a particular flow. It is a ratio between turbulent and laminar influence in a cell. The low Y+wall treatment is suitable only for low Reynolds turbulence models in which it is assumed that the viscous sub-layer is properly resolved. Turbulence flows are significantly affected by the presence of walls, where the viscous-affected regions have large gradients in the solution variable and accurate presentation of the near-wall determines successful prediction of wall bounded turbulence flow. Very close to the wall, viscous damping reduces the tangential velocity fluctuations, while kinematic blocking reduces the normal fluctuations. Towards the outer part of the near-wall region, however, the turbulence is rapidly augmented by the production of turbulent kinetic energy due to the large gradients in mean velocity . Figure 4.3 shows the heat fluxes of 20, 23, 30 and 40 kW/m2, it indicates that the heat fluxes increases as the Y+ wall profile moves away from the near – wall region, this is due to the effect of viscosity, buoyance, acceleration and the friction of the turbulence modification. Figure 4.3: Y+wall treatment at different heat fluxes. 5.0 CONCLUSION AND SUMMARY Analysis of heat transfer in circular tube with super critical CO 2 has been carried out using STAR- CCM+. The results of the simulation has been found to be able to reproduce the general features exhibited in the experimental data, The wall Y+ treatment simulation was performed for various heat fluxes (20, 23, 30, 40 kW/m2) and a constant mass flux of 314 kg/m2s. It was noted that as the heat flux was increased with constant mass flux, the predictability of the models was better in the case of the -ɛ models. The models were able to improve their predictability in terms of where the deterioration was observed experimentally. Results from the various simulations show that there is a laminar sub-layer or what is termed as a buffer, present. The effect of this buffer can be seen in the plot of velocity against the radial distance of the circular tube. Based on the results of simulations, it is noted that at low heat flux, the effect of deterioration is still evident in the -ɛ models even though this phenomena is not experienced in the experiment. This was attributed to the heat transfer correlation used in the development of the code. At high heat fluxes the deterioration observed in experiment was also evident during the simulation even though it was over predicted by the models used. It can be said that the effect of deterioration has been effectively modeled and accounted for. The effect of heat transfer enhancement was neither seen in any of the cases studied during the research nor was it seen in the experimental case used as the benchmark experiment for this research.