Allow turbulentDispersion model to use kinetic energy either from object registry or calculated from turbulence models.

This commit is contained in:
Thomas Lichtenegger
2020-11-13 10:13:37 +01:00
parent d3c10ba24d
commit 1a37e99dfd
2 changed files with 45 additions and 9 deletions

View File

@ -74,13 +74,16 @@ turbulentDispersion::turbulentDispersion
sm.mesh(), sm.mesh(),
dimensionedScalar("zero", dimensionSet(0,0,0,0,0,0,0), 0.0) dimensionedScalar("zero", dimensionSet(0,0,0,0,0,0,0), 0.0)
), ),
minTurbKinetcEnergy_(propsDict_.lookupOrDefault<scalar>("minTurbKinetcEnergy", 0.0)), minTurbKineticEnergy_(propsDict_.lookupOrDefault<scalar>("minTurbKineticEnergy", 0.0)),
turbKineticEnergyFieldName_(propsDict_.lookupOrDefault<word>("turbKineticEnergyFieldName","")),
turbKineticEnergy_(NULL),
existTurbKineticEnergyInObjReg_(false),
voidfractionFieldName_(propsDict_.lookupOrDefault<word>("voidfractionFieldName","voidfraction")), voidfractionFieldName_(propsDict_.lookupOrDefault<word>("voidfractionFieldName","voidfraction")),
voidfraction_(sm.mesh().lookupObject<volScalarField> (voidfractionFieldName_)), voidfraction_(sm.mesh().lookupObject<volScalarField> (voidfractionFieldName_)),
critVoidfraction_(propsDict_.lookupOrDefault<scalar>("critVoidfraction", 0.9)), critVoidfraction_(propsDict_.lookupOrDefault<scalar>("critVoidfraction", 0.9)),
ranGen_(clock::getTime()+pid()) ranGen_(clock::getTime()+pid())
{ {
if(ignoreCellsName_ != "none") if (ignoreCellsName_ != "none")
{ {
ignoreCells_.set(new cellSet(particleCloud_.mesh(),ignoreCellsName_)); ignoreCells_.set(new cellSet(particleCloud_.mesh(),ignoreCellsName_));
Info << type << ": ignoring fluctuations in cellSet " << ignoreCells_().name() << Info << type << ": ignoring fluctuations in cellSet " << ignoreCells_().name() <<
@ -90,7 +93,7 @@ turbulentDispersion::turbulentDispersion
// define a field to indicate if a cell is next to boundary // define a field to indicate if a cell is next to boundary
label cellI = -1; label cellI = -1;
forAll(mesh_.boundary(),patchI) forAll (mesh_.boundary(),patchI)
{ {
word patchName = mesh_.boundary()[patchI].name(); word patchName = mesh_.boundary()[patchI].name();
if (patchName.rfind("procB",0) == 0) continue; if (patchName.rfind("procB",0) == 0) continue;
@ -102,6 +105,29 @@ turbulentDispersion::turbulentDispersion
} }
} }
if (turbKineticEnergyFieldName_ != "")
{
existTurbKineticEnergyInObjReg_ = true;
volScalarField& k(const_cast<volScalarField&>(sm.mesh().lookupObject<volScalarField> (turbKineticEnergyFieldName_)));
turbKineticEnergy_ = &k;
}
else
{
turbKineticEnergy_ = new volScalarField
(
IOobject
(
"turbKinEnergy",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh_,
dimensionedScalar("zero", dimensionSet(0,2,-2,0,0), 0)
);
}
// make sure this is the last force model in list so that fluid velocity does not get overwritten // make sure this is the last force model in list so that fluid velocity does not get overwritten
label numLastForceModel = sm.nrForceModels(); label numLastForceModel = sm.nrForceModels();
word lastForceModel = sm.forceModels()[numLastForceModel-1]; word lastForceModel = sm.forceModels()[numLastForceModel-1];
@ -116,6 +142,7 @@ turbulentDispersion::turbulentDispersion
turbulentDispersion::~turbulentDispersion() turbulentDispersion::~turbulentDispersion()
{ {
if (!existTurbKineticEnergyInObjReg_) delete turbKineticEnergy_;
} }
// * * * * * * * * * * * * * * * private Member Functions * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * private Member Functions * * * * * * * * * * * * * //
@ -130,7 +157,10 @@ bool turbulentDispersion::ignoreCell(label cell) const
void turbulentDispersion::setForce() const void turbulentDispersion::setForce() const
{ {
const volScalarField turbKinetcEnergy(particleCloud_.turbulence().k()); if (!existTurbKineticEnergyInObjReg_)
{
*turbKineticEnergy_ = particleCloud_.turbulence().k()();
}
label cellI = -1; label cellI = -1;
label patchID = -1; label patchID = -1;
@ -142,7 +172,7 @@ void turbulentDispersion::setForce() const
vector position = vector::zero; vector position = vector::zero;
word patchName(""); word patchName("");
interpolationCellPoint<scalar> turbKinetcEnergyInterpolator_(turbKinetcEnergy); interpolationCellPoint<scalar> turbKineticEnergyInterpolator_(*turbKineticEnergy_);
for(int index = 0;index < particleCloud_.numberOfParticles(); ++index) for(int index = 0;index < particleCloud_.numberOfParticles(); ++index)
{ {
@ -156,14 +186,14 @@ void turbulentDispersion::setForce() const
if (interpolate_) if (interpolate_)
{ {
position = particleCloud_.position(index); position = particleCloud_.position(index);
k = turbKinetcEnergyInterpolator_.interpolate(position,cellI); k = turbKineticEnergyInterpolator_.interpolate(position,cellI);
} }
else else
{ {
k = turbKinetcEnergy[cellI]; k = (*turbKineticEnergy_)[cellI];
} }
if (k < minTurbKinetcEnergy_) k = minTurbKinetcEnergy_; if (k < minTurbKineticEnergy_) k = minTurbKineticEnergy_;
flucU=unitFlucDir()*Foam::sqrt(2.0*k); flucU=unitFlucDir()*Foam::sqrt(2.0*k);

View File

@ -60,7 +60,13 @@ protected:
mutable volScalarField wallIndicatorField_; mutable volScalarField wallIndicatorField_;
scalar minTurbKinetcEnergy_; scalar minTurbKineticEnergy_;
word turbKineticEnergyFieldName_;
volScalarField* turbKineticEnergy_;
bool existTurbKineticEnergyInObjReg_;
word voidfractionFieldName_; word voidfractionFieldName_;