effect of turbulence models on predicting convective heat transfer to hydrocarbon fuel at supercritical pressure

Tao Zhi,Cheng Zeyuan,Zhu Jianqin,Li Haiwang

National Key Laboratory of Science and Technology on Aero-Engine Aero-Thermodynamics,School of Energy and Power Engineering,Beihang University,Beijing 100083,China

effect of turbulence models on predicting convective heat transfer to hydrocarbon fuel at supercritical pressure

Tao Zhi,Cheng Zeyuan,Zhu Jianqin*,Li Haiwang

National Key Laboratory of Science and Technology on Aero-Engine Aero-Thermodynamics,School of Energy and Power Engineering,Beihang University,Beijing 100083,China

A variety of turbulence models were used to perform numerical simulations of heat transfer for hydrocarbon fuel flowing upward and downward through uniformly heated vertical pipes at supercritical pressure.Inlet temperatures varied from 373 K to 663 K,with heat flux ranging from 300 kW/m2to 550 kW/m2.Comparative analyses between predicted and experimental results were used to evaluate the ability of turbulence models to respond to variable thermophysical properties of hydrocarbon fuel at supercritical pressure.It was found that the prediction performance of turbulence models is mainly determined by the damping function,which enables the m to respond differently to local flow conditions.Although prediction accuracy for experimental results varied from condition to condition,the shear stress transport(SST)and launder and sharma models performed better than all other models used in the study.For very small buoyancy-influenced runs,the thermal-induced acceleration due to variations in density lead to the impairment of heat transfer occurring in the vicinity of pseudo-critical points,and heat transfer was enhanced at higher temperatures through the combined action of four thermophysical properties:density,viscosity,thermal conductivity and specific heat.For very large buoyancyinfluenced runs,the thermal-induced acceleration effect was over predicted by the LS and AB models.

1.Introduction

Supercritical fluids have been used in many fields since critical phenomena was observed by Thomas Andrews in 1869.1Regular fossil fuel power plants operated by using supercritical water to drive steam turbines in the 1960s and 1970s.In the brief period of time since the n,the development of emerging nuclear power stations such as supercritical water reactor (SCWR),2air-conditioners3and cooled cooling air (CCA) technology for advanced aircraft engines4has stimulated a renewed interest in the flow and heat transfer of supercritical fluids.

Particular focus has been placed upon its use in regenerative cooling technology to deal with heat management problems that occur due to aerodynamic heating of scramjet engines.5,6This cooling operates by the heat of the combustor first being absorbed by endothermic hydrocarbon fuel flowing through a cooling channel at supercritical pressure and the n released as the fuel is injected from nozzles into the combustor.7The most striking feature of supercritical fluids is that the y have no phase change under supercritical pressure conditions but undergo a dramatic alteration of physical properties at temperatures in the vicinity of critical point.The nonlinear relationship between thermal-physical properties and temperature leads to strong coupling between velocity fields and scalar functions.This may arouse distinct distortions of velocity fields,especially under the influence of buoyancy force and thermal-induced acceleration.Hence very significant phenomena,different from heat transfer under subcritical conditions,can be seen in supercritical heat transfer.

There is a wide variety of reviews about experimental studies and numerical analysis on heat transfer at supercritical pressure.8–11Moreover,there are some surveys on hydraulic resistance of fluids at supercritical pressure.12–15It should be noted that experimental data obtained in literature have mainly been limited to that related to heat transfer,such as wall temperature.With detailed predictive information lacking,numerical simulation might provide for what cannot be obtained from experimentation and help to give a better understanding of heat transfer mechanisms.Turbulence modeling has played an important role in the mathematical modeling for supercritical heat transfer using computational fluid dynamics (CFD) methods.Early studies usually used mixed length models to consider variations in physical properties.16Afterwards,more complicated models including one equation and two-equation models were used.17–21In recent years,in order to evaluate the performance of a range of turbulence models in terms of simulating heat transfer to fluids at supercritical pressure,studies comparing turbulence models adopted by CFD tools for numerical simulation have been performed.22–28It has been found that low-Reynolds number turbulence models can reproduce experimental data qualitatively in the heat transfer enhanced and deterioration regimes though vary obviously in terms of quantitative prediction.

Many comparison studies of turbulence models assisting in simulated heat transfer at supercritical pressure have been performed with water and CO2as working fluids,but very few have been done that study hydrocarbon fuels as working fluids.Compared to pure substances,hydrocarbon fuels show more complicated behaviors that result from their being composed of thousands of components that may undergo cracking and coking reactions as temperature increases.29,30Zhu et al.31have applied the SST k-ω two-equation turbulence model with enhanced wall treatment to the numerical simulation for the thermal cracking of n-decane flowing at supercritical pressure,while Goel and Boehman32used the laminar model in simulation of jet fuel degradation in flow reactors.To this author's knowledge,there is no detailed study on turbulence modelling of hydrocarbon fuels under supercritical pressure.

Numerical computations with commercial software(Fluent)have been performed that simulate the flow and heat transfer of hydrocarbon fuel at supercritical pressure,and a comparison study of computational results using a range of turbulence models with experimental data is presented in this paper.This is done in order to assess the ability of turbulence models to predict supercritical heat transfer of hydrocarbon fuel with a special focus on the response to variable thermophysical properties under supercritical pressure.

2.Simulated experiment

This study simulated and investigated the experiments by Zhang33where the characteristics of local convection heat transfer of a supercritical hydrocarbon fuel called RP-3 flowed downward through a vertical mini-tube.The test section was made of stainless steel with an inner diameter of 1.8 mm and wall thickness of 0.2 mm,comprising a 90 mm long adiabatic entrance section,and a uniformly-heated section length of 300 mm followed by an adiabatic exit section 60 mm long.The temperature,pressure and mass flow rate were measured by a K-type thermocouple welded onto the test tube;the pressure sensor mounted downstream of the test tube and mass flowmeter installed upstream of the test tube,respectively.Detailed information about the experiment can be found in Zhang.33Due to there being no obvious heat transfer deterioration data on RP-3,the experimental results of Liu34are selected,which study the heat transfer deterioration of ndecane at supercritical pressures flowing upward and downward in small vertical tubes.The test section,with an inner tube diameter of 2 mm and the wall thickness of 0.5 mm,was a 959 mm vertical tube including two 100 mm adiabatic sections before and after the heated section.

Table 1 shows the experimental conditions considered in present paper.At P=5 MPa and Tin=373 K,the deterioration of heat transfer could be found at the beginning of the heating section.At P=4 MPa,as Tinincreased,the process of heat transfer has four experimental regimes:the initial heating section,normal heat transfer,heat transfer enhancement,and heat transfer deterioration.At P=3 MPa,the buoyancy effect is very strong when flowing upward,but very weak when flowing downward.When computing all dimensionless numbers,the characteristic length and characteristic temperature are the inner diameter of tube and the average temperature of fluids in across section of the tube,respectively.The calculation of thermal physical properties of hydrocarbon fuel is described in Section 3.6.

3.Numerical modelling

The commercial CFD software,ANSYS Fluent 14.5,was adopted for numerical simulation in this paper.The following sections introduce the governing equations,turbulence modelling,boundary conditions,solution methods and the fluid properties to model supercritical convection heat transfer in a vertical tube,in sequence.

Table 1 Experimental conditions.

3.1.Governing equations

The governing equations employed are obtained via averaging of instantaneous equations over time,including mass conservation,momentum conservation and energy conservation,along with turbulence transport equations.It's worth mentioning that the irregular fluctuations over time have been included in the time-averaged form of pressure,temperature and velocity,excluding for thermal physical properties.The transient terms in all the governing equations have been removed,since only the steady-state processes of flow and heat transfer at supercritical pressure are considered in this study.The tensor forms of governing equations for mean flow and energy transport are:

Continuity equation:

where ρ is the density,uiis the velocity tensor and xiis the coordinate tensor.

Momentum equation:

in which

where P is the pressure,δijis the kronecker symbol,and μ is the laminar viscosity.

Reynolds stress based on the Boussinesq hypothesis is shown below:

Energy equation:

in which turbulent heat flux is modelled by using Boussinesq hypo the sis and Prtis the turbulent Prandtl number taken as 0.85 in the present study.k is the turbulent kinetic energy,μtthe turbulent viscosity,H the enthalpy,λ the thermal conductivity,˙Q the inner heat source and φ the dissipation term.

During the whole physical process there is no inner heat source if the pyrolysis reaction of fuel is not taken into consideration when the temperature of fluids is below the onset temperature of pyrolysis.According to Mohseni and Bazargen's derivation process,23the pressure and viscosity terms in the energy equation could be neglected in variable property flow.So,the energy equation could be simplified into the following form:

The convection heat transfer coefficient is calculated from the expression below:

where qwis the wall heat flux,Twis the wall temperature and Tbis the bulk temperature.

This is used to compute the Nusselt number:

where d is the tube internal diameter.

3.2.Turbulence modelling

Previous studies reported that low-Reynolds number k-ε type turbulence models respond well to the local conditions considered in the present study,and some recent models have been found to perform better in predicting supercritical heat transfer than some LRN models.Based on such recognition,five LRN models,one k-ω type model plus one four equation model were selected.The LRN models are those of Abid (ABID),35Abe,Kondoh and Nagano (AKN),36Lam and Bremhorst (LB),37Launder and Sharma (LS),38Yang and Shih (YS);39the k-ω type model is shear stress transport (SST) model from Menter40and one four equation model is the k-ε-υ2-f model(V2F)developed by Behnia et al.41

The principal differences between all those models are:(1)the damping effect on turbulence due to the wall;(2)the construction of additional terms in turbulent kinetic energy and the rate of turbulent dissipation transport equations;(3)the definition of boundary conditions on the wall in turbulent kinetic energy and rate of turbulent dissipation transport equations.

The constitutive equation and the generic form of transport equation for k-ε type models are shown below:

Constitutive equation:

Turbulent kinetic energy:

The rate of turbulent dissipation:

where Cμ,Cε1,Cε2, σkand σεare Constants,fμ,f1and f2are functions,ε is the Dissipation rate of k,Skand Sεare the source terms in k equation and ε equation,respectively.Gkand Pkare the turbulent production due to buoyancy and turbulent production due to shear stress,respectively.

The term Pkcan be computed by:

The buoyancy effect on flow and heat transfer can be divided into two types:the indirect(external)effect and the direct(structural)effect.42The indirect effect reflects the modification of the mean flow fields on turbulence by the body-force term in the momentum equation.The buoyancyinduced production of turbulent kinetic energy is regarded as the direct effect of the buoyancy on flow and heat transfer.In this study,the Simple Gradient Diffusion Hypothe sis (SGDH)28is adopted for modelling the direct effect of the buoyancy on turbulence,as expressed below:

Sharabi et al.24applied the SGDH to model the direct effect of buoyancy on turbulence and saw good performance,indicating the capacity of SGDH for assisting numerical simulation at supercritical pressure.

The SST model has a similar form to the k-ε type model,as shown below:

Turbulent kinetic energy:

Specific dissipation rate(identified as the ratio of ε to k):

Detailed descriptions of the constants and the coefficients can be found in the reference by Menter.40

The V2F model has two additional equation for υ2and f,as expressed below:

Turbulent velocity scale(υ2):

Production(f):

Turbulent viscosity μthas another definition aside from Eq.(9)as shown below.

where

In the V2F model,the definition of constants is shown below:

The summary of model constants,the function,additional terms and boundary conditions on the wall are listed in Tables 2–4.

Table 2 Summary of constants in the turbulence models.

Table 3 Summary of functions in the turbulence models.

Table 4 Summary of additional terms and boundary conditions on the wall in the turbulence models.

3.3.Computational mesh

A 2D structured mesh is used to discretize all of the computational domain(Fig.1),covering the adiabatic and heated sections along length direction,ranging from the fluid domain to the solid wall.The grid independent solution is inspected by trials and comparison in each run(Fig.2).At each trial-anderror process,the mesh is carefully adjusted to ensure to satisfy a basic principle that more nodes be placed near the inner wall and the heating region.The optimal mesh of grids,typically,(10+60)X(20+100+15), representing (radial solid + radial fluid)×(axial unheated + axial heated + axial unheated),is obtained in which the distance between two adjacent nodes grows by a factor of 1.05 with the increase of distance from the inner wall in the fluid domain.Meanwhile,the optimal mesh could ensure that the first node close to the walls has y+lt;0.5 in all runs in order to compute from the viscous sub-layer.

3.4.Boundary conditions

Based on the axisymmetric geometry and flow conditions present in this study,a two-dimensional cylindrical-coordinate system was employed to solve the governing equation,and an actual physical model was created for mathematical computation with half of the symmetric longitudinal section of test tube aimed at reducing total calculations.

Fig.1 Schematic of computational domain.

Fig.2 Grid independent solution.

The mass-flow inlet type was chosen for the inlet boundary condition,where the mass flow rate and total temperature are specified,and the pressure at the outlet of pipe is given for outlet boundary condition,i.e.pressure-outlet type.On the surface of pipe,the constant heat flux imposed on the uniformly-heated section served as wall boundary condition,meanwhile the heat flux is zero on the adiabatic wall of pipe.The symmetry axis of pipe uses the axis type as a boundary condition,implying that any quantities ϖ existing in one grid must meet the following expression:

3.5.Solution method

The pressure-based segregated algorithm was adopted to solve the discretization equations,and the SIMPLE scheme was used as the pressure-velocity coupling method.The Least Squares Cell-Based Gradient Evaluation was applied for modelling the gradient of the convection and diffusion terms in the flow conservation equations,and the STANDARD scheme was used to interpolate the pressure values at the face.The second order upwind scheme was used for computing the face values in the momentum and energy equations for the improvement of accuracy,and the first order upwind scheme was used for the turbulence equations to enhance the robustness.The absolute convergence criterion of residual for all of the equations was set to be less than 10-6.

The gravity option was enabled to account for buoyancy,and a transient calculation was performed to model buoyancy-driven flow on the premise that density could be obtained at known pressure and temperature.The calculation of fluid properties is described below.

3.6.Fluid properties

The fluids selected include RP-3 and n-decane,which is one of the components of RP-3,so the following description is based on the RP-3.The endothermic hydrocarbon fuel RP-3 includes alk ane,cyclanes,olefin and aromatic hydrocarbons,and many other substances.It is impossible to analyze every component in actual study,so surrogate fuels with a mixture of typical hydrocarbon compounds were adopted to simulate the thermal physical properties of the real fuel.

There are very few studies on the surrogate fuels for RP-3 from around the world,but Fan et al.in China have proposed two surrogate models for RP-3:the three-species surrogate model43and the ten-species surrogate model.44The numerical simulation code has been developed to compute the thermodynamic and transport properties of surrogate fuels by using extended corresponding state laws.A four-species surrogate model,consisting of 19.1% n-decane,36.5% n-dodecane,29.9% n-butylbenzene and 14.5% methyl-cyclohexane,has been found to give the best predictions for the experimental data from Zhang et al.,45–47maintaining the consistency of experimental data.Further,the predictions of the fourspecies surrogate model by this author's method agree well with those by SUPERTRAPP.48

Fig.3 shows the density,specific heat,viscosity and thermal conductivity of the four-species surrogate model computed using in-house code varying with temperature under different pressures.As can be seen from Fig.3,all of the thermophysical properties undergo a dramatic variation with increase of temperature at constant pressure in the vicinity of pseudo-critical temperature in which the specific heat is at maximum.The significant physical properties that vary with temperature at constant pressure near pseudo-critical temperature can weaken,or even disappear,when the pressure increases.

Considering the pressure drop through the pipe is very small,the thermal physical properties,as a function of temperature,are added to a User-Defined Database in a piecewiselinear manner.At a pressure of 4 MPa,the computed pseudo-critical temperature is 680 K;at the pressure of 5 MPa,the computed pseudo-critical temperature is 698 K.

4.Results and discussion

4.1.Different wall heat fluxes at 5 MPa

4.1.1.Wall temperature and Nusselt number

As mentioned earlier,the value of wall temperature can be measured conveniently from the experiments,so this variable can be compare results predicted by turbulence models with those from experimentation.In addition,it is a common practice to examine the Nu/Nufvarying with axial distance.Nufis calculated using the Dittus-Boelter correlation while considering the variations of the physical properties shown below49:

where Re is the Reynolds number and Pr is the molecular Prandtl number.ρwis the density based on the wall temperature,and ρbis the density based on the bulk temperature.Cpbis the specific heat based on the bulk temperature,andCpis computed as shown below:

Fig.3 Thermophysical properties of the four-species surrogate model varying with temperature under different pressures.

and n is obtained below:

Fig.4 shows the prediction of the inner wall temperature and the Nusselt number varying with the axial relative distance x/d using various turbulence models,together with the measured inner wall temperature for Runs 1–4.The axial relative distance x/d is the ratio of the distance from the inlet of heating section to the internal diameter of pipe.The experimental conditions for Runs 1–4 are under the same pressure of 5 MPa with different wall heat flux.

For Run 1,the experimental results show that near the start of heated section,heat transfer deterioration occurs,i.e.the inner wall temperature grows rapidly,the n,at x/d=40,the growth rate reduces and wall temperature increases steadily due to constant flux applied to the wall of the pipe.It can be seen from Fig.4(a)that all of the turbulence models were able to reproduce the general trend of the variation of inner wall temperature along the tube measured in the experiments qualitatively,but differ significantly in quantitative prediction.The LS model under-predicts the wall temperature while the AB,LB,LS,YS and AKN models over-predict it.A closer look shows that the SST and LS model produce a shorter distance for heat transfer deterioration than that exhibited in the experimental results,but the other models give a more accurate prediction.

As seen in the Fig.4(b),the Nu/Nufcalculated by most of the models decreased first and the n have an approximate constant as it moves downstream.The LS and SST model show greater values of Nu/Nufthan those of experiments,leading to lower wall temperature predictions by the LS and SST model when compared with the measured experimental values shown in Fig.4(a).All turbulence models(except the LS model)and the experimental results show the Nusselt number lower than that calculated using the Dittus-Boelter correlation considering the variations of the physical properties throughout the whole pipe.The Dittus-Boelter correlation used currently does not include the effects of buoyancy and thermal-induced acceleration,so the buoyancy and/or the thermal-induced acceleration causes the values of Nu/Nufless than unity in most computations,including the experimental data,based on the current heat transfer correlation.

The similar behavior of turbulence models in Run 2,where the buoyancy parameter is greater than that in Run 1,is presented in Fig.4(c)and 4(b).Compared with results in Fig.4(a)and 4(b),the main difference is that the AB and LB models predict the higher wall temperature and lower Nusselt number.

Fig.4 Wall temperature and Nusselt number along the tube for Run 1–4.

As shown in Fig.4(g),the wall temperature calculated using the AB model reaches beyond 1700 K,while the maximum wall temperature predicted by the LB model is about 900 K at x/d=60.The AB model over responds to the local varia-tion of physical properties in the current case and predicts the results of heat transfer deterioration.The poor performance of the AB model in predicting supercritical heat transfer can also be found in Refs.21,24It can be seen from Fig.4(f)that the AB and LS models predict a trend of variation for Nu/Nuf,which is inconsistent with the experimental data,and the predicted results of the V2F model agree very well with the experimental results of heat transfer characteristics for the first half of the heated section;better prediction precision for heat transfer characteristics is obtained by using the SST model for the second half of the heated section.

In Run 4,the heat flux of 550 kW/m2is loaded uniformly into the heated section,and the maximum value of the buoyancy parameter is obtained among all runs considered in present study.Fig.4(g)demonstrates that,in addition to the AB and LB models,the AKN model predicts a maximum value of wall temperature at the location of x/d=110,suggesting that the deterioration of heat transfer occurred,which is unobserved in the predicted results from YS,V2F,and SST models.According to the conclusion above,different models have different response speeds to the increase in heat flux,in other words,the increase of the buoyancy parameter.Kim et al.28divided turbulence models into two categories:turbulence models that can predict the deterioration of heat transfer due to the buoyancy well,including LS,CK,AKN,YS and V2F models,and those which are not able to do so,including CH,MK,HL and WI models.He et al.27also noted that the damping function of turbulence models that respond rather intensely to the local conditions were based on Ret=k2/ευ and those which did not respond to the local conditions were based on y+.Based on the results obtained in the Fig.4(g),the seven turbulence models can be grouped into two categories:Group 1,consisting of AB,LB and AKN models,which seriously over-predict the wall temperature,and Group 2,consisting of SST,V2F,YS and LS models,which predict the trends of wall temperature much better than Group 1.At first glance,there are discrepancies between current results and those obtained by Kim et al.,but this isn't the case.Though limited by this author's ability,some reasons are shown as below:(1)the evaluation criteria and working conditions for the current study are different from those used by Kim et al.The benchmark database of the study by Kim et al.was the DNS study conducted for turbulent mixed convection heat transfer to air flowing upward in a vertical tube with the assumption of ideal gas law at subcritical pressure and constant properties with the Boussinesq simplification.Kim et al.'s study examined the predictive ability for buoyancy-induced impairment of heat transfer,while there is no obvious deterioration of heat transfer due to the buoyancy occurring in the experiments conducted by Zhang,33which is analyzed later.(2)The conclusion that the turbulence models in which the damping functions have the Retas the feature variable respond to the local flow condition strongly,as obtained by He et al.,27is basically applicable to the cases in the present study.Group 1,in which the damping functions have the Retand Reyas the feature variables,respond to the local flow conditions too strongly,leading to the obvious deterioration of heat transfer though under very low buoyancy parameters.The YS and LS models where the damping function has the Retas the feature variable respond weaker to the local flow conditions than those in Group 1,so these two models are classified in Group 2,along with the SST and V2F models.

He et al.22found that a shorter distance is taken to approach a ''developedquot; state for heat transfer characteristic as the buoyancy parameter increases.However,there is no similar phenomenon as described by He et al.observed in present study;the turning point between the downward section and upward section of change curve of Nu/Nufoccurs in the location of x/d=40 in all runs,while the buoyancy parameter increases with the increase of heat flux applied in the heated section from Run 1 to Run 4.The comparative study of Figs.4(b),4(d),4(f)and 4(h)can reveal that the models in Group 1 predict very similar results,whereas the models in Group 2 predict different results in cases that have different buoyancy parameters,verifying the different response speeds to local conditions between Group 1 and Group 2.

Since the SST model performs relatively well compared to other models,the SST model is chosen for further analysis of the deterioration of heat transfer occurring in the start of heated section.

4.1.2.Detailed flow and turbulence fields

The radial distributions of temperature,density,axial-velocity and turbulent kinetic energy for Run 1 predicted by using the SST model at different axial locations along the pipe are shown in Fig.5.With the purpose of explaining the impairment of heat transfer in the start of heated section,four axial location are selected:x/d=5 representing the inlet of the heating section,x/d=10 representing the location during the impairment of heat transfer,x/d=30 representing the end of the impairment,x/d=80 in the normal region and x/d=120 also in the normal region but further downstream,as shown in Fig.4(a)and(b).

In Fig.5(a),at any cross section,the temperature of fluids varies dramatically in the neighborhood of the wall and remains basically unchanged in the rest region.Proceeding downstream,the temperature of fluids rises due to the heating effect by wall heat flux.As seen in Fig.5(b),the density varies rapidly near the wall,but within an order of magnitude,leading to a smaller local buoyancy effect.That,along with limited density change within the region very close to the wall,which can be considered as viscous sublayer,could not cause the distortion of velocity as shown in Fig.5(c).In addition,while proceeding downstream,the axial velocity sees only small increases,thus the influence of thermal-induced acceleration could be neglected.In Fig.5(d),at any cross section of pipe,the value of turbulent kinetic energy undergoes a sharp increase in the vicinity of the wall,reaches a peak and the n reduces linearly until it becomes stable at the centerline of the pipe.At x/d=5,the difference in temperature between the wall and the fluids very close to the wall is at maximum,so the effectiveness of heat transfer is fine,as seen from Fig.4(a)and(b),and a large peak value of turbulent kinetic energy can be observed.At x/d=10 and x/d=30 downstream,the turbulent kinetic energy is reduced,apparently,with the lower peak value,indicating that the turbulence weakens.On the other hand,the thermal boundary layer is at seedtime with increasing thickness and thermal conductivity decreasing with increasing temperature,leading to a large thermal resistance.Finally,the wall temperature increases and the Nusselt number decreases,owing to the above two canceling effects.Further downstream,at x/d=80 and x/d=120,the peak value of turbulent kinetic energy occurs in a location very close to the wall and returns to a higher level.At the same time,the thermal boundary layer continues to evolve from laminar flow to turbulent flow.So,the effectiveness of heat transfer gradually recovers.The impairment of heat transfer at the start of heated section may happen not only under supercritical pressure but also under subcritical pressure in this author's opinion;Hall considered the deterioration of heat transfer in the initial heated section the result of buoyancy force near the wall.50Differences in operating conditions and thermal physical properties are considered to be among the main reasons for such discrepancies.

Fig.5 Radial distribution at different axial locations along the pipefor Run 1.

4.2.Different inlet temperatures at 4 MPa

4.2.1.Wall temperature

Due to the lack of experimental data for the Nusselt number,only the comparison of wall temperatures between computational results and experimental data is conducted in current section.Fig.6 shows the prediction of the inner wall temperature varying with the axial relative distance x/d using various turbulence models together with the measured inner wall temperature obtained in the experiments for Runs 5–10.

In Fig.6(a),except for the AKN model,all other models exactly agree with the classification methods summarized in Section 4.1.1.In Group 2,the SST model predicts a change curve of wall temperature most close to the experimental curve among all models though,as expected,there is poor prediction accuracy.

Seen from Fig.6(b),beyond the initial heated section where heat transfer is impaired,all of the models used in current study show the orderly arrangement of wall temperature curves similar to the experimental curve,the difference being the distance of the impairment in the initial heated section and the initial temperature of the normal heat transfer region.It was also found that only the LS model predicts a lower temperature than the measured value.

As seen from Fig.6(c),for measured wall temperature,the initial heating section where the heat transfer is impaired disappears at the same time the gradient of the variation of measured wall temperature starts to reduce at the location where the wall temperature reaches to the pseudo-critical point.No turbulence models could predict such phenomenon qualitatively.As shown in Fig.6(d),as the inlet temperature is increases,the extent of excessive prediction of wall temperature increases,so the LS model predicts the wall temperature closest to the experimental result.The deterioration of heat transfer predicted by the AB,AKN,YS and V2F models occurs in different locations with varying degrees,possibly due to the buoyancy effect and thermal-induced acceleration,along with the variation in physical properties beyond the pseudo-critical temperature.

In Fig.6(e),except for the LS model,all other models predict the enhancement of heat transfer occurring in the vicinity of 790 K.The same phenomenon happens in Fig.6(f).The enhancement of heat transfer obtained by using the V2F,SST,AB,YS and AKN models ends near 800 K,the n the wall temperature increases basically along a straight line.The LS model responds very weak to variations of inlet temperature,indicating no obvious difference between cases.

Fig.6 Wall temperature along the tube for Runs 5–10.

Based on the above analyses,the prediction ability of the SST model is at an acceptable level for the flow and heat transfer of RP-3 both at the P=4 MPa and at P=5 MPa,and the subsequent analyses use the SST model to follow the above practice.

4.2.2.Buoyancy effect

It is easy to examine the buoyancy effect on the heat transfer of RP-3 by reversing the flow direction in the numerical simulation.Fig.7 gives the comparison between up flow and down flow of the Nusselt number varying with bulk temperature using the SST model for all runs at P=4 MPa.It is observed that the curve for down flow is in good agreement with the other curve for up flow,but shows slightly higher at the location of the pseudo-critical point if observed very closely.This subtle difference could be explained as follows:for buoyancy opposed-flow (downflow),the buoyancy strengthens the mixing flow,enhancing the turbulence.For buoyancy aided-flow (upflow),the buoyancy reduces the velocity gradient,flattening the velocity and,as a result,turbulence production is reduced and heat transfer is impaired.So,for the same bulk temperature,the effectiveness of heat transfer downflow is higher than that for upflow.

Fig.7 Nusselt number varying with the bulk temperature flowing both upward and downward.

Heat transfer is impaired at the location of pseudo-critical point and it is enhanced in the vicinity of 780 K(1.1 times of pseudo-critical temperature)for the bulk temperature.Due to buoyancy having little or no effect on the heat transfer in current study,these runs could be identified as forced convection flow under supercritical pressure.The following section gives further analysis of the enhancement and deterioration of heat transfer.

4.2.3.Further investigation of the deterioration and enhancement of heat transfer

Fig.8 gives the Nusselt number varying with the bulk temperature using different physical property models,including four models in which each of four thermophysical properties are fixed as constant at the inlet temperature in every run(CD,CT,CV and CS),along with the real physical property models(REAL)used in the above study.

For the same physical property model,the variation in the Nusselt number is discontinuous due to different thermophysical properties based the inlet temperature in every run.The results calculated using the CT and CS thermophysical models agree well with those calculated using the real model,indicating that thermal conductivity has little effect on the deterioration and enhancement of heat transfer occurring near Tpcand 1.1Tpc.The same can be said for specific heat.Although the results calculated using the CV the rmophysical model are able to show the trends of deterioration and enhancement of heat transfer,it underestimates the Nusselt number calculated using the real model in most of the temperature range.In the real physical property model,the viscosity decreases with the increasing temperature,which can enhance the heat transfer between the wall and the fluid,so the Nusselt number in the CV model is less than that in the real model,which is especially notable in the range below the pseudo-critical temperature where the viscosity changes drastically.When the calculation is conducted using the CD model,the results only shows the peak of the Nusselt number at the location ahead of 1.1Tpc,without the deterioration of heat transfer occurring in the vicinity of Tpc.This proves that the variation of density is the decisive factor among the four thermophysical parameters for the impairment of heat transfer near the Tpc.The thermalinduced acceleration due to the variation of density weakens the capability of heat transfer near Tpcand delays the onset of enhanced heat transfer.These four thermophysical properties vary with temperature according to the law shown in Fig.1,and reach an optimal combination for effectiveness of heat transfer at the temperature of 1.1Tpc.

Fig.8 Nusselt number varying with the bulk temperature using different physical property models.

4.3.Different inlet Reynolds and flow directions at 3 MPa

4.3.1.Wall temperature

The operating pressure is 3 MPa,very close to the critical pressure of the n-decane,2.11 MPa,resulting in the remarkable variations in physical properties,which may cause the buoyancy and thermal-induced acceleration effects.Fig.9 shows the prediction of the inner wall temperature varying with the axial relative distance x/d using various turbulence models together with the measured inner wall temperature obtained in the experiments for Run 11.

In Fig.9,all the turbulence models can qualitatively predict heat transfer deterioration occurring in the start of heated section,but over predict the wall temperature quantitatively in the region of heat transfer deterioration.The simulation using the LS model can obtain an accurate range for where heat transfer deterioration occurs,while the other models predict the location of the onset of heat transfer deterioration ahead of the experimental results.The range of heat transfer deterioration predicted by AB,LB,YS,AKN and Y2F models is wider than that of the experimental data.The trend of heat transfer recovery can be predicted by most of the models used in the current study,however,the LB model cannot predict that wall temperature decreases after heat transfer deterioration.Compared with the poor prediction accuracy of heat transfer deterioration in the first half of the heated section,most of the turbulence models show better prediction results akin to the experimental data.

According to the classification of turbulence models in Section 4.1.1,the models in Group 2 show better performance than the models in Group 1 for predicting the wall temperature in cases where the buoyancy parameter is large.In order to contrastively investigate the performance of the two groups in the large buoyancy effect cases,the AKN model in Group 1 and the LS model in Group 2 are selected,and the predictions of the inner wall temperature using the LS and AKN models together with the measured inner wall temperature for Runs 11–12 are shown in Fig.10.

Fig.9 Wall temperature along the tubefor Run 11.

Fig.10 Wall temperature along the tube for Runs 11–12.

While flowing upward,buoyancy can make the velocity profile near the wall become more flat and decrease the turbulent shear stress indicating weakened turbulence,which results in the occurrence of heat transfer deterioration.On the contrary,while flowing downward,the buoyancy effect increases turbulence,resulting in enhanced heat transfer,so the wall temperature flowing upward is higher than the wall temperature flowing downward.Both the LS model and AKN model can predict such phenomenon qualitatively,but show different performances quantitatively.The difference between the results flowing upward and downward can determine the occurrence of heat transfer deterioration due to the buoyancy effect.As seen in Fig.10,the LS model more accurately predicts the range where buoyancy effects work than the AKN model does.The results from using the LS model and the AKN model show that heat transfer deterioration occurs in cases flowing downward that cannot be seen in the experimental data,indicating that the obvious thermal-induced acceleration effect leads to heat transfer deterioration.The AKN model responds to the local flow condition much stronger than the LS model does,accounting for the former's higher temperatures.

Fig.11 Wall temperature along the tube for Runs 13–14.

Fig.12 Radial distribution at different axial locations along the pipe using LS for Runs 11–12.

Fig.11 shows the predictions of the inner wall temperature using the LS and AKN models together with the measured inner wall temperature for Runs 13–14.In Run 13,there are two wall temperature peaks in the experiments that can be predicted qualitatively by those two models.From the view of a quantitative analysis,similar results can be obtained as those above.The buoyancy parameter in Runs 13–14 is larger than that in Runs 11–12,but the wall temperature in Runs 13–14 is lower than that in Runs 11–12,meaning that the stronger thermal-induced acceleration effect in Runs 11–12 caused the pronounced heat transfer deterioration.Therefore,both the LS model and AKN model cannot predict the thermalinduced acceleration effect well in the cases considered.

4.3.2.Detailed flow and turbulence fields

Fig.12 shows the radial distributions of axial-velocity and turbulent kinetic energy for Runs 11–12 as predicted by using the LS model at different axial location along the pipe.The axial location x/d=15 represents the inlet of the heated section,with x/d=100 representing the location during the impairment of heat transfer,and x/d=215 in the normal region.

Seen from Fig.12(a),the velocity profiles at the locations of x/d=15,215 in both the upflow and downflow stay the same.As for x/d=100,the velocity profile in the downflow seems to be the typical trend of forced convection while the velocity profile in the upflow is distorted and the 'M' type appears,accompanied by very weak turbulent kinetic energy(Fig.12(b))due to the low velocity gradient near wall and the low turbulent shear stress.This further confirms that the buoyancy effect is the main factor causing the heat transfer deterioration occurring in Run 11.

5.Conclusions

In this study,the performance of turbulence models for predicting the flow and heat transfer of the hydrocarbon RP-3 under supercritical pressure was investigated by comparing computed results with experimental results.The conclusions are shown below:

(1)The performance contrast with experiments can vary markedly from model to model and for the same model with varying conditions.In most cases,turbulence models can reproduce the general trend of heat transfer,but no turbulence models used could predict it quantitatively.

(2)For different states of heat flux at 5 MPa,seven turbulence models are classified into two categories:the first,whose damping functions are based on Retand Rey,respond to the local flow condition too strongly (i.e.,Group 1 models,AB,LB and AKN),leading to the deterioration of heat transfer appearing obviously although under very low buoyancy parameters.The second is comprised of those whose damping functions are based on variables that respond to the local flow condition well(i.e.,Group 2 models,SST,V2F,YS and LS),resulting in a close prediction of the heat transfer exhibited in experiments.

(3)Different models have different response speeds to the variations in the buoyancy parameter.For different inlet temperature operating conditions at 4 MPa,due to the further decrease of buoyancy parameters,the AKN models break away from Group 1.At the same time,all other models agree well with the classification methods summarized above.

(4)The LS model in Group 2 can predict the accurate range of heat transfer deterioration caused by the buoyancy effect.Both the AB model in the Group 1 and the LS model over respond to the thermal-induced acceleration effect in the large buoyancy-influenced cases.

(5)The SST model performs well in predicting the heat transfer of hydrocarbons under supercritical pressure among all the turbulence models used in present study.The LS model performs the most stably in predicting the heat transfer to hydrocarbons under supercritical pressure among all the turbulence models used and has the most potential to be modified for use in simulating such heat transfer.

Acknowledgements

The authors gratefully acknowledge funding support from National Natural Science Foundation of China (No.51406005)and Defense Industrial Technology Development Program of China(No.B2120132006).

1.Andrews T.On the properties of matter in the gaseous and liquid states under various conditions of temperature and pressure.Philos Trans Roy Soc London Ser A Math Phys Eng Sci 1887;178:45–56.

2.Reiss T,Csom G,Fehe´r S,Czifrus S.The simplified supercritical water-cooled reactor (SSCWR),a new SCWR design.Prog Nucl Energy 2010;52(2):177–89.

3.Dang C,Hihara E.In-tube cooling heat transfer of supercritical carbon dioxide.Part 1.Experimental measurement.Int J Refrig 2004;27(7):736–47.

4.Qin J,Zhang S,Bao W,Zhou W,Yu D.Thermal management method of fuel in advanced aeroengines.Energy 2013;49(1):459–68.

5.Mcclinton CR,Hunt JL,Ricketts RH,Reukauf P,Peddie CL.Airbreathing hypersonic technology vision vehicles and development dreams.9th international space planes and hypersonic systems and technologies conference and 3rd weakly lonized gases workshop;1999 Nov 1–5;VA,USA.Reston:AIAA;1999.p.1–16.

6.Curran ET.Scramjet engines:the first forty years.J Propul Power 2001;17(6):1138–48.

7.Edwards T.USAF supercritical hydrocarbon fuels interests.31st aerospace science meeting;1993 Jan 11–14;NV,USA.Reston:AIAA;1993.p.1–11.

8.Jackson JD.Fluid flow and convective heat transfer to fluids at supercritical pressure.Nucl Eng Des 2013;264:24–40.

9.Jackson JD,Cotton MA,Axcell BP.Studies of mixed convection in vertical tubes.Int J Heat Fluid Flow 1989;10(1):2–15.

10.Pioro IL,Duffey RB.Experimental heat transfer in supercritical water flowing inside channels(survey).Nucl Eng Des 2005;235(22):2407–30.

11.Duffey RB,Pioro IL.Experimental heat transfer of supercritical carbon dioxide flowing inside channels(survey).Nucl Eng Des 2005;235(8):913–24.

12.Pioro IL,Duffey RB,Dumouchel TJ.Hydraulic resistance of fluids flowing in channels at supercritical pressures(survey).Nucl Eng Des 2004;231(2):187–97.

13.Kurganov VA,Zeigarnik YA,Maslakova IV.Heat transfer and hydraulic resistance of supercritical-pressure coolants.Part I:Specifics of thermophysical properties of supercritical pressure fluids and turbulent heat transfer under heating conditions in round tubes(state of the art).Int J Heat Mass Transfer 2012;55(11):3061–75.

14.Kurganov VA,Zeigarnik YA,Maslakova IV.Heat transfer and hydraulic resistance of supercritical-pressure coolants.Part II:Experimental data on hydraulic resistance and averaged turbulent flow structure of supercritical pressure fluids during heating in round tubes under normal and deteriorated heat transfer conditions.Int J Heat Mass Transfer 2013;58(1):152–67.

15.Kurganov VA,Zeigarnik YA,Maslakova IV.Heat transfer and hydraulic resistance of supercritical pressure coolants.Part III:Generalized description of SCP fluids normal heat transfer,empirical calculating correlations,integral method of theoretical calculations.Int J Heat Mass Transfer 2013;67:535–47.

16.Hess HL,Kunz HR.A study of forced convection heat transfer to supercritical hydrogen.J Heat Transfer 1965;87(1):41–6.

17.Forooghi P,Hooman K.Effect of buoyancy on turbulent convection heat transfer in corrugated channels–a nu merical study.Int J Heat Mass Transfer 2013;64(4):850–62.

18.Sharabi M,Ambrosini W.Discussion of heat transfer phenomena in fluids at supercritical pressure with the aid of CFD models.Ann Nucl Energy 2009;36(1):60–71.

19.He S,Kim WS,Jackson JD.A computational study of convective heat transfer to carbon dioxide at a pressure just above the critical value.Appl Therm Eng 2008;28(13):1662–75.

20.Masuda Y,Suzuki A,Ikushima Y.Calculation method of heat and fluid flow in a microreactor for supercritical water and its solution.Int J Heat Mass Transfer 2006;33(4):419–25.

21.Wen QL,Gu HY.Numerical simulation of heat transfer deterioration phenomenon in supercritical water through vertical tube.Ann Nucl Energy 2010;37(10):1272–80.

22.He S,Kim WS,Jiang PX,Jackson JD.Simulation of mixed convection heat transfer to carbon dioxide at super critical pressure.Proc Inst Mech Eng Part C-J Mech Eng Sci 2004;218(11):1281–96.

23.Mohseni M,Bazargan M.The effect of the low Reynolds number ke turbulence models on simulation of the enhanced and deteriorated convective heat transfer to the supercritical fluid flows.Heat Mass Transfer 2011;47(5):609–19.

24.Sharabi M,Ambrosini W,He S,Jackson JD.Prediction of turbulent convective heat transfer to a fluid at super critical pressure in square and triangular channels.Ann Nucl Energy 2008;35(6):993–1005.

25.Xu XY,Zeng M,Zhu HB,Wang QW.Numerical simulation and comparison of turbulent heat transfer in super critical and subcritical water.Prog Comput Fluid Dyn 2013;13(3):141–51.

26.Xu XY,Ma T,Zeng M,Wang QW.Numerical study of the effects of different buoyancy models on supercritical flow and heat transfer.Proceedings of the ASME 2013 heat transfer summer conference;2013 July 14–19;Minneapolis,USA.New York:ASME;2013.

27.He S,Kim WS,Bae JH.Assessment of performance of turbulence models in predicting supercritical pressure heat transfer in a vertical tube.Int J Heat Mass Transfer 2008;51(19):4659–75.

28.Kim WS,He S,Jackson JD.Assessment by comparison with DNS data of turbulence models used in simulations of mixed convection.Int J Heat Mass Transfer 2008;51(5):1293–312.

29.Xie W,Fang W,Li D,Xing Y,Guo Y,Lin R.Coking of model hydrocarbon fuels under supercritical condition.Energy Fuels 2009;23(6):2997–3001.

30.Jiang R,Liu G,Zhang X.Thermal cracking of hydrocarbon aviation fuels in regenerative cooling microchannels.Energy Fuels 2013;27(5):2563–77.

31.Zhu Y,Liu B,Jiang P.Experimental and numerical investigations on n-Decane thermal cracking at supercritical pressures in a vertical tube.Energy Fuels 2014;28(1):466–74.

32.Goel P,Boehman AL.Numerical simulation of jet fuel degradation in flow reactors.Energy Fuels 2000;14(5):953–62.

33.Zhang CB.Investigation of flow and heat transfer characteristics of hydrocarbon fuel at supercritical pressures[dissertation].Beijing:Beihang University;2011.

34.Liu B,Zhu YH,Yan JJ,Lei YT,Zhang B,Jiang PX.Experimental investigation of convection heat transfer of n-decane at supercritical pressures in small vertical tubes.Int J Heat Mass Transfer 2015;91:734–46.

35.Abid R.Evaluation of two-equation turbulence models for predicting transitional flows.Int J Eng Sci 1993;31(6):831–40.

36.Abe K,Kondoh T,Nagano Y.A new turbulence model for predicting fluid flow and heat transfer in separating and reattaching flows-I.Flow field calculations.Int J Heat Mass Transfer 1994;37(1):139–1351.

37.Lam CG,Bremhorst K.A modified form of the k-ε model for predicting wall turbulence.J Fluids Eng-Trans ASME 1981;103(3):456–60.

38.Launder BE,Sharma BI.Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc.Int Commun Heat Mass Transfer 1974;1(2):131–7.

39.Yang Z,Shih TH.New time scale based k-epsilon model for nearwall turbulence.AIAA J 1993;31(7):1191–8.

40.Menter FR,Zonal two equation k-turbulence models for aerodynamic flows.24th fluid dynamics conference;1993 Jul 6–9;Orlando,USA.Reston:AIAA;1993.p.1–21.

41.Behnia M,Parneix S,Durbin PA.Prediction of heat transfer in an axisymmetric turbulent jet impinging on a flat plate.Int J Heat Mass Transfer 1998;41(12):1845–55.

42.Petukhov BS,Polyakov AF,Launder BE.Heat transfer in turbulent mixed convection.New York:Hemisphere Publishing Corporation;1988.p.109–13.

43.Fan XJ,Yu G.Analysis of thermophysical properties of Daqing RP-3 aviation kerosene.J Propulsion Technol 2006;27(2):187–92[Chinese].

44.Zhong FQ,Fan XJ,Yu G,Li JG.Heat transfer of aviation kerosene at supercritical conditions. J Thermophys Heat Transfer 2009;23(3):543–50.

45.Deng HW,Zhang CB,Xu GQ,Tao Z,Zhang B,Liu GZ.Density measurements of endothermic hydrocarbon fuel at sub-and supercritical conditions.J Chem Eng Data 2011;56(6):2980–6.

46.Deng HW,Zhang CB,Xu GQ,Tao Z,Zhang B,Zhu K.Viscosity measurements of endothermic hydrocarbon fuel from(298 to 788)K under supercritical pressure conditions.J Chem Eng Data 2012;57(2):358–65.

47.Deng HW,Zhu K,Xu GQ,Tao Z,Zhang CB,Liu GZ.Isobaric specific heat capacity measurement for kerosene RP-3 in the nearcritical and supercritical regions.J Chem Eng Data 2011;57(2):263–8.

48.Huber ML.NIST thermophysical properties of hydrocarbon mixtures database (SUPERTRAPP).NIST Standard Reference Database;2003.

49.Jackson JD.Consideration of the heat transfer properties of supercritical pressure water in connection with the cooling of advanced nuclear reactors.Proceedings of the 13th pacific basin nuclear conference;2002 Jan 1;Shenzhen,China.Beijing:Chinese Nuclear Society:2002.p.21–5.

50.Hall WB.Heat transfer near the critical point.Adv Heat Transfer 1971;7:1–86.

Tao Zhiis a professor and Ph.D.supervisor at School of Energy and Power Engineering,Beihang University,Beijing,China.His current research interests are heat transfer under rotation,supercritical heat transfer,micro-scale heat transfer and heat transfer in porous media,etc.

Cheng Zeyuanis a Ph.D.student at School of Energy and Power Engineering,Beihang University.His area of research is supercritical heat transfer.

Zhu Jianqinreceived the B.S.and Ph.D.degrees in engineering thermophysics from Beihang University in 2005 and 2010 respectively,and the n became a teacher there.Her main research interests are engine thermal protection and cooling technology in aero-engine.

Li Haiwangis an associate professor and Ph.D.supervisor at School of Energy and Power Engineering,Beihang University,Beijing,China.His current research interests are engine thermal protection and cooling technology in aero-engine.

25 May 2015;revised 1 December 2015;accepted 13 May 2016

Available online 26 August 2016

Buoyancy effect;

Hydrocarbon fuel;

Supercritical pressure;

Turbulence models;

Variable properties

©2016 Chinese Society of Aeronautics and Astronautics.Production and hosting by Elsevier Ltd.This is an open access article under the CC BY-NC-ND license(http://creativecommons.org/licenses/by-nc-nd/4.0/).

*Corresponding author.Tel.:+86 10 82339181.

E-mail addresses:tao_zhi@buaa.edu.cn(Z.Tao),cheng_zeyuan@buaa.edu.cn(Z.Cheng),zhujianqin@buaa.edu.cn(J.Zhu),09620@buaa.edu.cn(H.Li).

Peer review under responsibility of Editorial Committee of CJA.