42 { equilibriumModelType::LANGMUIR,
"Langmuir" }
52 { kineticModelType::PseudoFirstOrder,
"PseudoFirstOrder" }
59Foam::speciesSorptionFvPatchScalarField::calcMoleFractions()
const
62 auto& Mole = tMole.ref();
71 const PtrList<volScalarField>&
Y =
thermo.composition().Y();
77 const label speciesId =
83 thermo.composition().W(speciesId)
90 const label
cellId = faceCells[i];
97 <<
"Thermo type is not 'rhoReactionThermo'. " <<
nl
98 <<
"This BC is designed to operate with a rho based thermo."
109 const word& fieldName,
110 const dimensionSet& dim
113 const fvMesh&
mesh = this->internalField().mesh();
147 zeroGradientFvPatchScalarField(
p, iF),
150 thicknessPtr_(nullptr),
168 zeroGradientFvPatchScalarField(
p, iF,
dict),
169 equilibriumModel_(equilibriumModelTypeNames.get(
"equilibriumModel",
dict)),
170 kinematicModel_(kinematicModelTypeNames.get(
"kinematicModel",
dict)),
175 rhoS_(
dict.get<scalar>(
"rhoS")),
176 pName_(
dict.getOrDefault<
word>(
"p",
"p")),
212 zeroGradientFvPatchScalarField(ptf,
p, iF, mapper),
213 equilibriumModel_(ptf.equilibriumModel_),
214 kinematicModel_(ptf.kinematicModel_),
215 thicknessPtr_(ptf.thicknessPtr_.clone(patch().patch())),
221 dfldp_(ptf.dfldp_, mapper),
222 mass_(ptf.mass_, mapper)
231 zeroGradientFvPatchScalarField(ptf),
232 equilibriumModel_(ptf.equilibriumModel_),
233 kinematicModel_(ptf.kinematicModel_),
234 thicknessPtr_(ptf.thicknessPtr_.clone(patch().patch())),
251 zeroGradientFvPatchScalarField(ptf, iF),
252 equilibriumModel_(ptf.equilibriumModel_),
253 kinematicModel_(ptf.kinematicModel_),
254 thicknessPtr_(ptf.thicknessPtr_.clone(patch().patch())),
272 zeroGradientFvPatchScalarField::autoMap(m);
279 thicknessPtr_->autoMap(m);
290 zeroGradientFvPatchScalarField::rmap(ptf, addr);
292 const auto& tiptf = refCast<const speciesSorptionFvPatchScalarField>(ptf);
294 dfldp_.rmap(tiptf.dfldp_, addr);
295 mass_.rmap(tiptf.mass_, addr);
299 thicknessPtr_->rmap(tiptf.thicknessPtr_(), addr);
312 const label speciesId =
313 thermo.composition().species()[this->internalField().name()];
315 const scalar Wi(
thermo.composition().W(speciesId));
317 const scalar t = db().time().timeOutputValue();
330 const label faceCelli = this->patch().faceCells()[facei];
331 Vol[facei] = this->internalField().mesh().V()[faceCelli];
340 Info<<
" Patch mass rate min/max [kg/m3/sec]: "
365 switch (equilibriumModel_)
367 case equilibriumModelType::LANGMUIR:
375 cEq = max_*(kl_*tco()*pp/(1 + kl_*tco()*pp));
385 switch (kinematicModel_)
387 case kineticModelType::PseudoFirstOrder:
389 dfldp_ = kabs_*(cEq - mass_);
396 const scalar dt = db().time().deltaTValue();
398 mass_ =
max(mass_, scalar(0));
403 "absorbedMass" + this->internalField().
name(),
405 ).boundaryFieldRef()[patch().index()];
411 Info<<
" Absorption rate min/max [mol/kg/sec]: "
415 zeroGradientFvPatchScalarField::updateCoeffs();
425 "equilibriumModel", equilibriumModelTypeNames[equilibriumModel_]
429 "kinematicModel", kinematicModelTypeNames[kinematicModel_]
433 thicknessPtr_->writeData(
os);
444 writeEntry(
"value",
os);
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Field< Type > & field() const
Return field.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
label size() const noexcept
The number of elements in the UList.
void size(const label n)
Older name for setAddressableSize.
static const word dictName
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
virtual bool write()
Write the output fields.
const Time & time() const
Return the top-level database.
A FieldMapper for finite-volume patch fields.
const objectRegistry & db() const
Return local objectRegistry.
virtual void operator=(const UList< Type > &)
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
friend Ostream & operator(Ostream &, const fvPatchField< scalar > &)
const fvPatch & patch() const
Return patch.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual const labelUList & faceCells() const
Return faceCells.
const Type & lookupObject(const word &name, const bool recursive=false) const
Type * getObjectPtr(const word &name, const bool recursive=false) const
This boundary condition provides a first-order zero-gradient condition for a given scalar field to mo...
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
static const Enum< kineticModelType > kinematicModelTypeNames
Names for kineticModelType.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual tmp< scalarField > patchSource() const
Source of cells next to the patch.
tmp< scalarField > mass() const
Access to mass.
equilibriumModelType
Options for the equilibrum model.
static const Enum< equilibriumModelType > equilibriumModelTypeNames
Names for equilibriumModelType.
kineticModelType
Options for the kinematic model.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
PtrList< volScalarField > & Y
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
static constexpr const zero Zero
Global zero (0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Type gMin(const FieldField< Field, Type > &f)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
errorManipArg< error, int > exit(error &err, const int errNo=1)
UList< label > labelUList
A UList of labels.
Type gMax(const FieldField< Field, Type > &f)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.