Heat Transfer and Flow Characteristics of Coal Slurries under the Temperature Difference between Inside and Outside of the Channel

Abstract

The pipeline transportation of coal slurries is always subject to a temperature difference between the outdoors environment and the fluid body. As slurries’ viscosity is typically temperature dependent, the flow is accompanied by the heat transfer. In this study, we used the CFD method to investigate temperature distributions and flow structures in straight and curved channels, which has not previously been investigated, according to our knowledge. First, the results demonstrate that the cooling process influences the flow structures along the stream. The fluid turns more sharply in the cooler fluid in the curved channel, the streamlines overlap at an earlier position within the bend, and the velocity maximum zone is wider. Cooling also has a significant impact on transverse flow. Because of the higher viscosity of the more cooled fluid and thus the difficulty of shearing the fluid in the stream-wise direction, the vorticity and strength of the vortex flow are greater. The fluid velocity at the central part of the channel points toward the inner wall at the beginning of the bend, resulting in an inner-wall biased temperature distribution, as the heat transfer is partially carried out by the fluid velocity. The central velocity points toward the outer wall at the end of the bend, resulting in the outer-wall biased temperature profile.

1. Introduction

Pipeline transportation has become the primary mode for transport of coal slurries, such as coal mud tailings or coal-water slurries. Because the rheological properties of these fluids, such as viscosity, are shear rate and temperature dependent due to the particle interaction, and components of clay [1,2,3] or polymers [4,5,6,7], their flow accompanied by the heat transfer process becomes rather complex when there is a temperature difference between the outdoor environment and the fluid body [8,9,10,11], which, in fact, is also encountered in some other geofluids [12,13,14]. Furthermore, the pipelines used for transportation in industries commonly consists of multiple straight tubes and bend sections. These bends induce a secondary vortex flow perpendicular to the main-stream direction, which is also complicated by the heat transfer processes.

In reality, even for Newtonian fluids in straight channels, the heat transfer and the accompanying flow characteristics in both laminar and turbulence regimes have been questioned for a long time [15,16,17,18,19,20,21,22,23,24,25,26,27,28]. Chandratilleke et al. carried our extensive studies on the flow and heat transfer in curved channels [29,30]. In one of their work, they used a helicity function and precisely described the three-dimensional secondary vortex structures in the rectangular channel bend [29]. Their simulations successfully determine the onset of hydrodynamic instability and the effects of the flow rates, rectangular aspect ratio and thermal fluxes on the flow characteristics and the heat transfer. In another work of Chandratilleke, they numerically simulated the flow in a rectangular bent duct with the wall heated, incorporating the buoyance force due to the heating in the formulas [30]. They found that under the wall heated, the Dean vortex in the curved channel is asymmetric compared to the symmetric vortex structure without wall heating, and the secondary vortex flow enhances the convective heat transfer. Yanase et al. used a spectral method and numerically calculated the flow in a curved rectangular duct, with the vertical outer wall heated but the inner wall cooled [31]. They found that at a smaller Grashof number, the flow turns from steady to periodic and then to chaotic pattern. Whereas, at a larger Grashof number, before the flow reaches the chaotic state, it experiences three steady-to-periodic cycles. At the periodic and chaotic states, the convective heat transfer is enhanced significantly by the secondary vortex flow. Zhang et al. numerically studied the flow and heat transfer in the curved wavy ducts [32]. They concluded that there exits an optimal combination of the wave amplitude and the wavy length to maximize the heat transfer. The flow resistance is more sensitive to the wave length compared to the amplitude because the wave length determines the frequency of the flow disturbance. Apart from these, the surface roughness has been proved to affect the heat transfer and flow behaviors. Wu et al. found that the efficiency of the heat transfer increases as the surface roughness or hydrophilicity increases, particularly at high Reynolds numbers [33]. However, as a cost, at these surface conditions, the pressure drop is larger.

For this research topic for non-Newtonian fluids, almost all studies are discussing the flow with heat transfer in geometries without bend. Dong et al. used a density-based topology optimization algorithm to calculate both Newtonian and non-Newtonian thermal fluids in rectangular channels, with the purpose of minimizing the pressure drop but maximizing the heat transfer [34]. The optimal design structures of the channels are different from Newtonian to non-Newtonian fluids, and this difference is maximized below Re = 0.01. Non-Newtonian fluids always perform better in heat transfer, but at the expense of a higher pressure drop. Shojaeian et al. studied the flow of non-Newtonian fluid in circular-cross-section channels under a thermal isoflux imposed on the channel wall, with the consideration that the viscosity and the thermal conductivity are functions of the temperature [35]. They gave an analytical solution of the fully developed flow and found its deviation from the one using the temperature-independent fluid properties that could reach 32%. Sayed-Ahmed et al. used a finite element method to study the flow of power-law fluids with different index numbers $n$ in an annular gap between two coaxial rotating cylinders [36]. The boundary condition of the simulation is an adiabatic and stationary outer wall and an inner wall with a constant temperature and angular velocity. They found that the increase of $n$ or the thermal conductivity causes the decrease of the average Nusselt number. For the shear thinning fluids (index smaller than one), the average Nusselt number is obviously affected by the viscosity parameter. With the elasticity considered further, Nóbrega et al. theoretically calculated the flow of a Phan-Thien–Tanner fluid in a channel with a temperature difference between the wall and the fluid from the inlet [37]. The results correct the expressions of the Nusselt number and the friction factor for the non-Newtonian fluids. The corrections are shown to be dependent less on the viscous dissipation but largely on the elasticity that is represented by the Weissenberg number. Barkhordari et al. used control volume finite difference method to solve the flow of non-Newtonian fluids described by Ostwald–de Waele power law in circular channels, but with the consideration of the wall slip [38]. The results demonstrate that the change on the Nusselt number due to the slip and no-slip is particularly pronouncing for higher power index. As the slip coefficient increases, the center-line velocity decreases but the Nusselt number increases.

The flow of non-Newtonian fluids in channel bends has been studied systematically, but with no heat transfer involved. Li et al. used lattice Boltzmann method to simulate the flow of Herschel-Bulkley fluids in 90° bent circular channels, with the power-law index, yield stress and the curvature radius of the channel varied [39]. They concluded that the flow resistance increases with the power index, yield stress and the curvature radius, while the additional pressure loss due to the bend reduces with the curvature radius. The bend introduces the Dean vortex secondary flow, whose intensity decreases as the yield stress or the power index increases, because the larger viscosity hinders the centrifugal force. Roberts et al. derived an analytical expression for the velocity profile of a Bingham fluid in a two-dimensional curved pipe [40]. As the Bingham number increases, the fluid velocity and stress at the inner wall decreases. For a fixed Bingham number, the width of the unyielded plug is smaller for the larger curvature. Marn et al. numerically solved the flow of shear thickening fluids in a pipe with curvature [41]. Since the flow in the bend part is not fully developed when the curvature radius is small, the numerical friction factor was demonstrated to be much greater than the classical empirical predictions given by White and Mishra-Gupta. They concluded a new experienced expression for the friction factor based on two shear thickening rheological models: the power-law model and Quadratic model, respectively.

In this study, we used CFD approach to simulate the flow of non-Newtonian coal slurries in straight and curved channels at different conveying rates and under different thermal isoflux on the wall. In the formulas employed in the simulations, the thermal conductivity of the fluid is only a function of temperature, while the viscosity is dependent on both shear rate and temperature. In straight channels, we studied the temperature distribution at the cross section, and the pressure drop that is affected by the temperature since the viscosity changes as the heat transfer processes. In the curved channels, we looked at the secondary flow structures in the bend and the impacts of heat transfer, which, to our knowledge, have not been specifically addressed in earlier studies.

2. Simulation Method

The simulations were carried out using the Ansys Fluent software [42]. There were three different types of circular cross-section channels used. One is a one-meter-long straight channel, and the other is 1.5-m-long. The third channel is a U-shaped one (Figure 1a), with two straight arms of 1 m length, a 0.25 m middle curved part of the outer-wall perimeter, and a 0.05 m wall–wall separation between the two straight parts. All of the channels have the same diameter of 55 mm. The structured mesh of the three channels has approximately 600,000 nodes (Figure 1b), and we tested a mesh with more nodes that would produce nearly identical results. The simulations’ boundary conditions are an inlet velocity $V$ and a constant temperature $T_0=300 \mathrm{K}$ for the incoming fluid, zero pressure and an open-to-the-air temperature $T_1$ at the outlet, and no-slip on the wall with a constant temperature $T_1$. Two inlet velocities and two temperatures $T_1=273/ 283 \mathrm{K}$ were imposed for each channel. The Bingham model describes slurry viscosity, and the Arrhenius model [43] reflects temperature influence:

$$\{\begin{array}{c}\eta =\left({\tau }_y/\dot{\gamma }+b\right)\exp \left(\frac{280.394}{T}\right), \left|\tau \right|>{\tau }_y\\ \dot{\gamma }=0, \left|\tau \right|<{\tau }_y\end{array}$$

(1)

Here, ${\tau }_y=7.1 \mathrm{Pa}$ is the yield stress, $b=2.3 \mathrm{Pa}·\mathrm{s}$ is the consistency index, and $T$ is the temperature whose unit is K. The values of ${\tau }_y$ and $b$ are referred from [44], in which the rheological properties of coal slurries are systematically measured using a rheometer. The thermal conductivity of the slurry is referred to [37]:

$$k=k_0\left(k_0^{*}+k_s^{*}T\right),$$

(2)

where $k_0=0.08 \mathrm{W}/\left(\mathrm{mK}\right)$, $k_0^{*}=0.7753$, and $k_s^{*}=0.001181/℃$. The temperature-dependent viscosity $\eta$ and thermal conductivity $k$ were written into Fluent through a User Defined Function. The Reynolds numbers $Re=\frac{D\rho V}{b}$, as defined by [39], are approximately 12 and 24, respectively, for the two flow rates, demonstrating that the turbulence can be neglected. Here, $\rho$ is the density of the fluid and is taken as $1250 \mathrm{kg}/{\mathrm{m}}^3$ in the simulations.

The equations to be solved for the incompressible flow of non-Newtonian fluid are the conservation of mass

$$\frac{\partial u_i}{\partial x_i}=0,$$

(3)

linear momentum

$$\frac{\partial \left(\rho u_i\right)}{\partial t}+\frac{\partial \left(\rho u_j{u}_i\right)}{\partial x_j}=-\frac{\partial P}{\partial x_i}+\frac{\partial {\tau }_{ij}}{\partial x_j},$$

(4)

and thermal energy

$$\frac{\partial \left(\rho c_{p}T\right)}{\partial t}+\frac{\partial \left(\rho c_p{u}_{j}T\right)}{\partial x_j}=-\frac{\partial }{\partial x_j}\left(k\frac{\partial T}{\partial x_j}\right)+{\tau }_{ij}\frac{\partial u_i}{\partial x_j}.$$

(5)

$c_p=4182 \mathrm{J}/\mathrm{kg}℃$ is the specific heat, $u_i$ is the velocity in the $i$ direction, $P$ is the fluid pressure, and ${\tau }_{ij}=\eta \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)$ is the stress tensor.

3. Results and Discussions

3.1. Straight Channels

Figure 2 depicts the pressure in the straight channels along the axial line. Because of the increased viscosity caused by cooling operations, the pressure values at the cooling temperature boundary conditions are frequently slightly higher than those without cooling. The pressure values at the larger temperature difference’s boundary conditions are slightly higher than those at the smaller temperature difference, with an exception of one case in the 1 m channel at $V=0.5 \mathrm{m}/\mathrm{s}$ and a wall temperature of 273 K, as shown by the green line in Figure 2a. Furthermore, as the channel length increases, the pressure rises slightly because the fluid is chilled further, increasing its viscosity. It appears that the fluid inlet velocity has no effect on the pressure readings. However, it has an effect on the temperature distributions, which will be discussed in more detail later.

Figure 3 depicts the temperature profiles at the channel outlet. The more cooled slurry has a slightly narrower central high-temperature core than the less cooled slurry. Again, increasing the channel length causes more heat transfer from the fluid body to the wall, resulting in a more-slender temperature maximum zone. These images also demonstrate the effect of flow rate on heat transfer, which is greater at lower flow rates and is represented by the wider low-temperature ring near the wall. Because the flow in straight channels is always along the axial direction, there is no transverse flow, and thus the heat flux is merely due to thermal conduction rather than the heat carried by the transverse fluid velocity. As a result, the heat transfer increases as the fluid stays for a longer time or at a lower flow rate.

3.2. U-Shaped Channel with a Larger Temperature Difference

In this section, we will demonstrate the results at a temperature difference of 300–273 K between the incoming fluid and the channel wall of the U-shaped channel. Figure 4 shows the pressure and the vorticity at the middle-height plane of the channel at uncooling and cooling temperature boundary conditions and at different flow rates, respectively. Comparing the cooling flow and the uncooled flow at the lower flow rate, there is no obvious difference. However, at the higher flow rate, both the pressure and the vorticity distributions change as the cooling processes. The pressure for the uncooled fluid decreases more gently along the channel, whereas the pressure of the cooled fluid drops more suddenly, especially behind the bend. From this, we can suspect that the cooled and quicker fluid dissipates the energy most at the bend, which can be partly evidenced by the fact that the vorticity of the cooled fluid is generally larger in space than that of the uncooled fluid, although the maximum vorticity of the cooled fluid is smaller. If we compare the data at the same temperature conditions but at different flow rates, the uncooled results demonstrate no obvious change, whereas the cooled fluid does. For the cooled fluid, as the flow rate increases, the pressure drops more sharply, and the vorticity is universally larger, the latter of which suggests that the secondary flow for the quicker fluid is stronger. The following results of the transverse flow structures will verify this.

To represent the flow details in the bend, Figure 5 shows the 3D streamlines, vorticity magnitude, velocity magnitude, and the velocity vector arrows on the middle-height plane. Again, let us first compare the results at the same flow rate but under different temperature conditions. In the cooled fluid, the maximum vorticity value seems to decrease, and the overlap of the streamlines occurs at an earlier position, which corresponds to the earlier twisted velocity vectors. As the fluid processes cool, the maximum velocity zones become wider. Secondly, let us compare the results at the same temperature conditions but at different flow rates. As the flow rate increases, the vorticity at the inner wall within the bend is less symmetric, the maximum velocity zone shrinks, and the turn of the flow is larger. Interestingly, at the higher flow rate, and especially in the cooled fluid, the flow deflects towards the outer wall behind the bend but then recovers slightly back to the inner wall, which is presumably because that the more viscous fluid has more motion in the transverse direction since it is more difficult to be sheared in the stream direction. This will be confirmed by the results of the transverse flow structures discussed in the following.

Figure 6 displays the distributions of the temperature at the five sections and the outlet plane of the U-shaped channel. From the left to the right, the order of the sections is from the inlet to the outlet. In the curved channel, the heat transfer differs from that in the straight channel since there is transverse flow in the curved one and thus the heat transfer also has the component carried by the fluid velocity. First, at the lower flow rate, it can already be observed that the temperature distribution is asymmetric about the inner and outer walls. As the fluid enters the bend, of which the result is shown in the left column, the temperature on the inner wall side is larger than on the outer wall side, since, at this section, the fluid velocity points to the inner wall (Figure 7, left column). At the 135° section (the fourth column), the temperature almost returns to the symmetry, which is pulled back by the outer-wall pointed fluid velocity at the central part of the cross-section. In the 180° section, the temperature is completely biased on the outer wall side. In contrast, when the flow rate is higher, the temperature distribution is more asymmetrical, whether it is toward the inner wall or the outer wall, and the symmetric distribution occurs earlier, at the 90° section, meaning that the outer-wall pointed central flow is stronger than that at the lower flow rate. Another difference is that at the outlet for the quicker flow, the temperature is still biased towards the outer wall due to the less decayed secondary flow.

Figure 7 demonstrates the velocity magnitude and the vector arrows for the five sections of the U-shaped channel. We can observe that, except for the 180° section, the maximum velocity zone is always closer to the inner wall, since at the inner wall side the fluid travels a shorter distance. Initially, let us compare the results at the same flow rate but with the different temperature boundary conditions. First, the cooled fluid shows a much stronger secondary vortex flow regardless of the lower and higher flow rates, which corresponds to the larger vorticity values observed in the cooled fluid (Figure 4). To understand this, the cooled fluid has a larger viscosity and thus the fluid is harder to shear in the main-stream direction. However, the fluid’s velocity has to be that large. As a result, more velocity is in the transverse direction and converts to the secondary flow. Second, we compare the results at the same temperature boundary conditions but different flow rates. Obviously, as the flow rate increases, the secondary flow is stronger, and the deflection of the velocity maximum zone is larger.

Figure 8 shows the vorticity profiles in the five sections of the channel. In all sections, the maximum vorticity is always biased towards the inner wall. At the beginning of the bend, the flow points to the inner wall, whereas in the rest sections, the central flow points to the outer wall. At the end of the bend and at the higher flow rate, the vortex flow is still there, whereas at the lower flow rate, the vortex has decayed out. Referring to the cooling, again, the same as the results shown in Figure 4, the vorticity values in the cooled fluid are larger, and the streamlines are denser, meaning the stronger vortex flow.

3.3. U-Shaped Channel with a Smaller Temperature Difference

This part describes the heat transfer and the flow characteristics in the U-shaped channel, at the temperature boundary conditions of 283 K on the wall and 300 K of the incoming fluid at the inlet. Figure 9 shows the pressure and the vorticity magnitude at the middle-height plane of the channel. The inlet–outlet pressure drop values show no change as the wall temperature changes, compared to Figure 4. However, the vorticity values change. At the lower flow rate, the vorticity shows a narrower minimum zone compared to the one in the more cooled fluid, meaning that the fluid in this case is less rigid since the fluid is less cooled and thus the viscosity is smaller. At the higher flow rate, the maximum vorticity is larger than that in which the fluid cools more, again showing that the fluid is less rigid. Another noticeable phenomenon is that, at both lower and higher flow rates, the vorticity is less symmetric than that of the more cooled fluid, which also reflected this. From the streamlines, vorticity, and velocity magnitudes/vectors displayed in Figure 10, the less cooled fluid shows a more asymmetric vorticity distribution at both flow rates. Moreover, the flow turns more gently and the maximum velocity zone is more oblique, together with the later overlap of the streamlines in the curved part, meaning that the deflection of the less cooled fluid is easier and smoother. Ass the flow is faster, it seems that the maximum velocity zone is even oblique. The disappearing velocity vectors at the strongly curved points indicate that there should be a strong cross-section secondary flow at these positions.

The temperature profiles at the five sections and the outlet plane of the U-shaped channel are shown in Figure 11. Overall, compared to the cases of the larger temperature difference, the distribution is more symmetrical because of fewer heat transfer processes, and at the higher flow rate, it is even more symmetrical. Again, at the beginning of the bend, the temperature at the inner wall side is higher due to the inner-wall-pointed velocities (Figure 12), and thus the heat is transferred from the outer wall to the inner wall. From the 0° section to the 180° section, the temperature profile biased to the inner wall converts to the profile biased to the outer wall due to the stronger and stronger outer-wall-pointed flow at the central part of the sections. At the higher flow rate, at the outlet plane of the channel, the profile is still biased to the outer wall, which is attributed to the vortex flow that has not decayed completely.

The velocity magnitudes, vectors, and the vorticity magnitudes and streamlines are shown in Figure 12 and Figure 13. As the fluid turns more gently and it starts to turn at an earlier point, the velocity maximum zone at the 0° section deflects less towards the inner wall, particularly at the higher flow rate, and similarly, the flow at the 145° section deflects less towards the outer wall. However, the maximum velocity values are not changed as the temperature changes, but the maximum velocity zone is obviously larger compared to the results in the more cooled fluid at the same flow rates. At the end of the bend, at the higher flow rate, the flow temporarily deflects towards the outer wall as an inertia, and this deflection is closer to the outer wall compared to the flow in the more cooled fluid at the same flow rate, due to the smaller viscosity and thus the easier twist. From the results of the vorticity and the streamlines, the vorticity on the cross sections in the less cooled fluid is smaller and the streamlines at the channel center are generally looser because the fluid’s viscosity is smaller and the fluid is easier to shear in the stream-wise direction.

3.4. Pressure Loss and Energy Dissipation Due to the Bend

In industries, it is important to accurately predict the pressure drop in the pipe bend for optimizing the design of the complex geometries of the pipes. For the fully developed flow in the bend with a large curvature radius, the friction coefficient has a ready-use expression. However, in the bend with a small curvature radius, such as the one in the present article, the flow is not fully developed. Moreover, our fluid is plastic and there is an unyielded plug at the center of the channel. All these factors may induce abnormal phenomenon, and the existing formula of the friction may not be applicable. Therefore, this section will analyze the pressure loss, friction coefficient, and the energy dissipation rate induced by the bend, to compare with the previous studies and to provide comprehensive understandings.

The extra pressure loss induced by the bend is the difference between the total pressure drop and that of a straight pipe with a same length, which is defined by [41]:

$$\mathsf{\Delta }{P}_{loss}=\mathsf{\Delta }{P}_t-\frac{dP}{dx}\left(L_s+L_b\right).$$

(6)

$\mathsf{\Delta }{P}_t$ is the total pressure drop of the U-shaped channel, $\frac{dP}{dx}$ is the pressure gradient of the straight channel with the same length, $L_s$ is the length of the straight part of the U-shaped channel, $L_b$ and is the length of the bend part. Here, we employ a dimensionless pressure loss coefficient $\xi =\mathsf{\Delta }{P}_{loss}/0.5\rho V^2$ to represent the influences of the flow rate and the temperature. The values of $\xi$ are summerized in Table. 1. When calculating $\frac{dP}{dx}$, we used the following relation between $\frac{dP}{dx}$ and $V$ for the Bingham fluid [45]:

$$V\left(r\right)=\{\begin{array}{c}\frac{1}{2}\left(\frac{\mathsf{\Delta }P}{2bL}\right){\left(R-R_p\right)}^2 r<r_p\\ \frac{1}{2}\left(\frac{\mathsf{\Delta }P}{2bL}\right){\left(R-R_p\right)}^2\left[1-{\left(\frac{r-R_p}{R-R_p}\right)}^2\right] otherwise \end{array}.$$

(7)

$R$ and $R_p=\frac{2L{\tau }_y}{\mathsf{\Delta }P}$ are the channel radius and the radius of the unyielded plug, respectively. $V=\frac{1}{R}\underset{0}{\overset{R}{{\displaystyle \int }}}V\left(r\right)$ is the average flow velocity. As

Table 1. Pressure loss coefficient in different cases.

$\mathit{\xi }$	0.5 m/s	1 m/s
300/300 (K)	39.29	11.18
283/300 (K)	44.18	13.03
273/300 (k)	47.23	14.18

shows, $\xi$ surprisingly decreases with the flow rate, and increases as the fluid cools more because of the increase in the viscosity. According to the previous studies of laminar flow of non-Newtonain fluids in a bend [32], it is reasonable to suppose that the pressure drop in the bend is almost totally associated to the energy dissipation rate $\epsilon$ due to the viscosity. Herein, we took $\epsilon$ along the horizontal diameter of the cross section at the middle of the channel as the average to calculate the total energy dissipation rate in the whole bend part:

$$\epsilon =\frac{Vo{l}_b}{D}{{\displaystyle \int }}_0^D\eta {\left|\overrightarrow{\omega }\right|}^{2}dl,$$

(8)

where $\overrightarrow{\omega }$ is the vorticity vector and $Vo{l}_b$ is the volume of the whole bend part of the channel. Moreover, we made $\epsilon$ dimensionless through dividing it by $Q\mathsf{\Delta }{P}_{loss}$, in which $Q$ is the volumetric flow rate. The values of $\epsilon /Q\mathsf{\Delta }{P}_{loss}$ in different cases are listed in

Table 2. Energy dissipation rate in different cases.

$\mathit{\epsilon }/\mathit{Q}\mathsf{\Delta }{\mathit{P}}_{\mathit{l}\mathit{o}\mathit{s}\mathit{s}}$	0.5 m/s	1 m/s
300/300 (K)	1.1080	1.0160
283/300 (K)	1.1187	1.0235
273/300 (k)	1.1254	1.0281

, and it shows the same tendency as $\xi$, which confirms our supposition that the viscosity mainly dissipates the energy in the bend. Actually, in the early studies, it was widely accepted that the energy loss due to the bend, which was sometimes represented by the friction coefficient $f_c$, increases as the Dean number $De$ that is propertional to the fluid velocity. For example, Singh and Mishra gave the expressions of the friction coefficient: $\frac{f_c}{f_p}={\left[1-{\left\{1-{\left(\frac{11.6}{De}\right)}^{0.45}\right\}}^{2.222}\right]}^{-1}$ and $\frac{f_c}{f_p}=1+0.033{\left(logDe\right)}^{4.0}$ for Newtonian fluids [46]. Marn concluded $\frac{f_c}{f_p}=1+K_1{\left(logDe\right)}^{K_2}$ for power-law fluids [41]. Here, $f_p$ is the friction coefficient of the straight pipe with a same length. $De=Re\sqrt{R/R_c}$, in which $R_c$ is the curvature radius. However, in the present study, we obtained the converse relation between $f_c$ and $De$, as observed more clearly in the plot of Figure 14. We speculate that this is due to two effects: one is the special rheology of our fluid, particularly the plasticity–yield stress; the other is the not fully developed flow in the strong curvature. At the lower flow rate, the yield stress fluid flows with an almost unyielded plug at the central part of the pipe, resulting in a nearly infinitely large viscosity [45]. Nonetheless, there is still a weak Dean vortex secondary flow inside the plug [47]. As the flow rate increases, the width of the plug shrinks, and in the places where the plug locates at the lower flow rate, the infinite viscosity sharply reduces to an intermediate value. In contrast, the magnitude of the squared vorticity increases by a smaller magnitude. Remember that in our case, the vorticity always weighs more on the inner-wall side of the channel (Figure 8 and Figure 13). Compared to the outer-wall-weighted vorticity in the fully developed flow in a smaller curvature, the vorticity in our case is relatively smaller and hindered due to the outer-wall pointed fluid velocity. Therefore, the decrease in the viscosity is greater than the increase in the squared vorticity, and the variation of the viscous energy dissipation rate $\epsilon ~\eta {\left|\overrightarrow{\omega }\right|}^2$ with the flow rate follows that of the viscosity, which decreases as the flow rate increases. At least, this is true in the range of the flow rates in the present study. It is uncertain to predict if the relation between the energy loss and the flow rate would return to the positive correlation at even higher flow rates, at which the decrease of the viscosity might be less than the increase in the squared vorticity magnitude. Because in engineering, the coal slurries are not conveyed at higher rates, that range is not discussed here. Moreover, the turbulence begins to arise at the higher rates, and the energy loss would be not only due to the viscosity, which would make the problem more complicated.

We also plot the friction coefficient $\frac{f_c}{f_p}$ as a function of the Nusselt number $Nu$ that represents the ratio of the convective heat transfer to the conductive transfer, shown in Figure 14. When calculating $Nu$, similar to $\epsilon$, we selected the horizontal radius near the outer-wall side of the cross section at the middle of the U-shaped channel as representative of the whole bend part:

$$Nu=\frac{|\rho C_p\mathsf{\Delta }T\overrightarrow{u|}}{\left|k\mathsf{\nabla }T\right|}.$$

(9)

$\overrightarrow{u}$ is the average velocity in the main heat transfer direction (positive x) outside the temperature-plug region (for example, Figure 6 and Figure 11). $\mathsf{\Delta }T$ is the temperature difference between the fluid body and the pipe wall. $\mathsf{\nabla }T$ is the temperature gradient in x direction. As shown in Figure 14, $\frac{f_c}{f_p}$ again decreases as $Nu$ increases, since $Nu$ is similar to $De$, which is proportional to the flow rate. From both plots in Figure 14, we notice that the deviation of $\frac{f_c}{f_p}$ from different thermal fluxes is reduced at larger $De$ or $Nu$, implying that the pressure loss is more affected by the temperature in the slower flow.

3.5. Inverse Thermal Flux from the Wall to the Fluid

In some real circumstances, the slurry’s temperature is lower than the environment due to the hot weather, and as a result, the thermal flux flows from the pipe wall to the inner slurry. We simulated such a case, in which the temperature on the pipe wall is constantly 300 K whereas the slurry coming from the inlet is 283 K. The inlet velocity is 0.5 m/s. Figure 15 shows the temperature distributions at the five slices in the bend and the outlet. Under such a thermal boundary condition, the heat flux is from the wall to the fluid. Therefore, opposite to the cases in the cooler pipe, the lower temperature is biased to the inner wall at the beginning of the bend, while it is biased to the outer wall at the end of the bend.

The velocity magnitudes and the vectors at the five cross sections within the bend are shown in Figure 16. Compared with the results in the case of the inverse temperature flux at the same flow rate (Figure 12, upper row), because of the lower temperature and thus the higher viscosity of the fluid in this case, the central velocity becomes smaller. The vorticity distributions (Figure 17) are not quite different from the case with the inverse thermal flux (Figure 13, upper row). There is only a minor deviation that the streamlines are denser near the channel edge in in this case with 300/283 K, especially at the 45° section. The same phenomenon is also observed, comparing the cases of 273/300 K and 300/300 K at the same flow rates (Figure 8). This is mainly due to the almost unyielded plug flow at the center of the cross section. The fluid’s temperature and the viscosities are lower in the case of 300/283 K compared with the case of 283/300 K. Thus, the plug is wider and the fluid within the plug is more rigid, resulting in a greater difficulty of creating a rotational flow inside, and therefore the most rotational flow is outside the plug and near the edge.

Comparing both cases of 283/300 K and 300/283 K with the case of 300/300 K (Figure 7, upper row), it is found that the deflection of the minimum or the maximum vorticity core is switched to different sides from the case with zero thermal flux to the cases with non-zero thermal flux, regardless of the dirction of the heat transfer. Under a temperature difference, at the 0° section, the minimum vorticity is always biased to the outer-wall side, whereas it is closer to the inner-wall side with no thermal flux. The same trend is observed at the 180° section. In conclusion, at both 0° and 180° sections, the maximum and minimum vorticity zones are always concentrated on the same side under zero thermal flux, or at the smaller fluid viscosity. On the other hand, with the temperature difference boundary conditions or the larger fluid viscosity, the maximum and minimum vorticity zones are less concentrated and thus distributed on the opposite sides. To understand this, when the fluid’s viscosity is small, its flow structures, including the vorticity maximum or the minimum zones, are easier to change within the space, and the characteristic length of the structure is smaller [29,30].

4. Conclusions

In this work, we used the CFD method to study the flow of slurries described by a Bingham model subjected to a temperature difference between the channel boundary wall and the inlet incoming fluid. Since the flow is accompanied by the heat transfer, its flow characteristics, including the flow structure, are obviously different from those without a temperature difference. As a conclusion, we obtained the following important points:

• As the fluid cools, the viscosity increases, and the inlet–outlet pressure drop increases. And no matter the shapes of the channels, as the length of the channel or the temperature difference increases, or as the flow rate decreases, the heat transfer increases, which is represented by a narrower maximum temperature zone at the channel outlet.
• In the straight channels, the heat transfer is merely the part because of the thermal conductivity, whereas in the curved channels, the heat transfer also has the component that the fluid velocity carries. As a result, at the beginning of the channel bend, the temperature biases towards the inner wall, since at this cross-section, the central velocity points towards the inner wall. In contrast, in the sections on the second half of the semicircle of the bend, the temperature biases toward the outer wall, which again comforts the direction of the fluid velocity.
• Temperature influences the flow structure in the stream-wise direction in this way: in the more cooled fluid, the maximum velocity zone is wider than that in the less cooled fluid because of the larger viscosity; the flow turns more sharply near the bend; the overlap of the streamlines occurs at the earlier position within the bend.
• The temperature influences the flow structure in the transverse direction in this way: the rotational flow in the more cooled fluid is much stronger, particularly near the edge of the channel, which is represented by the larger vorticity magnitude and the stronger strength of the vortex flow and the denser streamlines. Furthermore, when the fluid is in a higher temperature or has a smaller viscosity, the maximum and the minimum vortex cores are easily concentrated at the same side of the cross section. In contrast, under a lower temperature, the maximum and the minimum vorticity cores are separated to opposite sides. In summary, we attribute these effects to the wider unyielded central plug in the lower temperature and the larger characteristic length of the flow structures with a smaller viscosity.
• The pressure loss coefficient is found to decrease as the flow rate increases, and the energy dissipation rate due to the viscosity follows the same trend, implying that the pressure loss mainly results from the viscous energy dissipation in the laminar regime. This converse relation is against the previous studies because of the yield stress of our fluid and the not fully developed flow in the strongly curved bend, in which the vorticity weighs more on the inner-wall side and is hindered by the outer-wall pointed velocities. In this case, as the flow rate increases, the vorticity magnitude increases by a smaller magnitude compared to that in the fully developed flow, where the vorticity weights more on the outer-wall side.