Numerical investigation on convective heat transfer to aviation kerosene flowing in vertical tubes at supercritical pressures

International Journal of Heat and Mass Transfer 118 (2018) 857–871

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Numerical investigation on convective heat transfer to aviation kerosene flowing in vertical tubes at supercritical pressures Hang Pu, Sufen Li ⇑, Si Jiao, Ming Dong, Yan Shang School of Energy and Power Engineering, Key Laboratory of Ocean Energy Utilization and Energy Conservation of Ministry of Education, Dalian University of Technology, Dalian 116024, PR China

a r t i c l e

i n f o

Article history: Received 2 June 2017 Received in revised form 14 October 2017 Accepted 8 November 2017

Keywords: Supercritical pressure heat transfer Aviation kerosene Turbulence model Buoyancy effect Turbulent Prandtl number

a b s t r a c t Numerical simulations of convective heat transfer to aviation kerosene (China RP-3) flowing in vertical circular tubes at supercritical pressures are reported in this study. Firstly, performance of a variety of Reynolds-Averaged Navier-Stokes turbulence models in predicting the fluid-thermal behaviours under both forced and mixed convection conditions are evaluated. Under forced convection conditions, all models predict a gentler growth of wall temperature along the flow direction than experimental measurements. Under mixed convection conditions, the effect of buoyancy become significant and there are large discrepancies in the predicted wall temperature by different models. Only the low-Reynolds number k-e models are found to be able to qualitatively predict the flow laminarization and heat transfer deterioration. Profiles of thermal, flow and turbulence fields obtained using various models are studied to explain the differences in predictions. For mixed convection conditions, an examination on the turbulence production due to shear and density fluctuation indicates that the direct effect of buoyancy on the turbulence production is negligible compared with the indirect effect. Furthermore, the effect of turbulent Prandtl number on the predicted heat transfer is studied. It is found that turbulent Prandtl number has a significant influence on the simulation results. Under the conditions considered in the present study, the value of 1.0 for turbulent Prandtl number leads to a closest agreement with the experimental data. Ó 2017 Elsevier Ltd. All rights reserved.

1. Introduction Development of sustained hypersonic air-breathing propulsion (HAP) technologies have received global attentions since the 1960s [1]. At the hypersonic flight regime (Ma > 6), the specific impulse Isp of scramjet engine surpasses that of any other type of propulsion options, which makes it the key enabling technology for hypersonic flight applications within atmosphere [2,3]. Due to the drastic aerodynamic heating of vehicle body as well as the extreme high heat load released by the supersonic combustion process, thermal management remains to be a significant technical challenge in the designs of scramjet engines. For operations aimed at long-range flights, regenerative cooling that utilizing onboard fuel as the primary coolant has been considered to be an effective way to prevent the engine wall temperature from exceeding the material limit [4]. The choice of fuel for scramjet engine is another important issue that needs to be carefully considered. Hydrogen could provide the ⇑ Corresponding author. E-mail address: [email protected] (S. Li). https://doi.org/10.1016/j.ijheatmasstransfer.2017.11.029 0017-9310/Ó 2017 Elsevier Ltd. All rights reserved.

highest energy release and cooling capacity, however it has to face with logistic, cost, and safety problems. In comparison with hydrogen, hydrocarbons have larger densities, which enables a smaller volume and weight of the engine body. The unfavourable safety and operation problems associated with hydrogen could also be avoided [5]. At Ma < 8, hydrocarbon fuel is able to provide sufficient heat sink to meet the cooling requirements, which makes it a good fuel candidate for the lower range of hypersonic flight [6–8]. There have been several demonstrations on the viability of regeneratively cooled scramjet engine using hydrocarbon fuels [9–12]. Before entering the combustion chamber, the fuel flows along the cooling channels (in the order of millimetres) which surround the combustion chamber and takes away heat from the wall [13,14]. The typical aircraft fuel system pressure is higher than the critical pressure of most hydrocarbon fuels [15]. For fluids at supercritical pressures, the variation of thermodynamic and transport properties with temperature becomes significant in the heat transfer process. Particularly, the temperature at which the specific heat capacity achieves its peak value at a given pressure is known as the pseudo-critical temperature (Tpc) [16]. In the vicinity of Tpc, the rapid decrease of fluid density would bring in the effects of

858

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

Nomenclature A Bo⁄

cross sectional area of the tube non-dimensional buoyancy parameter, Bo⁄ = Gr/(Re3.425Pr0.8) cp isobaric specific heat capacity Ce1, Ce2, Ce3 constants in the e-equation Cl constant in the constitutive equation of eddy viscosity model D inner diameter of the tube; damping function term in the k-equation E damping function term in the e-equation f body force; elliptic relation parameter f1, f2 damping functions in the e-equation fl damping function in the constitutive equation of eddy viscosity model g gravitational acceleration G mass flux Gk buoyancy production term of turbulent kinetic energy Gr Grashof number, Gr = gbD4qw/(kv2) h convective heat transfer coefficient Isp specific impulse k turbulent kinetic energy Kv non-dimensional thermal acceleration parameter, Kv = 4bqw/(qcpUbRe) L turbulent length scale Ma Mach number Nu Nusselt number, Nu = hD/kf p pressure pc critical pressure shear production term of turbulent kinetic energy Pk Pr molecular Prandtl number qw wall heat flux r radial coordinate; radius of the tube Re Reynolds number, Re = UbD/v Ret turbulent Reynolds number, Ret = k2/(ve) Sh energy source term ST volumetric heat source t time T temperature; turbulent time scale

buoyancy and thermal acceleration, which would induce modifications of mean flow and turbulence fields. These effects, combined with the non-uniform distribution of fluid properties across the boundary layer would have a significant influence on the heat transfer process. In order to investigate the convective heat transfer characteristics of hydrocarbon fuel at supercritical pressures, a lot of experimental efforts have been carried out. Hitch and Karpuk [17] conducted experiments of heat transfer to supercritical JP-7 in vertical circular tubes. Heat transfer deterioration (HTD) along with large temperature and pressure oscillations were observed when the reduced pressure (p/pc) was below 1.5 and the wall temperature was higher than Tpc. Hu et al. [18], Zhong et al. [13], Zhang et al. [19], Liu et al. [20], and Li et al. [21] experimentally studied convective heat transfer of China RP-3 aviation kerosene flowing in circular tubes at supercritical pressures. In these studies, the effects of operating parameters (mass flux, wall heat flux, pressure, inlet temperature, and flow direction) on the heat transfer process were carefully investigated. Under large mass flux conditions, heat transfer enhancement (HTE) was observed when the wall temperature exceeded Tpc; under small mass flux conditions, dramatic increase of wall temperature was observed at the entrance region for both downward and upward flows. This kind of HTD was attributed to the slow development of thermal boundary layer and the

Tpc Twi U U+ u0j t 0 u0i u0j v

v2

x y y+

pseudo-critical temperature inner wall temperature velocity dimensionless velocity turbulent heat flux Reynolds stress kinematic viscosity variance of the normal component of turbulent velocity axial coordinate normal distance from the wall dimensionless distance from the wall

Greek symbols b isobaric thermal expansivity d Kronecker delta e dissipation rate of turbulent kinetic energy k thermal conductivity l dynamic viscosity lt turbulent viscosity q density rk turbulent Prandtl number for k re turbulent Prandtl number for e rt turbulent Prandtl number s viscous stress x dissipation per unit turbulent kinetic energy Subscripts b bulk f fluid; forced convection in inlet of the tube i, j spatial indices out outlet of the tube ref reference w wall s solid

buoyancy effect, respectively. Furthermore, a number of semiempirical correlations for predicting heat transfer to RP-3 at supercritical pressures have been proposed based on the experimental data, as shown in the technical note by Chen and Fang [22]. It should be noted that the data obtained in most of the experimental studies are limited to wall temperature, inlet/outlet fluid temperature, and pressure drop along the test tube. As an alternative approach, computational fluid dynamics (CFD) method could offer more detailed information of thermal, flow, and turbulence fields, which is necessary to gain deeper understanding of the underlying physics of the heat transfer process. Therefore, growing attentions have been paid on numerical simulations on flow and heat transfer of hydrocarbon fuel at supercritical pressures. In turbulent heat transfer simulations, the predicted flow field and heat transfer rate is highly dependent on the turbulence modelling method. Different types of Reynolds-Averaged Navier-Stokes (RANS) turbulence models were employed in the existing simulations: the shear stress transport (SST) k-x model (Meng et al. [23]), the standard k-x model (Zhu et al. [24]), and the renormalization group (RNG) k-e model with the enhanced wall functions (Zhong et al. [25], Liu et al. [26]). However, it is acknowledged that most of the turbulence models are developed for constant property fluids. Therefore, the applicability of these models in predicting heat transfer to fluids at supercritical pressures (associated with

859

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

severe property variations) is questionable and needs to be carefully investigated. Recently, some efforts have been made concerned with this issue. Jiang et al. [27] compared three turbulence models and concluded that using the RNG k-e model with the enhanced wall functions could obtain the best agreement with the experimental data. In the study by Tao et al. [28], seven types of turbulence models were evaluated and the results indicated that the performance of turbulence model was determined by their responses to the variation of local flow conditions. Due to the few existing studies by far, no firm conclusions on the evaluation of turbulence models could be made. Further investigations over this problem considering a wider range of flow conditions are still needed. Another important aspect in the simulation of supercritical heat transfer lies in the modelling of turbulent heat flux. In most of the numerical works, it is calculated based on the Reynolds analogy, in which the ratio of momentum diffusivity to thermal diffusivity is defined as the turbulent Prandtl number. The turbulent Prandtl number is set to a constant value of unity or close to it (0.9 or 0.85). A number of investigations on the effect of turbulent Prandtl number on convective heat transfer to supercritical fluids (carbon dioxide) have been performed [29,30]. However, there are some discrepancies between the results. To the authors’ knowledge, there has been no relevant work associated with this issue for hydrocarbon fuel. The effect of turbulent Prandtl number on predicted heat transfer will also be discussed in this work. In this paper, convective heat transfer to a specific hydrocarbon fuel China RP-3 aviation kerosene at supercritical pressures has been numerically investigated. The numerical simulation is performed using an in-house code developed using the open-source C++ library OpenFOAM [31,32]. The results obtained using a variety of Reynolds-Averaged Navier-Stokes (RANS) turbulence models are compared with previous conducted experimental studies to examine their accuracy. Thermal, flow, and turbulence fields across the turbulent boundary layer are examined to explain the difference in prediction by various models. For the mixed convection, the turbulent kinetic energy production due to shear and density fluctuation are examined to investigate the mechanism of buoyancy effect on the turbulence production. The effect of turbulent Prandtl number on predicted heat transfers rate is also studied by using different constant values of turbulent Prandtl number in the simulations. Some useful conclusions are derived from the results obtained in the present study.

pressure and temperature of RP-3 is 2.39 MPa and 645.5 K, respectively [33]. In the study, both forced and mixed convection conditions were considered. The buoyancy and thermal acceleration effect was estimated by calculating non-dimensional parameter Boin and Kv,in, respectively. The conditions that have been simulated in the study, together with the buoyancy and thermal acceleration parameters are summarized in Table 1. 2.2. Governing equations RP-3 flowing in a uniformly heated tube at supercritical pressures was considered in the simulation. The governing equations of fluid flow are the Reynolds-Averaged Navier-Stokes equations, which are written in the conservation form and cylindrical coordinates as follows: 2.2.1. Continuity

@ qf @ 1 @ ðrqf U r Þ ¼ 0 þ ðqf U x Þ þ @x r @r @t where qf is the density of fluid, t is time. 2.2.2. Momentum

@ @ 1 @ ðq U x Þ þ ðqf U x U x Þ þ ðr qf U x U r Þ @t f @x r @r @p @ sxx 1 @ @ ¼ þ ðrsxr Þ þ ðqf u0x u0x Þ þ @x r @r @x @x 1 @ þ ðr qf u0x u0r Þ þ qf f x r @r

2.1. Description of simulation conditions and experimental data used A range of simulation conditions representing different degrees of buoyancy effect on convective heat transfer to China RP-3 aviation kerosene flowing in vertical circular tubes were selected from the experimental studies by Zhang et al. [19] and Li et al. [21]. The accuracy of measured wall temperature is ±0.4%. The critical

ð2Þ

@ @ 1 @ ðq U r Þ þ ðqf U r U x Þ þ ðr qf U r U r Þ @t f @x r @r @p @ srx 1 @ @ ¼ þ ðrsrr Þ þ ðqf u0r u0x Þ þ @r r @r @x @x 1 @ ðr qf u0r u0r Þ þ qf f r þ r @r where U is velocity, p is pressure, f is body force,

ð3Þ

sij is the viscous

stress tensor, and qf u0i u0j is the Reynolds stress tensor. The viscous stress tensor is calculated by

sij ¼ lf 2. Methodology

ð1Þ


 @U i @U j þ @xj @xi

ð4Þ

where lf is the molecular viscosity of fluid. For linear eddy viscosity models, the Reynolds stress tensor is modelled using the Boussinesq hypothesis:


 2 @U i @U j qf u0i u0j ¼  qf kdij þ lt;f þ 3 @xj @xi

ð5Þ

where k is the turbulent kinetic energy, dij is the Kronecker delta,

lt;f is the turbulent viscosity.

Table 1 Simulation conditions of the cases. qw (kW m2)

qw/G

Rein

Boin

Kv,in

Dir.

(a) Forced convection A1 1572.7 A2 1572.7

300 500

0.19 0.32

14905.03 14905.03

5.54  109 9.24  109

2.09  108 3.49  108

Down Down

(b) Mixed convection B1 578 B2 481

250 250

0.43 0.52

3569.81 2970.73

4.34  107 8.14  107

3.38  107 4.88  107

Up Up

Case

G (kg m2 s1)

Notes: (a) D = 1.805 mm, p = 5 MPa, Tf,in = 473 K. Notes: (b) D = 2.3 mm, p = 3.5 MPa, Tf,in = 418 K.

860

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

@ @ 1 @ ðq kÞ þ ðqf U x kÞ þ ðr qf U r kÞ @t f @x r @r     

@ l @k 1 @ l @k þ r lf þ t þ Pk þ Gk lf þ t ¼ @x rk @x r @r rk @r

2.2.3. Energy

@ @ 1 @ ðq T f Þ þ ðqf U x T f Þ þ ðr qf U r T f Þ @t f
@x  r @r
 @ kf @T f 1 @ kf @T f @ ¼ r þ þ ðqf u0x T 0 Þ @x cp;f @x r @r cp;f @r @x 1 @ ðr qf u0r T 0 Þ þ r @r

 qf e þ qf D ð6Þ

where Tf is fluid temperature, kf is the thermal conductivity of fluid, cp;f is the isobaric specific heat capacity of fluid. qf u0i T 0 is the turbulent heat flux tensor, which is modelled by

lt;f @T f rt @xi

ð11Þ

@ @ 1 @ ðq eÞ þ ðqf U x eÞ þ ðr qf U r eÞ @t f @x r @r     

@ l @e 1 @ l @e e r lf þ t lf þ t þ þ f 1 C e1 P k ¼ @x re @x r @r re @r k

e2

e

þ C e1 C e3 Gk  f 2 C e2 qf þ qf E k k

ð12Þ

ð7Þ

where lt is the turbulent viscosity. Specially, in the LS k-e and CH ke model the transport equation for e ¼ e  D is solved instead of e. Pk

where rt is the turbulent Prandtl number. rt was set to be constant as 1.0 in the simulations. The convective heat transfer coefficient h is calculated by

is the shear production term of turbulent kinetic energy, which is modelled by

qf u0i T 0 ¼

h¼

qw T wi  T b

ð8Þ

where Twi is the inner wall temperature, the bulk temperature Tb is calculated by

R q U f cp;f T f dA Tb ¼ R f qf U f cp;f dA

ð9Þ

The heat conduction in the tube wall is governed by the Fourier’s law:

qs cp;s


 @T s @ @T s þ ST ks ¼ @xi @t @xi


 @U i @U j @U i þ @xj @xi @xj

2.3. Turbulence models Accurate prediction of the flow and thermal field in the nearwall region plays an important role in simulating convective heat transfer of fluids at supercritical pressures. In this study, five kinds of RANS turbulence models are considered: the low-Reynolds number (LRN) k-e model by Launder and Sharma (LS) [34] and by Chien (CH) [35], the modified version of v2-f model by Lien and Kalitzin [36], the shear stress transport (SST) k-x model by Menter [37], and the renormalization group (RNG) k-e model by Yakhot et al. [38]. These models use different near-wall treatments for wall bounded flows. In the LS k-e and CH k-e model, damping functions are added to the definition of turbulent eddy viscosity and in the transport equation of turbulent dissipation rate to consider the near-wall effects. The v2-f model uses the turbulent velocity scale instead of turbulent kinetic energy k to represent the suppression of the turbulent intensity at wall, thus the damping function are no longer required [39,40]. The SST k-x model switches to the original k-x form in the near-wall region [37]. The above four turbulence models use the integration method to resolve the boundary layer with a large number of cells close to the wall [41]. Another class of near-wall treatment is to use wall functions to predict fluid flow in the near-wall region. In this study, the RNG k-e model is solved with the low-Reynolds type wall functions (WF) for the turbulent kinetic energy k and the turbulent energy dissipation rate e. The general form of transport equations for the k-e family turbulence models are as follows:

ð13Þ

The buoyancy production term of turbulent kinetic energy Gk is modelled by the generalized gradient-diffusion hypothesis (GGDH) method.

Gk ¼ 


 C h b k @T @U i @U j gi lt þ T e @xi @xj @xi

ð14Þ

where b is the isobaric thermal expansivity, gi is the gravitational acceleration, Ch = 0.3. The turbulent eddy viscosity lt is calculated by

ð10Þ

where qs is the density of solid, cp,s is the heat capacity of solid, ST is the volumetric heat source.

v2

P k ¼ lt

lt ¼ f l C l qf

k

2

ð15Þ

e

where fl is the damping function, Cl is a constant. In the v2-f model, the transport equations for v 2 and the elliptic relaxation equation for f are shown as below:

@ @ 1 @ ðq v 2 Þ þ ðqf U x v 2 Þ þ ðr qf U r v 2 Þ @t f @x r @r "
"
#  2#  @ l @v 1 @ l @v 2 e ¼ þ þ kf  v 2 r lþ t lþ t rk @x rk @r k @x r @r " #
# @2f 1 @ @f 1 v2 2 Pk  f ¼ ðC 1  1Þ  C2 L r þ  @x2 r @r @r T k 3 k

ð16Þ

"

2

ð17Þ

where kf represents the unclosed redistribution of the turbulent energy from the streamwise component, C1 = 1.4, C2 = 0.3. The turbulent eddy viscosity is defined as

lt ¼ C l v 2 T

ð18Þ

The turbulent time scale T and turbulent length scale L are given by:

 T ¼ max

k

e

rffiffiffiffi

;6 "

L ¼ C L max

v




3=2

k

e

ð19Þ

e ; Cg

v3 e

1=4 # ð20Þ

where CL = 0.23, Cg = 70. The constants and damping functions used in the k-e family turbulence models are listed in Table 2. The transport equations for the SST k-x model are as follows:

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

model are: b⁄ = 0.09, rk1 = 0.85, rk2 = 1.0, rx1 = 0.5, rx2 = 0.856, a1 = 0.31, b1 = 0.075, b2 = 0.00828.

Table 2 Details of the k-e family turbulence models. Model

Cl

Ce1

(a) Model constants and additional terms LS k-e 0.09 1.44 CH k-e 0.09 1.35 2 v -f 0.22 1.40 0.0845 1.42 RNG k-e Model

Ce2

rk

re

1.92 1.8 1.9 1.68

1.0 1.0 1.0 0.71942

1.3 1.3 1.3 0.71942

2.4. Numerical methods

D

(b) Additional source terms LS k-e

E

CH k-e

2v yk2

v -f RNG k-e

 2 2 2v lt @@yU2   2v ye2 expð0:5yþ Þ

0 0

0 0

2

Model

2v



pffiffi2 @ k @y

fl

(c) Damping functions h LS k-e exp

3:4 ð1þRet =50Þ2

i

f1

f2

1.0

1  0:3 exp½ minðRe2t ; 50Þ h 2 i 1  0:22 exp  Re6 t

CH k-e

1  expð0:0115yþ Þ

1.0

v -f RNG k-e

– 1.0

– 1.0

2

pffiffiffiffiffiffiffiffiffiffiffiffi 2 Note: Ret ¼ k =ev , yþ ¼ y sw =q=v .

– 1.0

@ @ 1 @ ðq kÞ þ ðq U x kÞ þ ðr qf U r kÞ @t f  @x f  r @r   @ @k 1 @ @k þ ¼ ðl þ rk lt Þ rðl þ rk lt Þ @x @x r @r @r þ Pk þ Gk  qf b kx @ @ 1 @ ðq xÞ þ ðqf U x xÞ þ ðr qf U r xÞ @t f @x r @r     @ @x 1 @ @x ðl þ rx lt Þ rðl þ rx lt Þ þ þ Px þ Gx ¼ @x r @r @x @r

 rx2 @k @k @ x @ x þ þ 2ð1  F 1 Þ þ  qf bx2 q x f @x @r @x @r

ð21Þ

ð22Þ

where x is the dissipation per unit turbulence kinetic energy, F1 is the blending function calculated as follows:

F 1 ¼ tanhðarg41 Þ ! # pffiffiffi k 500v 4qf rx2 k arg1 ¼ min max  ; ; CDkx y2 b xy y2 x

ð23Þ

"


 rx2 @k @ x 10 ; 10 CDkx ¼ max 2qf x @xi @xi

ð24Þ

ð25Þ

The turbulent eddy viscosity in the SST k-x model is defined by:

lt ¼

a1 qf k maxða1 x; SF 2 Þ

ð26Þ

where S is the invariant measure of the strain rate, the blending function F2 is calculated as follows:

F 2 ¼ tanhðarg22 Þ arg2 ¼ max

! pffiffiffi 2 k 500v ; b xy y2 x

ð27Þ ð28Þ

All of the constants / in the model are computed by:

/ ¼ F 1 /1 þ ð1  F 1 Þ/2

861

Since the fluid flows vertically inside the tube, the flow was assumed to be two-dimensional and axisymmetric. The computational configuration is shown in Fig. 1. For the solid region, two insulated section were included at both upstream and downstream of the heated section. Volumetric heat source was placed in the heated section of the tube wall corresponding to the different wall heat flux conditions. The outside walls and both end of the tube were assumed to be adiabatic. For the flow field, uniform velocity, temperature and turbulent parameter profiles were specified at the inlet, and outflow boundary were set at the outlet. Zero gradient of pressure was set at the inlet, and constant pressure was specified at the outlet (according to the operating pressure of experiments). At the fluid-solid interface, the Robin-Dirichlet condition was implemented to deal with the coupled thermal interaction and no-slip condition was specified for velocity. Besides, specific boundary conditions at the fluid-solid interface for turbulence parameters are summarized in Table 3. It should be pointed out that a small value 1010 is used instead of 0 at wall to obtain numerical stability. The numerical simulation was conducted using an in-house code developed in the framework of OpenFOAM [31,32]. Coupled heat transfer between the fluid convection and solid conduction was solved using a partitioned approach. The detailed algorithm, coupling procedure, and validation of the code could be found in the authors’ previous work [42]. Basically, discretization and solution of the partial differential equations were carried out using the finite volume method (FVM). A modified version of PISO algorithm was used with one additional temperature predictor-correction step [43]. The QUICK scheme was used for approximation of the convection terms in the equations and the diffusion terms were approximated by the central differencing scheme. Thermodynamic and transport properties of fluid and solid were written in a lookup table and linked to the solver. The steady-state solution was obtained using the time-marching approach, it was considered to be converged if the following criteria were satisfied: (1) the normalized residual of all variables were smaller than 105. (2) The continuity of temperature and heat flux across the fluid-solid interface was satisfied. The computational domain, including both solid subdomain and fluid subdomain, was discretized into structured quadrilateral grids. The grids were refined in the radial direction towards the fluid-solid interface in both subdomains. In the fluid subdomain, the dimensionless wall distance y+ at the first cell next to the wall was ensured to be smaller than 0.5 for all the runs. The grid study was carried out for all the test cases with a single turbulence model (the SST k-x model) to ensure that the solutions were independent of the grid numbers. As an example, Table 4 shows the calculated outlet bulk fluid temperature by the SST k-x model using different grid numbers for Case B2. It is shown that the difference between Grid 2 and Grid 4 in the outlet bulk temperature is less than 1 K. Therefore, the grid with the first level of refinement (Grid 2) was used for Case B2 (also for Case B1). Following a similar procedure, a grid of (10 + 60)  (90 + 300 + 60), representing (radial solid + radial fluid)  (axial insulated section + axial heated section + axial insulated section) was used for Case A1 and Case A2. 2.5. Fluid properties

ð29Þ

where /1 are the constants of the original k-x model, /2 are the constants of the standard k-e model. The constants in the used

The working fluid China RP-3 aviation kerosene is a multicomponent hydrocarbon fuel, which is composed of alkanes, alkenes, benzenes, cycloalkanes, naphthalenes, and other components

862

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

Fig. 1. Computational configuration.

Table 3 Boundary conditions at fluid-solid interface for turbulence parameters. Turbulence model

e(e)

k 10

10

LS k-e CH k-e v2-f

10 1010 1010

10 1010

RNG k-e SST k-x

WF 1010

WF /

lf @ k qf @y2 2

Table 4 Comparison of calculated outlet bulk fluid temperature with different gird numbers for Case B2.

a

Grid

Solid subdomaina

Fluid subdomain

Tb,out (K)

1 2 3 4

20  600 20  1200 20  2400 20  4800

60  600 60  1200 60  2400 60  4800

624.07 637.71 637.95 638.21

Nodes in the radial direction  axial direction.

[44]. Surrogate mixtures which have generally the same properties as the fuel are often used in experimental and numerical studies for tractability and reproducibility [45]. In this paper, thermodynamic and transport properties of RP-3 were calculated using the NIST SUPERTRAPP program [46] with a 10-component surrogate model proposed by Zhong et al. [13]. Fig. 2 shows the density, isobaric specific heat capacity, dynamic viscosity, and thermal conductivity of RP-3 under different pressures. 3. Results and discussion 3.1. Turbulence model studies under forced convection conditions The results of convective heat transfer under forced convection conditions are firstly presented. According to the investigation by Liu et al. [47], the threshold of Bo⁄ for n-decane at supercritical pressures is about 2  107. This criterion is used to evaluate the effect of buoyancy in this study. The effect of thermal acceleration is evaluated using the criterion Kv = 9.5  107 proposed by Murphy et al. [48]. As shown in Table 1, both the effect of buoyancy and thermal acceleration are negligible for Case A1 and Case A2. Therefore, it is reasonable to consider that these two cases are under forced convection conditions. Fig. 3 shows the variation of local inner wall temperature along the tube predicted by various turbulence models together with the experimental data for Case A1 and Case A2. The location at x/d = 0 corresponds to the starting point of the heated section. Experimental data show that the inner wall temperature increase monotonically along the tube for both cases. The maximum wall

x

v2

f

/ / /

/ / 1010

/ / 1010

/

/ /

/ /

60v f 0:075y2

temperature is about 663 K, which is lower than the pseudocritical temperature (Tpc = 716 K at p = 5 MPa). It could be seen from Fig. 3 that all of the turbulence models predict a gentler growth of inner wall temperature along the flow direction than experimental measurements. The uncertainty on the fluid properties (calculated using the surrogate model) could play a relevant role in such differences. Specifically, the LS k-e, CH k-e, v2-f, and RNG k-e model (with wall functions) over-predict the inner wall temperature in the front part of the tube. Under the larger wall heat flux condition (Case A2), the difference between the predicted local wall temperature and experimental data could be as large as about 40 K when the RNG k-e model (with wall functions) is used. In general, the results obtained by the SST k-x model are in good agreement with the experimental data, except that the wall temperature is under-estimated in a short distance near the tube outlet. Heat transfer under forced convection conditions are further studied by considering the ratio of the Nusselt number calculated from the simulation to that calculated using the semi-empirical correlation. Only the numerical results obtained by the SST k-x model are used here since this model performs best in the prediction of inner wall temperature. The reference Nusselt number Nuref is calculated using the correlation proposed by Hu et al. [18], in which the effects of property variation, buoyancy or thermal acceleration are not considered:

Nuref ¼ 0:008Re0:873 Pr 0:451 b b

ð30Þ

The Nusselt number ratios of Case A1 and Case A2 are shown in Fig. 4. Six axial locations are selected to show the variation of Nu/ Nuref along the tube: x/d = 20, 50, 80, 110, 140, and 166. As seen from Fig. 4, Nu/Nuref increases monotonously as the fluid flows downstream. The SST k-x model predicts Nusselt numbers very close to that calculated using the Hu correlation in both cases (0.91 < Nu/Nuref < 1.05). Based on the discussion above it is concluded that heat transfer under forced convection conditions could be reasonably predicted with the SST k-x model. 3.2. Turbulence model studies under mixed convection conditions Comparing with the forced convection conditions, the effect of buoyancy brings in difficulties in simulating heat transfer under

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

863

(a) p=3.5MPa

(b) p=5MPa Fig. 2. Thermodynamic and transport properties of China RP-3 aviation kerosene under different pressures.

(a) Case A1

(b) Case A2

Fig. 3. Comparison of predicted inner wall temperature with experimental data under forced convection conditions.

mixed convection conditions, which will be discussed in this section. The buoyancy parameters Boin for Case B1 and Case B2 are larger than 2  107, which indicates that buoyancy force has significant influences on fluid flow and heat transfer. Again the effect of thermal acceleration is negligible in these two cases as shown in Table 1.

Fig. 5 shows the variation of inner wall temperature along the tube predicted by various turbulence models together with the experimental data for Case B1 and Case B2. It is found that only the low-Reynolds number (LRN) k-e models are able to qualitatively predict the deteriorated heat transfer at the entrance region. In Case B1, the LS k-e model under-predicts the wall temperature.

864

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

Fig. 4. Nusselt number ratios under forced convection conditions.

Besides, the predicted range of heat transfer deterioration (HTD) is smaller than that of experimental measurements. The CH k-e model predicts a higher wall temperature peak as well as a wider

(a) Case B1

range of HTD. Under the stronger buoyancy-influenced condition (Case B2), both models predict a wider region influenced by buoyancy compared with the experimental data: the inner wall temperature begins to increase again at about x/d = 400 in the prediction using the LS k-e model; and the predicted inner wall temperature by the CH k-e model show a continuous decrease as the fluid travels downstream of the wall temperature peak. According to Kim et al. [49], such discrepancies are due to the different damping function used by the models, which respond differently to changes of local flow conditions under the buoyancy effect. For the other three models that fails to predict the wall temperature peak near the tube inlet, the SST k-x model reproduces a closest wall temperature distribution compared with experimental data. Large discrepancies could be observed between the experimental data and numerical results using the v2-f and the RNG k-e model (with wall functions). This indicates that these two models predict unrealistic near-wall flow field and heat fluxes across the fluid-solid interface, therefore the temperature field could not be accurately predicted. The effect of buoyancy on convective heat transfer is studied by examining the ratio of Nusselt number from simulation to that calculated using a semi-empirical correlation which considers the effect of property variation, but not the buoyancy. The correlation proposed by Deng et al. [50] is used:

(b) Case B2

Fig. 5. Comparison of predicted inner wall temperature with experimental data under mixed convection conditions.

(a) Case B1

(b) Case B2

Fig. 6. Nusselt number ratios under mixed convection conditions.

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871 0:8777 Nuref ¼ 0:008151 Re0:95 Pr 0:4 ðlw =lb Þ0:6352 ðT b =T pc 6 0:9Þ b b ðc p;w =c p;b Þ

ð31Þ 1:447 Pr 0:4 ðlw =lb Þ0:4336 ðT b =T pc > 0:9Þ Nuref ¼ 0:002317 Re0:87 b b ðc p;w =c p;b Þ

ð32Þ The numerical results obtained by the LS k-e and SST k-x model are discussed here. The LS k-e model is selected as a representative of the LRN k-e models, and the SST k-x model represents the group of models that are unable to reproduce the buoyancy effect on flow and heat transfer. The Nusselt number ratios of Case B1 and Case B2 are shown in Fig. 6. Six axial locations are selected to show the variation of Nu/Nuref along the tube: x/d = 35, 100, 200, 300, 400, and 448. For Case B1, both models predict about 20%

(a) LS k-ε model

865

reduction in Nu at the entrance region, which is caused by the buoyancy effect. The Nusselt ratios calculated from the results using the SST k-x model are larger than that using the LS k-e model, therefore no local wall temperature peak is predicted by the SST k-x model as shown in Fig. 5(a). Case B2 represents a condition in which the buoyancy effect is stronger. It is worth noting that the LS k-e model predicts about 15% increase in Nu at the downstream part of the tube. This is due to the inaccurate prediction of buoyancy-influenced region as discussed before. The large discrepancies between the numerical results using various turbulence models and experimental data are due to the differences in response to the modification of flow and turbulence fields induced by the buoyancy effect. Fig. 7 shows the fluid temperature, isobaric specific heat capacity, and density profiles at several axial locations obtained using the LS k-e and SST k-x model in

(b) SST k-ω model

Fig. 7. Profiles of fluid temperature, isobaric specific heat capacity, and density at different axial cross sections (Case B1).

866

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

Case B1. The dimensionless wall distance y+ is calculated using the semi-local scaling form [51]:

yþ ¼

yðsw =qÞ1=2

v

ð33Þ

As shown in Fig. 7, no abnormal distributions of fluid temperature, isobaric specific heat capacity and density are predicted using both models. At various axial locations, the largest temperature gradient remains to be located inside the laminar sub-layer (y+ < 5). That is, the thickness of thermal boundary layer remains to be thin from the tube inlet to the outlet. The location where fluid property varies the most is near the wall, which is identical with convective heat transfer of normal fluids. Figs. 8 and 9 shows the mean axial velocity, dimensionless axial velocity and turbulent kinetic energy profiles at several axial locations in Case B1 using the LS k-e and SST k-x model, respectively. Using either model, the predicted mean velocity profile shows a typical one for turbulent boundary layer. At x/d = 100, a minor reduction of turbulent kinetic energy in the core region of flow is predicted by the LS k-e model, which corresponds to the small peak of local wall temperature as shown in Fig. 5(a). In comparison, the SST k-x model fails to predict the reduction of turbulent kinetic energy and hence the rapid increase of wall temperature. When qw/G is larger, more fluid is heated to exceed the pseudo critical temperature (Tpc = 706 K at p = 3.5 MPa), the drastic variation of fluid property leads to a stronger buoyancy effect. The differences in thermal, flow and turbulence fields predicted by the two models are more obvious. Fig. 10 shows the fluid temperature, isobaric specific heat capacity and density profiles at several axial locations obtained using the LS k-e and SST k-x models in Case B2. A large gradient in fluid temperature is predicted by the LS k-e model from x/d = 35 to x/d = 300, where the deteriorated heat transfer occurs. This large temperature gradient moves to a further distance from the wall as the fluid travels downstream. It could be also observed that the region within which the fluid temperature is

near Tpc is located in the buffer layer (5 < y+ < 30). According to He et al. [52], such is referred to as the large-property-variation (LPV) region. Large density gradient is felt at a location away from the wall as the LPV region enters the buffer layer, indicating the buoyancy begins to have significant influences on the flow and turbulence fields. This will be discussed later in this section. The predicted temperature gradient using the SST k-x model is smaller than that by the LS k-e model from location x/d = 35 to x/d = 300. The corresponding density gradient is also smaller and remains close the wall. As a result, the effect of buoyancy on the flow and turbulence fields is also weaker. Figs. 11 and 12 shows the mean axial velocity, dimensionless axial velocity and turbulent kinetic energy profiles at several axial locations in Case B2 using the LS k-e and SST k-x model, respectively. As can be seen from Fig. 11(1), the axial velocity profile predicted by the LS k-e model becomes flattened at location x/d = 100, indicating the laminarization of flow. At further downstream locations x/d = 200 and x/d = 300, the strong effect of buoyancy induces an overshoot of velocity close to the wall. The modification of mean velocity field is better illustrated in Fig. 11(2). The velocity overshoot first appears at y+10 and moves towards the core region of the flow as the fluid travels downstream. Due to the buoyancy effect, the turbulent kinetic energy continuously decreases along the flow direction and almost diminishes at location x/d = 200 and x/d = 300. As discussed before, the LPV region moves to the buffer layer at these axial locations. In the buffer layer, the production of turbulence kinetic energy reaches its peak [53], the effect of buoyancy strongly reduces the turbulence production and thus results a significant reduction of turbulent kinetic energy and thus deteriorated heat transfer. This corresponds to the rapid increase of wall temperature from x/d = 35 to x/d = 300 predicted by the LS k-e model as shown in Fig. 5(b). At x/d = 400 and x/d = 448, the velocity profile begins to return to the normal state. The turbulence intensity also begins to recover, corresponding to the decrease of wall temperature.

Fig. 8. Profiles of mean flow and turbulent kinetic energy at different axial cross sections (Case B1, LS k-e model).

Fig. 9. Profiles of mean flow and turbulent kinetic energy at different axial cross sections (Case B1, SST k-x model).

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

(a) LS k-ε model

(b) SST k-ω model

Fig. 10. Profiles of fluid temperature, isobaric specific heat capacity, and density at different axial cross sections (Case B2).

Fig. 11. Profiles of mean flow and turbulent kinetic energy at different axial cross sections (Case B2, LS k-e model).

867

868

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

Fig. 12. Profiles of mean flow and turbulent kinetic energy at different axial cross sections (Case B2, SST k-x model).

Fig. 13. Comparison of turbulent kinetic energy production due to shear and buoyancy (Case B2, LS k-e model).

(a) Shear production term

(b) Buoyancy production term

Fig. 14. Comparison of turbulent kinetic energy production due to shear and buoyancy (Case B2, SST k-x model).

In comparison, no velocity overshoot is predicted by the SST k-

3.3. Effect of buoyancy on turbulence production

x model as shown in Fig. 12. The SST k-x model also predicts that the turbulent kinetic energy decreases from location x/d = 35 to x/ d = 300. However, the reduction in turbulent kinetic energy is less severe comparing with the results obtained using the LS k-e model. This is in accordance with the fact that no obvious deteriorated heat transfer is predicted using the SST k-x model (see Fig. 5(b)).

As discussed in Section 3.2, the buoyancy effect has significant influences on the flow and turbulence fields, which may lead to impairment of heat transfer under mixed convection conditions. Particularly, the buoyancy force can affect the turbulence production in two ways: the indirect (external) effect and the direct

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

(a) LS k-ε model

869

(b) SST k-ω model

Fig. 15. Comparison of predicted inner wall temperature using various constant turbulent Prandtl number (Case A1).

(a) LS k-ε model

(b) SST k-ω model

Fig. 16. Comparison of predicted inner wall temperature using various constant turbulent Prandtl number (Case B1).

(structural) effect [52]. The indirect effect means that the turbulence production (by shear) is affected by the modification of mean velocity fields. Whereas the direct effect is caused by density fluctuations which leads to the buoyant production of turbulent kinetic energy [54]. In order to determine which effect is the primary cause of the reduction in turbulence production for deteriorated regime of heat transfer, Case B2 is studied in the following discussion. Figs. 13 and 14 show the comparison of turbulent kinetic energy production due to shear and buoyancy in Case B2 obtained using the LS k-e model and the SST k-x model, respectively. Both models predict negative values of Gk at most of the axial locations, indicating a weakening effect on turbulence production. However, it is observed that the magnitude of Gk is much small than the reduction in Pk. Therefore, it is concluded that the indirect effect of buoyancy on turbulence production is more important under the mixed convection heat transfer conditions. It could be also deduced that the predicted heat transfer is highly dependent on the performance of turbulence models in simulating the mean velocity fields. 3.4. Effect of the turbulent Prandtl number Apart from the turbulence models, the choice of turbulent Prandtl number may also influence the simulation results of convective heat transfer to fluids at supercritical pressures. To

investigate the effect of turbulent Prandtl number on predicted heat transfer, all of the simulation cases are revisited using various constant value of rt : 0.85, 0.9, and 1. In addition, the LS k-e model and SST k-x model are used to examine whether the effect of turbulent Prandtl number varies with different turbulence models. As examples, the predicted inner wall temperature in Case A1 and Case B1 is shown in Figs. 15 and 16, respectively. The results show that the discrepancies of predicted inner wall temperature using different turbulent Prandtl number are more dependent on the flow conditions other than the turbulence models. The predicted inner wall temperature becomes lower when a smaller rt is used. This could be explained that with the decrease of rt , the turbulent heat flux increases and therefore enhances the heat transfer from the wall to the core region of flow. Under all of the conditions considered in the current study, using a constant value of 1.0 for rt leads to a closest agreement with the experimental data. It is pointed out that several studies indicate that turbulent Prandtl number is likely a function of fluid-thermal variables and fluid properties [51]. Therefore, using constant turbulent Prandtl number in the simulations may not be the optimal choice. This is particularly true in mixed convection conditions, where the effect of fluid property variation on the turbulent Prandtl number should be carefully considered. More physical-based models for turbulent heat flux could help improve the accuracy of simulations.

870

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871

4. Conclusions In the study, numerical investigation of convective heat transfer to China RP-3 aviation kerosene at supercritical pressures is carried out using an in-house code developed using OpenFOAM. A variety of RANS turbulence models are used to simulate heat transfer under both forced and mixed convection conditions. The following conclusions may be derived based on the results from this study. (1) Under forced convection conditions, all turbulence models predict a gentler growth of wall temperature along the tube compared with the experimental data. The SST k-x model performs best in predicting the wall temperature distribution, and the corresponding predicted Nusselt number are close to that calculated using the correlation for forced convection proposed by Hu et al. (2) Under mixed convection conditions, turbulence models using different near-wall treatments respond differently to the buoyancy effect on flow and heat transfer, which leads to large differences in predicted thermal, velocity and turbulence fields. Only the LRN k-e models (LS k-e and CH k-e model) are able to qualitatively predict the heat transfer deterioration. In general, the LS k-e model performs better than the CH k-e model. The SST k-x, v2-f, and RNG k-e model (using wall functions) fails to predict the impairment of heat transfer at the entrance region. (3) In the deteriorated heat transfer regime, flow laminarization is mainly due to the modification of mean velocity fields, which results in significant reduction of turbulent kinetic energy. The direct effect of buoyancy on turbulence production (density fluctuations) is negligible compared with the indirect effect (modification of shear). Accurate simulation of the modification in velocity fields is important to obtain an accurate prediction of heat transfer under mixed convection conditions. (4) The value of turbulent Prandtl number needs to be carefully considered in conducting simulation of convective heat transfer to aviation kerosene at supercritical pressures. The predicted wall temperature becomes lower when smaller value of turbulent Prandtl number is used in simulation. Under the conditions considered in the current study, the value of 1.0 for turbulent Prandtl number leads to a closest agreement with the experimental data.

Acknowledgements This work was funded by National Natural Science Foundation of China (No. 51576027), the financial support is gratefully acknowledged. Conflict of interest The authors declare that they have no conflict of interests. References [1] W.H. Heiser, D.T. Pratt, Hypersonic Airbreathing Propulsion, American Institute of Aeronautics and Astronautics, Washington DC, 1994, pp. 2–26. [2] R.S. Fry, A century of ramjet propulsion technology evolution, J. Propul. Power 20 (1) (2004) 27–58, https://doi.org/10.2514/1.9178. [3] C. Segal, The Scramjet Engine: Processes and Characteristics, Cambridge University Press, 2009, p. 2. [4] H. Huang, L.J. Spadaccini, D.R. Sobel, Fuel-cooled thermal management for advanced aero engines, J. Eng. Gas Turb. Power 126 (2) (2002) 284–293, https://doi.org/10.1115/1.1689361.

[5] D.R. Sobel, L.J. Spadaccini, Hydrocarbon fuel cooling technologies for advanced propulsion, J. Eng. Gas Turb. Power 119 (2) (1997) 344–351, https://doi.org/ 10.1115/1.2815581. [6] H. Huang, L. Spadaccini, D. Sobel, Endothermic heat-sink of jet fuels for scramjet cooling, in: 38th AIAA/ASME/SAE/ASEE Joint Propulsion Conference & Exhibit, American Institute of Aeronautics and Astronautics, Indianapolis, Indiana, USA, 2002, 10.2514/6.2002-3871. [7] M. Bouchez, G. Cahuzac, S. Beyer, V. Avrashkov, PTAH-SOCAR fuel-cooled composite materials structure in 2003, in: 12th AIAA International Space Planes and Hypersonic Systems and Technologies, American Institute of Aeronautics and Astronautics, Norfolk, Virginia, USA, 2003, 10.2514/6.20036919. [8] T. Edwards, Liquid fuels and propellants for aerospace propulsion: 1903–2003, J. Propul. Power 19 (6) (2003) 1089–1107, https://doi.org/10.2514/2.6946. [9] I.W. Kay, W.T. Peschke, R.N. Guile, Hydrocarbon-fueled scramjet combustor investigation, J. Propul. Power 8 (2) (1992) 507–512, https://doi.org/10.2514/ 3.23505. [10] O.A. Powell, J.T. Edwards, R.B. Norris, K.E. Numbers, J.A. Pearce, Development of hydrocarbon-fueled scramjet engines: the hypersonic technology (Hytech) program, J. Propul. Power 17 (6) (2001) 1170–1176, https://doi.org/10.2514/ 2.5891. [11] A. Boudreau, Hypersonic air-breathing propulsion efforts in the air force research laboratory, in: AIAA/CIRA 13th International Space Planes and Hypersonics Systems and Technologies Conference, American Institute of Aeronautics and Astronautics, Reston, Virgina, USA, 10.2514/6.2005-3255. [12] A. Siebenhaar, F. Chen, M. Karpuk, B. Hitch, T. Edwards, Engineering scale titanium endothermic fuel reactor demonstration for hypersonic scramjet engine, in: 9th International Space Planes and Hypersonic Systems and Technologies Conference, American Institute of Aeronautics and Astronautics, Norfolk, Virgina, USA, 1999, 10.2514/6.1999-4909. [13] F. Zhong, X. Fan, G. Yu, J. Li, C.J. Sung, Heat transfer of aviation kerosene at supercritical conditions, J. Thermophys. Heat Transf. 23 (3) (2009) 543–550, https://doi.org/10.2514/1.41619. [14] J. Yan, Z. Liu, Q. Bi, Y. Guo, Z. Yang, Heat transfer of hydrocarbon fuel under steady states and pressure-transient states, J. Propul. Power 32 (1) (2015) 38– 45, https://doi.org/10.2514/1.B35694. [15] T. Edwards, S. Zabarnick, Supercritical fuel deposition mechanisms, Ind. Eng. Chem. Res. 32 (12) (1993) 3117–3122, https://doi.org/10.1021/ie00024a022. [16] S. He, W.S. Kim, J.H. Bae, Assessment of performance of turbulence models in predicting supercritical pressure heat transfer in a vertical tube, Int. J. Heat Mass Transf. 51 (19–20) (2008) 4659–4675, https://doi.org/10.1016/j. ijheatmasstransfer. 2007.12.028. [17] B. Hitch, M. Karpuk, Experimental investigation of heat transfer and flow instabilities in supercritical fuels, in: 33rd Joint Propulsion Conference and Exhibit, American Institute of Aeronautics and Astronautics, Seattle, Washington, USA, 1997, 10.2514/6.1997-3043. [18] Z. Hu, T. Chen, Y. Luo, J. Zheng, M. Tang, Experimental research on the heattransfer characteristics of kerosene under supercritical pressure, Hsi-An Chiao Tung Ta Hsueh/J. Xi’an Jiaotong Univ. 33 (9) (1999) 62–65, 70. [19] C. Zhang, G. Xu, L. Gao, Z. Tao, H. Deng, K. Zhu, Experimental investigation on heat transfer of a specific fuel (RP-3) flows through downward tubes at supercritical pressure, J. Supercrit. Fluids 72 (9) (2012) 90–99, https://doi.org/ 10.1016/j.supflu. 2012.07.011. [20] Z. Liu, Q. Bi, Y. Guo, J. Yan, Z. Yang, Convective heat transfer and pressure drop characteristics of near-critical-pressure hydrocarbon fuel in a minichannel, Appl. Therm. Eng. 51 (1–2) (2013) 1047–1054, https://doi.org/10.1016/j. applthermaleng. 2012.10.029. [21] W. Li, D. Huang, G.Q. Xu, Z. Tao, Z. Wu, H.T. Zhu, Heat transfer to aviation kerosene flowing upward in smooth tubes at supercritical pressures, Int. J. Heat Mass Transf. 85 (2015) 1084–1094, https://doi.org/10.1016/j. ijheatmasstransfer. 2015.01.079. [22] W. Chen, X. Fang, Modeling of convective heat transfer of RP-3 aviation kerosene in vertical miniature tubes under supercritical pressure, Int. J. Heat Mass Transf. 95 (2016) 272–277, https://doi.org/10.1016/j. ijheatmasstransfer.2015.12.018. [23] Y.X. Hua, Y.Z. Wang, H. Meng, A numerical study of supercritical forced convective heat transfer of n-heptane inside a horizontal miniature tube, J. Supercrit. Fluids 52 (1) (2010) 36–46, https://doi.org/10.1016/ j.supflu.2009.12.003. [24] J. Zhu, Z. Tao, H. Deng, K. Wang, X. Yu, Numerical investigation of heat transfer characteristics and flow resistance of kerosene RP-3 under supercritical pressure, Int. J. Heat Mass Transf. 91 (2015) 330–341, https://doi.org/ 10.1016/j.ijheatmasstransfer.2015.07.118. [25] G. Dang, F. Zhong, L. Chen, X. Chang, Numerical investigation on flow and convective heat transfer of aviation kerosene at supercritical conditions, Sci. China Technol. Sci. 56 (2) (2013) 416–422, https://doi.org/10.1007/s11431012-5075-3. [26] Z. Liu, J. Liang, Y. Pan, Numerical analysis of heat transfer deterioration of China RP-3 aviation kerosene in a circular tube at supercritical pressures, in: 11th AIAA/ASME Joint Thermophysics and Heat Transfer Conference, American Institute of Aeronautics and Astronautics, Atlanta, Georgia, U.S.A., 2014, 10.2514/6.2014-3358. [27] H. Jiang, J. Ervin, Z. West, S. Zabarnick, Turbulent flow, heat transfer deterioration, and thermal oxidation of jet fuel, J. Thermophys. Heat Transf. 27 (4) (2013) 668–678, https://doi.org/10.2514/1.T4103.

H. Pu et al. / International Journal of Heat and Mass Transfer 118 (2018) 857–871 [28] Z. Tao, Z. Cheng, J. Zhu, H. Li, Effect of turbulence models on predicting convective heat transfer to hydrocarbon fuel at supercritical pressure, Chin. J. Aeronaut. 29 (5) (2016) 1247–1261, https://doi.org/10.1016/j.cja.2016.08.007. [29] C. Dang, E. Hihara, In-tube cooling heat transfer of supercritical carbon dioxide. Part 2. Comparison of numerical calculation with different turbulence models, Int. J. Refrig 27 (7) (2004) 748–760, https://doi.org/10.1016/j. ijrefrig.2004.04.017. [30] M. Mohseni, M. Bazargan, Effect of turbulent Prandtl number on convective heat transfer to turbulent flow of a supercritical fluid in a vertical round tube, J. Heat Transf. 133 (7) (2011), https://doi.org/10.1115/1.4003570. [31] H.G. Weller, G. Tabor, H. Jasak, C. Fureby, A tensorial approach to computational continuum mechanics using object-oriented techniques, Comput. in Phys. 12 (6) (1998) 620–631, https://doi.org/10.1063/1.168744. [32] OpenFOAM Foundation Ltd., OpenFOAM user guide, 2015. [33] Q.M. Sun, Z.T. Mi, X.W. Zhang, Determination of critical properties (tc, pc) of endothermic hydrocarbon fuels-RP-3 and simulated JP-7, Ranliao Huaxue Xuebao/J. Fuel Chem. Technol. 34 (4) (2006) 466–470. [34] B.E. Launder, B.I. Sharma, Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc, Lett. Heat Mass Transf. 1 (2) (1974) 131–137, https://doi.org/10.1016/0094-4548(74)90150-7. [35] K.Y. Chien, Predictions of channel and boundary-layer flows with a lowReynolds-number turbulence model, AIAA J. 20 (1) (1982) 33–38, https://doi. org/10.2514/3.51043. [36] F.S. Lien, G. Kalitzin, Computations of transonic flow with the v2–f turbulence model, Int. J. Heat Fluid Flow 22 (1) (2001) 53–61, https://doi.org/10.1016/ S0142-727X(00)00073-4. [37] F.R. Menter, Two-equation eddy-viscosity turbulence models for engineering applications, AIAA J. 32 (8) (1994) 1598–1605, https://doi.org/10.2514/ 3.12149. [38] V. Yakhot, S.A. Orszag, S. Thangam, T.B. Gatski, C.G. Speziale, Development of turbulence models for shear flows by a double expansion technique, Phys. Fluids A 4 (7) (1992) 1510–1520, https://doi.org/10.1063/1.858424. [39] P.A. Durbin, Near-wall turbulence closure modeling without ‘‘damping functions”, Theoret. Comput. Fluid Dyn. 3 (1) (1991) 1–13, https://doi.org/ 10.1007/bf00271513. [40] P.A. Durbin, Application of a near-wall turbulence model to boundary layers and heat transfer, Int. J. Heat Fluid Flow 14 (4) (1993) 316–323, https://doi. org/10.1016/0142-727X(93)90004-7. [41] G. Kalitzin, G. Medic, G. Iaccarino, P. Durbin, Near-wall behavior of RANS turbulence models and implications for wall functions, J. Comput. Phys. 204 (1) (2005) 265–291, https://doi.org/10.1016/j.jcp.2004.10.018. [42] H. Pu, S. Li, M. Dong, S. Jiao, Y. Shang Development and validation of a conjugate heat transfer code for regenerative cooling structure with

[43]

[44]

[45]

[46]

[47]

[48]

[49]

[50]

[51]

[52]

[53] [54]

871

openfoam, in: 21st AIAA International Space Planes and Hypersonics Technologies Conference, American Institute of Aeronautics and Astronautics, Xiamen, China, 2017. 10.2514/6.2017-2103. R.I.I. Paulo, J. Oliveira, An improved PISO algorithm for the computation of buoyancy-driven flows, Numer. Heat Transf. Part B: Fundam. 40 (6) (2001) 473–493, https://doi.org/10.1080/104077901753306601. H.W. Deng, C.B. Zhang, G.Q. Xu, Z. Tao, B. Zhang, G.Z. Liu, Density measurements of endothermic hydrocarbon fuel at sub- and supercritical conditions, J. Chem. Eng. Data 56 (6) (2011) 2980–2986, https://doi.org/ 10.1021/je200258g. T. Edwards, L.Q. Maurice, Surrogate mixtures to represent complex aviation and rocket fuels, J. Propul. Power 17 (2) (2001) 461–466, https://doi.org/ 10.2514/2.5765. M.L. Huber, NIST Thermophysical Properties of Hydrocarbon Mixtures Database (SUPERTRAPP) Users’ Guide, National Institute of Standards and Technology, Gaithersburg, 2007. B. Liu, Y. Zhu, J.J. Yan, Y. Lei, B. Zhang, P.X. Jiang, Experimental investigation of convection heat transfer of n-decane at supercritical pressures in small vertical tubes, Int. J. Heat Mass Transf. 91 (2015) 734–746, https://doi.org/10.1016/j. ijheatmasstransfer.2015.07.006. H.D. Murphy, F.W. Chambers, D.M. Mceligot, Laterally converging flow. Part 1. Mean flow, J. Fluid Mech. 127 (-1) (1983) 379–401, https://doi.org/10.1017/ S0022112083002785. W.S. Kim, S. He, J.D. Jackson, Assessment by comparison with DNS data of turbulence models used in simulations of mixed convection, Int. J. Heat Mass Transf. 51 (5–6) (2008) 1293–1312, https://doi.org/10.1016/j. ijheatmasstransfer. 2007.12.002. H. Deng, K. Zhu, G. Xu, Z. Tao, J. Sun Heat transfer characteristics of RP-3 kerosene at supercritical pressure in a vertical circular tube, J. Enhanced Heat Transf. 19 (5) (2012) 409–421, https://doi.org/10.1615/JEnhHeatTransf. 2012004966. Y.Y. Bae, E.S. Kim, M. Kim, Numerical simulation of supercritical pressure fluids with property-dependent turbulent Prandtl number and variable damping function, Int. J. Heat Mass Transf. 101 (2016) 488–501, https://doi.org/ 10.1016/j.ijheatmasstransfer.2016.05.051. S. He, W.S. Kim, J.D. Jackson, A computational study of convective heat transfer to carbon dioxide at a pressure just above the critical value, Appl. Therm. Eng. 28 (13) (2008) 1662–1675, https://doi.org/10.1016/j. applthermaleng.2007.11.001. S.B. Pope, Turbulent Flows, Cambridge University Press, New York, 2000, p. 285. G.F. Hewitt, J.C. Vassilicos, Prediction of Turbulent Flows, Cambridge University Press, 2005, p. 65.