Velocity-vorticity correlation and turbulent diffusivity in buoyancy driven fluid dynamics

TL;DR

This paper explores how the velocity-vorticity correlation integral can be used to describe turbulent diffusivity in buoyancy-driven flows, supported by numerical simulations of various convection and mixing phenomena.

Contribution

It demonstrates that the velocity-vorticity correlation integral effectively characterizes turbulence in buoyancy-driven systems, linking invariant properties to turbulent diffusivity.

Findings

Distributed chaos governed by the invariant explains Rayleigh-Bénard convection.

The approach applies to stably stratified turbulence.

It models Rayleigh-Taylor mixing at Prandtl number around 1.

Abstract

The velocity-vorticity correlation integral (Chkhetiani invariant) is an invariant of a viscous Karman-Howarth equation. In the tree-dimensional space this invariant can naturally determine a turbulent diffusivity (viscosity). It is shown, using results of direct numerical simulations, that distributed chaos dominated by this integral can provide an adequate description of the turbulent Rayleigh-B\'{e}nard convection, stably stratified turbulence and Rayleigh-Taylor mixing at Prandtl number $Pr \sim 1$.

Equations

$$\mathcal{I}=\int ⟨\mathbf{u}(\mathbf{x},t)\cdot \mathbf{\omega }(\mathbf{x}+\mathbf{r},t)⟩d\mathbf{r}$$

$$⟨u_i(\mathbf{x})u_j(\mathbf{x}+\mathbf{r})⟩=A(r){\delta }_{ij}+B(r)r_i{r}_j+C(r){\epsilon }_{ijl}{r}_l$$

$$⟨u_i(\mathbf{x})u_k(\mathbf{x})u_l(\mathbf{x}+\mathbf{r})⟩=D(r){\delta }_{ik}{n}_l+E(r)({\delta }_{il}{n}_k+{\delta }_{kl}{n}_i)$$

$$+F(r)n_i{n}_k{n}_l+S(r)\left({\epsilon }_{ilj}{n}_j{n}_k+{\epsilon }_{klj}{n}_j{n}_i\right)$$

$$⟨\mathbf{u}(\mathbf{x},t)\cdot \mathbf{\omega }(\mathbf{x}+\mathbf{r},t)⟩=-6C(r)-2r\frac{\partial C(r)}{\partial r}$$

$$\frac{\partial }{\partial t}C(r)=\frac{2\nu }{r^4}\frac{\partial }{\partial r}\left(r^4\frac{\partial }{\partial r}C(r)\right)+\frac{2}{r^4}\frac{\partial }{\partial r}\left(r^{3}S(r)\right)$$

$$\frac{\partial\mathcal{I}(r)}{\partial t}=-2\left(6rS(r)+2r^{2}\frac{\partial S(r)}{\partial r}\right)\bigg{\rvert}_{0}^{\infty}$$

$$-2\nu\left(8r^{2}\frac{\partial C(r)}{\partial r}+2r^{3}\frac{\partial^{2}C(r)}{\partial r^{2}}\right)\bigg{\rvert}_{0}^{\infty}=0$$

$$\mathcal{D}=c|\mathcal{I}{|}^{1/2}$$

$$\frac{\partial ({\nabla }^2{u}_z)}{\partial t}={\nabla }^4{u}_z-{\mathbf{e}}_z\cdot \nabla \times \left[(\mathbf{\omega }\cdot \nabla )\mathbf{u}-(\mathbf{u}\cdot \nabla )\mathbf{\omega }\right]$$

$$+Ra\left(\frac{{\partial }^2\theta }{\partial x^2}+\frac{{\partial }^2\theta }{\partial y^2}\right),$$

$${\frac{\partial}{\partial t}}\omega_{z}=\nabla^{2}\omega_{z}+\left[(\mbox{\boldmath$\omega$}{\cdot}\nabla)v_{z}-({\bf u}{\cdot}\nabla)\omega_{z}\right],$$

$${\nabla }^2\theta =-u_z,\nabla \cdot \mathbf{u}=0,$$

$$E(k)=a\exp -(k/k_0)$$

$$E(k)=\int P(a,k_0)\exp -(k/k_0)dad{k}_0$$

$$E(k)\propto \int P(k_0)\exp -(k/k_0)d{k}_0$$

$$E(k)\propto \int P(k_0)\exp -(k/k_0)d{k}_0\propto \exp -(k/k_{\beta }{)}^{\beta }$$

$$u_0\propto k_0^{\alpha }$$

$$u_0\propto |\mathcal{I}{|}^{1/2}{k}_0\propto \mathcal{D}{k}_0$$

$$\mathcal{P}(u_0)d{u}_0\propto P(k_0)d{k}_0$$

$$P(k_0)\propto k_0^{\alpha -1}\mathcal{P}(u_0(k_0))$$

$$P(k_0)\propto k_0^{-1+\beta /[2(1-\beta )]}\exp (-b{k}_0^{\beta /(1-\beta )})$$

$$\beta =\frac{2\alpha }{1+2\alpha }$$

$$E(k)\propto \exp -(k/k_{\beta }{)}^{2/3}$$

$$\frac{\partial \mathbf{u}}{\partial t}+(\mathbf{u}\cdot \nabla )\mathbf{u}=-\nabla p+\theta {\mathbf{e}}_z+\sqrt{\frac{\mathrm{Pr}}{\mathrm{Ra}}}\nabla {}^2\mathbf{u}$$

$$\frac{\partial \theta }{\partial t}+(\mathbf{u}\cdot \nabla )\theta =\mathbf{u}\cdot {\mathbf{e}}_{\mathbf{z}}+\frac{1}{\sqrt{\mathrm{Ra}\mathrm{Pr}}}\nabla {}^2\theta$$

$$\nabla \cdot \mathbf{u}=0$$

$$\frac{\partial \theta }{\partial t}+(\mathbf{u}\cdot \nabla )\theta =-\mathbf{u}\cdot {\mathbf{e}}_{\mathbf{z}}+\frac{1}{\sqrt{\mathrm{Ra}\mathrm{Pr}}}\nabla {}^2\theta$$

$$\frac{\partial \mathbf{u}}{\partial t}+\mathbf{u}\cdot \nabla \mathbf{u}=-\frac{1}{{\rho }_0}\nabla p+N\theta {\mathbf{e}}_{\mathbf{g}}+\nu {\nabla }^2\mathbf{u}$$

$$\frac{\partial \theta }{\partial t}+\mathbf{u}\cdot \nabla \theta =N\mathbf{u}\cdot {\mathbf{e}}_{\mathbf{z}}+D{\nabla }^2\theta$$

Full text

I Viscous invariants and effective turbulent diffusivity

Let us, following the Ref. otto1 , consider invariance of the velocity-vorticity correlation integral:

[TABLE]

in freely decaying homogeneous incompressible fluid turbulence ($<...>$ denotes an ensemble average).

The asymmetric characteristics of the velocity correlation tensors: $C(r)$ and $S(r)$, are defined by the pair correlation tensor of velocities:

[TABLE]

and by the third-rank two-point correlation tensor of velocities

[TABLE]

[TABLE]

where ${\bf n}={\bf r}/|{\bf r}|$.

The velocity-vorticity correlation can be represented as

[TABLE]

Then it follows from an analogue of the Karman-Howarth equation otto2 :

[TABLE]

(the $\nu$ is viscosity) that

[TABLE]

[TABLE]

if $C(r)$ and $S(r)$ have the asymptotic behavior $O(r^{-2})$ at $r\rightarrow\infty$. This equation means dynamical invariance of the velocity-vorticity correlation integral $\mathcal{I}$ - Chkhetiani invariant.

One can see from the last term in the Eq. (6) that, unlike the energy and helicity, the invariance of the $\mathcal{I}$ takes place also for the molecular diffusive (viscous) dynamics. Moreover, in the three-dimensional space this invariant can determine an effective turbulent diffusivity (viscosity):

[TABLE]

where $c$ is a dimensionless constant.

These results can be generalized on the buoyancy driven turbulence (cf., for instance, Ref. davidson ). It should be noted that there is a whole family of the viscous invariants (the Loitsyanskii’s, Birkhoff-Saffman’s and many others, see for instance Refs. otto1 , davidson ,vas ,dok ). Domination of one of these integrals depends on fundamental symmetries of the flow (isotropy, homogeneity etc.) and initial conditions. For the buoyancy driven turbulence the type of turbulence behaviour is also strongly dependent on the boundary conditions. The physical relevance is crucial for application of these viscous invariants in each particular case.

II Distributed chaos in buoyancy driven turbulence

Chaotic or turbulent motion is characterized by a broadband spectrum. The particular case of the Rayleigh-Bénard (thermal) convection with the Prandtl number $\lim Pr\rightarrow 0$ exhibits such spectrum already at the onset paul . The equations describing this special case under Boussinesq approximation have the form spi ,t

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where ${\bf e}_{z}$ is the unit vector in vertical direction and $\theta$ is the temperature field deviation from the steady profile.

In the Ref. paul a direct numerical simulation (DNS) of these equations was performed with the perfectly conducting boundary conditions for the bottom (heated form below) and top horizontal plates. For the velocity field the freeslip boundary condition were used.

Figure 1 shows kinetic energy spectrum, computed at $Ra/Ra_{c}=1.1$ (where $Ra_{c}$ is a critical value of the Rayleigh number $Ra$ for the onset of the convection). The spectral data were taken from Fig. 1 of the Ref. paul and are shown in the log-log scales (here and in all other figures $\log k\equiv\log_{10}k$). The dashed curve indicates exponential spectrum

[TABLE]

with $k_{0}\simeq 0.96$, that indicates a large-scale nature of the exponential decay even for the smallest scales. The local interactions results in the scaling spectra my , whereas the non-local interactions cause the exponential spectral decay b2 .

If, in a more developed chaotic (turbulent) motion, there is an ensemble of the fields with the statistically varying parameters $a$ and $k_{0}$, then the ensemble averaged spectrum should be computed using equation

[TABLE]

where the $P(a,k_{0})$ is a joint probability distribution for the variables $k_{0}$ and $a$. If the parameters $a$ and $k_{0}$ are statistically independent, then

[TABLE]

where the $P(k_{0})$ is probability distribution of $k_{0}$.

A natural generalization of the exponential spectrum Eq. (12) is a stretched exponential spectrum

[TABLE]

Let us assume a scaling dependence of characteristic velocity $u_{0}$ on the scale $k_{0}$ at large values of the $k_{0}$

[TABLE]

Then one obtains from the dimensional considerations

[TABLE]

cf. Eq. (7).

Let us denote distribution of the characteristic velocity $u_{0}$ as $\mathcal{P}(u_{0})$. From the relationship

[TABLE]

and the Eq. (16) one obtains

[TABLE]

From the Eq. (15) one can estimate asymptote of $P(k_{0})$ at $k_{0}\rightarrow\infty$ jon

[TABLE]

where $b$ is a constant.

Then it follows from the Eqs. (16),(19) and (20) that for Gaussian probability distribution $\mathcal{P}(u_{0})$

[TABLE]

Substituting $\alpha=1$, from the Eq. (17), into Eq. (21) one obtains $\beta=2/3$, i.e.

[TABLE]

The spectrum Eq. (22) was suggested for the first time as a possible viscous (molecular diffusion) asymptote for isotropic homogeneous turbulence in Ref. b1 .

III Rayleigh-Bénard convection at small Prandtl numbers

Now we can consider turbululent Rayleigh-Bénard convection at non-zero, although small, Prandtl numbers. In a recent paper Ref. vvs results of a DNS for the Rayleigh-Bénard convection at $Pr=1$ and $Ra=10^{8}$ were reported. The dimensionless equations describing this case under Boussinesq approximation have the form

[TABLE]

[TABLE]

[TABLE]

where $p$ and $\theta$ are the pressure and temperature fluctuations from the conduction state correspondingly. The conducting and free-slip boundary conditions were used at the bottom and top horizontal plates, and periodic boundary condition were used at the side walls of the computational domain.

Figure 2 shows kinetic energy spectrum computed for this DNS at $Ra=10^{8}$ and $Pr=1$. The dashed curve indicates the stretched exponential spectrum Eq. (22) in the log-log scales. One can see that the spectrum Eq. (22) covers the spectral data both for the large and small scales. This phenomenon has two-fold origin. First is the chaotic nature of the spectrum (cf. Fig. 1 and position of the $k_{\beta}$), and second is the diffusive nature of the spectrum: cf. Eqs. (7) and (17).

IV STABLY STRATIFIED TURBULENCE

For stably stratified turbulence the temperature gradient is opposite to that of the Rayleigh-Bénard convection and the Eq. (24) should be replaced by equation

[TABLE]

The first term in the right hand side of this equation has its sign changed in comparison with the Eq. (24).

In paper Ref. kcv results of direct numerical simulations of the stably stratified turbulence with periodic boundary conditions in all directions (horizontal and vertical) were reported. In order to obtain a statistically steady state of the stably stratified turbulence a random forcing was applied in Eq. (23) at small wavenumbers according to the Ref. kh .

Figure 3 shows kinetic energy spectrum computed for this DNS at the $Pr=1$ and the Richardson number $Ri=0.01$. The Richardson number represents relation of the buoyancy to the nonlinear velocity term. The spectral data were taken from Fig. 2a of the Ref. kcv . The dashed curve in the figure indicates the stretched exponential spectrum Eq. (22) in the log-log scales.

V UNSTABLY STRATIFIED TURBULENCE and RAYLEIGH-TAYLOR MIXING ZONE

In paper Ref. bur the Rayleigh-Taylor mixing zone is considered as a unstably stratified turbulent flow and described by equations

[TABLE]

[TABLE]

[TABLE]

with periodic boundary conditions in all directions. In these equations the buoyancy field $\theta$ is rescaled as a velocity, $N=\sqrt{2gA(d{\bar{c}}/dz)}$, where $d{\bar{c}}/dz$ is the mean concentration z-gradient of the heavy fluid, $\rho_{0}=(\rho_{heavy}+\rho_{light})/2$ and $A=(\rho_{heavy}-\rho_{light})/\rho_{heavy}+\rho_{light})$ is the Atwood number representing the relative density contrast between the heavy and light mixed fluids.

Figure 4 shows kinetic energy spectrum computed in this DNS at N=4, Pr = 1 and the time of the mixing zone development $t=6$. The spectral data were taken from Fig. 6 of the Ref. bur . The dashed curve indicates the stretched exponential spectral law Eq. (22) and the dotted arrow is drawn in order to indicate a tuning of the distributed chaos at all scales to the large-scale coherent structures.

VI Acknowledgement

I thank M.K. Verma and R. Samuel, for sharing their data, and O.G. Chkhetiani for sending his paper and comments.