Dynamic Smagorinsky in Fluent: why is clipping needed?

Marvin · December 11, 2012, 06:08

Dear all,

A small question concerning Large Eddy Simulation and the implementation of the dynamic Smagorinsky subgrid-scale model in Fluent (Lilly's version).

In this procedure, the value of the Smagorinsky coefficient (C_s) is clipped to [0; 0.23].

For the lower bound, it's actually the square of C_s (and hence the subgrid-scale viscosity) which is kept positive to avoid numerical instabilities.

What about the upper bound? Why is it needed to keep C_s below a maximum value? Is it also due to numerical issues? And why this particular value (0.23)?

Any idea or comment is welcome. Thanks!

LuckyTran · December 11, 2012, 13:52

Good question!

0.23 is the value of the Smagorinsky coefficient used in standard method (non-dynamic) as derived by Lily for homogeneous isotropic turbulence. It is known (by observation) that the value of 0.23 is too high for flows with mean shear and other effects. Therefore in general, C_s must be less than 0.23 since that is the constant derived for the homogeneous isotropic case (all other flows must have a smaller C_s). Hence for semi-empirical reasons, the case of C_s > 0.23 is also non-physical.

Recall that C_s captures the influence of the subgrid scale. If C_s > 0.23 is needed then that is like saying the contribution of the subgrid scales is very significant, which means that subgrid scale model is not appropriate in the first place.

Hence, nearly all codes clip C_s at 0 and 0.23 since having values outside this range is either 1) non-physical or 2) not in the spirit of a subgrid scale model or both 1) and 2).

sbaffini · December 13, 2012, 15:06

I have to slightly disagree on the explanation.

First of all, the value of Lilly is slightly different. It is 0.17 but only for the assumption of a spectral cut-off filter and turbulence with an infinite inertial range (k^-5/3 all over). For different filters the constant is different (but never 0.23).

It is common practice to limit the constant to the value that, for the specific code and the static Smagorinsky model there implemented, can correctly reproduce the decay of homogeneous isotropic turbulence (a value which can be different for different implementations and/or codes). As a matter of fact, you can find the value 0.23 only in the Fluent manual and not in the AIAA paper 2004-2548 promoting the implementation of LES in Fluent.

However, there is no theoretical reason to assume the value for isotropic turbulence decay as the limiting factor (also because there is no theoretical reason for an eddy viscosity at all). Actually, the idea of a dynamic constant is to calibrate the dissipation with the estimated flux of energy toward the smallest scales (that's the actual meaning of the test filtering in the dynamic smagorinsky model); it is very easy to understand that if i limit the constant for an inertial range spectrum, than all the cases where the energy flux is higher are then underestimating the constant.

Also, while the procedure for computing the constant is dynamic, it is also very badly conditioned leading most of the times to very high values which, positive or not, have little sense (the fact that you minimize the error in the germano procedure doesn't mean that you make that 0... actually it is quite high).

This leads to the final part of the motivation, you can have strong spatial variations which do not help the solver (also because they are used in part as explicit).

The limiting of the constants is actually more an art than a science and slight variations in the clippings can also give enormous differences. Indeed, it is very likely that you will never read these values in any paper.

LuckyTran · December 13, 2012, 20:01

sbaffini made some very valid points. My explanation seemed to suggest that C_s < 0.2-0.23 is universally valid, this is not the case.

Quote:

Originally Posted by sbaffini

Also, while the procedure for computing the constant is dynamic, it is also very badly conditioned leading most of the times to very high values which, positive or not, have little sense (the fact that you minimize the error in the germano procedure doesn't mean that you make that 0... actually it is quite high).

The numerical procedure implemented for calculating C_s is numerically unstable and leads to values that we know (from observation & experience) do not make much sense. Hence some clipping is needed.
It is true that it is more of an art than a science. We should not be too worried about the exact clip (whether it is 0.2 or 0.23). Through the various explanations, it is possible to arrive at the conclusion that a C_s of 0.2 - 0.23 is WAY TOO HIGH and that a value of C_s much lower is necessary.

Quote:

Originally Posted by sbaffini

However, there is no theoretical reason to assume the value for isotropic turbulence decay as the limiting factor (also because there is no theoretical reason for an eddy viscosity at all).

I agree that in principle that there is no need to assume an eddy viscosity model. But in LES with SGS models (at least in the Smagorinsky model) the eddy viscosity model is invoked explicitly. The Smagorinsky constant (whether static of dynamic) is used to calculate the turbulent eddy viscosity.

The correct value of C_s is the value that predicts the correct eddy viscosity (regardless of sign or magnitude). But in the spirit of LES, the contribution of the filtered subgrid scales should be small and hence the value of C_s should be appropriately small (otherwise we are not doing good LES).

Marvin · December 14, 2012, 09:42

Thank you very much for your replies!
A few points to continue the discussion.

** I had the feeling that maintaining C_s below a maximum value was a numerical trick. Why: if too high values are not physical - let's assume this for now - then there is no reason to find such high values at the scales at which C_s is evaluated with the dynamic procedure (i.e. the smallest resolved scale, those between the filter and the test filter). And in this case there would be no need for clipping.

** However, while it seems rather logical that negative values of nu_sgs (the subgrid-scale viscosity) can lead to instabilities, high values of C_s - even if meaningless - lead to a high nu_sgs which in fact should somehow stabilize the flow, shouldn't it?

** On the C_s value for homogeneous isotropic turbulence:

I also had in mind the theoretical value of C_s=0.17.
Unfortunately Lilly's papers from the 60's are rather difficult to access. His results are well summarized in this very interesting article, regarding our topic:

J.W. Deardorff, On the magnitude of the subgrid scale eddy coefficient. Journal of Computational Physics 7, 120-133 (1971).

Briefly: 0.176 is indeed the value. Hypotheses: eddy-viscosity/Smagorinsky model; filter width in inertial subrange/local equilibrium between production and dissipation; spectral cut-off filter.
However, it also appears that this value has to be increased when the characteristic strain rate is evaluated with a finite difference method, leading to a value of 0.22. This last point is still unclear to me but the value of 0.23 would make sense regarding this point...

** On the relevance of using homogeneous isotropic turbulence as a limiting case:

The energy spectrum is affected by the presence of shear or buoyancy in the flow, which leads to a decrease in the value of C_s and explains why lower values than Lilly's are needed to accurately simulate "real" flows. (That's actually the same reason why damping functions are needed to decrease C_s close to walls.)
See for example:

V.M. Canuto and Y. Cheng, Determination of the Smagorinsky-Lilly constant C_s. Physics of Fluids 9, 1368 (1997).

** So if I try to summarize:

Very high C_s values are not physical but they can appear in the numerical flow due to the least-square method used to compute this coefficient. These high values "do not help the solver" so we get rid of them by imposing a maximum value. 0.23 is just an arbitrary number, rather close to the maximum values that can be encountered.

sbaffini, can you elaborate a bit more on why these high values are not good for the solver?

It would definitely be interesting to have Fluent people's point of view on this question...

Thanks again.

sbaffini · December 14, 2012, 19:42

Dear Marvin,

a high positive Cs (hence turbulent viscosity) can be positive as long as it is treated implicitly. However, this is rather hard in dynamic models (or any turbulent model, for that matter) and the usual cure is to use an estimate of it to be treated implicitly and a correction which is made explicitly (however, i don't know if it is really done in Fluent).

Now, imagine using 0.23 in Cs for the implicit estimate and then getting 5 (or even 50, in my case i usually get up to 200-300) for the real Cs. If you treat the correction explicitly, well, it is just a matter of numerical stability, it will blow up (imagine discretizing a diffusion equation explicitly, your time step has to satisfy some very stringent criteria, e.g., dt<K*rho*dx^2/(mu_sgs_exp) with K some constant, say 1/2).

This is not the only reason; you have to consider also the spatial variability of the constant. This can cause problems to the solver too if it is too severe; if you don't limit the constant, jumps of two order of magnitudes are likely to occur to your viscosity just from one cell to another.

The way the Smagorinsky constant is computed is very well explained in the book of Pope, in the related chapter.

For what concerns the theoretical maximum value, you can read the work of Yoshizawa: Subgridscale modeling with a variable length scale, Physics of Fluids A 1, 1293 (1989). There, a certain theoretical approach is shown to give:

Cs/Cs_0 = 1 + (Ca/S^2) * DS/Dt

where Cs is the constant, Cs_0 is 0.16 (classical value), Ca is 0.64 and S is the strain rate magnitude (or a very similar quantity, according to the paper). So, no matter what, this theory gives you Cs>Cs_0 as long as DS/Dt >0 (here D/Dt is the material derivative).

sbaffini · December 16, 2012, 07:42

Just to add an additional information.

For Code_Saturne, the open source code developed by EDF, the clipping of the constant for the Smagorinsky model is on 10*Cs_0 (page 156 of the user manual).

So, as you can see, different codes use different approaches... mostly because there is no definitive answer on this matter nor theoretical justification.

December 11, 2012, 06:08	Dynamic Smagorinsky in Fluent: why is clipping needed?	#1
Marvin New Member Join Date: Nov 2009 Posts: 10 Rep Power: 17	Dear all, A small question concerning Large Eddy Simulation and the implementation of the dynamic Smagorinsky subgrid-scale model in Fluent (Lilly's version). In this procedure, the value of the Smagorinsky coefficient (C_s) is clipped to [0; 0.23]. For the lower bound, it's actually the square of C_s (and hence the subgrid-scale viscosity) which is kept positive to avoid numerical instabilities. What about the upper bound? Why is it needed to keep C_s below a maximum value? Is it also due to numerical issues? And why this particular value (0.23)? Any idea or comment is welcome. Thanks!

December 11, 2012, 13:52		#2
LuckyTran Senior Member Lucky Join Date: Apr 2011 Location: Orlando, FL USA Posts: 5,754 Rep Power: 66	Good question! 0.23 is the value of the Smagorinsky coefficient used in standard method (non-dynamic) as derived by Lily for homogeneous isotropic turbulence. It is known (by observation) that the value of 0.23 is too high for flows with mean shear and other effects. Therefore in general, C_s must be less than 0.23 since that is the constant derived for the homogeneous isotropic case (all other flows must have a smaller C_s). Hence for semi-empirical reasons, the case of C_s > 0.23 is also non-physical. Recall that C_s captures the influence of the subgrid scale. If C_s > 0.23 is needed then that is like saying the contribution of the subgrid scales is very significant, which means that subgrid scale model is not appropriate in the first place. Hence, nearly all codes clip C_s at 0 and 0.23 since having values outside this range is either 1) non-physical or 2) not in the spirit of a subgrid scale model or both 1) and 2). Far, Ebrahim, ghost82 and 1 others like this.

December 13, 2012, 15:06		#3
sbaffini Senior Member Paolo Lampitella Join Date: Mar 2009 Location: Italy Posts: 2,192 Blog Entries: 29 Rep Power: 39	I have to slightly disagree on the explanation. First of all, the value of Lilly is slightly different. It is 0.17 but only for the assumption of a spectral cut-off filter and turbulence with an infinite inertial range (k^-5/3 all over). For different filters the constant is different (but never 0.23). It is common practice to limit the constant to the value that, for the specific code and the static Smagorinsky model there implemented, can correctly reproduce the decay of homogeneous isotropic turbulence (a value which can be different for different implementations and/or codes). As a matter of fact, you can find the value 0.23 only in the Fluent manual and not in the AIAA paper 2004-2548 promoting the implementation of LES in Fluent. However, there is no theoretical reason to assume the value for isotropic turbulence decay as the limiting factor (also because there is no theoretical reason for an eddy viscosity at all). Actually, the idea of a dynamic constant is to calibrate the dissipation with the estimated flux of energy toward the smallest scales (that's the actual meaning of the test filtering in the dynamic smagorinsky model); it is very easy to understand that if i limit the constant for an inertial range spectrum, than all the cases where the energy flux is higher are then underestimating the constant. Also, while the procedure for computing the constant is dynamic, it is also very badly conditioned leading most of the times to very high values which, positive or not, have little sense (the fact that you minimize the error in the germano procedure doesn't mean that you make that 0... actually it is quite high). This leads to the final part of the motivation, you can have strong spatial variations which do not help the solver (also because they are used in part as explicit). The limiting of the constants is actually more an art than a science and slight variations in the clippings can also give enormous differences. Indeed, it is very likely that you will never read these values in any paper. Far likes this.

December 14, 2012, 19:42		#6
sbaffini Senior Member Paolo Lampitella Join Date: Mar 2009 Location: Italy Posts: 2,192 Blog Entries: 29 Rep Power: 39	Dear Marvin, a high positive Cs (hence turbulent viscosity) can be positive as long as it is treated implicitly. However, this is rather hard in dynamic models (or any turbulent model, for that matter) and the usual cure is to use an estimate of it to be treated implicitly and a correction which is made explicitly (however, i don't know if it is really done in Fluent). Now, imagine using 0.23 in Cs for the implicit estimate and then getting 5 (or even 50, in my case i usually get up to 200-300) for the real Cs. If you treat the correction explicitly, well, it is just a matter of numerical stability, it will blow up (imagine discretizing a diffusion equation explicitly, your time step has to satisfy some very stringent criteria, e.g., dt<Krhodx^2/(mu_sgs_exp) with K some constant, say 1/2). This is not the only reason; you have to consider also the spatial variability of the constant. This can cause problems to the solver too if it is too severe; if you don't limit the constant, jumps of two order of magnitudes are likely to occur to your viscosity just from one cell to another. The way the Smagorinsky constant is computed is very well explained in the book of Pope, in the related chapter. For what concerns the theoretical maximum value, you can read the work of Yoshizawa: Subgridscale modeling with a variable length scale, Physics of Fluids A 1, 1293 (1989). There, a certain theoretical approach is shown to give: Cs/Cs_0 = 1 + (Ca/S^2) * DS/Dt where Cs is the constant, Cs_0 is 0.16 (classical value), Ca is 0.64 and S is the strain rate magnitude (or a very similar quantity, according to the paper). So, no matter what, this theory gives you Cs>Cs_0 as long as DS/Dt >0 (here D/Dt is the material derivative). Far likes this.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Guideline needed for Fluent 6.3	CCSC26	FLUENT	1	January 29, 2012 13:03
Dynamic Mesh in Fluent	haris05	Main CFD Forum	2	December 5, 2011 02:54
Dynamic Mesh - FLUENT 6.3	Lourival	FLUENT	2	February 20, 2008 13:47
Problem with Dynamic Mesh and Fluent	Claude	Main CFD Forum	0	December 12, 2006 07:03
Dynamic Mesh Fluent	gianluca	Main CFD Forum	3	December 13, 2004 12:09

December 14, 2012, 09:42		#5
Marvin New Member Join Date: Nov 2009 Posts: 10 Rep Power: 17	Thank you very much for your replies! A few points to continue the discussion. I had the feeling that maintaining C_s below a maximum value was a numerical trick. Why: if too high values are not physical - let's assume this for now - then there is no reason to find such high values at the scales at which C_s is evaluated with the dynamic procedure (i.e. the smallest resolved scale, those between the filter and the test filter). And in this case there would be no need for clipping. However, while it seems rather logical that negative values of nu_sgs (the subgrid-scale viscosity) can lead to instabilities, high values of C_s - even if meaningless - lead to a high nu_sgs which in fact should somehow stabilize the flow, shouldn't it? On the C_s value for homogeneous isotropic turbulence: I also had in mind the theoretical value of C_s=0.17. Unfortunately Lilly's papers from the 60's are rather difficult to access. His results are well summarized in this very interesting article, regarding our topic: J.W. Deardorff, On the magnitude of the subgrid scale eddy coefficient. Journal of Computational Physics 7, 120-133 (1971). Briefly: 0.176 is indeed the value. Hypotheses: eddy-viscosity/Smagorinsky model; filter width in inertial subrange/local equilibrium between production and dissipation; spectral cut-off filter. However, it also appears that this value has to be increased when the characteristic strain rate is evaluated with a finite difference method, leading to a value of 0.22. This last point is still unclear to me but the value of 0.23 would make sense regarding this point... On the relevance of using homogeneous isotropic turbulence as a limiting case: The energy spectrum is affected by the presence of shear or buoyancy in the flow, which leads to a decrease in the value of C_s and explains why lower values than Lilly's are needed to accurately simulate "real" flows. (That's actually the same reason why damping functions are needed to decrease C_s close to walls.) See for example: V.M. Canuto and Y. Cheng, Determination of the Smagorinsky-Lilly constant C_s. Physics of Fluids 9, 1368 (1997). ** So if I try to summarize: Very high C_s values are not physical but they can appear in the numerical flow due to the least-square method used to compute this coefficient. These high values "do not help the solver" so we get rid of them by imposing a maximum value. 0.23 is just an arbitrary number, rather close to the maximum values that can be encountered. sbaffini, can you elaborate a bit more on why these high values are not good for the solver? It would definitely be interesting to have Fluent people's point of view on this question... Thanks again.

December 16, 2012, 07:42		#7
sbaffini Senior Member Paolo Lampitella Join Date: Mar 2009 Location: Italy Posts: 2,192 Blog Entries: 29 Rep Power: 39	Just to add an additional information. For Code_Saturne, the open source code developed by EDF, the clipping of the constant for the Smagorinsky model is on 10*Cs_0 (page 156 of the user manual). So, as you can see, different codes use different approaches... mostly because there is no definitive answer on this matter nor theoretical justification.