Adding fvOptions to turbulence model
 

Does anyone know if it's possible to get fvOptions sources to work in the equations for turbulence quantities like k and epsilon? I've been trying to use a scalarExplicitSetValue 'source' for both epsilon and k using the standard kEpsilon model, with the subdict:

turb1
{
type scalarExplicitSetValue;
active true;
selectionMode cellZone;
cellZone Zone_1;

scalarExplicitSetValueCoeffs
{
injectionRate
{
epsilon 1;
k 1;
}
}
}

but what I get when I run is:

--> FOAM Warning :
From function void option::checkApplied() const
in file fvOptions/fvOption.C at line 368
Source turb1 defined for field epsilon but never used
--> FOAM Warning :
From function void option::checkApplied() const
in file fvOptions/fvOption.C at line 368
Source turb1 defined for field k but never used

I've had a look in the code for the kEpsilon model, and tried to recompile with the equations as follows (adding fvOptions source at the end):

// Dissipation equation
tmp<fvScalarMatrix> epsEqn
(
fvm::ddt(epsilon_)
+ fvm::div(phi_, epsilon_)
- fvm::laplacian(DepsilonEff(), epsilon_)
==
C1_*G*epsilon_/k_
- fvm::Sp(C2_*epsilon_/k_, epsilon_)
+ fvOptions(epsilon_)
);

and:

// Turbulent kinetic energy equation
tmp<fvScalarMatrix> kEqn
(
fvm::ddt(k_)
+ fvm::div(phi_, k_)
- fvm::laplacian(DkEff(), k_)
==
G
- fvm::Sp(epsilon_/k_, k_)
+ fvOptions(k_)
);

but I get the compilation error:

mykEpsilon.C: In member function ‘virtual void Foam::incompressible::RASModels::mykEpsilon::corre ct()’:
mykEpsilon.C:254:27: error: ‘fvOptions’ was not declared in this scope

I'm afraid that my very limited knowledge ends roughly here. The problem is especially irksome as the sample in the code for explicitSetValue.H is:

32 \verbatim
33 <Type>ExplicitSetValueCoeffs
34 {
35 injectionRate
36 {
37 k 30.7;
38 epsilon 1.5;
39 }
40 }
41 \endverbatim

Did this work in a previous OF version? I'm using 2.3.0.

Hi,

if you'd like to use fvOptions you have to create it. You add a property to class

and initialize it in constructor with something like

in constructor of your turbulence model.

Hi Alexey,

Thanks so much for getting back to me so quickly!

I added the lines you recommended, and with an #include "fvIOoptionList.H" in the mykEpsilon.H file the modified turbulence model compiled correctly and I was able to call it in the solver. Unfortunately it still doesn't seem to be working correctly.

I'm running a version of pimpleFoam that I've modified slightly to include the transport of a passive scalar, and I have fvOptions sources for both momentum and the scalar, vectorCodedSources applied to set time-varying values of U and scalarSemiImplicitSources to set zonal injection rates of the scalar. These both have been working as expected. To recap my objective, I'm looking to prescribe zonal values of k and epsilon with scalarExplicitSetValue sources. When I run pimpleFoam now, I get:

--> FOAM Warning :
From function void option::checkApplied() const
in file fvOptions/fvOption.C at line 368
Source turb1 defined for field epsilon but never used
--> FOAM Warning :
From function void option::checkApplied() const
in file fvOptions/fvOption.C at line 368
Source turb1 defined for field k but never used

before the Ux, Uy, Uz and scalar solution steps, and:

--> FOAM Warning :
From function void option::checkApplied() const
in file fvOptions/fvOption.C at line 368
Source codedSourceMom_1 defined for field U but never used
--> FOAM Warning :
From function void option::checkApplied() const
in file fvOptions/fvOption.C at line 368
Source smoke1 defined for field smoke but never used

before the k, epsilon solution steps, which I figure is all okay, as I wouldn't expect the turbulent sources to be used in the momentum/scalar equations and vice versa. The problem I'm having is that although the velocity and scalar are still being set correctly, the turbulence quantities are unchanged, whereas I was hoping for my fixed, specified levels (currently unity) to be set in those zones. Any idea what I'm missing?

Thanks again for your help,

Robin

Hi,

Can you post output of the simulation for one time step? I'd like to see the errors you've posted here in the context of solver output.

Hi Alexey,

Here is the output of the first two iterations. Seemingly I only get the warnings on the second time step?

Robin

Hi,

In case of fvOption in turbulence model, it is really applied only on the last PIMPLE iteration (as it is instructed by turbOnFinalIterOnly option in PIMPLE dictionary in fvSolution).

In case of codedSource fvOption, I guess, your code is in charge of setting applied state. As I don't know the source of Mom_x options there's nothing I can add.

Hi Alexey,

The momentum sources are currently just applying fixed values of velocity over the zones, although I intend to make them time/spatially varying in due course.

I think I've fixed the problem! I hadn't applied the functions fvOptions.constrain() and fvOptions.correct() to the turbulence quantities. Having done that, the relevant sections of mykEpsilon.H are something like:

tmp<fvScalarMatrix> epsEqn
(
fvm::ddt(epsilon_)
+ fvm::div(phi_, epsilon_)
- fvm::laplacian(DepsilonEff(), epsilon_)
==
C1_*G*epsilon_/k_
- fvm::Sp(C2_*epsilon_/k_, epsilon_)
+ fvOptions(epsilon_)
);

epsEqn().relax();

fvOptions.constrain(epsEqn());

epsEqn().boundaryManipulate(epsilon_.boundaryField ());

solve(epsEqn);

fvOptions.correct(epsilon_);

bound(epsilon_, epsilonMin_);

and the same for the k equation. Seems to be working now as the turbulence quantities are constant within the prescribed zones.

Thanks very much for your help on this, I don't think I would have been able to get started without your help.

Robin

Hi,

good to know you've solved the problem. Actually it was my first suggestion but then I've decided to read source code of fvOption class and found that operator () sets applied state. So I've decided that constrain and correct calls are in fact not necessary. It seems I was wrong ;)

slottedLemon,

I am looking to add sources to kEpsilon as well, but I can't get the code to compile. Would you be willing to share your code?

What I have so far is up at https://github.com/petebachant/kEpsilonSources

I got similar issue trying to make fvOptions working with k-e model by adding a source term 'scalarCodedSource'.

Did you guys find a solution about your issues?

here is the link of the thread:

https://www.cfd-online.com/Forums/op...k-e-model.html

Hi I just manged to get fvOptions codedScalarSource to work for my case, I did not use it for turbulence, but for a heat source. Maybe you can take some of it I don't know.



https://www.cfd-online.com/Forums/op...at-source.html

Thank you very much for your reply!

Is there a reason why your heSource[i] function is multiply by the volume V[i]?

Some day a go, I suspected that I had to multiply my function by V[i](I assumed that I had to do this since I didn't use the 'fvm' for the coding), and the simulation went well and converge but I was wondering if it was a simple luck.

Thank in advance for your feedback.

Regards,

Hey, I did not bother about that^^. But I mean it depends. It does not make really sense in that case to have the source be dependent on the cell volume, since he is a specific quantity. It would mean, that in smaller cell the specific source would be smaller what would not make sense. It was more of an example for me how to include different attributes of the simulation into fvOptions. But I could be totally wrong, because in the official example:


https://www.openfoam.com/documentati...ces-coded.html


they also included V[i], for a specific heat source. I have to think about that, and look for what the EEqn really solves.

I understand what do you mean. Since it a specific value, doesn't make sense to multiply by the volume.

The thing is, the term source that I try to implement is not a heat source but a source term apply to the dissipation equation (epsilon).

Like you saying, it depend on what OpenFoam do during the resolving process with the function in 'fvOptions' file. And I noticed that the V[i] was included in other example also, like you said.

For instance, we have to divide the term by 'rho', because we use 'nut' in the transport equation and I guess it multiply by 'rho' later on.

So, I was suspecting that my issue could be something similar to this.

P.S.: If you look at the thread I linked, the results calculated by the function in the fvOptions file is right(I validated the calculation with MathCad) but when it carry on in the transport equation to resolve it, the digit go really high.