The effect of Prandtl number on turbulent sheared thermal convection

In turbulent wall sheared thermal convection, there are three different flow regimes, depending on the relative relevance of thermal forcing and wall shear. In this paper we report the results of direct numerical simulations of such sheared Rayleigh-B\'enard convection, at fixed Rayleigh number $Ra=10^6$, varying the wall Reynolds number in the range $0 \leq Re_w \leq 4000$ and Prandtl number $0.22 \leq Pr \leq 4.6$, extending our prior work by Blass et al. (2020), where $Pr$ was kept constant at unity and the thermal forcing ($Ra$) varied. We cover a wide span of bulk Richardson numbers $0.014 \leq Ri \leq 100$ and show that the Prandtl number strongly influences the morphology and dynamics of the flow structures. In particular, at fixed $Ra$ and $Re_w$, a high Prandtl number causes stronger momentum transport from the walls and therefore yields a greater impact of the wall shear on the flow structures, resulting in an increased effect of $Re_w$ on the Nusselt number. Furthermore, we analyse the thermal and kinetic boundary layer thicknesses and relate their behaviour to the resulting flow regimes. For the largest shear rates and $Pr$ numbers, we observe the emergence of a Prandtl-von Karman log-layer, signalling the onset of turbulent dynamics in the boundary layer. Finally, our results allow to extend the Grossmann-Lohse theory for heat transport in Rayleigh-B\'enard convection to the sheared case, universally describing $Nu(Ra,Pr,Re_w)$.


Introduction
Buoyancy and shear are crucial processes in fluid dynamics and key for many flow related processes in nature and technology. A paradigmatic example of buoyancy driven flow is Rayleigh-Bénard (RB) convection, a system where the fluid is heated from below and cooled from above (Ahlers et al. 2009;Lohse & Xia 2010;Chilla & Schumacher 2012;Xia 2013). The flow is controlled by the Rayleigh number Ra = βgH 3 ∆/(κν), which quantifies the the non-dimensional temperature difference between the two horizontal plates. Here, H is their distance, β the thermal expansion coefficient of the fluid, g the gravitational acceleration, ∆ the temperature difference across the fluid layer, κ and ν the thermal diffusivity and kinematic viscosity, respectively. Furthermore, the Prandtl number is defined as P r = ν/κ, which is the ratio between momentum and thermal diffusivities. An important output of the flow is the heat transport between the plates, which can be non-dimensionally quantified by the Nusselt number N u = QH/(κ∆), with Q = w T − κ∂T /∂z the mean vertical heat flux (w and T are vertical velocity and temperature fluctuations, z is the vertical direction).
On the other hand, for flows driven by wall shear stress, a commonly used model problem is the Couette flow (Thurlow & Klewicki 2000;Barkley & Tuckerman 2005;Tuckerman & Barkley 2011). We adopt a geometry in which the bottom and top walls slide in opposite directions with a wall-tangential velocity u w and the forcing can be expressed non-dimensionally by the wall Reynolds number Re w = Hu w /ν. The relevant flow output is now the wall friction, quantified by the friction coefficient C f = 2τ w /(ρu 2 w ), with ρ the fluid density and τ w the surface-and time-averaged wall shear stress. Turbulent Couette flow is dominated by large-scale streaks (Lee & Kim 1991;Tsukahara et al. 2006;Kitoh & Umeki 2008;Pirozzoli et al. 2011Pirozzoli et al. , 2014Orlandi et al. 2015;Chantry et al. 2017). These remain correlated in the streamwise direction for a length up to about 160 times the distance between the plates (Lee & Moser 2018).
Combining both, buoyancy and wall shear forcings, yields a complex system that is relevant in many applications, especially for atmospheric and oceanic flows (Deardorff 1972;Moeng 1984;Khanna & Brasseur 1998). Also in sheared thermal convection large-scale structures emerge, as experiments have shown (Ingersoll 1966;Solomon & Gollub 1990). Investigations on channel flows with unstable stratification (Fukui & Nakajima 1985) revealed that temperature fluctuations in the bulk decrease while velocity fluctuations close to the wall increase for stronger unstable stratification.
Numerical simulations of wall sheared convection (Hathaway & Somerville 1986;Domaradzki & Metcalfe 1988) have revealed that adding shear to buoyancy increases the heat transport for low Ra, but causes also the large-scale structures to weaken thus decreasing the heat transport for Ra 150.000. Similar phenomena have been observed in Poiseuille-RB, where the wall parallel mean flow is driven by a pressure gradient rather than the wall shear: in this case the N u decrease was attributed to the disturbance of the longitudinal wind on the thermal plumes (Scagliarini et al. 2014(Scagliarini et al. , 2015Pirozzoli et al. 2017). This plume-sweeping mechanism, causing a Nusselt number drop, was also observed in Blass et al. (2020), who report very long, thin streaks, similar to those of the atmospheric boundary layer where these convection rolls are called cloud streets (Etling & Brown 1993;Kim et al. 2003;Jayaraman & Brasseur 2018).
In both flows, Couette-RB and Poiseuille-RB, the ratio between buoyancy and mechanical forcings can be best quantified by the bulk Richardson number which is a combination of the flow governing parameters Ra, Re w and P r. In the Couette-RB flow of Blass et al. (2020), Ri was in fact used to distinguish between three different flow regimes, namely thermal buoyancy dominated, transitional, and shear dominated, similarly to the case of stably stratified wall turbulence, where Zonta & Soldati (2018) distinguish between the buoyancy dominated, buoyancy affected and turbulence dominated regimes. Indeed, sheared stably or unstably stratified flows are present in many different situations involving both liquids and gases. Therefore the fluid properties, as reflected in the Prandtl number, play a major role (Chong et al. 2018). In the atmosphere it results P r = O(1) while in ocean dynamics P r = O(10). However, a much larger P r variation is found in industrial applications. E.g. P r ≈ O(10 −3 ) for liquid metals (Teimurazov & Frick 2017), which are for example in use for cooling applications in nuclear reactors Figure 1. Phase diagram of simulation runs. We show two panels to better illustrate our choice of simulation input parameters, which were determined based on Rew (left panel) and Ri (right panel). Rew = 2000, 3000, 4000 were chosen to be consistent with Blass et al. (2020) and to cover the shear dominated regime. The squared symbols show the datapoints for Rew = 0 for completeness and independently of the y-axis, since they cannot be directly included in the logarithmic scale. To have a sufficient amount of data in the thermal buoyancy dominated regime, we picked Ri = 100 as most thermal dominated case and then logarithmically spaced three more datapoints. (Usanov et al. 1999) or P r ≈ O(10 3 ) for molten salts or silicone oils (Vignarooban et al. 2015) for high-performance heat exchangers.
Despite this staggering range of Prandtl numbers encountered in real applications, the vast majority of studies on sheared, thermally stratified flows have been performed only at P r = O(1). To overcome this limitation, in this paper we extend the work of Blass et al. (2020) for P r = 1 by analysing the parameter space 0 Re w 4000 and 0.22 P r 4.6 while keeping the Rayleigh number constant at Ra = 10 6 (see figure 1 for the complete set of simulations).
The present study can be considered similar and complementary to that of Zhou et al. (2017) who carried out numerical simulations with a large P r variation for a stably stratified Couette flow.
The manuscript is divided in the following manner. Section 2 briefly reports the numerical method. Section 3 focusses on the global transport properties and section 4 on the boundary layers. The paper ends with conclusions (section 5).

Numerical method
The three-dimensional incompressible Navier-Stokes equations with the Boussinesq approximation are integrated numerically. Once non-dimensionalised, the equations read: with u and θ the velocity, normalized by √ gβ∆H, and temperature, normalized by ∆, respectively. t is the time normalized by H/(gβ∆) and P the pressure in units of gβ∆/H. Equations (2.1) and (2.2) are solved using the AFiD GPU package (Zhu et al. 2018) which bases on a second-order finite-difference scheme (van der Poel et al. 2015). The code has been validated and verified several times (Verzicco & Orlandi 1996;Verzicco & Camussi 1997Stevens et al. 2010Stevens et al. , 2011Ostilla-Mónico et al. 2014;Kooij et al. 2018). We use a uniform discretization in horizontal, periodic directions and a nonuniform mesh, with an error function-like node distribution in the wall-normal direction.
Following Blass et al. (2020), we performed our simulations in a 9πH × 4πH × H domain, which are the streamwise, spanwise and wall-normal directions, respectively. The grid resolutions are also based on Blass et al. (2020) and then further modified to account for the Prandtl number variation in this study.

Organization of turbulent structures
Using as guideline the description of Blass et al. (2020) we observe that also in the present case the flow can be classified in buoyancy dominated, transitional and shear dominated regimes (see figure 2 and Table 1 for a full overview). As shown in Blass et al. (2020), for P r = 1 and increasing Re w , we observe the thermal buoyancy dominated regime at Re w = 0 while already at Re w = 1000, 2000 the compact thermal structures elongate into streaks and evidence the transitional regime. Further increasing the wall shear causes the streaks to meander in the spanwise direction which indicates the shear dominated regime (Re w = 3000, 4000).
As P r = ν/κ exceeds unity, kinematic viscosity overtakes thermal diffusivity and the wall shear affects the flow structures in the bulk earlier. In fact, it can be observed that already for Re w = 1000 the flow shows the meandering behaviour of the shear dominated regime. For P r = 4.6 and Re w = 4000 the shear is strong enough to make the effect of the thermal forcing negligible, as confirmed by the flow structures similar to the plane Couette flow.
Conversely, for Prandtl numbers smaller than unity, the shear is less effective for a given Re w and the bulk flow is more dominated by the thermal structures. In the case of P r = 0.22, a wall shear of Re w = 1000 is not strong enough to fully disturb the plumes and only the next data point at Re w = 2000 shows signs of elongated streaks.
From the panels of figure 2 it is evident how P r changes the relative strength of momentum and thermal diffusivities: A higher Prandtl number increases the momentum transfer from the boundaries to the bulk and the transition to the shear dominated regime occurs at a lower Re w than for a corresponding low P r flow. Vice versa, for small Prandtl numbers, the thermal dominated regime is more persistent and the shear dominated flow features appear only at high Re w .

Heat transfer
The Nusselt number N u is plotted in figure 3 as function of Re w , showing a nonmonotonic behaviour. The common feature is that for increasing wall shear, N u first decreases and then increases as already observed in Blass et al. (2020) for P r = 1. In the present case, however, the specific values are strongly dependent on P r, as seen in figure 3c. The effect of P r is strongly dependent on the amount of shear added to  the system. For pure Rayleigh-Bénard convection (Re w = 0), N u increases with P r for P r < 1 and saturates to a constant value for 1 < P r < 4.6, see figure 3b, in agreement with the findings of van der Poel et al. (2013) and Stevens et al. (2013). For increasing Re w , Prandtl number effects on the heat transfer are more pronounced, because of the higher momentum transfer from the boundaries to the bulk. This is confirmed both by the initial N u decrease up to 20% of the RB value at P r = 4.6 and the subsequent strong increase by more than 50% for the highest Re w . In both cases the effects of the momentum transfer are enhanced by the high Prandtl number. We wish to stress that the non-monotonic behaviour of the Nusselt number observed here is a common feature of flows in which more than one parameter concur to determine the value of the heat transfer; for example similar dynamics are reported by Yang et al. (2020) and Wang et al. (2020) for thermal convection with rotation or Chong & Xia (2016) for severe lateral confinement.

Flow layering
The initial N u decrease can be understood upon considering that the added wall shear perturbs the thermal RB structures and produces a horizontal flow layering that weakens the vertical heat flux. Once the wall shear is strong enough, however, the flow undergoes a transition to a shear dominated regime and the vertical cross-stream motion generated by the elongated streaks makes up for the suppressed RB structures, thus starting the Nusselt number monotonic increase (Blass et al. 2020). To better understand the effect of the horizontal flow layering, we discuss the results of figure 4. In these 'side views' (i.e., streamwise cross-sections) of the temperature field snapshots and the corresponding top views of figure 2, we can observe how the flow changes from thermal plumes to straight thin streaks and then to meandering structures. As expected, the increase in wall shear causes the flow to become more turbulent. But the change in the large-scale structures is also very recognizable. Here, the transitional regime displays a more unexpected behaviour. In contrast to what is seen in panels 4(a,c), where the flow structures appear clearly divided into hot and cold columns, in panel 4b the structures are more complex. Due to the wall shear and the thereby imposed horizontal flow, the vertical structures are disturbed, the flow is not able to reach the opposite hot/cold wall, but is instead trapped in a warm/cool state in the bulk of the flow. The fluctuations in the flow are not strong enough to mix the bulk and therefore the heat gets insulated in a stably stratified layer in the middle of the flow. This layering causes the total heat transfer to decrease and is the reason for the drop in N u for low Re w in figure 3. Because of the heat entrapment in the bulk layer, relatively cold fluid comes very close to relatively warm fluid and the temperature gradients in wall-normal direction increase significantly. In the atmosphere, this phenomenon can be observed as cloud streets, which, similar to the here observed high-shear end of the transitional regime, manifests as long streaks of convection rolls (Etling & Brown 1993;Kim et al. 2003;Jayaraman & Brasseur 2018).

Boundary layer thicknesses
A complementary way to better understand the P r-dependence of the flow dynamics and the transport properties is to study the viscous and thermal boundary layer thicknesses λ u and λ θ , respectively. Here, we define both λ θ and λ u by extrapolating the linear slopes of the mean temperature and mean streamwise velocity close to the walls, similarly to Shishkina et al. (2010). The dependence of λ u and λ θ on Ri and P r is shown in figure  5. Here we use as abscissa the Richardson number. Given that Ra = 10 6 is constant, we have Ri ∝ (P rRe 2 w ) −1 . At every P r, for increasing Ri -and therefore decreasing shearλ θ initially grows, then reaches a plateau around Ri ≈ 1 and eventually decreases slowly to converge to the pure RB value (figure 5a). For comparison, we also plot N u(Ri) in figure 5c. Given that λ θ ∝ (N u) −1 to a good approximation, the behaviour of the thermal boundary layer thickness is consistent with the Nusselt number of figures 3 and 5c. The different flow regimes can be identified either from the different slopes of λ θ versus Ri or from those of N u(Ri). The slope is positive in the shear dominated region (small Ri), approximately zero in the transitional regime and then negative in the thermal buoyancy dominated regime.
As the Richardson number indicates the relative strength of buoyancy and shear, one might think that for increasing Ri there should be a monotonic λ θ decrease which, instead, is observed only for Ri 1. The reason for the counter-intuitive λ θ increase for Ri 1 is that in this region the thermal forcing is weak and the flow is mainly driven by the shear. In this case the thermal boundary layer is slaved to the viscous boundary layer which, according to the expectations, monotonically thickens as the wall shear weakens. From figure 5b we can see that indeed λ u monotonously increases with increasing Ri. Note that the viscous boundary layer thickness has a stronger dependence on P r than the thermal boundary layer thickness. Qualitatively, larger P r reflects stronger momentum diffusivity and therefore a thicker viscous boundary layer. In the shear dominated regime (high P r or low Ri), however, λ u grows faster than in the other regimes and this is especially true for the flows with higher P r. In fact, in these cases the thermal boundary layer is nested within the viscous one and the dynamics of the latter is not influenced by the former. This is not the case for small P r < 1 because then λ u evolves inside λ θ whose thinning with increasing Ri counteracts the thickening of the viscous boundary layer.
To further stress the importance of the relative thicknesses of the thermal and the viscous boundary layer, we show their ratio versus Ri in figure 5d. We can see that λ θ /λ u increases for decreasing P r at fixed Ri since the kinetic boundary layer thickness is driven by the momentum diffusivity. At fixed P r the behaviour of the boundary layer ratio is more complex: it always shows a decreasing trend in the high end of Ri which is due to the thinning of the thermal boundary layer. On the other hand, at the low end of Ri one can observe an increase only for P r > 1, which is due to the steep growth of λ θ with Ri observed in figure 5a.
Due to the limited amount of datapoints, we cannot show a more detailed behaviour in the extreme case of pure shear forcing. In contrast, in the limit of pure Rayleigh-Bénard convection we do observe the asymptotic trend for λ θ /λ u ; there the effect of the shear becomes very small (no imposed shear, all shear due to natural convection roll) and the ratio depends on P r only. This saturation occurs earlier for smaller P r, because the thermal forcing dominates over the shear forcing at smaller Ra.

Velocity and temperature wall profiles
For strong enough shear the boundary layers, which are first of laminar type, will eventually become turbulent, considerably enhancing the heat transport. However, for most of the values of the control parameters (Re w and P r) of this paper this is hardly the case. This can best be judged from the velocity profiles, which we show in figures 6(a-c) for three different values of P r and various Re w . Only in the high-P r range, towards the limit of plane Couette flow, we can see that u + evolves towards the well-known Prandtlvon Karman logarithmic behaviour u + (z + ) = κ −1 log z + + B for high Re w . Since the shear strongly affects the flow, the boundary layers can transition to turbulence earlier than without shear. But also the large P r enhance the shear. In fact, at P r = 4.6 already the flow at Re w = 3000 shows the onset of a log-law behaviour. The more P r is decreased, the harder it becomes for the wall shear to disturb the thermal plumes and, as a result, at Re w 4000 and P r 2.2, the log-scaling cannot be attained in our simulations.
Panels 6(d-f) show a similar behaviour for the mean temperature profiles. One can observe that the temperature profiles converge earlier towards some type of logarithmic behaviour. For P r = 1, we can see such behaviour for Re w = 4000, whereas at larger P r = 4.6, it already shows up even at Re w = 2000. From the shown temperature profiles, we can also identify the flow layering that was previously discussed in section 3.3. When the flow layering occurs, heat gets entrapped in the bulk flow. Since now an additional layer of warm and cool fluid exists in between of the cold and hot regions, T + shows a non-monotonic behaviour with a drop after the initial peak. This can most prominently be seen in figure 6d (P r = 0.22) for the strongest shear Re w = 4000.

Conclusion
In this manuscript we performed DNS of wall sheared thermal convection with 0 Re w 4000 and 0.22 P r 4.6 at constant Rayleigh number Ra = 10 6 . Similarly to Blass et al. (2020), who analysed the Ra-dependence of sall sheared thermal convection, we found three flow regimes and quantified them by using the bulk Richardson number and a visual analysis of two-dimensional cross-sectional snapshots. The flow undergoes a transition from the thermal buoyancy dominated to the transitional state when Ri 10. We found that the meandering streaks of the shear dominated regime start to emerge at Ri 0.1. Also the behaviour of the Nusselt number strongly depends on P r. For high Prandtl number, the momentum transfer from the walls to the flow is increased and therefore the flow can easier reach the shear dominated regime where the heat transfer is again increased. We analysed both the thermal and the kinetic boundary layer thicknesses to better understand the transitions of the flow between its different regimes. We found that the thermal boundary layer thickness λ θ shows a peak in the transitional regime and decreases for both lower and higher Ri. The kinetic boundary layer thickness λ u increases with increasing Ri and increasing P r. For very strong Re w and in particular large P r we notice the appearance of logarithmic boundary layer profiles, signalling the onset of turbulent boundary layer dynamics, leading to an enhanced heat transport.
Together with the results of Blass et al. (2020), we now have analysed two orthogonal cross-sections of the three-dimensional parameter space (Ra, P r, Re w ). More specifically, we have determined N u(Ra, P r, Re w ) for the two cross-sections N u(Ra, P r = 1, Re w ) in Blass et al. (2020) and N u(Ra = 10 6 , P r, Re w ) here. From standard RB without shear we of course know N u(Ra, P r, Re w = 0), which is perfectly described by the unifying theory of thermal convection by Grossmann & Lohse (2000, 2001 and Stevens et al. (2013). The knowledge of the two new cross-sections in parameter space may enable us to extend this unifying theory to sheared convection.