Go to the documentation of this file.
32 template<
class Thermo,
class OtherThermo>
45 IOobject::groupName(
"gamma1", pair.
name()),
46 pair.
phase1().mesh().time().timeName(),
56 IOobject::groupName(
"gamma2", pair.
name()),
57 pair.
phase1().mesh().time().timeName(),
68 if (this->speciesNames_.size() != 2)
71 <<
"NonRandomTwoLiquid model is suitable for two species only."
75 species1Name_ = this->speciesNames_[0];
76 species2Name_ = this->speciesNames_[1];
78 species1Index_ = this->thermo_.composition().species()[species1Name_];
79 species2Index_ = this->thermo_.composition().species()[species2Name_];
81 alpha12_.read(
"alpha",
dict.subDict(species1Name_));
82 alpha21_.read(
"alpha",
dict.subDict(species2Name_));
83 beta12_.read(
"beta",
dict.subDict(species1Name_));
84 beta21_.read(
"beta",
dict.subDict(species2Name_));
86 saturationModel12_.reset
90 dict.subDict(species1Name_).subDict(
"interaction"),
94 saturationModel21_.reset
98 dict.subDict(species2Name_).subDict(
"interaction"),
107 dict.subDict(species1Name_),
116 dict.subDict(species2Name_),
125 template<
class Thermo,
class OtherThermo>
133 template<
class Thermo,
class OtherThermo>
145 this->thermo_.composition().Y(species1Index_)
151 this->thermo_.composition().W(species1Index_)
157 this->thermo_.composition().Y(species2Index_)
163 this->thermo_.composition().W(species2Index_)
181 tau21*
sqr(G21)/
max(
sqr(X1 + X2*G21), SMALL)
182 + tau12*G12/
max(
sqr(X2 + X1*G12), SMALL)
190 tau12*
sqr(G12)/
max(
sqr(X2 + X1*G12), SMALL)
191 + tau21*G21/
max(
sqr(X1 + X2*G21), SMALL)
197 template<
class Thermo,
class OtherThermo>
201 const word& speciesName,
205 if (speciesName == species1Name_)
208 this->otherThermo_.composition().Y(speciesName)
209 *speciesModel1_->Yf(speciesName, Tf)
212 else if (speciesName == species2Name_)
215 this->otherThermo_.composition().Y(speciesName)
216 *speciesModel2_->Yf(speciesName, Tf)
222 this->thermo_.composition().Y(speciesName)
223 *(scalar(1) - Yf(species1Name_, Tf) - Yf(species2Name_, Tf));
228 template<
class Thermo,
class OtherThermo>
233 const word& speciesName,
237 if (speciesName == species1Name_)
240 this->otherThermo_.composition().Y(speciesName)
241 *speciesModel1_->YfPrime(speciesName, Tf)
244 else if (speciesName == species2Name_)
247 this->otherThermo_.composition().Y(speciesName)
248 *speciesModel2_->YfPrime(speciesName, Tf)
254 - this->thermo_.composition().Y(speciesName)
255 *(YfPrime(species1Name_, Tf) + YfPrime(species2Name_, Tf));
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Base class for interface composition models, templated on the two thermodynamic models either side of...
NonRandomTwoLiquid(const dictionary &dict, const phasePair &pair)
Construct from components.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
dimensionedScalar exp(const dimensionedScalar &ds)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void update(const volScalarField &Tf)
Update the composition.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual word name() const
Pair name.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
const phaseModel & phase1() const
virtual ~NonRandomTwoLiquid()
Destructor.