2.2.1. Pseudo-Homogeneous *λ*eff,r-*α*<sup>w</sup> Model

The pseudo-homogeneous two-dimensional plug flow heat transfer model under steady-state conditions is described by:

$$
\rho\_{\text{f}} \mathbf{c}\_{\text{p},\text{f}} \mu\_{\text{z}} \frac{\partial T}{\partial z} = \lambda\_{\text{eff},\text{r}} \left( \frac{\partial^2 T}{\partial r^2} + \frac{1}{r} \frac{\partial T}{\partial r} \right) + \lambda\_{\text{eff},\text{z}} \frac{\partial^2 T}{\partial z^2} \tag{1}
$$

whereas the following boundary conditions are used:

$$-\lambda\_{\rm eff,r} \frac{\partial T}{\partial r} = a\_{\rm w} (T - T\_{\rm w}) \tag{2}$$

$$\frac{\partial T}{\partial r} = 0 \tag{3}$$

$$\begin{cases} T = T\_0 \\ \dots \end{cases} \tag{4}$$

$$\frac{\partial T}{\partial z} = 0 \tag{5}$$

$$\text{at } z = L. \tag{5}$$

Here, *u*<sup>z</sup> = *u*0/*ε* is the constant interstitial velocity, *ρ*<sup>f</sup> the fluid density, and *c*p,f the specific heat of the fluid. The *λ*eff,r-*α*<sup>w</sup> model lumps all radial heat transfer mechanisms into a constant effective radial thermal conductivity *λ*eff,r. The steep temperature drop at the tube walls is modeled by the introduction of an artificial wall heat transfer coefficient *α*w, using Equation (2). The thermal conductivity in axial direction can be assumed to be equal to the stagnant effective thermal conductivity *λ*eff,z = *λ*<sup>0</sup> eff,r, or it can be neglected if the system is dominated by convective effects. The model itself has been critically discussed by many authors [31,45,46], but nevertheless widely spread due to its simplistic nature.

It was found by Yagi and Kunii [47] that the radial effective thermal conductivity can be expressed as:

$$\frac{\lambda\_{\rm eff,r}}{\lambda\_{\rm f}} = \frac{\lambda\_{\rm eff,r}^{0}}{\lambda\_{\rm f}} + \frac{1}{P\varepsilon\_{\rm f,r}(\infty)} PrR\varepsilon\_{\rm P}.\tag{6}$$

The first term on the right-hand side is the effective radial thermal conductivity of the stagnant bed. A huge number of correlations exists to determine *λ*<sup>0</sup> eff,r, which have been reviewed by van Antwerpen et al. [48]. Based on a unit cell approach, Zehner and Schlünder [49] derived the following correlation that is widely used:

$$\frac{\lambda\_{\text{eff},r}^{0}}{\lambda\_{\text{f}}} = \left(1 - 1\sqrt{1 - \varepsilon}\right) + \frac{2\sqrt{1 - \varepsilon}}{1 - \kappa^{-1}\overline{B}} \cdot \left(\frac{\left(1 - \kappa^{-1}\right)B}{\left(1 - \kappa^{-1}B\right)^{2}} \ln\left(\frac{1}{\kappa^{-1}B}\right) - \frac{B + 1}{2} - \frac{B - 1}{1 - \kappa^{-1}B}\right) \tag{7}$$

Here, *κ* is the ratio of solid to fluid thermal conductivity and *B* is the deformation parameter, which is related to the void fraction by *B* = 1.25((1 − *ε*)/*ε*) 10/9. The correlation can be further extended by incorporating secondary effects like radiative heat transfer or the effect of particle–particle contacts on the heat transfer. For a more detailed description, the interested reader is referred to the work of Tsotsas [50] and van Antwerpen et al. [48].

For the heat transfer coefficient at the wall, Yagi and Kunii [51] proposed the following correlation for *Nu*<sup>w</sup> = *α*w*d*p,v/*λ*f:

$$Nu\_{\rm W} = Nu\_{\rm w}^{0} + \frac{1}{(1/Nu\_{\rm w}^{\*}) + (1/Nu\_{\rm m})} \tag{8}$$

using:

$$Nu\_{\rm in} = 0.054 Pr Re\_{\rm P} \tag{9}$$

$$Nu\_{\rm w}^{\*} = 0.3 Pr^{1/3} Re\_{\rm P}^{3/4} \tag{10}$$

$$Nu\_{\rm W}^{0} = \left(1.3 + \frac{5}{N}\right) \frac{\lambda\_{\rm eff,r}^{0}}{\lambda\_{\rm f}}.\tag{11}$$

Nilles and Martin [52,53] developed the following correlation that is widely used:

$$Nu\_W = \left(1.3 + \frac{5}{N}\right) \frac{\lambda\_{\rm eff,r}^0}{\lambda\_{\rm f}} + 0.19 Pr^{1/3} Re\_{\rm p}^{3/4}.\tag{12}$$

According to Dixon [31] two methods are commonly used to determine *λ*eff,r and *α*w. The first option is a parameter estimation done by conducting an optimization study based on the *λ*eff,r-*α*<sup>w</sup> model, whereas the objective is to minimize the sum of squared error regarding the radial temperature profile at one or more axial positions. Alternatively, the method described by Wakau and Kaguei [54] can be used. This method is based on the approximate solution of the pseudo-homogeneous *λ*eff,r-*α*<sup>w</sup> model and allows determining *λ*eff,r and *α*<sup>w</sup> from the axial temperature profile in the core of the bed and the average outlet temperature. Both can easily be extracted from the particle-resolved simulation results.

The latter method was used in this study and presented in great detail by Wakao and Kaguei [54]. By neglecting the axial thermal conductivity, the analytical solution of Equation (1) is:

$$\frac{T\_\mathrm{w} - T}{T\_\mathrm{w} - T\_0} = 2 \sum\_{n=1}^\infty \left( \frac{J\_0(2a\_n r/D) \exp\left(-a\_n^2 y\right)}{a\_n \left(1 + \left(a\_n / Bi\right)^2\right) l\_1(a\_n)} \right). \tag{13}$$

Here, *r* is the radial position and *Bi* the Biot number *Bi* = *<sup>α</sup>*w*<sup>D</sup>* <sup>2</sup>*λ*eff,r . *an* is the n-th root of the following equations that include the Bessel function of the first kind and zeroth-order *J*<sup>0</sup> and the first kind and first-order *J*1:

$$
\delta B \mathfrak{i}\_0(a\_\mathfrak{n}) = a\_\mathfrak{n} \mathfrak{i}\_1(a\_\mathfrak{n}).\tag{14}
$$

The parameter *y* is expressed by:

$$y = \frac{\lambda\_{\text{eff},r}z}{\rho\_{\text{f}}u\_{\text{z}}c\_{\text{p},\text{f}}(D/2)^2},\tag{15}$$

whereas *z* is the axial position, *ρ*<sup>f</sup> the fluid density, and *c*p,f its specific heat. Deep in the bed, when *y* ≥ 0.2, the first term in the series of Equation (13) becomes predominant, leading to:

$$\frac{T\_\mathrm{W} - T}{T\_\mathrm{W} - T\_0} = \frac{2J\_0(2a\_1 r/D) \exp\left(-a\_1^2 y\right)}{a\_1 \left(1 + \left(a\_1 / Bi\right)^2\right) J\_1(a\_1)}.\tag{16}$$

with:

$$
\hat{a}I\!\_0(a\_1) = a\_1 I\_1(a\_1).\tag{17}
$$

In the center of the bed (*r* = 0 and *T* = *T*core), Equation (16) is reduced to:

$$\frac{T\_\mathrm{W} - T\_\mathrm{core}(z)}{T\_\mathrm{W} - T\_0} = \frac{2\exp\left(-a\_1^2 y\right)}{a\_1 \left(1 + \left(a\_1 / Bi\right)^2\right) I\_1(a\_1)}.\tag{18}$$

If Equation (18) is logarithmized, it gives:

$$\ln\left(\frac{T\_{\rm w} - T\_{\rm core}(z)}{T\_{\rm w} - T\_0}\right) = -a\_1^2 \left(\frac{\lambda\_{\rm eff,r}}{\rho\_1 u\_{\rm x} c\_{\rm p,t} \left(D/2\right)^2}\right) z + \ln\left(\frac{2}{a\_1 \left(1 + \left(a\_1/Bi\right)^2 f\_1(a\_1)\right)}\right). \tag{19}$$

It was shown by Wakao and Kaguei [54] that the following relationship for the average outlet temperature *T*m is valid for a reasonably large axial position:

$$\frac{T\_{\rm w} - T\_{\rm m}}{T\_{\rm W} - T\_{\rm core}} = \frac{2J\_1(a\_1)}{a\_1}.\tag{20}$$

From Equation (20), *a*<sup>1</sup> can be solved iteratively, and *λ*eff,r can be calculated from the slope of Equation (19), subsequently. The wall heat transfer coefficient can either be determined from the intercept of Equation (19) or from Equation (17). The latter method was promoted by Wakao and Kaguei [54], since the authors argued that *α*<sup>w</sup> is very sensitive to slight changes of the intercept.

Both methods were tested during this study. A sensitivity test was conducted based on the particle-resolved CFD results for a packing of spherical particles at *Re*<sup>p</sup> = 100. The sensitivity analysis of *α*<sup>w</sup> and *λ*eff,r towards the accounted temperature range was conducted, whereas the range of Θcore = (*T*core − *T*0)/(*T*<sup>w</sup> − *T*0) was varied as follows: {Θcore <sup>∈</sup> <sup>R</sup> : 0.05 <sup>≤</sup> <sup>Θ</sup>core <sup>≤</sup> 0.4 <sup>∧</sup> 0.6 <sup>≤</sup> <sup>Θ</sup>core <sup>≤</sup> 0.95}. It was found that *<sup>λ</sup>*eff,r had a relative standard deviation (RSD) of ±6 %. Then, from the intercept of Equation (19) calculated, the values of *α*<sup>w</sup> had a very low RSD of ±3 %, while the suggested method

of Wakao and Kaguei increased the RSD to ±15 %, which was in contrast to the authors' argumentation. Nevertheless, a huge discrepancy in *α*<sup>w</sup> was found: the values of the intercept-method were up to three times lower in comparison to the values determined from Equation (17). A comparison of particle-resolved CFD results against the results of the two-dimensional plug flow model in terms of axial and radial temperature profiles revealed that the temperature profiles were mispredicted if the intercept method was used, while the method promoted by Wakao and Kaguei led to reasonable results. Therefore, as Wakao and Kaguei did, we also highly recommend calculating *α*<sup>w</sup> from Equation (17) instead of the intercept of Equation (19). To evaluate *λ*eff,r and *α*w, the axial temperature profile was limited to 0.2 ≤ Θ ≤ 0.8 in this work.
