Hostname: page-component-69cd664f8f-fgph7 Total loading time: 0 Render date: 2025-03-13T05:14:42.608Z Has data issue: false hasContentIssue false

Sensitivity of three-dimensional boundary layer stability to intrinsic uncertainties of fluid properties: a study on supercritical CO2

Published online by Cambridge University Press:  12 March 2025

Jie Ren*
Affiliation:
State Key Laboratory of Explosion Science and Safety Protection, Beijing Institute of Technology, Beijing 100081, PR China Institute of Aerodynamics and Gas Dynamics, University of Stuttgart, Pfaffenwaldring 21, D-70569 Stuttgart, Germany Yangtze Delta Region Academy of Beijing Institute of Technology, Jiaxing 314003, PR China
Yongxiang Wu
Affiliation:
Institute of Aerodynamics and Gas Dynamics, University of Stuttgart, Pfaffenwaldring 21, D-70569 Stuttgart, Germany
Xuerui Mao
Affiliation:
State Key Laboratory of Explosion Science and Safety Protection, Beijing Institute of Technology, Beijing 100081, PR China Yangtze Delta Region Academy of Beijing Institute of Technology, Jiaxing 314003, PR China Advanced Research Institute of Multidisciplinary Sciences, Beijing Institute of Technology, Beijing 100081, PR China
Cheng Wang
Affiliation:
State Key Laboratory of Explosion Science and Safety Protection, Beijing Institute of Technology, Beijing 100081, PR China Yangtze Delta Region Academy of Beijing Institute of Technology, Jiaxing 314003, PR China
Markus Kloker
Affiliation:
Institute of Aerodynamics and Gas Dynamics, University of Stuttgart, Pfaffenwaldring 21, D-70569 Stuttgart, Germany
*
Corresponding author: Jie Ren, [email protected]

Abstract

The intrinsic uncertainty of fluid properties, including the equation-of-state, viscosity and thermal conductivity, on boundary layer stability has scarcely been addressed. When a fluid is operating in the vicinity of the Widom line (defined as the maximum of isobaric specific heat) in supercritical state, its properties exhibit highly non-ideal behavior, which is an ongoing research field leading to refined and more accurate fluid property databases. Upon crossing the Widom line, new mechanisms of flow instability emerge, feasibly leading to changes in dominating modes that yield turbulence. The present work investigates the sensitivity of three-dimensional boundary layer modal instability to these intrinsic uncertainties in fluid properties. The uncertainty, regardless of its source and the fluid regimes, gives rise to distortions of all profiles that constitute the inputs of the stability operator. The effect of these distortions on flow stability is measured by sensitivity coefficients, which are formulated with the adjoint operator and validated against linear modal stability analysis. The results are presented for carbon dioxide at a representative supercritical pressure of approximately 80 bar. The sensitivity to different inputs of the stability operator across various thermodynamic regimes shows an immense range of sensitivity amplitude. A balancing relationship between the density gradient and its perturbation leads to a quadratic effect across the Widom line, provoking significant sensitivity to distortions of the second derivative of the pressure with respect to the density, $\partial ^2 p/\partial \rho ^2$. From an application-oriented point of view, one important question is whether the correct baseflow profiles can be meaningfully analysed by the simplified ideal-fluid model. The integrated modal disturbance growth – the N factor calculated with different partly idealised models – indicates that the answer depends strongly on the thermodynamic regime investigated.

Type
JFM Papers
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution and reproduction, provided the original article is properly cited.
Copyright
© The Author(s), 2025. Published by Cambridge University Press

1. Introduction

Non-ideal compressible fluids are increasingly used in various applications for improved efficiency and reduced pollution (Guardone et al. Reference Guardone, Colonna, Pini and Spinelli2024). Examples include power generation (White et al. Reference White, Bianchi, Chai, Tassou and Sayma2021), heat exchangers (Chai & Tassou Reference Chai and Tassou2020) and fuel injections (Bellan Reference Bellan2020). ‘Non-ideal’ describes fluids that do not conform to the ideal gas equation-of-state and exhibit unique phenomena such as pseudo-boiling (Banuti Reference Banuti2015), heat transfer deterioration (Pizzarelli Reference Pizzarelli2018) and non-classical rarefaction shock waves (Alferez & Touber Reference Alferez and Touber2017). These non-ideal characteristics pose significant challenges to the traditional ideal gas framework used for predicting flows subject to distortions and laminar–turbulent transition (Li et al. Reference Li, Banuti, Ren, Lyu, Wang and Chu2024).

In addition to the complexity of a non-ideal fluid, the difficulty of predicting flow transition is essentially due to its high sensitivity and the multi-fold path from laminar to turbulence (Reshotko Reference Reshotko2008), which depends not only on the flow configuration, but also the form and amplitude of external perturbations present in the environment. Consequently, the dominant mechanisms are varied. In the linear regime, well-known examples are Tollmien–Schlichting waves due to eigenmodal growth of instabilities, the second mode in hypersonic boundary layer flows (Mack Reference Mack1984), cross-flow waves resulting from a three-dimensional swept flow (Saric, Reed & White Reference Saric, Reed and White2003), centrifugal instabilities due to the presence of concave surfaces (Saric Reference Saric1994) and streamwise velocity streaks following non-modal growth (Trefethen et al. Reference Trefethen, Trefethen, Reddy and Driscoll1993; Schmid & Henningson Reference Schmid and Henningson2001), among others.

In a typical linear stability analysis, the growth rate and dispersion relations are obtained for a pre-calculated laminar baseflow. However, actual flows are inevitably affected by numerous extraneous factors that are not thoroughly accounted for by the theoretical model. Examples include free stream turbulence (Hunt & Graham Reference Hunt and Graham1978), particles (Browne et al. Reference Browne, Hasnine, Sayed and Brehm2021), noise (Schneider Reference Schneider2001) and leading-edge contamination (Spalart Reference Spalart1989), to name a few. To connect to realistic configurations, a key question is how robust the analytical results are and to what extent the growth rate will change when certain distortions are present. Additionally, determining the appropriate distortion that leads to a desired transition promotion or delay is crucial for controling purposes.

The above requirement aligns with the operator perturbation theory (Kato Reference Kato2013), a well-developed field in mathematics (Bottaro, Corbett & Luchini Reference Bottaro, Corbett and Luchini2003). In the context of flow instability, seminal works were performed by Pralits et al. (Reference Pralits, Airiau, Hanifi and Henningson2000); Bottaro et al. (Reference Bottaro, Corbett and Luchini2003); Marquet, Sipp & Jacquin (Reference Marquet, Sipp and Jacquin2008); Bagheri, Brandt & Henningson (Reference Bagheri, Brandt and Henningson2009) and Brandt et al. (Reference Brandt, Sipp, Pralits and Marquet2011) for parabolised stability equations, local and global modal stability analyses, feedback control design, and non-modal growth, respectively. By adopting the adjoint equations (see the review by Luchini & Bottaro Reference Luchini and Bottaro2014), a measure of the system’s response to input variations is formulated. The computed sensitivity field indicates the regions where flow distortions most effectively modify the growth rate, thereby pointing to optimal control strategies.

Recent research on sensitivity analysis has extended to account for high-speed boundary layer flows. Park & Zaki (Reference Park and Zaki2019) investigated a Mach-4.5 flat-plate boundary layer, focusing on the sensitivity properties of the fast and slow modes, whose synchronisation gives rise to Mack’s second mode (Fedorov & Tumin Reference Fedorov and Tumin2011). Guo et al. (Reference Guo, Gao, Jiang and Lee2021) recognised two routes for sensitivity: one where distortion influences the baseflow, leading to variations in the linear stability operator, and another where the stability changes directly. Chen et al. (Reference Chen, Wan, Tu, Duan, Li and Chen2024) found that for an inclined blunt cone, the structural sensitive region is located on the windward side, just downstream of the inlet. Poulain et al. (Reference Poulain, Content, Rigas, Garnier and Sipp2024) formulated the sensitivity based on global instability and resolvent analysis. They identified the optimal locations for steady wall blow/suction and heating/cooling, some of which were shown to successfully damp Mack’s first/second modes and boundary layer streaks simultaneously.

In relation to the current study, Brynjell-Rahkola et al. (Reference Brynjell-Rahkola, Shahriari, Schlatter, Hanifi and Henningson2017) conducted a meaningful analysis on the sensitivity of a three-dimensional (3-D) Falkner–Skan–Cooke boundary layer flow to numerical details. Despite significant knowledge gained regarding the sensitivity of boundary layer stability, studies so far have been mostly limited to ideal gases. The idea here aligns with Juniper & Sujith (Reference Juniper and Sujith2018), who emphasised that ‘the systematic approach in adjoint methods requires an accurate thermoacoustic model ’. Recent efforts on supercritical fluids have discovered new inviscid instabilities (see Robinet & Gloerfelt Reference Robinet and Gloerfelt2019, for a short review) occurring during pseudo-boiling, where the fluid shifts from liquid-like to gas-like behaviour (Simeoni et al. Reference Simeoni, Bryk, Gorelli, Krisch, Ruocco, Santoro and Scopigno2010). Pseudo-boiling represents significant non-ideal thermodynamic regions of a supercritical fluid where the phase change vanishes and is replaced by substantial mutations in thermodynamic and transport properties. These highly non-ideal regions are recognised by the Widom line (Banuti Reference Banuti2015), typically defined as the maximum of the isobaric-specific heat ( $C_p$ ).

Linear stability analyses on canonical flows of supercritical fluids explored so far have demonstrated commonalities – the presence of a new inviscid instability that dominates. For example, the binary mixing layer was shown to be destabilised by a new thermodynamically induced instability (Ly & Ihme Reference Ly and Ihme2022), which can affect the performance of supercritical fuel injection systems. Plane Poiseuille and Couette flows can both become inviscidly unstable upon crossing the Widom line (Ren, Fu & Pecnik Reference Ren, Fu and Pecnik2019a ; Bugeat et al. Reference Bugeat, Boldini, Hasan and Pecnik2024). Under similar conditions, two-dimensional (2-D) boundary layers are subject to dual-mode instability, where the flow is dominated by new inviscid instability in addition to the conventional viscous Tollmien–Schlichting waves (Ren et al. Reference Ren, Marxen and Pecnik2019b ). Moreover, in 3-D boundary layers of accelerating flows with wall cooling, the dominating cross-flow (CF) modes are replaced by the 2-D inviscid mode, which has a growth rate significantly more prominent than that of the CF modes, despite the strong favourable pressure gradient (Ren & Kloker Reference Ren and Kloker2022a ).

The above correlative phenomena have motivated recent efforts to develop novel solvers and to understand the fundamental mechanisms. To the best of the authors’ knowledge, Boldini et al. (Reference Boldini, Hirai, Costa, Peeters and Pecnik2025) has recently introduced the first open-source high-order solver ‘CUBENS’ for single-phase non-ideal fluids in canonical geometries. Bugeat et al. (Reference Bugeat, Boldini, Hasan and Pecnik2024) showed that for stratified plane Couette flow, a minimum in the kinematic viscosity of the baseflow profile produces a generalised inflection point that fulfils Fjørtoft’s generalised inviscid instability criterion (Fjørtoft Reference Fjørtoft1950). They further extended Rayleigh’s criterion (Rayleigh Reference Rayleigh1880) to stratified flows, demonstrating that the excess of density-weighted vorticity (attributed to shear and inertial baroclinic effects) relative to its spatial thickness is responsible for the inviscid instability. Specifically, different fluid models featuring an extremum of the kinematic viscosity were devised to verify the generalised stability model.

In contrast to Couette flow, boundary layer flows are non-parallel, requiring further theoretical work to better understand the dual-mode instability (Bugeat, Boldini & Pecnik Reference Bugeat, Boldini and Pecnik2022). By analysing the boundary layer equation, Ren & Kloker (Reference Ren and Kloker2022a ) concluded that in a 3-D boundary layer, the tremendous negative near-wall viscosity gradient $\partial \mu /\partial y$ is responsible for the inflectional shape of the streamwise velocity profile. This viscosity gradient, upon wall cooling, can be mathematically expressed as

(1.1) \begin{equation} \left .\frac {\partial \mu }{\partial y}\right |_{\left (-\right )}=\left .\frac {\partial \mu }{\partial T}\right |_{\left (+\right )}\left .\frac {\partial T}{\partial y}\right |_{\left (+\right )}+\left .\frac {\partial \mu }{\partial \rho }\right |_{\left (+\right )}\left .\frac {\partial \rho }{\partial y}\right |_{\left (-,\:\mathrm {dominating}\right ).} \end{equation}

Along with the increase of $\partial \rho /\partial y$ in the pseudo-boiling regime, an inflection point is established (note that the first term on the right-hand side is positive for a gas or gas-like fluid, cf. figure 3). Meanwhile, the cross-flow component $w_s$ , responsible for the CF instability, is influenced by the boundary layer’s density distributions through the balancing relationship of the centripetal and centrifugal forces on a curved streamline.

A rational equation-of-state (EOS) and laws of transport properties are essential for both the baseflow and the stability operator, to capture the behaviour correctly. Numerically, look-up tables are used to obtain non-ideal fluid properties during the integration of flow equations. Widely used databases include the NIST (National Institute of Standards and Technology) Reference Fluid Thermodynamic and Transport Properties Database (RefProp) (Huber et al. Reference Huber, Lemmon, Bell and McLinden2022) and the open-source library CoolProp (Bell et al. Reference Bell, Wronski, Quoilin and Lemort2014). In both Refprop and Coolprop, the Helmholtz energy form for a fundamental EOS is used:

(1.2) \begin{equation} \alpha (\tau , \delta ) = \alpha ^{\text {id}} + \alpha ^{\text {r}} = \alpha ^{\text {id}} + \sum _{k} N_k \delta ^{{\rm d}_k} \tau ^{t_k} + \sum _{k} N_k \delta ^{{\rm d}_k} \tau ^{t_k} \exp (-\delta ^{l_k}), \end{equation}

where $\tau = T/T_c$ and $\delta = \rho /\rho _c$ are the reduced temperature and density (by critical values), respectively, $\alpha$ is the reduced molar Helmholtz energy, $\alpha ^{\text {id}}$ is the ideal gas contribution, and $\alpha ^{\text {r}}$ is the real-fluid contribution.Additionally, $N_k$ are coefficients obtained by fitting experimental data, and the exponents $d_k$ , $t_k$ and $l_k$ are also determined by regression. To recover the state equation $p=f(\rho , T)$ , the following thermodynamic relationship is used:

(1.3) \begin{equation} \frac {p}{\rho RT} = 1 + \delta \left ( \frac {\partial \alpha ^{\text {r}}}{\partial \delta } \right )_{\tau }. \end{equation}

However, it is important to note that to better represent properties in the critical region, additional terms are necessary in (1.2), which can make the expression extremely complex (sometimes with more than 50 terms). Due to the empirical nature of the equation, RefProp and CoolProp retain the equation-of-state in an implicit form, but one that closely aligns with physical values – achieving uncertainties that approach the level of the underlying experimental data, in line with the objectives of such libraries. Therefore, the constraints on EOS imposed by thermodynamic stability (see Section II.C of Menikoff & Plohr Reference Menikoff and Plohr1989) shall be satisfied. Accurately describing transport properties, especially in non-ideal regimes, is a challenging and ongoing task that depends heavily on precise experimental measurements. These libraries incorporate various empirical and theoretical (fluid-specific) models for transport properties (see Huber et al. Reference Huber, Sykioti, Assael and Perkins2016, for an example of thermal conductivity in CO2).

Figure 1 compares the gradient of thermal conductivity with respect to temperature, one of the inputs for the flow stability operator. Data were generated using RefProp (version 8.0.4) and CoolProp (version 6.4.1). As shown, panel (b) suffers from model imperfections and non-smooth behaviour near the critical temperature. Even far from the critical point, the values differ significantly between the two databases. The overall difference (see panel c) is generally of a similar order as the values, reflecting considerable uncertainties. Ren & Kloker (Reference Ren and Kloker2022b ) demonstrates that even slight adjustments to the stability operator can substantially influence modal growth, emphasising the importance of considering non-ideal gas behaviours. However, the root causes of this sensitivity variation are still unexplored, pointing to the need for a formalised sensitivity framework to gain deeper understanding. To date, the sensitivity characteristics of boundary layer stability in relation to these fluid properties remain unclear.

Figure 1. Comparison of the fluid property $\partial \kappa /\partial T$ (in dimension $\,\rm Wm^{-1}\,\rm K^-{^2}$ ) for carbon dioxide. The values are generated using (a) RefProp and (b) CoolProp. (c) Differences between the values generated by the two look-up tables. The red circle marks the critical point and the red line represents the isobar at 80.

We aim to characterise and quantify the behaviour by examining the sensitivity of linear stability to each of the inputs of the stability operator. This study will highlight the importance of accurately including the thermodynamic and transport properties of the non-ideal fluid, which are often missed in conventional hydrodynamic stability theory when a temperature gradient crosses the Widom line. Additionally, the possibility of simplifying the fluid model in other regimes will be discussed. The organisation of the paper is as follows: § 2 defines the problem, clarifies the coupling between fluid properties, the laminar baseflow and the linear instability, followed by the derivation of sensitivity profiles; § 3 discusses the results for sensitivity; and § 4 presents the conclusions.

Figure 2. Distribution of pressure (left axis) and pressure coefficient (right axis) for (a) wall cooling and (b) heating cases. (c) Pressure–temperature ( $P{-}T$ ) diagram of carbon dioxide. The contours represent the compressibility factor $\bar {z}=p/(\rho RT)$ . Above the critical point (red dot), the Widom line is plotted using a white dashed line. The fluid regimes are considered along the isobar of 80 ( $p/p_c=1.0844$ ). Four groups of cases are shown with yellow (wall-heating) and cyan (wall-cooling) lines. These lines characterise the distribution of flow temperature, with wall and free stream values denoted by $w$ and $\infty$ , respectively.

2. Problem definition and sensitivity

2.1. Governing equations and flow conditions

The flow satisfies the conservation laws for mass, momentum and energy for a generic fluid (Navier–Stokes equations). In Cartesian coordinates and dimensionless form, this reads

(2.1a) \begin{align} \frac {\partial \rho }{\partial t}+\frac {\partial \left (\rho u_{j}\right )}{\partial x_{j}} &= 0, \end{align}
(2.1b) \begin{align} \frac {\partial \left (\rho u_{i}\right )}{\partial t}+\frac {\partial \left (\rho u_{i}u_{j}-\sigma _{ij}\right )}{\partial x_{j}}&=0, \end{align}
(2.1c) \begin{align} \frac {\partial \left (\rho e\right )}{\partial t}+\frac {\partial \left (\rho eu_{j}+q_{j}\right )}{\partial x_{j}}-\sigma _{ij}\frac {\partial u_{i}}{\partial x_{j}} &=0. \end{align}

Here, the subscripts $i, j$ denote vector/matrix components in a 3-D space and follow Einstein’s summation rule. With this denotation, $(x_1,x_2,x_3) =(x,y,z)$ are coordinated along the streamwise, wall-normal and spanwise directions. Similarly, $(u_1,u_2,u_3) =(u,v,w)$ stand for velocity components along $(x,y,z)$ . The stress tensor $\sigma _{ij}$ and heat flux $q_j$ are given by

(2.2a) \begin{align} \sigma _{ij}=\frac {\mu }{Re}\left (\frac {\partial u_{i}}{\partial x_{j}}+\frac {\partial u_{j}}{\partial x_{i}}\right )+\frac {\lambda }{Re}\delta _{ij}\frac {\partial u_{k}}{\partial x_{k}}-p\delta _{ij}, \end{align}
(2.2b) \begin{align} q_{j}=-\frac {\kappa }{RePrEc}\frac {\partial T}{\partial x_{j}}, \end{align}

where $Re$ , $Pr$ and $\textit {Ec}$ are the Reynolds, Prandtl and Eckert numbers, $\delta _{ij}$ stands for the Kronecker delta, $\mu$ for viscosity, $\kappa$ for thermal conductivity, $p$ for pressure, and $\lambda$ is the Lamé constant (set to $-2\mu /3$ in this study). Equations (2.1) and (2.2) are closed with the relations for thermodynamic and transport properties. For a genetical fluid of pure substance, these relations are written as binary functions of density and temperature,

(2.3) \begin{equation} \left [p,e,\mu ,\kappa \right ]=\left [p,e,\mu ,\kappa \right ]\left (\rho ,T\right ). \end{equation}

We consider 3-D laminar boundary layers with a favourable pressure gradient. The pressure coefficient distribution matches the redesigned DLR experiment (Barth, Hein & Rosemann Reference Barth, Hein and Rosemann2018) on CF instability, leading to an established $p(x)$ , as shown in figure 2(a,b). In all cases, the static pressure $p_\infty$ is fixed at 80 bar and the pressure thus follows Bernoulli’s relation. We prescribe various free stream and wall temperatures $(T_\infty , T_w)$ to explore the physics of four representative regimes: liquid-like, pseudo-boiling, gas-like and ideal gas, as shown in figure 2(c). Specifically, relative to the pseudo-critical (pseudo-boiling) temperature $T_{pc}=307.7\,\rm K$ , the boundary temperatures are summarised in table 1. The temperature ratios are maintained at $T_w/T_\infty =16/15$ and $15/16$ for wall heating and cooling, respectively. An identical Reynolds number ( $Re=\rho _\infty U_\infty L_{\textrm {ref}}/\mu _\infty =1.4687\times 10^5$ ) as in the experiment and a low Mach number ( ${\textit {Ma}}=U_\infty /c_\infty =0.2$ ) have been used. The choice of these parameters ensures that the flow has exact comparability with typical cross-flow instabilities in the ideal-gas regime, in which the same dimensionless baseflow and neutral curve is obtained (Dörr & Kloker Reference Dörr and Kloker2017; Ren & Kloker Reference Ren and Kloker2022a ).

Table 1. A summary of flow cases investigated.

2.2. The laminar baseflow

Matching the pressure coefficient distribution of the DLR experiments results in a non-self-similar laminar baseflow where the shape factor (ratio of displacement to momentum-loss thickness) is not constant. This laminar baseflow is obtained by solving the parabolised Navier–Stokes (PNS) equations. In steady and non-separating boundary layer flows, the streamwise viscous gradient is notably smaller than the wall-normal component. Therefore, the PNS equations are derived from the full Navier–Stokes equations by neglecting the streamwise gradient in the viscous terms. The resulting equations in their 2-D form are given by

(2.4) \begin{equation} {\mathsf {A}}\dfrac {\partial \boldsymbol {Q}}{\partial x}+\boldsymbol {\mathsf {B}}\dfrac {\partial \boldsymbol {Q}}{\partial y}=\textrm {RHS}, \end{equation}
(2.5) \begin{equation} \boldsymbol {\mathsf {A}}=\left (\begin{array}{ccccc} u & \rho & 0 & 0 & 0\\ 0 & \rho u & -\dfrac {1}{Re}\dfrac {\partial \mu }{\partial y} & 0 & 0\\ 0 & -\dfrac {1}{Re}\dfrac {\partial \lambda }{\partial y} & \rho u & 0 & 0\\ 0 & 0 & 0 & \rho u & 0\\ \rho u\dfrac {\partial e}{\partial \rho } & p & 0 & 0 & \rho u\dfrac {\partial e}{\partial T} \end{array}\right ), \end{equation}
(2.6) \begin{equation} \boldsymbol {\mathsf {B}}=\left (\begin{array}{ccccc} v & 0 & \rho & 0 & 0\\ 0 & b_{2,2} & 0 & 0 & 0\\ \dfrac {\partial p}{\partial \rho } & 0 & b_{3,3} & 0 & \dfrac {\partial p}{\partial T}\\ 0 & 0 & 0 & b_{4,4} & 0\\ \rho v\dfrac {\partial e}{\partial \rho } & -\dfrac {\mu }{Re}\dfrac {\partial u}{\partial y} & p-\dfrac {2\mu +\lambda }{Re}\dfrac {\partial v}{\partial y} & -\dfrac {\mu }{Re}\dfrac {\partial w}{\partial y} & b_{5,5} \end{array}\right ), \end{equation}
(2.7) \begin{equation} \left .\begin{array}{ll} b_{2,2} = & \rho v-\dfrac {1}{Re}\dfrac {\partial \mu }{\partial y}-\dfrac {\mu }{Re}D\\ b_{3,3} = & \rho v-\dfrac {1}{Re}\dfrac {\partial \left (2\mu +\lambda \right )}{\partial y}-\dfrac {2\mu +\lambda }{Re}D\\ b_{4,4} = & \rho v-\dfrac {1}{Re}\dfrac {\partial \mu }{\partial y}-\dfrac {\mu }{Re}D\\ b_{5,5} = & \rho v\dfrac {\partial e}{\partial T}-\dfrac {1}{RePrEc}\dfrac {\partial \kappa }{\partial y}-\dfrac {\kappa }{RePrEc}D \end{array}\right \}, \end{equation}
(2.8) \begin{equation} \textrm {RHS}=(0,-{\rm d}p/{\rm d}x,0,0,0)^T. \end{equation}

The Prandtl ( $Pr={\mu _{\infty }C_{p\infty }}/{\kappa _\infty }$ ) and Eckert ( ${\textit {Ec}}={u_{\infty }^{2}}/{C_{p\infty }T_\infty }$ ) numbers are not independent and can be calculated based on $Re$ , $\textit {Ma}$ and the temperature conditions prescribed in table 1. The symbol $e$ denotes the internal energy. The operator $D$ in (2.7) stands for the wall-normal derivative. In numerically solving the PNS equations by an implicit Euler scheme, the system is linearised by ‘lagging’ the coefficients $\boldsymbol {\mathsf {A}}$ and $\boldsymbol {\mathsf {B}}$ relative to the solution vector $\boldsymbol {Q}= (\rho ,u,v,w,T )^{T}$ in an iterative procedure. Namely, sub-iterations are carried out to update $\boldsymbol {\mathsf {A}}$ and $\boldsymbol {\mathsf {B}}$ from $\boldsymbol {q}$ at each station of the streamwise marching to obtain the correct, fully nonlinear values. When external perturbations are not present, the boundary conditions are

(2.9a) \begin{align} y=y_e:\;\frac {\partial u}{\partial y}=\frac {\partial w}{\partial y}=0,\;\rho =\rho _{e}\left (x\right ),\;T=T_{e}\left (x\right ); \end{align}
(2.9b) \begin{align} y=0:\;u=v=w=0,\;T=T_{w}. \end{align}

We employ the subscript $e$ to denote local boundary layer edge values. At the upper edge of the boundary layer, both velocity components $u$ and $w$ are subject to Neumann conditions, while the gradient $\partial v/\partial y$ does not vanish due to the presence of streamwise pressure gradients. Thus, $v_e$ is calculated based on the continuity equation. The potential flow values $\rho _{e}(x)$ and $T_{e}(x)$ are given by the isentropic relations (where $S$ stands for entropy):

(2.10) \begin{equation} S\left (\rho _{e}\left (x\right ),p\left (x\right )\right )=S\left (\rho _{\infty },p_{\infty }\right ),\;S\left (T_{e}\left (x\right ),p\left (x\right )\right )=S\left (T_{\infty },p_{\infty }\right ). \end{equation}

At the wall, no-slip velocity and a specified wall temperature apply (see table 1). The wall density is not prescribed; instead, it is allowed to vary to ensure that the pressure gradient at the wall is zero, thereby satisfying the momentum equation in the wall-normal direction. The PNS equations are integrated downstream, starting from an initial profile at $x=x_0$ . In this study, we specify the streamwise and spanwise velocities $u(x_0,y)$ and $w(x_0,y)$ using the Falkner–Skan–Cooke (FSC) solution (Cooke Reference Cooke1950), with $v(x_0,y)=0$ . The thermodynamic variables ( $\rho$ , $T$ ) are either given as the potential-flow values (applying isentropic relations) or extrapolated from existing downstream data, ensuring that the influence of the initial profiles is insignificant.

2.3. Linear instability and the sensitivity framework

Considering flow instability, the flow variables $\tilde {\boldsymbol {q}}= (\tilde {\rho },\tilde {u},\tilde {v},\tilde {w},\tilde {T} )^{T}$ are decomposed into the laminar baseflow $\boldsymbol {Q}$ (steady state) and its perturbations $\boldsymbol {q}^{\prime }$ :

(2.11) \begin{equation} \tilde {\boldsymbol {q}}\left (x,y,z,t\right )=\boldsymbol {Q}\left (x,y\right )+\boldsymbol {q}^{\prime }\left (x,y,z,t\right ). \end{equation}

The stability equations are derived by subtracting the governing equations for $\tilde {\boldsymbol {q}}$ and $\boldsymbol {Q}$ , both of which satisfy the Navier–Stokes equations. We investigate the perturbations in Fourier space by introducing the ansatz:

(2.12) \begin{equation} \boldsymbol {q}^{\prime }\left (x,y,z,t\right )=\hat {\boldsymbol {q}}\left (y\right )\exp \left (i\alpha x+i\beta z-i\omega t\right )+c.c., \end{equation}

with $c.c.$ denoting the complex conjugate. The linearised stability equations are derived and written in a compact form:

(2.13) \begin{equation} \boldsymbol {\mathsf {L}}\left (\boldsymbol {Q},\alpha ,\beta ,\omega ,Re,{\textit {Ma}}\right )\hat {\boldsymbol {q}}=0. \end{equation}

Equation (2.13) constitutes an eigenvalue problem whose dimensions are $(5\times 5)$ before spatial discretisation. The detailed expressions are provided in Appendix A. We investigate the problem with $\boldsymbol {Q}$ , $\beta$ , $\omega$ , $Re$ and $\textit {Ma}$ as inputs, while $\alpha$ and $\hat {\boldsymbol {q}}$ are the eigenvalue and eigenfunction to be solved (spatial mode). We present the constituent elements of operator $\boldsymbol {\mathsf {L}}$ in (2.14). Compared with the scalar parameters $\beta$ , $\omega$ , $Re$ and $\textit {Ma}$ (not shown in (2.14)), $\boldsymbol {Q}$ includes not only the baseflow profiles ( $\rho$ , $u$ , $v$ , $T$ – including boundary conditions, $p$ ) but also the equation-of-state, viscosity, and thermal conductivity, which depend on these profiles. In other words, the intrinsic properties of the fluid, presented in $\boldsymbol {Q}$ , critically influence the stability and determine the type and outcome of the problem.

(2.14) \begin{equation} \boldsymbol {\mathsf {L}}=\underset {\begin{array}{l} {{[b]}}\textrm { - baseflow profiles: }\rho ,u,w,T\\[8pt] {[\textrm {EoS}]}\textrm { - equation of state: }\dfrac {\partial p}{\partial \rho },\dfrac {\partial p}{\partial T},\dfrac {\partial ^{2}p}{\partial \rho ^{2}},\dfrac {\partial ^{2}p}{\partial T^{2}},\dfrac {\partial ^{2}p}{\partial \rho \partial T},\dfrac {\partial e}{\partial \rho },\dfrac {\partial e}{\partial T}\\[8pt] {[\mu ]}\textrm { - viscosity:} {\mu },{\dfrac {\partial \mu }{\partial \rho }},{\dfrac {\partial \mu }{\partial T}},{\dfrac {\partial ^{2}\mu }{\partial \rho ^{2}}},{\dfrac {\partial ^{2}\mu }{\partial T^{2}}},{\dfrac {\partial ^{2}\mu }{\partial \rho \partial T}}\\[8pt] {[\kappa ]}\textrm { - thermal conductivity: }\kappa ,\dfrac {\partial \kappa }{\partial \rho },\dfrac {\partial \kappa }{\partial T},\dfrac {\partial ^{2}\kappa }{\partial \rho ^{2}},\dfrac {\partial ^{2}\kappa }{\partial T^{2}},\dfrac {\partial ^{2}\kappa }{\partial \rho \partial T} \end{array}}{\underbrace {{\left (\begin{array}{ccccc} {{[b]}} & {{[b]}} & {{[b]}} & {{[b]}} & 0\\ {{[b]}}{[\textrm {EoS}]}{[\mu ]} & {{[b]}}{[\mu ]} & {{[b]}}{[\mu ]} & {[\mu ]} & {{[b]}}{[\textrm {EoS}]}{[\mu ]}\\ {{[b]}}{[\textrm {EoS}]}{[\mu ]} & {[\mu ]} & {{[b]}}{[\mu ]} & {[\mu ]} & {{[b]}}{[\textrm {EoS}]}{[\mu ]}\\ {{[b]}}{[\textrm {EoS}]}{[\mu ]} & {[\mu ]} & {{[b]}}{[\mu ]} & {{[b]}}{[\mu ]} & {{[b]}}{[\textrm {EoS}]}{[\mu ]}\\ {{[b]}}{[\textrm {EoS}]}{[\mu ]}{[\kappa ]} & {{[b]}}{[\mu ]} & {{[b]}}{[\textrm {EoS}]}{[\mu ]} & {{[b]}}{[\mu ]} & {{[b]}}{[\textrm {EoS}]}{[\mu ]}{[\kappa ]} \end{array}\right )}}} \end{equation}

Figure 3. An overview of the inputs in the stability operator. Panels (a)–(d) show $[b]$ baseflow profiles, [EoS] equation-of-state, $[\mu ]$ viscosity and $[\kappa ]$ thermal conductivity, respectively. The pseudo-boiling and ideal gas regimes are plotted with solid and dashed lines, respectively (both under wall cooling). The red circle denotes the pseudo-critical point.

Symbols with square brackets have been used to distinguish the four groups of inputs, with the variables associated with each group listed at the bottom of (2.14). Specifically, $[b]$ represents baseflow profiles (density, velocities in the $x$ - and $z$ -directions, and temperature as functions of $y$ ), [EoS] represents the equation-of-state (profiles of thermodynamic derivatives of pressure and internal energy), and $[\mu ]$ and $[\kappa ]$ represent viscosity, thermal conductivity and their gradients with respect to temperature and density. Compared with the Orr–Sommerfeld equation (and also to the compressible perfect/ideal gas set-up), there is an increase in the number of inputs, specifically within the three groups [EoS], $[\mu ]$ and $[\kappa ]$ . To gain an overview of all the inputs, we compare the wall-normal profiles between pseudo-boiling (solid lines) and ideal gas (dashed lines) regimes in figure 3. According to the state postulate, the groups [EoS], $[\mu ]$ and $[\kappa ]$ are functions of $\rho$ and $T$ . In the stage of numerically obtaining the baseflow, $\rho$ and $T$ (and $u$ and $w$ ) in turn depend on the models for [EoS], $[\mu ]$ and $[\kappa ]$ (see § 2.2 and the discussion in § 2.4). From earlier works (see reviews in Guardone et al. Reference Guardone, Colonna, Pini and Spinelli2024), the non-ideal properties of a fluid can drive the baseflow to be inflectional, supporting a new inviscid mode that dominates the instability, which may significantly promote flow transition. However, it remains poorly understood how sensitive the flow stability is to the intrinsic properties of a fluid. Specifically, how important is each input, and to what degree does it influence the results? Is it feasible to ignore some of them or use idealised laws?

The answer will be pursued following the sensitivity framework. We consider a distortion $\delta \boldsymbol {Q}$ to the input $\boldsymbol {Q}$ , and thus (2.13) becomes

(2.15) \begin{equation} \boldsymbol {\mathsf {L}}\left (\boldsymbol {Q}+\delta \boldsymbol {Q},\alpha +\delta \alpha \right )(\hat {\boldsymbol {q}}+\delta \hat {\boldsymbol {q}})=0, \end{equation}

in which $\delta \alpha$ and $\delta \hat {\boldsymbol {q}}$ are corresponding changes in the eigenvalue and eigenfunction induced by $\delta \boldsymbol {Q}$ . For brevity, we have omitted the other fixed parameters ( $\beta$ , $\omega$ , $Re$ and $\textit {Ma}$ ). Taking Taylor’s expansion of operator $\boldsymbol {\mathsf {L}}$ yields

(2.16) \begin{equation} \boldsymbol {\mathsf {L}}\left (\boldsymbol {Q}+\delta \boldsymbol {Q},\alpha +\delta \alpha \right )= \boldsymbol {\mathsf {L}}\left (\boldsymbol {Q},\alpha \right )+\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}+\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \alpha }\delta \alpha +O\left (\delta ^{2}\right ), \end{equation}

where the gradient is defined as

(2.17) \begin{equation} \frac {\partial \boldsymbol {\mathsf {L}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}=\underset {s\rightarrow 0}{\lim }\frac {\boldsymbol {\mathsf {L}}\left (\boldsymbol {Q}+s\delta \boldsymbol {Q}\right )-\boldsymbol {\mathsf {L}}\left (\boldsymbol {Q}\right )}{s}. \end{equation}

We substitute (2.16) into (2.15) and subtract the equation for the undistorted state (2.13). The equation at $O(\delta )$ reads

(2.18) \begin{equation} \frac {\partial \boldsymbol {\mathsf {L}}}{\partial \alpha }\hat {\boldsymbol {q}}\delta \alpha +\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}\hat {\boldsymbol {q}}+\boldsymbol {\mathsf {L}}\delta \hat {\boldsymbol {q}}=0. \end{equation}

Equation (2.18) forms a tripartite relation between the distortion $\delta \boldsymbol {Q}$ and the induced $\delta \alpha$ and $\delta \hat {\boldsymbol {q}}$ . Further to this relation, we define the sensitivity coefficient $\boldsymbol {S}_{\boldsymbol {Q}}$ , which satisfies

(2.19) \begin{equation} \delta \alpha =\left \langle \boldsymbol {S}_{\boldsymbol {Q}},\delta \boldsymbol {Q}\right \rangle . \end{equation}

Here, $\left \langle \right \rangle$ stands for the inner product between two vectors: $\left \langle \boldsymbol {a},\boldsymbol {b}\right \rangle =\int _{0}^{\infty }\boldsymbol {a}^{H}\boldsymbol {b}\:{\rm d}y$ , $H$ implies Hermitian transpose. The sensitivity coefficient, through (2.19), measures the reactivity of the eigenvalue relative to the baseflow distortion $\delta \boldsymbol {Q}$ . To find the analytical expression for $\boldsymbol {S}_{\boldsymbol {Q}}$ , we seek the adjoint problem defined by the following relations:

(2.20) \begin{equation} \left \langle \hat {\boldsymbol {q}}^{\dagger },\boldsymbol {\mathsf {L}}\hat {\boldsymbol {q}}\right \rangle =\left \langle \boldsymbol {\mathsf {L}}^{\dagger }\hat {\boldsymbol {q}}^{\dagger },\hat {\boldsymbol {q}}\right \rangle =0. \end{equation}

The superscript $\dagger$ is used for adjoint variables and operators. The analytical process deriving the adjoint equations has been provided in Appendix B. The inner product of $\hat {\boldsymbol {q}}^{\dagger }$ with (2.18) eliminates its last term, giving

(2.21) \begin{equation} \delta \alpha =-\frac {\left \langle \hat {\boldsymbol {q}}^{\dagger },\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}\hat {\boldsymbol {q}}\right \rangle }{\left \langle \hat {\boldsymbol {q}}^{\dagger },\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \alpha }\hat {\boldsymbol {q}}\right \rangle }. \end{equation}

Upon normalisation of the denominator (see Appendix C),

(2.22) \begin{equation} \left \langle \hat {\boldsymbol {q}}^{\dagger },\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \alpha }\hat {\boldsymbol {q}}\right \rangle =1. \end{equation}

Equation (2.21) is simplified as

(2.23) \begin{equation} \delta \alpha =-\left \langle \hat {\boldsymbol {q}}^{\dagger },\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}\hat {\boldsymbol {q}}\right \rangle . \end{equation}

Comparing (2.23) and (2.19), the sensitivity coefficients $\boldsymbol {S}_{\boldsymbol {Q}}$ are obtained by integrating (2.23) by parts. We provide their specific expressions in Appendix D (the superscript $*$ indicates the complex conjugate of a variable) for all the inputs discussed in figure 3.

2.4. Distortions of the baseflow and coupling to the linear stability

Recall § 2.2, where the governing equations of the baseflow are nonlinear and the coefficient matrix in (2.4) is a function of the solution vector. Here, we rewrite the equation in a more specific form:

(2.24) \begin{equation} \boldsymbol {\mathsf {A}}\left (\boldsymbol {Q},[\mu ], [\kappa ], [\textrm {EoS}]\right )\dfrac {\partial \boldsymbol {Q}}{\partial x}+\boldsymbol {\mathsf {B}}\left (\boldsymbol {Q},[\mu ], [\kappa ], [\textrm {EoS}]\right )\dfrac {\partial \boldsymbol {Q}}{\partial y}=\textrm {RHS}. \end{equation}

We consider uncertainties in the fluid properties, such as distortions in viscosity $\delta \mu$ (as well as in other terms of $[\mu ]$ , $[\kappa ]$ and [EoS]). We solve (2.4) twice (with and without distorted fluid properties) to determine $\delta \boldsymbol {Q}$ , which does not require small distortion amplitudes. Additionally, to understand the influence of $\delta \mu$ on the other components, take Taylor’s expansion of (2.24),

(2.25) \begin{equation} \begin{aligned} \left (\boldsymbol {\mathsf {A}}+\frac {\partial \boldsymbol {\mathsf {A}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}+\frac {\partial \boldsymbol {\mathsf {A}}}{\partial \mu }\delta \mu +O\left (\delta ^{2}\right )\right )\frac {\partial \left (\boldsymbol {Q}+\delta \boldsymbol {Q}\right )}{\partial x}&+\\ \left (\boldsymbol {\mathsf {B}}+\frac {\partial \boldsymbol {\mathsf {B}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}+\frac {\partial \boldsymbol {\mathsf {B}}}{\partial \mu }\delta \mu +O\left (\delta ^{2}\right )\right )\frac {\partial \left (\boldsymbol {Q}+\delta \boldsymbol {Q}\right )}{\partial y}&=\boldsymbol {\textrm {RHS}}. \end{aligned} \end{equation}

In the linear regime (when $\delta$ is small), by subtracting (2.25) and (2.24), keeping terms of order $O (\delta )$ :

(2.26) \begin{equation} \begin{aligned} \boldsymbol {\mathsf {A}}\frac {\partial \delta \boldsymbol {Q}}{\partial x}+\boldsymbol {\mathsf {B}}\frac {\partial \delta \boldsymbol {Q}}{\partial y}=\mathrm {\delta {}_{RHS}} \end{aligned}, \end{equation}

where

(2.27) \begin{equation} \mathrm {\delta {}_{RHS}}=-\left (\boldsymbol {\mathsf {A}}\left (\delta \boldsymbol {Q}\right )+\boldsymbol {\mathsf {A}}\left (\delta \mu \right )\right )\frac {\partial \boldsymbol {Q}}{\partial x}-\left (\boldsymbol {\mathsf {B}}\left (\delta \boldsymbol {Q}\right )+\boldsymbol {\mathsf {B}}\left (\delta \mu \right )\right )\frac {\partial \boldsymbol {Q}}{\partial y}, \end{equation}
(2.28) \begin{equation} \boldsymbol {\mathsf {A}}\left (\delta \boldsymbol {Q}\right )=\dfrac {\partial \boldsymbol {\mathsf {A}}}{\partial \boldsymbol {Q}}\delta \boldsymbol {Q}=\underset {s\rightarrow 0}{\lim }\frac {\boldsymbol {\mathsf {A}}\left (\boldsymbol {Q}+s\delta \boldsymbol {Q}\right )-\boldsymbol {\mathsf {A}}\left (\boldsymbol {Q}\right )}{s}, \end{equation}
(2.29) \begin{equation} \boldsymbol {\mathsf {A}}\left (\delta \mu \right )=\frac {\partial \boldsymbol {\mathsf {A}}}{\partial \mu }\delta \mu . \end{equation}

In this linear regime, $\delta \boldsymbol {Q}$ can be obtained iteratively, as in PNS, by ‘lagging’ the coefficients of $\boldsymbol {\mathsf {A}} (\delta \boldsymbol {Q} )$ and $\boldsymbol {\mathsf {B}} (\delta \boldsymbol {Q} )$ . Scrutinising (2.26), one recognises that a scalar deviation $\delta \mu$ can induce distortions in all components of the baseflow $\delta \boldsymbol {Q}= (\delta \rho , \delta u, \delta v, \delta w, \delta T )^{T}$ through the coupling of the distorted matrices ( $\boldsymbol {\mathsf {A}} (\delta \boldsymbol {Q} )$ , $\boldsymbol {\mathsf {A}} (\delta \mu )$ , $\boldsymbol {\mathsf {B}} (\delta \boldsymbol {Q} )$ , $\boldsymbol {\mathsf {B}} (\delta \mu )$ , and the corresponding raw state $\boldsymbol {Q}$ . Figure 4 provides an overview of the baseflow distortions induced by viscosity alterations: $\delta \mu =\sigma \mu (\rho ,T)$ . Here, $\sigma =\pm 10\,\%$ and $\pm 20\,\%$ , leading to a bulk alteration in the $\mu (y)$ profile as given by the fluid property database. These profiles show that distortions remain within the boundary layer thickness and their amplitudes stay largely constant moving downstream. Specifically, $\delta \mu = \pm 20\,\%\mu$ gives rise to $(\delta u)_{\max }\approx 5\,\% U_\infty$ . The stability diagram showing the imaginary part of the eigenvalue $\alpha _i$ indicates the steady cross-flow mode grows from around $x = 0.5$ , reaching a maximum around $x = 1.0$ and continues growing downstream.

Figure 4. (a) Distortions of the streamwise velocity, $\delta u/u_e$ , induced by viscosity alterations ( $\delta \mu = \pm 10\,\%\mu , \pm 20\,\%\mu$ ). (b) The raw state baseflow profiles of $u$ and $w$ and (c) the corresponding stability diagram of the steady ( $\omega =0$ ) cross-flow instability. In panel (b), the dashed line stands for the boundary-layer thickness based on $0.99u/U_\infty$ , and the dotted lines denote $u/U_\infty =w/U_\infty =1$ . The flow is in the supercritical regime, subject to wall heating, with gas-like fluid properties.

Figure 5. Baseflow distortions induced by uncertainties in the viscosity model. The flow is in the supercritical regime. Viscosity distortion was applied in proportion to the original model, ranging from 0 % to 50 %, as shown on the $x$ -axis of each panel.

The influence of uncertainties in viscosity on the other terms of the baseflow is presented in figure 5. Solid and dashed lines represent the actual distortion and the linear approximation, respectively. The viscosity distortion was applied in proportion to the original model: $\delta \mu = \sigma \mu (\rho ,T)$ with $\sigma \in [0, 50\,\%]$ . Linear behaviour is evident when $\sigma \leq 10\,\%$ . The influence is measured with the normalised 2-norm, with some results scaled according to the line legend (for better presentation). Here, the normalised 2-norm,

(2.30) \begin{equation} ||\boldsymbol {Q}||_2=\sqrt {\sideset {\frac {1}{N}}{_{i=1}^{N}}\sum Q_{i}^{2}}, \end{equation}

is a measure of the amplitude of the distortions accounting for the overall distortions. One may also choose to exclude the influence of the number of grid points ( $N$ ) by using the $\infty$ -norm, which considers only the maximum distortion over the wall-normal coordinates. Since the presented uncertainty is driven by viscosity, the term $\delta (\partial \mu /\partial T) \approx \sigma (\partial \mu /\partial T)$ , though large, stays within the linear range. The uncertainty in viscosity induces distortions in all components of the baseflow, with considerable magnitudes in $\delta u$ , $\delta w$ and terms related to pressure.

Figure 6. Distortion magnitude of the baseflow components $\log (||\delta \boldsymbol {Q}||_2)$ . The flow in the subcritical, transcritical, supercritical and ideal regimes are compared. In each regime, the uncertainty is driven by $[\mu ]$ , $[\kappa ]$ and [EoS], corresponding to columns 1–3 (with wall heating) and 4–6 (with wall cooling). The driving terms are highlighted with a rectangle of dotted lines.

It is meaningful to compare fluid regimes, wall temperature and sources of uncertainties in figure 6. The distortion magnitude $\log (||\delta \boldsymbol {Q}||_2)$ of the baseflow components is presented as a heat map. The driving uncertainty is prescribed as $\delta \mu =\sigma \mu (\rho ,T)$ , $\delta \kappa =\sigma \kappa (\rho ,T)$ and $\delta p=\sigma p(\rho ,T)$ , with $\sigma =10\,\%$ . Although the figure contains a wealth of data, key knowledge is obtained as follows.

  1. 1. Comparing different fluid regimes, distortions in the ideal fluid are the smallest, followed by the supercritical, subcritical and transcritical cases, indicating an enhancing effect of non-ideal fluid properties on the induced distortion.

  2. 2. The driving uncertainty (see dotted rectangles in each column) leads to distortions in all other components of the baseflow. When the uncertainty comes from the equations of state, the driving terms remain the largest (compared with their induced distortions) across different flow regimes. Both viscosity and thermal conductivity uncertainties give rise to significant distortions in the pressure terms. The influence of thermal conductivity uncertainties is relatively smallest.

  3. 3. Uncertainties in the viscosity lead to the strongest distortion of the primary baseflow profile $[b]$ , regardless of the flow regime and wall temperature. Velocities ( $u$ and $w$ ) are more distorted by any of the studied uncertainties compared with temperature $T$ and density $\rho$ .

  4. 4. Wall heating, compared with wall cooling, does not make an essential difference except for the transcritical case, where wall cooling shows larger distortions than the heating counterpart.

Going through (2.14) and Appendix A, one notes that the fluid properties are directly inputs to the stability operator (Guo et al. Reference Guo, Gao, Jiang and Lee2021; Poulain et al. Reference Poulain, Content, Rigas, Garnier and Sipp2024). The fluid properties, the baseflow and the stability operator form a coupled system as summarised in figure 7. The intrinsic uncertainty influences the system in an integrated manner: the primary baseflow profiles $[b]$ are distorted, which, together with the distorted fluid properties, form inputs to the stability operator. These distortions in the eigenvalue problem are measured by the sensitivity coefficients, which will be discussed next.

Figure 7. Relational diagram depicting the uncertainties in intrinsic fluid properties, the baseflow distortions and their influences on the eigenvalue problem.

3. Results and discussions

3.1. Input distortions and the induced eigenvalue shifts

Since latent uncertainty can lead to distortions of various shapes, we begin by comparing three types of distortions (listed in table 2). First, white noise is assigned to assess the robustness of stability when exposed to a realistic environment, such as a laboratory experiment or a flight test. The second signal corresponds to structural sensitivity, where some inputs are biased due to an inaccurate fluid model. Third, we consider the distortion generating the maximum growth rate shift, which represents the border of sensitivity at a specific amplitude and serves as a measure of sensitivity. Figure 8 provides a profile of $\partial \mu /\partial \rho$ in the pseudo-boiling regime with wall cooling. Unlike an ideal gas, the term $\partial \mu /\partial \rho$ is non-zero. Figure 8 illustrates the shapes of these distortions, with their amplitudes intentionally assigned large for illustration purposes.

Table 2. Distortions of different types.

Figure 8. A portray of different distortions on $\partial \mu /\partial \rho$ .

Figure 9. Illustration of eigenvalue shift validating the sensitivity framework. (a) The eigenvalue trajectory with distorted $\partial \mu /\partial \rho$ ; (b) error of eigenvalue prediction as a function of the departure parameter $\epsilon$ according to (3.1).

The sensitivity coefficients derived in Appendix D have been validated by comparing the distorted eigenvalue using direct estimation (2.13) and the sensitivity framework (2.19). Figure 9(a) illustrates this comparison. The case investigated involves pseudo-boiling with wall cooling. We evaluate the sensitivity of the eigenvalue at $\omega =40$ , $\beta =100$ and $x=1.0$ , corresponding to a typical inviscid instability mode found in a recent investigation (see figure 12 b of Ren & Kloker Reference Ren and Kloker2022a ). The distortion

(3.1) \begin{equation} \delta \boldsymbol {Q}_{\textrm {property}} |_{\textrm {non-zero component}} = -\epsilon \partial \mu /\partial \rho \end{equation}

is applied to $\partial \mu /\partial \rho$ , where $0 \leq \epsilon \leq 1$ is the parameter controlling the amplitude. As $\epsilon$ increases from 0 to 1, the input $\partial \mu /\partial \rho$ is incrementally deformed by an ideal gas model until $\partial \mu /\partial \rho =0$ (fully ideal). As shown in figure 9(a), the sensitivity coefficients provide an accurate prediction in the linear range, and the error (see panel b) becomes noticeable from $\epsilon =0.12$ . We note that this value does however not provide a constant limit criterion for the linear regime, as the distortions can act on other input variables.

3.2. The sensitivity coefficients and a scalar measure

Figure 10. Assembly of sensitivity coefficients in the pseudo-boiling regime with wall cooling. Solid and dashed lines represent the real and imaginary parts, respectively.

For the pseudo-boiling regime with wall cooling, we present profiles of all the sensitivity coefficients in figure 10. Solid lines represent the real parts, while dashed lines represent the imaginary parts. The wave number and frequency values correspond to the case in figure 9. Arrows indicate the location of the Widom line and the boundary layer thickness $\delta _0$ ( $u/u_\infty =0.99$ ). All panels are plotted from the wall to the height of $\delta _0$ . Without exception, each profile shows shears near the Widom line and decays to zero towards the wall and upper boundaries. Distinctive differences between the profiles lie in the amplitude, spanning from $O(10^{-1})$ to $O(10^9)$ . Given such a range of magnitudes on a logarithmic scale, the dominating terms are recognised. For clarity, we have applied a yellow background for values of $O(10^4)$ and higher. Before discussing different thermodynamic regimes, we propose a scalar measure for the degree of sensitivity.

Since the distortion shape can be arbitrary, the focus is placed on the distortion that leads to the maximum shift of the eigenvalue (Giannetti & Luchini Reference Giannetti and Luchini2007). According to (2.19), and accounting for the non-imaginary nature of physical distortions, the inner product reaches its maximum when $\delta \boldsymbol {Q}$ equals the real or imaginary part of $\boldsymbol {S}_{\boldsymbol {Q}}$ . In other words, $\delta \boldsymbol {Q} = \pm \textrm {imag}(\boldsymbol {S}_{\boldsymbol {Q}})$ results in the highest growth rate shift $\alpha _i$ , while $\delta \boldsymbol {Q} = \pm \textrm {real}(\boldsymbol {S}_{\boldsymbol {Q}})$ leads to the maximum shift of the chordwise wavenumber $\alpha _r$ . To avoid ambiguity, we focus on the growth rate and define $\delta \boldsymbol {Q}_{\textrm {max}}$ as $\pm \textrm {imag}(\boldsymbol {S}_{\boldsymbol {Q}})$ .

To better understand the eigenvalue shift induced by $\delta \boldsymbol {Q}_{\textrm {max}}$ relative to $\delta \boldsymbol {Q}_{\textrm {noise}}$ of the same amplitude, we present the comparison in figure 11. As seen in panel (a), the distortion

(3.2) \begin{equation} \delta \boldsymbol {Q} = \textrm {real}(\boldsymbol {S}_{\boldsymbol {Q}} \exp {i\theta }), \end{equation}

with $0\leq \theta \lt 2\pi$ , leads to an elliptic trajectory of the eigenvalue. The eigenvalues corresponding to $\delta \boldsymbol {Q}$ ( $\theta =0,\pi /2,\pi ,3\pi /2$ ) are highlighted with red dots. When the coherence of the distortions is lost and replaced with random noise, the eigenvalue shifts are significantly narrowed, as shown in panel (b), which corresponds to the central rectangle highlighted in panel (a). In both cases, we have maintained an amplitude of

(3.3) \begin{equation} \left \Vert \delta \boldsymbol {Q}\right \Vert _{2}=10^{-3}/\sqrt {N},\,\mathrm {with}\,N=200. \end{equation}

In contrast to panel (a), where the eigenvalue shift is ‘optimised’ along an elliptical trajectory, the introduction of random distortion causes the eigenvalue to remain closer to its undistorted value, with the displacement appearing random as well. This difference aligns with the sensitivity coefficient defined in (2.19). Consequently, within the linear sensitivity regime, a random distortion is less likely to cause significant deviations in the eigenvalue.

Figure 11. (a) The eigenvalue shift due to $\delta \boldsymbol {Q}_{\textrm {max}}$ (red dots), the term distorted is $\partial \mu /\partial \rho$ ; (b) similar to panel (a), but with $\delta \boldsymbol {Q}_{\textrm {noise}}$ at the same amplitude (measured with 2-norm).

So far, we have examined the eigenvalue shift at a fixed chordwise position (corresponding to $x=1$ ). To account for the integral effects and gain an initial understanding of the sensitivity to different inputs, we plot the N factor (integral of the growth rate) as a function of the streamwise coordinate in figure 12. Twenty-three terms corresponding to (2.14) have been calculated. The distortion corresponds to $\delta \boldsymbol {Q}_{\textrm {max}}$ with a 2-norm of $10^{-3}/\sqrt {200}$ . Recalling the inviscid nature of the new mode (connected to the main flow), one notices that the growth rate is significantly larger than that of conventional (also largely inviscid) cross-flow modes (connected to the true cross-flow $w_s$ ) and reaches an N factor of 33.3 at $x=5$ . Temperature distortions lead to a dramatic N factor change, followed by density and $\partial p/\partial \rho$ , which agrees with the local amplitude seen in figure 10. We have tabulated the top seven terms in table 3. If the temperature and density are correctly obtained in the baseflow stage, ensuring an accurate EOS is essential for a rational transition prediction.

Figure 12. N factor influenced by $\delta \boldsymbol {Q}_{\textrm {max}}$ . The amplitude of an individual distortion is $\Vert \delta \boldsymbol {Q}\Vert _{2}=10^{-3}/\sqrt {200}$ . The red and green curves show positive and negative shifts, respectively. The numbers $1\cdots 23$ correspond to input terms (see (2.14)) sorted by the impact on the N factor in descending order (see table 3). The instability is in the pseudo-boiling regime with dominating inviscid instability.

The result has also been tested as a weak function of the spanwise wavenumber $\beta$ , as shown in Appendix E. Based on the above discussions, we propose a scalar measure of sensitivity:

(3.4) \begin{equation} M=\frac {\left |\delta \alpha _{i}\right |}{\left \Vert \delta \boldsymbol {Q}_{\mathrm {max}}\right \Vert _2 }. \end{equation}

This measure demonstrates the maximum possible response of the growth rate to distortions scaled by its normalised norm, thereby indicating the sensitivity to certain inputs. Note that the choice of the norm is not exclusive. See also the discussion for baseflow distortions (2.30). The normalised norm provides an intuitive indicator of the sensitivity amplitude which does not influence the sensitivity coefficients and the induced eigenvalue distortions. Instead of concentrating on the shape of the sensitivity coefficients, $M$ helps to compare across different cases and regimes.

Table 3. $\delta N/N$ at four uniformly distributed observation points.

3.3. Sensitivity crossing different thermodynamic regimes

Figure 13. A comparison of sensitivity across different regimes and between different input variables. The numbers on the heatmap specify the log value of the sensitivity measure: $\log (M)$ .

The non-ideal gas properties primarily depend on the thermodynamic regimes of the fluid. We compare the stability diagram and the sensitivity measure $M$ in figure 13 for regimes introduced in table 1, grouped into wall heating and cooling. To allow for comparison across different regimes and without loss of generality, the parameters ( $x=1.0$ , $\omega =15$ , $\beta =80$ ) have been chosen for the investigations below (see figure 2 of Ren & Kloker Reference Ren and Kloker2022b ). According their panels (b) and (c), the neutral curves and growth rates are only quantitatively different, except for the transcritical regime with wall cooling, where the instability is dominated by a new inviscid mode. In figure 13, the term ‘distorted’ is listed on the figure’s left and sorted according to the sensitivity measure in the transcritical regime. Considering a representative growth rate of $O(1)$ , a sensitivity measure of $M =1$ ( $\log M =0$ ) will give $\delta \alpha _i=0.001$ at a distortion amplitude of $10^{-3}$ , which starts to become visible. Note that the shape of the distortion was assumed to be the best case (recall figure 11). With this estimation, the sensitivity to distortions with negative values of $\log (M)$ is minimal.

Going through the sensitivity measure of all eight cases presented in figure 13 implies an analogous ranking of various terms, showing that different flow regimes only moderately influence the relative position of each input of the stability operator. The temperature and density profiles are the most consequential. Subsequently, the term $\partial ^2p/\partial \rho ^2$ stands out only for the pseudo-boiling case with wall cooling. For the ideal regime, baseflow profiles remain more influential than other thermodynamic and transport properties. Some terms are zero (e.g. $\partial \mu /\partial \rho$ , $\partial \kappa /\partial \rho$ ) per the ideal assumption, and random or structural distortions with moderate amplitude will not invoke a significant shift in the growth rate, indicating their less active roles in the system. We notice that the transport properties ( $[\mu ]$ and $[\kappa ]$ profiles) are less influential for all the regimes. This implies that the ideal case will be robust provided the six key terms, $T,\rho ,u,w,\partial p/\partial T, \partial p/\partial \rho$ , are not distorted. In contrast, the pseudo-boiling regimes are sensitive to more terms: $\partial ^2p/\partial \rho ^2, \partial ^2p/\partial \rho \partial T$ for wall heating and, in addition, $\partial \mu /\partial \rho , \mu , \partial ^2\mu /\partial \rho ^2, \kappa$ for wall cooling.

3.4. Magnitude range of sensitivities to various inputs

The sensitivity profiles (figure 10) and their measure (figure 13) demonstrate a magnitude range of approximately $O(10^{10})$ . The reason is explored here by analysing different sensitivity profiles and their mathematical constituents. We examine the pseudo-boiling case with wall cooling. From (D4), the dominating terms are written as

(3.5) \begin{equation} \begin{aligned} S_{T}= & \underset {S_{T}\left (1\right )}{\underbrace {\frac {\partial }{\partial y}\left (\frac {\partial ^{2}p}{\partial \rho \partial T}\hat {\rho }^{*}\hat {v}^{\dagger }+\frac {\partial ^{2}p}{\partial T^{2}}\hat {T}^{*}\hat {v}^{\dagger }\right )}}+\textrm {visc.}\\ = & \underset {S_{T}\left (1a\right )}{\underbrace {\frac {\partial }{\partial y}\left (\frac {\partial ^{2}p}{\partial \rho \partial T}\hat {\rho }^{*}\hat {v}^{\dagger }\right )}}+\underset {S_{T}\left (1b\right )}{\underbrace {\frac {\partial }{\partial y}\left (\frac {\partial ^{2}p}{\partial T^{2}}\hat {T}^{*}\hat {v}^{\dagger }\right )}}+\textrm {visc.}\\ = & \underset {S_{T}\left (1a_{1}\right )}{\underbrace {\frac {\partial }{\partial y}\left (\frac {\partial ^{2}p}{\partial \rho \partial T}\right )\hat {\rho }^{*}\hat {v}^{\dagger }}}+\underset {S_{T}\left (1a_{2}\right )}{\underbrace {\frac {\partial ^{2}p}{\partial \rho \partial T}\frac {\partial \hat {\rho }^{*}}{\partial y}\hat {v}^{\dagger }}}+\underset {S_{T}\left (1a_{3}\right )}{\underbrace {\frac {\partial ^{2}p}{\partial \rho \partial T}\hat {\rho }^{*}\frac {\partial \hat {v}^{\dagger }}{\partial y}}}+S_T(1b)+\textrm {visc.}, \end{aligned} \end{equation}

with ‘visc.’ standing for viscous terms. The distribution of key terms are plotted in figure 14. As inferred from figure 14(ac), the terms $S_{T}(1)$ significantly dominate the sum of the rest terms, which are viscous (scaled by $Re=1.4687\times 10^5$ ). Among $S_{T}(1)$ , $S_{T}(1a)$ is much larger than $S_{T}(1b)$ . We plot in panels (d)–(f) the physical terms related to $S_{T}(1)$ . As can be seen, the gradients of the terms $\hat {\rho }^*$ , $\hat {v}^\dagger$ and $\partial ^2 p/\partial \rho \partial T$ synchronise around the Widom line, leading to the large amplitude of $S_T$ .

Figure 14. Illustration of key terms for the sensitivity profile $S_T$ corresponding to (3.5): (a) $S_T$ and $S_T(1)$ ; (b) $S_T(1a)$ and $S_T(1b)$ ; (c) $S_T(1a_1)$ , $S_T(1a_2)$ and $S_T(1a_3)$ ; (df) the components composing $S_T(1)$ .

Figure 15. Sensitivity to distortions of $\rho$ , $u$ and $w$ (panels ac). We show dominating terms (the real parts) of each profile.

Are viscous terms always negligible? The composition of $S_\rho$ , $S_u$ and $S_w$ is plotted in figure 15. Similar to $S_T$ , $S_\rho$ is dominated by an inviscid term

(3.6) \begin{equation} S_{\rho }= \underset {S_{\rho }\left (1\right )}{\underbrace {\frac {\partial }{\partial y}\left (\frac {\partial ^{2}p}{\partial \rho ^{2}}\hat {\rho }^{*}\hat {v}^{\dagger }+\frac {\partial ^{2}p}{\partial \rho \partial T}\hat {T}^{*}\hat {v}^{\dagger }\right )}}+\cdots \end{equation}

Compared with (D5), the rest of the terms include viscous and inviscid terms. However, a single dominating term does not exist for $S_u$ and $S_w$ . We re-write these two profiles in (3.7) and (3.8). As seen in figure 15, the viscous terms influence the profile in the leading order (only $S_u(3)$ and $S_w(3)$ can be ignored), showing the significance of viscous terms though the Reynolds number scales them.

(3.7) \begin{equation} \begin{aligned} S_{u}= & \underset {S_{u}\left (1\right )}{\underbrace {\frac {\partial }{\partial y}\left (\rho \hat {v}^{*}\hat {u}^{\dagger }\right )}}+\underset {S_{u}\left (2\right )}{\underbrace {i\alpha ^{*}\left (\hat {\rho }^{*}\hat {\rho }^{\dagger }+\rho \hat {u}^{*}\hat {u}^{\dagger }+\rho \hat {v}^{*}\hat {v}^{\dagger }+\rho \hat {w}^{*}\hat {w}^{\dagger }+\rho \frac {\partial e}{\partial \rho }\hat {\rho }^{*}\hat {T}^{\dagger }+\rho \frac {\partial e}{\partial T}\hat {T}^{*}\hat {T}^{\dagger }\right )}}\\ & +\underset {S_{u}\left (3\right )}{\underbrace {\frac {1}{Re}\frac {\partial }{\partial y}\left [2\mu \left (i\alpha ^{*}\hat {v}^{*}-\frac {\partial \hat {u}^{*}}{\partial y}\right )\hat {T}^{\dagger }+\frac {\partial \mu }{\partial T}\hat {T}^{*}\left (\frac {\partial \hat {u}^{\dagger }}{\partial y}-2\frac {\partial u}{\partial y}\hat {T}^{\dagger }+i\alpha ^{*}\hat {v}^{\dagger }\right )\right ]}}\\ & +\underset {S_{u}\left (4\right )}{\underbrace {\frac {1}{Re}\frac {\partial }{\partial y}\left [\frac {\partial \mu }{\partial \rho }\hat {\rho }^{*}\left (\frac {\partial \hat {u}^{\dagger }}{\partial y}-2\frac {\partial u}{\partial y}\hat {T}^{\dagger }+i\alpha ^{*}\hat {v}^{\dagger }\right )\right ]}}, \end{aligned} \end{equation}
(3.8) \begin{equation} \begin{aligned} &S_{w}= \underset {S_{w}\left (1\right )}{\underbrace {\frac {\partial }{\partial y}\left (\rho \hat {v}^{*}\hat {w}^{\dagger }\right )}}+\underset {S_{w}\left (2\right )}{\underbrace {i\beta \left (\hat {\rho }^{*}\hat {\rho }^{\dagger }+\rho \hat {u}^{*}\hat {u}^{\dagger }+\rho \hat {v}^{*}\hat {v}^{\dagger }+\rho \hat {w}^{*}\hat {w}^{\dagger }+\rho \frac {\partial e}{\partial \rho }\hat {\rho }^{*}\hat {T}^{\dagger }+\rho \frac {\partial e}{\partial T}\hat {T}^{*}\hat {T}^{\dagger }\right )}}\\ & +\underset {S_{w}\left (3\right )}{\underbrace {\frac {i\beta }{Re}\frac {\partial }{\partial y}\left (2\mu \hat {v}^{*}\hat {T}^{\dagger }+\frac {\partial \mu }{\partial \rho }\hat {\rho }^{*}\hat {v}^{\dagger }+\frac {\partial \mu }{\partial T}\hat {T}^{*}\hat {v}^{\dagger }\right )}}\\ & +\underset {S_{w}\left (4\right )}{\underbrace {\frac {1}{Re}\frac {\partial }{\partial y}\left [\frac {\partial \mu }{\partial T}\hat {T}^{*}D\hat {w}^{\dagger }+\frac {\partial \mu }{\partial \rho }\hat {\rho }^{*}D\hat {w}^{\dagger }-2\left (\mu D\hat {w}^{*}\hat {T}^{\dagger }+\frac {\partial \mu }{\partial \rho }\frac {\partial w}{\partial y}\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial \mu }{\partial T}\frac {\partial w}{\partial y}\hat {T}^{*}\hat {T}^{\dagger }\right )\right ]}} \end{aligned}. \end{equation}

Take a column-wise comparison on figure 13, the sensitivity to $\partial ^2 p/\partial \rho ^2$ stands out in the pseudo-boiling regime with wall cooling. Recall (D7), the sensitvity coefficient amounts to a single term $-{\partial \rho }/{\partial y}(\hat {\rho }^{*}\hat {v}^{\dagger })$ . In figure 16, we compare the real parts of these terms across different regimes. From panel (a), we observe that the terms in the pseudo-boiling case is four orders of magnitude larger than in other regimes, in accordance with figure 13. Panels (b)–(d) unveil that the differences are due to $\partial \rho /\partial y$ and $\hat {\rho }$ . Both terms are largely two orders larger than the rest of the cases. In fact, according to the continuity equation,

(3.9)

the term $D\rho$ is balanced by terms $\rho$ and $\hat {\rho }$ .

As seen from figures 16(e)–16(h), all three terms are in the same order of magnitude regardless of the thermodynamic regime. Simultaneously, juxtaposing the four panels reveals that the pseudo-boiling regime has a significantly larger amplitude of the three balancing terms according to (3.8), clarifying the differences observed in panel (a). The non-ideal gas behaviour constitutes a quadratic effect, magnifying the influence of $\partial \rho /\partial y$ on sensitivities into its square (through the product with term $\hat {\rho }$ ). At the same time, we notice that the thermodynamic derivatives (e.g. $\partial ^2p/(\partial \rho \partial T)$ and $\partial ^2p/\partial T^2$ ) and their wall-normal gradients are significantly larger in the pseudo-boiling regime (see figure 14 f). Multiplication with $\hat {\rho }$ leads to tremendous sensitivity to temperature and density distortions.

Figure 16. (a) Sensitivity (the real part) to $\partial ^2 p/\partial \rho ^2$ with wall cooling. The sensitivity equals the derivatives of the products of three terms shown in panels (b)–(d). We illustrate the balance of term $\hat {\rho }$ , $\rho$ and $D{\rho }$ according to the continuity equation (3.9) in panels (e)–(h) for the four regimes considered (with wall cooling).

As another example, the sensitivity to $\partial \mu /\partial \rho$ is studied according to (D9). The sensitivity coefficients contain 15 terms, which are presented in figure 17, comparing across different thermodynamic regimes. Since the amplitude for the pseudo-boiling case is more prominent in the wall-cooling case, we have scaled their amplitudes with a factor of 10, as indicated on the axis of panel (b). With wall cooling, we find terms 1 and 10 stand out in the transcritical case, say, $({\partial u}/{\partial y})D\hat {\rho }^{*}\hat {u}^{\dagger }$ and $ ({\partial w}/{\partial y})D\hat {\rho }^{*}\hat {w}^{\dagger }$ , both terms are led by $D\hat {\rho }^{*}$ . This again attributes to the variation of the density near the Widom line. We also noticed that some terms 4–9 and 12, 14 and 15 remain small for all the cases.

Figure 17. Exposition of sensitivity terms for $\partial \mu /\partial \rho$ with (a) wall heating and (b) cooling. The terms corresponding to (D9) have been labelled as 1–15 on the plot. The height of each pillar stands for the 2-norm of a term.

Taking an overview of the sensitivity coefficients, it is unambiguous that the sensitivity can be uniformly written as

(3.10)

where $q$ , $\hat {q}^*$ and $\hat {q}^\dagger$ stand for the baseflow, the eigenfunction (complex conjugate) and the adjoint eigenfunction inclusive of their spatial gradients. The sensitivity is fully viscous for the viscosity and thermal conductivity, while mixed for the baseflow and EOS profiles, inducing a significant decrease of sensitivity to viscosity and thermal conductivity terms shown in figure 13.

This subsection examined the sensitivity both term wise and case wise. The impact of all inputs differs significantly, with a magnitude range of approximately $O(10^{10})$ . In the pseudo-boiling regime, the wall-normal gradient of baseflow profiles (e.g. $\rho$ , $\partial ^2p/\partial \rho \partial T$ and $\partial ^2p/\partial T^2$ ) and the density perturbation $\hat {\rho }$ become large simultaneously following the balancing relation of the continuity equation. The overall sensitivity thus increases with a boost due to a quadratic effect in certain terms (e.g. sensitivity to $\partial ^2p/\partial \rho ^2$ ) that remain small in other regimes. Special care must be taken in modelling, experiments or applications, as a slight deviation may lead to a significant shift in the transition location.

3.5. Fluid model relaxation

Based on the above discussions and towards an effective engineering prediction, one naturally raises the question of how accurate it is to use an empirical law for specific inputs in different thermodynamic regimes. Considering the dependence of the baseflow and stability operator on fluid properties (recall figure 7), the outcome depends on two factors: the departure of the input from its actual value and the eigenvalue sensitivity to it. The error is thus proportional to the non-ideality of the fluid regime. For example, if one uses ideal gas laws for all the fluid properties, the non-ideal gas behaviour will be totally ignored, and the difference will be significant in the highly non-ideal regime. An extreme example is in the transcritical regime for wall cooling, where the calculated mode will be a cross-flow mode rather than the 2-D inviscid mode, whose growth rate is much larger.

Figure 18. Comparison of N factors using non-ideal and ideal models (for [EoS], $[\mu ]$ and $[\kappa ]$ ). The thermodynamic regime (gas-like, liquid-like, pseudo-boiling) has been signified for relevant panels. All cases are subject to $\beta =80$ , $\omega =0$ (steady CF mode), except the last panel with $\beta =80$ , $\omega =20$ (inviscid TS mode).

However, if the baseflow has been correctly solved numerically or measured in experiments, while in the stage of stability analysis, the inputs related to the EOS or transport properties are supplied with ideal-gas models, how reliable are the results? This is the classic issue of using a ‘standard’ simplified stability solver on complex-physics baseflow profiles. This also means that no distortion amplitudes need to be specified as before. We supplement the discussion by analysing the N-factors based on stability calculations. The biased inputs are subject to ideal models for EOS ( $p=\rho RT$ ), viscosity (Sutherland’s law) and thermal conductivity (Sutherland’s law), as shown in figure 18. Recall the sensitivity measure presented in figure 13. Provided the primary baseflow profiles (group $[b]$ in (2.14)) are correct, in the gas-like regime, the sensitivity is only noteworthy for terms $\partial p/\partial T$ and $\partial p/\partial \rho$ . Both inputs can be satisfactorily characterised using the ideal gas equation-of-state, not leading to a visible difference in the N-factor. Furthermore, the liquid-like regime shows concern for $[\mu ]$ terms, while the pseudo-boiling regime requires that [EoS] and $[\mu ]$ terms be adequately accounted for. In particular, the inviscid TS mode additionally requires $[\kappa ]$ terms, as the N-factor and the absolute error are much more considerable.

Table 4 summarises the idealisable fluid models for stability analysis, indicating, for each fluid regime, which terms can be safely taken from the simple ideal model when correct baseflow profiles are supplied. This means that the induced distortions and the corresponding sensitivity are weak enough not to induce a notable change in the N factor. Table 4 provides a first estimation of standard stability solvers in predicting flow transition in different fluid regimes, where only the supercritical regime (gas-like), approximately similar to an ideal gas, is applicable to ideal gas models. This assessment is qualitative and assumes a low-Mach-number configuration; the error can increase or decrease when the flow temperature range approaches or recedes from the Widom line, where viscous heating effects may become significant at higher Mach numbers.

Table 4. A summary of idealisable fluid models for linear stability analysis, provided correct baseflow profiles ( $u$ , $w$ , $\rho$ , $T$ ) are supplied.

4. Conclusions

The investigation examines the sensitivity of three-dimensional boundary layer modal instability to intrinsic uncertainties in fluid properties, including the equation-of-state (EoS) and transport properties. We adopt a representative fluid, supercritical CO2 at 80 bar (static pressure), which exhibits significant non-ideal gas behaviour crossing the Widom line. The flow matches the redesigned DLR experiment with identical pressure coefficient distribution (strong favourable pressure gradient), Reynolds number and a low Mach number. Distinguishable thermodynamic regimes – liquid-like, pseudo-boiling, gas-like and ideal gas – emerge by prescribing respective temperature boundary conditions. The mechanisms for linear instability found recently (Ren & Kloker Reference Ren and Kloker2022a ) demonstrate a changeover of the dominant mode from conventional cross-flow (CF) instability to inviscid Tollmien–Schlichting (TS) type in the pseudo-boiling regime with wall cooling, despite the strong favourable pressure gradient.

The sensitivity of linear stability is a response of a coupled system, including the acquisition of the laminar baseflow and the eigenvalue problem. Uncertainties in any one or combination of terms of the EoS, viscosity and thermal conductivity give rise to distortions in primary baseflow profiles ( $u$ , $w$ , $\rho$ , $T$ ) and all profiles related to fluid properties needed for the stability analysis. The shape and amplitude of these distortions are discussed and found to be dependent on the source of uncertainty and fluid regimes.

Further, these distortions are measured by the sensitivity coefficients, whose product indicates the alterations of the eigenvalue. The sensitivity coefficients are theoretically formulated for each input of the stability operator using the adjoint equations, indicating the eigenvalue shift when performing the inner product with corresponding physical distortions. The adjoint equations have been validated using the bi-orthogonal relationship, and all sensitivity coefficients have been validated against linear stability calculations with distorted inputs. Considering that an actual distortion can be arbitrary (e.g. random noise) or structural (e.g. model deviation), we propose a unified measure for the degree of sensitivity that describes the maximum possible eigenvalue shift. It is worth noting that the adjoint-equation-based linear sensitivity provides details to understand the system’s response to the bias of specific inputs, thereby directing the most effective control strategies. The method also efficiently predicts growth rate shifts due to different distortions. However, one may circumvent the details and obtain a stability prediction directly with distorted inputs, which is less effective. However, the amplitude of the distortion can be considerable in this case.

The sensitivity has been investigated for all inputs of the stability operator and compared horizontally across different thermodynamic regimes. A range of its amplitude spanning the order of $10^{10}$ for different inputs is found and attributed to the scaling of the Reynolds number for viscosity-related terms and different wall-normal gradients of the baseflow and perturbation profiles. The sensitivity is significantly larger (by 1–2 orders of magnitude) in the pseudo-boiling regime. In particular, a quadratic effect has been uncovered. The boost of $\partial \rho /\partial y$ results in an equivalent rise for $\hat {\rho }$ through a balancing relation of the continuity equation. The product of $\partial \rho /\partial y$ and $\hat {\rho }$ gives rise to a quadratic sensitivity increase to $\partial ^2 p/\partial \rho ^2$ across the Widom line. Meanwhile, the tremendous sensitivity to temperature and density terms – especially regarding the respective direct baseflow profiles – is also due to large thermodynamic derivatives (e.g. $\partial ^2 p/(\partial \rho \partial T)$ and $\partial ^2 p/\partial T^2$ ) and their wall-normal gradients.

Towards practical applications for transition prediction, the N-factor deviation has been investigated when certain groups of terms are given following an ideal gas assumption, based on the use of the correct, non-ideal baseflow profiles. In other words, a standard instability solver is applied to the fully non-ideal baseflow. In particular, we show that in the gas-like regime, using the idealised equation-of-state (EOS), viscosity and thermal conductivity do not provoke noticeable differences. The liquid-like regime requires viscosity to be correctly modelled, while the pseudo-boiling regime needs all elements to be correct except for the thermal conductivity for the CF mode. Accurately predicting the inviscid TS mode’s growth rate requires proper modelling of all the above elements – employing an ideal-fluid solver to the correct baseflow profiles gives a meaningless instability result.

The present research focused on the sensitivity of flow stability to intrinsic uncertainties in fluid properties. This question arises from the ever-growing industrial demand and the ongoing upgrades and improvements of fluid property databases. The different sensitivity behaviours found in this research indicate that an accurate model describing the fluid is essential in highly non-ideal regimes, particularly for the accurate description of the new dominating inviscid mode.

Acknowledgements

J.R. acknowledges helpful discussions with Drs. J. Park (Coventry University), Y. Xi, S. Fu (Tsinghua University), P.C. Boldini and R. Pecnik (TU Delft).

Funding

This work was supported by National Natural Science Foundation of China (Grant Nos.12372215 and 92152109), the Aeronautical Science Fund (Grant No. 2024M005072001) and the Alexander von Humboldt foundation.

Declaration of interests

The authors report no conflict of interest.

Data availability statement

The data that support the findings of this study are available from the first and corresponding author upon reasonable request.

Appendix A. The stability operator $\boldsymbol {\mathsf {L}}$

The specific expressions of the linear stability operator introduced in (2.14) are as follows.

(A1) \begin{equation} \boldsymbol {\mathsf {L}}_{1,1}=-i\alpha u+i\omega -i\beta w,\; \boldsymbol {\mathsf {L}}_{1,2}=-i\alpha \rho ,\; \boldsymbol {\mathsf {L}}_{1,3}=-\rho D-\frac {\partial \rho }{\partial y},\; \boldsymbol {\mathsf {L}}_{1,4}=-i\beta \rho . \end{equation}
(A2) \begin{equation} \boldsymbol {\mathsf {L}}_{2,1}=-i\alpha \frac {\partial p}{\partial \rho }+\frac {1}{Re}\frac {\partial \mu }{\partial \rho }\frac {\partial u}{\partial y}D+\frac {1}{Re}\frac {\partial \mu }{\partial \rho }\frac {\partial ^{2}u}{\partial y^{2}}+\frac {1}{Re}\frac {\partial u}{\partial y}\left (\frac {\partial ^{2}\mu }{\partial \rho ^{2}}\frac {\partial \rho }{\partial y}+\frac {\partial ^{2}\mu }{\partial \rho \partial T}\frac {\partial T}{\partial y}\right ). \end{equation}
(A3) \begin{equation} \boldsymbol {\mathsf {L}}_{2,2}=-\frac {\alpha ^{2}\left (2\mu +\lambda \right )}{Re}-i\alpha \rho u+i\omega \rho +\frac {1}{Re}\frac {\partial \mu }{\partial y}D-i\beta \rho w+\frac {\mu }{Re}D^{2}-\frac {\beta ^{2}\mu }{Re}. \end{equation}
(A4) \begin{equation} \boldsymbol {\mathsf {L}}_{2,3}=\frac {i\alpha \left (\mu +\lambda \right )D}{Re}+\frac {i\alpha }{Re}\frac {\partial \mu }{\partial y}-\rho \frac {\partial u}{\partial y},\; \boldsymbol {\mathsf {L}}_{2,4}=-\frac {\alpha \beta \left (\mu +\lambda \right )}{Re}. \end{equation}
(A5) \begin{equation} \boldsymbol {\mathsf {L}}_{2,5}=-i\alpha \frac {\partial p}{\partial T}+\frac {1}{Re}\frac {\partial \mu }{\partial T}\frac {\partial u}{\partial y}D+\frac {1}{Re}\frac {\partial \mu }{\partial T}\frac {\partial ^{2}u}{\partial y^{2}}+\frac {1}{Re}\frac {\partial u}{\partial y}\left (\frac {\partial ^{2}\mu }{\partial T^{2}}\frac {\partial T}{\partial y}+\frac {\partial ^{2}\mu }{\partial T\partial \rho }\frac {\partial \rho }{\partial y}\right ). \end{equation}
(A6) \begin{equation} \boldsymbol {\mathsf {L}}_{3,1}=\frac {i\alpha }{Re}\frac {\partial \mu }{\partial \rho }\frac {\partial u}{\partial y}-\frac {\partial p}{\partial \rho }D+i\beta \frac {1}{Re}\frac {\partial \mu }{\partial \rho }\frac {\partial w}{\partial y}-\frac {\partial ^{2}p}{\partial \rho ^{2}}\frac {\partial \rho }{\partial y}-\frac {\partial ^{2}p}{\partial \rho \partial T}\frac {\partial T}{\partial y}. \end{equation}
(A7) \begin{equation} \boldsymbol {\mathsf {L}}_{3,2}=\frac {i\alpha \left (\mu +\lambda \right )D}{Re}+\frac {i\alpha }{Re}\frac {\partial \lambda }{\partial y}. \end{equation}
(A8) \begin{equation} \boldsymbol {\mathsf {L}}_{3,3}=-\frac {\alpha ^{2}\mu }{Re}-i\alpha \rho u+i\omega \rho +\frac {1}{Re}\frac {\partial \left (2\mu +\lambda \right )}{\partial y}D-i\beta \rho w+\frac {2\mu +\lambda }{Re}D^{2}-\frac {\beta ^{2}\mu }{Re}. \end{equation}
(A9) \begin{equation} \boldsymbol {\mathsf {L}}_{3,4}=i\beta \frac {1}{Re}\frac {\partial \lambda }{\partial y}+\frac {i\beta \left (\mu +\lambda \right )D}{Re}. \end{equation}
(A10) \begin{equation} \boldsymbol {\mathsf {L}}_{3,5}=\frac {i\alpha }{Re}\frac {\partial \mu }{\partial T}\frac {\partial u}{\partial y}-\frac {\partial p}{\partial T}D+i\beta \frac {1}{Re}\frac {\partial \mu }{\partial T}\frac {\partial w}{\partial y}-\frac {\partial ^{2}p}{\partial T^{2}}\frac {\partial T}{\partial y}-\frac {\partial ^{2}p}{\partial \rho \partial T}\frac {\partial \rho }{\partial y}. \end{equation}
(A11) \begin{equation} \boldsymbol {\mathsf {L}}_{4,1}=\frac {1}{Re}\frac {\partial \mu }{\partial \rho }\frac {\partial w}{\partial y}D-i\beta \frac {\partial p}{\partial \rho }+\frac {1}{Re}\frac {\partial \mu }{\partial \rho }\frac {\partial ^{2}w}{\partial y^{2}}+\frac {1}{Re}\frac {\partial w}{\partial y}\left (\frac {\partial ^{2}\mu }{\partial \rho ^{2}}\frac {\partial \rho }{\partial y}+\frac {\partial ^{2}\mu }{\partial \rho \partial T}\frac {\partial T}{\partial y}\right ). \end{equation}
(A12) \begin{equation} \boldsymbol {\mathsf {L}}_{4,2}=-\frac {\alpha \beta \left (\mu +\lambda \right )}{Re},\; \boldsymbol {\mathsf {L}}_{4,3}=i\beta \frac {1}{Re}\frac {\partial \mu }{\partial y}+\frac {i\beta \left (\mu +\lambda \right )D}{Re}-\rho \frac {\partial w}{\partial y}. \end{equation}
(A13) \begin{equation} \boldsymbol {\mathsf {L}}_{4,4}=-\frac {\alpha ^{2}\mu }{Re}-i\alpha \rho u+i\omega \rho +\frac {1}{Re}\frac {\partial \mu }{\partial y}D-i\beta \rho w+\frac {\mu }{Re}D^{2}-\frac {\beta ^{2}\left (2\mu +\lambda \right )}{Re}. \end{equation}
(A14) \begin{equation} \boldsymbol {\mathsf {L}}_{4,5}=\frac {1}{Re}\frac {\partial \mu }{\partial T}\frac {\partial w}{\partial y}D-i\beta \frac {\partial p}{\partial T}+\frac {1}{Re}\frac {\partial \mu }{\partial T}\frac {\partial ^{2}w}{\partial y^{2}}+\frac {1}{Re}\frac {\partial w}{\partial y}\left (\frac {\partial ^{2}\mu }{\partial T^{2}}\frac {\partial T}{\partial y}+\frac {\partial ^{2}\mu }{\partial T\partial \rho }\frac {\partial \rho }{\partial y}\right ). \end{equation}
(A15) \begin{equation} \begin{aligned}\boldsymbol {\mathsf {L}}_{5,1}= & -i\alpha \rho u\frac {\partial e}{\partial \rho }+i\omega \rho \frac {\partial e}{\partial \rho }+\frac {1}{RePrEc}\frac {\partial \kappa }{\partial \rho }\frac {\partial T}{\partial y}D-i\beta \rho w\frac {\partial e}{\partial \rho }\\ & +\frac {1}{RePrEc}\left (\frac {\partial \kappa }{\partial \rho }\frac {\partial ^{2}T}{\partial y^{2}}+\frac {\partial ^{2}\kappa }{\partial \rho ^{2}}\frac {\partial \rho }{\partial y}\frac {\partial T}{\partial y}+\frac {\partial ^{2}\kappa }{\partial \rho \partial T}\left (\frac {\partial T}{\partial y}\right )^{2}\right )\\ & +\frac {1}{Re}\frac {\partial \mu }{\partial \rho }\left (\left (\frac {\partial u}{\partial y}\right )^{2}+\left (\frac {\partial w}{\partial y}\right )^{2}\right ). \end{aligned} \end{equation}
(A16) \begin{equation} \boldsymbol {\mathsf {L}}_{5,2}=-i\alpha p+\frac {2\mu }{Re}\frac {\partial u}{\partial y}D. \end{equation}
(A17) \begin{equation} \boldsymbol {\mathsf {L}}_{5,3}=\frac {2i\alpha \mu }{Re}\frac {\partial u}{\partial y}-pD+i\beta \frac {2\mu }{Re}\frac {\partial w}{\partial y}-\rho \frac {\partial e}{\partial y}. \end{equation}
(A18) \begin{equation} \boldsymbol {\mathsf {L}}_{5,4}=\frac {2\mu }{Re}\frac {\partial w}{\partial y}D-i\beta p. \end{equation}
(A19) \begin{equation} \begin{aligned} \boldsymbol {\mathsf {L}}_{5,5}= &-\frac {\alpha ^{2}\kappa }{RePrEc}-i\alpha \rho u\frac {\partial e}{\partial T}+i\omega \rho \frac {\partial e}{\partial T}+\frac {1}{RePrEc}\left (\frac {\partial \kappa }{\partial y}+\frac {\partial \kappa }{\partial T}\frac {\partial T}{\partial y}\right )D\\ &-i\beta \rho w\frac {\partial e}{\partial T}+\frac {\kappa }{RePrEc}D^{2}-\frac {\kappa \beta ^{2}}{RePrEc}\\ &+\frac {1}{RePrEc}\left (\frac {\partial \kappa }{\partial T}\frac {\partial ^{2}T}{\partial y^{2}}+\frac {\partial ^{2}\kappa }{\partial T^{2}}\left (\frac {\partial T}{\partial y}\right )^{2}+\frac {\partial ^{2}\kappa }{\partial \rho \partial T}\frac {\partial T}{\partial y}\frac {\partial \rho }{\partial y}\right )\\ &+\frac {1}{Re}\frac {\partial \mu }{\partial T}\left (\left (\frac {\partial u}{\partial y}\right )^{2}+\left (\frac {\partial w}{\partial y}\right )^{2}\right ). \end{aligned} \end{equation}

Appendix B. Derivation of the adjoint equations

We re-write the stability operator based on the order of wall-normal derivative:

(B1) \begin{equation} \boldsymbol {\mathsf {L}}=\boldsymbol {\mathsf {L}}_{0}+\boldsymbol {\mathsf {L}}_{1}D+\boldsymbol {\mathsf {L}}_{2}D^{2}. \end{equation}

Upon integration by parts,

(B2) \begin{equation} \left \langle \hat {\boldsymbol {q}}^{\dagger },\boldsymbol {\mathsf {L}}\hat {\boldsymbol {q}}\right \rangle = \left \langle \left (\boldsymbol {\mathsf {L}}_{0}^{H}-\boldsymbol {\mathsf {L}}_{1}^{H}D-\frac {\partial \boldsymbol {\mathsf {L}}_{1}^{H}}{\partial y}+\boldsymbol {\mathsf {L}}_{2}^{H}D^{2}+2\frac {\partial \boldsymbol {\mathsf {L}}_{2}^{H}}{\partial y}D+\frac {\partial ^{2}\boldsymbol {\mathsf {L}}_{2}^{H}}{\partial y^{2}}\right )\hat {\boldsymbol {q}}^{\dagger },\hat {\boldsymbol {q}}\right \rangle +\textrm {B.T.}, \end{equation}

where $D=\partial /\partial y$ , B.T. stands for boundary terms produced during integration:

(B3) \begin{equation} \textrm {B.T.}=\underset {=0}{\underbrace {\left .\hat {\boldsymbol {q}}^{\dagger H}\boldsymbol {\mathsf {L}}_{1}\hat {\boldsymbol {q}}\right |_{0}^{\infty }}}+\underset {\textrm {need }\hat {u}^{\dagger }=\hat {v}^{\dagger }=\hat {w}^{\dagger }=\hat {T}^{\dagger }=0}{\underbrace {\left .\hat {\boldsymbol {q}}^{\dagger H}\boldsymbol {\mathsf {L}}_{2}D\hat {\boldsymbol {q}}\right |_{0}^{\infty }}}\underset {=0}{\underbrace {-\left .D\hat {\boldsymbol {q}}^{\dagger H}\boldsymbol {\mathsf {L}}_{2}\hat {\boldsymbol {q}}\right |_{0}^{\infty }-\left .\hat {\boldsymbol {q}}^{\dagger H}D\boldsymbol {\mathsf {L}}_{2}\hat {\boldsymbol {q}}\right |_{0}^{\infty }}}. \end{equation}

To make sure the B.T. vanishes, we specify the boundary conditions for adjoint vectors:

(B4) \begin{equation} \hat {u}^{\dagger }=\hat {v}^{\dagger }=\hat {w}^{\dagger }=\hat {T}^{\dagger }=0,\qquad (y=0,\;\infty ). \end{equation}

According to the definition of the adjoint equations (2.20), the following relation holds:

(B5) \begin{equation} \boldsymbol {\mathsf {L}}^{\dagger }\left (\alpha ^{\dagger }\right )\hat {\boldsymbol {q}}^{\dagger }=\boldsymbol {\mathsf {L}}^{\dagger }\left (\alpha ^{*}\right )\hat {\boldsymbol {q}}^{\dagger }\Leftrightarrow \alpha ^{\dagger }=\alpha ^{*}. \end{equation}

In addition, the bi-orthogonal relationship is obtained from the definition of the adjoint equations. Examine $\hat {\boldsymbol {q}}_{i}$ and $\hat {\boldsymbol {q}}_{j}^{\dagger }$ ,

(B6) \begin{equation} \left \langle \boldsymbol {\mathsf {L}}^{\dagger }\left (\alpha _{j}^{\dagger }\right )\hat {\boldsymbol {q}}_{j}^{\dagger },\hat {\boldsymbol {q}}_{i}\right \rangle =\left \langle \hat {\boldsymbol {q}}_{j}^{\dagger },\boldsymbol {\mathsf {L}}\left (\alpha _{j}\right )\hat {\boldsymbol {q}}_{i}\right \rangle =\left \langle \hat {\boldsymbol {q}}_{j}^{\dagger },\boldsymbol {\mathsf {L}}\left (\alpha _{i}\right )\hat {\boldsymbol {q}}_{i}\right \rangle =0 \end{equation}

gives

(B7) \begin{equation} \left \langle \hat {\boldsymbol {q}}_{j}^{\dagger },\left [\boldsymbol {\mathsf {L}}\left (\alpha _{i}\right )-\boldsymbol {\mathsf {L}}\left (\alpha _{j}\right )\right ]\hat {\boldsymbol {q}}_{i}\right \rangle =0. \end{equation}

Equations (B5) and (B7) have been used to validate the adjoint equations derived.

Appendix C. Normalisation of the eigenvector and its adjoint

Following discussions in § 2.3, we normalise $\hat {\boldsymbol {q}}$ and $\hat {\boldsymbol {q}}^{\dagger }$ such that the following term is equal to unity:

(C1) \begin{align} \begin{split} \left \langle \hat {\boldsymbol {q}}^{\dagger },\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \alpha }\hat {\boldsymbol {q}}\right \rangle = & \int _{0}^{\infty }\hat {\boldsymbol {q}}^{\dagger H}\frac {\partial \boldsymbol {\mathsf {L}}}{\partial \alpha }\hat {\boldsymbol {q}}\:{\rm d}y\\[5pt] =\int _{0}^{\infty } & \hat {\rho }^{\dagger *}\left (-iu\hat {\rho }-i\rho \hat {u}\right )\:{\rm d}y\\[5pt] +\int _{0}^{\infty } & \hat {u}^{\dagger *}\left (-i\frac {\partial p}{\partial \rho }\hat {\rho }+i\rho u\hat {u}-i\frac {\partial p}{\partial T}\hat {T}\right )\:{\rm d}y\\[5pt] +\int _{0}^{\infty } & \hat {u}^{\dagger *}\left (-\frac {2\alpha \left (2\mu +\lambda \right )}{Re}\hat {u}+\frac {i\left (\mu +\lambda \right )D\hat {v}}{Re}+\frac {i}{Re}\frac {\partial \mu }{\partial y}\hat {v}-\frac {\beta \left (\mu +\lambda \right )}{Re}\hat {w}\right )\:{\rm d}y\\[5pt] +\int _{0}^{\infty } & \hat {v}^{\dagger *}\left (\frac {i}{Re}\frac {\partial \mu }{\partial \rho }\frac {\partial u}{\partial y}\hat {\rho }+\frac {i\left (\mu +\lambda \right )D\hat {u}}{Re}+\frac {i}{Re}\frac {\partial \lambda }{\partial y}\hat {u}\right )\:{\rm d}y\\[5pt] +\int _{0}^{\infty } & \hat {v}^{\dagger *}\left (-\frac {2\alpha \mu }{Re}\hat {v}-i\rho u\hat {v}+\frac {i}{Re}\frac {\partial \mu }{\partial T}\frac {\partial u}{\partial y}\hat {T}\right )\:{\rm d}y\\[5pt] +\int _{0}^{\infty } & \hat {w}^{\dagger *}\left (-\frac {\beta \left (\mu +\lambda \right )}{Re}\hat {u}-\frac {2\alpha \mu }{Re}\hat {w}-i\rho u\hat {w}\right )\:{\rm d}y\\[5pt] +\int _{0}^{\infty } & \hat {T}^{\dagger *}\left (-i\rho u\frac {\partial e}{\partial \rho }\hat {\rho }-ip\hat {u}+\frac {2i\mu }{Re}\frac {\partial u}{\partial y}\hat {v}-\frac {2\alpha \kappa }{RePrEc}\hat {T}-i\rho u\frac {\partial e}{\partial T}\hat {T}\right )\:{\rm d}y \end{split}. \end{align}

Appendix D. Sensitivity coefficients

(D1) \begin{align} \begin{split} \left \langle \boldsymbol {S}_{\boldsymbol {Q}},\delta \boldsymbol {Q}\right \rangle = & \left (\begin{array}{c} S_{u}\\ S_{w}\\ S_{\rho }\\ S_{T}\\ S_{p} \end{array}\right )^{H}\left (\begin{array}{c} \delta _{u}\\ \delta _{w}\\ \delta _{\rho }\\ \delta _{T}\\ \delta _{p} \end{array}\right )+\left (\begin{array}{c} S_{\frac {\partial p}{\partial \rho }}\\ S_{\frac {\partial p}{\partial T}}\\ S_{\frac {\partial ^{2}p}{\partial \rho ^{2}}}\\ S_{\frac {\partial ^{2}p}{\partial T^{2}}}\\ S_{\frac {\partial ^{2}p}{\partial \rho \partial T}} \end{array}\right )^{H}\left (\begin{array}{c} \delta _{\frac {\partial p}{\partial \rho }}\\ \delta _{\frac {\partial p}{\partial T}}\\ \delta _{\frac {\partial ^{2}p}{\partial \rho ^{2}}}\\ \delta _{\frac {\partial ^{2}p}{\partial T^{2}}}\\ \delta _{\frac {\partial ^{2}p}{\partial \rho \partial T}} \end{array}\right )+\left (\begin{array}{c} S_{\frac {\partial e}{\partial \rho }}\\ S_{\frac {\partial e}{\partial T}} \end{array}\right )^{H}\left (\begin{array}{c} \delta _{\frac {\partial e}{\partial \rho }}\\ \delta _{\frac {\partial e}{\partial T}} \end{array}\right )+\\ & \left (\begin{array}{c} S_{\mu }\\ S_{\frac {\partial \mu }{\partial \rho }}\\ S_{\frac {\partial \mu }{\partial T}}\\ S_{\frac {\partial ^{2}\mu }{\partial \rho ^{2}}}\\ S_{\frac {\partial ^{2}\mu }{\partial T^{2}}}\\ S_{\frac {\partial ^{2}\mu }{\partial \rho \partial T}} \end{array}\right )^{H}\left (\begin{array}{c} \delta _{\mu }\\ \delta _{\frac {\partial \mu }{\partial \rho }}\\ \delta _{\frac {\partial \mu }{\partial T}}\\ \delta _{\frac {\partial ^{2}\mu }{\partial \rho ^{2}}}\\ \delta _{\frac {\partial ^{2}\mu }{\partial T^{2}}}\\ \delta _{\frac {\partial ^{2}\mu }{\partial \rho \partial T}} \end{array}\right )+\left (\begin{array}{c} S_{\kappa }\\ S_{\frac {\partial \kappa }{\partial \rho }}\\ S_{\frac {\partial \kappa }{\partial T}}\\ S_{\frac {\partial ^{2}\kappa }{\partial \rho ^{2}}}\\ S_{\frac {\partial ^{2}\kappa }{\partial T^{2}}}\\ S_{\frac {\partial ^{2}\kappa }{\partial \rho \partial T}} \end{array}\right )^{H}\left (\begin{array}{c} \delta _{\kappa }\\ \delta _{\frac {\partial \kappa }{\partial \rho }}\\ \delta _{\frac {\partial \kappa }{\partial T}}\\ \delta _{\frac {\partial ^{2}\kappa }{\partial \rho ^{2}}}\\ \delta _{\frac {\partial ^{2}\kappa }{\partial T^{2}}}\\ \delta _{\frac {\partial ^{2}\kappa }{\partial \rho \partial T}} \end{array}\right ) \end{split} \end{align}

(D2) \begin{align} S_{u}= & \frac {\partial }{\partial y}\left (\rho \hat {v}^{*}\hat {u}^{\dagger }\right )\nonumber\\ & +i\alpha ^{*}\left (\hat {\rho }^{*}\hat {\rho }^{\dagger }+\rho \hat {u}^{*}\hat {u}^{\dagger }+\rho \hat {v}^{*}\hat {v}^{\dagger }+\rho \hat {w}^{*}\hat {w}^{\dagger }+\rho \frac {\partial e}{\partial \rho }\hat {\rho }^{*}\hat {T}^{\dagger }+\rho \frac {\partial e}{\partial T}\hat {T}^{*}\hat {T}^{\dagger }\right ) \end{align}
\begin{align*} & +\frac {2}{Re}\frac {\partial }{\partial y}\left [\mu \left (i\alpha ^{*}\hat {v}^{*}-\frac {\partial \hat {u}^{*}}{\partial y}\right )\hat {T}^{\dagger }\right ]\nonumber\\ & +\frac {1}{Re}\frac {\partial }{\partial y}\left [\frac {\partial \mu }{\partial T}\hat {T}^{*}\left (\frac {\partial \hat {u}^{\dagger }}{\partial y}-2\frac {\partial u}{\partial y}\hat {T}^{\dagger }+i\alpha ^{*}\hat {v}^{\dagger }\right )\right ]\nonumber\\ & +\frac {1}{Re}\frac {\partial }{\partial y}\left [\frac {\partial \mu }{\partial \rho }\hat {\rho }^{*}\left (\frac {\partial \hat {u}^{\dagger }}{\partial y}-2\frac {\partial u}{\partial y}\hat {T}^{\dagger }+i\alpha ^{*}\hat {v}^{\dagger }\right )\right ] \nonumber \end{align*}
(D3) \begin{equation} \begin{aligned} S_{w}= & \frac {\partial }{\partial y}\left (\rho \hat {v}^{*}\hat {w}^{\dagger }\right )+i\beta \left (\hat {\rho }^{*}\hat {\rho }^{\dagger }+\rho \hat {u}^{*}\hat {u}^{\dagger }+\rho \hat {v}^{*}\hat {v}^{\dagger }+\rho \hat {w}^{*}\hat {w}^{\dagger }+\rho \frac {\partial e}{\partial \rho }\hat {\rho }^{*}\hat {T}^{\dagger }+\rho \frac {\partial e}{\partial T}\hat {T}^{*}\hat {T}^{\dagger }\right )\\[10pt] & +\frac {i\beta }{Re}\frac {\partial }{\partial y}\left (2\mu \hat {v}^{*}\hat {T}^{\dagger }+\frac {\partial \mu }{\partial \rho }\hat {\rho }^{*}\hat {v}^{\dagger }+\frac {\partial \mu }{\partial T}\hat {T}^{*}\hat {v}^{\dagger }\right )\\[10pt] & +\frac {1}{Re}\frac {\partial }{\partial y}\left (\frac {\partial \mu }{\partial T}\hat {T}^{*}D\hat {w}^{\dagger }+\frac {\partial \mu }{\partial \rho }\hat {\rho }^{*}D\hat {w}^{\dagger }\right )\\[10pt] & -\frac {2}{Re}\frac {\partial }{\partial y}\left (\mu D\hat {w}^{*}\hat {T}^{\dagger }+\frac {\partial \mu }{\partial \rho }\frac {\partial w}{\partial y}\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial \mu }{\partial T}\frac {\partial w}{\partial y}\hat {T}^{*}\hat {T}^{\dagger }\right ) \end{aligned} \end{equation}
(D4) \begin{equation} \begin{aligned} S_{T}= & \frac {\partial }{\partial y}\left (\frac {\partial ^{2}p}{\partial \rho \partial T}\hat {\rho }^{*}\hat {v}^{\dagger }+\frac {\partial ^{2}p}{\partial T^{2}}\hat {T}^{*}\hat {v}^{\dagger }\right )\\[10pt] & -\frac {1}{Re}\frac {\partial }{\partial y}\left (\frac {\partial u}{\partial y}\frac {\partial ^{2}\mu }{\partial \rho \partial T}\hat {\rho }^{*}\hat {u}^{\dagger }+\frac {\partial u}{\partial y}\frac {\partial ^{2}\mu }{\partial T^{2}}\hat {T}^{*}\hat {u}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial ^{2}\mu }{\partial \rho \partial T}\hat {\rho }^{*}\hat {w}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial ^{2}\mu }{\partial T^{2}}\hat {T}^{*}\hat {w}^{\dagger }\right )\\[10pt] & -\frac {1}{RePrEc}\frac {\partial }{\partial y}\left (\frac {\partial \kappa }{\partial \rho }D\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial ^{2}\kappa }{\partial \rho ^{2}}\frac {\partial \rho }{\partial y}\hat {\rho }^{*}\hat {T}^{\dagger }+2\frac {\partial ^{2}\kappa }{\partial \rho \partial T}\frac {\partial T}{\partial y}\hat {\rho }^{*}\hat {T}^{\dagger }\right )\\[10pt] & -\frac {1}{RePrEc}\frac {\partial }{\partial y}\left (\frac {\partial \kappa }{\partial T}D\hat {T}^{*}\hat {T}^{\dagger }+\frac {\partial ^{2}\kappa }{\partial \rho \partial T}\frac {\partial \rho }{\partial y}\hat {T}^{*}\hat {T}^{\dagger }+2\frac {\partial ^{2}\kappa }{\partial T^{2}}\frac {\partial T}{\partial y}\hat {T}^{*}\hat {T}^{\dagger }\right )\\[10pt] & +\frac {1}{RePrEc}\frac {\partial ^{2}}{\partial y^{2}}\left (\frac {\partial \kappa }{\partial \rho }\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial \kappa }{\partial T}\hat {T}^{*}\hat {T}^{\dagger }\right ) \end{aligned}\end{equation}
(D5) \begin{align} S_{\rho }= & i\alpha ^{*}\hat {u}^{*}\hat {\rho }^{\dagger }-\left (D\hat {v}^{*}\right )\hat {\rho }^{\dagger }+\frac {\partial }{\partial y}\left (\hat {v}^{*}\hat {\rho }^{\dagger }\right )+i\beta \hat {w}^{*}\hat {\rho }^{\dagger }\nonumber\\ & -\frac {\partial u}{\partial y}\hat {v}^{*}\hat {u}^{\dagger }-\frac {\partial w}{\partial y}\hat {v}^{*}\hat {w}^{\dagger }-\frac {\partial e}{\partial y}\hat {v}^{*}\hat {T}^{\dagger }\nonumber\\ & +\frac {\partial }{\partial y}\left (\frac {\partial ^{2}p}{\partial \rho ^{2}}\hat {\rho }^{*}\hat {v}^{\dagger }+\frac {\partial ^{2}p}{\partial \rho \partial T}\hat {T}^{*}\hat {v}^{\dagger }\right )\nonumber\\ & +\left (i\alpha ^{*}u-i\omega +i\beta w\right )\left (\hat {u}^{*}\hat {u}^{\dagger }+\hat {v}^{*}\hat {v}^{\dagger }+\hat {w}^{*}\hat {w}^{\dagger }+\frac {\partial e}{\partial \rho }\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial e}{\partial T}\hat {T}^{*}\hat {T}^{\dagger }\right )+ \\ & -\frac {1}{Re}\frac {\partial }{\partial y}\left (\frac {\partial u}{\partial y}\frac {\partial ^{2}\mu }{\partial \rho ^{2}}\hat {\rho }^{*}\hat {u}^{\dagger }+\frac {\partial u}{\partial y}\frac {\partial ^{2}\mu }{\partial T\partial \rho }\hat {T}^{*}\hat {u}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial ^{2}\mu }{\partial \rho ^{2}}\hat {\rho }^{*}\hat {w}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial ^{2}\mu }{\partial T\partial \rho }\hat {T}^{*}\hat {w}^{\dagger }\right )\nonumber\\ & -\frac {1}{RePrEc}\frac {\partial }{\partial y}\left (\frac {\partial ^{2}\kappa }{\partial \rho ^{2}}\frac {\partial T}{\partial y}\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial ^{2}\kappa }{\partial \rho \partial T}\frac {\partial T}{\partial y}\hat {T}^{*}\hat {T}^{\dagger }\right ) \nonumber \end{align}
(D6) \begin{equation} S_{p}=\left (i\alpha ^{*}\hat {u}^{*}-D\hat {v}^{*}+i\beta \hat {w}^{*}\right )\hat {T}^{\dagger } \end{equation}
(D7) \begin{equation} \left . \begin{array}{l} \displaystyle S_{\frac {\partial p}{\partial \rho }}=i\alpha ^{*}\hat {\rho }^{*}\hat {u}^{\dagger }-D\hat {\rho }^{*}\hat {v}^{\dagger }+i\beta \hat {\rho }^{*}\hat {w}^{\dagger } \\[10pt] \displaystyle S_{\frac {\partial p}{\partial T}}=i\alpha ^{*}\hat {T}^{*}\hat {u}^{\dagger }-D\hat {T}\hat {v}^{\dagger }+i\beta \hat {T}^{*}\hat {w}^{\dagger } \\[10pt] \displaystyle S_{\frac {\partial ^{2}p}{\partial \rho ^{2}}}=-\frac {\partial \rho }{\partial y}\hat {\rho }^{*}\hat {v}^{\dagger }\\[10pt] \displaystyle S_{\frac {\partial ^{2}p}{\partial T^{2}}}=-\frac {\partial T}{\partial y}\hat {T}^{*}\hat {v}^{\dagger }\\[10pt] \displaystyle S_{\frac {\partial ^{2}p}{\partial \rho \partial T}}=-\left (\frac {\partial T}{\partial y}\hat {\rho }^{*}+\frac {\partial \rho }{\partial y}\hat {T}^{*}\right )\hat {v}^{\dagger }\\[10pt] \displaystyle S_{\frac {\partial e}{\partial \rho }}=-\rho \frac {\partial \rho }{\partial y}\hat {v}^{*}\hat {T}^{\dagger }+\left (i\alpha ^{*}\rho u-i\omega \rho +i\beta \rho w\right )\hat {\rho }^{*}\hat {T}^{\dagger }\\[10pt] \displaystyle S_{\frac {\partial e}{\partial T}}=-\rho \frac {\partial T}{\partial y}\hat {v}^{*}\hat {T}^{\dagger }+\left (i\alpha ^{*}\rho u-i\omega \rho +i\beta \rho w\right )\hat {T}^{*}\hat {T}^{\dagger } \end{array} \right \} \end{equation}
(D8) \begin{equation} \begin{aligned} S_{\mu }= & \frac {1}{Re}\left (-\frac {4}{3}\alpha ^{*2}\hat {u}^{*}\hat {u}^{\dagger }+D^{2}\hat {u}^{*}\hat {u}^{\dagger }-\beta ^{2}\hat {u}^{*}\hat {u}^{\dagger }-\frac {1}{3}i\alpha ^{*}D\hat {v}^{*}\hat {u}^{\dagger }-\frac {1}{3}\alpha ^{*}\beta \hat {w}^{*}\hat {u}^{\dagger }\right )\\ & +\frac {1}{Re}\left (-\frac {1}{3}i\alpha ^{*}D\hat {u}^{*}\hat {v}^{\dagger }-\alpha ^{*2}\hat {v}^{*}\hat {v}^{\dagger }+\frac {4}{3}D^{2}\hat {v}^{*}\hat {v}^{\dagger }-\beta ^{2}\hat {v}^{*}\hat {v}^{\dagger }-\frac {1}{3}i\beta D\hat {w}^{*}\hat {v}^{\dagger }\right )\\ & +\frac {1}{Re}\left (-\frac {1}{3}\alpha ^{*}\beta \hat {u}^{*}\hat {w}^{\dagger }-\frac {1}{3}i\beta D\hat {v}^{*}\hat {w}^{\dagger }-\alpha ^{*2}\hat {w}^{*}\hat {w}^{\dagger }+D^{2}\hat {w}^{*}\hat {w}^{\dagger }-\frac {4}{3}\beta ^{2}\hat {w}^{*}\hat {w}^{\dagger }\right )\\ & +\frac {2}{Re}\left (\frac {\partial u}{\partial y}D\hat {u}^{*}\hat {T}^{\dagger }-i\alpha ^{*}\frac {\partial u}{\partial y}\hat {v}^{*}\hat {T}^{\dagger }-i\beta \frac {\partial w}{\partial y}\hat {v}^{*}\hat {T}^{\dagger }+\frac {\partial w}{\partial y}D\hat {w}^{*}\hat {T}^{\dagger }\right ) \end{aligned} \end{equation}
(D9) \begin{align} S_{\frac {\partial \mu }{\partial \rho }}= & \frac {1}{Re}\!\!\left (\underset {1}{\underbrace {\frac {\partial u}{\partial y}D\hat {\rho }^{*}\hat {u}^{\dagger }}}+\underset {2}{\underbrace {\frac {\partial ^{2}u}{\partial y^{2}}\hat {\rho }^{*}\hat {u}^{\dagger }}}+\underset {3}{\underbrace {\frac {\partial \rho }{\partial y}D\hat {u}^{*}\hat {u}^{\dagger }}}+\underset {4}{\underbrace {-i\alpha ^{*}\frac {\partial \rho }{\partial y}\hat {v}^{*}\hat {u}^{\dagger }}}+\underset {5}{\underbrace {-i\alpha ^{*}\frac {\partial u}{\partial y}\hat {\rho }^{*}\hat {v}^{\dagger }}}\right )\nonumber\\ + &\frac {1}{Re}\!\!\left (\underset {6}{\underbrace {-i\beta \frac {\partial w}{\partial y}\hat {\rho }^{*}\hat {v}^{\dagger }}}+\underset {7}{\underbrace {\frac {2}{3}i\alpha ^{*}\frac {\partial \rho }{\partial y}\hat {u}^{*}\hat {v}^{\dagger }}}+\underset {8}{\underbrace {\frac {4}{3}\frac {\partial \rho }{\partial y}D\hat {v}^{*}\hat {v}^{\dagger }}}+\underset {9}{\underbrace {\frac {2}{3}i\beta \frac {\partial \rho }{\partial y}\hat {w}^{*}\hat {v}^{\dagger }}}+\underset {10}{\underbrace {\frac {\partial w}{\partial y}D\hat {\rho }^{*}\hat {w}^{\dagger }}}\right )\nonumber\\ +&\frac {1}{Re}\!\!\left (\underset {11}{\underbrace {\frac {\partial ^{2}w}{\partial y^{2}}\hat {\rho }^{*}\hat {w}^{\dagger }}}+\underset {12}{\underbrace {-i\beta \frac {\partial \rho }{\partial y}\hat {v}^{*}\hat {w}^{\dagger }}}+\underset {13}{\underbrace {\frac {\partial \rho }{\partial y}D\hat {w}^{*}\hat {w}^{\dagger }}}+\underset {14}{\underbrace {\left (\frac {\partial u}{\partial y}\right )^{2}\hat {\rho }^{*}\hat {T}^{\dagger }}}+\underset {15}{\underbrace {\left (\frac {\partial w}{\partial y}\right )^{2}\hat {\rho }^{*}\hat {T}^{\dagger }}}\right ) \end{align}

(D10) \begin{equation} \begin{aligned} S_{\frac {\partial \mu }{\partial T}}= & \frac {1}{Re}\left (\frac {\partial T}{\partial y}D\hat {u}^{*}\hat {u}^{\dagger }-i\alpha ^{*}\frac {\partial T}{\partial y}\hat {v}^{*}\hat {u}^{\dagger }+\frac {\partial u}{\partial y}D\hat {T}^{*}\hat {u}^{\dagger }+\frac {\partial ^{2}u}{\partial y^{2}}\hat {T}^{*}\hat {u}^{\dagger }\right )\\[10pt] & +\frac {1}{Re}\left (\frac {2}{3}i\alpha ^{*}\frac {\partial T}{\partial y}\hat {u}^{*}\hat {v}^{\dagger }+\frac {4}{3}\frac {\partial T}{\partial y}D\hat {v}^{*}\hat {v}^{\dagger }+\frac {2}{3}i\beta \frac {\partial T}{\partial y}\hat {w}^{*}\hat {v}^{\dagger }\right )\\[10pt] & +\frac {1}{Re}\left (-i\alpha ^{*}\frac {\partial u}{\partial y}\hat {T}^{*}\hat {v}^{\dagger }-i\beta \frac {\partial w}{\partial y}\hat {T}^{*}\hat {v}^{\dagger }-i\beta \frac {\partial T}{\partial y}\hat {v}^{*}\hat {w}^{\dagger }+\frac {\partial T}{\partial y}D\hat {w}^{*}\hat {w}^{\dagger }\right )\\[10pt] & +\frac {1}{Re}\left (+\frac {\partial w}{\partial y}D\hat {T}^{*}\hat {w}^{\dagger }+\frac {\partial ^{2}w}{\partial y^{2}}\hat {T}^{*}\hat {w}^{\dagger }+\left (\frac {\partial u}{\partial y}\right )^{2}\hat {T}^{*}\hat {T}^{\dagger }+\left (\frac {\partial w}{\partial y}\right )^{2}\hat {T}^{*}\hat {T}^{\dagger }\right ) \end{aligned} \end{equation}
(D11) \begin{equation} \begin{aligned} S_{\frac {\partial ^{2}\mu }{\partial \rho ^{2}}}=\frac {1}{Re}\left (\frac {\partial u}{\partial y}\frac {\partial \rho }{\partial y}\hat {\rho }^{*}\hat {u}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial \rho }{\partial y}\hat {\rho }^{*}\hat {w}^{\dagger }\right ) \end{aligned} \end{equation}
(D12) \begin{equation} \begin{aligned} S_{\frac {\partial ^{2}\mu }{\partial T^{2}}}=\frac {1}{Re}\left (\frac {\partial u}{\partial y}\frac {\partial T}{\partial y}\hat {T}^{*}\hat {u}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial T}{\partial y}\hat {T}^{*}\hat {w}^{\dagger }\right ) \end{aligned} \end{equation}
(D13) \begin{equation} \begin{aligned} S_{\frac {\partial ^{2}\mu }{\partial \rho \partial T}}=\frac {1}{Re}\left (\frac {\partial u}{\partial y}\frac {\partial T}{\partial y}\hat {\rho }^{*}\hat {u}^{\dagger }+\frac {\partial u}{\partial y}\frac {\partial \rho }{\partial y}\hat {T}^{*}\hat {u}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial T}{\partial y}\hat {\rho }^{*}\hat {w}^{\dagger }+\frac {\partial w}{\partial y}\frac {\partial \rho }{\partial y}\hat {T}^{*}\hat {w}^{\dagger }\right ) \end{aligned} \end{equation}
(D14) \begin{equation} \left . \begin{array}{l} \displaystyle S_{\kappa }=\frac {1}{RePrEc}\left (-\alpha ^{*2}+D^{2}-\beta ^{2}\right )\hat {T}^{*}\hat {T}^{\dagger } \\[10pt] \displaystyle S_{\frac {\partial \kappa }{\partial \rho }}=\frac {1}{RePrEc}\left (\frac {\partial T}{\partial y}D\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial ^{2}T}{\partial y^{2}}\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial \rho }{\partial y}D\hat {T}^{*}\hat {T}^{\dagger }\right ) \\[10pt] \displaystyle S_{\frac {\partial \kappa }{\partial T}}=\frac {1}{RePrEc}\left (2\frac {\partial T}{\partial y}D\hat {T}^{*}\hat {T}^{\dagger }+\frac {\partial ^{2}T}{\partial y^{2}}\hat {T}^{*}\hat {T}^{\dagger }\right )\\[10pt] \displaystyle S_{\frac {\partial ^{2}\kappa }{\partial \rho ^{2}}}=\frac {1}{RePrEc}\frac {\partial \rho }{\partial y}\frac {\partial T}{\partial y}\hat {\rho }^{*}\hat {T}^{\dagger }\\[10pt] \displaystyle S_{\frac {\partial ^{2}\kappa }{\partial T^{2}}}=\frac {1}{RePrEc}\left (\frac {\partial T}{\partial y}\right )^{2}\hat {T}^{*}\hat {T}^{\dagger }\\[10pt] \displaystyle S_{\frac {\partial ^{2}\kappa }{\partial \rho \partial T}}=\frac {1}{RePrEc}\left (\left (\frac {\partial T}{\partial y}\right )^{2}\hat {\rho }^{*}\hat {T}^{\dagger }+\frac {\partial T}{\partial y}\frac {\partial \rho }{\partial y}\hat {T}^{*}\hat {T}^{\dagger }\right ) \end{array} \right \} \end{equation}

Appendix E. Influence of the wavenumber on sensitivity

Figure 19. (a) Sensitivity measure $\log (M)$ as a function of the wavenumber $\beta$ . (b) Growth rate versus $\beta$ . The fluid is in the pseudo-boiling regime with $x=1$ , $\omega =40$ (inviscid TS mode).

We provide an example of the influence of the spanwise wavenumber on the sensitivity. Figure 19 shows that the sensitivity measure is only mildly influenced by $\beta$ , with the rank of input terms largely unaffected. The fluid is in the pseudo-boiling regime.

References

Alferez, N. & Touber, E. 2017 One-dimensional refraction properties of compression shocks in non-ideal gases. J. Fluid Mech. 814, 185221.CrossRefGoogle Scholar
Bagheri, S., Brandt, L. & Henningson, D.S. 2009 Input–output analysis, model reduction and control of the flat-plate boundary layer. J. Fluid Mech. 620, 263298.CrossRefGoogle Scholar
Banuti, D. 2015 Crossing the widom-line–supercritical pseudo-boiling. J. Supercritical Fluids 98, 1216.CrossRefGoogle Scholar
Barth, H. P., Hein, S. & Rosemann, H. 2018 Redesigned swept flat-plate experiment for crossflow-induced transition studies. In New Results in Numerical and Experimental Fluid Mechanics XI, (eds. A. Dilmann et al.), pp. 155165. Springer.CrossRefGoogle Scholar
Bell, I.H., Wronski, J., Quoilin, S. & Lemort, V. 2014 Pure and pseudo-pure fluid thermophysical property evaluation and the open-source thermophysical property library coolprop. Ind. Engng Chem. Res. 53 (6), 24982508.CrossRefGoogle ScholarPubMed
Bellan, J. 2020 High-Pressure Flows for Propulsion Applications. American Institute of Aeronautics and Astronautics, Inc.CrossRefGoogle Scholar
Boldini, P.C., Hirai, R., Costa, P., Peeters, J.W. & Pecnik, R. 2025 CUBENS: a GPU-accelerated high-order solver for wall-bounded flows with non-ideal fluids. Comput. Phys. Commun. 309, 109507.CrossRefGoogle Scholar
Bottaro, A., Corbett, P. & Luchini, P. 2003 The effect of base flow variation on flow stability. J. Fluid Mech. 476, 293302.CrossRefGoogle Scholar
Brandt, L., Sipp, D., Pralits, J.O. & Marquet, O. 2011 Effect of base-flow variation in noise amplifiers: the flat-plate boundary layer. J. Fluid Mech. 687, 503528.CrossRefGoogle Scholar
Browne, O.M., Hasnine, A., Sayed, M. & Brehm, C. 2021 Numerical method for particulate-induced high-speed boundary-layer transition simulations. AIAA J. 59 (4), 11961213.CrossRefGoogle Scholar
Brynjell-Rahkola, M., Shahriari, N., Schlatter, P., Hanifi, A. & Henningson, D.S. 2017 Stability and sensitivity of a cross-flow-dominated Falkner–Skan–Cooke boundary layer with discrete surface roughness. J. Fluid Mech. 826, 830850.CrossRefGoogle Scholar
Bugeat, B., Boldini, P., Hasan, A. & Pecnik, R. 2024 Instability in strongly stratified plane couette flow with application to supercritical fluids. J. Fluid Mech. 984, A31.CrossRefGoogle ScholarPubMed
Bugeat, B., Boldini, P.C. & Pecnik, R. 2022 On the new unstable mode in the boundary layer flow of supercritical fluids. In Proceedings of the 12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP-12). Available at: http://www.tsfp-conference.org/proceedings/proceedings-of-tsfp-12-2022-osaka.html.Google Scholar
Chai, L. & Tassou, S. A. 2020 A review of printed circuit heat exchangers for helium and supercritical CO2 Brayton cycles. Therm. Sci. Engng Prog. 18, 100543.CrossRefGoogle Scholar
Chen, X., Wan, B., Tu, G., Duan, M., Li, X. & Chen, J. 2024 Non-modal growth analysis of high-speed flows over an inclined cone. arXiv: 2406.18803.Google Scholar
Cooke, J. 1950 The boundary layer of a class of infinite yawed cylinders. In Mathematical Proceedings of the Cambridge Philosophical Society, vol. 46, pp. 645648. Cambridge University Press.CrossRefGoogle Scholar
Dörr, P.C. & Kloker, M.J. 2017 Crossflow transition control by upstream flow deformation using plasma actuators. J. Appl. Phys. 121 (6), 063303.CrossRefGoogle Scholar
Fedorov, A. & Tumin, A. 2011 High-speed boundary-layer instability: old terminology and a new framework. AIAA J. 49 (8), 16471657.CrossRefGoogle Scholar
Fjørtoft, R. 1950 Application of integral theorems in deriving criteria of stability for laminar flows and for the baroclinic circular vortex. Geophys. Publ. 17 (6), 152.Google Scholar
Giannetti, F. & Luchini, P. 2007 Structural sensitivity of the first instability of the cylinder wake. J. Fluid Mech. 581, 167197.CrossRefGoogle Scholar
Guardone, A., Colonna, P., Pini, M. & Spinelli, A. 2024 Nonideal compressible fluid dynamics of dense vapors and supercritical fluids. Annu. Rev. Fluid Mech. 56 (2024), 241269.CrossRefGoogle Scholar
Guo, P., Gao, Z., Jiang, C. & Lee, C.-H. 2021 Sensitivity analysis on supersonic-boundary-layer stability subject to perturbation of flow parameters. Phys. Fluids 33 (8), 084111.CrossRefGoogle Scholar
Huber, M.L., Lemmon, E.W., Bell, I.H. & McLinden, M.O. 2022 The NIST refprop database for highly accurate properties of industrially important fluids. Ind. Engng Chem. Res. 61 (42), 1544915472.CrossRefGoogle ScholarPubMed
Huber, M.L., Sykioti, E., Assael, M.J. & Perkins, R.A. 2016 Reference correlation of the thermal conductivity of carbon dioxide from the triple point to 1100K and up to 200 MPa. J. Phys. Chem. Ref. Data 45 (1), 013102.CrossRefGoogle Scholar
Hunt, J. & Graham, J. 1978 Free-stream turbulence near plane boundaries. J. Fluid Mech. 84 (2), 209235.CrossRefGoogle Scholar
Juniper, M. P. & Sujith, R. I. 2018 Sensitivity and nonlinearity of thermoacoustic oscillations. Annu. Rev. Fluid Mech. 50 (1), 661689.CrossRefGoogle Scholar
Kato, T. 2013 Perturbation Theory for Linear Operators. Vol. 132. Springer Science & Business Media.Google Scholar
Li, Z., Banuti, D.T., Ren, J., Lyu, J., Wang, H. & Chu, X. 2024 Microscale physics and macroscale convective heat transfer in supercritical fluids. arXiv: 2406.04119.Google Scholar
Luchini, P. & Bottaro, A. 2014 Adjoint equations in stability analysis. Annu. Rev. Fluid Mech. 46 (1), 493517.CrossRefGoogle Scholar
Ly, N. & Ihme, M. 2022 Destabilization of binary mixing layer in supercritical conditions. J. Fluid Mech. 945, R2.CrossRefGoogle Scholar
Mack, L. M. 1984 Boundary-layer linear stability theory. Agard Rep 709 (3), 13.Google Scholar
Marquet, O., Sipp, D. & Jacquin, L. 2008 Sensitivity analysis and passive control of cylinder flow. J. Fluid Mech. 615, 221252.CrossRefGoogle Scholar
Menikoff, R. & Plohr, B J. 1989 The Riemann problem for fluid flow of real materials. Rev. Mod. Phys. 61 (1), 75130.CrossRefGoogle Scholar
Park, J. & Zaki, T.A. 2019 Sensitivity of high-speed boundary-layer stability to base-flow distortion. J. Fluid Mech. 859, 476515.CrossRefGoogle Scholar
Pizzarelli, M. 2018 The status of the research on the heat transfer deterioration in supercritical fluids: a review. Intl Commun. Heat Mass 95, 132138.CrossRefGoogle Scholar
Poulain, A., Content, C., Rigas, G., Garnier, E. & Sipp, D. 2024 Adjoint-based linear sensitivity of a supersonic boundary layer to steady wall blowing–suction/heating–cooling. J. Fluid Mech. 978, A16.CrossRefGoogle Scholar
Pralits, J.O., Airiau, C., Hanifi, A. & Henningson, D.S. 2000 Sensitivity analysis using adjoint parabolized stability equations for compressible flows. Flow Turbul. Combust. 65 (3/4), 321346.CrossRefGoogle Scholar
Rayleigh, L. 1880 On the stability, or instability, of certain fluid motions. Proc. Lond. Math. Soc. 9, 5770.Google Scholar
Ren, J., Fu, S. & Pecnik, R. 2019 a Linear instability of Poiseuille flows with highly non-ideal fluids. J. Fluid Mech. 859, 89125.CrossRefGoogle Scholar
Ren, J. & Kloker, M. 2022 a Instabilities in three-dimensional boundary-layer flows with a highly non-ideal fluid. J. Fluid Mech. 950, A9.CrossRefGoogle Scholar
Ren, J. & Kloker, M. 2022 b Non-ideal gas behavior matters in hydrodynamic instability. Sci. Rep. 12 (1), 22089.CrossRefGoogle ScholarPubMed
Ren, J., Marxen, O. & Pecnik, R. 2019 b Boundary-layer stability of supercritical fluids in the vicinity of the widom line. J. Fluid Mech. 871, 831864.CrossRefGoogle Scholar
Reshotko, E. 2008 Transition issues for atmospheric entry. J. Spacecr. Rockets 45 (2), 161164.CrossRefGoogle Scholar
Robinet, J.-C. & Gloerfelt, X. 2019 Instabilities in non-ideal fluids. J. Fluid Mech. 880, 14.CrossRefGoogle Scholar
Saric, W.S. 1994 Görtler vortices. Annu. Rev. Fluid Mech. 26 (1), 379409.CrossRefGoogle Scholar
Saric, W.S., Reed, H.L. & White, E.B. 2003 Stability and transition of three-dimensional boundary layers. Annu. Rev. Fluid Mech. 35 (1), 413440.CrossRefGoogle Scholar
Schmid, P.J. & Henningson, D.S. 2001 Stability and Transition in Shear Flows. Springer.CrossRefGoogle Scholar
Schneider, S.P. 2001 Effects of high-speed tunnel noise on laminar-turbulent transition. J. Spacecr. Rockets 38 (3), 323333.CrossRefGoogle Scholar
Simeoni, G., Bryk, T., Gorelli, F., Krisch, M., Ruocco, G., Santoro, M. & Scopigno, T. 2010 The widom line as the crossover between liquid-like and gas-like behaviour in supercritical fluids. Nat. Phys. 6 (7), 503507.CrossRefGoogle Scholar
Spalart, P. 1989 Direct numerical study of leading-edge contamination. In AGARD Conference Proceedings No. 438: Fluid Dynamics of Three-dimensional Turbulent Shear Flows and Transition, pp. 5–1.Google Scholar
Trefethen, L.N., Trefethen, A.E., Reddy, S.C. & Driscoll, T.A. 1993 Hydrodynamic stability without eigenvalues. Science 261 (5121), 578584.CrossRefGoogle ScholarPubMed
White, M.T., Bianchi, G., Chai, L., Tassou, S.A. & Sayma, A.I. 2021 Review of supercritical CO2 technologies and systems for power generation. Appl. Therm. Engng 185, 116447.CrossRefGoogle Scholar
Figure 0

Figure 1. Comparison of the fluid property $\partial \kappa /\partial T$ (in dimension $\,\rm Wm^{-1}\,\rm K^-{^2}$) for carbon dioxide. The values are generated using (a) RefProp and (b) CoolProp. (c) Differences between the values generated by the two look-up tables. The red circle marks the critical point and the red line represents the isobar at 80.

Figure 1

Figure 2. Distribution of pressure (left axis) and pressure coefficient (right axis) for (a) wall cooling and (b) heating cases. (c) Pressure–temperature ($P{-}T$) diagram of carbon dioxide. The contours represent the compressibility factor $\bar {z}=p/(\rho RT)$. Above the critical point (red dot), the Widom line is plotted using a white dashed line. The fluid regimes are considered along the isobar of 80 ($p/p_c=1.0844$). Four groups of cases are shown with yellow (wall-heating) and cyan (wall-cooling) lines. These lines characterise the distribution of flow temperature, with wall and free stream values denoted by $w$ and $\infty$, respectively.

Figure 2

Table 1. A summary of flow cases investigated.

Figure 3

Figure 3. An overview of the inputs in the stability operator. Panels (a)–(d) show $[b]$ baseflow profiles, [EoS] equation-of-state, $[\mu ]$ viscosity and $[\kappa ]$ thermal conductivity, respectively. The pseudo-boiling and ideal gas regimes are plotted with solid and dashed lines, respectively (both under wall cooling). The red circle denotes the pseudo-critical point.

Figure 4

Figure 4. (a) Distortions of the streamwise velocity, $\delta u/u_e$, induced by viscosity alterations ($\delta \mu = \pm 10\,\%\mu , \pm 20\,\%\mu$). (b) The raw state baseflow profiles of $u$ and $w$ and (c) the corresponding stability diagram of the steady ($\omega =0$) cross-flow instability. In panel (b), the dashed line stands for the boundary-layer thickness based on $0.99u/U_\infty$, and the dotted lines denote $u/U_\infty =w/U_\infty =1$. The flow is in the supercritical regime, subject to wall heating, with gas-like fluid properties.

Figure 5

Figure 5. Baseflow distortions induced by uncertainties in the viscosity model. The flow is in the supercritical regime. Viscosity distortion was applied in proportion to the original model, ranging from 0 % to 50 %, as shown on the $x$-axis of each panel.

Figure 6

Figure 6. Distortion magnitude of the baseflow components $\log (||\delta \boldsymbol {Q}||_2)$. The flow in the subcritical, transcritical, supercritical and ideal regimes are compared. In each regime, the uncertainty is driven by $[\mu ]$, $[\kappa ]$ and [EoS], corresponding to columns 1–3 (with wall heating) and 4–6 (with wall cooling). The driving terms are highlighted with a rectangle of dotted lines.

Figure 7

Figure 7. Relational diagram depicting the uncertainties in intrinsic fluid properties, the baseflow distortions and their influences on the eigenvalue problem.

Figure 8

Table 2. Distortions of different types.

Figure 9

Figure 8. A portray of different distortions on $\partial \mu /\partial \rho$.

Figure 10

Figure 9. Illustration of eigenvalue shift validating the sensitivity framework. (a) The eigenvalue trajectory with distorted $\partial \mu /\partial \rho$; (b) error of eigenvalue prediction as a function of the departure parameter $\epsilon$ according to (3.1).

Figure 11

Figure 10. Assembly of sensitivity coefficients in the pseudo-boiling regime with wall cooling. Solid and dashed lines represent the real and imaginary parts, respectively.

Figure 12

Figure 11. (a) The eigenvalue shift due to $\delta \boldsymbol {Q}_{\textrm {max}}$ (red dots), the term distorted is $\partial \mu /\partial \rho$; (b) similar to panel (a), but with $\delta \boldsymbol {Q}_{\textrm {noise}}$ at the same amplitude (measured with 2-norm).

Figure 13

Figure 12. N factor influenced by $\delta \boldsymbol {Q}_{\textrm {max}}$. The amplitude of an individual distortion is $\Vert \delta \boldsymbol {Q}\Vert _{2}=10^{-3}/\sqrt {200}$. The red and green curves show positive and negative shifts, respectively. The numbers $1\cdots 23$ correspond to input terms (see (2.14)) sorted by the impact on the N factor in descending order (see table 3). The instability is in the pseudo-boiling regime with dominating inviscid instability.

Figure 14

Table 3. $\delta N/N$ at four uniformly distributed observation points.

Figure 15

Figure 13. A comparison of sensitivity across different regimes and between different input variables. The numbers on the heatmap specify the log value of the sensitivity measure: $\log (M)$.

Figure 16

Figure 14. Illustration of key terms for the sensitivity profile $S_T$ corresponding to (3.5): (a) $S_T$ and $S_T(1)$; (b) $S_T(1a)$ and $S_T(1b)$; (c) $S_T(1a_1)$, $S_T(1a_2)$ and $S_T(1a_3)$; (df) the components composing $S_T(1)$.

Figure 17

Figure 15. Sensitivity to distortions of $\rho$, $u$ and $w$ (panels ac). We show dominating terms (the real parts) of each profile.

Figure 18

Figure 16. (a) Sensitivity (the real part) to $\partial ^2 p/\partial \rho ^2$ with wall cooling. The sensitivity equals the derivatives of the products of three terms shown in panels (b)–(d). We illustrate the balance of term $\hat {\rho }$, $\rho$ and $D{\rho }$ according to the continuity equation (3.9) in panels (e)–(h) for the four regimes considered (with wall cooling).

Figure 19

Figure 17. Exposition of sensitivity terms for $\partial \mu /\partial \rho$ with (a) wall heating and (b) cooling. The terms corresponding to (D9) have been labelled as 1–15 on the plot. The height of each pillar stands for the 2-norm of a term.

Figure 20

Figure 18. Comparison of N factors using non-ideal and ideal models (for [EoS], $[\mu ]$ and $[\kappa ]$). The thermodynamic regime (gas-like, liquid-like, pseudo-boiling) has been signified for relevant panels. All cases are subject to $\beta =80$, $\omega =0$ (steady CF mode), except the last panel with $\beta =80$, $\omega =20$ (inviscid TS mode).

Figure 21

Table 4. A summary of idealisable fluid models for linear stability analysis, provided correct baseflow profiles ($u$, $w$, $\rho$, $T$) are supplied.

Figure 22

Figure 19. (a) Sensitivity measure $\log (M)$ as a function of the wavenumber $\beta$. (b) Growth rate versus $\beta$. The fluid is in the pseudo-boiling regime with $x=1$, $\omega =40$ (inviscid TS mode).