LESRegion computed wrongly when using kOmegaSSTIDDES turbulence model
Functionality to add/problem to solve
The LESRegion()
function in the kOmegaSSTIDDES
turbulence model is not implemented.
Instead the implementation from the base class kOmegaSSTDES
is used:
template<class BasicTurbulenceModel>
tmp<volScalarField> kOmegaSSTDES<BasicTurbulenceModel>::LESRegion() const
{
const volScalarField& k = this->k_;
const volScalarField& omega = this->omega_;
const volVectorField& U = this->U_;
const volScalarField CDkOmega
(
(2*this->alphaOmega2_)*(fvc::grad(k) & fvc::grad(omega))/omega
);
const volScalarField F1(this->F1(CDkOmega));
return tmp<volScalarField>::New
(
IOobject::scopedName("DES", "LESRegion"),
neg(dTilda(mag(fvc::grad(U)), CDES(F1)) - lengthScaleRAS())
);
}
leading to LES region indicator computed with different equation than the one actually used in the model:
template<class BasicTurbulenceModel>
tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::dTilda
(
const volScalarField& magGradU,
const volScalarField& CDES
) const
{
const volScalarField lRAS(this->lengthScaleRAS());
const volScalarField lLES(this->lengthScaleLES(CDES));
const volScalarField alpha(this->alpha());
const volScalarField expTerm(exp(sqr(alpha)));
tmp<volScalarField> fB = min(2*pow(expTerm, -9.0), scalar(1));
const volScalarField fdTilda(max(1 - fdt(magGradU), fB));
if (fe_)
{
/*
...
*/
}
// Simplified formulation from Gritskevich et al. paper (2011) where fe = 0
return max
(
fdTilda*lRAS + (1 - fdTilda)*lLES,
dimensionedScalar(dimLength, SMALL)
);
}
While writing this issue, I became aware that outputting fd
is also possible now, as of v2212.
It seems that for IDDES model it is actually equivalent to the way fdTilda
is computed in the code,
and 1-fdTilda
multiplies the LES length scale in both versions of blending, so it serves the same purpose (as far as I understand it).
Perhaps one option to resolve this would be implementing a virtual LESRegion()
that throws FatalError
telling the user to
output fd
instead and subtract it from 1?
Target audience
Target audience are the users of hybrid RANS/LES turbulence models. Proper indication in which region the model works in a LES mode is important for debugging and setting up cases with those models properly.
Proposal
I have made a quick solution for myself in an out-of-repository copy of kOmegaSSTIDDES.
- Added a protected member function, to refactor computation of
fdTilda
. Header:
//- Return the fdTilda function (moved to a function for
// LESRegions function)
tmp<volScalarField> fdTilda
(
const volScalarField &magGradU,
const volScalarField &expTerm
) const;
Implementation:
template <class BasicTurbulenceModel>
tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::fdTilda
(
const volScalarField &magGradU,
const volScalarField &expTerm
) const
{
tmp<volScalarField> fB = min(2 * pow(expTerm, -9.0), scalar(1));
return max(1 - fdt(magGradU), fB);
}
- Added the
LESRegion()
: Header:
//- Return the LES field indicator
virtual tmp<volScalarField> LESRegion() const;
Implementation:
template <class BasicTurbulenceModel>
tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::LESRegion() const
{
const volScalarField magGradU(mag(fvc::grad(this->U_)));
const volScalarField expTerm(exp(sqr(this->alpha())));
tmp<volScalarField> tLESRegion
(
new volScalarField
(
IOobject::scopedName("IDDES", "LESRegion"),
1 - this->fdTilda(magGradU, expTerm)
)
);
// Solver is in RANS mode next to the wall.
// i have no clue if this can cause problems
tLESRegion.ref().boundaryFieldRef() == 0.0;
return tLESRegion;
}
What does success look like, and how can we measure that?
The LES region computed with unchanged code from a sample test case:
The LES region computed with the code above, using fdTilda
as it is computed in dTilda
function:
Links / references
Reference from KOmegaSSTIDDES model: Gritskevich, M.S., Garbaruk, A.V., Schuetze, J., Menter, F.R. (2011) Development of DDES and IDDES Formulations for the k-omega Shear Stress Transport Model, Flow, Turbulence and Combustion, pp. 1-19