2.1 Governing equations

The governing equations of mass, momentum and total energy for a fully compressible flow are solved in conservative form which reads

(2.1a ) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70C}}{\unicode[STIX]{x2202}t}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70C}u_{j}}{\unicode[STIX]{x2202}x_{j}}=0, & \displaystyle\end{eqnarray}$$

(2.1b ) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70C}u_{i}}{\unicode[STIX]{x2202}t}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70C}u_{i}u_{j}}{\unicode[STIX]{x2202}x_{j}}=-\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}x_{j}}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70F}_{ij}}{\unicode[STIX]{x2202}x_{j}}, & \displaystyle\end{eqnarray}$$

(2.1c ) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D70C}E}{\unicode[STIX]{x2202}t}+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x_{j}}[u_{j}(\unicode[STIX]{x1D70C}E+p)]=\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x_{j}}(u_{i}\unicode[STIX]{x1D70F}_{ij}-q_{j}), & \displaystyle\end{eqnarray}$$

$x_{1}$

$x_{2}$

$x_{3}$

$x$

$y$

$z$

$u_{i}$

$i$

$t$

$\unicode[STIX]{x1D70C}$

$p$

$E$

The viscous and conductive heat fluxes in (2.1b ) and (2.1c ) are, respectively,

(2.2a ) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D70F}_{ij}=2\unicode[STIX]{x1D707}\left[\unicode[STIX]{x1D61A}_{ij}-\frac{1}{3}\frac{\unicode[STIX]{x2202}u_{k}}{\unicode[STIX]{x2202}x_{k}}\unicode[STIX]{x1D6FF}_{ij}\right], & \displaystyle\end{eqnarray}$$

(2.2b ) $$\begin{eqnarray}\displaystyle & \displaystyle q_{j}=-\unicode[STIX]{x1D706}\frac{\unicode[STIX]{x2202}T}{\unicode[STIX]{x2202}x_{j}}=-\frac{c_{p}\unicode[STIX]{x1D707}}{Pr}\frac{\unicode[STIX]{x2202}T}{\unicode[STIX]{x2202}x_{j}}, & \displaystyle\end{eqnarray}$$

$\unicode[STIX]{x1D707}$

$\unicode[STIX]{x1D61A}_{ij}$

$\unicode[STIX]{x1D61A}_{ij}=(\unicode[STIX]{x2202}u_{j}/\unicode[STIX]{x2202}x_{i}+\unicode[STIX]{x2202}u_{i}/\unicode[STIX]{x2202}x_{j})/2$

$\unicode[STIX]{x1D706}$

$c_{p}$

$Pr$

$T$

2.2 Real fluid model

The Peng–Robinson (PR) equation of state (EoS) (Peng & Robinson Reference Peng and Robinson1976) is used to model the working fluid of choice for this study, R-134a (1,1,1,2-tetrafluoroethane, $\text{CH}_{2}\text{FCF}_{3}$ ), which benefits from experimentally accessible critical pressures and temperatures of $p_{cr}=40.59$  bar and $T_{cr}=374.26$  K and is widely used in turbomachinery and heat exchangers as a cooling fluid because of its non-toxic and non-flammable characteristics. Departure functions guaranteeing full thermodynamic consistency with the chosen EoS have been derived following Ewing & Peters (Reference Ewing and Peters2000). Transport properties such as the dynamic viscosity and thermal conductivity are estimated via Chung’s method (Chung et al. Reference Chung, Ajlan, Lee and Starling1988) which predicts experimental values within 5 % error (Poling, Prausnitz & O’Connell Reference Poling, Prausnitz and O’Connell2001). The choice of an accurate and simple EoS such as the PR EoS provides a consistent thermodynamic model, computationally less expensive than interpolating tabulated values. In order to prove adequacy and accuracy of implementation of the equations, detailed derivations and comparisons against the NIST database (Lemmon, McLinden & Friend Reference Lemmon, McLinden and Friend2016) (figure 27) are included in appendix A.

2.3 Computational set-up

The proposed numerical simulations have been carried out with Hybrid, a fully compressible Navier–Stokes solver originally written by Johan Larsson. This code utilizes a finite central difference scheme with a fourth-order accuracy by summation-by-parts operators for the inviscid terms and a second-order accuracy for the viscous terms. The time advancement is achieved by a fourth-order-accurate Runge–Kutta method. Hybrid has been used in several canonical numerical investigations such as those involving shock–vortex interaction, compressible homogeneous isotropic turbulence (Larsson, Lele & Moin Reference Larsson, Lele and Moin2007) and shock–turbulence interaction (Larsson & Lele Reference Larsson and Lele2009; Larsson, Bermejo-Moreno & Lele Reference Larsson, Bermejo-Moreno and Lele2013). The code solves single-component fluid, which is a suitable modelling approach for a supercritical pressure flow since surface tension becomes negligible for supercritical pressures, $p>p_{cr}$ , and numerical techniques for multiphase simulations such as interface tracking or reconstruction are not required. New features that have been added to the code include parallel HDF5 (The HDF Group 1998) input/output capabilities and a generic EoS.

The computational setup is a compressible turbulent channel flow (figure 2) kept at a nominal bulk pressure of $p_{b}\simeq 1.1p_{cr}$ corresponding to a pseudoboiling temperature of $T_{pb}=379.1$  K defined based on the maximum isobaric thermal expansion coefficient, $\unicode[STIX]{x1D6FC}_{p}=-(\unicode[STIX]{x2202}\unicode[STIX]{x1D70C}/\unicode[STIX]{x2202}T)_{p}/\unicode[STIX]{x1D70C}$ (figure 1 a). The assigned isothermal top- and bottom-wall boundary conditions bracket the pseudoboiling temperature ( $T_{top/bot}=T_{pb}\pm \unicode[STIX]{x0394}T/2$ ) maintaining transcritical temperature conditions (figure1 b).

Figure 1.

Phase diagram for R-134a showing the critical point ( $p_{cr}=40.59$  bar,  $T_{cr}=374.26$  K) (▫), the PBL (- - -) and the isolines of isobaric thermal expansion coefficient, $\unicode[STIX]{x1D6FC}_{p}$ (——, $\text{K}^{-1}$ ) (a). Density and isobaric heat capacity versus temperature for $p=1.1p_{cr}$ with the pseudoboiling point (●) and top-to-bottom temperature differences, $\unicode[STIX]{x0394}T$ , bracketing $T_{pb}=379.1$  K (b).

Figure 2.

Computational set-up for supercritical compressible channel flow simulations in transcritical temperature conditions. Simulation parameters are given in tables 1 and 2.

Table 1.

Simulation parameters achieving transcritical temperature conditions for R-134a. Bulk parameters are indicated with a subscript ‘ $b$ ’, while pseudoboiling or pseudo(phase)transitioning are indicated with ‘ $pb$ ’.

Top- to bottom-wall temperature differences investigated are $\unicode[STIX]{x0394}T=T_{top}-T_{bot}=5$ , 10 and 20 K, with bulk density set to $\unicode[STIX]{x1D70C}_{b}=450$ , $474$ and $520~\text{kg}~\text{m}^{-3}$ , respectively, determined via trial and error to obtain the desired bulk pressure for all cases (see tables 1 and 2). Periodic boundary conditions are applied to the streamwise and spanwise directions and the grid is stretched in the wall-normal direction with a hyperbolic tangent law. To guarantee feasibility of the simulations on the finest grid and the highest temperature difference considered where the time step is acoustically limited to $\unicode[STIX]{x0394}t=1.4\times 10^{-8}$ , the bulk velocity has been set for all cases to the relatively high value (for typical heat transfer applications) of $U_{b}=36~\text{m}~\text{s}^{-1}$ corresponding to a Mach number in the low-subsonic range of $M_{b}=0.26$ with a range of turbulent Mach number $M_{t}=0.015$ (centre region) to 0.051 (near-wall peak). Given the large density variation near the pseudocritical point, buoyancy effects may be important in the mean as well as in the turbulent quantities. In this study, however, the buoyancy effects are neglected in order to focus on structural changes in compressible channel flow turbulence due to wall heat transfer in the presence of real fluid effects.

Table 2.

Friction Reynolds number and grid resolution in wall units $(u_{\unicode[STIX]{x1D70F}}/\unicode[STIX]{x1D708})^{-1}$ for the bottom and top portion of the channel evaluated with respective wall quantities. See also table 1.

To assess the sensitivity of the flow to the thermodynamic gradients near the critical point, a reference simulation at twice the critical pressure ( $p_{b}=2p_{cr}$ ) is carried out with the same working fluid (R-134a) and simulation parameters below. Given the higher pressure, the gradients at the pseudoboiling point are much weaker compared to the near-critical conditions $p_{b}=1.1p_{cr}$ chosen for the rest of the runs (see comparison in figure 27 in appendix A). It should be noted that some real fluid effects are still present in the higher-pressure case; the thermophysical and thermodynamic variations as a function of temperature are in fact still present, albeit much weaker. For example, for the same $\unicode[STIX]{x0394}T=20$  K, the relative density difference between the top and bottom walls is only 25 % at $p_{b}=2p_{cr}$ ; it is 62 % at $p_{b}=1.1p_{cr}$ (see table 3). The parameters of this reference simulation are $\unicode[STIX]{x1D70C}_{b}=596~\text{kg}~\text{m}^{-3}$ , $p_{b}=2p_{cr}=81.18$  bar, $\unicode[STIX]{x0394}T=20$  K ( $T_{bot}=394.9$  K and $T_{top}=414.9$  K where the pseudoboiling temperature is $T_{pb}=404.9$  K at the given pressure) and $U_{b}=36~\text{m}~\text{s}^{-1}$ .

To ensure the proper spatial resolution of all relevant hydrodynamic and thermodynamic scales, a systematic grid refinement study has been carried out (see appendix B and table 2); this is especially important in simulations of supercritical flows in the near-critical or pseudophase transitioning conditions (see Introduction). The relevant metric of the spectral broadening level for channel flow turbulence is the friction Reynolds number,

(2.3) $$\begin{eqnarray}\displaystyle Re_{\unicode[STIX]{x1D70F}}=\frac{u_{\unicode[STIX]{x1D70F}}h}{\unicode[STIX]{x1D708}_{w}}, & & \displaystyle\end{eqnarray}$$

based on the friction velocity, $u_{\unicode[STIX]{x1D70F}}$ , the channel half-height, $h$ , and the kinematic viscosity at the wall, $\unicode[STIX]{x1D708}_{w}$ , of the fluid. It can be viewed as the channel half-height normalized by the viscous length scale, $\unicode[STIX]{x1D708}_{w}/u_{\unicode[STIX]{x1D70F}}=\unicode[STIX]{x1D708}_{w}/(\unicode[STIX]{x2202}u/\unicode[STIX]{x2202}x_{2})_{x_{2}=0}$ , hence $Re_{\unicode[STIX]{x1D70F}}=h^{+}$ . Therefore, $Re_{\unicode[STIX]{x1D70F}}$ is the ratio of an integral length scale, ${\sim}h$ , to a viscous scale evaluated at the wall. Typical practice in DNS is to adopt relatively low values of friction Reynolds number to enable full resolution of the relevant scales. For the present simulations, this is achieved by augmenting dynamic viscosity and thermal conductivity by the same multiplicative factor (figure 3) resulting in $Re_{\unicode[STIX]{x1D70F}}$ in the range of 342–394 (table 2). This choice leaves the Prandtl number unaltered and reproduces the correct trend of transport properties in the transcritical regime (see appendix A). The reference simulation has been carried out at $Re_{\unicode[STIX]{x1D70F}}=366$ , $\unicode[STIX]{x0394}x^{+}=10.98$ , $\unicode[STIX]{x0394}y^{+}=0.43$ –5.93, $\unicode[STIX]{x0394}z^{+}=5.63$ for the bottom wall and $Re_{\unicode[STIX]{x1D70F}}=386$ , $\unicode[STIX]{x0394}x^{+}=11.58$ , $\unicode[STIX]{x0394}y^{+}=0.46$ –6.25, $\unicode[STIX]{x0394}z^{+}=5.94$ for the top wall.

Figure 3.

Dynamic viscosity $\unicode[STIX]{x1D707}$ , thermal conductivity $\unicode[STIX]{x1D706}$ and Prandtl number $Pr$ for R-134a taken from Chung’s model (——) (see appendix A); scaled dynamic viscosity and conductivity (- - -), augmented by a factor of 60, used in the computations, yielding the same Prandtl number.

Table 3.

Top- and bottom-wall values of mean density and compressibility factor and average location of pseudophase transition $y_{pb}$ for various temperature conditions. With the exception of $\unicode[STIX]{x0394}T$ , all values reported are a result of the calculations. The first three rows are related to the $p_{b}=1.1p_{cr}$ cases, the last row to the reference high-pressure case $p_{b}=2p_{cr}$ .

3.1 Mean flow quantities

Figure 4.

Reynolds-averaged density (a), temperature (b) and compressibility factor (c) for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ) and reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (●) (see § 2.3). Average location of pseudotransition for $\unicode[STIX]{x0394}T=5~\text{K}$ (○), 10 K (▵) and 20 K (▫).

Figure 4 shows Reynolds-averaged profiles of density, temperature and compressibility factor $Z=p/(\unicode[STIX]{x1D70C}R_{gas}T)$ , where $R_{gas}=81.49~\text{J}~\text{kg}^{-1}~\text{K}^{-1}$ is the gas constant for R-134a. The top-to-bottom density difference (table 3) of $\unicode[STIX]{x0394}\unicode[STIX]{x1D70C}=447.5~\text{kg}~\text{m}^{-3}$ achieved under transcritical conditions at $1.1p_{cr}$ for $\unicode[STIX]{x0394}T=20~\text{K}$ is about 2.67 times higher than $\unicode[STIX]{x0394}\unicode[STIX]{x1D70C}=167.4~\text{kg}~\text{m}^{-3}$ , obtained at twice the critical pressure for the same $\unicode[STIX]{x0394}T$ . Remarkably, as shown later, density fluctuation intensities for $\unicode[STIX]{x0394}T=20$  K and $2p_{cr}$ are very similar to those observed for $\unicode[STIX]{x0394}T=5$  K and $1.1p_{cr}$ ; in the latter case, however, higher temperature fluctuations are obtained. As expected for higher pressure conditions (see figure 27), a gradual change is also observed in the compressibility factor for the reference data analogous to the density.

As $\unicode[STIX]{x0394}T$ is increased, the average location of pseudotransition $y_{pb}$ , where the real fluid effects are expected to be the most accentuated, moves from a near-centreplane location to the upper wall. The same holds when keeping $\unicode[STIX]{x0394}T$ constant and lowering the base pressure. It is important to recall that the isothermal wall conditions are selected to be exactly $\unicode[STIX]{x0394}T/2$ warmer (top) and colder (bottom) than the pseudoboiling temperature; yet the location of the pseudotransition in the channel is not located at the centreplane and it is an output of the calculation. The shift of the pseudotransition location is related to the highly nonlinear and asymmetric thermophysical properties of the fluid about the pseudoboiling point, especially the specific heat capacity. In the transcritical regime, the pseudophase change is accompanied by a finite peak in the specific heat capacity which acts as a thermal barrier. We make an approximation that the enthalpy difference from the cold wall to the pseudoboiling point, $\unicode[STIX]{x0394}h_{pl}$ (pseudoliquid), is approximately equal to the enthalpy difference from the pseudoboiling point to the warm wall, $\unicode[STIX]{x0394}h_{pg}$ (pseudogas). In other words, we assume: $\unicode[STIX]{x0394}h_{pl}\approx \unicode[STIX]{x0394}h_{pg}$ . Given the small temperature difference between both walls, the error of this approximation is between 3.5 % (for the $\unicode[STIX]{x0394}5~\text{K}$ case) and 4.6 % (for the $\unicode[STIX]{x0394}20~\text{K}$ case) based on the NIST data. Using this simplifying assumption, we can relate the total enthalpy of the pseudogas and pseudoliquid domain of the flow by the following relationship:

(3.1) $$\begin{eqnarray}\displaystyle {\mathcal{V}}_{pg}\int _{T_{pb}-\unicode[STIX]{x0394}T/2}^{T_{pb}}\unicode[STIX]{x1D70C}c_{p}\,\text{d}T={\mathcal{V}}_{pl}\int _{T_{pb}}^{T_{pb}+\unicode[STIX]{x0394}T/2}\unicode[STIX]{x1D70C}c_{p}\,\text{d}T, & & \displaystyle\end{eqnarray}$$

where the terms ${\mathcal{V}}_{pg}$ and ${\mathcal{V}}_{pl}$ respectively denote the volume of the pseudogas and pseudoliquid phases. By computing this relation based on the tabulated properties, we can estimate the relative volume of pseudogas and pseudoliquid in the simulation. It is clear that under perfect gas conditions and/or modest density variation (e.g. away from the critical point), the volume of fluid below and above the average temperature would be equal. In fact, our reference simulation at $p_{b}=2p_{cr}$ , which shows a very modest variation in thermophysical properties compared to the other runs at $p_{b}=1.1p_{cr}$ , contains about 55 % volume of pseudoliquid based on the above relation. This compares favourably to the simulation results which show the pseudotransition point at the middle of the domain (figure 4 b). When considering the cases at $p_{b}\approx 1.1p_{cr}$ , we find the pseudoliquid volume takes up 37.0 %, 69.6 % and 76.1 % for the $\unicode[STIX]{x0394}T=5$ , 10 and 20 K cases, respectively, based on this simple model. This trend corresponds favourably to the average pseudotransition height computed from the DNS which is 38.5 %, 77.5 % and 94.5 %. These results support the idea that high nonlinearity of the thermophysics near the critical point will dictate the location of the pseudoboiling point relative to the solid walls. As shown here, this is in fact a first-order real fluid effect.

In all cases, the transition from a seemingly fully thermally mixed region in channel core (i.e. $\overline{T}(y)$ is relatively uniform and close to the pseudoboiling value) is more defined than in the reference simulation. Such steep mean flow gradients near the walls sustain significant density and enthalpy fluctuations, up to $\unicode[STIX]{x1D70C}_{rms,max}=44.1~\text{kg}~\text{m}^{-3}$ and $h_{rms,max}=8.9~\text{kJ}~\text{kg}^{-1}$ , respectively (as discussed later in figures 11 and 12) for the $\unicode[STIX]{x0394}T=20$  K case. The very high heat capacity of the fluid undergoing pseudophase transition, on the other hand, limits the temperature fluctuations to $T_{rms,max}<2~\text{K}$ .

Figure 5.

Reynolds-averaged streamwise velocity component (a) and its wall-normal gradient (b) for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ) and reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (●) (see § 2.3). Average location of pseudotransition for $\unicode[STIX]{x0394}T=5~\text{K}$ (○), 10 K (▵) and 20 K (▫).

Figure 6.

Wall-normal gradient of Reynolds-averaged mean density (a) and temperature (b) for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ) and reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (●) (see § 2.3). Average location of pseudotransition for $\unicode[STIX]{x0394}T=5~\text{K}$ (○), 10 K (▵) and 20 K (▫).

The mean turbulent streamwise velocity profile (figure 5 a) becomes more asymmetric (with a slight acceleration of the pseudogaseous layer) with increasing $\unicode[STIX]{x0394}T$ , with an upwards shift in the maximum velocity location, $y/h=0.06$ for $\unicode[STIX]{x0394}T=5$  K, 0.11 for $\unicode[STIX]{x0394}T=10$  K and 0.17 for $\unicode[STIX]{x0394}T=20$  K (see inset in figure 5 b), following the same trend of the pseudotransition location, $y_{pb}$ . As a result, a larger velocity gradient magnitude is found near the top wall (the magnitude ratio of top-to-bottom velocity gradient is 1.24 for $\unicode[STIX]{x0394}T=5$  K, 1.32 for $\unicode[STIX]{x0394}T=10$  K, 1.44 for $\unicode[STIX]{x0394}T=20$  K in figure 5 b). In figure 6, while top-down asymmetries in the temperature gradient are confined to the sublayer regions, the mean density gradient profile is more visibly affected by the location of pseudotransition. A logarithmic increment of the centreplane of the temperature gradient is observed as $\unicode[STIX]{x0394}T$ is also increased logarithmically (i.e. $\text{d}(\unicode[STIX]{x0394}T)/\unicode[STIX]{x0394}T=\text{const.}$ ), suggesting a linear relation between the overall top-to-bottom equilibrium heat flux and $\unicode[STIX]{x0394}T$ . The latter is a surprising result given that the level of hydrodynamic and thermodynamic nonlinearity of the problem. These results also suggest that transcritical heat flux rates are amenable to straightforward dimensionless scaling in similar canonical set-ups. While the velocity gradient increase (decrease) in the pseudogaseous (pseudoliquid) region as $\unicode[STIX]{x0394}T$ is increased is not as significant as the corresponding variations in density and temperature gradients, the real fluid effects are very apparent when attempting to scale the mean velocity profiles with commonly used scaling laws.

Table 4.

Semi-local scaling factors where $\overline{u}_{\unicode[STIX]{x1D70F}}^{\ast }(y)=\sqrt{\overline{\unicode[STIX]{x1D70F}}_{w}/\overline{\unicode[STIX]{x1D70C}}(y)}$ .

For all $\unicode[STIX]{x0394}T$ values, the mean streamwise velocity profiles are scaled following the recently proposed approach by Trettel & Larsson (Reference Trettel and Larsson2016) which accounts for the wall heat transfer effects, the van Driest transformation (van Driest Reference van Driest1951) and the semi-local scaling (Huang et al. Reference Huang, Coleman and Bradshaw1995) (figure 7). The expressions of the three transformations considered are reported here for convenience and completeness.

Figure 7.

Mean streamwise velocity versus wall-normal coordinate in wall units scaled based on the conventional van Driest transformation plotted against wall-normal distance in classic wall units (a) and semi-locally scaled (Huang et al. Reference Huang, Coleman and Bradshaw1995) (b), transformed following Trettel & Larsson (Reference Trettel and Larsson2016) (c) for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——, thickened), 10 K (- - -) and 20 K ( $\cdots \cdots$ ); reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (circles); bottom wall (blue, ●) and top wall (red, ○). Profiles of the law of the wall ( $\overline{u}^{+}=y^{+}$ for the viscous sublayer; $\overline{u}^{+}=1/\unicode[STIX]{x1D705}\ln y^{+}+C$ , where $\unicode[STIX]{x1D705}=0.41$ and $C=5.2$ for the log-law region) are shown with a thin solid black line for reference.

The van Driest transformation (van Driest Reference van Driest1951) is given by

(3.2) $$\begin{eqnarray}\displaystyle \overline{u}_{VD}^{+}=\int _{0}^{\overline{u}^{+}}\left(\frac{\overline{\unicode[STIX]{x1D70C}}(y)}{\overline{\unicode[STIX]{x1D70C}}_{w}}\right)^{1/2}\,\text{d}\overline{u}^{+}, & & \displaystyle\end{eqnarray}$$

where $\overline{u}^{+}=\overline{u}(y)/\overline{u}_{\unicode[STIX]{x1D70F}}$ and the conventional set of scaling parameters reads

(3.3) $$\begin{eqnarray}\displaystyle y^{+}=\frac{y}{\unicode[STIX]{x1D6FF}_{v}}=\frac{y}{\overline{\unicode[STIX]{x1D707}}_{w}/(\overline{\unicode[STIX]{x1D70C}}_{w}\overline{u}_{\unicode[STIX]{x1D70F}})},\quad \overline{u}_{\unicode[STIX]{x1D70F}}=\sqrt{\overline{\unicode[STIX]{x1D70F}}_{w}/\overline{\unicode[STIX]{x1D70C}}_{w}}, & & \displaystyle\end{eqnarray}$$

whereas, for the semi-local scaling (Huang et al. Reference Huang, Coleman and Bradshaw1995), it reads

(3.4) $$\begin{eqnarray}\displaystyle y^{\ast }=\frac{y}{\unicode[STIX]{x1D6FF}_{v}^{\ast }}=\frac{y}{\overline{\unicode[STIX]{x1D707}}(y)/(\overline{\unicode[STIX]{x1D70C}}(y)\overline{u}_{\unicode[STIX]{x1D70F}}^{\ast }(y))},\quad \overline{u}_{\unicode[STIX]{x1D70F}}^{\ast }(y)=\sqrt{\overline{\unicode[STIX]{x1D70F}}_{w}/\overline{\unicode[STIX]{x1D70C}}(y)}, & & \displaystyle\end{eqnarray}$$

and its factors for various quantities are shown in table 4.

Figure 8.

Semi-local friction Reynolds number at the bottom (a) and top (b) wall for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ).

Finally, the transformation by Trettel & Larsson (Reference Trettel and Larsson2016) reads

(3.5) $$\begin{eqnarray}\displaystyle \overline{u}_{TL}^{+}=\int _{0}^{\overline{u}^{+}}\left(\frac{\overline{\unicode[STIX]{x1D70C}}(y)}{\overline{\unicode[STIX]{x1D70C}}_{w}}\right)^{1/2}\left[1+\frac{1}{2}\frac{1}{\overline{\unicode[STIX]{x1D70C}}(y)}\frac{\text{d}\overline{\unicode[STIX]{x1D70C}}(y)}{\text{d}y}y-\frac{1}{\overline{\unicode[STIX]{x1D707}}(y)}\frac{\text{d}\overline{\unicode[STIX]{x1D707}}(y)}{\text{d}y}y\right]\text{d}\overline{u}^{+}. & & \displaystyle\end{eqnarray}$$

While the reference results at $p_{b}=2p_{cr}$ collapse profiles from both walls in the log-law region with the Trettel & Larsson (Reference Trettel and Larsson2016) transformation, the widest spread is observed with the semi-local scaling. The transformed top- and bottom-wall streamwise velocity profiles at transcritical temperature conditions result in higher intercepts than the classic incompressible log-law as the pressure increases. In the recent publication by Ma et al. (Reference Ma, Yang and Ihme2018) replicating our same set-up, very large values of the transformed velocity $\overline{u}_{TL}^{+}$ in the log region were obtained by erroneously adopting a semi-locally scaled velocity differential in (3.5), suggesting inadequacy of this transformation for this flow. The present results show, instead, an acceptable collapse using the Trettel & Larsson (Reference Trettel and Larsson2016) transformation.

Effects of varying $\unicode[STIX]{x0394}T$ are visible (hence not collapsed perfectly) in all adopted transformations. Increasing $\unicode[STIX]{x0394}T$ results in an enhancement of real fluid effects (at the present conditions) yielding variations of the state of turbulence in the wall-normal direction, analysed below via extraction of the semi-local friction Reynolds number and density fluctuation intensity profiles.

Figure 8 shows the semi-local friction Reynolds number

(3.6) $$\begin{eqnarray}\displaystyle \overline{Re}_{\unicode[STIX]{x1D70F}}^{\ast }=\overline{Re}_{\unicode[STIX]{x1D70F}}\sqrt{\overline{\unicode[STIX]{x1D70C}}(y)/\overline{\unicode[STIX]{x1D70C}}_{w}}/(\overline{\unicode[STIX]{x1D707}}(y)/\overline{\unicode[STIX]{x1D707}}_{w}), & & \displaystyle\end{eqnarray}$$

where $\overline{Re}_{\unicode[STIX]{x1D70F}}=\overline{\unicode[STIX]{x1D70C}}_{w}\overline{u}_{\unicode[STIX]{x1D70F}}h/\overline{\unicode[STIX]{x1D707}}_{w}$ and the semi-locally scaled wall-normal coordinate, $y^{\ast }$ , is here intended as the relative distance from either the top or the bottom wall. The values of $\overline{Re}_{\unicode[STIX]{x1D70F}}^{\ast }$ in the bottom-wall viscous sublayer are lower than those near the top wall; the opposite occurs in the respective log-law regions. However, values of $\overline{Re}_{\unicode[STIX]{x1D70F},bot}^{\ast }$ in the log-law region are comparable across the different $\unicode[STIX]{x0394}T$ considered, while $\overline{Re}_{\unicode[STIX]{x1D70F},top}^{\ast }$ systematically decreases in the respective log-law region as $\unicode[STIX]{x0394}T$ increases (and as the pseudophase transitioning region of the flow approaches the top-wall buffer layer). The overall higher sensitivity of $\overline{Re}_{\unicode[STIX]{x1D70F}}^{\ast }$ to $\unicode[STIX]{x0394}T$ on the heated top wall is manifest in the van Driest transformed velocity (figure 7 a,b) showing the systematic increase in value with $\unicode[STIX]{x0394}T$ for the top wall more than the bottom wall. Also, as noted by Patel, Boersma & Pecnik (Reference Patel, Boersma and Pecnik2016), the van Driest transformed velocity profiles plotted as a function of $y^{\ast }$ in figure 7(b) for the bottom-wall region resemble the profiles of $\overline{Re}_{\unicode[STIX]{x1D70F},bot}^{\ast }$ and the inverse of $\overline{Re}_{\unicode[STIX]{x1D70F},top}^{\ast }$ (vice versa for the top wall). Changes in mean velocity profiles, as well as in the turbulence, are related to wall-normal variations of $\overline{Re}_{\unicode[STIX]{x1D70F}}^{\ast }$ . On the other hand, the Trettel & Larsson transformed profiles collapse the data across the various $\unicode[STIX]{x0394}T$ cases considered, in equal manners on both walls, despite the real fluid effects being more pronounced at the top wall (as also discussed later and illustrated in figures 12 and 13).

Although the Trettel & Larsson transform results in a collapse between the top and bottom walls for a given simulation (or given $\unicode[STIX]{x0394}T$ ), we do note a non-negligible spread among the simulations. More specifically, the same slope in the log law is obtained, but the intercept value increases with the wall-to-wall temperature difference. Such spread in the intercept values (compared to a universal profile) can arise either due to the error in the derivation of the stress balance condition or due the log-law assumptions; both are investigated here.

In regard to the closure of the stress balance, all the above wall scaling models rely on Morkovin’s hypothesis, which assumes that the compressibility effects on turbulence are only related to the mean density variations, not the density fluctuations. Further, all the scalings neglect the correlation between the viscosity and fluctuating velocity derivative as a component of the total shear. These hypotheses are evoked to simplify the near-wall stress balance, yielding

(3.7) $$\begin{eqnarray}\displaystyle \overline{\unicode[STIX]{x1D707}\frac{\unicode[STIX]{x2202}u}{\unicode[STIX]{x2202}y}}-\overline{\unicode[STIX]{x1D70C}u^{\prime }v^{\prime }}\approx \overline{\unicode[STIX]{x1D707}}\frac{\unicode[STIX]{x2202}\overline{u}}{\unicode[STIX]{x2202}y}-\overline{\unicode[STIX]{x1D70C}}\overline{u^{\prime }v^{\prime }}+\underbrace{\overline{\unicode[STIX]{x1D707}^{\prime }\frac{\unicode[STIX]{x2202}u^{\prime }}{\unicode[STIX]{x2202}y}}-\overline{\unicode[STIX]{x1D70C}^{\prime }u^{\prime }v^{\prime }}}_{}=\unicode[STIX]{x1D70F}_{w}. & & \displaystyle\end{eqnarray}$$

For ideal-gas flows in the subsonic and low-supersonic regime, the terms highlighted by the underbrace are negligible. In the transcritical regime, where the fluctuations of viscosity and density are significant, the error due to neglecting $\overline{\unicode[STIX]{x1D707}^{\prime }(\unicode[STIX]{x2202}u^{\prime }/\unicode[STIX]{x2202}y)}$ and $\overline{\unicode[STIX]{x1D70C}^{\prime }u^{\prime }v^{\prime }}$ may be significant. In the simulations presented herein, they are, respectively, up to 2 % and 5 % of the total shear (for $\unicode[STIX]{x0394}T=20$  K). This makes Morkovin’s hypothesis and the assumption of negligible viscosity fluctuation correlations questionable for the low-speed transcritical flows analysed herein.

We now move on to considering the assumptions underlying the log law. Its intercept value corresponds to the integration constant in the scaled velocity profile (Marusic et al. Reference Marusic, Monty, Hultmark and Smits2013), and is known to vary for wall roughness (Bradshaw Reference Bradshaw1994), heated (Lee et al. Reference Lee, Jung, Sung and Zaki2013) or superhydrophobic (Min & Kim Reference Min and Kim2004) boundary-layer flows. Therefore, the intercept value is highly dependent on the near-wall conditions. Using dimensional scaling arguments, Bradshaw (Reference Bradshaw1994) argues that the velocity derivative in the log layer is proportional to a velocity scale (friction velocity) and inversely to the wall distance:

(3.8) $$\begin{eqnarray}\displaystyle \frac{\text{d}u}{\text{d}y}=\frac{1}{\unicode[STIX]{x1D705}}\frac{1}{y}u_{\unicode[STIX]{x1D70F}}. & & \displaystyle\end{eqnarray}$$

In the compressible regime, a heuristic argument (Bradshaw Reference Bradshaw1994) is evoked in which the friction factor is replaced by a locally varying velocity scale (Trettel & Larsson Reference Trettel and Larsson2016):

(3.9) $$\begin{eqnarray}\displaystyle \frac{\text{d}u}{\text{d}y}=\frac{1}{\unicode[STIX]{x1D705}}\frac{1}{y}\left(\frac{\unicode[STIX]{x1D70F}_{w}}{\overline{\unicode[STIX]{x1D70C}}}\right)^{1/2}. & & \displaystyle\end{eqnarray}$$

The above equation, along with the stress-balance condition, is central to the Trettel & Larsson transform. When assessing the validity of (3.9) in the transcritical flows analysed herein, we note a significant relative error in the log region when applying the above relation. The error is above 20 % in the log layer in the bottom portion of the channel for $\unicode[STIX]{x0394}T=20$  K; the error on the top half of the channel is even more significant. In other words, this observation informs us that the selection of the characteristic velocity scales within the van Driest and Trettel & Larsson derivation of the log-law conditions is not straightforwardly extendable to transcritical flows.

3.2 Turbulent fluctuation intensities

Figure 9.

Root-mean-square of streamwise (a), wall-normal (b) and spanwise (c) Favre fluctuating velocity component and Reynolds shear stress (d) for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ) and reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (●) (see § 2.3). Average location of pseudotransition for $\unicode[STIX]{x0394}T=5~\text{K}$ (○), 10 K (▵) and 20 K (▫).

Table 5.

Top-to-bottom difference in root-mean-square peak values of streamwise, wall-normal and spanwise velocity components, density, temperature and pressure as percentage of the bottom peak root-mean-square value.

Figure 10.

Semi-local scaled root-mean-square of streamwise, wall-normal and spanwise Favre fluctuating velocity component (a,c) and Reynolds shear stress (b,d) at the bottom (a,b) and top (c,d) wall for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ) and reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (●) (see § 2.3). Semi-local scaling factors are given in table 4.

Other real fluid effects associated with transcritical thermal conditions are observable in the variance of the hydrodynamic turbulent fluctuations as shown in figure 9. As $\unicode[STIX]{x0394}T$ increases, the asymmetries with respect to the channel centreplane grow and the peak fluctuation intensity values at the top wall (pseudogaseous region towards which the pseudotransition location migrates) are attenuated corresponding to the values in the pseudoliquid flow as shown in table 5. This suggests the occurrence of damping of hydrodynamic turbulence due to the proximity to the region of pseudophase change. Such attenuation is noted in all Reynolds stress terms, but is strongest in the wall-normal velocity component, directly involved in the turbulent heat and mass transport working against the steep mean temperature and density gradient. Following Morinishi, Tamano & Nakabayashi (Reference Morinishi, Tamano and Nakabayashi2004), in figure 10 we show the semi-local scaled profiles of turbulent velocity fluctuations (see table 4 for the semi-local scaling factors). The semi-local scaling collapses well the high-pressure data $2p_{cr}$ onto the near-critical $1.1p_{cr}$ cases, with the exception of the streamwise velocity fluctuations; a much better collapse is observed across the various $\unicode[STIX]{x0394}T$ conditions for the same near-critical case in terms of the peak value level and its location in the semi-local wall unit, $y^{\ast }$ (see table 6).

Table 6.

Peak ranges of the root-mean-square of streamwise, wall-normal and spanwise Favre fluctuating velocity component and Reynolds shear stress at the bottom and top wall and their wall-normal location using the semi-local scaling.

Figure 11.

Root-mean-square of Reynolds fluctuations (a–c) and their normalized quantities (d–f) with respect to the local mean values for density (a,d), temperature (b,e) and pressure (c,f) for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ) and reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (●) (see § 2.3). Average location of pseudotransition for $\unicode[STIX]{x0394}T=5~\text{K}$ (○), 10 K (▵) and 20 K (▫).

Figure 12.

Root-mean-square of Favre fluctuations for enthalpy (a) and wall-normal turbulent enthalpy flux (b) for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ). Average location of pseudotransition for $\unicode[STIX]{x0394}T=5~\text{K}$ (○), 10 K (▵) and 20 K (▫).

Figure 13.

Semi-local scaled root-mean-square of Reynolds fluctuating density at the bottom (a) and top (b) wall for $p_{b}=1.1p_{cr}$ and $\unicode[STIX]{x0394}T=5$  K (——), 10 K (- - -) and 20 K ( $\cdots \cdots$ ) and reference data for $p_{b}=2p_{cr}$ and $\unicode[STIX]{x0394}T=20$  K (●) (see § 2.3). Semi-local scaling factors are shown in table 4.

On the contrary, the proximity to the pseudophase change location $y_{pb}$ (locally) and the increasing top-to-bottom temperature difference $\unicode[STIX]{x0394}T$ (globally) enhance the intensity of all thermodynamic fluctuations (figures 11 and 12 a). In spite of damping in the wall-normal velocity fluctuations, the wall-normal turbulent enthalpy flux is enhanced (figure 12 b) by increasing $\unicode[STIX]{x0394}T$ , as expected by the statistical steadiness of the flow, implying equilibrium conditions for the turbulent heat transfer.

For any given $\unicode[STIX]{x0394}T$ , the root-mean-square peak of density, temperature and enthalpy closer to the location of pseudophase transition, $y_{pb}$ , has a higher value than that farther away. This is not surprising as the thermophysical gradients are the largest at $y\sim y_{pb}$ (see figure 1), and therefore pseudoboiling conditions are expected to enhance the local thermodynamic fluctuations, concomitantly enhancing the asymmetry in the fluctuating profiles; in fact, based on the transport equation for density variance (Chassaing et al. Reference Chassaing, Antonia, Anselmet, Joly and Sarkar2013), the mean density gradient is responsible for production of density fluctuations, consistently with the data presented in this paper. Moreover, as $y_{pb}$ moves upwards for increasing $\unicode[STIX]{x0394}T$ , it approaches the peak of the shear Reynolds stress and enthalpy flux (figure 12 b), only significantly increasing the latter roughly proportionally to $\unicode[STIX]{x0394}T$ .

The density fluctuation intensity normalized by its mean local value shown in figure 11 is significantly higher than the temperature fluctuation (for $\unicode[STIX]{x0394}T=20$  K, $(\unicode[STIX]{x1D70C}_{rms}^{\prime }/\overline{\unicode[STIX]{x1D70C}})_{max}\simeq 26\times (T_{rms}^{\prime }/\overline{T})_{max}$ ). These real fluid effects have a direct impact on the structure of near-wall turbulence and are attenuated with increasing pressure.

Previous calculations involving dense gases in a supersonic channel flow by Sciacovelli et al. (Reference Sciacovelli, Cinnella and Gloerfelt2017) observed the maximum density fluctuation intensity located in the viscous sublayer and, as such, it was argued that it did not alter the turbulence structure significantly while satisfying Morkovin’s hypothesis. In the present results, the $\unicode[STIX]{x1D70C}_{rms}^{\prime }$ peak is located in the buffer layer, and much more pronounced at the top wall (where the real fluid effects become more concentrated as $\unicode[STIX]{x0394}T$ increases) therein contributing to the inhomogeneity of the state of turbulence, as discussed earlier in the analysis of the semi-local friction Reynolds number (figure 8). The bottom and top near-wall peaks of $\unicode[STIX]{x1D70C}_{rms}^{\prime }$ based on the semi-local scaling for the $\unicode[STIX]{x0394}T$ conditions in this study are located in $y^{\ast }=7$ –13 approximately corresponding to the buffer layer and the location moves slightly towards the channel centreplane with increasing $\unicode[STIX]{x0394}T$ (see figure 13).