LES ​

Of course we can resolve eddies bigger than a mesh using velocity vectors.
But how about eddies smaller than a mesh?
The below picture shows the concept of the problem. ^[1]
We will going to use sub-grid model that represent eddies smaller than a single mesh size

Wave number indicates the size of the eddy. Actually, it is the reciprocal of the size of the eddy.

Wavenumber

IDK why the wavenumber is defined like this, but it is how it is.

The below diagram shows the relation between the size of the eddy(wavenumber) and the kinetic energy density.

The area under the curve is the turbulent kinetic energy, TKE.

As we can't shrink the size of the mash to infinitesimal size, there's a certain point that we cannot resolve eddies using cells.
If the blue area is smaller than 20%, we call it Good LES.
IDK know why, but this is kind of convention.

So we will have something called 'integral length scale which is representative of all the eddies at a location.
Because it is easier to look at a single value than the spectrum.

integral length scale is simply size of the averaged energy density eddy.
The area devided by integral length scale is the same.
We can calculate it by mathematical expression.

But the above mathematical expression is a bit absurd,as we don't know the function of energy density.
So how can we calculate it?
We can calculate from a precursor RANS calculation using either or model.

We said the good LES should resolve more than 80% of energy density.
In order to to this, the size of a cell should be smaller than one fifth of the integral length scale.
So that we can resolve more than 80% of the turbulent kinetic energy.

So if we want to evaluate if the mesh is good or not, we can define a new function f with variable integral length scale and cell volume

and is from the relations above.

As we discussed above, good LES should resolve more than 80% of the TKE.
So our goal is to resolve more than 80% of the TKE.

First, we will going to calculate the mean velocity of the flow.
As the CFD Code computes the instantaneous velocity , we will going to time average and get mean velocity

We will average the velocity after the trasient phase.

We will do almost the same process as RANS

As we all know, kinetic energy per mass is ,
So we can multiply fluctuating velocity components together, and it will lead us to some kinetic energy term.

We have three veolcity components, in each directions, so we will have 9 possible combinations,
, , , , , , , ,
There are instaneous Reynolds-Stresses that are resolved by the mesh
But, only normal components are used to calcaulte the resolved TKE.(, , )

And then, we time-average all those Reynolds-Stresses, and get Reynolds Stress tensor per unit density.

Anyways, we get Reynolds STress tensor per unit density.

From the diagonal components, we get resolved TKE,

just using normal components.

The reason why it is called resolved TKE is, we resolved turbulent kinetic energy only by using out mesh.
This is the best we can do here.
We can't resolve TKE smaller then a cell right now, but we will do this later.

The image shows the rewolved TKE

So our next goal is to get sub-grid scale TKE.

is the TKE of the eddies smaller than the mesh size.
But how do we do this?

depends on the method you choose to determine it.
There are several methods like smagorinky, WALE, etc.
We will cover this later on.

But the easiest one is solving kinetic energy transport model.
What about other methods?
One way is to calculate from sub-grid legnth scale
has the same concept as integral length scale
$l_{sgs} represents the size of eddies within the cell.

is defined as

If this is Smagornisky coefficient.

But, if we look at some of CFD mannuals, we can see that formula is slightly modified.
This is because near the walls, eddies are damped, so we modify the function.
If we are so close to the wall, we will have smaller eddies than we calculated from above.
So, modified formula is

This also represents the effect of high aspect ratio.
This will have better approximation.

Now, we have calculated and
Below two picture shows and

If we want to check the quality of our CFD, we can calculate the ratio and we want >0.8 in the entire domain.

If we don't have good quality, we can refine the mesh and increase so that we can resolve tubulent kinetic energy.

What can we do to remove eddies from the mesh?

We can do this by increasing the turbulence dissipation rate.

Turbulence Dissipation rate, is the rate that turbulence is converted to thermal energy.
It has units of turbulent kinetic energy per unit time.
We have large , eddies will quickly be dissipated into heat.
It's quite straightforward.

Then what's the mathmatical expression?
In real turbulence,

But this is just mathamteical deifintion.
How do we calculate it?
In RANS, we solve transport equation for , and this was shown in model
But this is different in LES.

If we look at below picture, we can see that velocity gradient gets larger as eddies get smaller.
This means, at first, as velocity grdient is not large, turbulenece dissipation rate is also not large and it takes a while for the eddies to be dissipated.
However, as time goes by, eddies get smaller and velocity gradient also gets smaller and turbulence dissipation rate will accelerate.

We know that eddies just bigger than the cell size are not being dissipated, and it means velocity gradient is not decreasing anymore, leading to constant turbulence dissipation rate.
But our goal is to dissipate eddies.
So, we can increase molecular viscosity , and it will lead to larger turbulence dissipation rate.
If the turbulence dissipate rate is high enough, eddies will be dissipated.
So the mathmatical expression of turbulence dissipation rate will be

If we select adequte to dissipate eddies just larger than the cell size, our goal is achieved.
But we should remind that we're not solving turbulence dissipation rate in LES.
How can we reflect this concept to the real CFD?

We looked at the concept of how to dissipate eddies in LES.
But it was only a concept.
How do we increase molecular viscoity in LES simulation?
We do this by modyfing N-S equation

We can see that extra stress term is added to the N-S equation.
viscosity stress term does the same role as molecular viscosity.
As we increased to , turbulence dissipation will also be increased, and it will dissipate eddies.
This can be derived by filtering the N-S equation, but we'll going to skip this.

It's important to know the meaning of sub-grid stress, which is here.

Molecular viscosity can be represented as Sub-Grid Stress in N-S equation.
We're going to look at the meaning of Sub-Grid Stress here.
If we think of eddies smaller than the size of cells, they can be considered as stresses.
Because eddies exert stresses to fluid particles and wall.
So we're substituing the concept of eddies into stresses instead of considering the real velocity of the eddies.

Sub-Grid Stress can be modeled as

And this is from eddy viscosity model.
You can find this on RANS page
If we see the equation, as the eddies get smaller, velocity gradient increases and it makes Stress tensor bigger.
But term depends on the size of the eddy, which means we don't have control for it.
Instead, we can control .
By controlling we can dissipate the eddy just bigger than the cell size.
This also means if we have larger cell size, also have to be bigger.

Because of the cell size - dependence, meshing and LES is incorporated.
In RANS, we can change the mesh and use the same turbulence model.
But in LES, as we change the mesh, molecular viscosity also changes.
This is the characteristic of the LES.

We can express sub-gird kinematic viscosity model as velocity multiplied by length according to the units of kinematic viscosity.

This simple approach was proposed by Smagorinky. So,

As eddies are isotropic, we only need a length to categorise their shape.

So our goal is to determine and to get

We are going to take velocity difference as a velocity scale here.
Velocity difference will be,

But what if we draw a horizontal line across the eddy?

Velocity difference will be

We can express this term in strain rate tensor

But, the problem is, we need a scalar value for the velocity scale.
What if we just take magnitude of ?
If we look at the simple case(Couette Flow), we can know that

So we have to multiply to the velocity gradient.

Smagorinsky Model doesn't take wall effect into consideration.
Wall makes eddy small

There are several ways.
But we're only going to cover modifying

we can change , , or
Simplest way is to modify .
We can do this by using RANS equation - eddy viscosity model

from and relation, and especially in logarithmic region, we can get velocity gradient.
And in log region, reynolds stress is constant which is wall shear stress.
From these relations, we can quantify the size of the eddy with the varient
The result will be shown in next section.

In logarithmic region, the mean velocity profile is modelled by

The dots are given from DNS simulation.
Also, in logarithmic region, Reynolds shear stress is relatively constant, hence,

So, if we plug those equations to Eddy Viscosity Equation , we get

The Turbulent Viscosity can be written as the product of a length scale and a velocity scale, so,

Turbulent Viscosity has the dimension of product of length and velocity.
So what will be the adequate choice for and to represent eddies?
will be the mixing length and will be the reynolds stress

Anyways, if we massage the equation, we get

So what does it mean?
As we get close to the wall, will decrease, and it will make smaller.
So this is the behavior that we wanted.
Taking wall effect into account.
Now we will return to where we start from, modeling sub-grid scale molecular viscosity which is the control of the dssipation of the eddy.

Now we have more accurate model for the sub-grid scale eddies.
So when modeling sub-grid length scale , instead of just using Smagorinsky Coefficient , we can use function. So,

But our was just for logarithmic region.
How can we consider other region?
We can do this by using Van Driest Solution.
We can get general equation that rerpresent relation between and

So the final Sub-Grid Length Scale will be