Influence of vibrational non-equilibrium on the dynamics of velocity gradients in compressible flows

Abstract

Analysing the dynamics of velocity gradients is useful for understanding various nonlinear turbulence processes. This work focuses on how the vibrational non-equilibrium of the constituent molecules in a gaseous medium affects the dynamics of the velocity gradient and the pressure Hessian tensors. We first derive the exact evolution equation of the pressure-Hessian tensor in the presence of the vibrational non-equilibrium process. Subsequently, we perform several direct numerical simulations of compressible isotropic turbulence, including the vibrational relaxation process therein. Using flow fields extracted from these simulations, we conduct several parametric studies over a range of the Damköhler number (ratio of the relevant fluid time scale to that of the mean vibrational relaxation process) and the initial ratio of the vibrational temperature to the mean local temperature. We find that a variation in the initial Damköhler number does influence the evolution of the pressure-Hessian and the velocity gradient tensors. As the vibrational relaxation process becomes more rapid (an increase in the value of the initial Damköhler number), it causes a decrease in the strength of the pressure-Hessian tensor and simultaneous suppression of dilatational fluctuations in the flow field. On the other hand, a variation in the initial value of the ratio of the vibrational temperature to the local temperature does not seem to affect the pressure-Hessian or the velocity gradient tensor. These findings are expected to aid in the development of closure models for the pressure-Hessian tensor in compressible flows under vibrational non-equilibrium conditions.

1. Introduction

The velocity gradients and their temporal evolution in a turbulent flow hold the key toward deeper insights and improved understanding of several nonlinear turbulence phenomena like cascade (Suman & Girimaji 2011), material element deformation (Batchelor 1952; Girimaji & Pope 1990), intermittency (Li & Meneveau 2005, 2006), energy dissipation and scalar mixing (O'Neill & Soria 2005; Meneveau 2011). Experimental measurements (Tsinober, Kit & Dracos 1992) of velocity gradients are challenging. Experimentally, performing temporal tracking of the gradients is even more complex and rare (Xu, Pumir & Bodenschatz 2011). On the other hand, direct numerical simulations (DNS) can provide a realization of the raw Eulerian velocity field, which can then be used to assess the velocity gradient tensor. Still, developing simple (ordinary differential equation) dynamical models of the velocity gradient tensor is desirable. These models are computationally inexpensive compared with DNS, which require enormous computing resources and pose challenges in handling (post-processing) massive datasets. Further, these simple dynamical models provide a more natural framework in which the temporal evolution of the velocity gradients can be tracked following individual fluid elements (Lagrangian evolution), allowing for more insightful and physics-based investigations and analyses. Moreover, such dynamical models, whenever available, can directly serve as closure models for alternate methods of turbulence computations, such as the Lagrangian probability density function (PDF) methods (Pope 1985).

Given the advantages of such dynamical models of velocity gradients, much research has been directed toward their development and validation. The pioneering work in this regard is that of Vieillefosse (1982), who ignored the viscous effects and the unclosed, non-local anisotropic portion of the pressure-Hessian tensor. He presented the first approximate autonomous dynamical equation of the velocity gradients. This equation is referred to as the restricted Euler equation (REE). Based solely on an REE-based analysis, the author conjectured that the gradient steepening observed in turbulent flows must involve the tendency of the vorticity tensor to align well with the positive eigenvalues of the strain-rate tensor. Subsequently, Ohkitani & Kishiba (1995) used the REE model to provide insightful explanations about the alignment tendencies of the pressure-Hessian eigenvectors with the vorticity vector. With the subsequent availability of the DNS database, these conjectures have indeed been confirmed to be quite accurate (Ashurst et al. 1987; Kalelkar 2006).

Over the last three decades, the REE model has undergone several modifications. The missing viscous effects have been included using the linear diffusion model of Martın, Dopazo & Valiño (1998), and later, using the linear Lagrangian diffusion model of Jeong & Girimaji (2003), which could control the finite-time singularity problem but led to an over-amplified, exponential increase in the viscous effects. Recently Chevillard & Meneveau (2006) proposed their recent fluid deformation closure by which they encompass several essential physics, namely (i) inclusion of the anisotropic part of the pressure-Hessian tensor and (ii) a significantly improved representation of the viscous process. This not only addressed the finite-time singularity issue but also allowed the modified model to capture various statistical behaviours related to the velocity gradients and the pressure-Hessian tensor. Subsequently, Parashar, Suman & Srinivasan (2019) suggested alternate modelling for some parts of the viscous processes applicable to non-stationary flow fields.

Despite its current advanced form and maturity, the REE model is inherently unsuitable for compressible flows because of its dependence on Poisson's pressure equation. Following a thermodynamic-based approach that uses the continuity, energy and the state equations of a compressible medium Suman & Girimaji (2009, 2011), and later Danish, Suman & Srinivasan (2014), have proposed models for the dynamics of velocity gradients in compressible turbulence as well. However, despite these attempts, there is a vast scope for improvement of such models, especially for hypersonic flow fields wherein the impact of thermal non-equilibrium is also expected to be significant, over and above the effects of compressibility (Anderson 2000).

Flow fields past hypersonic vehicles, owing to the high-temperature, vibrational and chemical non-equilibrium (thermal non-equilibrium) of the flow, add further complexity to the turbulent flow field (Anderson 2000; Smits, Martin & Girimaji 2009). For air, vibrational excitations become significant beyond 800 K, and disassociation of oxygen and nitrogen molecules is probable beyond 2000 and 4000 K, respectively. As a consequence of these phenomena, the specific heat capacity values of air do not remain constant. How these affect turbulence is indeed of pertinence.

Recently, some fundamental DNS-based studies have been performed to understand the influence of thermal non-equilibrium on turbulence. Neville et al. (2014) studied the effects of thermal non-equilibrium processes on decaying isotropic turbulence, wherein they examined the influence of the Damköhler number on the flow field. Donzis & Maqui (2016) and Khurshid & Donzis (2019) performed a comparative study of isotropic turbulence in thermal equilibrium and non-equilibrium modes. In another recent study, Zheng et al. (2021) explained the transfer of internal energy fluctuation in isotropic turbulence under thermal non-equilibrium conditions. Even though these studies have contributed towards an improved understanding of the influence of thermal non-equilibrium on turbulence, to the best of the authors’ knowledge, no previous attempt has been made to investigate in detail the influence of thermal non-equilibrium directly on the evolution of the velocity gradient and the pressure-Hessian tensor. Such is the objective of this paper.

Specifically, we focus on understanding the unclosed pressure-Hessian tensor under the influence of vibrational non-equilibrium in the velocity gradient equation. Toward this goal, we identify the following tasks for this study:

1. (i) Derivation of the exact pressure-Hessian evolution equation in a vibrational non-equilibrium state of a gaseous medium.

2. (ii) Identification of new unclosed processes appearing in this evolution equation.

3. (iii) Performing DNS of homogeneous turbulence, including the phenomenon of vibrational non-equilibrium over a range of Damköhler number (ratio of the relevant fluid time scale to the time scale of vibrational relaxation) and the ratio of initial mean temperature to the mean vibrational temperature in the flow field.

4. (iv) Employing the generated DNS database to identify, isolate and understand the influence of the vibrational non-equilibrium phenomenon on the velocity gradient and the pressure Hessian fields.

This paper is organized as follows. In § 2, we derive the exact evolution equation of the pressure-Hessian tensor in the presence of vibrational non-equilibrium in the flow field. Section 3 identifies the closure problem associated with the pressure Hessian equation. In § 4, we present details of the DNS employed in this work. Section 5 discusses the vibrational relaxation process and the plan of study. Section 6 examines the influence of vibrational non-equilibrium on various aspects of the velocity gradient and the pressure-Hessian tensors. Section 7 examines the behaviour of the vibrational mechanism itself along with the influence of the initial Damköhler number and the ratio of vibrational temperature to local temperature. Finally, § 8 concludes the paper with a summary.

2. Velocity gradient dynamics in the presence of vibrational non-equilibrium

In this section, we derive the evolution equation of the velocity gradient ($\boldsymbol {A}$) and the pressure-Hessian ($\boldsymbol {P}$) tensors under vibrational non-equilibrium conditions. Mathematically, these tensors are defined as

(2.1a,b)\begin{equation} A_{ij}=\frac{\partial u_i}{\partial x_j} \quad \mathrm{and} \quad P_{ij}=\frac{\partial}{\partial x_j}\left(\frac{1}{\rho} \frac{\partial p}{\partial x_i}\right), \end{equation}

where $u_i$ and $x_i$ represent the $i$th velocity component and $i$th spatial coordinate in the Eulerian description of the flow field. The symbol $\rho$ represents density, and $p$ is pressure.

We start with the governing equations of mass and momentum for a viscous compressible medium

(2.2)$$\begin{gather} \frac{\partial\rho}{\partial t}+u_j\frac{\partial\rho}{\partial x_j}={-}\rho\frac{\partial u_j}{\partial x_j}, \end{gather}$$

(2.3)$$\begin{gather}\frac{\partial u_i}{\partial t}+u_j\frac{\partial u_i}{\partial x_j}={-}\frac{1}{\rho}\frac{\partial p}{\partial x_i}+\frac{1}{\rho}\frac{\partial\tau_{ji}}{\partial x_j}. \end{gather}$$

The viscous stress $(\tau _{ji})$ is expressed using the following relationship:

(2.4)\begin{equation} \tau_{ji}= 2 \mu\left[\frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}+ \frac{\partial u_j}{\partial x_i}\right)\right]+\delta_{ji}\zeta\frac{\partial u_k}{\partial x_k}, \end{equation}

where, $\delta _{ji}$ is the Kronecker delta symbol, $\mu$ denotes dynamic viscosity and $\zeta$ is the bulk viscosity of the fluid medium. Bulk viscosity ($\zeta$) is taken to be equal to $-(2/3)\mu$ (White 1979).

The governing equation of the total energy per unit mass $(e)$ is

(2.5)\begin{equation} \frac{\partial e}{\partial t}+ u_i\frac{\partial e}{\partial x_i} ={-}\frac{1}{\rho} \frac{\partial}{\partial x_i}(u_i p)+\frac{1}{\rho} \frac{\partial}{\partial x_j}(u_i\tau_{ij})-\frac{1}{\rho} \frac{\partial q_i}{\partial x_i}-\frac{1}{\rho}\frac{\partial q_{v,i}}{\partial x_i}, \end{equation}

where $e$ represents the total energy (per unit mass), which equals the sum of the kinetic energy, translational–rotational energy and vibrational energy ($e_v$) per unit mass. It is expressed as

(2.6)\begin{equation} e=\underbrace{\frac{u_iu_i}{2}}_{Kinetic\ energy}+\underbrace{\frac{5}{2}RT}_{Translational\text{-}rotational\ energy}+\underbrace{e_v}_{Vibrational\ energy}. \end{equation}

The symbol $T$ is the temperature of the gas, and $R$ is the gas constant (J kg$^{-1}$ K$^{-1}$). The symbols $q_i$ and $q_{v,i}$ are the thermal and vibrational heat fluxes (J s$^{-1}$ m$^{-2}$) in the $i$th direction.

In order to understand the influence of vibrational non-equilibrium on compressible flows, an additional equation of the evolution of vibrational energy per unit mass is introduced

(2.7)\begin{equation} \frac{\partial e_v}{\partial t}+u_i\frac{\partial e_v}{\partial x_i} ={-}\frac{1}{\rho} \frac{\partial q_{v,i}}{\partial x_i}+\frac{w_v}{\rho}. \end{equation}

This equation uses a source term ($w_v$) that governs the energy transfer between translational–rotational and vibrational modes. We employ the Landau–Teller model for the source term $(w_v)$ expressed as (Vincenti & Kruger 1982)

(2.8)\begin{equation} w_v ={-} \rho \frac{(e_{v}-e_{v}^*)}{\tau}, \end{equation}

where $e_{v}^*$ is the vibrational energy per unit mass in equilibrium with local temperature $(T)$

(2.9)\begin{equation} e_{v}^* = \frac{\dfrac{h\nu}{k_BT}}{\exp\left({\dfrac{h\nu}{k_B T}}\right)-1}RT. \end{equation}

The symbol $h$ represents Plank's constant, $\nu$ is the vibrational frequency of molecules of the gas and $k_B$ is the Boltzmann constant. The quantity $\tau$ appearing in (2.8) is the vibrational relaxation time (Park 1993)

(2.10)\begin{equation} \tau = \frac{1}{p} \exp[A (T^{{-}1/3} - B) - 18.42], \end{equation}

where $p$ is in atm. Coefficients $A$ and $B$ are calculated using (2.11) and (2.12), respectively,

(2.11)$$\begin{gather} A = 1.16 \times 10^{{-}3} \sqrt{M} \theta^{4/3}, \end{gather}$$

(2.12)$$\begin{gather}B = 0.015 M^{1/4}. \end{gather}$$

The symbol $M$ is the molecular weight of the gas in gram/mol, and $\theta$ is the characteristic vibrational temperature of the gas in $K$. The values of $A$ and $B$ for nitrogen $(N_2)$ are 221.46 and 0.029, respectively (Park 1993). In line with (2.9), we can express $e_{v}$ also as a function of new variable $T_v$, the so-called local vibrational temperature

(2.13)\begin{equation} e_{v} = \frac{\dfrac{h\nu}{k_BT_v}}{\exp\left({\dfrac{h\nu}{k_B T_v}}\right)-1}RT_v. \end{equation}

Finally, the thermodynamic variables like pressure ($p$), density ($\rho$) and temperature ($T$) are related algebraically through the state equation of a perfect gas

(2.14)\begin{equation} p =\rho RT. \end{equation}

Equation set (2.2)–(2.14) describes a compressible flow field under the vibrational non-equilibrium condition with density ($\rho$), pressure ($p$), temperature ($T$), velocity ($u_i$) and vibrational energy ($e_v$) as seven dependent variables.

To derive the exact evolution equation of the velocity gradient tensor $\boldsymbol {A}$ under a vibrational non-equilibrium condition, we take the gradient of (2.3) to arrive at the following equation:

(2.15)\begin{equation} \frac{{\rm D} A_{ij}}{{\rm D} t}={-A_{ik}A_{kj}}-{P_{ij}}+{\varGamma_{ij}}. \end{equation}

The symbol, ${{\rm D}}/{{\rm D} t}={\partial }/{\partial t}+u_{k}({\partial }/{\partial x_{k}})$ in (2.15) is the material derivative representing the change in $A_{ij}$ following the motion of the local small fluid element. The first term on the right-hand side is the self-stretching term, representing the inertial action on the local state of the velocity gradient tensor. The second term is the pressure-Hessian tensor ($\boldsymbol {P}$). The pressure-Hessian tensor represents the process by which density and pressure gradients directly influence the evolution of velocity gradients in compressible flows. The third term is the viscous transport term

(2.16)\begin{equation} \varGamma_{ij}=\frac{\partial}{\partial x_j}\left\{\frac{1}{\rho} \frac{\partial}{\partial x_k}\left[\mu\left(A_{ik} +A_{ki}- \frac{2}{3}A_{pp}\delta_{ki}\right)\right]\right\}. \end{equation}

Our interest is in having a simple dynamical equation set to represent the dynamics of velocity gradients. While (2.15) describes the exact evolution of the velocity gradient tensor ($\boldsymbol {A}$), it is mathematically unclosed. The unclosed terms are $P_{ij}$ and $\varGamma _{ij}$. In incompressible flows, pressure is completely governed by Poisson's pressure equation. However, in compressible flows, the acoustic time scale responsible for pressure evolution may not be significantly smaller than the fluid time scale, and thus pressure behaves like a $bona fide$ thermodynamic variable. Therefore, in compressible flows, it is desirable to have an evolution equation of the pressure-Hessian tensor ($\boldsymbol {P}$) as well. Suman & Girimaji (2011) derived the exact evolution equation of the pressure-Hessian tensor using the continuity, energy and the state equations for compressible flows with no reference to the vibrational excitation of the constituent molecules. Consequently, their evolution equation of the $\boldsymbol {P}$ tensor does not account for the effects of vibrational non-equilibrium in a flow field. In contrast, in this work, we aim to derive the exact evolution equation of the pressure Hessian for a compressible flow under the influence of the vibrational non-equilibrium process.

To derive such an evolution equation of the pressure Hessian, we first express the material derivative of the pressure Hessian as

(2.17)\begin{align} \frac{{\rm D} P_{ij}}{{\rm D} t} &= \frac{{\rm D}}{{\rm D} t} \left\{\frac{\partial}{\partial x_j}\left(\frac{1}{\rho}\frac{\partial p}{\partial x_{i}}\right)\right\} \nonumber\\ &= \frac{2}{\rho^3}\frac{{\rm D}\rho}{{\rm D} t} \frac{\partial\rho}{\partial x_j} \frac{\partial p}{\partial x_i} - \frac{1}{\rho^2} \frac{\partial p}{\partial x_i}\frac{{\rm D}}{{\rm D} t} \left(\frac{\partial\rho}{\partial x_j}\right)-\frac{1}{\rho^2} \frac{\partial\rho}{\partial x_j}\frac{{\rm D}}{{\rm D} t}\left(\frac{\partial p}{\partial x_i}\right) \nonumber\\ &\quad -\frac{1}{\rho^2} \frac{{\rm D}\rho}{{\rm D} t}\left(\frac{\partial^2p}{\partial x_j\partial x_i}\right) + \frac{1}{\rho}\frac{{\rm D}}{{\rm D} t}\left(\frac{\partial^2p}{\partial x_{j}\partial {x_i}}\right). \end{align}

Using (2.14) we write

(2.18)\begin{equation} \frac{{\rm D} p}{{\rm D} t}=\frac{{\rm D}\rho}{{\rm D} t}RT+\rho R\frac{{\rm D} T}{{\rm D} t}. \end{equation}

Next, we derive the material derivative of temperature. Equation (2.5) is the evolution equation of total energy $[e={u_iu_i}/{2}+\frac {5}{2}RT+e_v]$. The governing equation for kinetic energy is obtained by multiplying (2.3) with $u_i$. Subsequently, subtracting the evolution equation of kinetic energy and that of vibrational energy (2.7) from the total energy equation (2.5) leads us to the evolution equation of the temperature $(T)$

(2.19)\begin{equation} \frac{\partial T}{\partial t}+u_i \frac{\partial T}{\partial x_i}={-}\frac{2}{5}T \frac{\partial u_i}{\partial x_i}+\frac{2}{5}\frac{\tau_{ij}}{\rho R} \frac{\partial u_i}{\partial x_j}+\frac{2}{5}\frac{k}{\rho R} \frac{\partial}{\partial x_k}\left(\frac{\partial T}{\partial x_k}\right)-\frac{2}{5}\frac{w_v}{\rho R}, \end{equation}

where $k$ is the thermal conductivity of the medium (J s$^{-1}$ m$^{-1}$ K$^{-1}$), and it appears when the thermal flux $q_i$ is represented using Fourier's law of heat conduction: ($q_i=-k({\partial {T}}/{\partial {x}_i})$). Now, using (2.19) in (2.18) along with (2.2) leads to the evolution equation of pressure under the vibrational non-equilibrium conditions

(2.20)\begin{equation} \frac{\partial p}{\partial t}+u_i \frac{\partial p}{\partial x_i}={-}\frac{7}{5}p \frac{\partial u_i}{\partial x_i}+\frac{2}{5}\tau_{ij}\frac{\partial u_i}{\partial x_j}+ \frac{2}{5}\frac{k}{R}\frac{\partial}{\partial x_k}\left\{\frac{\partial}{\partial x_k} \left(\frac{p}{\rho}\right)\right\}-\frac{2}{5} w_v. \end{equation}

We subject (2.20) to the gradient operator successively to arrive at the equation $({{\rm D}}/{{\rm D} t})({\partial p}/{\partial x_i})$ and $({{\rm D}}/{{\rm D} t})({\partial ^2p}/{\partial x_i\partial x_j})$. Subsequently, employing (2.17), we finally arrive at the exact evolution equation of the pressure-Hessian tensor in the presence of vibrational non-equilibrium

(2.21)\begin{align} \frac{{\rm D} P_{ij}}{{\rm D} t} &= {-A_{kj}P_{ik}-A_{ki}P_{kj}-\frac{2}{5}A_{kk}P_{ij}} \nonumber\\ &\quad -\frac{7}{5\rho}\frac{\partial A_{kk}}{\partial x_j}\frac{\partial p}{\partial x_i}- \frac{7}{5\rho}\frac{\partial A_{kk}}{\partial x_i}\frac{\partial p}{\partial x_j}-\frac{1}{\rho} \frac{\partial A_{ki}}{\partial x_j}\frac{\partial p}{\partial x_k}+\frac{1}{\rho} \frac{\partial A_{kk}}{\partial x_j}\frac{\partial p}{\partial x_i} \nonumber\\ &\quad -\frac{7}{5}\frac{p}{\rho}\left(\frac{\partial^2A_{kk}}{\partial x_i\partial x_j}-\frac{1}{\rho} \frac{\partial A_{kk}}{\partial x_i}\frac{\partial\rho}{\partial x_j}\right) +\left(\frac{2}{5R}\right)\frac{\partial}{\partial x_j}\left\{\frac{1}{\rho} \frac{\partial^2}{\partial x_i\partial x_k}\left[k \frac{\partial}{\partial x_k} \left(\frac{p}{\rho}\right)\right]\right\} \nonumber\\ &\quad +\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho}\frac{\partial\tau_{km}}{\partial x_i}A_{km}+\frac{1}{\rho}\tau_{km}\frac{\partial A_{km}}{\partial x_i}\right] - \frac{2}{5} \frac{\partial}{\partial{x}_j}\left(\frac{1}{\rho}\frac{\partial{w}_v}{\partial{x}_i}\right). \end{align}

To make further progress, we must have a closure model for $w_v$. Employing the Landau–Teller model in (2.21), the pressure-Hessian equation modifies to

(2.22)\begin{align} \frac{{\rm D} P_{ij}}{{\rm D} t} &= \underbrace{-A_{kj}P_{ik}-A_{ki}P_{kj}-\frac{2}{5}A_{kk}P_{ij}}_{(I)} \nonumber\\ &\quad \underbrace{-\frac{7}{5\rho}\frac{\partial A_{kk}}{\partial x_j} \frac{\partial p}{\partial x_i}-\frac{7}{5\rho}\frac{\partial A_{kk}}{\partial x_i} \frac{\partial p}{\partial x_j}-\frac{1}{\rho}\frac{\partial A_{ki}}{\partial x_j} \frac{\partial p}{\partial x_k}+\frac{1}{\rho}\frac{\partial A_{kk}}{\partial x_j} \frac{\partial p}{\partial x_i}}_{(II)} \nonumber\\ &\quad \underbrace{-\frac{7}{5}\frac{p}{\rho}\left(\frac{\partial^2A_{kk}}{\partial x_i\partial x_j}- \frac{1}{\rho}\frac{\partial A_{kk}}{\partial x_i}\frac{\partial\rho}{\partial x_j}\right)}_{(III)} \underbrace{+\left(\frac{2}{5R}\right)\frac{\partial}{\partial x_j}\left\{\frac{1}{\rho} \frac{\partial^2}{\partial x_i\partial x_k}\left[k\frac{\partial}{\partial x_k} \left(\frac{p}{\rho}\right)\right]\right\}}_{(IV)} \nonumber\\ &\quad +\underbrace{\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho} \frac{\partial\tau_{km}}{\partial x_i}A_{km}+\frac{1}{\rho}\tau_{km} \frac{\partial A_{km}}{\partial x_i}\right]}_{(V)} \nonumber\\ &\quad +\underbrace{\frac{2}{5}\frac{\partial}{\partial x_j} \left[\frac{1}{\rho}\frac{\partial}{\partial x_i}\left(\frac{\rho{e}_v}{\tau}\right)\right]}_{(VI)} - \underbrace{\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho} \frac{\partial}{\partial x_i}\left(\frac{\rho{e}_v^{*}}{\tau}\right)\right]}_{(VII)}, \end{align}

where $e_v^*$ and $e_v$ are defined in (2.9) and (2.13), respectively.

To make this equation mathematically as well as physically more amenable to our understanding, we have gathered various terms as groups $(I-VII)$. We refer to these groups as mechanisms that influence the evolution of the pressure-Hessian tensor (in line with Suman & Girimaji 2011).

The first mechanism on the right-hand side involves interactions between our primary unknowns, the $\boldsymbol {A}$ and $\boldsymbol {P}$ tensors only. In line with the nomenclature for the $-A_{ik}A_{kj}$ term in (2.15), this mechanism is called the local stretching of the pressure Hessian field by the local velocity gradient field (Suman & Girimaji 2011).

The second mechanism involves higher-order gradients of the velocity gradient components and the first gradient of the pressure field. In contrast, mechanism $III$ involves terms in which higher-order gradients of $\boldsymbol {A}$ interact with the raw thermodynamic field $p/\rho$ itself. Suman & Girimaji (2011) argued that the presence of the raw thermodynamic variable $p/\rho$ (rather than its gradients) makes the underlying physics dependent on an acoustic time scale (and, accordingly, a Mach number). Indeed, the authors subsequently modelled the essential physics incumbent in this mechanism to recover the well-known Mach number effects from their enhanced homogenized Euler equation model of the compressible velocity gradient tensor. Mechanism $IV$ involves thermal conductivity $(k)$. This clearly represents the influence of thermal conduction on the pressure-Hessian evolution. On the other hand, mechanism $V$ involves the effect of viscous heating on the evolution of the $\boldsymbol {P}$ tensor. Mechanisms $I$ to $V$ were discussed by Suman & Girimaji (2011) as well. However, the authors examined a compressible flow field in the absence of any vibrational excitation. Their total energy equation involved changes in the kinetic energy and translational–rotational energy only. In contrast to mechanisms $I$ to $V$, mechanisms $VI$ and $VII$ are the ‘new’ mechanisms appearing in the evolution equation in the presence of vibrational non-equilibrium. These two ‘new’ mechanisms show the explicit influence of vibrational energy on the evolution of the dynamics of the pressure-Hessian field and, in turn, that on the velocity gradient field, as well. These mechanisms are attributable to the employment of the Landau–Teller model (2.8) for the interconversion of translational–rotational and vibrational energy modes. In subsequent sections, we wish to understand how these mechanisms affect the pressure-Hessian and the velocity gradient tensors.

3. The closure problem

We now examine (2.22) from the point of view of mathematical closure. This examination is motivated by the eventual goal of having a computationally viable set of ordinary differential equations for the $\boldsymbol {A}$ and $\boldsymbol {P}$ tensors. Among mechanisms $I$ to $V$, all are mathematically unclosed except for mechanism $I$. Suman & Girimaji (2011) and Danish et al. (2014) have already addressed this closure problem and proposed some models. In this paper, we focus on the ‘new’ mechanisms $VI$ and $VII$. These mechanisms are also mathematically unclosed for multiple reasons: the algebraic involvement of higher-order gradients of the velocity gradients, density, pressure, temperature and the gradients of vibrational temperature.

Danish et al. (2014) has earlier employed a DNS database of isotropic compressible turbulence to examine the relative importance of mechanisms $I$ to $V$ of the evolution equation of the $\boldsymbol {P}$ tensor for a compressible medium without any vibrational excitation. To identify the dominant mechanisms, they contracted the evolution equation of $P_{ij}$ and arrived at the evolution equation of $({{\rm D}}/{{\rm D} t})(P_{ij}P_{ij})$. All the original mechanisms were still present on the right-hand side of that equation; however, they appeared in their contracted scalar forms. The relative significance of these scalar forms was then compared by plotting the PDFs of the relative magnitudes of these mechanisms using a DNS database. In line with this approach, we, too, first contract equation (2.22) with ${Pij}$ and arrive at the scalar version of (2.22)

(3.1) \begin{align} \frac{1}{2} \frac{{\rm D}}{{\rm D} t}(P_{ij}P_{ij}) &= \underbrace{\left({-}A_{kj}P_{ik}-A_{ki}P_{kj}-\frac{2}{5}A_{kk}P_{ij}\right)P_{ij}}_{(I')} \nonumber\\ &\quad -\underbrace{\left(\frac{7}{5\rho}\frac{\partial A_{kk}}{\partial x_j} \frac{\partial p}{\partial x_i}+\frac{7}{5\rho}\frac{\partial A_{kk}}{\partial x_i} \frac{\partial p}{\partial x_j}+\frac{1}{\rho}\frac{\partial A_{ki}}{\partial x_j} \frac{\partial p}{\partial x_k}-\frac{1}{\rho}\frac{\partial A_{kk}}{\partial x_j} \frac{\partial p}{\partial x_i}\right)P_{ij}}_{(II')} \nonumber\\ &\quad -\underbrace{\frac{7}{5}\frac{p}{\rho}\left(\frac{\partial^2A_{kk}}{\partial x_i\partial x_j}- \frac{1}{\rho}\frac{\partial A_{kk}}{\partial x_i} \frac{\partial\rho}{\partial x_j}\right)P_{ij}}_{(III')}\nonumber\\ &\quad + \underbrace{\left(\frac{2}{5R}\right)\ \frac{\partial}{\partial x_j} \left[\frac{1}{\rho}\frac{\partial^2}{\partial x_i\partial x_k} \left\{k\frac{\partial}{\partial x_k}\left(\frac{p}{\rho}\right)\right\}\right]P_{ij}}_{(IV')} \nonumber\\ &\quad + \underbrace{\left(\frac{2}{5}\ \frac{\partial}{\partial x_j}\left[\frac{1}{\rho} \frac{\partial\tau_{km}}{\partial x_i}A_{km}+\frac{1}{\rho}\tau_{km} \frac{\partial A_{km}}{\partial x_i}\right]\right)P_{ij}}_{(V')} \nonumber\\ &\quad + \underbrace{\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho} \frac{\partial}{\partial x_i}\left(\frac{\rho{e}_v}{\tau}\right)\right]P_{ij}}_{(VI')} - \underbrace{\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho} \frac{\partial}{\partial x_i}\left(\frac{\rho{e}_v^{*}}{\tau}\right)\right]P_{ij}}_{(VII')}. \end{align}

In this equation, mechanisms $VI'$ and $VII'$ are the scalar counterparts of mechanisms $VI$ and $VII$ of (2.22). On the other hand, mechanisms $I'$ to $V'$ are the scalar counterparts of the non-vibration mechanisms $I$ to $V$ of (2.22). We intend to examine the relative importance of the collective action of the vibrational mechanisms $VI'$–$VII'$ to that of the inertial mechanisms $I', II', III'$. At this point, the physics related to the viscous processes ($IV'$ and $V'$) are not of interest. Thus, for a direct and unambiguous study, we define a fraction $f$

(3.2)\begin{equation} f = \frac{|VI' - VII'|}{| I' + II' + III'| + | VI' - VII'|}. \end{equation}

Examining $f$ can be useful to quantify the extent to which the evolution of the pressure Hessian is affected by the vibrational non-equilibrium process, if at all. The natural algebraic range of $f$ is $[0,1]$. A small (near-zero value) of $f$ would imply that the contribution of the vibrational mechanism is negligible in the governing equation of the pressure Hessian. On the other hand, a high magnitude of $f$ would imply a substantial contribution of the vibrational mechanism in governing and shaping the pressure-Hessian tensor in a compressible flow field.

4. Direct numerical simulations of isotropic turbulence with vibrational non-equilibrium

We have employed the hy2Foam solver of Casseau et al. (2016) and Casseau (2017) to perform DNS for our study. hy2foam is an open-source (Git 2022) solver that has been developed on the OpenFOAM platform and has the capability to handle both vibrational and chemical non-equilibrium. Even though the solver offers various turbulence models to perform Favre-averaged simulations, no turbulence model has been used to generate a well-resolved direct numerical simulation database. In this case, the solver solves the instantaneous continuity equation (2.2), velocity equation (2.3) and total energy equation (2.5). The state equation of a perfect gas is used as given by (2.14). Further, the evolution equation of vibrational energy ($e_v$) (2.7) is solved using the Landau–Teller model as presented in (2.8). The vibrational relaxation time scale appearing in (2.8) is calculated using the Millikan–White equation, as detailed in (2.10).

The computational domain is a cubical box with each side of length $2{\rm \pi}$. The computational grid is uniform, with the number of nodes being $192^3$/$256^3$/$512^3$. Periodic boundary conditions are imposed on the opposite sides of the cubic domain. A random initial velocity field is generated with zero mean and zero divergence. The initial spectrum of turbulent kinetic energy $E(\kappa )$ is

(4.1)\begin{equation} E(\kappa,t=0) = A_0\kappa^4\exp({-}2\kappa^2/\kappa^2_0), \end{equation}

where $A_0$ represents the spectrum constants. The symbol $\kappa$ is the magnitude of the wavenumber vector, and $\kappa _0$ is the $\kappa$ where $E(k)$ is maximum at $t=0$.

In homogeneous isotropic decaying turbulence, the Reynolds number based on the Taylor microscale $(Re_{\lambda })$ and initial turbulent Mach number ($M_t$) are defined as

(4.2a,b)\begin{equation} Re_{\lambda}=\sqrt{\frac{20}{3\eta\epsilon}}K \quad \text{and}\quad M_t=\sqrt{\frac{ 2K }{n R\langle T \rangle}}, \end{equation}

where $\eta$, $\epsilon$ and $K$ represent kinematic viscosity, dissipation rate and turbulent kinetic energy, respectively, and $\langle T \rangle$ is the volume-averaged temperature. The symbol $n$ represents the ratio of the specific heat capacity at constant pressure to that at constant volume.

In addition to $M_t$ and $Re_{\lambda }$, two more non-dimensional numbers characterize a simulation in the presence of vibrational non-equilibrium: (i) the initial Damköhler number $(D)$, and the initial ratio of the mean vibrational temperature to mean temperature $(\xi )$. These numbers are defined as in (4.3) and (4.4), respectively,

(4.3)$$\begin{gather} D = \frac{\tau_{LE}}{\tau}, \end{gather}$$

(4.4)$$\begin{gather}\xi=\frac{\langle{T_v}\rangle}{\langle{T}\rangle}, \end{gather}$$

where $\tau _{LE}$ is the large eddy-turnover time defined as $\tau _{LE} = \lambda /u^{\prime }$. The root mean squared value of the fluctuating velocity, denoted as $u'$, is defined as

(4.5)\begin{equation} u' = \sqrt{\left\langle \frac{u_1^2 + u_2^2 + u_3^2}{3} \right\rangle}, \end{equation}

where $\langle \text { }\rangle$ represents the volume-averaged value of the quantity of interest (equivalently, the mean of the quantity), and $\lambda$ is the Taylor microscale, equal to (Honein & Moin 2004)

(4.6)\begin{equation} \lambda = \sqrt{\frac{u'^2}{\langle(\partial u_1/\partial x_1)^2 \rangle}}. \end{equation}

The definition of the Damköhler number ($D$) used in this study is defined in line with the similar dimensionless number used in chemically reacting flows to relate the characteristic time scale of chemical reactions to that of the advective transport phenomenon (or the flow time scale) rate occurring therein (Fogler 2006). As the initial Damköhler number becomes large, the disparity between the fluid time scale and the vibrational relaxation time scale increases, and an initial vibrational non-equilibrium state will more rapidly tend to approach its equilibrium state. In the limit of the Damköhler number being zero, the vibrational state of the flow will remain effectively frozen even though the flow itself evolves. In this study, our objective is to examine the influence of this parameter having a value in the moderate range of 0.25 to 4, wherein both turbulence and vibrational relaxation processes are expected to evolve substantially.

The range of $\xi$ used in our simulation set spans from 0.125 to 0.75. In hypersonic aerospace applications, a shock upstream of a solid obstruction typically causes a substantial disparity in the vibrational temperature values upstream and downstream of the shock front. The translational energy rapidly attains its Maxwellian equilibrium state immediately downstream of the shock. This increases the temperature immediately downstream of the shock front. However, the vibrational temperature after the shock might remain at its pre-shock value. Downstream of such a shock front, the ratio of vibrational temperature to the local temperature is less than unity. Indeed, our chosen range for $\xi$ reflects such typical and realistic scenarios that initiate a vibrational relaxation process alongside the evolution of the background flow field.

In table 1, we present the list of simulation cases employed in our study. Simulation A has identical initial conditions as employed by Honein & Moin (2004) and is used to validate the hy2FOAM solver against the results of Honein & Moin (2004). Like Honein & Moin (2004), in simulation A, the viscosity dependence on temperature is modelled using a power law with the exponent being $0.76$. This simulation does not involve any vibrational non-equilibrium process. On the other hand, simulation cases B1–B4 and C–F all have the vibrational relaxation process included in the simulation. Simulation cases B1–B4 have been performed for the most ‘severe’ combination of the initial parameters employed in our study. For all these simulations, the turbulent Mach number, the Reynolds number, the Damköhler number and $\xi$ have been set as 0.5, 60, 4 and 0.125, respectively. While B1–B3 differ only in terms of the grid size ($192^3, 256^3$ and $512^3$), case B4 differs from B2 only in terms of the Courant–Friedrich–Lewy (CFL) number.

Table 1. The DNS employed in this study.

Cases B2, C and D are employed to evaluate the influence of the initial Damköhler number on the physics of interest since these sets of simulations differ only in terms of that parameter. The initial turbulent Mach number and the Taylor microscale-based Reynolds number are 0.5 and 60, respectively, for these simulations. The initial value of $\xi$ equals 0.5. On the other hand, cases B2, E and F differ in terms of $\xi$. Simulation case G has its settings identical to B2, but has the vibrational relaxation processes switched off. This simulation is used as the baseline case with which the results from the simulation cases involving vibrational relaxation processes are compared.

In all these simulations (B1–B4 and C–G) the gaseous medium is chosen to be nitrogen $(N_2)$. In these simulations, we have kept the coefficient of viscosity ($\mu$) a constant (as determined by the initial Reynolds number). This is a deliberate choice because we wish to isolate and identify the direct influence of the variation in the Damköhler number and $\xi$ on the velocity gradient and the pressure-Hessian dynamics. We wish to mention that our simulations are different from those of Khurshid & Donzis (2019), who not only use a temperature-dependent viscosity coefficient, but also allow the flow field to decay from an initial state which has been extracted from an already fully developed state of forced turbulence. We plan to include a variable viscosity coefficient in our future work.

Subject to the availability of the library of discretization schemes in the hy2Foam package, we have employed the second-order accurate central-upwind scheme developed by Kurganov, Noelle and Petrova (Kurganov, Noelle & Petrova 2001; Greenshields et al. 2010) for the advection term, whereas the diffusive fluxes utilize the central difference method. For temporal discretization, a second-order accurate, backward implicit formulation is used. Uniform time stepping is used in all the simulations. The interpolation schemes employed are of second-order and unbounded. The value of the time step is chosen to be such that, initially, the CFL number is in the range of 0.00625–0.025 (case-specific values are listed in table 1).

Before discussing the results related to the vibrational non-equilibrium phenomenon, we first assess the capability of the hy2Foam solver in simulating standard compressible isotropic turbulence (case A). We compare the results from this simulation against the results of Honein & Moin (2004). In figure 1(a–c), we present comparisons in terms of the temporal evolution of appropriately normalized turbulence kinetic energy ($K$), root-mean-square (r.m.s.) pressure fluctuations ($p'$) and the r.m.s. of the temperature ($T'$). The normalizing quantities in figure 1(a–c) are the initial value of the turbulent kinetic energy, the initial values of $(np_oM_t^2)$ and the initial values of $(n-1)(T_oM_t^2)$, respectively. The subscript $o$ denotes the initial value of the respective quantity. In each case, the time axis is normalized by the eddy-turnover time ($\tau _{LE} = {\lambda }/{u'}$) calculated based on the initial flow field. The initial condition set for these simulations corresponds to Honein & Moin (2004), wherein all thermodynamic fluctuations are zero initially. This choice explains the steep variations in the thermodynamic quantities during the first few time steps of the simulations. In figure 1(d), we present the evolution of skewness ($S_u$) from case A

(4.7)\begin{equation} S_u=\frac{\left\langle \left[ \left( \dfrac{\partial u_1}{\partial x_1} \right)^3 + \left(\dfrac{\partial u_2}{\partial{x}_2} \right)^3 + \left(\dfrac{\partial u_3}{\partial{x}_3} \right)^3 \right] / 3 \right\rangle}{\left\langle \left[ \left( \dfrac{\partial u_1}{\partial{x}_1} \right)^2 + \left(\dfrac{\partial u_2}{\partial{x}_2} \right)^2 + \left(\dfrac{\partial u_3}{\partial{x}_3} \right)^2 \right] / 3 \right\rangle^{3/2}}. \end{equation}

The included reference results in figure 1(a–d) are the digitized versions of the data available in Honein & Moin (2004) and Chen et al. (2020) (for $S_u$). Evidently, the results of the case A simulation agree quite well with the reference DNS results.

Figure 1. Simulation case A. Temporal evolution of (a) turbulent kinetic energy, (b) pressure fluctuations, (c) temperature fluctuations and (d) skewness. Symbols represent the results of Honein & Moin (2004) and Chen et al. (2020).

In figures 2(a–c) and 3(a–c) we present results from simulation cases B1–B3. Figure 2(a–c) shows the temporal evolution of the turbulent kinetic energy and the thermodynamic fluctuations, whereas figure 3(a–c) shows the temporal evolution involving different aspects of the velocity gradient field: skewness ($S_u$), r.m.s. of vorticity ($\omega ^{'}$) and r.m.s. of dilatation rate ($\theta ^{'}$)

(4.8)$$\begin{gather} \omega^{'}={\langle \omega_i \omega_i \rangle}^{1/2}, \end{gather}$$

(4.9)$$\begin{gather}\theta'=\left\langle \left( \frac{\partial u_i}{\partial x_i} \right)^2 \right\rangle^{1/2}, \end{gather}$$

where $\omega _i$ represents the $i$th component of the instantaneous vorticity vector. All these time-dependent statistics have been normalized by their respective initial values. The difference between the results obtained on the $256^3$ grid and those obtained on the $512^3$ grid is small. The maximum error percentage of $\omega ^{'}$ between the $256^3$ grid and the $512^3$ grid is less than $3.5\,\%$ over the entire duration of $4$ eddy-turnover times. Similarly, the maximum error percentage of $\theta ^{'}$ between the $256^3$ grid and the $512^3$ grid is less than $3.0\,\%$ over the entire duration of $4$ eddy-turnover times. Based on these results, we deem the results of simulation B2 ($256^3$) to be adequately grid independent. Thus, other simulations (cases C–G) have also been performed on a $256^3$ grid. We have performed a similar grid-independence study corresponding to the initial conditions of the simulation case G, as well. Even though those results are not included in the paper, the $256^3$ grid is found to be adequate there too.

Figure 2. Grid independence. Temporal evolution of (a) turbulent kinetic energy, (b) r.m.s. of density, (c) r.m.s. of temperature fluctuations in simulation cases B1, B2 and B3.

Figure 3. Grid independence: (a) skewness, (b) vorticity fluctuations and (c) dilatational fluctuations in simulation cases B1, B2 and B3.

We have performed a time-step-independence study using the results of the simulation cases B2 and B4. In simulation B4, the employed time step is twice that of the time step used in B2. The corresponding difference in the results (turbulent kinetic energy, $p^{'}$, $T^{'}$, $S_u$, $\omega ^{'}$ and $\theta ^{'}$) from these two simulations is found to be insignificant (results not included here to avoid repetition). Thus, the time step employed in all other simulations of table 1 is chosen to be the same as employed in the simulation case B2.

In figure 4, we present the evolution of skewness for simulations B2, C, D, E, F and G (4.7). We observe that, beyond $t/\tau _{LE}>1.5$, in every case, the skewness stabilizes to a value close to $-$0.45, indicating that, beyond that time, the nonlinear turbulent processes are in full effect. In figure 5 we present the decay of $Re_{\lambda }$ and turbulent Mach number ($M_t$, 4.2a,b) in various simulations that differ in terms of initial Damköhler number and the initial value of $\xi$. Evidently, there is no influence of these initial parameters on the decay of $Re_\lambda$ and $M_t$.

Figure 4. Temporal evolution of skewness for simulation cases (a) B2, C, D and G and (b) B2, E, F and G.

Figure 5. Temporal evolution of (a) Taylor micro-scale Reynolds number and (b) turbulent Mach number in simulation cases B2, C, D, E and F.

5. The vibrational relaxation process and plan of study

To visualize the vibrational relaxation process in our simulations, we first present the temporal evolution of the statistics of temperature and the vibrational temperature. In figure 6(a) we present the temporal evolutions of mean temperatures ($\langle T \rangle$ and $\langle T_v \rangle$) in the sets of simulations B2, C and D, which differ in terms of the initial Damköhler number. In figure 6(b) we present the temporal evolutions of the r.m.s. of the two temperature quantities ($\sqrt {\langle T^{'}T^{'} \rangle } \text { and } \sqrt {\langle T_v^{'}T_v^{'} \rangle }$) in the same set of simulations. In figure 7(a) we present the temporal evolutions of mean temperatures ($\langle T \rangle$ and $\langle T_v \rangle$) in the sets of simulations B2, E and F, which differ in terms of the initial $\xi$. In figure 7(b) we present the temporal evolutions of r.m.s. of fluctuations in the two temperature quantities ($\sqrt {\langle T^{'}T^{'} \rangle } \text { and } \sqrt {\langle T_v^{'}T_v^{'} \rangle }$) in the same set of simulations. The mean temperature profiles have been non-dimensionalized with $T_o$, whereas the r.m.s. values of fluctuations have been normalized with the factor $(n-1)T_o{M}_t^2$.

Figure 6. (a) Evolution of mean temperature (with symbols) and mean vibrational temperature (without symbols) in simulation cases B2, C and D. (b) Evolution of r.m.s. of temperature (with symbols) and that of the vibrational temperature (without symbols) in simulation cases B2, C and D.

Figure 7. (a) Evolution of mean temperature (with symbols) and mean vibrational temperature (without symbols) in simulation cases B2, E and F. (b) Evolution of r.m.s. of temperature (with symbols) and that of the vibrational temperature (without symbols) in simulation cases B2, E and F.

In all these simulations (figures 6(a) and 7(a)), it is evident that the mean translational–rotational energy (represented by $\langle T \rangle$) of the flow is diverted to the mean vibrational energy (represented by $\langle T_v \rangle$). Consequently, the mean vibrational temperature increases with a simultaneous decrease in the mean temperature. In simulation B2, ($D=4$), by the end of 1.5 eddy-turnover times, $\langle T \rangle$ and $\langle T_v \rangle$ become almost equal, indicating the achievement of a vibrational equilibrium state for the mean flow field. However, this does not necessarily mean that the fluctuating flow field has also attained vibrational non-equilibrium. In figures 6(b) and 7(b), we observe that the r.m.s. values of the two temperatures tend to differ substantially from each other even at later times (up to 5 eddy-turnover times).

The new vibrational mechanisms appearing in (2.22) can be expressed in terms of the instantaneous and fluctuating energies. To understand the significance of contributions to these mechanisms from the mean and the fluctuating flow fields separately, we perform the Reynolds decomposition of $e_v$ and $e_v^{*}$

(5.1a,b)\begin{equation} e_v=\langle e_v \rangle + e_v^{'} \quad \text{and}\quad e_v^{*}=\langle e_v^{*}\rangle + e_v^{*'}. \end{equation}

Substituting these in the vibrational mechanism $VI^{'}$–$VII^{'}$ in (2.22) gives

(5.2)\begin{align} & \frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho}\frac{\partial}{\partial x_i} \left( \frac{\rho}{\tau}e_v\right)\right]-\frac{2}{5}\frac{\partial}{\partial x_j} \left[\frac{1}{\rho}\frac{\partial}{\partial x_i}\left( \frac{\rho}{\tau}e_v^{*}\right)\right] \nonumber\\ &\quad =\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho}\frac{\partial}{\partial x_i} \left( \frac{\rho}{\tau}\langle{e_v}\rangle+\frac{\rho}{\tau}e_v^{'}\right)\right] -\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho}\frac{\partial}{\partial x_i} \left(\frac{\rho}{\tau}\langle{e_v^{*}}\rangle+\frac{\rho}{\tau}e_v^{*'}\right)\right], \end{align}

which, upon regrouping, leads to

(5.3)\begin{align} & \frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho}\frac{\partial}{\partial x_i} \left(\frac{\rho}{\tau}e_v\right)\right]-\frac{2}{5}\frac{\partial}{\partial x_j} \left[\frac{1}{\rho}\frac{\partial}{\partial x_i}\left( \frac{\rho}{\tau}e_v^{*}\right)\right] \nonumber\\ &\quad =\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho}\frac{\partial}{\partial x_i} \left( \frac{\rho}{\tau}\langle{e_v}\rangle-\frac{\rho}{\tau} \langle{e}_v^{*}\rangle\right)\right]+ \frac{2}{5} \frac{\partial}{\partial x_j}\left[\frac{1}{\rho}\frac{\partial}{\partial x_i} \left(\frac{\rho}{\tau}{e_v^{'}}-\frac{\rho}{\tau}e_v^{*'}\right)\right]. \end{align}

Equation (5.3) shows that the vibrational mechanisms take two kinds of contributions: (i) from the disparity of $\langle {e_v}\rangle -\langle {e_v^{*}}\rangle$, and (ii) from the disparity of ${e_v^{'}}-{e_v^{*'}}$. To gauge the temporal variation in the relative importance of the two contributions, we define, with reference to (5.3), the following fraction ($\chi$):

(5.4)\begin{equation} \chi = \frac{| F_{ij}F_{ij}|}{|F_{ij}F_{ij}| + | M_{ij}M_{ij}|}, \end{equation}

where

(5.5)\begin{equation} M_{ij}=\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho} \frac{\partial}{\partial x_i}\left(\frac{\rho}{\tau}\langle{e_v}\rangle- \frac{\rho}{\tau}\langle{e}_v^{*}\rangle\right)\right], \end{equation}

and

(5.6)\begin{equation} F_{ij}=\frac{2}{5}\frac{\partial}{\partial x_j}\left[\frac{1}{\rho} \frac{\partial}{\partial x_i}\left( \frac{\rho}{\tau}{e_v^{'}}-\frac{\rho}{\tau}e_v^{*'}\right)\right]. \end{equation}

Clearly, at a given location in the flow field, the quantity $M_{ij}M_{ij}$ represents a scalar measure of vibrational disparity in the mean flow, whereas $F_{ij}F_{ij}$ represents a scalar measure of vibrational disparity in the fluctuating flow in the context of the vibrational mechanisms $VI$–$VII$. With the initial thermodynamic fluctuations in the flow field being zero in our simulations, the value of $\chi$ is not well defined at that instant. However, as the flow field evolves, $\chi$ starts assuming meaningful values. An increasing value of $\chi$ at a location in the flow field would indicate the reducing significance of the disparity in the mean flow. In figures 8 we present the evolution of the volume-averaged (or the mean) value of $\chi$. In figure 8(a) we present the results from simulations B2, C and D, which differ only in terms of the initial Damköhler number. In figure 8(b) we present the results from simulations B2, E and F, which differ only in terms of the initial value of $\xi$. We observe that, in both sets of simulations, $\langle \chi \rangle$ grows and tends to reach the asymptotic value of unity. In the first set of simulations (B2, C and D), the approach towards unity is most rapid for the case with $D=4$ and is the slowest for the case with $D=0.25$. For the $D=4$ case, $\langle \chi \rangle$ reaches its asymptotic state within one eddy-turnover time, which coincides with the mean vibrational temperature equilibrating with the mean temperature (figure 6a). Similarly, for other simulations, the time instant when $\langle \chi \rangle$ reaches its asymptotic state coincides with the mean vibrational temperature equilibrating with the mean temperature. Furthermore, we observe that, in figure 8(a), $\langle \chi \rangle$ crosses the value of $0.5$ quite early (within two eddy-turnover times) in all the simulations. This is clear evidence that, beyond two eddy-turnover times, in all three simulations (B2, C and D), the disparity in the fluctuating flow field (represented by $F_{ij}F_{ij}$) is actually the major contributor to the vibrational mechanisms ($VI$–$VII$). The same conclusion emerges from figure 8(b), wherein we have included the results of the simulations B2, E and F, which differ only in terms of the initial value of $\xi$. In all these cases, $\langle \chi \rangle$ reaches its asymptotic state at the same time (in all cases, the initial Damköhler number is identical). Evidently, beyond one eddy-turnover time, the disparity in the fluctuating flow field is almost the sole contributor to the vibrational mechanisms ($VI$–$VII$) in the pressure Hessian equation (2.22). These results demonstrate that the disparity in the fluctuating field persists for a long time in decaying turbulence (long after the mean vibrational temperature has equilibrated with the mean temperature in decaying turbulence). This, in turn, means that the vibrational mechanism would remain active throughout the simulation time even when the kinetic energy has decayed to a much lower level compared with its initial state.

Figure 8. Temporal evolution of $\chi$. (a) Influence of Damköhler number. (b) Influence of $\xi$.

Our goal is to examine if the presence of the vibrational non-equilibrium phenomenon affects the velocity gradients and the pressure-Hessian tensors at all in a compressible flow field. Accordingly, in §§ 6 and 7, we examine the influence of the Damköhler number on these tensors. Specifically, we examine (i) the principal strain-rate components (eigenvalues of the strain-rate tensor), (ii) the magnitude of the velocity gradient tensor and (iii) the magnitude of the pressure-Hessian tensor relative to that of the velocity gradient tensor with varying Damköhler number. Flow fields from simulation cases B2, C, D and G are employed for this study.

Subsequently, in § 7, to gain deeper insight into the exact mechanisms by which the vibrational non-equilibrium phenomenon influences the evolution of the pressure Hessian tensor (and in turn the velocity gradient tensor itself), we perform a detailed parametric study of the ‘new’ mechanisms $VI$–$VII$ appearing in the governing equation of the pressure Hessian tensor (2.22). Here, we examine the influence of the initial Damköhler number $\xi$ on these mechanisms. Further, we also discuss the influence of the sign of the local dilatational state on the vibrational mechanism.

6. Influence of vibrational non-equilibrium on the velocity gradient and pressure-Hessian tensors

In compressible turbulence, the scalar aspects of the fluctuating velocity gradient tensor are often examined in terms of the so-called solenoidal dissipation ($\epsilon _s$) and the dilatational dissipation ($\epsilon _d$) of the turbulent kinetic energy

(6.1a,b)\begin{equation} \epsilon_s=\langle\mu\omega_i\omega_i\rangle,\quad \epsilon_d = \left\langle \frac{4\mu}{3} \left(\frac{\partial u_i}{\partial x_i} \right)^2 \right\rangle, \end{equation}

where $\omega _i$ represents the $i$th component of the instantaneous vorticity vector and $({\partial u_i}/{\partial x_i} )$ is the velocity divergence. While, in general, the variation in $\mu$ could also alter the temporal variations of $\epsilon _s$ and $\epsilon _d$, for our simulations (B2, C, D, E, F and G) it is held constant. Thus, the temporal changes in $\epsilon _s$ and $\epsilon _d$ must be directly attributable to the changes in $\omega ^{'2}$ and $\theta ^{'2}$. Accordingly, we directly examine the temporal evolutions of $\omega ^{'}$ and $\theta ^{'}$ in these simulation cases to identify and understand the influence of the vibrational parameters $D$ and $\xi$.

In figure 9(a,b), we present the temporal evolutions of $\theta ^{\prime ^2}$ and $\omega ^{\prime 2}$ from the simulation cases B2, C and D (which differ in terms of the initial Damköhler number alone). In figure 9 (c,d), we present the temporal evolutions of $\theta ^{\prime ^2}$ and $\omega ^{\prime 2}$, from the simulation cases B2, E and F (which differ in terms of the initial ratio of vibrational temperature to local temperature ($\xi$) alone). For reference, we also present the results from simulation case G, which has the same initial Mach number and Reynolds number but no vibrational relaxation process present therein. For each curve, the instantaneous variable has been normalized by the initial value of $\omega ^{\prime 2}$.

Figure 9. Temporal evolutions in simulation cases B2, C, D and G; (a) $\theta ^{'2}$ and (b) $\omega ^{'2}$. Temporal evolutions in simulation cases B2, E, F and G;(c) $\theta ^{'2}$ and (d) $\omega ^{'2}$.

We observe that the influence of the Damköhler number (over the range of this parameter considered in this work) on $\omega ^{\prime 2}$ is negligible. In contrast, the initial Damköhler number does affect the evolution of $\theta ^{\prime 2}$. We observe that, as the initial Damköhler number increases, the magnitude of $\theta ^{\prime ^2}$ decreases, and this difference tends to become larger at later times. In contrast, figure 9(c,d) shows that changes in the initial value of $\xi$ influence neither the statistics of dilatational fluctuations ($\theta ^{\prime ^2}$) nor the vortical fluctuations ($\omega ^{\prime ^2}$).

It is well known that, in an isotropic turbulent flow field, the probability density function (PDF) of $A_{ii}$ shows a single peak distribution centred at zero. However, the distribution becomes wider in a simulation with a higher initial turbulent Mach number (Suman & Girimaji 2013). This tendency indicates the existence of a significant number of fluid elements having non-zero $A_{ii}$. While ${A_{ii}}<0$ indicates contracting fluid elements, $A_{ii}>0$ indicates expanding fluid elements. We wish to examine if such a distribution of $A_{ii}$ is affected by the presence of vibrational non-equilibrium processes. In figure 10, we present the PDFs of $A_{ii}$ extracted from simulations B2, C, D and G at a representative time of $t/\tau _{LE}=3.5$. As the vibrational process becomes more rapid (represented by higher $D$), the distribution of $A_{ii}$ is indeed affected. The distributions become narrower with an increase in $D$. This shrinking of the distribution function explains the reduction in $\theta ^{\prime 2}$ evident in figure 9(a).

Figure 10. The PDF of $A_{ii}$. Influence of initial Damköhler number ($D$). Cases B2, C, D and G.

To further examine the influence of the Damköhler number on the velocity gradient field, we turn our attention to the strain-rate tensor. We specifically examine the locally normalized strain-rate tensor ($s_{ij}$), where

(6.2a,b)\begin{equation} s_{ij} = \frac{1}{2} (a_{ij} + a_{ji})\quad \text{and}\quad a_{ij} = \frac{A_{ij}}{\sqrt{A_{mn} A_{mn}}}. \end{equation}

The eigenvalues of the locally normalized strain-rate tensor may not give direct information about the absolute magnitudes. Still, they are an excellent indicator of the relative stretching/compression of a small fluid element in the principal coordinate system of $s_{ij}$. We express these self-normalized eigenvalues of the $\boldsymbol {s}$ tensor ((6.2a,b), Suman & Girimaji 2009)

(6.3a–c)\begin{equation} \alpha^* = \frac{\alpha}{\sqrt{\alpha^2 + \beta^2 + \gamma^2}}, \quad \beta^* = \frac{\beta}{\sqrt{\alpha^2 + \beta^2 + \gamma^2}}, \quad \gamma^* = \frac{\gamma}{\sqrt{\alpha^2 + \beta^2 + \gamma^2}}, \end{equation}

where the eigenvalues $(\alpha, \beta, \gamma )$ satisfy the condition

(6.4)\begin{equation} \alpha \geq \beta \geq \gamma. \end{equation}

In incompressible turbulence, it is a well-known behaviour that the ratio $\alpha ^*:\beta ^*:\gamma ^*$ attains a most-probable state of 3:1:–4 (Ashurst et al. 1987). This presents the physical picture of a small fluid element being severely stretched and severely compressed in two mutually perpendicular directions. In contrast, the third direction still has a moderate stretching to ensure that the volume of the fluid element is preserved ($\alpha ^*+\beta ^*+\gamma ^*=0$). In compressible turbulence, however, there would be an obvious departure from these proportions because the velocity field is no longer divergence free $(\alpha ^*+\beta ^*+\gamma ^*\neq 0)$. Our intent here is to identify whether the presence of vibrational non-equilibrium in a compressible flow field further influences the proportionality of the three eigenvalues or not. In line with this objective, in figure 11(a–c) we present the PDFs of $\alpha ^{*}$, $\beta ^{*}$ and $\gamma ^{*}$ extracted from simulations which differ in terms of the initial Damköhler number at a representative time $t/\tau _{LE}=3.5$. Figure 11(a) includes PDFs of $\alpha ^{*}$, figure 11(b) includes PDFs of $\beta ^{*}$ and figure 11(c) includes PDFs of $\gamma ^{*}$. All the simulations represented here have identical Reynolds number $(Re_{\lambda }=60)$ and Mach number ($M_t = 0.5$). Considering the no-vibration case (simulation case G) to be our datum for comparison at the specified Mach number and Reynolds number, we observe that the Damköhler number does have an influence on the distributions of normalized eigenvalues. The overall trend is that, at increased Damköhler number, the extreme eigenvalues ($\alpha$ and $\gamma$) tend to decrease in terms of their intensity, but they still retain predominantly positive and negative signs, respectively. On the other hand, the intermediate eigenvalue PDFs seem to be mostly unaffected as the Damköhler number increases.

Figure 11. The PDFs of normalized strain-rate eigenvalues ($\alpha ^* \geq \beta ^* \geq \gamma ^*$). Cases B2, C, D and G.

To further understand how a change in the Damköhler number results in the changes observed in figure 9, we examine the evolution equation of the instantaneous rotation-rate component $W_{ij}$ and the trace of the velocity gradient tensor $A_{ii}$

(6.5)$$\begin{gather} \frac{{\rm D} W_{ij}}{{\rm D} t} ={-}S_{ik}W_{kj} - W_{ik}S_{kj} - 0.5(P_{ij} - P_{ji})+ 0.5(\varGamma_{ij} - \varGamma_{ji}), \end{gather}$$

(6.6)$$\begin{gather}\frac{{\rm D} A_{ii}}{{\rm D} t}={-A_{ik}A_{ki}}-{P_{ii}}+{\varGamma_{ii}}, \end{gather}$$

where

(6.7a,b)\begin{equation} S_{ij} = \tfrac{1}{2} (A_{ij} + A_{ji})\quad \text{and} \quad W_{ij} = \tfrac{1}{2} (A_{ij} - A_{ji}). \end{equation}

In (6.5), $S_{ij}$ is the strain-rate tensor, $W_{ij}$ is the rotation-rate tensor, $P_{ij}$ is the pressure-Hessian tensor and $\varGamma _{ij}$ is the viscous transport term. While some aspects of the dynamics of $\theta ^{\prime 2}$ can be understood with reference to the evolution equation of $A_{ii}$, the evolution equation of the components of the rotation-rate tensor ($W_{ij}$) can prove to be helpful in understanding the behaviour of $\omega ^{\prime 2}$. Using (6.5) and (6.6), we write the evolution equation for the magnitude of the vorticity vector (Wilczek & Meneveau 2014) and that for $A_{ii}A_{mm}$ (neglecting the viscous processes on the right-hand sides)

(6.8)$$\begin{gather} \frac{{\rm D}}{{\rm D} t} (\omega_i \omega_i) = S_{ik} \omega_k \omega_i - 2 \varepsilon_{ijk} \tilde{P}_{jk} \omega_i, \end{gather}$$

(6.9)$$\begin{gather}\frac{{\rm D}}{{\rm D} t} (A_{ii}A_{mm})={-A_{ik}A_{kj}A_{mm}}-{P_{ii}A_{mm}}, \end{gather}$$

where $\tilde {P}_{ij} = \frac {1}{2} (P_{ij} - P_{ji}),\ \omega _i = -\varepsilon _{ijk} W_{jk}, \text {and}\ \varepsilon _{ijk}$ is the Levi-Civita tensor.

We have neglected the viscous processes in (6.8) and (6.9) as our immediate focus in this study to examine the roles of the inertial and pressure-related processes only. On the right-hand side of these equations, we identify the presence of two types of mechanisms: inertial (the respective first terms in the two equations) and pressure related (the respective second terms in the two equations). To better understand how the initial Damköhler number influences the evolution of $\theta ^{^\prime 2}$ in figure 9(a) but does not influence the evolution of $\omega ^{\prime 2}$ in figure 9(b), we examine the behaviour of these inertial and pressure mechanisms in the simulation cases B2, C, D and G.

In figure 12(a,b), we plot the temporal variations in the mean values of the inertial ($\langle |S_{ik}\omega _k\omega _i|\rangle$) and pressure mechanisms ($\langle |-2\varepsilon _{ijk} \tilde {P}_{jk} \omega _i|\rangle$) in (6.8). The symbol $|\text { }|$ denotes the absolute value of the scalar quantity in context. Every plot in these figures has been normalized with the initial value of $\omega ^{^\prime 3}$. In figure 12(a), the inertial mechanism of (6.8) is not much affected by the changed Damköhler number. On the other hand, in figure 12(b), even though the pressure mechanism of (6.8) seems to be affected by the initial value of the Damköhler number, its strength is evidently too small to have any consequence compared with the inertial mechanism of ((6.8), figure 12a). The strength of the pressure mechanism is small because it takes a contribution from the antisymmetric part of the pressure-Hessian tensor (6.5). Danish et al. (2014) demonstrate that, in compressible decaying turbulence, the generation of the antisymmetric portion of the pressure-Hessian tensor is directly influenced by the occurrence of shocklets. With the Mach number of our current simulations being small, the antisymmetric part of the pressure Hessian must be quite small, and thus (even though it might be influenced by the initial Damköhler number) is not large enough to have any consequence for the evolution of the strength of the vorticity vector (figure 9b).

Figure 12. Temporal variation in (a) the inertial mechanism and (b) the pressure mechanism in the evolution equation of (6.8). Simulation cases B2, C, D and G.

In figure 13(a,b) we plot the temporal variations in the mean values of the corresponding inertial $(\langle |-A_{ik}A_{kj}A_{mm}|\rangle )$ and pressure mechanisms $(\langle |-P_{ii}A_{mm}|\rangle )$ of (6.9). We observe that a change in the initial Damköhler number has negligible influence on the inertial mechanism of (6.9), as well (like figure 12a). However, the pressure mechanism does show substantial influence of the initial Damköhler number. Further, unlike the behaviour observed in figure 12(a,b), the order of magnitude of the pressure mechanism in figure 13(b) is indeed comparable to that of the inertial mechanism in figure 13(a). Based on these observations, we conclude that the perceivable influence of the initial Damköhler number on the magnitude of the dilatation fluctuations in figure 9(a) is attributable to the changes in the pressure-Hessian tensor induced by the initial Damköhler number.

Figure 13. Temporal variation in (a) the inertial mechanism and (b) the pressure mechanism in the evolution equation of (6.9). Simulation cases B2, C, D and G.

Suman & Girimaji (2011) have examined the statistics of an aptly normalized form of the algebraic sum of the inertial and pressure mechanisms in the evolution equation of $A_{ii}$ (6.6). This algebraic sum does represent the rate of change in the dilatation rate ($A_{ii}$) following a local small fluid element (while neglecting the role of viscous processes). The authors reported that, as the initial turbulent Mach number increases in their simulation (even with no vibrational physics included therein), the PDF of this quantity spreads, implying the presence of higher rates of changes in the dilatation rates of fluid elements. Following a similar line of investigation as followed by those authors, in this work, we define a quantity, $\delta$

(6.10)\begin{equation} \delta={-}A_{ik}A_{ki}-P_{ii}. \end{equation}

In figure 14(a), we present the temporal evolution of the mean of $\delta ^2$ from the simulation cases B2, C, D and G. All these curves are normalized by the initial value of $\omega ^{\prime 4}$. We observe that, in figure 14(a), as the initial Damköhler number increases, there is indeed a suppression in the magnitude of $\delta$. This suppression does provide an explanation for the suppression in $\theta ^{\prime 2}$ observed earlier in figure 9(a). In figure 14(b) we present the temporal evolution of the mean of $\delta ^2$ in simulation cases B2, E, F (which differ in terms of initial $\xi$ alone) and G. Evidently, the initial value of $\xi$ has no consequence for the evolution of the mean of $\delta ^2$. This observation is consistent with the observation made in figure 9(c) wherein no influence of $\xi$ has been observed on the evolution of $\theta ^{\prime 2}$.

Figure 14. Temporal evolution of $\langle \delta ^2\rangle$ in simulation cases (a) B2, C, D and G and (b) B2, E, F and G.

Suman & Girimaji (2011) have demonstrated that, in the context of the velocity gradient dynamics (2.15), normalizing the pressure-Hessian tensor ($P_{ij}$) with the square of the velocity gradient tensor's magnitude is appropriate. This normalized form delineates the relative importance of the pressure-Hessian process in comparison with the self-stretching terms in (2.15). This relative magnitude is presented as $\phi =\sqrt {P_{mn}P_{mn}}/(A_{ij}A_{ij})$, which is a dimensionless quantity. In figure 15(a), we present the PDF of $\ln \phi$ extracted from simulations B2, C, D and G at the time instant $t/\tau _{LE}=3.5$. All these PDFs are single-peaked distributions. However, we do observe some perceivable effect on these PDFs as the Damköhler number changes. With an increase in the value of the initial Damköhler number, the peak of the distribution shifts to the left, indicating a general decrease in the strength of the pressure-Hessian tensor compared with that of the velocity gradient tensor. In figure 15(b), we present the PDF of $\ln \phi$ extracted from simulations B2, E, F and G at the time instant $t/\tau _{LE}=3.5$. Evidently, the influence of $\xi$ on these PDFs is negligible.

Figure 15. The PDF of $\ln {{\sqrt {P_{mn}P_{mn}}}/{A_{ij}A_{ij}}}$: (a) simulation cases B2, C, D and G and (b) simulation cases B2, E, F and G.

Based on the results of this section, we conclude that, even within the moderate range of the Damköhler number considered in this work, there is a considerable influence of the vibrational non-equilibrium process (represented by varying $D$) on the pressure-Hessian tensor, which in turn, affects the velocity gradient statistics. Overall, an increased Damköhler number suppresses the dilatational fluctuations, whereas the initial value of $\xi$ does not seem to have any significant effect. To gain further insights, in the next section, we directly examine the temporal behaviour of the ‘new’ vibrational mechanisms instrumental in the evolution equation of the pressure-Hessian tensor (3.1).

7. Behaviour of vibrational mechanisms $VI^{'}$ and $VII^{'}$

In this section, we focus on the vibrational mechanisms that have emerged in the evolution equation of the pressure-Hessian tensor (3.1). Accordingly, we examine the behaviour of $f$ (3.2) and attempt to understand how the initial Damköhler number affects the vibrational mechanisms.

To understand the influence of the initial Damköhler number on $f$, we present the temporal evolution of $\langle |\,f|\rangle$ from simulation cases B2, C and D in figure 16(a). We observe that, in each case, $\langle |\,f| \rangle$ undergoes a two-phase evolution. In the first phase of evolution, $\langle |\,f| \rangle$ increases rapidly, reaches a peak value and then decays monotonically. The peak value of $\langle |\,f| \rangle$ is higher in a simulation with a higher initial Damköhler number. Indeed, over the entire duration of the simulations (up to 5 eddy-turnover times), there is an apparent incremental shift in $\langle |\,f| \rangle$ with an increase in the initial Damköhler number at all times. We observe the mean value of $f$ is predominantly negative in all simulations at all times. This implies that $VI < VII$ in (2.22) at most locations of the flow domain. This bias in the distribution is plausibly due to our choice of the initial state of the flow field wherein $T_v< T$ (or equivalently, $e_v< e$), which means $\xi < 1$. Figure 16(a) shows that such bias continues to exist in the flow field at later times too.

Figure 16. Temporal evolution of $|\,f|$. (a) Simulation cases B2, C and D and (b) simulation cases B2, E and F.

We have examined the influence of the initial ratio $\xi = T_v/T$, as well, on the vibrational mechanisms $VI'$ and $VII'$. The temporal evolution of $|\,f|$ from simulations B2, E and F are presented in figure 16(b). Evidently, the influence of initial $\xi$ on $|\,f|$ is significantly smaller than the influence of the initial Damköhler number (figure 16a). This contrasting behaviour does provide an explanation for the observations made earlier in the previous section wherein the statistics of the pressure-Hessian tensor and the velocity gradient tensor seem to be affected by the initial Damköhler number but not much by the initial value of $\xi$.

In incompressible fluid flows, the term $a_{ii}$ is consistently zero for all fluid elements. However, in compressible flows, $a_{ii}$ can assume positive, negative or zero values. A value of $a_{ii} > 0$ signifies an expanding fluid element, $a_{ii} < 0$ indicates contraction and $a_{ii} = 0$ implies an instantaneously volume-preserving fluid element.

Suman & Girimaji (2011) demonstrated that some of the flow statistics (related to velocity gradients) in a compressible flow field may not show any variation by changing a global compressibility parameter like $M_t$. However, when examining the statistics of the same quantities conditioned upon discrete levels of $a_{ii}$ or even upon the sign of $a_{ii}$, significant variations are observed. In our study, it is pertinent to examine whether $f$ exhibits any sensitivity to the local compressibility parameter $a_{ii}$. To examine this, we focus on the simulation case B2 ($M_t=0.5$). In figure 17, we present the PDF of the mean of $f$ conditioned upon the sign of instantaneous $a_{ii}$ at the time instant $t/\tau _{LE}=2.5$. We include two PDFs, separately conditioned upon (a) the contracting and (b) the expanding fluid elements. Clearly, these evolutions do not show much difference compared with the unconditional evolution (included in figure 17 for easy comparison). We have examined the conditioned PDFs of $a_{ii}$ at two more time instants ($t/\tau _{LE}=1.5$, and $t/\tau _{LE}=3.5$) (results not included to avoid repetition). However, the PDFs of $f$ remain identical for all time instants. Based on these results, we conclude that the importance of the vibrational mechanism in (3.1) seems to be similar for the contracting and the expanding fluid elements. It is plausible to attribute these observations to the fact that there is no explicit appearance of $A_{ii}$ and, for that matter, any component of $A_{ij}$, in vibrational mechanisms $VI$ and $VII$ in (2.22).

Figure 17. Conditional PDF of $f$ (3.2). Influence of normalized dilatation rate ($a_{ii}$) (case B2).

8. Conclusions

The overarching motivation behind this work is to enhance our understanding of the velocity gradient dynamics in compressible turbulent flows. Specifically, we are interested in flows wherein a non-equilibrium exists between the vibrational and the translational–rotational energy modes of the constituent molecules. Such flow fields typically exist past aerospace vehicles like missiles and spacecraft, which move through planetary atmospheres at hypersonic speeds. Specifically, our interest is to understand how the vibrational relaxation process influences the pressure-Hessian tensor and the velocity gradient field.

In this work, we first derive the exact evolution equation of the velocity gradient and the pressure-Hessian tensors in a vibrationally excited flow field. We identify ‘new’ unclosed vibrational mechanisms that emerge in the evolution equation of the pressure-Hessian tensor. To further understand the physics of these new mechanisms, we perform an extensive set of DNS of compressible decaying turbulence. Besides the continuity, momentum and energy equations, these simulations solve an evolution equation for the vibrational energy, as well. The Landau–Teller model is used to simulate the exchange of energy between the vibrational mode and the translational–rotational modes. These simulations are performed over a range of Damköhler number and the initial ratio of the mean vibrational temperature to the mean temperature. The simulations employ the hy2Foam solver, which is based on the OpenFOAM platform.

The major findings of this study are:

1. (i) Among the two different non-dimensional parameters considered in this study, the Damköhler number is found to be the major influencer on the ‘new’ vibrational mechanisms of the pressure-Hessian tensor identified in this work. In contrast, the other parameter, the initial ratio of the mean vibrational temperature to the mean temperature of the flow field, does not seem to have any significant impact.

2. (ii) Even though the mean vibrational energy tends to reach an equilibrium with the mean translational–rotational energy over a few (3–5) eddy-turnover times, the difference in the fluctuating vibrational energy and the fluctuating translational–rotational energy tends to persist in the flow field for a long time. This keeps the vibrational mechanism active over a long duration, even when the turbulent kinetic energy has diminished significantly.

3. (iii) The vibrational mechanisms are found to be unaffected by the sign of the dilatational state of the local fluid element (expanding/contracting fluid elements).

4. (iv) As the vibrational relaxation becomes more rapid (in other words, as the initial Damköhler number increases in a simulation), the vibrational mechanisms of the pressure-Hessian equation become more pronounced. However, this pronouncement has a mitigating effect on the strength of the pressure-Hessian tensor.

5. (v) An increase in the Damköhler number causes a significant reduction in the strength of the dilatational fluctuations of the flow field. However, the vortical fluctuations are not affected by a change in the initial Damköhler number.

   The findings of this study are expected to aid in the development of closure models for the pressure-Hessian tensor in compressible flows under the influence of vibrational non-equilibrium.