40template<
class BasePhaseSystem>
78template<
class BasePhaseSystem>
103 if (massTransferModels_.found(keyik))
106 massTransferModels_[keyik];
108 word speciesName = interfacePtr->transferSpecie();
110 const word species(speciesName.substr(0, speciesName.find(
'.')));
112 L -=
neg(dmdtNetki)*interfacePtr->L(species,
T);
115 if (massTransferModels_.found(keyki))
118 massTransferModels_[keyki];
120 word speciesName = interfacePtr->transferSpecie();
122 const word species(speciesName.substr(0, speciesName.find(
'.')));
124 L -=
pos(dmdtNetki)*interfacePtr->L(species,
T);
133template<
class BasePhaseSystem>
151 auto& dmdt = tdmdt.ref();
153 if (dmdt_.found(key))
162template<
class BasePhaseSystem>
170 auto& eqn = teqn.ref();
178 for (++iterk; iterk != this->phaseModels_.end(); ++iterk)
180 if (iteri()().
name() != iterk()().
name())
202 auto& dmdtNetki = tdmdtNetki.ref();
215 auto&
Sp = tSp.ref();
228 auto&
Su = tSu.ref();
231 if (massTransferModels_.found(keyik))
234 massTransferModels_[keyik];
236 dmdtNetki -= *dmdt_[keyik];
262 if (massTransferModels_.found(keyki))
265 massTransferModels_[keyki];
267 dmdtNetki += *dmdt_[keyki];
304template<
class BasePhaseSystem>
312 auto& eqn = teqn.ref();
325 auto&
Sp = tSp.ref();
338 auto&
Su = tSu.ref();
354 if (massTransferModels_.found(key12))
357 massTransferModels_[key12];
404 if (massTransferModels_.found(key21))
407 massTransferModels_[key21];
454template<
class BasePhaseSystem>
466 for (++iterk; iterk != this->phaseModels_.end(); ++iterk)
468 if (iteri()().
name() != iterk()().
name())
478 if (massTransferModels_.found(keyik))
481 massTransferModels_[keyik];
485 *dmdt_[keyik] = Kexp.
ref();
489 if (massTransferModels_.found(keyki))
492 massTransferModels_[keyki];
497 *dmdt_[keyki] = Kexp.
ref();
505template<
class BasePhaseSystem>
514 bool includeDivU(
true);
543 if (massTransferModels_.found(key12))
546 massTransferModels_[key12];
556 includeDivU = interfacePtr->includeDivU();
571 if (massTransferModels_.found(key21))
574 massTransferModels_[key21];
584 includeDivU = interfacePtr->includeDivU();
617 scalar dmdt21 = dmdtNet[celli];
618 scalar coeffs12Cell = coeffs12[celli];
620 scalar alpha1Limited =
max(
min(
alpha1[celli], 1.0), 0.0);
623 SuPhase1[celli] += coeffs1[celli]*dmdt21;
629 if (coeffs12Cell > 0)
632 SpPhase1[celli] -= dmdt21*coeffs12Cell;
634 else if (coeffs12Cell < 0)
638 dmdt21*coeffs12Cell*alpha1Limited;
643 if (coeffs12Cell > 0)
647 dmdt21*coeffs12Cell*alpha1Limited;
649 else if (coeffs12Cell < 0)
652 SpPhase1[celli] -= dmdt21*coeffs12Cell;
660 scalar dmdt12 = -dmdtNet[celli];
661 scalar coeffs21Cell = -coeffs12[celli];
663 scalar alpha2Limited =
max(
min(
alpha2[celli], 1.0), 0.0);
666 SuPhase2[celli] += coeffs2[celli]*dmdt12;
672 if (coeffs21Cell > 0)
675 SpPhase2[celli] -= dmdt12*coeffs21Cell;
677 else if (coeffs21Cell < 0)
681 dmdt12*coeffs21Cell*alpha2Limited;
686 if (coeffs21Cell > 0)
690 coeffs21Cell*dmdt12*alpha2Limited;
692 else if (coeffs21Cell < 0)
695 SpPhase2[celli] -= dmdt12*coeffs21Cell;
704 max(
gMax((dmdt21*coeffs1)()),
gMax((dmdt12*coeffs2)()));
709template<
class BasePhaseSystem>
715 const word speciesName
721 if (iter()->transferSpecie() == speciesName)
732template<
class BasePhaseSystem>
735 bool includeVolChange(
true);
738 if (!iter()->includeVolChange())
740 includeVolChange =
false;
743 return includeVolChange;
const volScalarField & alpha1
const volScalarField & alpha2
tmp< volScalarField > dmdt() const
Return the blended mass transfer rate.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
bool set(const Key &key, T *ptr)
Assign a new entry, overwriting existing entries.
bool found(const Key &key) const
Return true if hashed entry is found in table.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Class for mass transfer between phases.
virtual void correctMassSources(const volScalarField &T)
Correct/calculates mass sources dmdt for phases.
virtual void massSpeciesTransfer(const Foam::phaseModel &phase, volScalarField::Internal &Su, volScalarField::Internal &Sp, const word speciesName)
Calculate mass transfer for species.
dmdtTable dmdt_
Overall inter-phase mass transfer rates [Kg/s].
virtual tmp< fvScalarMatrix > volTransfer(const volScalarField &p)
Return the volumetric rate transfer matrix.
virtual bool includeVolChange()
Add volume change in pEq.
massTransferModelTable massTransferModels_
Mass transfer models.
virtual void alphaTransfer(SuSpTable &Su, SuSpTable &Sp)
Calculate mass transfer for alpha's.
tmp< volScalarField > calculateL(const volScalarField &dmdtNetki, const phasePairKey &keyik, const phasePairKey &keyki, const volScalarField &T) const
Calculate L between phases.
const HeatTransferModel< ThermoCloud< CloudType > > & heatTransfer() const
Return reference to heat transfer model.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Mesh data needed to do the Finite Volume discretisation.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
const word & name() const
An ordered or unorder pair of phase names. Typically specified as follows.
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
const multiphaseInter::phaseModel & phase1() const
const multiphaseInter::phaseModel & phase2() const
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
A class for managing temporary objects.
bool valid() const noexcept
Identical to good(), or bool operator.
A class for handling words, derived from Foam::string.
Fundamental dimensioned constants.
Calculate the divergence of the given field.
Calculate the finiteVolume matrix for implicit and explicit sources.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
const dimensionSet dimPressure
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
const dimensionSet dimless
Dimensionless.
const dimensionSet dimEnergy
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
dimensionedScalar neg(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVolume(pow3(dimLength))
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const dimensionSet dimDensity
Type gMax(const FieldField< Field, Type > &f)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
#define forAll(list, i)
Loop across all elements in list.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
const vector L(dict.get< vector >("L"))