5.1 Computation of PTHF

First we describe the computation of the PTHF from our PR-DNS set-up, and then the computation of the corresponding transport term. In our thermal fully developed gas–solid flow, the fluid velocity is statistically homogeneous whereas the fluid temperature is inhomogeneous along the axial location $x_{\Vert }$ . Therefore, any average involving fluid temperature has to be computed over a cross-sectional plane at a given axial location. Here the PTHF is computed from PR-DNS data using cross-sectional averages over $M$ realizations as

(5.2) $$\begin{eqnarray}\langle I_{f}u_{i}^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle (x_{\Vert })\approx \frac{1}{M}\mathop{\sum }_{{\it\omega}=1}^{M}\left\{\frac{1}{A}\int _{A}I_{f}u_{i}^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}(\boldsymbol{x};{\it\omega})\,\text{d}A\right\},\end{eqnarray}$$

where the fluid velocity fluctuation for each realization is defined as

(5.3) $$\begin{eqnarray}u_{i}^{\prime \prime (f)}(\boldsymbol{x};{\it\omega})=u_{i}(\boldsymbol{x};{\it\omega})-\{u_{i}^{(f)}\}_{V}({\it\omega}),\end{eqnarray}$$

and where $\{u_{i}^{(f)}\}_{V}$ is the volumetric mean fluid velocity that is computed as

(5.4) $$\begin{eqnarray}\{u_{i}^{(f)}\}_{V}({\it\omega})=\frac{\displaystyle \frac{1}{V}\int _{V}I_{f}(\boldsymbol{x};{\it\omega})u_{i}(\boldsymbol{x};{\it\omega})\,\text{d}V}{\displaystyle \frac{1}{V}\int _{V}I_{f}(\boldsymbol{x};{\it\omega})\,\text{d}V}=\frac{1}{V_{f}}\int _{V}I_{f}u_{i}\,\text{d}V.\end{eqnarray}$$

The non-dimensional temperature fluctuation ${\it\phi}^{\prime \prime }(\boldsymbol{x};{\it\omega})$ for each realization is defined as

(5.5) $$\begin{eqnarray}{\it\phi}^{\prime \prime }(\boldsymbol{x};{\it\omega})={\it\phi}(\boldsymbol{x};{\it\omega})-\{{\it\phi}^{(f)}\}_{cs}(x_{\Vert };{\it\omega}),\end{eqnarray}$$

where $\{{\it\phi}^{(f)}\}_{cs}$ is the cross-sectional average of the non-dimensional temperature along the axial location (see (4.13)). Note that due to periodicity in the $y$ and $z$ directions only the PTHF along the axial coordinate $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle$ is non-zero. Therefore, we only discuss the axial component of the PTHF in the following.

Since the PTHF is computed using the cross-sectional average in (5.2) (unlike the Nusselt number, which is computed using a volume average since it is statistically homogeneous), it is susceptible to higher statistical variability than the Nusselt number. Figure 4 shows the axial variation of the ensemble-averaged PTHF $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle$ for $Re_{m}=100$ and ${\it\varepsilon}_{s}=0.4$ . The square symbols are the ensemble average from 5 multiple independent simulations (MIS) while the downward triangles are the ensemble average from 50 MIS. Both averages are very close, indicating convergence. However, as expected, the one-sided error bars from 5 MIS (denoted below the square symbols) are larger than the error bars from 50 MIS (above the downward triangles). Since the ensemble-averaged PTHF obtained from 50 MIS is within the range of the error bars of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle$ obtained from 5 MIS, we use 5 MIS to quantify the PTHF over a range of mean slip Reynolds numbers and solid volume fractions.

Figure 4. Variation of the ensemble-averaged PTHF normalized by the magnitude of mean slip velocity $|\langle \boldsymbol{W}\rangle |$ along axial location $x_{\Vert }$ over 5 and 50 MIS at $Re_{m}=100$ and ${\it\varepsilon}_{s}=0.4$ . The square and triangle symbols represent the PTHF obtained using 5 and 50 MIS, respectively. One-sided error bars indicate 95 % confidence intervals: the error bars below the squares correspond to 5 MIS while the error bars above the triangles correspond to 50 MIS.

Figure 4 also shows that the PTHF decays along the axial coordinate. In the hydrodynamic solution of this gas–solid flow, the PTKE $\langle I_{f}u_{i}^{\prime \prime (f)}u_{i}^{\prime \prime (f)}\rangle /2$ does not change with axial location since the velocity field is statistically homogeneous. However, in the heat transfer problem, the fluid temperature variance $\langle I_{f}\boldsymbol{{\it\phi}}^{\prime \prime (f)}\boldsymbol{{\it\phi}}^{\prime \prime (f)}\rangle$ decays (result not shown here) because the mean temperature gradient decays along the axial coordinate. Correspondingly, the covariance of temperature and velocity also decays. This corresponds to the decay of $\langle {\it\phi}_{m}\rangle$ (see figure 2 a) that is shown to be the sole source of inhomogeneity in the PTHF (see (5.1)).

As shown in (5.1), another way to compute the PTHF is to obtain $\langle I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\rangle$ , and then multiply it by the average non-dimensional bulk fluid temperature $\langle {\it\phi}_{m}\rangle$ . Since the scaled temperature field ${\it\theta}$ is homogeneous, we also expect $\langle I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\rangle$ to be statistically homogeneous along the axial coordinate. Figure 5 shows the ensemble-averaged axial value of $\langle I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\rangle$ (only $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ is non-zero due to periodicity in the $y$ and $z$ directions), which is computed by

(5.6) $$\begin{eqnarray}\langle I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\rangle \approx \frac{1}{M}\mathop{\sum }_{{\it\omega}=1}^{M}\{I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\}(x_{\Vert };{\it\omega}),\end{eqnarray}$$

for 5 and 50 MIS. Again the two ensemble-averaged values are reasonably close to each other, with higher MIS yielding smaller error bars, as expected. The ensemble average using 50 MIS clearly shows that $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ is statistically homogeneous. Therefore, $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ can be computed using a volume average. For $Re_{m}=100$ and ${\it\varepsilon}_{s}=0.4$ , the volume-averaged value of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ (i.e. $\overline{\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle }=(1/L)\int _{0}^{L}\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle (x_{\Vert })\,\text{d}x_{\Vert }$ , for convenience we drop the overbar later) is approximately 0.22 from five MIS and 0.20 from 50 MIS. This indicates that the volume-averaged value with fewer realizations is close to the one from 50 MIS. Therefore, $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ can be calculated as a volume average from five realizations and only depends on mean slip Reynolds number and solid volume fraction.

Figure 5. Variation of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ normalized by mean slip velocity $|\langle \boldsymbol{W}\rangle |$ along axial location $x_{\Vert }$ at $Re_{m}=100$ and ${\it\varepsilon}_{s}=0.4$ . The red and blue symbols represent the PTHF obtained using 5 and 50 realizations, respectively. One-sided error bars indicate 95 % confidence intervals: the error bars above the squares correspond to 5 MIS and the error bars below the triangles correspond to 50 MIS.

In order to develop a model for the PTHF, we characterize the dependence of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ on mean slip Reynolds number and solid volume fraction as shown in figure 6. Since ${\it\theta}$ is non-dimensional, we expect $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ to scale with $|\langle \boldsymbol{W}\rangle |$ and therefore increase with mean slip Reynolds number ( $Re_{m}=(1-{\it\varepsilon}_{s})|\langle \boldsymbol{W}\rangle |D/{\it\nu}_{f}$ ) for a fixed solid volume fraction. For a fixed solid volume fraction (see figure 6 a), the volume average of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ normalized by the magnitude of the mean slip velocity $|\langle \boldsymbol{W}\rangle |$ is not constant but decreases slightly with increasing mean slip Reynolds number, indicating that the dependence of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ on $|\langle \boldsymbol{W}\rangle |$ is not exactly linear. The mean value of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ is not very sensitive to the mean slip Reynolds number but lies in the range 0.2–0.34 for $1\leqslant Re_{m}\leqslant 100$ . Figure 6(b) shows that for a fixed Reynolds number $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle /|\langle \boldsymbol{W}\rangle |$ first increases with increasing solid volume fraction up to ${\it\varepsilon}_{s}=0.2$ , and then decreases for ${\it\varepsilon}_{s}>0.2$ .

Figure 6. Dependence of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ on (a) mean slip Reynolds number at ${\it\varepsilon}_{s}=0.1{-}0.5$ and (b) solid volume fraction at $Re_{m}=1$ , $50$ and 100. The symbols are $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ from PR-DNS data and the lines are the correlation by fitting PR-DNS data in (5.7). Error bars indicate 95 % confidence intervals using 5 MIS.

In order to develop a PTHF model in § 5.3, a correlation for $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ is given by fitting PR-DNS data:

(5.7) $$\begin{eqnarray}\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle =(1-{\it\varepsilon}_{s})(0.2+1.2{\it\varepsilon}_{s}-1.24{\it\varepsilon}_{s}^{2})\exp (-0.002Re_{m})|\langle \boldsymbol{W}\rangle |.\end{eqnarray}$$

This correlation, shown by the lines in figure 6, fits the data with an average deviation of 8 %. It is valid in the range of $0.1\leqslant {\it\varepsilon}_{s}\leqslant 0.5$ and $1\leqslant Re_{m}\leqslant 100$ . Note that the homogeneity of $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ in (5.1) also implies that, in thermally fully developed homogeneous flow, the source of inhomogeneity in PTHF arises solely from the average bulk fluid temperature. More importantly, inhomogeneity of PTHF implies that the transport term corresponding to the PTHF in (1.1) is non-zero.

In the following, we first quantify the magnitude of PTHF relative to convective mean flux in § 5.2 and then propose a model for it in § 5.3. In § 5.4 the Péclet number dependence of the effective thermal diffusivity is explained using a wake scaling analysis. In § 5.5 we quantify the magnitude of the transport term involving the PTHF relative to the average gas–solid heat transfer term in (1.1).

5.2 Relative importance of the PTHF to convective mean flux

In order to verify the importance of PTHF, we compare the PTHF with the convective mean flux ${\it\varepsilon}_{f}\langle u_{i}^{(f)}\rangle \langle {\it\phi}^{(f)}\rangle$ that appears in (1.1). Based on the expression for the PTHF in (5.1) and the relation for $\langle {\it\phi}^{(f)}\rangle$ in (4.12), the ratio of the PTHF to the convective mean flux can be written as

(5.8) $$\begin{eqnarray}\frac{\langle I_{f}u_{i}^{\prime \prime (f)}{\it\phi}^{\prime \prime }\rangle (x_{\Vert })}{{\it\varepsilon}_{f}\langle u_{i}^{(f)}\rangle \langle {\it\phi}^{(f)}\rangle (x_{\Vert })}=\frac{\langle I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\rangle \langle {\it\phi}_{m}\rangle }{{\it\varepsilon}_{f}\langle u_{i}^{(f)}\rangle \langle {\it\phi}^{(f)}\rangle }=\frac{\langle I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\rangle }{{\it\varepsilon}_{f}\langle u_{i}^{(f)}\rangle \langle {\it\theta}^{(f)}\rangle }.\end{eqnarray}$$

Since all the terms in the last expression of the above equation are homogeneous, the ratio of the PTHF to the convective mean flux is independent of axial location. Since $\langle I_{f}u_{i}^{\prime \prime (f)}{\it\theta}\rangle$ , $\langle u_{i}^{(f)}\rangle$ and $\langle {\it\theta}^{(f)}\rangle$ are functions of mean slip Reynolds number and solid volume fraction, this ratio also depends only on mean slip Reynolds number and solid volume fraction.

It can be shown that the ratio of PTHF to the convective mean flux is independent of axial location. The reason for this homogeneity lies in the important property of the PTHF that we deduced in the paragraph preceding equation (5.1). There we used the relation ${\it\phi}=\langle {\it\phi}_{m}\rangle {\it\theta}$ to show in (5.1) that the spatial inhomogeneity in the PTHF $\langle I_{f}\boldsymbol{u}^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle (x_{\Vert })$ arises solely from the inhomogeneity of the non-dimensional bulk temperature $\langle {\it\phi}_{m}\rangle$ . When the PTHF is divided by the convective mean flux, the only spatially inhomogeneous terms appear as a ratio $\langle {\it\phi}_{m}\rangle /\langle {\it\phi}^{(f)}\rangle$ . However, by virtue of (4.12) this ratio is simply $1/\langle {\it\theta}^{(f)}\rangle$ , which is statistically homogeneous because the locally scaled excess fluid temperature field ${\it\theta}$ is statistically homogeneous. In order to verify this independence of the ratio of the PTHF to the convective mean flux on axial location, figure 7(a) shows the axial variation of this ratio using 5 and 50 realizations. Although the ratio obtained from 5 realizations has more statistical variability than the one from 50 MIS, neither shows any systematic dependence on axial location. This finding confirms the statistical homogeneity of this ratio in (5.8).

Figure 7. (a) Axial variation of the ratio of the PTHF to the convective mean flux at $Re_{m}=100$ and ${\it\varepsilon}_{s}=0.4$ . The open circles and the squares represent the ratio obtained from 5 MIS and 50 MIS, respectively. Error bars indicate 95 % confidence intervals from 5 MIS (blue) and 50 MIS (red). (b) Comparison of the PTHF with convective mean flux ${\it\varepsilon}_{f}\langle u_{\Vert }^{(f)}\rangle \langle {\it\phi}^{(f)}\rangle$ in the range $1\leqslant Re_{m}\leqslant 100$ and $0.1\leqslant {\it\varepsilon}_{s}\leqslant 0.5$ . The open symbols represent the ratio of the PTHF and ${\it\varepsilon}_{f}\langle u_{\Vert }^{(f)}\rangle \langle {\it\phi}^{(f)}\rangle$ obtained from PR-DNS data. Error bars indicate 95 % confidence intervals from 5 MIS.

Figure 7(b) shows a comparison of the axial PTHF with the convective mean flux ${\it\varepsilon}_{f}\langle u_{\Vert }^{(f)}\rangle \langle {\it\phi}^{(f)}\rangle$ over a range of mean slip Reynolds numbers and solid volume fractions. Note that again due to the periodicity in the $y$ and $z$ directions, we only compare the axial value $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ with ${\it\varepsilon}_{f}\langle u_{\Vert }^{(f)}\rangle \langle {\it\phi}^{(f)}\rangle$ . The symbols denote the ratio of the PTHF to the convective mean flux. For a fixed Reynolds number, the ratio increases with increasing solid volume fraction in most cases. For $Re_{m}=100$ , the PTHF is approximately 40 % of the convective mean flux at ${\it\varepsilon}_{s}=0.1$ but approximately 70 % of the convective mean flux at ${\it\varepsilon}_{s}=0.5$ . The increase in the magnitude of the PTHF results from higher fluctuations at higher solid volume fraction. For a fixed solid volume fraction, the ratio of the PTHF to the convective mean flux tends to decrease slightly with increasing mean slip Reynolds number. Overall, the magnitude of PTHF is in the range 40–100 % of the convective mean flux. Therefore, the PTHF is a significant fraction of the total convective flux even at low solid volume fraction.

5.3 Model for pseudo-turbulent heat flux

Since the PTHF is found to be significant for gas–solid heat transfer, the transport term involving the PTHF needs to be modelled in CFD simulations based on the two-fluid model (see (1.1)). In order to develop a model for the PTHF, we introduce a gradient-diffusion model by analogy with turbulent scalar flux models in single-phase flow (Fox Reference Fox2003):

(5.9) $$\begin{eqnarray}R_{\boldsymbol{u}{\it\phi}}=\frac{\langle I_{f}u_{j}^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle (x_{\Vert })}{\langle I_{f}\rangle }=-{\it\alpha}_{jk,PT}\frac{\partial \langle {\it\phi}^{(f)}\rangle }{\partial x_{k}},\end{eqnarray}$$

where ${\it\alpha}_{jk,PT}$ is the pseudo-turbulent thermal diffusivity. Note that in general the thermal diffusivity is a tensor rather than a scalar. However, in our gas–solid heat transfer problem, the only non-zero component of the PTHF is the axial component which is aligned with the gradient of the mean fluid temperature. Therefore, we can only deduce one component ${\it\alpha}_{\Vert ,\Vert }$ of the pseudo-turbulent thermal diffusivity tensor from the PR-DNS data, as follows:

(5.10) $$\begin{eqnarray}R_{u_{\Vert }{\it\phi}}=\frac{\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle (x_{\Vert })}{\langle I_{f}\rangle }=-{\it\alpha}_{PT}\frac{\partial \langle {\it\phi}^{(f)}\rangle }{\partial x_{\Vert }},\end{eqnarray}$$

where ${\it\alpha}_{PT}={\it\alpha}_{\Vert ,\Vert }$ .

Once the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ is computed, the transport term involving the PTHF at a given axial location can be obtained in terms of the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ and average non-dimensional fluid temperature $\langle {\it\phi}^{(f)}\rangle$ as:

(5.11) $$\begin{eqnarray}\frac{\partial }{\partial x_{\Vert }}\{{\it\rho}_{f}c_{pf}\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle (x_{\Vert })\}=\frac{\partial }{\partial x_{\Vert }}\left(-{\it\varepsilon}_{f}{\it\rho}_{f}c_{pf}{\it\alpha}_{PT}\frac{\partial \langle {\it\phi}^{(f)}\rangle }{\partial x_{\Vert }}\right).\end{eqnarray}$$

A model for the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ can be derived by substituting the relation $\langle {\it\phi}^{(f)}\rangle =\langle {\it\theta}^{(f)}\rangle \langle {\it\phi}_{m}\rangle$ in (4.12) into (5.10), then substituting (5.1) and replacing $\langle {\it\phi}_{m}\rangle$ with the exponential decay model $\langle {\it\phi}_{m}\rangle =\exp \left(-{\it\lambda}x_{\Vert }/D\right)$ (see (4.8)) to obtain:

(5.12) $$\begin{eqnarray}{\it\alpha}_{PT}=-\frac{\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle \exp (-{\it\lambda}x_{\Vert }/D)}{-{\it\varepsilon}_{f}\langle {\it\theta}^{(f)}\rangle {\displaystyle \frac{{\it\lambda}}{D}}\exp (-{\it\lambda}x_{\Vert }/D)}=\frac{D}{{\it\lambda}}\frac{\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle }{(1-{\it\varepsilon}_{s})\langle {\it\theta}^{(f)}\rangle }.\end{eqnarray}$$

Using the correlation for $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ given in (5.7), the expressions for ${\it\lambda}$ given in (4.9) and the average scaled fluid temperature $\langle {\it\theta}^{(f)}\rangle$ from Sun et al. (Reference Sun, Tenneti and Subramaniam2015), the final expression for ${\it\alpha}_{PT}$ is

(5.13) $$\begin{eqnarray}{\it\alpha}_{PT}=\frac{4D(Re_{m}+1.4)Pr}{6{\rm\pi}{\it\varepsilon}_{s}\langle Nu\rangle }\frac{(0.2+1.2{\it\varepsilon}_{s}-1.24{\it\varepsilon}_{s}^{2})\exp (-0.002Re_{m})|\langle \boldsymbol{W}\rangle |}{[1-1.6{\it\varepsilon}_{s}(1-{\it\varepsilon}_{s})-3{\it\varepsilon}_{s}(1-{\it\varepsilon}_{s})^{4}\exp (-Re_{m}^{0.4}{\it\varepsilon}_{s})]}.\end{eqnarray}$$

It is noteworthy that in this model ${\it\alpha}_{PT}$ scales as $D|\langle \boldsymbol{W}\rangle |$ . Also as expected ${\it\alpha}_{PT}$ is independent of axial location and depends only on $Re_{m}$ , ${\it\varepsilon}_{s}$ , and $Pr$ (note that $\langle Nu\rangle$ is also a function of $Re_{m}$ , ${\it\varepsilon}_{s}$ , and $Pr$ (cf. equation (27) in Sun et al. (Reference Sun, Tenneti and Subramaniam2015)).

In order to evaluate the performance of this model for the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ , we compare ${\it\alpha}_{PT}$ obtained from direct quantification of the PTHF $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle$ using PR-DNS data in (5.10) with the model expression given by (5.13). Figure 8 shows the dependence of the ratio ${\it\alpha}_{PT}/{\it\alpha}_{f}$ ( ${\it\alpha}_{f}$ is the constant molecular thermal diffusivity that is equal to ${\it\nu}_{f}/Pr$ ) on mean slip Reynolds number and solid volume fraction. The symbols denote the values of ${\it\alpha}_{PT}/{\it\alpha}_{f}$ extracted from PR-DNS data which show that the pseudo-turbulent thermal diffusivity is two orders of magnitude larger than its molecular counterpart. In figure 8 the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ increases with increasing mean slip Reynolds number for a fixed solid volume fraction. This increase is due to the increase in the magnitude of $\boldsymbol{u}^{\prime \prime (f)}$ with increasing Reynolds number.

Figure 8. Dependence of the pseudo-turbulent thermal diffusivity normalized by the molecular thermal diffusivity in the fluid phase ${\it\alpha}_{f}$ for gas–solid flow on mean slip Reynolds number and solid volume fraction. The symbols represent the average values from PR-DNS data using 5 MIS. The lines represent the model for the pseudo-turbulent thermal diffusivity for ${\it\varepsilon}_{s}=0.1$ , 0.3 and 0.5.

For a fixed Reynolds number, figure 8 shows that as the solid volume fraction increases the pseudo-turbulent thermal diffusivity decreases. Since the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ can be conceived as arising from the product of a velocity scale $\boldsymbol{u}^{\prime \prime (f)}$ (or equivalently the velocity scale $|\langle \boldsymbol{W}\rangle |$ , since the two are related by a correlation for $k_{f}$ given in Mehrabadi et al. (Reference Mehrabadi, Tenneti and Subramaniam2015)), and a length scale $\ell$ , the dependence for fixed Reynolds number must arise from a change in the length scale associated with ${\it\alpha}_{PT}$ . Looking at the expression for ${\it\alpha}_{PT}$ in (5.12), we can see that such a length scale dependence on solid volume fraction can arise from $D/{\it\lambda}(1-{\it\varepsilon}_{s})$ , $\langle {\it\theta}^{(f)}\rangle$ or $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle$ . Since the velocity field and the scaled temperature field ${\it\theta}$ are statistically homogeneous, the Eulerian two-point correlation corresponding to $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ is:

(5.14) $$\begin{eqnarray}{\it\rho}_{u_{\Vert }{\it\theta}}(\boldsymbol{r})=\frac{\langle I_{f}(\boldsymbol{x}){\it\theta}^{\prime \prime (f)}(\boldsymbol{x})\boldsymbol{\cdot }I_{f}(\boldsymbol{x}+\boldsymbol{r})u_{\Vert }^{\prime \prime (f)}(\boldsymbol{x}+\boldsymbol{r})\rangle }{\langle I_{f}(\boldsymbol{x}){\it\theta}^{\prime \prime (f)}(\boldsymbol{x})\boldsymbol{\cdot }I_{f}(\boldsymbol{x})u_{\Vert }^{\prime \prime (f)}(\boldsymbol{x})\rangle },\end{eqnarray}$$

and it can be computed to infer a length scale associated with $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle$ . Note that in the above equation $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}^{\prime \prime (f)}\rangle$ is equal to $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ because ${\it\theta}=\langle {\it\theta}^{(f)}\rangle +{\it\theta}^{\prime \prime (f)}$ and $\langle I_{f}u_{\Vert }^{\prime \prime (f)}\rangle =0$ .

Figure 9 shows the Eulerian two-point cross-correlation corresponding to scaled temperature–velocity for solid volume fractions of 0.1 and 0.4. The decay of the cross-correlation to zero within the computational domain establishes the adequacy of the domain size. In dispersion without interphase heat transfer, the temperature fluctuations are only driven by velocity fluctuations. Since the temperature–velocity correlation would be arising from the velocity–velocity correlation, the velocity–velocity correlation length (which is related to the Brinkman length for Stokes flow) is the important length scale for hydrodynamic dispersion (Koch & Brady Reference Koch and Brady1985). However, for the present study with interphase heat transfer, the temperature fluctuations arise from both the velocity fluctuation and the interphase heat transfer. The fluctuations from the interphase heat transfer may not scale with the velocity fluctuation. Therefore, this temperature–velocity cross-correlation is the appropriate correlation for our problem rather than the velocity–velocity correlation. The cross-correlation curves in figure 9 for the two volume fractions have comparable length scales. The length scale ( $L_{u_{\Vert }{\it\theta}}=\int _{0}^{\infty }{\it\rho}_{u_{\Vert }{\it\theta}}(\boldsymbol{r})\,\text{d}\boldsymbol{r}$ ) for the case with a solid volume fraction of 0.1 is 0.114, while for a solid volume fraction of 0.4 it is 0.078. While this is a 46 % increase, it alone cannot explain the 230 % increase in ${\it\alpha}_{PT}$ that is seen in figure 9. This implies that the length scale in the ${\it\theta}$ field is only weakly sensitive to solid volume fraction. Thus, $\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\theta}\rangle$ is not solely responsible for the change in length scale with solid volume fractions that is observed in ${\it\alpha}_{PT}$ . Also the scaled fluid temperature $\langle {\it\theta}^{(f)}\rangle$ varies only slightly with Reynolds number and solid volume fraction (see Sun et al. (Reference Sun, Tenneti and Subramaniam2015), Figure 8). Clearly, only the length scale $D/{\it\lambda}(1-{\it\varepsilon}_{s})$ is important in determining the magnitude of the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ . According to the expression for ${\it\lambda}$ (see (4.7)), with increasing solid volume fraction, the length scale $D/{\it\lambda}(1-{\it\varepsilon}_{s})$ decreases. This explains the decrease of the pseudo-turbulent thermal diffusivity ${\it\alpha}_{PT}$ with solid volume fraction in figure 8.

Figure 9. Decay of the scaled fluid temperature–velocity fluctuation cross-correlation functions with separation distance $\boldsymbol{r}$ obtained from steady flow past a random configuration of spheres at a solid volume fraction of 0.1 and 0.4, and mean slip Reynolds numbers of 100. The box length is $L=7.5D$ for solid volume fraction of 0.1 and $L=5D$ for for solid volume fraction of 0.4, respectively.

Figure 8 also compares the model for the pseudo-turbulent thermal diffusivity (see (5.13)) with the PR-DNS data. The lines represent the model for ${\it\alpha}_{PT}$ given by (5.13) at selected solid volume fractions. This figure shows that the model has a good agreement with PR-DNS data with an average difference of 18 %. Therefore, this model for the pseudo-turbulent thermal diffusivity can be used to compute the transport term involving the PTHF in the average fluid temperature equation (1.1) if we assume a simpler isotropic form of the pseudo-turbulent diffusivity tensor given in (5.9).

5.4 Scaling of pseudo-turbulent thermal diffusivity with Péclet number

Theoretical studies on hydrodynamic dispersion in fixed beds for Stokes/low Reynolds number flow predict the dependence of the effective diffusivity on the Péclet number (Brenner Reference Brenner1980; Carbonell & Whitaker Reference Carbonell and Whitaker1983; Eidsath et al. Reference Eidsath, Carbonell, Whitaker and Herrmann1983). Koch & Brady (Reference Koch and Brady1985) derived linear and $Pe\ln (Pe)$ dependencies in the effective diffusivity by solving the convection–diffusion equation for mass transfer with no source or sink of mass within the particles in Stokes flow using an asymptotic analysis that is valid in the dilute limit (low solid volume fraction). The analysis of Koch & Brady (Reference Koch and Brady1985) is for the case with no source or sink of mass within the particles, which is an important distinction from the present work. Here the Péclet number is defined as $Pe=Ua/D_{f}$ , where $U$ is the average velocity through the bed, $a$ is the particle radius and $D_{f}$ is the molecular diffusivity of the scalar. The linear dependence is attributed to the mechanical dispersion mechanism while the $Pe\ln (Pe)$ dependency is attributed to a non-mechanical dispersion mechanism that arises from the no-slip boundary condition (obtained from a boundary layer analysis).

While the Koch & Brady (Reference Koch and Brady1985) analysis is valid for Stokes flow at dilute solid volume fraction, the results of the present study span a range of Reynolds number from 1 to 100 and solid volume fraction values from 0.1 to 0.5. Furthermore, since in our heat transfer problem we impose the isothermal boundary condition on particle surfaces, we cannot compare directly with the results of Koch & Brady (Reference Koch and Brady1985). Acrivos et al. (Reference Acrivos, Hinch and Jeffrey1980) analysed Stokes flow past a fixed bed with heat transfer at low Péclet number when the local mean temperature profile is approximately linear, rather than exponential, corresponding to the decay length of the mean temperature being larger than the Brinkman screening length. However, the cases studied here correspond to $Pe\geqslant 1$ . Here we deduce the scaling of the effective thermal diffusivity with Péclet number in two ways. The first is based on correlations developed from the PR-DNS data for the average bulk fluid temperature and Nusselt number. We then present a scaling analysis similar to that described in Koch (Reference Koch1993) which is appropriate for $Re_{a}\gg 1$ , where $Re_{a}=U_{\Vert }a/{\it\nu}_{f}$ is the Reynolds number based on the radii of particle, $U_{\Vert }$ is the mean fluid velocity (which is the mean slip velocity for fixed particles considered in this study).

In appendix C we derive a model for the effective thermal diffusivity based on the correlations developed for the average bulk fluid temperature and Nusselt number. This reveals the scaling of the effective thermal diffusivity with Péclet number as:

(5.15) $$\begin{eqnarray}\frac{{\it\alpha}_{PT}+{\it\alpha}_{f}}{{\it\alpha}_{f}}=\frac{C_{1}C_{3}}{C_{2}(C_{4}+C_{5}Re_{m}^{0.7}Pr^{1/3})}Pe_{D}^{2}+1,\end{eqnarray}$$

where the coefficients $C_{1}$ through $C_{5}$ are functions of only the solid volume fraction. In figure 10, we compare this model evaluated at $Pr=0.7$ for a fixed solid volume fraction ( ${\it\varepsilon}_{s}=0.1$ ) with PR-DNS data. This derived model (represented by the red solid line) is very close to our PR-DNS data that are obtained for cases with different Prandtl number values. For these values the factor preceding the square of the Péclet number is approximately constant. This good agreement with the PR-DNS data shows that the effective thermal diffusivity has a $Pe_{D}^{2}$ scaling (see also the match with the blue dashed line representing $1+0.25Pe_{D}^{2}$ ).

Figure 10. Variation of $({\it\alpha}_{PT}+{\it\alpha}_{f})/{\it\alpha}_{f}$ with Péclet number $Pe_{D}=|\langle \boldsymbol{W}\rangle |D/{\it\alpha}_{f}$ at $Re_{m}=1$ (up-triangles), $Re_{m}=100$ (down-triangles) with $Pr=0.01$ , 0.1, 0.7 and 1, and $Re_{m}=10{-}50$ at $Pr=0.7$ (squares) for solid volume fraction of $0.1$ . The dashed line represents the $1+0.25Pe_{D}^{2}$ scaling and the solid line represents the model in (5.15). The dotted line represents the $0.065Pe_{D}^{2}[\ln (1/Pr)+1]+1$ scaling at $Re_{m}=100$ for different Prandtl numbers ( $0.01\leqslant Pr\leqslant 0.7$ ) in (D 12).

The dependence of effective thermal diffusivity on Péclet number can also be explained on the basis of scaling arguments in the hydrodynamic and thermal wakes behind a particle. The wake structure can be deduced from the conditionally averaged fluid velocity $\langle I_{f}U_{\Vert }\rangle _{c}(\boldsymbol{r}=\boldsymbol{X}-\boldsymbol{X}_{p}|\boldsymbol{X}_{p})$ and conditionally averaged scaled fluid temperature $\langle I_{f}(T-T_{s})/(\langle T_{m}\rangle -T_{s})\rangle _{c}(\boldsymbol{r}=\boldsymbol{X}-\boldsymbol{X}_{p}|\boldsymbol{X}_{p})$ , where $\boldsymbol{X}_{p}$ and $\boldsymbol{r}$ are the particle position and the relative separation between the particle and field point in the fluid, respectively. These conditional averages correspond to averaging the fluid velocity or temperature field over members of an ensemble where each particle’s centre has been translated to the origin.

Figure 11 shows that for dilute flow there exists a distinct hydrodynamic wake at $Re_{m}=100$ , and a distinct thermal wake behind the particle at both large and small Péclet number. It is observed in figure 11(b) and (c) that the thermal wake is longer and thinner at higher Péclet (or higher Prandtl number) than at lower Péclet number. The distance over which wake can diffuse due to the viscous diffusion in the cross-stream direction over the time it takes for the fluid to convect a distance $x_{\Vert }$ in the streamwise direction is $\sqrt{({\it\nu}_{f}x_{\Vert }/U_{\Vert })}$ . We can identify the width of the hydrodynamic wake $r_{WM}$ in the near-wake and far-wake regions as follows. For $x_{\Vert }<aRe_{a}$ , the diffusion of momentum in the near-wake region occurs over a smaller distance than the $O(a)$ size of the region disturbed by the particle. For $x_{\Vert }>aRe_{a}$ in the far-wake region the wake thickness is larger than the particle size leading to:

(5.16) $$\begin{eqnarray}\displaystyle r_{WM}=\left\{\begin{array}{@{}ll@{}}{\sim}O(a),\quad & x_{\Vert }<aRe_{a}\\ \displaystyle \left(\frac{{\it\nu}_{f}x_{\Vert }}{U_{\Vert }}\right)^{1/2}=\displaystyle a\left(\frac{x_{\Vert }}{aRe_{a}}\right)^{1/2},\quad & x_{\Vert }>aRe_{a}.\end{array}\right. & & \displaystyle\end{eqnarray}$$

Figure 11. Contour plot of (a) the conditionally averaged fluid velocity $\langle I_{f}U_{\Vert }\rangle _{c}/|\langle \boldsymbol{W}\rangle |$ , (b,c) the conditionally averaged scaled fluid temperature $\langle I_{f}(T-T_{s})/(\langle T_{m}\rangle -T_{s})\rangle _{c}$ based on $\langle T_{m}\rangle$ for solid volume fraction of 0.1 and mean slip Reynolds number of 100; (b) $Pr=0.01$ , (c) $Pr=1$ . The conditional average is obtained from 5 MIS.

The velocity fluctuation can be derived from the momentum balance equation ${\rm\pi}r_{WM}^{2}{\it\rho}_{f}U_{\Vert }u_{\Vert }^{\prime \prime }=F$ leading to

(5.17) $$\begin{eqnarray}\displaystyle \frac{u_{\Vert }^{\prime \prime }}{U_{\Vert }}=C_{D}\frac{a^{2}}{r_{WM}^{2}}=C_{D}\,\frac{aRe_{a}}{x_{\Vert }}=\left\{\begin{array}{@{}ll@{}}{\sim}O(C_{D}),\quad & x_{\Vert }<aRe_{a}\\ \displaystyle C_{D}\frac{aRe_{a}}{x_{\Vert }},\quad & x_{\Vert }>aRe_{a},\end{array}\right. & & \displaystyle\end{eqnarray}$$

where $F=C_{D}{\it\rho}_{f}U_{\Vert }^{2}{\rm\pi}a^{2}$ is the drag force and $C_{D}$ is the drag coefficient corresponding to a fixed particle bed. Essentially this says that in the near-wake region the fluid velocity does not vary much $(u_{\Vert }^{\prime \prime }=u-U_{\Vert }=O(U_{\Vert }))$ whereas in the far-wake region the fluctuation is far less than the mean velocity ( $u_{\Vert }^{\prime \prime }=u-U_{\Vert }\ll U_{\Vert }$ ).

Similarly, the width of the thermal wake $r_{WH}$ can be estimated on the basis of thermal diffusivity as

(5.18) $$\begin{eqnarray}\displaystyle r_{WH}=\left\{\begin{array}{@{}ll@{}}{\sim}O(a),\quad & x_{\Vert }<aPe_{a}\\ \displaystyle a\left(\frac{x_{\Vert }}{aPe_{a}}\right)^{1/2},\quad & x_{\Vert }>aPe_{a},\end{array}\right. & & \displaystyle\end{eqnarray}$$

and the temperature fluctuation can be derived from the energy balance equation ${\rm\pi}r_{WH}^{2}{\it\rho}_{f}c_{pf}U_{\Vert }T^{\prime \prime }=Q_{pf}=4{\rm\pi}a^{2}h(T_{s}-\langle T_{m}\rangle )$ as

(5.19) $$\begin{eqnarray}\displaystyle \frac{T^{\prime \prime }}{(T_{s}-\langle T_{m}\rangle )}=\frac{4h}{{\it\rho}_{f}c_{pf}U_{\Vert }}\frac{a^{2}}{r_{WH}^{2}}=\frac{4h}{{\it\rho}_{f}c_{pf}U_{\Vert }}\frac{aPe_{a}}{x_{\Vert }}=\left\{\begin{array}{@{}ll@{}}{\sim}O\displaystyle \left(\frac{4h}{{\it\rho}_{f}c_{pf}U_{\Vert }}\right),\quad & x_{\Vert }<aPe_{a}\\ \displaystyle \frac{4h}{{\it\rho}_{f}c_{pf}U_{\Vert }}\frac{aPe_{a}}{x_{\Vert }},\quad & x_{\Vert }>aPe_{a},\end{array}\right. & & \displaystyle\end{eqnarray}$$

where $Pe_{a}=U_{\Vert }a/{\it\alpha}_{f}$ is the Péclet number based on the radius of the particle, and $r_{WH}$ is the width of the thermal wake. The thermal wake for $x_{\Vert }<aPe_{a}$ depends on whether the temperature field is disturbed throughout on $O(a)$ region at the back of the particle or only in a thinner region where the thermal boundary layer separates from the particle. We assume that it is an $O(a)$ region and this assumption is verified by the thermal wakes in figure 11(b) and (c).

The unconditional ensemble-averaged PTHF $\langle I_{f}u_{\Vert }^{\prime \prime (f)}T^{\prime \prime (f)}\rangle$ is calculated from the wake scaling analysis as the particle number density $n_{p}$ times an integral over the probability density function (p.d.f.) of the conditionally averaged particle position $f$ :

(5.20) $$\begin{eqnarray}\langle I_{f}u_{\Vert }^{\prime \prime (f)}T^{\prime \prime (f)}\rangle =n_{p}\int _{0}^{L_{w}}\int \int f\,\text{d}x_{\Vert }\,\text{d}y\,\text{d}z,\end{eqnarray}$$

where $n_{p}$ is the particle number density defined as the ratio of the average number of particles to the volume of the domain, and $L_{w}$ is the length of the wake that represents the velocity contour surrounding the particle where the value of the conditionally ensemble-averaged velocity reaches $|\langle \boldsymbol{W}\rangle |$ (note that since the particles are stationary, the mean slip velocity is equal to the unconditionally averaged fluid velocity). Note that the full length of the far wake is not attained in the computational domain as shown in figure 11(a) due to hydrodynamic interactions with neighbour particles (note that the two-point velocity correlation has decayed to zero within the computational domain, indicating that the domain is large enough for this to not be an artefact of periodicity). In this study we have $Pr\leqslant 1$ , and the thermal far-wake and hydrodynamic near-wake region overlap in the interval $aPe_{a}<x_{\Vert }<aRe_{a}$ . By inserting (5.16)–(5.19), and integrating over the near-wake, intermediate and far-wake regions, the PTHF yields

(5.21) $$\begin{eqnarray}\langle I_{f}u_{\Vert }^{\prime \prime (f)}T^{\prime \prime (f)}\rangle =Pe_{a}\left[k_{2}\ln \left(\frac{1}{Pr}\right)+k_{3}\ln \left(\frac{L_{w}}{aRe_{a}}\right)+k_{1}\right],\end{eqnarray}$$

where $k_{1}$ , $k_{2}$ and $k_{3}$ are undetermined coefficients arising from the scaling estimates and uncertainty in the limits of the integral (see appendix D). In the above expression, the $\ln (1/Pr)$ term comes from the intermediate region and the constant term comes from the near wake. Note that since hydrodynamic interactions with neighbour particles cause the velocity to decay before achieving a far-wake behaviour, $\ln (L_{w}/aRe_{a})$ is not present in practice. The detailed derivation can be found in appendix D. Substituting this expression for the PTHF into the expression for the pseudo-turbulent thermal diffusivity:

(5.22) $$\begin{eqnarray}{\it\alpha}_{PT}=-\frac{\langle I_{f}u_{\Vert }^{\prime \prime (f)}T^{\prime \prime (f)}\rangle (x_{\Vert })}{{\displaystyle \frac{\partial \langle I_{f}T\rangle }{\partial x_{\Vert }}}}\sim \frac{\langle I_{f}u_{\Vert }^{\prime \prime (f)}T^{\prime \prime (f)}\rangle (x_{\Vert })}{{\displaystyle \frac{\partial (T_{s}-\langle T_{m}\rangle )}{\partial x_{\Vert }}}},\end{eqnarray}$$

and using the decay length scale of the bulk and mean fluid temperature $D/{\it\lambda}$ to write the gradient as $(T_{s}-\langle T_{m}\rangle )/(D/{\it\lambda})$ results in

(5.23) $$\begin{eqnarray}\frac{{\it\alpha}_{PT}+{\it\alpha}_{f}}{{\it\alpha}_{f}}=\frac{C_{D}Pe_{D}^{2}}{{\rm\pi}^{2}}\left[B_{2}\ln \left(\frac{1}{Pr}\right)+B_{3}\ln \left(\frac{L_{w}}{aRe_{a}}\right)+B_{1}\right]+1,\end{eqnarray}$$

where $B_{1}$ – $B_{3}$ are again undetermined coefficients. Note that using the correct length scale based on the mean temperature gradient is crucial to recovering the scaling observed in PR–DNS.

The wake analysis of the scaling of the effective thermal diffusivity with Péclet number is compared with PR-DNS data in figure 10. The results obtained from the wake scaling analysis (the dotted line) agree well with the PR-DNS data (the symbols) at $Re_{m}=100$ for $Pr<1$ . The $Pe_{D}^{2}$ scaling itself comes from there being a wake and from realizing that the decay length of the mean fluid temperature is the correct scaling to use for the mean temperature gradient. Therefore, this analysis of the hydrodynamic and thermal wakes behind the particle gives a physical explanation for the existence of a $Pe_{D}^{2}$ scaling in effective thermal diffusivity in the regime of high Reynolds number and low Prandtl number.

5.5 Relative importance of the PTHF in gas–solid heat transfer

We have found that the PTHF is significant when compared with the convective mean flux, especially for high solid volume fraction. In order to quantify the importance of the transport term involving the PTHF

(5.24) $$\begin{eqnarray}\langle T_{\boldsymbol{u}{\it\phi}}\rangle (x_{\Vert })\equiv \boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\{{\it\rho}_{f}c_{pf}\langle I_{f}\boldsymbol{u}^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle \},\end{eqnarray}$$

we need to compute the streamwise derivative of $\langle I_{f}\boldsymbol{u}^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle$ since the PTHF is only statistically inhomogeneous along the axial coordinate. However, as shown in figure 4, given that the ensemble-averaged statistical estimate in (5.2) has statistical variability, this can yield noisy results. In order to circumvent this difficulty, we integrate the transport term over the computational domain to express the mean value of the transport term involving the PTHF in the domain in terms of boundary values of the PTHF as follows:

(5.25) $$\begin{eqnarray}\displaystyle \overline{\langle T_{\boldsymbol{u}{\it\phi}}\rangle }=\overline{\langle T_{u_{\Vert }{\it\phi}}\rangle } & = & \displaystyle \frac{1}{L}\int _{L}\frac{\partial }{\partial x_{\Vert }}\{{\it\rho}_{f}c_{pf}\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle \}(x_{\Vert })\,\text{d}x_{\Vert }\nonumber\\ \displaystyle & = & \displaystyle \frac{1}{L}([{\it\rho}_{f}c_{pf}\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle ]_{out}-[{\it\rho}_{f}c_{pf}\langle I_{f}u_{\Vert }^{\prime \prime (f)}{\it\phi}^{\prime \prime (f)}\rangle ]_{in}),\end{eqnarray}$$

where $L$ is the length of the domain and $[\cdot ]_{in}$ and $[\cdot ]_{out}$ are obtained at the inlet and outlet of the computational domain, respectively. Note that due to periodic boundary conditions in the $y$ and $z$ directions the flux term in those directions is zero.

In the two-fluid equation (see (1.1)), since the gas–solid heat transfer is considered as a significant term, the importance of the transport term involving the PTHF can be quantified by comparing it with the average gas–solid heat transfer term (see (6.5)). Figure 12 shows a comparison of the transport term involving the PTHF scaled by the average gas–solid heat transfer over a range of mean slip Reynolds numbers and solid volume fractions. The colour symbols represent the ratio of the transport term involving the PTHF to the average gas–solid heat transfer. For a fixed Reynolds number beyond $Re_{m}=10$ , the transport term involving the PTHF is approximately 50 % of the average gas–solid heat transfer at ${\it\varepsilon}_{s}=0.5$ , and it drops to about 30 % of the average gas–solid heat transfer at ${\it\varepsilon}_{s}=0.1$ . This increase in the transport term involving the PTHF with increasing solid volume fraction is similar to the findings of Tenneti (Reference Tenneti2013) and Mehrabadi et al. (Reference Mehrabadi, Tenneti and Subramaniam2015) for PTKE in a fixed particle assembly. For a fixed solid volume fraction, the ratio of transport term involving the PTHF to the average gas–solid heat transfer does not vary significantly beyond $Re_{m}=10$ , but reduces to 15–20 % at $Re_{m}=1$ . Therefore, we conclude that the transport term involving the PTHF is important compared to average gas–solid heat transfer for high solid volume fraction, and it cannot be neglected in CFD simulations based on the two-fluid model.

Figure 12. Comparison of transport term involving the PTHF (see (5.25)) with the average gas–solid heat transfer (see (A 11)) in the range $1\leqslant Re_{m}\leqslant 100$ and $0.1\leqslant {\it\varepsilon}_{s}\leqslant 0.5$ . The symbols represent the transport term involving the PTHF obtained from PR-DNS data. Error bars indicate 95 % confidence intervals from 5 MIS. For clarity, only one-sided error bars are shown in this figure.

Note that for low Reynolds number ( $Re_{m}=1$ ) figure 12 shows that the ratio of PTHF to gas–solid heat transfer reduces. This results from the fact that the total convective term (mean convective heat flux and PTHF) needs to balance axial conduction in the fluid phase and average gas–solid heat transfer in the average temperature equation in the steady flow (see (1.1)). A budget analysis of these terms that is described later in this paper illustrates this point.

6.1 Verification of the fluid-phase axial conduction model

Since average axial conduction in the fluid phase is modelled in terms of the second derivative of $\langle {\it\phi}^{(f)}\rangle$ , a model for axial conduction in the fluid phase can be developed by using the expression for $\langle {\it\phi}^{(f)}\rangle$ ( $\langle {\it\phi}^{(f)}\rangle \sim \exp (-{\it\lambda}x_{\Vert }/D)$ ) and (4.12) given in § 4. Using these relations, the model for axial conduction in the fluid phase is:

(6.3) $$\begin{eqnarray}\frac{\partial }{\partial x_{\Vert }}\langle I_{f}q_{\Vert }^{{\it\phi}}\rangle =-k_{f}{\it\varepsilon}_{f}\frac{\partial ^{2}\langle {\it\phi}^{(f)}\rangle }{\partial x_{\Vert }\partial x_{\Vert }}\approx -k_{f}{\it\varepsilon}_{f}\langle {\it\theta}^{(f)}\rangle ({\it\lambda}/D)^{2}\exp (-{\it\lambda}x_{\Vert }/D).\end{eqnarray}$$

A comparison of the average axial conduction in the fluid phase from PR-DNS data and the above model is shown in figure 13. For the case of solid volume fraction of 0.1 in figure 13(a), the normalized axial conduction in the fluid phase (normalized by the reference scale $k_{f}/D^{2}$ ) obtained from PR-DNS data shows some scatter about the model prediction in (6.3), although the average trend is captured by the model. The scatter in the PR-DNS data is because of the finite number of realizations and statistical variability in $\langle {\it\phi}^{(f)}\rangle$ , and should reduce with more realizations. Computational resources limit these results to five realizations of the particle configuration. Nevertheless, the model does a fairly good job of capturing the trend in the PR-DNS data. For the case with solid volume fraction ${\it\varepsilon}_{s}=0.4$ shown in figure 13(b), the PR-DNS data show scatter within the length $L=2D$ , but the model still captures the trend of axial conduction in the fluid phase. This results from the fact that at the same Reynolds number, the fluid temperature at high solid volume fraction (0.4) decays faster to approach the particle surface temperature than the one at low solid volume fraction (0.1) (see figure 2).

Figure 14. Contours of the magnitude of the heat flux vector $|I_{f}\boldsymbol{q}^{{\it\phi}}|$ normalized by the reference scale $k_{f}/D^{2}$ at $Re_{m}=5$ for solid volume fraction of ${\it\varepsilon}_{s}=0.1$ (a) and ${\it\varepsilon}_{s}=0.4$  (b).

6.2 Relative importance of fluid-phase axial conduction in average gas–solid heat transfer

We now quantify the relative importance of the fluid-phase axial conduction $\partial \langle I_{f}q_{\Vert }^{{\it\phi}}\rangle /\partial x_{\Vert }$ with respect to average gas–solid heat transfer $\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle$ (see (A 11)) over the range of solid volume fraction and mean slip Reynolds number considered in this work. Since both of terms are spatially inhomogeneous and vary with axial location $x_{\Vert }$ , it is convenient to define volumetric averages of these quantities. We quantify the average volumetric axial conduction in the fluid phase by spatially averaging axial conduction in the fluid phase over the domain length $L$ to obtain:

(6.4) $$\begin{eqnarray}\overline{\langle q_{cond}^{\prime \prime \prime }\rangle }=\frac{1}{L}\int _{0}^{L}\frac{\partial }{\partial x_{\Vert }}\langle I_{f}q_{\Vert }^{{\it\phi}}\rangle (x_{\Vert })\,\text{d}x_{\Vert }.\end{eqnarray}$$

In order to validate the assumption of neglecting axial conduction in the fluid phase in 1-D models that are used to interpret experimental data, we compare this term with average gas–solid heat transfer that is

(6.5) $$\begin{eqnarray}\overline{\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle }=\frac{1}{L}\int _{0}^{L}\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle (x_{\Vert })\,\text{d}x_{\Vert },\end{eqnarray}$$

where $\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle$ is the local average interphase heat transfer rate (see (A 11)). Figure 15 compares the average volumetric axial conduction in the fluid phase $\overline{\langle q_{cond}^{\prime \prime \prime }\rangle }$ with the average gas–solid heat transfer $\overline{\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle }$ at selected volume fractions over a range of Reynolds number values. For a fixed solid volume fraction, the ratio $\overline{\langle q_{cond}^{\prime \prime \prime }\rangle }/\overline{\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle }$ decreases rapidly with increasing mean slip Reynolds number and goes to almost zero at high mean slip Reynolds number of 100. The scaled average volumetric axial conduction also increases with solid volume fraction at each Reynolds number. This results from higher temperature gradients (and heat flux) due to increase in the proximity of particle surfaces at high solid volume fraction. For the case of solid volume fraction of $0.4$ , the average volumetric axial conduction in the fluid phase $\overline{\langle q_{cond}^{\prime \prime \prime }\rangle }$ is approximately 84 % of the average gas–solid heat transfer $\overline{\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle }$ at $Re_{m}=1$ but only 3 % at $Re_{m}=20$ .

Figure 15. Dependence of the ratio of average volumetric axial conduction in the fluid phase to average gas–solid heat transfer $\overline{\langle q_{cond}^{\prime \prime \prime }\rangle }/\overline{\langle q_{{\it\phi}}^{\prime \prime \prime }\rangle }$ on mean slip Reynolds number at solid volume fraction ${\it\varepsilon}_{s}=0.1$ and $0.4$ . The symbols are the values of the ratio from PR-DNS data and error bars indicate 95 % confidence intervals from 5 MIS.

These findings imply that in the low Reynolds number regime, there are high gradients of heat flux in the fluid phase. It is clear that the average volumetric axial conduction in the fluid phase $\overline{\langle q_{cond}^{\prime \prime \prime }\rangle }$ is important only for $Re_{m}<20$ (when compared to average gas–solid heat transfer). Therefore, the neglect of axial conduction in 1-D models that are used to infer the Nusselt number corresponding to average gas–solid heat transfer from inlet/outlet temperature measurements is justified for $Re_{m}>20$ . In the low Reynolds number regime $Re_{m}<20$ (or low Péclet number since in gas–solid flow Prandtl number can be less than the order of one) this assumption is not justified. Now in order to understand the balance of various terms in (1.1) we perform a budget analysis in the following section.