Go to the documentation of this file.
30 #include "interfaceCompositionModel.H"
36 template<
class BasePhaseSystem>
47 this->
phases()[key.first()],
57 if (interfaceCompositionModels_.found(pair))
64 interfaceCompositionModels_[pair]->species(),
68 const word& member = *memberIter;
73 IOobject::groupName(member, otherPhase.
name())
79 *(*iDmdtSu_[pair])[member]
80 + *(*iDmdtSp_[pair])[member]*
phase.Y(member)
93 template<
class BasePhaseSystem>
100 BasePhaseSystem(
mesh),
101 nInterfaceCorrectors_
103 this->
template getOrDefault<label>(
"nInterfaceCorrectors", 1)
106 this->generatePairsAndSubModels
108 "interfaceComposition",
109 interfaceCompositionModels_
112 this->generatePairsAndSubModels
123 interfaceCompositionModels_,
124 interfaceCompositionModelIter
128 this->phasePairs_[interfaceCompositionModelIter.key()];
135 <<
"An interfacial composition model is specified for the "
136 <<
"unordered " << pair <<
" pair. Composition models only "
137 <<
"apply to ordered pairs. An entry for a "
138 <<
phasePairKey(
"A",
"B",
true) <<
" pair means a model for "
139 <<
"the A side of the A-B interface; i.e., \"A in the presence "
147 if (!this->phasePairs_.found(key))
150 <<
"A mass transfer model the " << key <<
" pair is not "
151 <<
"specified. This is required by the corresponding interface "
152 <<
"composition model."
156 const phasePair& uoPair = this->phasePairs_[key];
158 if (!massTransferModels_[uoPair][uoPair.
index(
phase)].valid())
161 <<
"A mass transfer model for the " << pair.
phase1().
name()
162 <<
" side of the " << uoPair <<
" pair is not "
163 <<
"specified. This is required by the corresponding interface "
164 <<
"composition model."
172 massTransferModelIter
176 this->phasePairs_[massTransferModelIter.key()];
178 if (!this->heatTransferModels_.found(pair))
181 <<
"A heat transfer model for " << pair <<
" pair is not "
182 <<
"specified. This is required by the corresponding species "
192 interfaceCompositionModels_,
193 interfaceCompositionModelIter
197 interfaceCompositionModelIter();
200 this->phasePairs_[interfaceCompositionModelIter.key()];
207 const word& member = *memberIter;
216 IOobject::groupName(
"iDmdtSu", pair.
name()),
217 this->
mesh().time().timeName(),
232 IOobject::groupName(
"iDmdtSp", pair.
name()),
233 this->
mesh().time().timeName(),
247 template<
class BasePhaseSystem>
255 template<
class BasePhaseSystem>
262 return BasePhaseSystem::dmdt(key) + this->iDmdt(key);
266 template<
class BasePhaseSystem>
275 interfaceCompositionModels_,
276 interfaceCompositionModelIter
280 interfaceCompositionModelIter();
283 this->phasePairs_[interfaceCompositionModelIter.key()];
289 const word& member = *memberIter;
294 IOobject::groupName(member, otherPhase.
name())
299 *(*iDmdtSu_[pair])[member]
300 + *(*iDmdtSp_[pair])[member]*
phase.Y(member)
303 this->addField(
phase,
"dmdt", iDmdt, dmdts);
304 this->addField(otherPhase,
"dmdt", - iDmdt, dmdts);
312 template<
class BasePhaseSystem>
326 interfaceCompositionModels_,
327 interfaceCompositionModelIter
331 interfaceCompositionModelIter();
334 this->phasePairs_[interfaceCompositionModelIter.key()];
344 massTransferModels_[unorderedPair][unorderedPair.index(
phase)]->K()
349 const word& member = *memberIter;
354 IOobject::groupName(member, otherPhase.
name())
361 *(*iDmdtSu_[pair])[member] =
phase.
rho()*KD*Yf;
362 *(*iDmdtSp_[pair])[member] = -
phase.
rho()*KD;
366 *(*iDmdtSu_[pair])[member]
372 *(*iDmdtSu_[pair])[member]
373 + *(*iDmdtSp_[pair])[member]*
phase.Y(member)
380 if (eqns.found(otherName))
382 *eqns[otherName] -= iDmdt;
391 template<
class BasePhaseSystem>
409 typename BasePhaseSystem::heatTransferModelTable,
410 this->heatTransferModels_,
411 heatTransferModelIter
415 this->phasePairs_[heatTransferModelIter.key()];
426 for (label i = 0; i < nInterfaceCorrectors_; ++ i)
452 if (this->interfaceCompositionModels_.found(key12))
454 this->interfaceCompositionModels_[key12]->addMDotL
456 massTransferModels_[pair].first()->
K(),
462 if (this->interfaceCompositionModels_.found(key21))
464 this->interfaceCompositionModels_[key21]->addMDotL
466 massTransferModels_[pair].second()->
K(),
482 max(H1 + H2 + mDotLPrime, HSmall)
494 if (this->interfaceCompositionModels_.found(key12))
496 this->interfaceCompositionModels_[key12]->update(Tf);
498 if (this->interfaceCompositionModels_.found(key21))
500 this->interfaceCompositionModels_[key21]->update(Tf);
507 template<
class BasePhaseSystem>
const T & second() const noexcept
Return second element, which is also the last element.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
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...
virtual bool read()
Read base phaseProperties dictionary.
A class for handling words, derived from Foam::string.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
A class for managing temporary objects.
const dimensionSet dimEnergy
const dimensionSet dimDensity
virtual tmp< volScalarField > D(const word &speciesName) const =0
Mass diffusivity.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
const phaseModel & otherPhase() const
Return the other phase in this two-phase system.
label index(const phaseModel &phase) const
Return the index of the given phase. Generates a FatalError if.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual const rhoThermo & thermo() const =0
Return the thermophysical model.
virtual void correctInterfaceThermo()
Correct the interface temperatures.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object.
Generic base class for interface composition models. These models describe the composition in phase 1...
const hashedWordList & species() const
Return the transferring species names.
CGAL::Exact_predicates_exact_constructions_kernel K
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< GeometricField< Type, faPatchField, areaMesh > > average(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Area-weighted average a edgeField creating a areaField.
A wordList with hashed named lookup, which can be faster in some situations than using the normal lis...
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the mass transfer rate for a pair.
Mesh data needed to do the Finite Volume discretisation.
bool ordered() const
Return the ordered flag.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const word & name() const
InterfaceCompositionPhaseChangePhaseSystem(const fvMesh &)
Construct from fvMesh.
void correctBoundaryConditions()
Correct boundary field.
const word & name() const
Return the name of this phase.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
An ordered pair of two objects of type <T> with first() and second() elements.
virtual const volScalarField & T() const
Temperature [K].
A HashTable of pointers to objects of type <T>.
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const =0
Interface mass fraction.
multiphaseSystem::phaseModelList & phases
virtual word name() const
Pair name.
virtual ~InterfaceCompositionPhaseChangePhaseSystem()
Destructor.
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
const phaseModel & phase2() const
Return phase 2.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
virtual tmp< volScalarField > iDmdt(const phasePairKey &key) const
Return the interfacial mass transfer rate for a pair for a pair.
const dimensionSet dimVolume(pow3(dimLength))
phaseSystem::massTransferTable & massTransfer(massTransferPtr())
const phaseModel & phase1() const
Return phase 1.
const dimensionedScalar & rho() const
Return const-access to phase1 density.