1. Introduction
Natural convection boundary layer (NCBL) flow is ubiquitous in a vast variety of industrial and geophysical applications where the temperature difference between a vertical surface and the ambient fluid induces a spontaneous fluid motion. The majority of past efforts have been dedicated to investigating such convective flows with small temperature differences where the fluid property variations are considered negligible and the buoyancy force is modelled by a linear variation with local temperature difference (e.g. Versteegh & Nieuwstadt Reference Versteegh and Nieuwstadt1999; Abedin, Tsuji & Hattori Reference Abedin, Tsuji and Hattori2009; Nakao, Hattori & Suto Reference Nakao, Hattori and Suto2017). However, many natural convection flows operate at large temperature differences – for instance, the temperature difference between the surface of a solar thermal receiver and the ambient air usually exceeds $650\,{\rm K}$ (Ho & Iverson Reference Ho and Iverson2014), in which cases, the local fluid density and properties would vary significantly. In an early study, Siebers, Moffatt & Schwind (Reference Siebers, Moffatt and Schwind1985) experimentally investigated the heat transfer characteristics of the NCBL along a vertical isothermal plate up to a temperature difference of $500\,{\rm K}$ for air. Based on their measurements, they found that the critical Grashof number which marks the onset of the laminar–turbulent transition, when evaluated at the film temperature $T_f = ( T_w+ T_\infty )/2$, decreases with increasing temperature ratio $T_w/T_\infty$ (where $T_w$ is the isothermal wall temperature and $T_\infty$ is the ambient temperature). Such an observation implies that the stability results obtained from the Oberbeck–Boussinesq (OB) limit may not be directly applied to the NCBL with large temperature differences using the film temperature.
While the linear stability for NCBLs under OB conditions has been studied extensively over decades (see, e.g. Xin & Le Quéré Reference Xin and Le Quéré2001; Tumin Reference Tumin2003; Xin & Le Quéré Reference Xin and Le Quéré2012), the non-Oberbeck–Boussinesq (NOB) effects on the stability for NCBL flows have received less attention. Carey & Mollendorf (Reference Carey and Mollendorf1978, Reference Carey and Mollendorf1980) analytically demonstrated the applicability of similarity analysis to NCBLs in liquids wherein the viscosity is assumed to vary linearly with temperature while treating all other properties as constants. Although their first-order viscosity variation approximation is only valid when the temperature difference is small (the viscosity variation itself being small enough), their results suggest that the temperature-dependent viscosity would have significant effects on the stability and transition of the NCBL flows by adjusting the mean velocity and temperature profiles. Carey & Mollendorf (Reference Carey and Mollendorf1978) further pointed out that the flow similarity exists only when the temperature distribution of the flow is fixed (i.e. a self-similar temperature profile). Sabhapathy & Cheng (Reference Sabhapathy and Cheng1986) numerically investigated the effects of temperature-dependent viscosity and coefficient of thermal expansion on the linear stability properties of NCBLs along a vertical isothermal wall for liquids with Prandtl numbers $7\leq {{Pr}}\leq 10$. In their study, the variation of viscosity is approximated by a linearised Taylor series expansion about the film temperature. Using Orr–Sommerfeld eigenvalue formulation, they found that the variable viscosity stabilises the liquid flow near heated walls and destabilises the flow adjacent to cooled walls, which appears consistent with the findings of Piau (Reference Piau1974); while the variable coefficient of thermal expansion would stabilise the flow near a heated wall locally but destabilise the flow farther downstream. A similar approach was also adopted by Jiin-Yuh & Mollendorf (Reference Jiin-Yuh and Mollendorf1988) who numerically studied the linear stability of a vertical plate NCBL in ethylene glycol (${{Pr}}=100$ at film temperature) with a temperature difference of $23\,{\rm K}$. In their analysis, the linear variation of viscosity with temperature is incorporated into the Orr–Sommerfeld equation based on self-similar velocity and temperature base flow profiles. Their numerical results confirmed the findings of Sabhapathy & Cheng (Reference Sabhapathy and Cheng1986), and showed that the glycol flow is more stable when heated than the cooled case with same temperature difference and film temperature. Chenoweth & Paolucci (Reference Chenoweth and Paolucci1985) investigated the vertical NCBL in differentially heated enclosed slots for gases. Using Sutherland's law, they analytically obtained the steady-state exact solution to the laminar flow up to $\epsilon =1$ with both viscosity and thermal conductivity considered temperature-dependent, where $\epsilon = 2(T_h-T_c)/(T_h+T_c)$ is a dimensionless temperature difference given by the differentially heated walls $T_h$ and $T_c$. In their study, they showed that although the variable property only has minimal effect on the wall heat transfer, the temperature and velocity distributions are sensitive to the property variations. Their study was extended by Chenoweth & Paolucci (Reference Chenoweth and Paolucci1986) who numerically investigated the laminar vertical NCBL in cavities. Using a weakly compressible formulation, Chenoweth & Paolucci (Reference Chenoweth and Paolucci1986) showed that flow instability may follow different mechanisms with increasing temperature difference. Suslov & Paolucci (Reference Suslov and Paolucci1995b) further investigated NOB effects on the linear stability. In their analysis, the authors demonstrated that the cavity NCBL becomes unstable to a shear instability at relatively low temperature difference and the critical Rayleigh number which marks the onset of flow instability increases with increasing temperature difference (i.e. stabilising the flow) until a critical temperature difference $\epsilon ^*$. For greater temperature differences $\epsilon >\epsilon ^*$, the flow instability is switched to a buoyancy-driven one, resulting in a sudden reduction in the critical Rayleigh number (i.e. destabilising the flow) and the most unstable wavenumber. A similar behaviour is also seen by Suslov & Paolucci (Reference Suslov and Paolucci1995a) for mixed convection flows in tall channel. More recently, Rajamanickam, Coenen & Sánchez (Reference Rajamanickam, Coenen and Sánchez2017) studied the NOB effects on the linear stability of an inclined NCBL using a simplified low-Mach-number formulation and the fluid properties are approximated using a power-law. Their numerical results show that the NOB air flow can have travelling waves over a significantly wider inclination angle than those predicted by the Boussinesq limit. Liu et al. (Reference Liu, Xia, Yan, Wan and Sun2018) investigated the linear stability in Rayleigh–Bénard convection, where the NOB effects are found to stabilise (destabilise) the flow for large (small) aspect ratios. For natural convection in cavities (including Rayleigh–Bénard systems) and differentially heated channels, the strong influences of the NOB effects on the flow instability are found mainly due to the breakdown of symmetry when heated and cooled (e.g. Ahlers Reference Ahlers1980; Suslov & Paolucci Reference Suslov and Paolucci1995a,Reference Suslov and Paoluccib; Liu et al. Reference Liu, Xia, Yan, Wan and Sun2018); whereas such a symmetry does not exist for NCBL flows along a vertical plate. Despite the existing efforts, a fundamental understanding of the NOB effects, induced by the bulk temperature difference, on the onset of such flow instability remains unclear.
In this paper, the linear stability of a two-dimensional laminar, periodic vertical NCBL in air is investigated by employing the linearised disturbance equations as an initial value problem. The periodic parallel flow is practically important to those start-up transient flows with large downstream distances where the one-dimensional conductive regime could exist for a significant amount of time before the arrival of the leading edge signal. The existence of such a one-dimensional conductive flow regime has been confirmed both experimentally and numerically (e.g. Siegel Reference Siegel1958; Miyamoto Reference Miyamoto1977; Gebhart Reference Gebhart1988; Sammakia, Gebhart & Qureshi Reference Sammakia, Gebhart and Qureshi1982) with analytical solutions made available for the OB case (Illingworth Reference Illingworth1950; Schetz & Eichhorn Reference Schetz and Eichhorn1962). In particular, Goldstein & Briggs (Reference Goldstein and Briggs1964) and Joshi & Gebhart (Reference Joshi and Gebhart1987) showed that the flow could become unstable and therefore transition to turbulence before the arrival of the leading edge signal. The linear stability of such parallel transient flow has been investigated under OB conditions for NCBL in cavities (Brooker et al. Reference Brooker, Patterson, Graham and Schöpf2000) and along a vertical flat plate (Joshi & Gebhart Reference Joshi and Gebhart1987; Krane & Gebhart Reference Krane and Gebhart1993; Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019). However, with increasing characteristic temperature difference, the stability results obtained from the OB studies become largely unreliable in predicting the flow instability and transition in real-world applications. Using a weakly compressible formulation (Batchelor Reference Batchelor1953; Paolucci Reference Paolucci1982), for the first time, we are able to obtain the marginal stability curve of such a periodic NCBL flow with large temperature differences (up to $100\,{\rm K}$); and provide a detailed insight of the NOB effects on its instability mechanisms.
The rest of paper is organised as follows. An overview and the mathematical formulation for this initial-value problem is given in § 2. Based on the linear stability results, the marginal stability curves for various temperature differences are obtained in § 3.1. The instability dynamics are further explored in § 3.2 by inspecting the energy budgets for two representative marginally unstable modes. The linear stability results are re-normalised using the film temperature in § 3.3 to allow comparisons with observations in the literature. In § 3.4, a direct numerical simulation is carried out to incorporate the unsteady and nonlinear effects – which is then compared with the linear results in § 3.1 to determine the linear stability regime. Section 4 briefly summarises the findings in this paper.
2. Numerical formulation
The problem under consideration is a vertical natural convection flow adjacent to an isothermally heated plate in air. With the use of Squire's theorem (Squire Reference Squire1933), the transient stability problem is reduced to a two-dimensional one. In the analysis that follow, equations and results are made dimensionless using the intrinsic length scale $L_{s,\infty }={\xi ^*}^{2/3}/(g^*\beta ^*\Delta T)^{1/3}$ and velocity scale $U_{s,\infty }=(\xi ^*g^*\beta ^*\Delta T)^{1/3}$ since the temporally developing parallel flow has no nature length scale. Here, $\xi$ represents the thermal diffusivity, $g$ is the gravitational acceleration, $\Delta T$ is the temperature difference given by the heated wall $T_w$ and the ambient $T_\infty$, and $^*$ indicates the dimensional quantities that are evaluated at the ambient temperature $T_\infty$. Fluid density, coefficient of thermal expansion and specific heats are normalised using the ambient temperature quantities so that $\rho^* = \rho \rho_\infty^*$, $\beta^* = \beta \beta_\infty^*$ and $C_p^* = C_p C_{p\infty}^*$.
The NCBL flow is governed by the conservation equations with low-Mach-number approximation in which the sound waves are filtered while allowing arbitrary density variations (Batchelor Reference Batchelor1953; Paolucci Reference Paolucci1982). The dimenionless form of these governing equations read
where the subscript $i,j=1,2$ denotes the streamwise ($x$) and wall-normal ($y$) directions, and $u_1, u_2 = u, v$ are the streamwise and wall-normal velocity, respectively; $\kappa$ is the thermal conductivity, $\mu$ is the dynamic viscosity, $p$ represents the hydrodynamic pressure; $n_i$ is the unit gravitational vector pointing in the $-x$ direction. The viscous stress tensor $\tau _{ij}$ is given by
and $\delta _{ij}$ is the Kronecker delta. The thermal field is made dimensionless using $\theta = (T-T_\infty ) /\Delta T$, whereas the bulk temperature difference is given by $\epsilon =\Delta T/T_\infty$; $\varpi$ represents the spatially uniform thermostatic pressure, and $\varGamma = (\gamma _\infty -1)/\gamma _\infty$ where $\gamma _\infty = C_p/C_v$ is the ratio of the specific heats under constant pressure and volume at $T_\infty$. In the present study, both dynamic viscosity and thermal conductivity of the air flow are assumed to vary with the local temperature $\theta$ following Sutherland's law (Sutherland Reference Sutherland1893), which in dimensionless form reads
where ${{Pr}}_\infty \equiv C_p\mu _\infty /\kappa _\infty$ is the Prandtl number evaluated at $T_\infty$, and $S_\mu =111\,{\rm K}$ and $S_\kappa =194\,{\rm K}$ are Sutherland constants for air (Hilsenrath Reference Hilsenrath1955; White Reference White1988). The use of Sutherland's law was shown to be able to accurately predict both viscosity and conductivity for air up to $\epsilon =0.6$, beyond which its accuracy degrades rapidly (Chenoweth & Paolucci Reference Chenoweth and Paolucci1985). For $T_\infty = 288\,{\rm K}$, this limit corresponds to a maximum $\Delta T = 172\,{\rm K}$. All other properties are considered constant and evaluated at $T_\infty$ (Bergman et al. Reference Bergman, Lavine, Incropera and DeWitt2011, e.g. the specific heat $C_p$ varies by ${\sim }1.5\,\%$ for a temperature difference of $\Delta T =172\,{\rm K}$ from $T_\infty =288\,{\rm K}$) – similar assumptions have also been invoked in early studies (e.g. Chenoweth & Paolucci Reference Chenoweth and Paolucci1985, Reference Chenoweth and Paolucci1986; Emery & Lee Reference Emery and Lee1999, for NCBL in vertical slots and cavities).
Since the laminar NCBL flow under consideration is perfectly parallel to the isothermal wall and homogeneous in the streamwise direction (see figure 1), the wall-normal velocity and streamwise gradients must vanish $v=\partial /\partial x=0$, with which system (2.1) is simplified to (cf. Suslov & Paolucci Reference Suslov and Paolucci1995a,Reference Suslov and Paoluccib)
where $\mathcal {V}$ is the computational volume; $U_b$, $\theta _b$ and $\rho _b$ are the streamwise velocity, temperature and density profiles of the laminar flow. For the purpose of linear stability analysis, an artificial temperature perturbation $\tilde {\theta }$ is superimposed on the laminar base flow ($U_b, \theta _b$) at a designated Grashof number $Gr_\infty$ – at which the base flow is ‘frozen’ (treated as quasi-steady) since the temporal evolution of the base flow is considered much slower than that of the perturbations (Otto Reference Otto1993; Brooker et al. Reference Brooker, Patterson, Graham and Schöpf2000; Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019). The validity of quasi-steady treatment for the transient one-dimensional periodic NCBL flows has been extensively discussed in the limit of both short and long waves under OB conditions (Daniels & Patterson Reference Daniels and Patterson1997, Reference Daniels and Patterson2001), where they showed that the quasi-steady approximation is most accurate for short waves, and the accuracy for the long-wave disturbances will increase with the amplification rate (or, with increase Reynolds number, Tromans Reference Tromans1978; Dwoyer & Hussaini Reference Dwoyer and Hussaini1987). The temporal evolution of the perturbations is then examined in a similar sense to Ke et al. (Reference Ke, Williamson, Armfield, McBain and Norris2019) as an initial value problem.
Here, the ambient Grashof number $Gr_\infty$ and the integral velocity boundary layer thickness is defined by
where $U_m$ is the maximum streamwise velocity of the laminar flow. For a parallel periodic NCBL, both $\delta$ and $Gr_\infty$ are functions of time $t$ only as the flow is invariant in the streamwise ($x$) direction but unsteady in time. Notably, the exact solutions to the boundary layer in (2.4) are analytically difficult to obtain since the flow fields lose their self-similarity due to the temperature-dependent thermal conductivity (Carey & Mollendorf Reference Carey and Mollendorf1978). Although the laminar streamwise velocity profile of the NCBL under different temperature differences do not appear self-similar, the use of the ambient Grashof number (2.5) allows a direct comparison of these flows with the same dimensionless boundary layer thickness $\delta$. In this study, the laminar base flow profiles are numerically obtained using a precursor simulation which solves (2.4) simultaneously with the boundary conditions $U_b=0$, $\theta _b=1$ at $y=0$ and $\partial U_b/\partial y=\theta _b=0$ at $y=\infty$ (see Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019, Reference Ke, Williamson, Armfield, Norris and Komiya2020, for numerical details of the precursor DNS), allowing the flow to develop from quiescence until a designated $Gr_\infty$ is reached. We also note that for larger temperature difference, a shorter development time $t$ is required to reach a given $Gr_\infty$ than that of the smaller $\epsilon$ case, resulting in a lower maximum streamwise velocity $U_m$. For the cases considered in the present study, the maximum streamwise velocities at $Gr_\infty =600$ are $U_m=3.67, 3.62,3.38, 3.14$ for $\epsilon =0$ (Oberbeck–Boussinesq case), $0.035, 0.174, 0.347$. The streamwise velocity magnitude difference is not as obvious as for the NCBL flows in cavities and vertical channels: for example, in an open cavity, Juárez et al. (Reference Juárez, Hinojosa, Xamán and Tello2011) showed the maximum velocity is decreased by ${\sim }30\,\%$ and the boundary layer thickness $\delta$ is increased by $\sim 70\,\%$ for a $0.3$ increase in the temperature difference $\epsilon$ at a given cavity height Rayleigh number ($Ra_H$) in the steady state; in a fully enclosed cavity, Zhong, Yang & Lloyd (Reference Zhong, Yang and Lloyd1985) showed that there is a ${\sim }20\,\%$ decrease in the $U_m$ for $\epsilon =0.1$ (although the fluid density is considered constant in their study) when compared with the OB flows. Notably, the steady-state boundary layer thickness decreases for the NCBL in cavities at larger $Ra_H$ since the geometric length scale $H$ and cavity aspect ratio limit the development time to reach steady state for different $\epsilon$ values (Juárez et al. Reference Juárez, Hinojosa, Xamán and Tello2011); whereas in our temporally developing flow, the boundary layer thickness increases with the development time. Typical laminar base flow profiles for our periodic NCBL are shown in figure 2 at $Gr_\infty =600$, where the OB results are shown to agree well with the low-temperature-difference case ($\Delta T = 10\,{\rm K}$, $\epsilon =0.035$) indicating the NOB effects on the base flow profiles are negligible at this temperature difference. With increasing $\epsilon$ (equivalently, $\Delta T$ in dimensional arguments), the laminar base flow profiles deviate from the OB approximation: as depicted by figure 2(a), with increasing $\epsilon$, the velocity profile shifts away from the wall (see inset). This trend is consistent with the observations of Zhong et al. (Reference Zhong, Yang and Lloyd1985) for the natural convection in a square cavity that the variable-fluid-property effect would reduce the velocity in the near-wall region but increase it in the outer region. Temperature profiles, shown in figure 2(b), exhibit a slight difference from the OB case, where the local temperature is increased with increasing $\epsilon$ in the boundary layer (see inset).
The temporal response of the two-dimensional (as per Squire's theorem) artificial perturbations, however, is characterised by the linearised disturbance equations, which can be obtained by substituting $\phi = \phi _b+\tilde {\phi }$ (where $\phi = u,v,\theta,p, \rho, \mu, \kappa$) into (2.1), dropping out the higher-order nonlinear terms and neglecting the base flow variation in time
where the Laplacian operator $\triangledown ^2=\partial ^2/\partial x^2+\partial ^2/\partial y^2$. The temperature perturbation, which takes the form
is fed into the base flow at time $t=0$ and interacts linearly with the base flow following (2.6). Here, $n$ is the mode number – a real positive integer indicating the wavelength $L_x/n$ (the dimensionless wavenumber evaluated at ambient temperature $T_\infty$ is therefore $k_\infty =2{\rm \pi} n /L_x$); $A_n=10^{-3}$ is the initial magnitude of the perturbation and $L_x=L_x^*/L_{s,\infty }$ the dimensionless streamwise domain length. Perturbations must vanish at the non-slip wall and decay towards zero in the far-field, giving the perturbation boundary conditions $\tilde {u}=\tilde {v}=\tilde {\theta } = 0$ at both wall ($y=0$) and far-field ($y=\infty$). No velocity perturbation is added at $t=0$ in the current study as the temperature perturbation will directly disturb the velocity field via the buoyancy term in (2.6b). The linearised system (2.6) is solved numerically for a given wavenumber $k_\infty$ (or mode number $n$) with a base flow at $Gr_\infty$ and temperature difference $\epsilon$ as an initial-value problem. In the present study, the NOB effects on the stability properties are explored for $400\leq Gr _\infty \leq 2000$, $4\leq n\leq 32$ ($0.02 \leq k_\infty \leq 0.19$) and $10K\leq \Delta T\leq 100\,{\rm K}$ (equivalently, $0.035\leq \epsilon \leq 0.347$). The ambient temperature is set to $T_\infty =288\,{\rm K}$ with an ambient Prandtl number of ${{Pr}}_\infty = 0.71$ (based on Bergman et al. Reference Bergman, Lavine, Incropera and DeWitt2011) for all cases considered. The calculations are carried out in a computational domain of size $L_x^* \times L_y^* = 1035L_{s,\infty } \times 210L_{s,\infty }$, using a fractional step method (Xia et al. Reference Xia, Wan, Liu, Wang and Sun2016; Demou, Frantzis & Grigoriadis Reference Demou, Frantzis and Grigoriadis2018) with a collocated finite volume grid ($N_x\times N_y =1024 \times 256$). The domain is made sufficiently large so that both gradients and magnitudes of the flow variables are reduced to their ambient values at the far-field $y=L_y$. A uniform mesh is employed in the streamwise direction while the grids in the wall-normal direction are stretched using a Gamma function with a maximum stretching rate of $1.24\,\%$.
3. Results and discussion
3.1. Linearised disturbance equations
With the laminar base flow at a given ambient Grashof number $Gr_\infty$, the flow is unstable (or stable) to the perturbation $k$ if the magnitude of its temporal response grows (or decreases) with time. Notably, due to the parallel nature of the unsteady periodic flow, perturbations have no streamwise dependence as would have been seen in the spatially developing flows (Huerre & Monkewitz Reference Huerre and Monkewitz1990). As a result, the perturbations would become unstable (or stable) simultaneously at the same amplification (or damping) rate in the flow anywhere, in which sense the convective instability (or stability) of the flow is global rather than local (see, e.g. later in figures 6 and 8). In this study, the temporal response to the initial perturbation $k$ is tracked by the time trace of the velocity and temperature signals ($\tilde {u}$ and $\tilde {\theta }$) at an arbitrary location $(x_p,y_p)$ within the boundary layer (see also figure 6 of Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019). Figure 3 shows the typical temporal responses of the perturbation signals for $k_\infty =0.0486$ ($n=8$) and $k_\infty =0.0971$ ($n=16$) at $Gr_\infty =600$ and $\epsilon = 0.174$. As seen in figure 3, at $Gr_\infty =600$, the oscillatory temporal signals quickly adjust their magnitude to a fully periodic state within $1\sim 2$ temporal periods (i.e. receptivity as the flow receives the external perturbation at $t=0$). In this fully periodic state, the envelope of the temporal signals follow an exponential growth (for $k_\infty =0.0486$) or decay (for $k_\infty =0.0971$) with time. The growth rate of the perturbation signals, termed the amplification rate, is then given by
where $A_t$ is the instantaneous magnitude of the perturbation signal at time $t$ (i.e. the envelope magnitude in figure 3), which is obtained by taking the fast Fourier transform in the streamwise ($x$) direction at $y=y_p$. For a given combination of $k_\infty$ and $Gr_\infty$, a positive $\sigma _k$ indicates that the temporal response is amplified by the laminar base flow while a negative $\sigma _k$ represents a damped signal. The neutral curve, which characterises the minimal Grashof numbers at which the flow is linearly unstable to a given perturbation $k_\infty$ (i.e. $\sigma _k=0$), can be traced by bisectioning the amplification rates in the $k$–$Gr$ plane.
Figure 4 compares the neutral curves ($\sigma _k=0$) for $\epsilon =0.035,0.174$ and $0.347$ with the OB case (Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019). While the neutral curves of all the NOB cases show apparent qualitative similarity to the OB case, these curves demonstrate a clear $\epsilon$-dependence: with increasing $\epsilon$, the neutral curve shifts towards higher Grashof numbers at both lower branch ($k_\infty < O(0.1)$) and upper branch ($k_\infty >O(0.1)$). As seen in figure 4, although the critical mode $k_\infty =0.049$ in the lower branch has minimal $\epsilon$-dependency, the critical Grashof number is increased to $Gr_\infty = 473, 536$ and $622$ for $\epsilon = 0.035,0.174$ and $\epsilon =0.347$, respectively, when compared with $Gr_{\infty }=454$ for the OB case. At larger wavenumbers in the upper branch ($k_\infty >O(0.1)$), such a dependency on the temperature difference $\epsilon$ is considerably less than that of the lower branch, where the marginally unstable Grashof number of $k_\infty =0.12$ is increased to $Gr_\infty =948, 959$ and $980$ for $\epsilon = 0.035,0.174$ and $\epsilon =0.347$, respectively. A similar $\epsilon$ dependency of the critical Grashof number has also been reported by Suslov & Paolucci (Reference Suslov and Paolucci1995b) for an NCBL in a differentially heated cavity when using a reference temperature obtained by averaging the hot and cold walls. Based on their eigenvalue analysis, Suslov & Paolucci (Reference Suslov and Paolucci1995b) found that the critical point (i.e. critical $Gr$ and $k$) for the NCBL in a cavity initially locates in the upper branch for relatively low $\epsilon$ and the critical Grashof number would increase with increasing $\epsilon$ – until $\epsilon$ is sufficiently large ($\epsilon >0.54$) such that the critical wavenumber would then shift to the lower branch. However, in the present study, such a drastic change in the critical wavenumber is not observed – we will show later in § 3.2 that this is due to a different instability mechanism in the cavity.
Since the perturbation wavelength $\lambda$ is fully determined by the mode number $n$ and the streamwise domain length $L_x$, the temporal response signal in figure 3 allows further examination of the wave speed of the perturbation travelling wave $v_p$,
where $\omega$ is the oscillation frequency and $\Delta \tau$ is the dimensionless time interval between two neighbouring peaks of the response signal, as shown in figure 3(b). Based on linear theory, the interaction between the perturbation and its base flow is linear so that $\Delta \tau$ is constant for a given choice of $Gr_\infty$, $k_\infty$ and $\epsilon$. Figure 5(a,b) compare the maximum base flow velocity $U_m$ with the wave velocity obtained by (3.2) at $Gr_\infty =600$ and $Gr_{\infty }=1000$, respectively. For all temperature differences and base flows considered, the wave velocity decreases with increasing wavenumber and temperature difference $\epsilon$; and the NOB effects on the phase velocity appear more apparent for higher wavenumbers. Notably, Daniels & Patterson (Reference Daniels and Patterson1997, Reference Daniels and Patterson2001) showed that for vertical NCBL flows, the ‘less accurate’ long-wave perturbations will always have a phase speed greater than the maximum base flow velocity ($v_p/U_m\approx 2$); whereas the phase velocity for the short-waves are always less than the maximum base flow velocity ($v_p< U_m$). In this sense, perturbations considered in our study are mostly short-wave disturbances, with a small portion of the ‘less accurate’ intermediate–long waves (since $v_p/U_m\approx 1.2 < 2$ in figure 5a,b) whose accuracy is improved with increasing Grashof (Reynolds) number. At $Gr_\infty =600$, the critical mode $k_\infty =0.049$ (cf. figure 4) is found to have a wave speed close to the maximum base flow velocity $U_m$. This observation is consistent with the findings of Armfield & Patterson (Reference Armfield and Patterson1992) in which the authors showed that the travelling waves amplified by the OB base flow have a wave speed close to the maximum base flow velocity; and with the eigenvalue analysis for the OB case (Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019) from which the phase velocity of the critical point is found to be equal to $U_m$. Such a result is seen more clearly in figure 5(c), where the amplification rate $\sigma _k$ is plotted against the wave speed $v_p$. From figure 5(c), only a small band of wavenumbers are amplified by the base flow ($\sigma _k>0$) at $Gr_\infty =600$ – all of which are shown to have a wave speed close to the velocity maximum in the base flow $U_m$, regardless of the temperature differences. Figure 5(c) also shows apparent NOB effects on the amplification rates: the base flow is unstable to $k_\infty =0.049$ ($v_p/U_m=0.98$) and $k_\infty =0.061$ ($v_p/U_m=0.89$) for $\varepsilon \leq 0.174$, but appears stable for $\varepsilon =0.347$ (see inset). At $Gr_{\infty }=1000$, the base flow becomes unstable to a broader band of wavenumbers at the lower phase velocity range at $v_p/U_m\approx 0.55$ (i.e. shorter waves with higher wavenumbers). As shown in figure 5(d), while the critical mode $k_\infty =0.049$ remains the most amplified mode at $Gr_{\infty }=1000$, its phase velocity is decreased to $v_p/U_m=0.88$.
3.2. Energy budgets for the buoyancy-driven and shear-driven modes
In figure 4, the lower wavenumbers demonstrate a greater sensitivity to the temperature difference $\epsilon$ compared with the upper branch (i.e. larger wavenumbers $k_\infty >O(0.1)$), indicating potential differences in their instability dynamics. In this subsection, we examine the instantaneous perturbation field as well as the perturbation kinetic energy budget of the representative modes in the lower and upper branches to identify and discern the underlying driving mechanisms for these two branches.
Figure 6 depicts the instantaneous contours of the temperature and velocity perturbations of the critical mode $k_\infty =0.049$ at $Gr_\infty =600$, which is representative of the lower branch instabilities. For brevity, only $\epsilon =0.174$ is shown here since all other temperature differences show qualitatively similar structures. From figure 6, it is seen that the velocity fluctuations extend farther out from the wall than the temperature fluctuations, and both instantaneous fluctuation fields have their peaks located beyond the maximum base flow velocity location – indicating that the instability is strongest in the outer part of the boundary layer. There also exists a secondary emerging peak in the velocity perturbation field at $y\approx 16$, which can be seen more clearly in figure 6(c), where the Reynolds normal stress $\overline {\tilde {u}\tilde {u}}$ profile is normalised by its instantaneous maximum (denoted by subscript $_m$). Here, $\bar {\phi }$ denotes the spatial average of flow variable $\phi$ over the homogeneous direction. Although it is difficult to make a direct quantitative comparison due to flow differences, a similar velocity perturbation structure has also been observed in earlier OB stability results for spatially developing natural convection in air (Versteegh & Nieuwstadt Reference Versteegh and Nieuwstadt1998; Aberra et al. Reference Aberra, Armfield, Behnia and McBain2012) and water (Brooker, Patterson & Armfield Reference Brooker, Patterson and Armfield1997; Brooker et al. Reference Brooker, Patterson, Graham and Schöpf2000).
To further understand the energy contributions to the perturbation pattern shown in figure 6, we inspect the two-dimensional perturbation energy $E=(\tilde {u}^2+\tilde {v}^2)/2$ budget for this critical mode at $Gr_{\infty }=600$. With streamwise homogeneity, the energy budget equation can be obtained by summing (2.6b) multiplied by $\tilde {u}$ and (2.6c) multiplied by $\tilde {v}$. The resulting budget equation for our temporally developing NCBL reads
where $\mathcal {P}$ is the shear production, $\mathcal {B}$ is the buoyancy production; $\mathcal {D}_b$ is the base flow diffusion; $\varepsilon _b$ is the base flow viscous pseudo dissipation; $\varPi$ is the pressure transport; $\varSigma$ represents the effects due to spatial variation in fluid density and properties. These terms read
Here, $\mathcal {G}_b$ and $\mathcal {G}_\mu$ are the viscosity gradient effects due to the base flow and local instantaneous viscosity fluctuation, respectively; and $\mathcal {D}_\mu$ is diffusion due to viscosity fluctuation:
Figure 7(a) demonstrates the individual contributors to the energy budget (3.3) for the critical mode $k_\infty =0.049$ at $Gr_\infty =600$. In linear analysis, the perturbations are allowed to grow/decay indefinitely, resulting in a time-dependent scale for the energy budget. However, we note that it is the relative contributions of the terms in (3.4) to the energy budget, rather than the scales, that is of crucial importance to the understanding of the instability mechanisms. In the present study, the energy budgets are normalised in a way such that the base flow viscous diffusion $\mathcal {D}_b$ at the wall ($y=0$) is one. As shown in figure 7(a), the buoyancy production $\mathcal {B}$, being the dominant term in driving the instability, correlates well with the peak of the $\overline {\tilde {u}\tilde {u}}$ seen in figure 6(c); whereas the shear production, $\mathcal {P}$, makes only a small contribution at the peak location and acts as a sink for $y>7$. The production terms shown here behave differently to what is observed in fully turbulent flow where the buoyancy production is usually negligible and negative, and the turbulence is predominantly maintained by shear production (for ${{Pr}}<1$, see Janssen & Armfield Reference Janssen and Armfield1996; Versteegh & Nieuwstadt Reference Versteegh and Nieuwstadt1998; Ke et al. Reference Ke, Williamson, Armfield and Komiya2023). The strong buoyancy production is often seen in the onset of the buoyant instabilities of the unstable atmospheric surface layers (e.g. Leclerc et al. Reference Leclerc, Beissner, Shaw, Den Hartog and Neumann1990, although the gravitational vector aligns differently to the mean flow direction). The high level buoyancy production in figure 7(a) suggests that the instability of the critical mode $k_\infty =0.049$ at $Gr_\infty =600$ is buoyancy driven. However, the effects of local viscosity gradient and fluctuations $\varSigma$ only have minimal effect on the balance. The contributors to $\varSigma$, given by (3.5), are shown in figure 7(b). It is seen that the viscosity fluctuation gradient effect $\mathcal {G}_\mu$ is negligible; while other terms are two orders of magnitude smaller than those in figure 7(a). The buoyancy production is then balanced by base flow viscous diffusion $\mathcal {D}_{b}$ and pressure transport $\varPi$, redistributing the production to the near-wall region and the outer free-shear layer. The sum of the contributors, giving the temporal growth of the perturbation kinetic energy $E$, is shown in figure 7(c). The net energy gain peaks beyond the velocity maximum location $\delta _m$ with a secondary peak at $y\approx 16$. These energy peaks correlate well with those of the instantaneous $\overline {\tilde {u}\tilde {u}}$ profile in figure 6(c).
At $Gr_\infty =1000$ and $k_\infty =0.12$, which is representative of the upper branch instabilities, a fundamentally different perturbation field is seen in figure 8. An additional peak in the velocity contour is apparently observed in the near-wall shear layer ($y\approx 2$), as shown in figure 8(b,c). The strong oscillations in the velocity perturbation correlate well with the critical layer locations (marked in violet in figure 8), where the local base flow velocity $U_b\approx 0.55U_m$ (cf. figure 5d) matches the phase velocity $v_p$ of $k_\infty =0.12$ (see also in, e.g. Mack Reference Mack1984; Wall & Wilson Reference Wall and Wilson1996). The existence of such a near-wall peak gives rise to a stronger shear stress and marks the onset of a shear instability for the wall-bounded NCBL (see also figure 8 of Janssen & Armfield Reference Janssen and Armfield1996), rather than the buoyant instability as seen in the lower branch. This can be more clearly seen in the perturbation energy budget, depicted by figure 9(a), where the kinetic energy is shown driven predominantly by shear and buoyancy production near the outer critical layer at $y\approx 8$. Unlike the buoyant instability mode shown in figure 7, for $Gr_{\infty }=1000$ and $k_\infty =0.12$, the shear production grows in strength and demonstrates a similar peak magnitude to that of the buoyancy production $\mathcal {B}$. In the region close to the base flow velocity maximum location $\delta _m$ (green dotted lines in figure 9), $\mathcal {B}$ and $\mathcal {P}$ are slightly negative – the near-wall energy is mainly sustained by the pressure transport from the outer critical layer. Near the production peaks (i.e. outer critical layer), both the base flow viscous dissipation and diffusion are relatively small; and the overall production is balanced by the pressure transport, which redistributes the excessive kinetic energy towards the near-wall and outer shear layers. A similar distribution is also seen by Versteegh & Nieuwstadt (Reference Versteegh and Nieuwstadt1998) for a vertical natural convection in a differentially heated channel, where the near-wall peak of the redistributive pressure term is shown to increase the flow anisotropy near the non-slip wall. The strong pressure transport in the near-wall region (cf. figure 9a) also results in an additional peak for the temporal evolution $\partial E/\partial t$ around the near-wall critical layer, shown in figure 9(c). Figure 9(b) demonstrates the terms involving spatial fluid property variations. Similar to the buoyant instability in figure 7(b), these contributors to $\varSigma$ are at least two orders of magnitude smaller than those of the dominant terms in the kinetic energy budgets. While $\mathcal {G}_\mu$ remains negligible, an extra local maximum is developed for $\mathcal {D}_\mu$ and $\mathcal {G}_b$ close to the critical layer in the near-wall region, which is also suggestive of a change in the energy generation in the near-wall region.
With the kinetic energy budgets, one is able to discern the driving mechanisms for arbitrary modes by comparing the relative contribution (defined by the $\mathcal {R}$) of shear production $\mathcal {P}$ to the buoyancy production $\mathcal {B}$ across the boundary layer,
The production ratio $\mathcal {R}$ conveniently distinguishes the flow regimes by the relative strength that drives the instabilities on the neutral curve: $\mathcal {R}>1$ for shear-driven instability and $\mathcal {R}<1$ for buoyancy-driven instability. As shown in figure 10(a), the shear-driven instability occurs at relatively large wavenumbers in the upper branch and it occupies increasingly greater area in the $k$–$Gr$ plane with increasing $Gr_\infty$ as the base flow gains more momentum. The buoyant instability, however, mainly happens at relatively low $Gr_{\infty }$ and small wavenumbers $k_\infty$ in the lower branch. This result is consistent with the findings in a differentially heated cavity (Suslov & Paolucci Reference Suslov and Paolucci1995b). However, we shall point out that the onset of shear and buoyant modes observed by Suslov & Paolucci (Reference Suslov and Paolucci1995b) may be due to a different mechanism – Suslov & Paolucci (Reference Suslov and Paolucci1995b) observed that the flow is first unstable to the shear-driven modes in the upper branch due to the breakdown of odd symmetry at the mid-span of the differentially heated cavity. For NCBLs in differentially heated cavities and channels, the shear-driven instability is typically stationary (with a phase velocity of $v_p=0$), while the buoyancy-driven modes are oscillatory in the Boussinesq limit (Suslov & Paolucci Reference Suslov and Paolucci1995b; McBain & Armfield Reference McBain and Armfield2003). With NOB effects, Suslov & Paolucci (Reference Suslov and Paolucci1995b) showed that the resulting shear-driven mode has a negative wave speed $v_p<0$ (i.e. propagating against the base flow) for cavity flows; whereas in the present study, perturbation waves are travelling in the direction of the base flow ($v_p>0$, cf. figure 5). Notably, for our NCBL, there also exists a region at low wavenumbers where the shear production is negative (resulting in $\mathcal {R}<0$). In this region, the instability is driven entirely by buoyancy (see, e.g. the critical mode budget in figure 7a, where $\mathcal {P}$ has a peak value an order of magnitude smaller than that of $\mathcal {B}$ and acts as a sink in the outer shear layer). Figure 10(b) depicts the temperature difference effects on the $\mathcal {R}$ contours. Evidently, with increasing temperature difference $\epsilon$, the $\mathcal {R}$ contours are shifted to higher $Gr_{\infty }$ and larger $k_\infty$, effectively expanding the buoyancy-driven region in the $k$–$Gr$ plane. Such a result is not surprising as the buoyancy effect becomes increasingly strong with greater temperature differences and therefore results in a broadened buoyancy-driven instability region in figure 10. The growing strength of buoyant instability with increasing $\epsilon$ is also seen by Suslov & Paolucci (Reference Suslov and Paolucci1995b), where the critical point (minimum $Gr_{\infty }$ at which the flow becomes marginally unstable to any perturbation) changes from shear driven to buoyancy driven when $\epsilon$ is large enough. However, for the temperature differences considered here up to $\epsilon =0.347$, we have not observed such a sudden change of the critical mode since the unsteady NCBL is unstable to the buoyancy-driven mode first and larger $\epsilon$ would only strengthen the buoyant instability – we conjecture that such a change in critical mode is more likely to occur in NCBLs with higher Prandtl number where the flow might become unstable to the upper shear-driven branch first.
3.3. Neutral curves re-normalised using the film temperature
In the context of NCBL literature, flow field characterisation parameters, such as Grashof number $Gr$ and Rayleigh number $Ra$, are often evaluated at the film temperature ${T_f=(T_w+T_\infty )/2}$ for practical purposes. Earlier investigations on the NOB NCBL revealed that greater temperature differences would destabilise the gaseous boundary layer flows with a lower transition Grashof number when evaluated at $T_f$ (e.g. Chenoweth & Paolucci Reference Chenoweth and Paolucci1986; Sabhapathy & Cheng Reference Sabhapathy and Cheng1986; Jiin-Yuh & Mollendorf Reference Jiin-Yuh and Mollendorf1988), suggesting an earlier onset of instabilities at larger $\epsilon$ for air. As shown in figure 4(a), our linear stability results show that a higher $Gr_\infty$ is needed to trigger the instability for greater temperature differences $\epsilon$, i.e. increasing $\epsilon$ would stabilise the NCBL. However, these results do not imply that our calculations contradict earlier investigations as the neutral curves in figure 4 are normalised using the ambient temperature. Figure 4(a) can be re-normalised using fluid properties at $T_f$, as shown in figure 11. Here, $k_f$ and $Gr_f$ are the dimensionless wavenumber and the film Grashof number,
As depicted by figure 11, while the neutral curves for the $\varepsilon$-range reported still share a similar shape, it is evident that with increasing $\epsilon$, the critical Grashof number decreases – indicating the air NCBL flow, when evaluated by the film temperature $T_f$, is destabilised when a larger temperature difference is applied to the base flow. Unlike figure 4(a), when normalised by $T_f$, the shear-driven instability demonstrates a greater dependence on the temperature difference $\epsilon$. For instance, the marginally unstable Grashof number for $k_f\approx 0.12$ is decreased from $Gr_f=892$ at $\epsilon =0.035$ to $Gr_{f}=564$ at $\epsilon =0.347$; however, the critical Grashof number for the buoyancy-driven instability in the lower branch decreases slightly to $Gr_{f}=355$ at $\epsilon =0.347$ from $Gr_{f}=445$ at $\epsilon =0.035$. Figure 11 indicates when evaluating at $T_f$, the NOB effects due to heating are to destabilise the air flow, prominent for the shear-driven modes but to a lesser extent for the buoyancy-driven modes. These results are consistent with the earlier linear stability investigations for NCBLs along vertical plates, e.g. Sabhapathy & Cheng (Reference Sabhapathy and Cheng1986) and Jiin-Yuh & Mollendorf (Reference Jiin-Yuh and Mollendorf1988) reported heating (decreasing viscosity and conductivity) would stabilise the liquid flows ($7\leq {{Pr}}\leq 100$). Note the viscosity and conductivity for liquid flows decrease with increasing temperature, which is different from the gaseous flow described by Sutherland's law.
3.4. Direct stability analysis
With the linear stability results obtained in § 3.1, the extent of the linear regime where a perturbation of choice, $k_\infty$, interacts linearly with the base flow can be identified by comparing the direct numerical solution to the perturbed system (2.1) in three-dimensional form. By imposing the boundary conditions as impermeable isothermal wall at $y=0$, and the quiescent ambient at $y=L_y$,
With periodic boundary condition applied in the spanwise $z$ direction, the three-dimensional direct numerical simulation (DNS) is carried out in an $L_x\times L_y \times L_z = 1035L_{s,\infty }\times 207 L_{s,\infty } \times 129 L_{s,\infty }$ computational domain with an $N_x \times N_y \times N_z = 1024\times 256 \times 128$ Cartesian grid (stretched in the $y$ direction but uniform in the homogeneous directions $x$ and $z$, with the same stretching strategy as used in the linear stability runs in § 2). The temperature perturbation, given by (2.7), is fed into the flow at $Gr_\infty =250$ with $n=8$ (so that $k_\infty =0.049$ is obtained) and $A_n=10^{-3}$. The ambient temperature is set to $T_\infty =288\,{\rm K}$ and the bulk temperature difference is set to $\Delta T=100\,{\rm K}$ to enable a direct comparison with the linear stability results for $\varepsilon =0.347$. For the DNS, the perturbation field of a flow variable $\tilde {\phi }$ is obtained by subtracting the homogeneous spatial average of a flow variable $\bar {\phi }$ from its instantaneous value $\phi$, i.e. $\tilde {\phi } = \phi -\bar {\phi }$. Figure 12(a) shows the instantaneous and the fluctuation signal of the temperature at $(x_p,y_p,z_p)=(781.6,0.25,64.5)$. Although the temperature fluctuation demonstrates a similar oscillatory behaviour as shown in the linearised solutions (e.g. figure 3), nonlinear developments become visible at $t\approx 120$ (which corresponds to $Gr_\infty \approx 2.2\times 10^4$). This can be clearly seen in figure 12(b), where the instantaneous amplitude of the temperature fluctuation signal $A_t$ is plotted against the Grashof number $Gr_\infty$. From figure 12(b), the temporal signal follows an adjustment immediately after the temperature perturbation is fed into the flow. Similar adjustments are also seen in the linearised solutions (cf. figure 3) and typically end at $t\approx 50$. During this process, the temperature fluctuation is tuned to an appropriate magnitude and phase condition in a way that it does not follow the linear amplification or decay as the boundary layer receives external artificial perturbation (also known as the receptivity process, see Reed, Saric & Arnal Reference Reed, Saric and Arnal1996; Saric, Reed & Kerschen Reference Saric, Reed and Kerschen2002; Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019). For the unsteady DNS, a development time of $t=50$ for the base flow corresponds to $Gr_{\infty }\approx 7\times 10^3$, as shown in figure 12(b). Beyond this inevitable initial receptivity, the temperature fluctuation is amplified almost linearly in the log–log plot as the mean flow continues to develop – up until the onset of nonlinearity at $Gr_\infty \approx 1.7\times 10^4$.
This limit of the linear range can be seen in figure 13, where the amplification rate, $\sigma _k$, of the LNS is compared with the DNS results. With increasing $\epsilon$, the NOB LNS solutions are shown to have a slightly lower amplification rate than that of the OB case for $Gr_{\infty }<5\times 10^3$; however, at greater $Gr_{\infty }$, the perturbation is amplified at a much faster rate. This indicates that a larger bulk characteristic temperature difference would result in the infinitesimal perturbations being amplified more quickly as the flow continues to evolve in time, and potentially would lead to an earlier onset of nonlinearities and laminar–turbulent transition. The DNS data, however, show good agreement with the linear stability results after its initial receptivity process at $Gr_{\infty }\approx 7\times 10^3$. Notably, there are instantaneous fluctuations showing up in the DNS amplification data as there have been weakly nonlinear effects due to the interactions between the perturbation field and property/density fluctuations, all of which, being the higher order nonlinear terms, have been dropped in the linear disturbance modelling (2.6). However, this nonlinear effect is considered negligibly weak since it only has minimal impact on the amplification rates obtained by DNS. By comparing the amplification rates obtained from the DNS and the linear stability results, the linear regime for $\varepsilon =0.347$ is identified to be up until $Gr_\infty \approx 1.7\times 10^4$ (or, equivalently $Gr_{f}\approx 9.7\times 10^3$ when evaluated using $T_f$) – beyond this point, the linear stability result diverges from the DNS data since nonlinearity takes over the perturbation amplification mechanism. The onset of nonlinearity, when evaluated at film temperature, is found to be earlier than that of the OB case, whose linear range extends up to $Gr=1.3\times 10^4$ for an initial perturbation amplitude of $A_n=10^{-3}$ (Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019). Such a result is indeed consistent with the observation of Siebers et al. (Reference Siebers, Moffatt and Schwind1985), where they experimentally found that for a fixed $T_\infty$, the film Grashof number $Gr_f$ at which transition occurs decreases with increasing wall temperature $T_w$ (and thus the temperature difference $\Delta T$).
4. Conclusions
In this paper, the linear stability of a vertical natural convection boundary layer was investigated using both direct numerical simulation and the linearised disturbance equations. The non-Oberbeck–Boussinesq effects are accounted for by using a weakly compressible formulation (Batchelor Reference Batchelor1953; Paolucci Reference Paolucci1982) with the temperature-dependent fluid viscosity and conductivity modelled by Sutherland's law. The linear stability results show that with increasing temperature difference $\epsilon$, the onset of flow instability occurs at a higher ambient Grashof number $Gr_{\infty }$. However, such a neutral curve is sensitive to the reference temperature used to normalise the flow: when evaluating at the film temperature $T_f = (T_w+T_\infty )/2$, a lower critical film Grashof number $Gr_f$ is found for greater $\epsilon$, suggesting a greater temperature difference destabilises the air flow. This observation at $T_f$ is fully consistent with earlier investigations for liquid NCBL flows (Sabhapathy & Cheng Reference Sabhapathy and Cheng1986; Jiin-Yuh & Mollendorf Reference Jiin-Yuh and Mollendorf1988).
The instability dynamics of the perturbations are further explored by investigating the energy budgets for the marginally unstable modes at $Gr_\infty =600$, $k_\infty =0.049$ and $Gr_\infty =1000$, $k_\infty =0.12$, representing the lower and upper branches of the neutral curve, respectively. By comparing the relative contributions of the shear and buoyancy production to the kinetic energy, we show that the marginally unstable mode near the critical point in the lower branch ($Gr_\infty =600$, $k_\infty =0.049$) is driven by buoyancy; whereas the instability at the upper branch ($Gr_\infty =1000$, $k_\infty =0.12$) is shear driven. With increasing temperature difference, the buoyancy-driven instability region occupies an increasingly greater region in the stability diagram due to stronger buoyancy; while the shear-driven instability region expands with increasing Grashof number as the base flow gains more momentum. It is also shown that the buoyancy-driven modes are more sensitive to the temperature difference when normalised by $T_\infty$; while the shear-driven modes are more sensitive to the temperature difference when normalised by $T_f$.
Finally, the extent of the linear regime is identified by comparing the linear stability results with a three-dimensional DNS for $\epsilon =0.347$ and $k_\infty =0.049$. Although the amplification rate obtained from DNS shows slight instantaneous fluctuations due to higher order nonlinear interactions between the fluid properties and density variations, it still demonstrates good agreement with the linear stability results after its initial receptivity, suggesting the linear regime extends up to $Gr_\infty \approx 1.7\times 10^4$. When evaluating at the film temperature, this corresponds to $Gr_{f}=9.7\times 10^3$ and the onset of nonlinearities is found to be earlier than that of the OB case (Ke et al. Reference Ke, Williamson, Armfield, McBain and Norris2019), which is consistent with earlier observations (e.g. Siebers et al. Reference Siebers, Moffatt and Schwind1985, the film Grashof number $Gr_f$ at which transition occurs decreases with increasing wall temperature $T_w$).
Funding
The authors gratefully acknowledge the support provided by the Australian Research Council under award DP220103209, and the computational resources provided by the Australian National Computational Infrastructure (NCI) and the University of Sydney Informatics Hub.
Declaration of interests
The authors report no conflict of interest.