Understanding turbulence

The subject of turbulence modelling conjures up images of page long equations that can be intimidating to the average reader (at least it is to me). The importance of understanding turbulence modelling was brought up to me in a recent project where the boundary layer takes up a critical part of the flow. I decided to take this opportunity to try and gain an intuitive understanding of turbulence.

Reynolds stresses from a statistical point of view

Most undergraduate courses introduce turbulence from the Reynolds Averaged Navier-Stokes (RANS) equations and highlight the additional $\left<u_i'u_j'\right>$ terms brought about by the Reynolds decomposition process. They then go on to explain the difficulty in closing the equations and introduce the Boussinesq approximation as way of solving the problem in an engineering way. This is a nice way of explaining turbulence in a nutshell, but what is turbulence?

The Reynolds stress tensor can actually be thought of as a covariance matrix between velocities in the various components. Recall that covariance between two random quantities $x$ and $y$ is taken to be:

$\begin{align*} cov\left(x,y\right) = \frac{1}{N} \sum_n^N \left(x_n-\left<x\right>\right)\left(y_n-\left<y\right>\right) = \left<\left(x-\left<x\right>\right)\left(y-\left<y\right>\right)\right> \end{align*}$

Also recall the Reynolds decomposition:

$\begin{align*} u_i &= \left<u_i\right> + u_i' \\ \Rightarrow \left<u_i'u_j'\right> &= \left<\left(u_i-\left<u_i\right>\right)\left(u_j-\left<u_j\right>\right)\right> = cov(u_i,u_j) \end{align*}$

Hence we see how the Reynolds stress tensor can be thought of as a covariance matrix. So what does this seemingly statistical quantity have to do with stresses being applied on a fluid? First, we have to imagine a boundary on an infinitesimal cube of fluid.

Take the top surface for example. At an instance of time, the instantaneous flux of fluid crossing the boundary into the cube injects momentum in the x-direction ( $M_1$ ):

$\begin{align*} \Delta M_1 = \rho \left(-\vec{u}\cdot \hat n\right)u_1 \Delta A = -\left(\rho u_2\right)u_1 \Delta A \end{align*}$

Where the terms in parenthesis represent the mass flux going into the cube. Likewise, flux leaving the cube can take positive/negative momentum along with it, reducing/increasing the momentum left in the cube. In turbulent flows, these fluxes fluctuate. Averaged over time:

$\begin{align*} \left<\Delta M_1\right> = -\rho \left<\left(\left<u_1\right>+u_1'\right)\left(\left<u_2\right>+u_2'\right)\right> \Delta A = -\rho \left( \left<u_1\right>\left<u_2\right> + \left<u_1'u_2'\right> \right) \Delta A \end{align*}$

Where the $\left<u_i\right>u_j'$ terms cancel in the averaging process. We see that the covariance terms add an additional momentum change on average. In more concrete terms, imagine we have a scatter of velocities on the top surface of the cube as shown on the left plot:

The $u_1$ and $u_2$ velocities have mean values of 0 and have the following covariance matrix:

$\begin{align*} \begin{bmatrix} 1 & 0.75 \\ 0.75 & 1 \end{bmatrix} \end{align*}$

The positive correlation between $u_1$ and $u_2$ tells us that there is a tendency for the fluid to have a negative x-velocity when the fluid enters the cube (i.e. $u_1<0$ when $u_2<0$ ), and a positive x-velocity when fluid is leaving the cube (i.e. $u_1>0$ when $u_2>0$ ). On average, fluid crossing the top boundary brings low momentum in and takes high momentum out, reducing the momentum in the cube.

The covariance measures how strong that tendency is and hence how strongly the momentum is affected. This is exemplified in the plot on the right which shows the $\Delta M_1$ term for each point on the scatter plot. The red line shows that the mean values of $-u_1u_2$ is $-0.74$ in this case. Theoretically this value should take on the off diagonal values in the covariance matrix, which is $-0.75$ .

Hence we can think of turbulence as the preference of the fluid to move in one direction in tandem with another on average. Because of the fluid’s preference to fluctuate in a certain way, it causes an average momentum flux which is non zero in general. The covariance of the fluid varies throughout the flow domain and between different flows, and is a difficult quantity to predict (if at all).

Decomposition and rotation of the Reynolds stress tensor

As one rotates the axes (i.e. rotating the cube of fluid), the Reynolds stress tensor will be transformed into the new axes. However, there are certain properties of the tensor which will remain invariant. One such property is the trace of a tensor, which works out to be twice the turbulent kinetic energy. Hence one can subtract $\frac{2}{3}k\delta_{ij}$ from the Reynolds stress tensor to extract the Reynolds stress anisotropy $a_{ij}$ :

$\begin{align*} a_{ij} = \left<u_i'u_j'\right> - \frac{2}{3}k\delta_{ij} \end{align*}$

One can think of $\frac{2}{3}k\delta_{ij}$ as an analog to the hydrostatic pressure faced by an infinitesimal cube of fluid. This hydrostatic pressure always acts in the normal direction with the same magnitude regardless of how one rotates the cube.

At certain orientations, one can obtain a Reynolds stress tensor which is diagonal. Linking this back to our interpretation of the tensor as a covariance matrix, this is equivalent to performing a principal component decomposition! The eigenvalues then tell us how the velocities are scattered.

For example, the velocities may vary greatly in one principal component but not much in the other 2, resulting in 1D turbulence. The velocity may also be scattered in an ellipse on a 2D plane with little variation in the third component (2D turbulence). Finally, the velocities may be totally uncorrelated, resulting in a scatter within a sphere (3D isotropic turbulence).

Lumley triangle

A more rigorous method of describing turbulence like we did above is through the turbulence triangle. First, we nondimensionalize the anisotropy tensor:

$\begin{align*} b_{ij} = \frac{a_{ij}}{2k} = \frac{\left<u_i'u_j'\right>}{2k} - \frac{1}{3}\delta_{ij} \end{align*}$

For a general tensor, there are 3 quantities which are invariant with rotation. These quantities are simply referred to as the invariants and they are some times labelled as $I_1$ , $I_2$ and $I_3$ . The trace and determinant makes up 2 of the invariants. All invariants can be calculated from the eigenvalues. More information can be found on Wikipedia.

Because the trace of $b_{ij}$ is zero, we really only have 2 independent invariants. This means that we are able to describe the anisotropy on a 2D plot independent of the tensor coordinate axes. Pulling a plot from a very useful paper where the shape of the velocity scatter is shown as a surface:

This tells us the type of turbulence that is not based on any particular frame of coordinate system. Sometimes instead of $I_2$ and $I_3$ , we use $\xi$ and $\eta$ , which is calculated slightly differently from the eigenvalues to stretch out the triangle for easier visualisation. As an example, let’s look at the variation in turbulence throughout a turbulent boundary layer. The data was taken here from DNS simulations.

On the left plot above, we see the variation of the various terms of the Reynolds stress tensor in the boundary layer. The $\left<u'u'\right>$ component rises quickly and reaches a peak around $y+ \approx 7$ above all other components. Far from the wall, all components on the diagonal of the Reynolds stress tensor converge.

The variation in turbulence is plotted on the Lumley triangle on the right in red. Close to the wall, the turbulence is highly 2D (the fluctuations are dominated by $\left<u'u'\right>$ and $\left<w'w'\right>$ ), it moves towards the 1D corner as $\left<u'u'\right>$ reaches its peak at $y+ \approx 7$ . As we move further from the wall, the red dots move towards the origin, indicating that turbulence approaches isotropic as we move far away enough from the wall. There are many other ways to visualize the type of turbulence on a flow domain, one of which is quite cool.

Limits of the Boussinesq approximation

According to the Boussinesq hypothesis, the Reynolds stress anisotropy can be closed with:

$\begin{align*} a_{ij} = -\nu_T\left(\frac{\partial \left<u_i\right>}{\partial x_j}+\frac{\partial \left<u_j\right>}{\partial x_i}\right) = -\nu_T s_{ij} \end{align*}$

There are 2 important assumptions made: (1) the turbulence is described by the local mean flow gradients and (2) the anisotropy is aligned with the mean rate-of-strain tensor. Both these assumptions can be shown to be plain wrong in certain cases. Take the flow in such a channel for example.

In this channel, the flow gradients are zero in the parallel flow sections. In the varying section:

$\begin{align*} \left<u\right> = U_0 \quad \left<v\right> = ay \quad \left<w\right> = -az \quad \end{align*}$

Where $U_0$ and $a$ are constants. This flow was investigated experimentally in this paper and the results summarized in Wilcox’s book which is a good resource (along with Pope’s book).

The experimental and computational results ( $k-\omega^2$ model) are shown in markers and lines respectively above. With the Boussinesq approximation, the model was able to obtain correct trends in anisotropy in the varying section and immediately returns to isotropy in the parallel section afterwards ( $x>2.3$ ). This mirrors the mean rate-of-strain tensor. However in experiments, the flow takes time to return to isotropy. Hence the anisotropy of turbulence is affected by flow history, which is not modeled with the Boussinesq approximation. This is an example of when assumption (1) is not correct.

A separate shear flow example taken from Pope’s book shows when assumption (2) is wrong. Consider a flow with the following profile below.

In this example, the only non-zero gradient is $\partial \left<u\right>/\partial y$ so:

$\begin{align*} s_{ij} = \begin{bmatrix} 0 & \frac{\partial\left<u\right>}{\partial y} & 0 \\ \frac{\partial\left<u\right>}{\partial y} & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} \end{align*}$

Since $\left<u_i'u_j'\right> = -\nu_T s_{ij} + \frac{2}{3}k\delta_{ij}$ , we expect all the diagonal terms of the Reynolds stress tensor to be equal. However, both experimental and DNS results show that the Reynolds stress tensor have non equal diagonal terms. For example, from the DNS results:

$\begin{align*} \frac{\left<u_i'u_j'\right>}{k} = \begin{bmatrix} 1.06 & -0.33 & 0 \\ -0.33 & 0.32 & 0 \\ 0 & 0 & 0.62 \end{bmatrix} \end{align*}$

The diagonal terms would have been calculated to be $2/3$ using the turbulence viscosity model. Hence we see where turbulence models based on the Boussinesq approximation is not expected to perform well.

Understanding turbulence

Leave a reply Cancel reply