CFD Online Discussion Forums - Smagorinsky model details

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)

- - Smagorinsky model details (https://www.cfd-online.com/Forums/openfoam-solving/65838-smagorinsky-model-details.html)

lakeat

June 28, 2009 08:55

Smagorinsky model details

Sorry, I do not understand, I saw in "Smagorinsky.H",

Code:

tmp<volScalarField> k(const tmp<volTensorField>& gradU) const

{

    return (2.0*ck_/ce_)*sqr(delta())*magSqr(dev(symm(gradU)));

}

As I remember:

$\begin{array}{l} {\nu _{SGS}} = {\left( {{C_S}\Delta } \right)^2}\left| {{\bf{\bar S}}} \right| \\ K = {\left( {{C_I}\Delta } \right)^2}{\left| {{\bf{\bar S}}} \right|^2} \\ \left| {{\bf{\bar S}}} \right| = {\left( {{\bf{\bar S}}{\rm{:}}{\bf{\bar S}}} \right)^{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-\nulldelimiterspace} 2}}} \\ \end{array}$

Question 1: Why using magSqr(dev(symm(gradU))) instead of symm(gradU) && symm(gradU) to get ${{\bf{\bar S}}{\rm{:}}{\bf{\bar S}}}$ ????

Question 2: If magSqr(dev(symm(gradU))) = symm(gradU) && symm(gradU) = ${{\bf{\bar S}}{\rm{:}}{\bf{\bar S}}}$ , then

$K = \frac{{2{C_K}}}{{{C_\varepsilon }}}{\Delta ^2}{\left| {{\bf{\bar S}}} \right|^2}$

But I saw in "Smagorinsky.C"

Code:

nuSgs_ = ck_*delta()*sqrt(k(gradU));

Which means

${\nu _{SGS}} = {C_K}\Delta \sqrt K$

Then, replace K with $K = \frac{{2{C_K}}}{{{C_\varepsilon }}}{\Delta ^2}{\left| {{\bf{\bar S}}} \right|^2}$

${\nu _{SGS}} = {C_K}\Delta \sqrt K = {C_K}\Delta \sqrt {\frac{{2{C_K}}}{{{C_\varepsilon }}}{\Delta ^2}{{\left| {{\bf{\bar S}}} \right|}^2}} = {C_K}\sqrt {\frac{{2{C_K}}}{{{C_\varepsilon }}}} {\Delta ^2}\left| {{\bf{\bar S}}} \right|$

Compare with ${\nu _{SGS}} = {\left( {{C_S}\Delta } \right)^2}\left| {{\bf{\bar S}}} \right|$

We'll get

${\left( {{C_S}} \right)^2} = {C_K}\sqrt {\frac{{2{C_K}}}{{{C_\varepsilon }}}}$

But I heard somone said ${\left( {{C_S}} \right)^2} = {C_K}\sqrt {\frac{{{C_K}}}{{{C_\varepsilon }}}}$

So, I'm puzzled, I wonder if it was a mistake, that k should be written as

Code:

tmp<volScalarField> k(const tmp<volTensorField>& gradU) const

{

    return (ck_/ce_)*sqr(delta())*magSqr(dev(symm(gradU)));

}

Thank you

gaby	August 19, 2009 12:42

Hi Daniel,
Did you find any answer to this question?. in the paper:
"A SUBGRID-SCALE MODEL FOR LARGE-EDDY SIMULATION OF
PLANETARY BOUNDARY-LAYER FLOWS
PETER E SULLIVAN, JAMES C. McWILLIAMS, and CHIN-HOH MOENG" 1994, they defined Cs as:

Cs=(Ck*(Ck/Ce)^0.5)^0.5

So, I think you are right...

Now, I'm confused, why it is defined in Smagorinsky.H like 2*Ck/Ce ??

Is it because of the symm(grad(U)) definition???

It would be great if you could share your opinion...

Gaby

lakeat

August 20, 2009 08:35

Sorry, you see, no one come and help.

And you have noticed, I have done a detailed deduction in my top post, I still don't know why they use

2.0*ck_/ce_
magSqr(dev(symm(gradU))) instead of symm(gradU) && symm(gradU)

sandy

August 20, 2009 19:59

Quote:

Originally Posted by gaby (Post 226862)

Hi Gaby, can you send a copy of the paper to me? I am finding the expression to estimate the k inlet boundary condition in oneEqEddy model. Thanks.

Sandy
sandy.lee37@gmail.com

andrea

October 6, 2010 05:00

tensor norm definition

Hallo,
I was also trying to understand the implementation of Smagorinsky model.
As it is said here
http://www.cfd-online.com/Forums/ope...rain-rate.html
the definition of the norm of a tensor differs from what is computed in OF, so in my opinion the 2 before ck_/ce_ is exaclty the missing sqr(2) in the definition.
so, if
|S| = sqrt(2 S:S)
then
2 magSqr(S) = 2 sqrt(S:S)^2 = sqrt(2*S:S)^2 = |S|^2

Am I right?

PS: default Smagorinsky constant should be Cs =0.1677 (in Pope's book is said to be around 0.17) given ck=0.094, ce= 1.048

MaximeIST

March 23, 2011 15:31

Hello

I am also trying to understand how the Smagorinsky model is coded, for the incompressible version and also for the compressible version.

And It seems that for the incompressible Smagorinsky model, the default constant Cs is equal to
Cs=sqrt(ck*sqrt(2*ck/ce))

If I define Cs such that the eddy viscosity is equal to
nuSgs=( Cs *delta)^2 * ||D||

if ck=0.094 and ce=1.048 then Cs=0.1995.. ~ 0.2

Same question as Andrea : am I right?

maysmech

April 22, 2011 08:13

Hi,

I have same questions too.

alberto

April 22, 2011 23:22

Reference: http://pof.aip.org/resource/1/phfle6/v9/i5/p1416_s1

maysmech

May 1, 2011 11:53

Thanks Alberto for the useful reference.
As MaximeIST and lakeat told if we compare what is stated in your reference and openFoam smagorinsky with Pope's book we reach to Cs=0.2.
As stated in Pope's book this constant can be 0.1~0.2 and using 0.2 can be cause of high diffusivity. Please tell me if this is true.
My cyclone simulation with Smagorinsky has high diffusivity and i want use 0.1 for Cs, How can i do that? Which one of Ck or Ce should be changed? I mean are they be used elsewhere or not?

Regards,

alberto

May 1, 2011 18:27

I'd just use the dynamic Smagorinsky model, so that you do not have to fiddle with the coefficient, and you do not need to play with dumping functions.

yingkun

May 2, 2011 21:59

Hi,Andrea:
I think what you have said is right,the reason is just the different expression of Vsgs

maysmech

May 12, 2011 05:48

So, what is Cs value in OpenFOAM's Smagorinsky? 0.2 or 0.167?

yingkun

May 12, 2011 06:50

hi,Maysam:
since Cs=sqrt(ck*sqrt(*ck/ce), in smagorinsky.c Ck=0.094 , in GenEddyVisc.c Ce=1.048, then Cs=0.167 as we all know ,am I right?

maysmech

May 12, 2011 08:30

I am in doubt because of first post of this thread that it becomes 0.2.
I don't know. Any idea will be appreciated.

andrea

May 12, 2011 08:33

It's Cs \approx 0.17. You need just to follow the way it's computed. The question is why the simple Smagorinsky model is implemented in a such confusing way? We re-wrote our Smag. with just Cs required.

MaximeIST

May 12, 2011 09:12

Quote:

Originally Posted by yingkun (Post 307334)

hi,Maysam:
since Cs=sqrt(ck*sqrt(*ck/ce), in smagorinsky.c Ck=0.094 , in GenEddyVisc.c Ce=1.048, then Cs=0.167 as we all know ,am I right?

Hello

I may keep on confusing people, but the way it is coded, if I am not doing mistake is
Cs=sqrt(ck*sqrt(2*ck/ce))
in the incompressible Smagorinsky.H line 114.
There is a factor 2 added in the root-mean squared.

And in the case where Ce=1.048 and Ck=0.094, and with this factor 2, we obtain Cs=0.1995.

May be I miss something?

Maxime

yingkun

May 13, 2011 08:23

Quote:

Originally Posted by MaximeIST (Post 307357)

Hi,Maxime,
In the incompressible Smagorinsky.H line 114,it's the expression of K,not Ck,it means k=2Ck/Ce*(delta^2)*(Sof^2),may be my blog can help you: http://blog.sina.com.cn/s/blog_6d9c27ab0100u9ez.html
it's just what I think, I'm not sure.

MaximeIST

May 13, 2011 11:10

Hello Yingkun

I agree with you that it is the expression of k which is in line 114 of the Smagorinsky.H.

But the factor 2 is still there and as lakeat shown in the first message of this thread, we still get
Cs=sqrt(ck*sqrt(2*ck/ce))

I have been in your blog, and I think you lost this factor 2 in the passage of equation (7) to the last equation of the page.

Maxime

yingkun

May 14, 2011 01:01

Quote:

Originally Posted by MaximeIST (Post 307553)

hi,
I think you don't understand exactly what I mean,the difference is just caused by the different expressions of S between turbulent therory and Openfoam,there is a factor of sqrt(2) difference,lakeat just regards them the same

MaximeIST

May 16, 2011 06:13

Quote:

Originally Posted by yingkun (Post 307591)

OK! that's what I missed!
It's true I don't understand what is written in your blog, chinese is not easy for me :D!

So you said that the sqrt(2) is a consequence of the definition of S in OpenFoam.
In this case, the default constant in the incompressible Smagorinsky is 0.167.
good to know!

Bernhard

July 12, 2011 11:31

How about dynamic Smagorinsky?

Hi all,

I have some questions about this same issue with respect to the dynamic Smagorinsky model. I think I am making a mistake somewhere, but please correct me if I am wrong. I've tried to look at other locations for the same issue, but couldn't find it.

I compare the dynamic model in OpenFOAM (not caring about the domain averages coefficient) with the dynamic model (described in Lilly and Pope).

I think we agree on the following (D defined in OpenFOAM, S common definition)
$\begin{array}{rcl} \overline{S}_{ij} &=& \frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right) \\ \overline{D}_{ij} &=& \overline{S}_{ij} - \frac{1}{3} \delta_{ij}\overline{S}_{kk} = \overline{S}_{ij} (incompressible) \\ |\overline{D}| &=& |\overline{S}|/\sqrt{2} \end{array}$

Now, looking at the code (dynSmagorinsky 1.7, homogeneousDynSmagorinsky 2.0), in the .C file at line 57, 62 respectively:
$M_{ij} = \Delta^2 \left[ \widehat{|\overline{D}|\overline{D}_{ij}} - 4 |\widehat{\overline{D}}_{ij}\widehat{\overline{D}}_{ij}|\right]$
as compared to the original of Lilly (in Pope this expression
$M_{ij} = \widehat{\Delta}^2\widehat{|\overline{S}|}\widehat{\overline{S}}_{ij}-\Delta^2\widehat{|\overline{S}|\overline{S}_{ij}}$

If D_ij=S_ij, but |S|=sqrt(2) |D|, then the second term in the OF implementation is off by a factor sqrt(2). (Where the factor 4 is from the double width of the test filter, and the minus sine is also present in L)

Somewhere I must have skipped a step, I hope one of you can point it out for me. Thanks in advance!

Bernhard

July 15, 2011 02:35

Anybody that can give a comment on this?

Bernhard

July 18, 2011 03:47

Mistake is this:

Quote:

Originally Posted by Bernhard (Post 315854)

$M_{ij} = \Delta^2 \left[ \widehat{|\overline{D}|\overline{D}_{ij}} - 4 |\widehat{\overline{D}}_{ij}\widehat{\overline{D}}_{ij}|\right]$

$M_{ij} = \Delta^2 \left[ \widehat{|\overline{D}|\overline{D}_{ij}} - 4 |\widehat{\overline{D}}|\widehat{\overline{D}}_{ij}\right]$

gregor

August 26, 2011 06:39

Quote:

Originally Posted by lakeat (Post 220718)

Question 1: Why using magSqr(dev(symm(gradU))) instead of symm(gradU) && symm(gradU) to get ${{\bf{\bar S}}{\rm{:}}{\bf{\bar S}}}$ ????
Thank you

I think the reason is to have a consistent formulation for incompressible an compressible LES models :

let symm(gradU) be S, then dev(S) = S - 1/3 trace(S)I

however in a incompressible case 1/3 trace(S)I = 0, since trace(S) is the continuity eq. . Therefore in an incrompessible case it doesnt matter whether you take dev(S) or not , but consider a compressible case then 1/3 trace(S) doesn't vanish.

In a compressible or variable density case the solver calls divDevRhoBeff to compute the source term due to SGS stress B = 2/3k I - 2 nu_t S_D. (See Fureby http://pof.aip.org/resource/1/phfle6/v9/i5/p1416_s1 Eq. 3) There you have the deviatoric part of D. But i guess openFoam uses nu_t = c Delta^2 ||S_D|| and B = 2/3k I - 2 nu_t S. So it takes S_D for the turbulent viscosity and S for the SGS stress tensor B.

hz283

January 24, 2013 10:40

Quote:

Originally Posted by MaximeIST (Post 307357)

Hi All,

For the compressible Smagorinsky model, the parameters for ck=0.02, ce=1.048

Following the following line:
Cs=sqrt(ck*sqrt(ck/ce))
=> Cs=0.0525...

Does anybody know the references for these specification of the ck and ce for compressible Smagorinsky model?

Sunxing

March 11, 2013 21:25

Hi Yingkun,

As you mentioned, in the incompressible solvers Cs=sqrt(Ck*sqrt(Ck/Ce)). So if I want to set Cs=1, Do I just need to modify Ck and Ce in the LESProperties？Or are there other rules I must obey?

Code:

SmagorinskyCoeffs

{

ce 1.05;

ck 0.0472;

}

regards

itchy

February 28, 2014 03:03

Hi Bernhard and everybody,

I compare the dynamic model in OpenFOAM (homoDynSmag OF 2.2.2) with the dynamic model (described in Lilly and Pope) and I have the same Problem. To reanimate the discussion:

In Pope (F(.) means filtered):
[1] nu_SGS = cS * delta^2 * sqrt(2 * S_ij S_ij)
[2] cs = (M_ij L_ij)/(M_kl M_kl)
where
[3] S_ij = 0.5 (ui,j + uj,i)
[4] M_ij = 2 * delta^2 * (F(sqrt(2 * S_ij S_ij) S_ij) - F(sqrt(2 * S_ij S_ij)) F(S_ij))
[5] L_ij = F(ui uj) - F(ui) F(uj)

In OF 2.2.2 (homogeneousDySmagorinsky, <.> means averaged):
[6] nu_SGS = cD * delta^2 * sqrt(S_ij S_ij)
[7] cD = 0.5 (<L_ij M_ij>)/(<M_kl M_kl>)
[8] S_ij = D_ij = S_ij Pope
[9] M_ij = delta^2 * (F(sqrt(S_kl S_kl) Sij) - 4 * sqrt(<S_kl> <S_kl>) <S_ij>)
[10] L_ij = F(ui uj) - F(ui) F(uj)

I marked the differences of these models.

- First difference is the factor 0.5 in Eq[7] in comparison to Eq[2]. This comes from the factor 2 in Eq[4]. If we put this in Eq[2] we get 0.5 (ok)

- Second difference is the different filtering in M_ij. What effect does this have???? (x)

- Third difference is the factor 4 in Eq[9] in comparison to Eq[4]. Bernhard:
here the factor 4 is from the double width of the test filter
can you please explain that?? Why we don`t get a fector in the other filtered terms??? (x)

- Fourth difference is the factor 2 in mag(S_ij). (x)

In my opinion these models are different or I`m not able to bring the OF-model in the form of pope-model. If it is possible, anyone can please give some advice?

kind regards
Florian

wenxu

October 28, 2014 01:39

Quote:

Originally Posted by lakeat (Post 220718)

Sorry, I do not understand, I saw in "Smagorinsky.H",

Code:

tmp<volScalarField> k(const tmp<volTensorField>& gradU) const

{

    return (2.0*ck_/ce_)*sqr(delta())*magSqr(dev(symm(gradU)));

}

Code:

nuSgs_ = ck_*delta()*sqrt(k(gradU));

Code:

tmp<volScalarField> k(const tmp<volTensorField>& gradU) const

{

    return (ck_/ce_)*sqr(delta())*magSqr(dev(symm(gradU)));

}

Thank you

$\begin{array}{l} {\nu _{SGS}} = {\left( {{C_S}\Delta } \right)^2}\left| {{\bf{\bar S}}} \right| \\ K = {\left( {{C_I}\Delta } \right)^2}{\left| {{\bf{\bar S}}} \right|^2} \\ \left| {{\bf{\bar S}}} \right| = {\left( {{\bf{\bar S}}{\rm{:}}{\bf{\bar S}}} \right)^{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-\nulldelimiterspace} 2}}} \\ \end{array}$

the third one shoub be sqrt(2*Sij*Sij):o

babakflame

December 11, 2014 05:40

How to change ck and ce in Smagorinsky approach

Greetings All

I have performed a Smagorinsky-based compressible LES simulation with the Coefficients as follows:

HTML Code:

{

  ce = 1.048;

  ck= 0.02;

}

According to the relation:

${\left( {{C_S}} \right)^2} = {C_K}\sqrt {\frac{{{C_K}}}{{{C_\varepsilon}}}}$

I will get

.
My case is a reacting non-premixed combustion with a bluff-body separating fuel and oxidizer streams.

I want to change the

value into 0.13. which value between

and $C_\varepsilon$ should be modified to retain the nature of the problem?

Best,
Bobi

openfoammaofnepo

March 21, 2015 18:27

Actually in this paper, we cannot find the information about how the model constants c_k=0.094 and c_{\epsilon}=1.048 come out. So which one is correct reference when I use these two constants? Thanks.

Quote:

Originally Posted by alberto (Post 304757)

Reference: http://pof.aip.org/resource/1/phfle6/v9/i5/p1416_s1

lyne00

September 21, 2015 10:11

Openfoam

Hi Lakeat,now I want to write the Scalar SikSkj,but how to write that in openFoam, Thanks.

Elyas_Mosibat

January 15, 2016 12:56

Change the smagorinsky coeff

hi friends
How can I change the smagorinsky Coeff (Cs) in OF?I want to use Cs=0.1...
Regard;)

Jingxue Wang

December 26, 2017 21:16

I read that all the previous posts and reply.

May I understand that Cs=sqrt(Ck*sqrt(2Ck/Ce)) is correct in OpenFoam?

Thanks a lot

manafaero

March 26, 2019 05:15

dear madam
As you know the smagorinsky constant Cs in OpenFOAM is defined as a function of Ce and Ck. But when i do the simulation for various combinations of Ck and Ce, that gives the same Cs, I am getting different results. So could you please tell Me where else Ck and Ce is used other than in smagorinsky model. My problem is incompressible LES with standard smagorinsky model and vandriest delta.

please help

my mailId is manafaero@gmail.com

wyldckat

March 31, 2019 15:24

Quote:

Originally Posted by manafaero (Post 728865)

But when i do the simulation for various combinations of Ck and Ce, that gives the same Cs, I am getting different results. So could you please tell Me where else Ck and Ce is used other than in smagorinsky model.

"Ce" is loaded here: https://github.com/OpenFOAM/OpenFOAM...ddyViscosity.C
"Ck" is loaded and used here; "Ce" is also used here: https://github.com/OpenFOAM/OpenFOAM.../Smagorinsky.C

The nut, k and epsilon fields are created by these classes, but I 'm not familiar with LES modelling. Either way, if you look for "Ce_" and "Ck_" in those files, you will see how they are used.

manafaero

April 3, 2019 08:54

OpenFoam

I would like to share some experience of mine regarding the discussion. I also verified the equations given in the discussion. But when i run the simulation for various ck and ce combinations I am getting different results. the nu sgs is not the same in all cases. so eventhough the formula is correct the openFoam is not calculating as per the discussion. Note that I had removed the dependence of k in calculating shear stress B. but still the results differ. So be careful before using this. Also I am a beginner to OpenFOAM and C++. If any body is interested and any comments please post. my mail ID is manafaero@gmail.com.

jorkolino

September 16, 2020 11:30

did you check this page with OF implementation https://caefn.com/openfoam/smagorinsky-sgs-model

All times are GMT -4. The time now is 06:57.