Non-ideal gas behavior matters in hydrodynamic instability

Hydrodynamic instability, the foundation for flow’s laminar-turbulent transition and various predicting models, has been helping to understand the physics and shape the design of aerodynamic devices. While for hypersonic flow it is clear that thermodynamic/-chemical effects need be accounted for due to the high temperatures occurring, this letter unveils that also for low-speed flow at ambient temperatures non-ideal, i.e. real-gas effects can play a strong role—a feature missed by the classic theory for Newtonian fluids. By considering a three-dimensional low-speed boundary-layer flow in different thermodynamic regimes—subcritical, supercritical and transcritical—we show the importance of coupling thermodynamics by sensitivity studies of the perturbation growth rate to various inputs of the full stability equations. High sensitivities are found, and not only the transition-onset location but also the transition mechanism may be concerned.

When studying flow instability it is broadly acknowledged that a flow with an inflectional velocity profile generally bears inviscid instability whose growth rate is strong, paving the way towards expeditious laminar-turbulent transition. The theoretical support dates back to Rayleigh's inviscid theory 1 , and Fjørtoft's improved criterion 2 . Both assume inviscid perturbations in the stability equation (inviscid Orr-Sommerfeld equation) and allow the baseflows to be viscous but incompressible. The inflection-point criterion, although serving formally only as a necessary condition for inviscid instability, has led to a substantial and intuitive understanding simply by looking at the base-flow velocity profiles, constituting the cornerstone of classic monographs on hydrodynamic stability 3,4 . The success of the inflection-point criterion extends to the perception of the mechanisms of various controlling strategies, e.g., pressure gradient, wall suction & blowing, and wall heating & cooling from the alteration of the base-flow velocity profiles 5 . Specifically, the general non-dimensional boundary-layer equation for a 2D baseflow reads, derived from the streamwise momentum for a fluid with non-constant viscosity and density: In this article, u, v, w are used for streamwise (x), wall-normal (y) and spanwise (z) velocity components. p and µ denote pressure and dynamic viscosity. All variables have been non-dimensionalized with fixed reference values and for the pressure (ρu 2 ) ref is used. Re is the fixed Reynolds number based on the reference values, and µ and ρ are identical to one if there is no alteration of the viscosity and density with respect to their reference values which is the standard case for incompressible flow. Subscript w is used for wall parameters. Using (1), one concludes that for a standard Blasius (flat-plate) boundary layer (no pressure gradient, controlling nor heat transfer), ∂ 2 u/∂y 2 | w = 0 , indicating that the inflection point falls strictly on the wall; hence inviscid instability is absent. Further it is concluded that, upon controlling, a favorable pressure gradient ( dp/dx < 0 ), wall suction ( v w < 0 ), or cooling wall for ideal gases ( ∂µ/∂y| w > 0 ) shift the inflection point inside the wall and shape the velocity profiles fuller. The flow is consequently linearly more stable (and vice versa). When friction is included in the stability consideration the viscous Orr-Sommerfeld (O-S) equation results (For given baseflow profiles u(y) and w(y), at given x and Re one solves for the spatial growth rate α i and streamwise wavenumber α r of a fixed-frequency disturbance growing exponentially in x-direction, see the description of Eqs. (5) and (6). Viscous instability adds, does not need an inflection point in the velocity profiles, is much weaker than inviscid instability, and the above findings equally hold. Though widely used by researchers and engineers, we shall note that the O-S framework essentially misses physics in the aspect of thermodynamics brought by the fluid itself. Therefore, the completeness and truthfulness of the conventional hydrodynamic stability theory remains unclear.
(1) www.nature.com/scientificreports/ The present work investigates Newtonian fluids , i.e. fluids where the viscosity does not depend on the shear in the flow but on temperature, density or pressure of the fluid. We note in passing that non-Newtonian, viscoelastic flow instabilities 6 amount to an own subject that witnesses essential applications in daily life and industry (e.g. blood flow, healthcare products, engineering inks). One of the most striking attributes is that the Reynolds number supporting the instabilities can be significantly lower than in Newtonian flows. A shred of earlier evidence for this is the work by Forbes 7 , who proposed that even weakly non-Newtonian effects can give rise to instabilities.

Results
Overview of the general framework for flow instability with non-ideal fluids. Figure 1a encapsulates the constituent inputs toward the linear stability analysis of flows. The problem is defined both by the flow and the running fluid. In the conventional O-S framework, fluids are assumed to be incompressible, ideal, and isothermal such that the energy equation is dropped, and the physical role of the fluid thus gets minimized. Interested readers are referred to Eq. (6) of "Methods" section. Assuming incompressible flows and a perfect gas, the classical stability theory solves a simplified eigenvalue problem. In a general framework, in addition to the velocity profiles, the state-variable relation (equation of state) & transport profiles (viscosity and thermal conductivity), as well as their thermodynamic gradients, become essential inputs in the stability problem. Moreover, the velocity profiles are altered by the fluids' properties (as the classical framework accounts for). All the key parameters, e.g., Reynolds, Mach, Prandtl, and Eckert numbers, correspondingly depend on fluid properties. Carry this general view with Eq. (1), the involvement of pressure gradient and wall suction & blowing in (1) stays within the O-S framework. However, the explanation of the effect of wall cooling/heating using the inflectionpoint criterion is ill-defined. The reason is that wall cooling/heating not only influences the velocity profile but also calls for integrating the intrinsic physics of the fluid, especially in the stability equations. As an example we consider the stability of a three-dimensional boundary-layer flow because this is the most representative for the flow on rotating blades of energy conversion machines 8 . Two fundamentally different instabilities can occur-cross-flow (CF) and Tollmien-Schlichting-type (TS) instability. An alteration in the primary stability www.nature.com/scientificreports/ properties can not only lead to a transition-location shift but also to a different transition mechanism and thus control measure 9 . The baseflow corresponds to an (ideal-air) experiment 10,11 and mimics a low-speed swept flow over a wing surface subject to a favorable pressure gradient, where CF modes (comprehensive reviews are available through 12,13 ) are responsible for the laminar-turbulent transition. Researchers and engineers rely on the integral of growth rates-the N factor-in determining the transition position 14 . Figure 1b compares the growth rates and N factors of steady CF modes between the full stability equations for an ideal gas (solid lines) and the standard O-S setting with constant density and viscosity (dashed lines). The N factor amounts to the integral of the spatial growth rate starting from the neutral position of a considered disturbance: The flow is subject to wall cooling with T w /T ∞ = 15/16 . We show the isolines of N-factors with values of 2, 4, 6, and 8, where the uncontrolled isothermal reference case (dotted lines) is included to indicate the effect of wall cooling. As can be seen, if a wind tunnel measured N = 8 to indicate transition onset, the isothermal case has a transition position x iso = 4.34 . Considering wall cooling, the full equations give the transition shift x full = 0.33 while the classical O-S setting predicts x O−S = 0.52 which is 58.8% more than the correct transition shift. Using the classical O-S framework thus considerably overestimates the impact of wall cooling because the thermodynamic & transport property effects in the stability equations partly compensate the effect of the altered base-flow profiles.
As an essential step in building a more physical system, Lees and Lin 15 derived the generalized-inflectionpoint criterion for compressible flows, incorporating the contribution of flow temperature. Mack 16 advanced the viscous theory for compressible flows. He found a change in the mathematical nature of the stability equations when the Mach number increases, leading to the well-known Mack's higher modes (often dominate flow transition in hypersonic flows). However, in the compressible scheme, the physics of the fluid largely stayed in the ideal gas regime due to the erstwhile aerodynamic demand. On the other hand, the O-S framework has been extended to account for temperature-dependent viscosity 17,18 , coupling with various controlling methods 19 . In line with the emphasis on viscosity, viscosity-stratified flows amount to a particular type that attracted significant attention 20 . For example, their instabilities have been revealed for flows with horizontal heat advection 21 , mixing layer of miscible fluids 22,23 . Along different paths towards turbulence 24 , viscosity stratification plays a minimal role in algebraic instability 25 . Contrarily, the edge state is considerably modified by viscosity gradients leading to a change of threshold for transition 26 . Seen from Fig. 1a, viscosity is only one part of the real fluid physics 27,28 . An other is the state-variable relation, provoking the effect of non-ideal fluid properties, also called real-gas effects, if the pressure is high and intermolecular forces have a direct impact on the fluid behavior. Will the stability framework with ideal-gas assumption, provided with the non-ideal base-flow profiles, give the right prediction?
We aim to answer this question through this article, which shall call for attention to integrate thermodynamics into hydrodynamic stability theory and contribute to a better understanding of flows with substantive fluids, e.g., blood flows 29 , geophysical flows 30,31 , superfluid turbulence 32 and dense flows 33 .
Three-dimensional boundary-layer instability of supercritical CO 2 . Without loss of generality we consider supercritical fluids near the critical point. These fluids are recognized as important working mediums in the energy sector 8 , while their highly non-ideal thermodynamic & transport properties challenge the conventional framework. Figure 2a shows the temperature-density ( T − ρ ) diagram of CO 2 colored by Prandtl number ( Pr = µC p /κ ). C p and κ stand for isobaric specific heat capacity and thermal conductivity, respectively. We will focus on the isobar of 80 along which strong non-ideality occurs on crossing the pseudo-critical point (coinciding with the critical point at the critical pressure), also known as the Widom line 34 that characterizes the physics of pseudo-boiling, defined as the isobaric maximum of C p . The ideal regime reported in Fig. 1a corresponds to the high-temperature regime (above 600 K). As can be seen, the Prandtl number only behaves as a constant in the ideal regime while showing strong variations near the critical point, highly dependent on the thermodynamic state. With this fluid and a look-up table for its thermodynamic properties 35 , we obtain the three-dimensional boundary-layer stability diagram in the ω − β frame for x = 1 in Fig. 2b and c for wall cooling and heating, respectively. A clear comparison indicates that the effect of wall cooling/heating is highly dependent on the fluid's thermodynamic regime. For example, the flow supports another mode with tremendous growth rates and instability band in the transcritical regime. This mode is of inviscid Tollmien-Schlichting(/ Kelvin-Helmholtz) type, i.e. with the wave phase speed about oriented in main-flow direction for highest amplification, which may be expected with adverse but not at all with favorable pressure gradient 9 . Instead of CF-type transition, with involvement of steady CF vortices originating from minute surface roughness, a TS-type flow transition will set in here, despite the favorable pressure gradient. The dramatic difference suggests that new physics is present when integrating thermodynamics. Besides, the effect of wall heating gets strengthened in the supercritical regimes but less pronounced in the sub-and transcritical regimes. What has caused such large differences, the baseflow velocity or other items shown in Fig. 1a? The role of non-ideality. Figure 3a and b compares a representative modal state from the stability diagram, β = 80 and ω = 15 , for wall cooling and heating respectively. The spatial growth rate −α i and wavenumber α r are displayed. The color style is matched to Fig. 2b and c and throughout the article, indicating different thermodynamic regimes. A zoom-in of the ideal regime indicates how different simplifications bias the results. For example, the O-S framework leads to a growth rate error of 7.6%/− 1.5% for wall cooling/heating. Solely freezing viscosity or density & temperature profiles respectively also leads to noticeable eigenvalue changes, showing that www.nature.com/scientificreports/ even with an ideal gas, an accurate N-factor can not be obtained with a simplified framework once a temperature gradient exists. Considering that a conventional stability framework takes empirical laws for transport properties, we show that it is inadequate to use this framework with the non-ideal baseflow as input. In particular, the growth rate demonstrates significant errors by replacing the viscosity ( µ and related gradients) with ideal laws (e.g., Sutherland laws). We shall also note that the influence of thermal conductivity ( κ and related gradients) is much smaller and only appears crucial in the transcritical regime with cooling. Besides transport properties, Fig. 3c compares the base-flow profiles for w, ρ , ∂p/∂T , ∂µ/∂ρ . Noticeably, the velocity profiles (w(y)) alter in different thermodynamic regimes. In particular, the transcritical case with wall cooling shows an inflection point alongside strong density alteration near the wall, causing inviscid TS instability. To understand a fluid's thermodynamic role comprehensively, we present the sensitivity ǫ of the growth rates to each input of the stability equations, defined as Here χ stands for any physical input, e.g., T and ∂µ/∂ρ . Thus the sensitivity measures how the stability depends on the non-ideal inputs relative to conventional ideal assumptions. Figure 4 shows a summary of each term at different thermodynamic regimes. A first conclusion is drawn for all the non-ideal regimes: The velocity and thermodynamic profiles & gradients exert the most critical influence (top 10 rows of the table except www.nature.com/scientificreports/  Fig. 2b and c, β = 80, ω = 15 . Panel (a) is for Fig. 2b, cooling wall, and panel (b) for Fig. 2c, heating wall. We consider the influence of freezing certain inputs and O-S equations in the ideal regime. The results are denoted with small solid circles with the error of growth rate labeled. In the three non-ideal regimes (red, blue, magenta), we present the influence using ideal viscosity and thermal conductivity inputs. (c) Distribution of w, ρ , ∂p/∂T , and ∂µ/∂ρ as functions of y for the wall-cooling case. They belong to the velocity and density profiles and thermodynamic gradients originating from the state-variable relation and the transport properties (see Fig. 1a). www.nature.com/scientificreports/ the pressure), whereas the transport properties, as discussed with Fig. 3a and b, are less influential, albeit they already show a significant influence (viscosity related terms). Therefore it is essential to account for all the thermodynamical physics of the fluid (including transport properties) to obtain base profiles correctly (e.g., ρ , u, T, ∂p/∂T ). Examing Fig. 4 column by column, one concludes that the transcritical regime bears the highest sensitivities. This is due to its most substantial non-ideality on crossing the Widom line. It is also noteworthy that terms like ∂µ/∂ρ and ∂κ/∂ρ are usually ignored in conventional studies. Therefore even with correct base-flow profiles, the hydrodynamic stability framework must be upgraded to incorporate all the related profiles unless the ideal-gas assumption is met in a practical situation.

Discussion
Our study highlights the importance of accurately including the thermodynamic properties of the (non-ideal) fluid missed in conventional hydrodynamic stability theory when a temperature gradient in the investigated shear layer exists. More generally, when a fluid is operated near its critical point such that highly non-ideal thermodynamic effects are present, the conventional framework fails to give the correct growth rate, dominating mode, and consequently, prediction of the transition physics, which may render a predesigned control measure pointless. A sensitivity study shows that thermodynamics, viz. the state-variable relation (equation of state) and the viscosity (including its thermodynamic gradients), but not the thermal conductivity, must be adequately accounted for in the base-flow profiles and hydrodynamic stability framework. This also means that using an p Figure 4. Sensitivity ǫ of the spatial growth rate to different inputs of the stability operator. These inputs correspond to different categories shown in Fig. 1a. The six columns are for supercritical, subcritical, and transcritical regimes with wall cooling and heating conditions, respectively (see bottom of the figure).
In steady boundary-layer flows without separation, the streamwise viscous gradient is much smaller than the other component in wall-normal direction. The PNS are thus derived from the complete Navier-Stokes equations by dropping the streamwise gradient in the viscous terms 36 . We employ subscript e for boundary-layer edge values. Here ρ e (x) and T e (x) are potential/irrotational-flow values given by the isentropic relations: where S stands for entropy, and p(x) is the measured pressure with dp/dx < 0 . The PNS are integrated downstream using an implicit Euler scheme, starting from an initial profile at x = x 0 . In this study, we prescribe the stream-and spanwise velocities u(x 0 , y) , w(x 0 , y) using the Falkner-Skan-Cooke (FSC) solution and v(x 0 , y) = 0 . The thermodynamic variables ( ρ, T ) are given either as the potential-flow values (applying isentropic relations) or extrapolated from existing downstream data. The validation of the PNS results with different initial profiles was made, ensuring that the influence of the initial profiles is insignificant.
Linear stability with fully included thermodynamics. We consider linear modal instability of the laminar baseflow. The theory was derived based on the state postulate that a simple compressible system is defined by two independent thermodynamics properties. The perturbation vector is defined as q = ρ ′ , u ′ , v ′ , w ′ , T ′ T . Therefore, perturbations in the other transport and thermodynamic properties (e.g., e ′ , κ ′ ) are formulated as a function of (ρ ′ , T ′ ) through two-dimensional Taylor expansion 37 . The stability equations are obtained from the N-S equations by subtracting the governing equations of the unperturbed baseflow, formulated as The definition of the matrices in Eq. (5) is provided in ref 27 . The normal-mode form, q x, y, z, t =q y exp (iαx + iβz − iωt) + c.c. is used, and we seek the solution in the spatial mode, where the spanwise wavenumber β and angular frequency ω are given as real input and α is the complex eigenvalue to be determined by the resulting nonlinear eigenvalue problem: Here L is a 5-by-5 matrix. Its detailed expressions are given as: (6) Lq = 0.   www.nature.com/scientificreports/ Here D stands for the differential matrix (using the Chebyshev spectral method) with Dû = dû/dy . The perturbations are subject to Dirichlet boundary conditions ( u ′ = v ′ = w ′ = T ′ = 0 ) at the wall and in the freestream. Compared to the conventional LST framework, the current method has been verified to reproduce the classic results in the ideal-gas regime 38 . As discussed in the previous section, when the fluid operates near its critical point, a perturbation's growth rate experiences enormous sensitivity to the input of the eigenvalue problem. Therefore, the thermodynamic and transport models are essential to obtaining physically accurate results. A comparison using different equations of state (van der Waals, Redlich-Kwong, Peng-Robinson) 27 shows no convergence of the result even if the other models are kept the same. Therefore, an accurate look-up table has been used to generate inputs as outlined in the left column of figure 4. The eigenvalue problem (6) is discretized with Chebyshev collocation points and numerically solved with 300 grid points along the wall-normal direction. The results have been tested to be grid-independent.

Data availability
The data that support the findings of this study are available from the corresponding author upon reasonable request. www.nature.com/scientificreports/ Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.