[matejmuller]

Hello!

I'm having some trouble with limiting the maximal velocities inside the whole domain. I've tried with a loop, something like in this thread http://www.cfd-online.com/Forums/ope...terfoam-2.html

When i try

    
forAll(U,celli)
{
    if (U[celli].component(vector::X) > 1){
        U[celli].component(vector::X) = 1;}
 
    if (U[celli].component(vector::Y) > 1){
        U[celli].component(vector::Y) = 1;}
 
    if (U[celli].component(vector::Z) > 1){
        U[celli].component(vector::Z) = 1;}
}

and call it right after the Ueqn.h, it compiles without errors, but the velocities from the results are still greater than 1. I'm pretty new to openfoam code...I've read that one should do this without a loop... has anyone any idea how to do that?

Thanks in advance!

Matej

[kmooney]

You probably want to double check that you're performing the limit after the last time U is updated. With piso this is often at the very end of the timestep loop.

See line #128:
https://github.com/OpenFOAM/OpenFOAM...oam/pisoFoam.C

[matejmuller]

Hi!

Thank you for your answer. That was actually the problem...

However, I also want to bound the max k and epsilon values when using a k-epsilon model and k and omega values when using a k-omega model. But when I use

forAll(k, celli)

in interfoam.c, at compiling I get the message:

error: ‘k’ was not declared in this scope

Do I have to compile a new turbulence model? Perhaps there is a simpler way to do this? I would imagine openfoam has this already implemented but where can i define the max values?

Regards, Matej

[kmooney]

Originally Posted by matejmuller

Hi!

Thank you for your answer. That was actually the problem...

However, I also want to bound the max k and epsilon values when using a k-epsilon model and k and omega values when using a k-omega model. But when I use

forAll(k, celli)

in interfoam.c, at compiling I get the message:

error: ‘k’ was not declared in this scope

Do I have to compile a new turbulence model? Perhaps there is a simpler way to do this? I would imagine openfoam has this already implemented but where can i define the max values?

Regards, Matej

Hi Matej,

You just need to lookup the field of interest before trying to operate on it at the top level of the solver.

This should explain how:
https://openfoamwiki.net/index.php/Contrib/IOReferencer

Cheers!
Kyle

[matejmuller]

Hi Kyle!

So I've looked in the code kEpsilon.H, and k nad epsilon are IObjects. But when compiling with:

const volScalarField& k = db().lookupObject<volScalarField>("k");

I get:

error: ‘db’ was not declared in this scope

What does " // The "db()" above may vary" in your provided link mean?

Thanks in advance!

Matej

[kmooney]

Instead of

db()

try:

mesh

Which should work in the top level solver (for most solvers).

[matejmuller]

Hi!

Thanks for the response. The mesh part did the trick. However, I couldn't access the k and epsilon values, because they are somehow write-protected. I have found it easier to compile a new turbulence model according to http://www.tfd.chalmers.se/~hani/kur...lenceModel.pdf and manipulated the k and epsilon values in it.

Regards, Matej

[vitor.geraldes@ist.utl.pt]

I managed to ensure that the velocity can not increase beyond a given threshold value by introducing an artificial hydrodynamic resistance that increases polynomially as the velocity magnitude approaches that limit.
The code is this one

forAll(Rdamp,cellI)
{
Rdamp[cellI] = Rref.value()*Foam:ow((mag(U[cellI])/maxVelocity.value()),10);
}
Rdamp.correctBoundaryConditions();
Rdamp.relax();

fvVectorMatrix UEqn
(
fvm::ddt(rho, U)
+ fvm::div(rhoPhi, U)
- fvm::laplacian(mu, U)
- fvc::div(mu*fvc::grad(U)().T())
+ fvm::Sp(Rdamp, U)
);

In my case, the units of Rdamp are [1 -3 -1 0 0 0 0].

With this technique, the pressure correction equation remains conservative. I imposed an high power of 10 to ensure that the dampening resistance is negligible when the velocity is slightly lower than the threshold value. Underelaxation of Rdamp on the order of 0.5 - 0.8 appears to have good effects on the numerical stability of the algorithm.

I would be glad to know if this numerical trick works for other cases too.

[akidess]

For completeness, OpenFOAM now includes fvOptions to do this: limitVelocity (dev), velocityDampingConstraint (plus).