2.1. Mathematical modelling and governing equations

Figure 1 shows the schematic of RB convection with and without bounded shear flows in viscoelastic fluids. The upper and lower plates extend infinitely in x (streamwise) and z (spanwise) directions, and maintain constant lower and higher temperatures, respectively. Three baseflow configurations are considered, i.e. hydrostatic state, Poiseuille flow and Couette flow. Regarding the viscoelastic effects, we consider the Oldroyd-B fluid since it focuses on the interaction between polymers and hydrodynamic fluctuations (Zhang et al. Reference Zhang, Zhang, Wang, Li, Yu and Li2022) and can capture the essential instability mechanism in the flow (e.g. centre mode instability in Garg et al. (Reference Garg, Chaudhary, Khalid, Shankar and Subramanian2018)).

Figure 1.

Schematic description of the problem.

It is assumed that the flow is incompressible and the Oberbeck–Boussinesq approximation (Chandrasekhar Reference Chandrasekhar1961) holds. Therefore, the governing equations are

(2.1a) \begin{align} &\nabla ^{*}\cdot \boldsymbol{U}^{*}=0, \end{align}

(2.1b) \begin{align} \rho ^{*}\frac{\partial \boldsymbol{U}^{*}}{\partial t^{*}}+\rho ^{*}\left(\boldsymbol{U}^{*}\cdot \nabla ^{*}\right)\boldsymbol{U}^{*}&=-\nabla ^{*}P^{*}+\nabla ^{*}\cdot \boldsymbol{\tau }^{*}+\rho ^{*}g^{*}\zeta ^{*}\left(\Theta ^{*}-\Theta _{1}^{*}\right)\mathbf{1}_{y}\,, \end{align}

(2.1c) \begin{align} \left(\frac{\partial \Theta ^{*}}{\partial t^{*}}+\boldsymbol{U}^{*}\cdot \nabla ^{*}\Theta ^{*}\right)&=\nabla ^{*}\cdot \left(\kappa ^{*}\nabla ^{*}\Theta ^{*}\right),\end{align}

(2.1d) \begin{align} \boldsymbol{\tau }^{*}&=\boldsymbol{\tau }_{s}^{*}+\boldsymbol{\tau }_{p}^{*}\,, \end{align}

(2.1e) \begin{align} \boldsymbol{\tau }_{s}^{*}&=\mu _{s}^{*}\big[\nabla ^{*}\boldsymbol{U}^{*}+(\nabla ^{*}\boldsymbol{U}^{*})^{T}\big], \end{align}

(2.1f) \begin{align} \boldsymbol{\tau }_{p}^{*}+\lambda ^{*}(\boldsymbol{\tau }_{p}^{*})^{\nabla }&=\mu _{p}^{*}\big[\nabla ^{*}\boldsymbol{U}^{*}+(\nabla ^{*}\boldsymbol{U}^{*})^{T}\big], \end{align}

where U * is the velocity field, ρ* the density, P * the pressure, τ * the total stress tensor in viscoelastic fluid, g* the acceleration of gravity, ζ* the thermal expansion coefficient, $ \boldsymbol{\unicode{x03C4}}^{\ast}_{{s}}$ the solvent contribution in τ *, Θ* the temperature, κ* the thermal diffusivity, $ \boldsymbol{\unicode{x03C4}}^{\boldsymbol{\ast}}_{\boldsymbol{p}}$ the polymer contribution in τ *, $ \mu^{\ast}_{{s}}$ the viscosity of the solvent, $ \mu^{\ast}_{{p}}$ the viscosity of the polymer, λ* the relaxation time of the polymer molecules and superscript ^T denotes the transpose of tensor ∇* U *.

Here ( $ \boldsymbol{\unicode{x03C4}}^{\boldsymbol{\ast}}_{\boldsymbol{p}}$ )^∇ is the upper-convected time derivative, defined by Larson (Reference Larson1988),

(2.2) \begin{equation} (\boldsymbol{\tau }_{p}^{*})^{\nabla }=\frac{\partial \boldsymbol{\tau }_{p}^{*}}{\partial t^{*}}+\nabla ^{*}\cdot (\boldsymbol{U}^{*}\boldsymbol{\tau }_{p}^{*})-(\nabla ^{*}\boldsymbol{U}^{\boldsymbol{*}})^{T}\cdot \boldsymbol{\tau }_{p}^{*}-\boldsymbol{\tau }_{p}^{*}\cdot \nabla ^{*}\boldsymbol{U}^{\boldsymbol{*}}. \end{equation}

Substituting (2.2) into (2.1f ), we have

(2.3) \begin{equation} \boldsymbol{\tau }_{p}^{*}\!+\!\lambda ^{*}\left(\frac{\partial \boldsymbol{\tau }_{p}^{*}}{\partial t^{*}}\!+\!\nabla ^{*}\cdot (\boldsymbol{U}^{*}\boldsymbol{\tau }_{p}^{*})-(\nabla ^{*}\boldsymbol{U}^{\boldsymbol{*}})^{T}\cdot \boldsymbol{\tau }_{p}^{*}-\boldsymbol{\tau }_{p}^{*}\cdot \nabla ^{*}\boldsymbol{U}^{\boldsymbol{*}}\right)\!=\!\mu _{p}^{*}\left[\nabla ^{*}\boldsymbol{U}^{*}+(\nabla ^{*}\boldsymbol{U}^{*})^{T}\right]. \end{equation}

According to Bird et al. (Reference Bird, Curtiss, Armstrong and Hassager1987) and Zhang et al. (Reference Zhang, Lashgari, Zaki and Brandt2013), we have

(2.4) \begin{equation} \boldsymbol{\tau }_{p}^{*}=\frac{\mu _{p}^{*}h}{\lambda ^{*}k_{B}\Theta ^{*}}(\boldsymbol{C}^{*}-\boldsymbol{C}^{*eq}),\quad \boldsymbol{C}^{*eq}=\frac{k_{B}\Theta ^{*}}{h},\quad\boldsymbol{\tau }_{p}^{*}=\frac{\mu _{p}^{*}}{\lambda ^{*}}\left(\boldsymbol{C}-\boldsymbol{I}\right), \end{equation}

where h is the spring constant of the elastic dumbbell, k _B the Boltzmann constant and C * the conformation tensor.

2.2. Non-dimensionalized governing equations

To non-dimensionalize the governing equations above, we introduce the following scalings:

(2.5) \begin{align} \left(x,y,z\right)&=\frac{\left(x^{*},y^{*},z^{*}\right)}{H^{*}},\boldsymbol{U}=\frac{\boldsymbol{U}^{*}}{U_{0}^{*}},t=t^{*}/\left(\frac{H^{*}}{U_{0}^{*}}\right),P=\frac{P^{*}}{\rho ^{*}U_{0}^{*^{2}}},\\ \nonumber \boldsymbol{\tau }_{p}&=\boldsymbol{\tau }_{p}^{*}/\left(\frac{\mu _{p}^{*}U_{0}^{*}}{H^{*}}\right),\boldsymbol{C}=\boldsymbol{C}^{*}/\left(\frac{k_{B}\Theta ^{*}}{h}\right),\Theta =\frac{\Theta ^{*}-\Theta _{1}^{*}}{\Theta _{0}^{*}-\Theta _{1}^{*}}, \end{align}

where $ U^{\ast}_{0}$ = (g*ζ*H* $\varDelta$ Θ*)^0.5 for the hydrostatic case, and $ U^{\ast}_{0}$ = U _max for the sheared cases. The use of these two velocity scalings does not affect the value of instability threshold and transient growth.

After non-dimensionalization, (2.1)–(2.4) read

(2.6a) \begin{align} &\nabla \cdot \boldsymbol{U}=0, \end{align}

(2.6b) \begin{align} \frac{\partial \boldsymbol{U}}{\partial t}+\left(\boldsymbol{U}\cdot \nabla \right)\boldsymbol{U}=-\nabla P&+\frac{\beta }{Re}\nabla ^{2}\boldsymbol{U}+\frac{1-\beta }{Re^{2}E}\nabla \cdot \boldsymbol{C}\mathbf{+}\frac{Ra}{Re^{2}Pr}\Theta \mathbf{1}_{\boldsymbol{y}}\,, \end{align}

(2.6c) \begin{align} \frac{\partial \Theta }{\partial t}+\boldsymbol{U}\cdot \nabla \Theta & =\frac{1}{RePr}\nabla ^{2}\Theta,\end{align}

(2.6d) \begin{align} \frac{\partial \boldsymbol{C}}{\partial t}+\boldsymbol{U}\cdot \nabla \boldsymbol{C}-\boldsymbol{C}\cdot \nabla \boldsymbol{U}&-\left(\nabla \boldsymbol{U}\right)^{T}\cdot \boldsymbol{C}=-\frac{\boldsymbol{C}-\boldsymbol{I}}{ReE},\end{align}

where

(2.7) \begin{equation} Re=\frac{\rho ^{*}U_{0}^{*}H^{*}}{\mu ^{*}}, Ra=\frac{g^{*}\zeta ^{*}(\Theta _{0}^{*}-\Theta _{1}^{*})H^{{*^{3}}}}{\kappa ^{*}\nu ^{*}}, Pr=\frac{\nu ^{*}}{\kappa ^{*}}, \beta =\frac{\mu _{s}^{*}}{\mu ^{*}}, E=\frac{\lambda ^{*}\mu ^{*}}{\rho ^{*}H^{{*^{2}}}}. \end{equation}

Here, Re is Reynolds number which measures the relative importance of the inertia force to the viscous force. The Rayleigh number Ra represents the ratio between buoyancy force and viscous force. The Prandtl number Pr is the ratio between the kinematic viscosity and the thermal diffusion coefficient. Here β is viscosity ratio, defined as the ratio of the solvent viscosity to the total viscosity. Here E is elasticity number, which measures the magnitude of elasticity in a viscoelastic fluid. Besides, the kinematic viscosity is defined as ν* = μ*/ρ*, and the total dynamic viscosity is defined as μ* = $ \mu^{\ast}_{{s}} + \mu^{\ast}_{{p}}$ .

The dimensionless boundary conditions are Θ(−1) = 1, Θ(1) = 0, U(±1) = 0 for the hydrostatic case and Poiseuille flow, and U(1) = 1, U(−1) = -1 for Couette flow, where the vertical distance y $\in$ [−1,1] after non-dimensionalization.

2.3. Linear stability problem

We decompose the flow variable into the sum of the base state and perturbation, i.e. $\boldsymbol{U}=\overline{\boldsymbol{U}}+u, P=\overline{P}+p, \Theta =\overline{\Theta}+\theta$ and $\boldsymbol{C}=\overline{\boldsymbol{C}}+c$ , where the uppercase letters with a horizontal above represent the base state quantity and the lowercase letters are perturbations. After substituting these decompositions into (2.6), the linear system reads

(2.8a) \begin{equation} \nabla \cdot \boldsymbol{u}=0, \end{equation}

(2.8b) \begin{equation} \frac{\partial \boldsymbol{u}}{\partial t}+\left(\boldsymbol{u}\cdot \nabla \right)\overline{\boldsymbol{U}}+\left(\overline{\boldsymbol{U}}\cdot \nabla \right)\boldsymbol{u}=-\nabla p+\frac{\beta }{Re}\nabla ^{2}\boldsymbol{u}+\frac{1-\beta }{Re^{2}E}\nabla \cdot \boldsymbol{c}+\frac{Ra}{Re^{2}Pr}\theta \mathbf{1}_{\boldsymbol{y}}\,, \end{equation}

(2.8c) \begin{equation} \frac{\partial \theta }{\partial t}+\overline{\boldsymbol{U}}\cdot \nabla \theta +\boldsymbol{u}\cdot \nabla \overline{\Theta}=\frac{1}{RePr}\nabla ^{2}\theta, \end{equation}

(2.8d) \begin{equation} \frac{\partial \boldsymbol{c}}{\partial t}+\overline{\boldsymbol{U}}\cdot \nabla \boldsymbol{c}+\boldsymbol{u}\cdot \nabla \overline{\boldsymbol{C}}-\overline{\boldsymbol{C}}\cdot \nabla \boldsymbol{u}-\boldsymbol{c}\cdot \nabla \overline{\boldsymbol{U}}-\left(\nabla \boldsymbol{u}\right)^{T}\cdot \boldsymbol{C}-\left(\nabla \overline{\boldsymbol{U}}\right)^{T}\cdot \boldsymbol{c}=-\frac{\boldsymbol{c}}{ReE}, \end{equation}

with the accompanying boundary conditions being u(±1) = 0 and θ(±1) = 0.

Based on (2.6), the base-state velocity field is $\overline{\boldsymbol{U}}=0$ (hydrostatic state), $\overline{\boldsymbol{U}}=1-y^{2}$ (Poiseuille flow) and $\overline{\boldsymbol{U}}=y$ (Couette flow), and the base temperature field is $\overline{\Theta}=\mathrm{0}.\mathrm{5}(1-y)$ . Moreover, the base conformation tensor is

(2.9) \begin{equation} \overline{C}=\left(\begin{array}{c@{\quad}c@{\quad}c} 2\left(ReE\right)^{2}\overline{U}^{\prime 2}+1 & ReE\overline{U}' & 0\\ ReE\overline{U}' & 1 & 0\\ 0 & 0 & 1 \end{array}\right). \end{equation}

Since the flow is homogeneous in the streamwise and spanwise directions, it is assumed that the perturbations are wave-like (periodic in these two directions), i.e.

(2.10) \begin{equation} \boldsymbol{f}\left(x,y,z,t\right)=\,\tilde{\!\boldsymbol{f}}\left(y,t\right)e^{\left(i\alpha x+i\gamma z\right)}=\,\hat{\!\boldsymbol{f}}\left(y\right)e^{-i\omega t}e^{\left(i\alpha x+i\gamma z\right)}, \end{equation}

where f can be an arbitrary perturbation, the real-valued α and γ are streamwise and spanwise wavenumbers, respectively, $\,\tilde{\!\boldsymbol{f}}$ and $\,\hat{\!\boldsymbol{f}}$ are the shape functions and ω is the complex circular frequency of disturbance, with its imaginary part ω _i representing the growth rate of linear disturbance (if ω _i > 0, the fluid system is linearly unstable, and ω _i < 0 indicates the stable flow). Using the wall-normal vorticity η = ∂u/∂z – ∂w/∂x to eliminate the pressure term, we get

(2.11) \begin{equation} \boldsymbol{A}\frac{\partial }{\partial t}\tilde{\boldsymbol{\chi }}=\boldsymbol{B}\tilde{\boldsymbol{\chi }}, \end{equation}

where $\tilde{\boldsymbol{\chi }}= (\tilde{{v}},\tilde{{\eta }},\tilde{\theta },\tilde{{c}}_{11},\tilde{{c}}_{22},\tilde{{c}}_{33},\tilde{{c}}_{12},\tilde{{c}}_{13},\tilde{{c}}_{23})$ ^T, and the elements in matrices A and B can be easily deduced from (A1)–(A9) in Appendix A.

For modal analysis, (2.11) can be further simplified as

(2.12) \begin{equation} -i\omega \hat{\boldsymbol{\chi }}=\boldsymbol{L}\hat{\boldsymbol{\chi }}, \end{equation}

where L = A ⁻¹ B , since matrix A is invertible, and the linear modal question becomes an eigenvalue problem.

2.4. Energy norm and input–output non-modal analysis

To predict the transient evolution of perturbation, we transform the eigenvalue problem (2.12) into an initial value problem, i.e.

(2.13) \begin{equation} \frac{\partial \tilde{\boldsymbol{\chi }}}{\partial t}=\boldsymbol{L}\tilde{\boldsymbol{\chi }}. \end{equation}

In this study, we focus on the non-modal instability related to the inertial turbulence instead of the elastic turbulence. Therefore, we only consider the temporal evolution of the perturbed kinetic energy and thermal energy but neglect the elastic energy stored in polymers. The case that polymer disturbances are input directly and elastic potential energy is defined is discussed in Appendix D. Following Hoda et al. (Reference Hoda, Jovanović and Kumar2008), Zhang et al. (Reference Zhang, Lashgari, Zaki and Brandt2013) and Ibarra & Park (Reference Ibarra and Park2023), at the initial time t = 0⁺, the polymer is not perturbed directly but indirectly forced through Hamiltonian operator ∇. In this case, the energy norm is defined as

(2.14) \begin{equation} \int _{{\Omega} }\varepsilon {\mathrm{d}}V=\frac{1}{2}\times \frac{1}{2}\int _{{\Omega} }(A_{1}(\left| \tilde{u}\right| ^{2}+\left| \tilde{{v}}\right| ^{2}+| \tilde{w}| ^{2})+A_{2}| \tilde{\theta }| ^{2}){\mathrm{d}}y=\int _{{\Omega} }\tilde{\boldsymbol{\chi }}^{\dagger }\boldsymbol{M}\tilde{\boldsymbol{\chi }}{\mathrm{d}}y, \end{equation}

where $\unicode{x1D6FA}$ = 2ab, a = 2π/α, b = 2π/γ, and the superscript † stands for the conjugate transpose. In addition, A ₁ = 1, and A ₂ = Ra/(Re ² Pr) (Sameen & Govindarajan Reference Sameen and Govindarajan2007; Jerome et al. Reference Jerome, Chomaz and Huerre2012).

The details of the energy weight matrix M are given as

(2.15) \begin{equation} \boldsymbol{M}=\frac{1}{4}diag\left(\boldsymbol{I}+\frac{1}{k^{2}}\boldsymbol{D}_{1}^{\dagger }\boldsymbol{D}_{1},\frac{1}{k^{2}}\boldsymbol{I},\frac{Ra}{Re^{2}Pr}\boldsymbol{I},\boldsymbol{Z},\boldsymbol{Z},\boldsymbol{Z},\boldsymbol{Z},\boldsymbol{Z,Z}\right), \end{equation}

where I denotes a N× N identity matrix, Z a N× N zero matrix, k ² = α ² + γ ² and D ₁ the first-derivation matrix with respect to the vertical direction y.

Applying a Cholesky decomposition to the weight matrix M , i.e.

(2.16) \begin{equation} \boldsymbol{M}=\boldsymbol{F}^{\dagger }\boldsymbol{F} \end{equation}

and defining $\tilde{\boldsymbol{\xi }}=\boldsymbol{F}\tilde{\boldsymbol{\chi }}$ , we arrive at

(2.17) \begin{equation} \int _{{\Omega} }\varepsilon {\mathrm{d}}V=\int _{{\Omega} }\tilde{\boldsymbol{\chi }}^{\dagger }\boldsymbol{M}\tilde{\boldsymbol{\chi }}{\mathrm{d}}y=\int _{\Omega }\tilde{\boldsymbol{\chi }}^{\dagger }\boldsymbol{F}^{\dagger }\boldsymbol{F}\tilde{\boldsymbol{\chi }}{\mathrm{d}}y=\int _{\Omega }\tilde{\boldsymbol{\xi }}^{\dagger }\tilde{\boldsymbol{\xi }}{\mathrm{d}}y=\left\| \tilde{\boldsymbol{\xi }}\right\| _{2}^{2}\,, \end{equation}

where || · ||₂ represents the L ₂-norm.

To perform the input–output analysis, the input filter B _in and the output filter C _out should be defined. The input initial condition is $\tilde{\boldsymbol{\chi }}$ _in = B _in $\tilde{\boldsymbol{\chi }}^0_{in}$ (t = 0) and the output is $\tilde{\boldsymbol{\chi }}$ _out = C _out $\tilde{\boldsymbol{\chi }}^0_{out}$ , where $\tilde{\boldsymbol{\chi }}^0_{out}$ = e ^{t L} $\tilde{\boldsymbol{\chi }}$ _in = e ^{t L} B _in $\tilde{\boldsymbol{\chi }}^0_{in}$ (t = 0). Therefore, the output becomes $\tilde{\boldsymbol{\chi }}$ _out = C _out $\tilde{\boldsymbol{\chi }}^0_{out}$ = C _out e ^{t L} B _in $\tilde{\boldsymbol{\chi }}^0_{in}$ (t = 0).

We redefine the output energy weight matrix M _out = C _out $\boldsymbol{MC}^{T}_{out}$ and the input energy weight matrix M _in = $\boldsymbol{B}^{T}_{in}$ MB _in . Reapplying a Cholesky decomposition to M _in and M _out : M _in = $\boldsymbol{F}^{\dagger}_{in}$ F _in , M _out = $\boldsymbol{F}^{\dagger}_{out}$ F _out , the maximum energy growth over all initial conditions with an input–output matrix is

(2.18) \begin{equation} \begin{array}{ll} G\left(t\right) & =\max\limits_{\tilde{\boldsymbol{\xi }}_{in}^{0}(t\,=\,0)}\frac{\left\| \tilde{\boldsymbol{\xi }}_{out}\left(t\right)\right\| _{2}^{2}}{\left\| \tilde{\boldsymbol{\xi }}_{in}^{0}(t\,=\,0)\right\| _{2}^{2}}=\max\limits _{\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)}\frac{\left\| \boldsymbol{F}_{out}\tilde{\boldsymbol{\chi }}_{out}\left(t\right)\right\| _{2}^{2}}{\left\| \boldsymbol{F}_{in}\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)\right\| _{2}^{2}}\\ &=\max\limits _{\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)}\frac{\left\| \boldsymbol{F}_{out}\boldsymbol{C}_{out}e^{t\boldsymbol{L}}\boldsymbol{B}_{in}\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)\right\| _{2}^{2}}{\left\| \boldsymbol{F}_{in}\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)\right\| _{2}^{2}}\\ & =\max\limits _{\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)}\frac{\left\| \boldsymbol{F}_{out}\boldsymbol{C}_{out}e^{t\boldsymbol{L}}\boldsymbol{B}_{in}\boldsymbol{F}_{in}^{-1}\boldsymbol{F}_{in}\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)\right\| _{2}^{2}}{\left\| \boldsymbol{F}_{in}\tilde{\boldsymbol{\chi }}_{in}^{0}(t\,=\,0)\right\| _{2}^{2}}=\left\| \boldsymbol{F}_{out}\boldsymbol{C}_{out}e^{t\boldsymbol{L}}\boldsymbol{B}_{in}\boldsymbol{F}_{in}^{-1}\right\| _{2}^{2}\, . \end{array} \end{equation}

The maximum value on the G(t) curve is defined as G _max , and the input initial condition that reaches this peak G _max is the optimal initial condition. The transient growths (of all disturbances) corresponding to different initial conditions are defined in table 1, while the transient growths of various energy components are defined in table 2. It is easy to know the expression of various B _in and C _out , see, for instance, Zhang et al. (Reference Zhang, Martinelli, Wu, Schmid and Quadrio2015).

Table 1.

The definition of input filters and the definition of its resulting transient growth, with C _out fixed as C _out4.

Table 2.

The definition of output filters and the definition of its resulting transient growth, with B _in fixed as B _in4.

According to Schmid & Henningson (Reference Schmid and Henningson2001), using singular value decomposition, the first column of T is extracted as the optimal initial condition and the first column of R is the optimal response, i.e.

(2.19) \begin{equation}\boldsymbol{F}_{out}\boldsymbol{C}_{out}e^{{t_{\textit{max} }}\boldsymbol{L}}\boldsymbol{B}_{in}\boldsymbol{F}_{in}^{-1}=\boldsymbol{R}\sum \boldsymbol{T}^{\dagger}\rightarrow \boldsymbol{F}_{out}\boldsymbol{C}_{out}e^{{t_{\textit{max} }}\boldsymbol{L}}\boldsymbol{B}_{in}\boldsymbol{F}_{in}^{-1}\boldsymbol{T}=\boldsymbol{R}\sum \end{equation}

where t _max is the time to reach G _max .

The numerical methods and code validation can be found in Appendix C.

4.1. Low Re cases

In the presence of shear, it is fundamental to identify the orientation of disturbances that achieves the largest energy growth. To this end, figure 11 shows the contours of G _max in the α–γ plane, from which it is seen that the presence of shear breaks the symmetry of G _max with respect to α and γ. Specifically, regardless of whether it is Couette flow or Poiseuille flow imposed, streamwise-uniform disturbances always can achieve the largest energy amplification in both weakly and strongly elastic regimes, despite of different α dependence of G _max in these two regimes. It should be noted that in the RB–Couette flow of Newtonian fluid, it is nearly longitudinal rolls that have the largest transient growth (Jerome et al. Reference Jerome, Chomaz and Huerre2012). The difference between our results and those of Jerome et al. (Reference Jerome, Chomaz and Huerre2012) is due to the existence of viscoelastic effects. It is well known that in RB–Couette/Poiseuille flow of Newtonian fluids, the mechanism of transient energy growth is identical to pure shear flows. Therefore, similar to the pure Couette flow, in RB–Couette flow of Newtonian fluid the nearly longitudinal rolls can achieve the largest transient growth due to the lift-up mechanism (Jerome et al. Reference Jerome, Chomaz and Huerre2012). However, RB in viscoelastic fluids can generate transient growth even in the absence of shear, as discussed in the last section. By comparing figure 11(a,d) and figure 11(b,c,e,f), it is found that the presence of a shear flow can inhibit the transient growth of viscoelastic RB convection with α> 0. In this case, longitudinal rolls can have the largest transient growth caused by the viscoelasticity effects. On the other hand, in pure Couette flow of both Newtonian fluids (Trefethen 1993; Schmid & Henningson Reference Schmid and Henningson2001) and viscoelastic fluids (Hoda et al. Reference Hoda, Jovanović and Kumar2008), nearly longitudinal rolls have the largest transient growth when α is close to zero. The superposition of transient growth by viscoelasticity and that by the lift-up mechanism makes longitudinal rolls instead of nearly longitudinal rolls have the largest transient growth in RB–Couette mixed flow of viscoelastic fluids.

Figure 11.

Contours of G _max in α-γ plane for (a) E = 4, hydrostatic state, (b) E = 4, Poiseuille flow, (c) E = 4, Couette flow, (d) E = 12, hydrostatic state, (e) E = 12, Poiseuille flow, (f) E = 12, Couette flow; and the other parameters are β = 0.75, Pr = 10, Re = 0.5 in (b), (c), (e) and (f), Ra = 213 for E = 4 and Ra = 196.846 for E = 12.

Figure 12.

In viscoelastic RB–Couette flow, (a) the effect of Re on different transient growths, (b) VT term and PD term in energy budget analysis as functions of Re (note that VT _rel = VT/(VT + PD) and PD _rel = PD/(VT + PD)) and (c) the transition Reynolds number Re _T as a function of E. Other parameters are β = 0.75, Pr = 10, α = 0, Ra = Ra _c and γ = γ _c .

In order to interpret how temperature gradient, viscoelasticity and shear affect the transient growth of viscoelastic RB-Couette flow, we perform an input–output analysis and transient energy analysis, see figure 12. Let us first have a look at the output. As shown in figure 12(a), there is a transition regarding G _{o_η} . Specifically, once Re exceeds 1.986, G _{o_η} will supersede G _{o_θ} (which is always unchanged during the increase of Re) to be the dominant component of disturbance energy in the optimal response. As a result, the overall energy output G is increased with larger values of Re. It is interesting that in the presence of shear, the optimal initial condition is nearly unchanged, which is deduced by inspecting the curve marked by G _{i_θ} shown in figure 12(a). Comparing G _{i_θ} and G, it is found that the temperature disturbance is always dominant in the optimal initial condition regardless of the values of Re. This means that the imposing of shear into a viscoelastic RB convection does not affect how the non-modal growth should be triggered but significantly increases the magnitude of this non-modal growth. Figure 12(a) illustrates how this increase is achieved from the aspect of output energy, but the energy source to amplify disturbances is still unknown. Figure 12(b) gives the answer by showing the results of energy budget analysis. In order to highlight the transition, we only show VT and PD which are normalized as VT _rel = VT/(VT + PD) and PD _rel = PD/(VT + PD). It is surprising that near the transition it is PD instead of Prod (Prod is close to 0, not shown) that supersedes VT to be the dominant energy source of transient amplification. Apparently, Prod being nearly zero is because Re is too small here. In the high-Re regime examined in the next section, the energy production from the background shear is definitely non-negligible, see, for instance, figure 17. We postpone to the next section the explanation of how and why polymeric stresses can make considerable positive work to the amplification of disturbances even if they are not perturbed initially. Here, it is noted that the transition showing in figure 12(b) occurs at a lower Re than that in figure 12(a). This hysteresis is because G _max is dependent on the integration of energy input (here it is VT and PD), rather than the energy input itself. Figure 17(a) is a perfect annotation of this hysteresis.

We use Re _T at which G _{o_η} exceeds G _{o_θ} to denote the transition between different regimes. It is seen from figure 12(c) that this transitional Re, i.e. Re _T , is a decreasing function of E. Specifically, there is no sudden jump of Re _T at E _c and after the Hopf bifurcation, Re _T asymptotes to a constant.

Figure 13 illustrates the transient growth of the sheared viscoelastic RB with α = 0, which shares some analogies with the unsheared cases. As shown in figure 13(a), G evolves smoothly in stationary instability regime (E = 4) but periodically in oscillatory instability regime (E = 12), regardless it is Poiseuille flow or Couette flow imposed, and whether Re < Re _T (i.e. Re = 0.5) or Re > Re _T (i.e. Re = 10). On the other hand, as shown in figure 13(c), the transient growth is still maximized at E = E _c . It is worth noting that while only shears of weak intensity are imposed, e.g. Re = 10, the transient growth at E = E _c is surprisingly enhanced by at least an order of magnitude. For instance, at parameter settings considered in figure 13(c), G _max = 580 in pure RB convection, G _max = 4570 in RB-Poiseuille flow and G _max = 12 750 in RB-Couette flow. It is clear that the transient growth here also benefits significantly from the viscoelastic effects and thermal effects, since the pure shear flows in this Re range cannot generate any transient growth, see figure 13(b). In addition, it is noted that different from the cases dominated by the viscoelastic effects at Re = 0.5, here t _max at Re = 10 is not a decreasing function of E at E > E _c . As shown in figure 13(d), for sufficiently large E, t _max in fact increases with larger values of E, despite that this increase has no influence on the maximum energy gain in the optimal response, i.e. G _max .

Figure 14 further studies the Re dependence of the transient growth of the viscoelastic mixed convection. As expected, in both cases of Poiseuille flow and Couette flow imposed, G _max and t _max are increasing functions of Re. Furthermore, it is surprising that except for the very small Re, i.e. the regime around Re = Re _T , the scaling laws G _max ∼Re ² and t _max ∼Re hold, which were originally found in pure bounded shear flows of high Re (Schmid & Henningson Reference Schmid and Henningson2001), are held here. Specifically, regarding the Re dependence of G _max , the slope at oscillatory instability regime (E = 12) is much larger than that at E = 4 (stationary instability regime); note that G _max at E = 4 is close to G _max at E = 12 in the hydrostatic case. Regarding the Re dependence of t _max , the slope is nearly independent of the type of shear flows and also independent of the elasticity, while t _max is affected by the type of shear flow in purely bounded Newtonian shear flows (Schmid & Henningson Reference Schmid and Henningson2001). Therefore, t _max here is only affected by the shear intensity, as shown in figure 14(c,d).

Figure 13.

Transient energy growth G of streamwise-uniform disturbances for (a) Ra = 213 (for E = 4), Ra = 196 (for E = 12) and (b) Ra = 0; effects of elasticity number E on (c) G _max and (d) t _max . Other parameters are α = 0, β = 0.75, Ra = Ra _c (for panels (c, d)), γ = γ _c , Pr = 10 and Re = 0.5, 10. Note that, for E < E _c , ra = 213 is used in (d) to ensure the accuracy of t _max .

Figure 14.

Effects of Re on (a) G _max , E = 4, (b) G _max , E = 12, (c) t _max , E = 4, (d) t _max , E = 12; and the other parameters are α = 0, β = 0.75, Ra = Ra _c , γ = γ _c and Pr = 10.

One remaining notable phenomenon here concerns the differences between the Poiseuille flow and Couette flow. As shown in figure 14(a,b), whether it is a strongly elastic fluid or a weakly elastic fluid, G _max in viscoelastic RB-Couette flow is always greater than that in viscoelastic RB-Poiseuille flow, and this difference becomes more significant as Re increases. Similar results can be found in the plane Poiseuille/Couette flow of viscoelastic fluids studied by Hoda et al. (Reference Hoda, Jovanović and Kumar2008) who assumed that this difference may depend on the correlation between the vortex stretching term and vertical scale y. Specifically, for Newtonian plane Poiseuille/Couette flow, it has been proved that the vortex stretching term is proportional to $\overline{U}^{\prime}\partial _{z}{v}$ , so the vortex stretching term depends linearly on y in Poiseuille flow, i.e. $\overline{U}^{\prime}=-2y$ and it is independent of y in Couette flow, i.e. $\overline{U}^{\prime}=1$ (Landahl Reference Landahl1975; Butler & Farrell Reference Butler and Farrell1992; Hoda et al. Reference Hoda, Jovanović and Kumar2008). Due to the strongly non-trivial coupling between thermal-hydrodynamic equations and polymer molecular constitutive equations, the exact reasons of aforementioned difference between viscoelastic RB–Couette flow and viscoelastic RB-Poiseuille flow cannot be determined at present, but the phenomenon in figure 14(a,b) is indeed consistent with those results of pure shear flows (Landahl Reference Landahl1975; Butler & Farrell Reference Butler and Farrell1992; Hoda et al. Reference Hoda, Jovanović and Kumar2008).

Figure 15.

Effects of Ra and E on transient energy growth of Couette flow for (a) viscoelastic fluid at Ra = 0, (b) Newtonian fluid and (c) viscoelastic fluid at Ra = 100, (d) viscoelastic fluid at E = 12 and (e) effects of β on transient growth at Ra = 213 and E = 1. Other parameters for are Re = 2000, β = 0.75, Pr = 10, α = 0, γ = γ _c .

4.2. High Re cases

Now we turn to investigate the transient energy growth of mixed convection at relatively high Re, i.e. Re∼O(10³). It should be noted that the Poiseuille flow can be linearly unstable to the Tollmien–Schlichting waves at Re = 5772.22 (Orszag Reference Orszag1971) for Newtonian fluid and particularly, to a newly found instability, i.e. the centre mode (Khalid et al. Reference Khalid, Chaudhary, Garg, Shankar and Subramanian2021), at much lower Re for viscoelastic fluid. Therefore, to highlight the non-modal instability of the sheared RB convection of viscoelastic fluids in high-Re regime, we only consider the case of Couette flow imposed, which, in both Newtonian (Schmid & Henningson Reference Schmid and Henningson2001) and viscoelastic fluids (Chaudhary et al. Reference Chaudhary, Garg, Shankar and Subramanian2019), is always linearly stable.

In the case of high Re, pure shear flow can cause a large transient energy growth in both Newtonian (Schmid & Henningson Reference Schmid and Henningson2001; Schmid Reference Schmid2007) and viscoelastic fluids (Hoda et al. Reference Hoda, Jovanović and Kumar2008; Zhang et al. Reference Zhang, Lashgari, Zaki and Brandt2013). In figure 15 we show how the viscoelasticity and temperature gradient affect the transient growth by the lift-up mechanism. Specifically, as shown in figure 15(a), viscoelasticity has an enhancing effect on the transient growth by the lift-up mechanism and with the increase of the elasticity, this enhancement is more significant. Similarly, an unstable temperature gradient can also enhance the transient growth by the lift-up mechanism, as shown in figure 15(b). The difference between these two cases lies in that viscoelasticity can significantly increase G _max with t _max being nearly unchanged, while heat can significantly extend t _max and cause the generation of the second peak which is higher than the first one for large Ra, see figure 15(b). This implies that the ways of viscoelasticity and heat to enhance the transient growth of high-Re shear flow are different. This point will be clarified later. As shown in figure 15(c,d), once the viscoelasticity and an unstable temperature gradient are simultaneously included, the second peak of the G(t) curve will be further increased and the disturbance energy will be most amplified in this peak instead of in the first peak. In this case, G _max is always an increasing function of both E and Ra, which is different from the non-monotonic variation of G _max with E in figure 5. It should be noted that for a concentrated solution, e.g. β = 0.15, the first peak is always the largest (see figure 15 e), indicating that the significant enhancement of the second peak is only effective for a dilute solution.

Figure 16.

(a) Input analysis and (b) output analysis at E = 12 and ra = 100, and the other parameters are the same as in figure 15(c).

Figure 17.

Transient energy analysis at high Re = 2000 for (a) pure Couette flow of Newtonian fluid, (b) RB–Couette flow of Newtonian fluid, (c) pure Couette flow of viscoelastic fluid and (d) RB–Couette flow of viscoelastic fluid. Other parameters are β = 1 (panel (a)) or 0.75 (panels (b), (c) and (d)), E = 12, Pr = 10, α = 0 and γ = γ _c = 1.558 (Newtonian) or 1.66 (viscoelastic).

Here we are interested in the case of dilute solution where the enhancement due to viscoelasticity and heat is related to the generation for the second peak of the G(t) curve. Figures 16 and 17 give an insight into the second peak. From the input analysis shown in figure 16(a), one sees that when wall-normal velocity ${v}$ is perturbed initially, the disturbance energy will be maximized in the first peak, while the disturbance energy will be maximized in the second peak when only temperature is perturbed initially, which indicates the dependence of these two peaks of G(t) on different input disturbances, i.e. the non-normality of the linear operator forming the first peak is mainly related to the hydrodynamics, while the non-normality of the linear operator forming the second peak is mainly related to the perturbed temperature field. Regarding the output analysis shown in figure 16(b), it is seen that, the main energy output of these two peaks is the kinetic energy of perturbed wall-normal vorticity (similar to the transient growth of pure Couette flow) while the thermal energy output is very small. Therefore, the shape of the optimal initial condition here is identical to the unsheared viscoelastic RB convection while in the optimal response, one expects the formation of streamwise streaks as with the pure Couette flow, which is evidenced by figure 18(a,b,c).

Figure 18.

(a) Cross-stream view of perturbed velocity vectors and temperature perturbation of the optimal initial condition (note that this colourbar is used to illustrate temperature perturbation rather than velocity perturbation, and wavelength λ = 2π/γ); (b) cross-stream view of perturbed streamwise velocity $\tilde{{u}}$ at the first peak (above) / the second peak (i.e. optimal response, below); (c) perturbed velocity of optimal response; (d) perturbed conformation tensor at the time of optimal growth and (e) PD components in transient energy analysis. Other parameters are the same as those in figure 17(d).

In the transient energy analysis shown in figure 17, we compare the short-term energy budget of four typical cases aforementioned, i.e. the pure Couette flow of Newtonian fluid shown in figure 17(a), the RB–Couette flow of Newtonian fluid shown in figure 17(b), the pure Couette flow of viscoelastic fluid shown in figure 17(c) and the RB–Couette flow of viscoelastic fluid shown in figure 17(d). As expected, in the first case, the energy from the mean shear is the only energy source to amplify the disturbance and the viscous dissipation is the only energy loss to damp the non-modal instability. Interestingly, it is found that in the second case shown in figure 17(b) the presence of an unstable temperature gradient does not change the dominant energy balance mentioned above, especially VT ≈ 0. However, the presence of the temperature gradient promotes the mean shear to release more energy into the perturbed hydrodynamic field, which is evidenced by the slower decay of Prod in figure 17(b) compared with that in figure 17(a). This result is consistent with the previous idea that the temperature gradient can enhance the classical lift-up mechanism by extending its duration (Jerome et al. Reference Jerome, Chomaz and Huerre2012). Consequently, the action of an unstable temperature gradient results in the second growth (peak) of the disturbance energy, as shown in figure 17(b). The case of the pure Couette flow of viscoelastic fluid shown in figure 17(c) is similar to that of Newtonian fluid shown in figure 17(a), with the presence of polymeric stresses favouring the energy production by the base flow. This point has been well understood in previous studies, e.g. Zhang et al. (Reference Zhang, Lashgari, Zaki and Brandt2013). When polymeric stresses and temperature gradient are simultaneously included, one expects an interesting energy balance shown in figure 17(d). As expected, the amplification of disturbances corresponding to the first peak of G(t) curve results from the energy production by the basis shear. Different from the temperature field whose energy contribution is always negligible as before, here polymeric stresses field can absorb non-negligible energy from the perturbed hydrodynamic field. However, it is surprising that as the second peak is achieved, PD becomes the dominant energy source while at the same time, Prod becomes negligible. This result is rather counterintuitive since the polymeric stresses field is not perturbed initially, i.e., it is impossible for polymeric stresses field to store any energy initially. Therefore, it is necessary to figure out how the disturbance energy is transferred into the polymeric stresses field and then into the perturbed hydrodynamic field. We have known from figure 8 that polymeric stresses field can absorb disturbance energy from the temperature field. Since the temperature field and the polymeric stresses field are not coupled directly, temperature field cannot transfer energy to polymeric stresses field directly. As mentioned before, this energy transfer includes two processes: (i) the energy transfer from temperature field to hydrodynamic field and (ii) energy transfer from hydrodynamic field to polymeric stresses field. Physically, temperature disturbance perturbs the hydrodynamic field, which in turn perturbs the polymeric stresses field. On the other hand, polymeric stresses field can also obtain energy from the base flow when a background shear is present. Specifically, as pointed out by Buza, Page & Kerswell (Reference Buza, Page and Kerswell2022), the polymer energy has four sources, i.e. the energies from the base flow, the perturbed hydrodynamic field (here, to reiterate, this kind of energy partially or mainly results from the temperature field), the polymeric relaxation and polymeric diffusion, see (4.8) in Buza et al. (Reference Buza, Page and Kerswell2022). In the present system, after the energy from the base flow and from the temperature gradient enters the polymeric stresses field, it is partially transferred into the perturbed hydrodynamic field (the remaining part causes the amplification of elastic energy, see figure 18 d), thus driving the amplification of the perturbed wall-normal vorticity. This process is essentially different from the classical lift-up mechanism. Specifically, as shown in figure 18(e), among various components in PD, $\tilde{{c}}$ ₁₂ and $\tilde{{c}}$ ₁₃ have significant influence on the energy transfer between the perturbed polymeric stresses field and velocity fluctuation field. Here $\tilde{{c}}$ ₁₂ causes the energy transfer from the perturbed polymeric stresses field to velocity fluctuation field, which exceeds the opposite energy transfer caused by $\tilde{{c}}$ ₁₃. Therefore, one expects large, positive values of PD in energy budget analysis shown in figure 17(d). Our results are consistent with Hoda et al. (Reference Hoda, Jovanović and Kumar2008, Reference Hoda, Jovanović and Kumar2009) who also found that $\tilde{{c}}$ ₁₂ and $\tilde{{c}}$ ₁₃ can enhance and suppress transient growth of velocity disturbances in channel flows with very high fluid elasticity, respectively.