Decaying compressible turbulence with thermal non-equilibrium

Turbulence is ubiquitous in a wide range of natural and engineering phenomena, such as formation of galaxies,^1,2 mixing of fuel and oxidizer in engines,³ and numerous processes in the atmosphere.⁴ In engineering applications, modeling and control of turbulence is of paramount interest. For example, development of planetary re-entry systems and future supersonic aircraft depend critically on the understanding of hypersonic flows especially under turbulent conditions.^5,6 In many of these systems, turbulence can be studied under the so-called continuum assumption in thermal equilibrium, where the molecular structure of the gas is not incorporated directly into the formulation. This is an appropriate assumption for low-speed flows and has been invoked widely to study the fundamental aspects of turbulence processes. However, high-speed flows are characterized by very high temperatures⁷ at which molecular energy modes, such as vibrational, rotational, or electronic modes, may not be in thermodynamic equilibrium and need to be accounted for to define the thermodynamic state of the gas.^7–11 It is therefore important to account for and characterize these processes in high-speed compressible flows. A substantial body of literature exists for laminar flows under those conditions.^7,8,12 However, studies in the case of turbulent flows are limited and thus studying turbulent flows is the objective of the present work.

Rapid changes in the thermodynamic state of fluid elements may be induced by natural means such as shocks⁷ and turbulent fluctuations^13–15 or induced deliberately in the laboratory with, for example, laser excitation.¹⁶ If these changes are fast enough, one or more modes may be out of equilibrium, in which case the gas is said to be in thermal non-equilibrium (TNE). Relaxation toward equilibrium of different modes is achieved through molecular collisions with varying degree of efficiency. Rotational modes only require order 10-100 collisions while vibrational modes may require more than three orders of magnitude more collisions to equilibrate.¹² In other words, rotational modes relax toward equilibrium faster than vibrational modes. In fact, in many situations rotational modes can be assumed to be in equilibrium, but vibrational modes relax at timescales comparable to hydrodynamic processes in the flow and can thus have profound effects on the flow dynamics.

There are several studies that highlight this strong interaction between hydrodynamic and TNE processes. Bertolotti,¹⁷ for example, studied the effect of thermal non-equilibrium on the stability of boundary layers and found that both rotational and vibrational non-equilibrium have an effect. Rotational non-equilibrium was observed to dampen high-frequency instabilities and vibrational non-equilibrium affected the growth of disturbances. However this, as most other studies, investigated laminar flows.^7,8 There are only a few studies on the interaction between TNE and turbulence.^13–16 It is thus not surprising that the general understanding of these flows is poor in terms of governing parameters, energy dynamics, and conditions and regimes in which this interaction could happen. This is a major thrust in the present work.

It has recently been shown¹³ that stationary turbulence can significantly alter the distribution of energy across vibrational, rotational, and translational modes. Furthermore, complete relaxation of TNE is never achieved in turbulent flows though the degree of non-equilibrium may be small depending on conditions. In the case of decaying flows, however, in addition to the effects studied in Ref. 13, one has to account for the time evolution of hydrodynamic and thermodynamic quantities. Liao, Peng, and Luo¹⁸ studied the effect of rotational non-equilibrium on decaying isotropic turbulence. They observed that relaxation of rotational non-equilibrium weakens compressibility after a short initial time. Rotational non-equilibrium was also observed to affect dissipation. The authors noted that the effect on compressibility is more significant for flows with high Mach numbers, which are generally characterized by stronger compressibility levels. Experimentally, it has also been found that vibrational non-equilibrium has a clear effect on the decay of turbulence. In particular, Fuller et al.¹⁶ found that the decay is faster when the degree of initial non-equilibrium is increased. This seems to be different from the conclusion arrived at by Neville et al.¹⁹ who found a negligible effect on the decay of isotropic turbulence in direct numerical simulations. They also reported a damping of thermodynamic fluctuations in the presence of vibrational non-equilibrium. As we show below, however, the observations and disagreement from both studies can be explained by the results presented here.

In this article, we study decaying isotropic turbulence in the presence of vibrational TNE to characterize and understand the observed changes in decay rates. For this, we use well-resolved direct numerical simulations (DNS) at a range of parameters characterizing the TNE, in particular, the degree of initial non-equilibrium and the characteristic relaxation time both of which depend on the excitation mechanism and the molecular system considered.

To study the two-way coupling between turbulence and TNE, we performed and analyzed a large database of direct numerical simulations (DNS), which resolve all dynamically relevant spatial and temporal turbulent scales as well as molecular relaxation processes. The governing equations expressing conservation of mass and momentum are given by

where u_i is the flow velocity, ρ is the fluid density, p is the pressure, and σ_ij is the viscous stress tensor, which corresponds to a Newtonian fluid, that is, $σij=μ∂ui∂xj+∂uj∂xi−23δij∂uk∂xk$⁠. Viscosity, as commonly done in compressible flows,^20,21 follows a power law with temperature of the form T^a with a = 0.5.

Since translational and rotational energies are assumed to be in equilibrium,¹³ we will not differentiate between them and use e to denote the sum of these two contributions. In the case of a perfect gas of diatomic molecules, as assumed here, we have e = (5/2)RT, where R is the gas constant. Conservation of translational–rotational energy can then be written as

where κ is the thermal conductivity and $sij=12∂ui∂xj+∂uj∂xi$ is the strain rate tensor. The Prandtl number is assumed to be constant, Pr = 0.72, which implies κ to follow the same temperature dependence as viscosity. The last term in Eq. (2), $Qv$⁠, accounts for the transfers of energy between different modes of colliding molecules, in this case between translational–rotational modes and vibrational modes. Correspondingly, the evolution equation for the vibrational energy takes the familiar advection–diffusion form,²² which is coupled to the system Eqs. (1) and (2),

where $ev$ is the vibrational energy per unit mass and $κv$ is the corresponding diffusion coefficient for vibrational energy.²² The appearance of $Qv$ in both Eqs. (2) and (3) with opposite signs highlights the exchange nature of this process, which does not alter the total energy in the system, but simply distributes it across internal molecular modes.

In a wide range of applications, the exchange term $Qv$ can be well approximated by the so-called Landau–Teller approximation,^23,24 which can be written as

where $τv$ is the characteristic relaxation timescale of the vibrational mode, and $ev*$ is the equilibrium vibrational energy. If we use the energy per unit volume instead of per unit mass (⁠$Ev=ρev$⁠, $Ev*=ρev*$⁠, E = ρe, etc.), we can simply write $Qv=(Ev−Ev*)/τv$⁠. For harmonic oscillators at a translational temperature T, quantum mechanics⁸ yields an explicit expression for $ev*$⁠,

where $θv$ is the characteristic vibrational temperature, which depends on the molecular system considered, and R is the gas constant. The vibrational relaxation timescale, $τv$⁠, can also be derived from first principles,²² 

where $c1τv$ and $c2τv$ are constants that depend on the molecular system. Some slightly different generalized forms with different number of adjustable constants have been put forth,^24,25 which appear to describe more accurately the observed vibrational relaxation of a wider range of molecular systems. The main conclusions presented here, however, do not depend on these differences. As in previous investigations, the system of equations is closed with a perfect gas equation of state.^13,17

The governing equations are solved in a triply periodic domain using a hybrid OpenMP/MPI massively parallel implementation, which has been extensively validated^13,21,26 and has been shown to scale up to very large number of cores.²⁷ The code is based on tenth-order compact differentiation in space coupled with a low-storage third-order Runge–Kutta in time. Adequate resolution of quantities of interest is assured by grid convergence studies, the details of which can be found in Jagannathan and Donzis.²¹ 

Numerical stability is maintained by limiting the time step size. Typically, this is imposed with a Courant–Friedrichs–Levy (CFL) condition based on the advective term. In the presence of TNE, the additional timescale $τv$ could provide stricter temporal resolution constraints. This will depend on the gas and thermodynamic state of the flow. In particular, to resolve TNE relaxation, we require the time step size Δt to also satisfy $Δt≪τvmin$⁠, where $τvmin$ is the smallest $τv$ in the domain at a given time step. Thus, the code dynamically selects the smallest time step size determined by the standard hydrodynamic CFLs (convective and diffusive) and the TNE CFL (≡Δt/$τv$⁠). From convergence studies, we found that proper capture of the relaxation process is achieved with a TNE CFL of order 0.01, a value that has been used in all simulations presented here. The small value of this CFL highlights the additional cost when resolving all temporal scales in these kinds of simulations.

The introduction of vibrational non-equilibrium is done in a manner consistent with the photoexcitation of vibrational modes in the experiments of Ref. 16 to facilitate comparisons. The initial condition for the decaying simulations correspond to fully developed turbulent fields stochastically forced at the largest scales, with Taylor Reynolds number (R_λ) of 60 and turbulent Mach number (M_t) of 0.1 and 0.4.²⁰ At t = 0 forcing is removed and turbulence is allowed to decay with the simultaneous introduction of additional vibrational energy. This process increases vibrational energy uniformly at different predetermined levels above its equilibrium value in an analogous way to the photoexcitation with lasers or RF plasma¹⁶ in experiments. This approach is also similar to that used to study rotational non-equilibrium in gas kinetic simulations.¹⁸ In this study, the initial vibrational energy, $⟨Ev0⟩$ (subscript zero is used throughout to indicate initial conditions, and angular brackets ⟨·⟩ are used to indicate volume averages) is higher than that at equilibrium, i.e., $⟨Ev0⟩>⟨Ev0*⟩$⁠, giving rise to an initially vibrationally hot state. The parameters characterizing the initial state of TNE are the ratio of vibrational and translational–rotational energy (ξ₀) and the initial degree of non-equilibrium (⁠$ξv0$⁠) where

and the characteristic vibrational relaxation time, typically compared with flow timescales, e.g.,

where τ_E is the eddy turnover time and τ_η is the Kolmogorov timescale. K_τ has been used to characterize turbulence-TNE interactions in statistically steady states.¹³ 

A large DNS database has been created over a range of flow conditions, which are summarized in Table I. As we will see below, the specific TNE parameters of the simulations allow us to understand the different regimes in which the interaction can be and explain the discrepancies observed between previous simulations¹⁹ and experiments.¹⁶ 

The decay of turbulence is accompanied by energy exchanges among different modes. Energy is transferred between the translational–rotational mode e and the vibrational mode $ev$⁠, through molecular collisions governed by Eqs. (2) and (3). In our simulations, the vibrational mode is initialized with excess energy. Thus, initially energy is transferred to the translational–rotational mode. For vibrationally cold flows, energy is transferred from the translational–rotational mode to the vibrational mode.

The third type of energy present here is the turbulent kinetic energy (⁠$K≡⟨ρuiui⟩/2$ with summation convention implied), which also interacts with the translational–rotational mode. The evolution equation for $K$ can be derived by taking a dot product of the momentum equation [Eq. (1)] with the velocity vector and averaging. In our flow configuration, homogeneity results, upon averaging, in vanishing transport terms leading to

where ϵ = u_j∂σ_ij/∂x_i is the viscous dissipation, which transfers irreversibly kinetic energy, $K$ to translational–rotational energy, ⟨E⟩. Because there are no terms in Eq. (9) that result in direct energy transfers between $K$ and ⟨$Ev$⟩, the dynamics of vibrational energy can only affect $K$ though an indirect path through the translational–rotational mode. This effect, as we show later, is indeed an important contributor to changes in the dynamics of the decay under certain conditions.

Turbulent kinetic energy also interacts with the translational–rotational mode by the so-called pressure–dilatation, the second term on the right hand side of Eq. (9). This mechanism, which can transfer energy reversibly between translational–rotational and kinetic energy, has been found to have a negligible effect in many instances.^28,29 We have also found this to be the case for simulations presented in this study. A schematic of the energy exchanges is provided in Fig. 1, which corresponds to the complete system

When turbulence decays, the exchange terms on the RHS of Eqs. (10) and (11) will act to decrease the degree of TNE and restore thermodynamic equilibrium by transferring energy from vibrational mode to translational–rotational mode through molecular collisions. The result would be an increase in ⟨T⟩ and a corresponding decrease in ⟨$Tv$⟩ (using its standard definition as the temperature that in Eq. (5) corresponds to a given value of $ev$⁠). The evolution of both ⟨T⟩ and ⟨$Tv$⟩ is shown in Fig. 2 for some selected cases. Both temperatures approach the same asymptotic value at TNE timescales (t ∼ O(⟨$τv0$⟩)) as the vibrational mode approaches equilibrium. Obviously longer times are needed to achieve equilibrium in flows with slow vibrational relaxation timescales.

In order to understand the role of turbulence in TNE decay, we compare the TNE decay process in a simplified model where no fluctuations are present, that is, an equivalent laminar model with the same mean initial conditions and properties. Since the mean velocity is zero, and no fluctuations exist, then TNE relaxation evolves according to

The above set of equations are highly nonlinear but can be solved numerically using an explicit forth-order Runge–Kutta method with variable time stepping for error control. The initial conditions are same as to those in the turbulent DNS simulations, i.e., $Ev0=⟨Ev0⟩,Ev0*=⟨Ev0*⟩,τv0=⟨τv0⟩$⁠.

The decay of vibrational non-equilibrium is plotted in Fig. 3 for turbulent flows (solid lines) at different conditions and the corresponding laminar model (dashed lines). $Ev^$ is the instantaneous degree of non-equilibrium as a fraction of initial non-equilibrium, i.e.,

Clearly TNE relaxation occurs at timescales $O∼(⟨τv0⟩)$ in both turbulent flows and the laminar model. We also see that the decay of TNE is faster for larger initial non-equilibrium energy, ξ₀ = ⟨$Ev0$⟩/⟨E₀⟩. This is not surprising given that for large values of ξ₀, a large amount of energy is transferred to ⟨E⟩, which creates a large increase in ⟨T⟩ and thus, a decrease in vibrational relaxation timescale ⟨$τv$⟩ [cf. Eq. (6)]. This effect is clearly seen in Fig. 6(b).

While the observations presented in Sec. III A apply to both laminar and turbulent flows, there are significant differences between them near equilibrium, that is, as $Ev^→0$⁠. In particular, an interesting observation is that initially vibrationally hot turbulent flows may become vibrationally cold $(Ev^<0)$ near the equilibrium as seen in the inset of Fig. 3(b). Note that this does not necessarily correspond to an overshoot of the equilibrium value at a spatially local point [i.e., $Ev(x,t)<Ev*(x,t)$] during the relaxation process. The negative values of $Ev^$⁠, instead, indicate that vibrationally cold regions become more dominant than vibrationally hot regions [i.e., $⟨Ev(t)⟩<⟨Ev*(t)⟩$] as we show momentarily. This effect interacts with turbulence, which near equilibrium at longer times is the dominant process as discussed by Donzis and Maqui.¹³ This regime is characterized with fluctuations in thermodynamic variables,²⁰ which locally create vibrationally hot or cold regions. The specific degree of TNE would then depend on the relative timescales of turbulent fluctuations driving TNE and molecular TNE relaxation, that is, K_τ. Indeed for TNE timescales comparable or slower than flow timescales, we observe overshoot of mean vibrational energy equilibrium energy.

This feature can be explained using the simple model represented by Eqs. (13) and (14), when spatial fluctuations and averages are accounted for in the evolution. Consider two realizations of this system of equations, one vibrationally cold and another vibrationally hot such that $⟨Ev0*⟩<⟨Ev0⟩$⁠, where the mean is trivially defined as the average of quantities at hot and cold realizations [e.g., $⟨Ev0⟩=(Ev0hot+Ev0cold)/2$]. The initial translational–rotational energy is identical for both.

There are in general two contributing effects in the evolution of $Ev$⁠: one due to changes in $Ev*$ and another one due to changes in $τv$ both on the right-hand side of Eq. (13). For the hot realization, $Ev$ decays toward $Ev*$ and in doing so, T increases due to the energy transfer. From Eq. (5), we see that $Ev*$ will increase and $Ev$ will decrease, a combination which tends to reduce the energy transfer rate in Eq. (14). On the other hand, as T increases, τ_v decreases [Eq. (6)], an effect that will tend to increase the energy transfer rate. From DNS data, we observe that it is the evolution of $τv$ that has a stronger effect in the energy transfer rate in vibrationally hot flows which, thus, results in a stronger $Qv$⁠. For the vibrationally cold realization, on the other hand, the relaxation would transfer energy from translational–rotational to the vibrational mode and T decreases. Thus, both effects [i.e., a decrease in $(Ev−Ev*)$ and increase in $τv$] contribute to decrease the rate at which energy is transferred. The overall effect is thus seen to be an asymmetry between hot and cold evolutions: relaxation in vibrationally cold zones becomes slower as the decay proceeds, in contrast to hot zones where relaxation becomes faster as the decay proceeds. The result of this dynamics is then a lower instantaneous average and a negative value for $⟨Ev⟩−⟨Ev*⟩$⁠. The evolution of the mean $Ev^$ for the simple model is presented in Fig. 4 for different values of $δ≡(Ev,hot−Ev*)/|Ev,cold−Ev*|$⁠, which is a measure of the degree of initial non-equilibrium. Decreasing δ represents an increased prevalence of vibrationally cold zones, which are slower to relax. Consistent with the discussion above, we see, in Fig. 4, an increasing overshoot at lower values of δ.

The applicability of this general argument for actual turbulent flows would obviously depend on whether the flow develops cold regions. This is indeed the case as can be seen in Fig. 5, where we show the probability density function (PDF) of $Y=Ev−Ev*$ at several representative times during the decay. At t = 0, the PDF presents nonzero values only for Y > 0, indicating the existence of vibrationally hot regions only. As time proceeds, however, we see that the PDF moves to the left. Beyond t/$τv$ ∼ O(1), we see the appearance of a negative tail. At this point, most of the initial TNE has decayed and the remaining TNE is sustained by turbulence.¹³ The PDF is then close to symmetric around positive and negative values, though the mean, while small in this case, is not zero.¹³ 

Fluctuations in thermodynamic variables due to turbulence also affect important TNE parameters such as $⟨Ev*⟩$ and ⟨$τv$⟩. Since both of these quantities are involved in the decay of TNE, it is of interest to assess turbulent effects especially for turbulent modeling purposes. Donzis and Maqui¹³ showed that turbulent fluctuations produce an increased storage of energy in vibrational mode. It is clear from Eqs. (5) and (6) that ⟨$τv$⟩ ≠ $τv$(⟨T⟩, ⟨ρ⟩) and $⟨Ev*⟩≠Ev*(⟨T⟩,⟨ρ⟩)$ due to their nonlinear dependence on their arguments. An expression for $⟨Ev*⟩$ can be derived¹³ using a Taylor expansion of $Ev*$ about ⟨T⟩, ⟨ρ⟩ as

with

where K_T = θ/⟨T⟩, h_T and h_ρT depend on molecular structure and mean temperature only. T^† ≡ (T − ⟨T⟩)/⟨T⟩ is the normalized temperature, such that $⟨T†2⟩=⟨(T−⟨T⟩)2⟩/⟨T⟩2$⁠. Note that the equilibrium vibrational energy also depends on density–temperature correlation (⟨ρT⟩), which depends on the dynamics of the flow, and ⟨ρ^†T^†⟩ is the corresponding mean normalized density–temperature correlation. The derived expression agrees very well with the present DNS data as shown in Fig. 6(a).

Similarly, one can obtain

where

depend on molecular structure (⁠$c2τv$⁠) and mean flow temperature. The fluctuations in density, $⟨ρ†2⟩$⁠, directly affect ⟨$τv$⟩ unlike $⟨Ev*⟩$ where h_ρ = 0. The role of turbulent fluctuations is again obvious through the temperature–density correlation where g_ρT is also a function of molecular structure and mean temperature only. The good agreement of this expression and the present DNS data is observed in Fig. 6(b).

The inset in both figures show directly the effect of turbulence as we plot the ratio of the average quantity to the equivalent laminar evaluation at the same mean conditions. Fluctuations contribute more significantly to the equilibrium vibrational energy than the vibrational relaxation time. The effect on the former is seen to be up to 7.5% for the conditions in the figure. The non-equilibrium observed when turbulence fluctuations dominate is due to the vibrational mode ⟨$Ev$⟩ lagging behind $⟨Ev*⟩$⁠. This effect and its consequences are discussed by Donzis and Maqui.¹³ It should also be noted that $⟨T†2⟩$ and ⟨ρ^†T^†⟩ grow strongly with the turbulent Mach number,^13,20 M_t. Thus, these corrections are expected also to grow in importance as compressibility levels increase.

In Sec. III, we investigated the relaxation of vibrational energy toward equilibrium as turbulence decays. Due to the energy exchanges present in this situation, however, it is also expected that the decay of turbulence itself be affected by the additional energy present in the system in the form of an initial vibrational excitation.

In fact, we observe that TNE does have a significant effect on the turbulence decay rate. In Fig. 7, we show the relative departure of kinetic energy in the presence of initial TNE from a decay without TNE for different levels of initial TNE measured by $ξ0=⟨Ev0⟩/⟨E0⟩$⁠. Note that all simulations start from identical initial conditions for velocity, density, and translational–rotational energy from a forced statistical steady state auxiliary simulation (Sec. II). We clearly see that higher levels of initial vibrational energy lead to faster decay of $K$⁠. This is consistent with the experimental observations in Ref. 16. We also see that while vibrational energy relaxes in timescales of O(⁠$τv$⁠), the effects of this process are apparent for very long times.

From the energy exchanges shown in Fig. 1, one sees that the decay of kinetic energy can be affected by two mechanisms, namely, dissipation and pressure dilatation. It has been observed that in a range of situations the latter is small.^21,28,29 Indeed we have verified this to be consistent with all the simulations presented here. For example, we computed the time integral of the ratio of average pressure–dilatation to dissipation over timescales of order $τv0$ and found this be only a few percentage points. The same result was obtained when we compute ensemble averages over a range of initial conditions with different levels of pressure-dilatation. Thus, it is expected that changes in dissipation are responsible for the distinct decays at different initial degrees of TNE. This is indeed supported by DNS below.

The behavior of dissipation is shown in Fig. 8 for varying values of ξ₀ [part (a)] and K_τ [part (b)]. From both figures, we observe a larger dissipation relative to a case without TNE shown as dashed lines. Furthermore, dissipation is higher for both large initial vibrational energies and fast vibrational relaxation timescales. We also see that the increased dissipation happens at different normalized times $t/τv0$⁠. As we will see later, the peak occurs at timescales comparable to the vibrational relaxation time (∼O(⁠$τv$⁠)), which can be faster or slower than hydrodynamic scales like the Kolmogorov timescale (τ_η) or the eddy turnover time (τ_E) depending of the molecular system under consideration and the flow conditions.

To understand this increase, we note that in fully developed turbulent flows, it is well-known that at high Reynolds numbers, dissipation becomes independent of viscosity,³⁰ a phenomenon called dissipative anomaly. More recently, this has also found support in compressible flows.²¹ However, dissipative anomaly is expected to hold when turbulence is in a state of statistical equilibrium where dissipation at the small scales balances turbulence production at the largest scales. If either one of these is quickly altered externally, the flow will in general evolve toward a new equilibrium where the balance between dissipation and production is restored. During this transition, dissipative anomaly is not expected to hold.

This is indeed the case here. The evolution of the flow configuration is such that energy deposition in the vibrational mode triggers, through the Landau–Teller term in Eqs. (10) and (11), an energy transfer to the translational–rotational mode in a timescale of order $τv$⁠. If $τv$ is relatively short compared to turbulence timescales, the translational–rotation temperature ⟨T⟩ will also increase in a timescale of order $τv$ while turbulence is essentially frozen. Thus, the temperature-dependent viscosity will increase too. Because turbulence (velocity gradients) has not time to adjust to the new viscosity in $τv$⁠, the dissipation ⟨ϵ⟩ = μ(2⟨s_ijs_ij⟩ − 2/3⟨(∂u_i/∂x_i)(∂u_j/∂x_j)⟩) will effectively follow μ as well. In other words, the decay of TNE increases ⟨T⟩ and thus viscosity μ which in turn increases dissipation. This is what we see in Fig. 9, where we show the evolution of viscosity normalized by its initial value for different vibrational energies. Viscosity increases in TNE timescales and becomes larger with larger initial vibrational energy. Further evidence of this mechanism is seen in the inset of Fig. 8(a), where we show the evolution of ⟨ϵ⟩/⟨μ⟩. If dissipation increases because of an increase in viscosity, this normalized quantity will not exhibit the peak seen in the main figure. This is indeed what is observed.

We note that the increase in dissipation shown in Fig. 8 is different from the increase in dissipation observed when simulations are started from an initial nonturbulent flow field (see, e.g., Fig. 12 in Ref. 29). The increase in dissipation in such simulations is due to the development of new scales in the flow as it approaches a turbulent state. Our simulations, on the other hand, are at a fully developed turbulent state at t = 0 (when TNE is introduced) where turbulent scales are fully developed. The increase in dissipation, as noted, is due to increase in the temperature-dependent viscosity caused by the decay of TNE.

The asymptotic value of viscosity can be calculated using global energy conservation which in the present system implies that the sum of all modes is constant during the decay, that is, $⟨E⟩+⟨Ev⟩+K=⟨E0⟩+⟨Ev0⟩+K0$⁠. At t → ∞, turbulent kinetic energy will vanish, and vibrational energy will tend to its equilibrium state at the final temperature, ⟨T_∞⟩. Correspondingly, the final translational–rotational energy is denoted by ⟨E_∞⟩. Thus, we can write

which simply states that the total change in temperature is due to both vibrational relaxation and decay of turbulent kinetic energy. Since ⟨E_∞⟩ and $⟨Ev∞*⟩$ depend on temperature, this expression is an implicit equation for the final temperature ⟨T_∞⟩. This is so because at t → ∞, there are no fluctuations of temperature or density, which implies ⟨E_∞⟩ = ⟨ρ_∞⟩⟨e_∞⟩ = ⟨ρ_∞⟩(R/γ − 1)⟨T_∞⟩ where averages, though unnecessary due to the absence of fluctuations, are kept for consistency. Similarly we can write $⟨Ev∞*⟩=⟨ρ∞⟩Rθv/(eθv/⟨T∞⟩−1)$⁠. Equation (23) can be readily solved numerically for ⟨T_∞⟩. Note that at very high temperatures if vibration is fully excited (i.e., ⟨$Ev$⟩ ≈ R⟨T⟩) both at t = 0 and t → ∞, then an analytical expression can be obtained for the final temperature, namely, $⟨T∞⟩≈(5/7)⟨T0⟩+2(⟨Ev0*⟩+K0)/7R$⁠. Finally, since μ ∝ T^1/2, at t → ∞ the increase of viscosity is then

Now consider the case of very rapid vibrational relaxation, that is, K_τ ≪ 1. In this case, turbulence can be considered frozen for times of the order of $τv$ in which the decay of $K$ is negligible. All excess vibrational energy is then transferred to translational–rotational energy. Conservation of energy can then be approximated as

where the prime in the subscript ∞′ stands for the intermediate asymptotic state t/$τv$ ≫ 1 but t/τ_E ≪ 1. This equation can again be solved numerically (or analytically if vibration is fully excited during the entire process) for $⟨T∞′⟩$⁠, which then yields the relative increase of viscosity as $⟨μ∞′⟩/⟨μ0⟩≈⟨T∞′⟩/⟨T0⟩1/2$⁠. Since $K$ is not accounted for here, this intermediate asymptotic temperature is lower than the t → ∞ asymptotic value, that is, $⟨T∞′⟩<⟨T∞⟩$⁠, which also implies $⟨μ∞′⟩<⟨μ∞⟩$⁠.

The two estimates for viscosity, namely, ⟨μ_∞⟩ and $⟨μ∞′⟩$⁠, are also are shown in Fig. 9 as dashed-dotted lines and dashed lines respectively. We see that for times of order $τv$⁠, viscosity approaches $⟨μ∞′⟩$ especially for small K_τ. For slow TNE relaxation (large K_τ), viscosity is expected to approach ⟨μ_∞⟩. Indeed, this is the case as shown in Fig. 9 for ξ₀ = 1.39 with K_τ = 0.07 (red) and K_τ = 94.8 (yellow). For the former, the increase in viscosity at early times is purely due to TNE decay and approaches the TNE limit, $⟨μ∞′⟩$⁠. For the latter, the increase in viscosity approaches ⟨μ_∞⟩, which results from the combined TNE-turbulence decay. The figure also contains an inset with the case ξ₀ = 0.1 and K_τ = 0.07, which exhibits an intermediate asymptotic state $⟨μ∞′⟩$ in TNE timescales and long-term asymptotic behavior in eddy turnover timescales, τ_E0 (inset).

In summary, consistent with experimental observations,¹⁶ TNE relaxation can affect the decay of turbulence. The main mechanism is found to be an increase in viscosity due to the energy transfer from vibrational to translational-rotational modes. This transfer depends on the relative timescales of TNE and turbulence as well as the initial degree of TNE. The governing parameters and the effect on the decay are discussed next, where we are also able to classify regimes in which turbulence is and is not affected by initial TNE.

In Sec. IV A, we have shown that the increase in viscosity due to TNE decay results in an increase in dissipation. We have also shown that excess vibrational energy $(⟨Ev0⟩−⟨Ev∞′*⟩)$ defines an asymptotic limit on this increase in viscosity $⟨μ∞′⟩$ and therefore in dissipation $⟨ϵ∞′⟩$⁠. This asymptotic increase in dissipation can now be related to the peak in dissipation ⟨ϵ_p⟩ at t = τ_p in Fig. 8. For small K_τ, TNE decay is much faster than turbulent processes in which case we expect all excess energy in vibration to quickly equilibrate yielding $⟨ϵp⟩≈⟨ϵ∞′⟩$⁠. For larger values of K_τ, TNE relaxation occurs simultaneously with turbulence decay. In such a case, we expect $⟨ϵp⟩<⟨ϵ∞′⟩$⁠. For t ≫ τ_p, the degree of non-equilibrium decreases and a classical turbulent decay dynamics becomes the dominant mechanism in which dissipation, as described earlier, decays over time. For very large values of K_τ, one may expect a very slow TNE relaxation relative to the turbulent energy dissipation rate decay—that is, the vibrational mode is frozen—in which case turbulence is only weakly affected and a peak may not even appear.

The conditions under which this occurs can be established by considering the behavior of dissipation at short times. For very small t, one can write a Taylor expansion of dissipation as

where b characterizes the decay rate of dissipation in time. For turbulence without TNE, dissipation decays monotonically and b > 0.

If TNE decay is very fast, the energy balance in Eq. (23) can be rewritten as

which, when divided by ⟨E₀⟩, reads

The second line was obtained by noting that $μ∼T∼E$⁠. This equation relates the initial excess of the energy in the vibrational mode to the increase in viscosity.

At the same time, Eq. (26) can be rewritten as

As explained above, dissipative anomaly does not hold for rapid changes, in which case, we have instead ⟨ϵ⟩/⟨ϵ₀⟩ ≈ ⟨μ⟩/⟨μ₀⟩. We can then evaluate Eq. (30) at t = t_p, and equate with Eq. (29), linearize the resulting equation, and solve for b to obtain

where we see that b < 0 for initially vibrationally hot flows (⁠$⟨Ev0⟩>⟨Ev∞′*⟩$⁠). This leads to an increase in dissipation due to TNE decay as explained. In this short-time analysis, this increase would clearly depend only on the nondimensional group observed in Eq. (31), which for convenience we define as

We see that β is comprised of two parts: an energy ratio, which is related to the maximum energy available that can be transferred to the translational–rotational mode, and a ratio of timescales, K_E = $τv$/τ_E, which determines the relative timescales at which this transfer occurs. The value of β depends solely on initial conditions and can be thus computed a priori. Consistent with the observations made for Fig. 8, the increase in dissipation depends not only on how much energy is available in the vibrational mode but also on how fast this energy is released.

We are thus interested in assessing whether β is indeed the parameter governing the increase in dissipation. In Fig. 10(a), we show the increase in dissipation (⟨ϵ_p⟩/⟨ϵ₀⟩) for the entire DNS database as a function of β. We observe peaks in dissipation (⟨ϵ_p⟩ > ⟨ϵ₀⟩) only when β > 1 reaching an asymptotic value at β > 200 (shown as horizontal dashed lines), which depends on the initial thermodynamic state (ξ₀) and excess energy in the vibrational mode (ξ_v0) according to Eq. (25). For small values of β, the vibrational energy available for increase in dissipation is small or K_E ≫ 1, which implies that turbulence is faster than TNE timescales and energy transfer from vibration to translational–rotational mode is slow. Then turbulence will largely remain unaffected and we observe no peak for β < 1 in Fig. 10(a). This does not however mean that TNE does not play a role. For β close to 1, for example, one still observes an increase in dissipation, which increases with initial vibrational energy. This increase, however, occurs over flow timescales and is not large enough to exceed the value of dissipation at initial conditions.

In Fig. 10(b), we show $ϵ^≡(⟨ϵp⟩−⟨ϵ0⟩)/(⟨ϵ∞′⟩−⟨ϵ0⟩)$⁠, the ratio of the difference between the peak dissipation from its initial value and the difference between its asymptotic peak dissipation $(⟨ϵ∞′⟩)$ and the initial dissipation for the entire database. If the energy exchange due to TNE relaxation is not strong enough, then dissipation will not increase sufficiently fast to compensate for the concurrent turbulent decay and dissipation is maximum at t = 0, which implies $ϵ^=0$⁠. If, on the other hand, all the excess vibrational energy is transferred to translational–rotational energy before turbulence evolves appreciably then $ϵ^=1$ as the peak in dissipation is equal to the asymptotic value, predicted by Eq. (29) with an $⟨ϵ∞′⟩/⟨ϵ0⟩≈⟨μ∞′⟩/⟨μ0⟩$⁠.

The substantial degree of universality observed in terms of characteristics of the diverse set of initial conditions of the flow (e.g., ξ₀, K_τ, M_t) suggests that β is indeed the governing parameter for the increased dissipation for all TNE and flow conditions. Figure 10(b) also identifies the asymptotic regime beyond β > 200 as mentioned earlier.

Of interest is also the time of occurrence of the peak (τ_p). This is shown in Fig. 10(c) relative to the hydrodynamic timescale τ_E. Because the peak in dissipation, if one exists, occurs at TNE timescales, τ_p/⟨τ_E⟩ is seen to scale approximately as K_E. For K_E > 1, the initial value of dissipation is the maximum, i.e., τ_p = 0, and turbulent processes largely govern the decay of dissipation in time. The instantaneous dissipation is still higher than cases without TNE. In this regime, TNE is slower than turbulent processes and leaves them largely unaffected (β < 1). Energy stays longer in vibrational modes and is released slowly to the translational–rotational mode. τ_p is seen to differ somewhat from $τv0$ due to latter becoming smaller as TNE relaxation proceeds. The additional scatter seen in the data at intermediate values of K_E seems to stem from second-order effects related to turbulent fluctuations.

The scaling proposed above is useful to understand other results in the literature. For example, Neville et al.¹⁹ conducted decaying simulations where the initial conditions correspond to complete thermal equilibrium. After a short initial period, the flow developed TNE due to difference in flow and TNE timescales, which will depend on mean and fluctuation levels of translational temperature.¹³ As turbulence decays, the mean temperature increases faster than the relaxation of TNE, which results in a vibrationally cold flow. They observed no increase in dissipation and thus no change in the decay of turbulent kinetic energy. Using the data available in Ref. 19, we can estimate that after a very short period of time, the simulations would have $|β|$ ∼ 10⁻³ at most. From Figs. 10(a) and 10(b), it is clear that at these conditions one would not observe any significant increase in dissipation consistent with their reported behavior.¹⁹ 

The experiments by Fuller et al.,¹⁶ on the other hand, consisted of a flow with initially vibrationally hot conditions, similar to the present study. They observed a faster decrease of $K$ when vibrational excitation was increased. In particular, they reported an increase in n in a classical power-law decay $K∼t−n$⁠, which can be used to estimate the additional dissipation produced by the introduction of vibrational TNE. Since $dK/dt≈−⟨ϵ⟩$⁠, it is readily shown that the ratio of dissipation for two different conditions at t = τ_p ≈ $τv$ [see Fig. 10(c)] is the ratio of exponents n. For their 300 W settings case, temperatures were T = 295 K and $Tv$ = 1550 K (corresponding to ξ₀ ≈ 0.5) and they measured n = 1.25 compared to n = 1 for the case without vibrational excitation. Thus, we estimate ⟨ϵ₃₀₀⟩/⟨ϵ₀⟩ ≈ n₃₀₀/n₀ ≈ 1.25, i.e., a 25% increase in dissipation. Using the initial vibrational and rotational–translational energy reported, we can use Eq. (29) to estimate precisely a 25% asymptotic increase in dissipation. However, the excellent quantitative agreement between the experimental measurement and our prediction is considered somewhat fortuitous given all the estimations and uncertainties involved in the comparison. Furthermore, using reported data, we estimate β ∼ O(1), which is not in the asymptotic regime. Thus, we consider the good comparison observed, though correct, only qualitative in nature.

We have investigated the interaction of turbulence with thermal non-equilibrium (TNE) using well-resolved direct numerical simulations (DNS). In particular, we focused on the two-way coupling between vibrational non-equilibrium and decaying turbulence, its governing parameters, and possible regimes of this interaction. For this, we generated a large DNS database with a wide range of parameters characterizing the TNE.

A general conclusion from this work is that the strength of the interaction and the degree in which TNE can affect the dynamics of turbulence depend on the initial degree of non-equilibrium and its relaxation time relative to hydrodynamic timescales. Therefore, the interaction is determined by the state of turbulence, its thermodynamic state, as well as the molecular structure of the fluid.

We showed that TNE decays at timescale of the order of its initial relaxation time. Larger initial energy in the vibrational mode leads to a faster effective decay rate of TNE. This relaxation toward equilibrium leads, through molecular collisions, to an increase of translational–rotational energy. Because of the temperature dependence of molecular transport coefficients, in particular viscosity, dissipation may increase temporarily leading to faster turbulence decay. From global energy conservation and classical scaling laws, we found the asymptotic increase in dissipation if all excess vibrational energy would relax relatively quickly.

We further showed that the governing parameter (β) contains a combination of the initial degree of non-equilibrium and the ratio of TNE and turbulence timescales. The effect of TNE on turbulence is determined by this parameter. When β is below unity (too little initial excess energy in vibrational modes or too long relaxation times), dissipation does not exhibit a peak. For β > 200, the asymptotic energy conversion is observed. These two cases would correspond, respectively, to a turbulence frozen case (fast TNE relaxation) and frozen TNE (slow TNE relaxation). The cases in between are seen to be characterized by a universal curve in β when normalized appropriately.

The identification of the governing parameter and possible regimes of the interaction helped explain the disagreements between simulations and experiments in the literature about the effect of TNE on decay rates: the specific conditions in these studies correspond to different values of β though some other parameters were matched.

Turbulence was also found to modify TNE relaxation. At long times compared to TNE, turbulence fluctuations were found to be the main driver of non-equilibrium consistent with the results for steady flows in Ref. 13. At short times, it was found that although the flow starts vibrationally hot everywhere, the flow can become vibrationally cold in the mean during a transient period. This new effect was explained to be due to the strongly nonlinear dependence of the relaxation timescale $τv$ on temperature, which leads to instantaneous and local relaxation rates that can be very different at different locations in the flow. The mechanism for this effect was illustrated by a simple system, which reproduces asymmetry observed between relaxation from vibrationally cold and hot regions.

We close by mentioning a potential application of the effects observed here. It is clear that vibrational excitation can lead to non-negligible changes in the fundamental turbulent cascade of energy. The rapid decay of turbulent fluctuations may, for example, be leveraged to reduce drag or heat exchange over surfaces in high-speed flows. The specifics of such a control mechanism by TNE is, in fact, part of ongoing research.

The authors gratefully acknowledge support from AFOSR (No. FA9550-17-1-0107) and the Hagler Institute for Advanced Study at Texas A&M University. The authors also thank R. Bowersox and S. North for helpful discussions.