Explanation bounding k and bounding epsilon
 

Hey Foamers,

I have no idea what is bounding k and bounding epsilon supposed to mean. Can anybody explaing, please? Druing my simulation I am getting such messages with even negative values.

Cheers in advance

Hi Claus,

as i understood it, bounding means, that the actual calculated values are limited not to get too huge. bounding then should mean that the calculated value was the bounded amount higher than it is being used for the next iteration/timestep.
hope this is right, if not please anyone correct me :).

bye,
mo

Bound the given scalar field if it has gone unbounded below a lower bound.

Having really large values for k/epsilon may not be great, but avoiding divide-by-zero is the real issue:
In any case, avoiding zero or near zero values for k/epsilon is required, which is what bound() does.

The source code is found under src/finiteVolume/cfdTools/general/bound/bound.C

Thanks for your replies!

When the minimum bound value for k is e.g. -136 do I have to worry about something? I think the value should never fall below 0.

Cheers

The k/epsilon are bounded with SMALL (1e-15 for double) so that will stop the divide-by-zero issues. If your problem always produces negative values throughout the solution, it could indicate that something else needs attention - eg, system needs more relaxation or something.

It was obvious to me that something is wrong! Okay, I used Gauss linear for all convection terms, because I wanted to keep accuracy as best as possible. I changed to limitedLinear(V) 1; know the negative min values increases toward to zero and, eventually, they vanished.

Cheers

Now, there is no occurrence of bounding k or epsilon anymore. However, it crossed my mind: can I use a better scheme instead of limitedLinear 1 for the approximation of epsilon and k, respectively?

To explain your problem, we need to look at the boundness of the differencing scheme of the convection term. Upwind differencing scheme usually gives bounded results, but introduces a numerical diffusion which is not good for the accuracy. However, the second-order central differencing scheme, i.e. Gauss linear, does not assure the boundness of the result. That is why you got a unbounded result in your earlier calculation.

The limitedLinear scheme you are using now is a blended one of the upwind and central differencing, I guess. (I am now considering the same problem to increase the accuracy of the calculation) So you can obtain a bounded result with it. The limitedLinear scheme may have a accuracy between first and second order???

I am looking forward to a discussion about this.

Thanks for reading.

Hallo JinBiao!

Actually, I am not really sure about the usage of lmitedLinear 1; I spotted it in a tutorial. Basically, I tend to use schemes which are both of high order and bounded. However, as we both know there a no higher order bounded schemes, therefore, I would like to know which schemes gives the best results? Which e.g. TVD scheme would you prefere? I simulation two-phase flows (water-air) through bendend pipe with a maximal velocity of 9.6 m/s. Additionally, in this pipes there holes where the water can flow out.

Cheers in advance

Hi there,

I was also struggling to find out a scheme for the convection term. But now I think no method will solve all the problems. One solution will be refining your mesh at the location where you are expecting a peak value or a large gradient. It improves the boundness of the result, I think. The reason is, from my understanding, that the face value is more inclined to fall into the range between the value of the two neighbor points, if your mesh is refined.

The scheme I am now using is SFCD and QUICK. I am not sure it is good for you. I have no experience with two phase flow simulation, so I cannot provide more information.

Good luck.

Hello,

just a quick clarification. The limitedLinear scheme is the linear scheme with a limiter (if I remember correctly, it uses Sweby's limiter). It's not a blended version of the "linear" scheme and the "upwind" scheme, which would be "linearUpwind" in OpenFOAM. So limitedLinear is a TVD scheme (take a look at ../OpenFOAM-1.6.x/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/limitedLinear), as also written in the user's guide.

What schemes do you usually like to use?

For your application, you can either use the default works fine. You can try linearUpwind / linearUpwindV with the cellLimited option to have what many would call "second order upwind".

Best,

Do you notice significant improvements compared to limitedLinear or linearUpwind with these schemes?

In my experience, QUICK doesn't offer any particular advantage (it's either second or first order on general grids, more diffusive than limitedLinear).

Best,

In my case, QUICK converges better than TVD and NVA schemes. I am not quite sure about the reason. One of the explanations can be NVA and TVD has artificial switch which cause non-smooth change of the coefficient matrix.

This can be particular in my case because I am simulating a backward-facing step flow. There is large gradient at the edge of the step.

Regards,

Interesting, because QUICK is known to have some trouble due to the possibility of having coefficients that change sign.

Did you try to estimate the influence of the numerical diffusion (QUICK is quite diffusive for what I know). This could explain the improved stability.

Hi Alberto,

Up to now, I have not analysed the numerical diffusion because it is difficult to find out a scheme to compare with. I think an ideal counterpart for comparison should be CDS, with which it is hard to obtain the bounded result. Also, I think TVD schemes, such as limitedLinear, have a better boundness because they introduce numerical diffusion. So it does not make sense to compare with them.

It is also good to know there are some defects in QUICK. I will be more careful when I am using it. I also tried to find some comments on QUICK in the textbooks. Mainly, people think QUICK has about the same accuracy with CDS because the integration on the surface is usually second-order in accuracy.

I have another question about TVD schemes. Does the prefix "limited" mean that Sweby's limit function is applied to the scheme indicated by the second half of the name. For example, limitedLinear is a limited CDS scheme?

Thanks for your discussion.

Regards,

Correct. However also upwinded schemes (QUICK, for instance) are diffusive by nature ;-)

To compare with CDS you don't have to run a huge practical case. Just compare in a realistic but small test case, where you don't have problems to use the appropriate discretization. After all it's not so hard to obtain a solution with CDS, being it one of the schemes used in LES due to its central nature (=> no dissipation)

I know a lot of people uses QUICK, also for LES, claiming it is "proven to give good results". In practice, if you analyze the dissipative error, you notice that it is significant, meaning the scheme is quite diffusive.
QUICK might have the same formal order of accuracy (not always), but CDS have zero dissipation error (only the dispersion error is present in central schemes).

Yes, limitedLinear is a limited CDS, as the name says. It's not a blend of upwind and linear.

About the implementation of the limiter, OF allows the implementation of generic limiter functions, so take a look at

~/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes

where you find the implementation of the abstract limited class and of the schemes.

Best,

Hallo alberto!

Actually, I really like to use central differencing, but I am facing ocillations when I am using it. Therefore, I used linearLimited(V) 1 successfully, however, there sill a question: Is there yet a better scheme?

Since limiters switches between UD and CD, I would like to use a limiter which switches in that moment to UD when it is really necessary. As far as I understand, the limiter switches depending on "r" and I think that limiters switches often to early. Is this correct?

Cheers

Claus

Thanks a lot for your patient explanation.

I am sorry but I really got confused about the words, diffusion, dispersion and dissipation.

I have not read the code quite in detail. It just come out in my mind that the limit function introduces nonlinearity which can be the reason for poor convergence of TVD scheme. How do you think?

But I still think the limitedLinear scheme falls into the blend of upwind and linear. Please see Hrvoje's doctor thesis page 98 eq.(3.66). The limit function performs as the blend factor. However, the blend factor goes larger than one, according to Sweby's constraint.

Thanks.

It's bit technical to explain in detail, so please refer to the literature. In short, the terms dissipation error/diffusion are used with to indicate the same thing, meaning the numerical diffusion/viscosity added by the scheme. This error is due to the presence of non-zero imaginary part in the modified wavenumber. The real part of the wavenumber is responsible of what is called dispersion error. Upwinded schemes are affected by both, while central schemes are only affected by the dispersion error.

As far as I remember, the limitedLinear scheme was introduced in OpenFOAM later, and Hrv worked on the Gamma scheme (NVD).

Best,

Hi Foamers,

as i also try to tune the solvers in the fvSchemes at the moment, i'd be great if you might answer if i the system works like this:

in the divSchemes dictionary i defined like that:

div(phi,U) Gauss upwind Gauss limitedLinear 0.7;

is it right that is is defined like that:

Gauss <interpolationScheme>
Gauss upwind <gradScheme>
Gauss upwind Gauss <interpolationScheme>
Gauss upwind Gauss limitedLinear 0.5


<interpolationScheme>
-> upwind <gradScheme>

<gradScheme>
-> Gauss <interpolationScheme>
<interpolationScheme>
-> limitedLinear 0.5



Another question: how can scalars be strictly bounded? As i understood the manual, the specific interpolation schemes are not defined in the interpolationSchemes dictionary, but in the divSchemes dictionary.
I dont understand how to limit scalar values, when divergence is calculated.
Would it maybe work like hits?

div(phi,k) Gauss limitedLinear Gauss linear 1e-10 30;

All the best,
Moritz

The gradient scheme is used only with some of the schemes, and OpenFOAM will return an error if the syntax is not appropriate for the scheme you are using.

The upwind schemes does not require anything:

div(phi, U) Gauss upwind;

The gradient scheme is used, for example, by the linearUpwind scheme:

div(phi, U) Gauss linearUpwind Gauss linear;

Depending on what you need to do, you can either use one of the bounded schemes between 0 and 1, if that's the range you need, or do as explained on the user's guide, using the prefix "limited" and specifying the range.

Best,

hi alberto,
thanks a lot for your reply.
i want to limit k and epsilon and dont know how to make that.

e.g. for k, would it work like that?

div(phi, k) Gauss linearUpwind Gauss limitedLinear 1e-10 30;

if not, could you give me an example?

thanks a lot,
Moritz

Are their values actually limited between the bounds you specify? Meaning, is your model including that limitation? If yes, it is correct to use a limiter on the scheme to keep the quantity bounded whenever it goes out of those bounds due to numerical reasons.

If the equations do not include that physical limit, I would not use a limiter, but a bound, as done in the codes already.

The syntax, explained at page U-113 of the user's guide, is, for example:

div(phi, U) Gauss limitedVanLeer a b;

where a and b are the bounds of your range [a,b]. Schemes that support the limitation in this way are listed in the same page.

Best,

I want to find out a drag coefficient with the LamBremhorst-model, and i have the problem that epsilon tends to explode without any obvious reason inside the volume mesh (not close to walls). So if i want to limit epsilon to a value between 1e-30 and 30, do i have to do that in both the divScheme and the laplacianScheme?

What i dont understand is, that in these Schemes you calculate divergence and divgrad, but not the scalar values themselves, so can you limit the values themselves there?



Is it right, that the bounding limits the change of a value during one interation? That means in a few timesteps it can diverge although values are bounded?

in this example, would i limit the magnitude of the velocity?
what parameters have to be given, when you use e.g. limitedVanLeerV?

Thank you so much for your time!!

You clearly have a problem in the setup of your case if some variable tends to explode. Limiters to impose arbitrary ranges won't be of much help if you do not understand the reason why this is happening.

Let's try to understand what is the cause. What solver are you using, what turbulence model, and what boundary condition setup?

Best,

Hi Alberto, I have a similar problem with my case using kOmegaSST model. However, I ran many different but similar meshes of this external flow and they all converge. Only this case diverges on k. All the BCs and other parameters are identical.

Also, the case runs fine in laminar, although it doesn't converge...

Could my mesh be the only cause of the divergence?


Thanks,


-Louis

Let me guess: is the mesh tetrahedral (maybe generated with Salomé)?

Please run checkMesh, and if you have very skewed cells and/or a high degree of non-orthogonality, yes, the mesh can be the reason of the problem.
In this case you can easily see where this happens by saving solutions right before the crash and looking for peaks of values in the diverging variables.

Best,

Hi Alberto,

it is a hexmesh, from snappyHexMesh. Yes, it does have skewed faces but those skewed faces didn't seem to create a problem on my other meshes (well at least the bound on k worked for those other meshes). Thanks for your advice, I will try to slightly modify my snappyHexMesh parameters to get a slightly different mesh!

Cheers,

-Louis

OK. If you used snappyHexMesh, it is probably easier to fix than in tet meshes. Check in particular around surfaces, because sometime snappy adds some very thin cell there :)

Best,

I got it to converge by modifying the mesh on the triSurface file.stl I gave to snappyHexMesh.. Not the best solution but works as a fix :)

Thx,

-Louis

MRFSimpleFoam bounding k & epsilon
 

Hi Foamers,
I have a great problem using MRFSimpleFoam with OF 1.6.
I'm studing (fluent Vs OpenFOAM) a 3D fan (mesh imported from GAMBIT) using the k-epsilon model and upwind schemes for convective fluxes, the under relaxation factors are 0.3, 0.7, 0.5, 0.5 for pressure, velocity, k and epsilon.
After 5600 iterations my OF residuals are very low and the message "bounding epsilon..." is printed very few times but the solution is far from the one obtained by Fluent using UDS. Trying to relaunch OF case from the 5600th iteration I faced the k and epsilon diverging and I haven't changed settings. I'm struggling to know what's going wrong.Can you please help me?
Is there anybody who know about bugs in MRFSimpleFoam in OF 1.6?
How to let k-epsilon converge?
Thanks
Aldo

Hi,

we do not have enough details to answer.

P.S. Please, when you have a question, open a new thread. Thanks! :-)

Resume:

MRFSimpleFoam 1st order upwind-> convergence of u, p, k, epsilon untill 5800 iterations-> stop and restart-> residuals step-up (3 order) and divergence of k and epsilon (under zero) without any change in settings (I thought to get crazy but it's true).

MRFSimpleFoam 2nd order upwind-> divergence of k and epsilon (under zero), It seem there is not the possibility using other discretization schemes.

How to set up a 2nd order disc. scheme without any problem?
Is there anybody who know have tested MRFSimpleFoam comparing it to fluent? How arethe resulting fields?

Thanks
aldo

How did you set a "second order upwind" scheme in OpenFOAM?

You might want to try

• For gradients: cellLimited Gauss linear 1;
• For divergence: Gauss linearUpwind cellLimited Gauss linear 1; (use linearUpwindV for div(phi,U) and other vector fields, with cellMDLimited).
• Limit also laplacians and surface normal gradients.

P.S. What k-epsilon are you using, and what boundary conditions did you specify?

K-epsilon standard, zerogradient on the outlet/inlet where the pressure is fixed to zero.

Sorry, I'm quite new in OF, can you please correct my fvSchemes dict (expecially the laplacian schemes)?
__________________________________________________ ______

ddtSchemes
{
default steadyState;
}

gradSchemes
{
default cellLimited Gauss linear 1;
grad(p) cellLimited Gauss linear 1;
grad(U) cellLimited Gauss linear 1;
}

divSchemes
{
default none;
div(phi,U) Gauss linearUpwindV cellMDLimited Gauss linear 1;
div(phi,k) Gauss linearUpwind cellLimited Gauss linear 1;
div(phi,epsilon) Gauss linearUpwind cellLimited Gauss linear 1;
div((nuEff*dev(grad(U).T()))) Gauss linearUpwind cellLimited Gauss linear 1;//Gauss linear;
}

laplacianSchemes
{
default none;
laplacian(nuEff,U) Gauss linear corrected;
laplacian((1|A(U)),p) Gauss linear corrected;
laplacian(DkEff,k) Gauss linear corrected;
laplacian(DepsilonEff,epsilon) Gauss linear corrected;
}

interpolationSchemes
{
default linear;
interpolate(U) linear;
}

snGradSchemes
{
default corrected;
}

fluxRequired
{
default no;
p ;
}
__________________________________________________ ___


Thanks

MRF Divergence
 

Hi Alberto,
I tried to follow your advices but It diverges after 4000 iter. and I don't understand why, can you please have a look at the attached log and fvSchemes files in order to correct my mistakes?
Thanks
Aldo

Are there zones of the system where you can observe some unphysical behaviour (velocity too high, for example?
Additionally, what under-relaxation factors are you using? Did you try to reduce them for k and epsilon?

P.S. You wrote you use OpenFOAM 1.6. Did you update it to 1.6.x?

Best,

How can I bound own Scalar value in equations?
 

Thanks all for information about bounded schemes
I need a method to bound own volScalarField that is solved in an equation has div , laplacian ,...

for example C1 in the concentration equation must be between [0,1]
anybody can tell me how can i add this limiter to the equation and force C1 to be in [0,1]?
thank you.

So just to be clear,
for example say the bounding value for epsilon is 3 and the calculated was to be 2000, so the bounding epsilon would show up to be 1997? and so in the next iteratio they will use epsilons value as 3 instead of the calculated 2000?

accuracy of simpleFoam - bounding epsilon/k
 

Hello,

I am working on a simulation of wind on buildings with a complex Mesh.
I am using a RASModel kEspilon with the simpleFoam solver.
I am applying ABL conditions (atmospheric boundary layer).

I am facing problems with bounding K and bounding epsilon.
the problem of bounding k and epsilon is caused by skew faces in my mesh.
I searched for schemes that can fix the problem, and found a modification that improved my simulation:
changing "laplacianSchemes" from "Gauss linear limited 1" to Gauss linear limited 0.333
changing "snGradSchemes" from "limited 1" to limited 0.333

simpleFoam ran for 800 steps before the "bounding k" / "bounding epsilon" warning appears.
the following is an output from the last steps of simpleFoam:

I am only interested in the results of p and U, and not interested in epsilon and k.
simpleFoam is stable and converging.

my question is:
-is simpleFoam converging to a correct solution (for p and U)?
-what is the accuracy of p and U results, does the bounding k and epsilon problem influence the results too much?

thank you