46 void Foam::fv::massSourceBase::readCoeffs()
65 const basicThermo&
thermo =
70 heName_ =
thermo.he().name();
71 TName_ =
thermo.T().name();
77 void Foam::fv::massSourceBase::addGeneralSupType
85 const scalar massFlowRate = this->massFlowRate();
90 fieldValues_[fieldName]->value<Type>(mesh().time().userTimeValue());
94 eqn.source()[
cells[i]] -=
95 mesh().V()[
cells[i]]/set_.V()*massFlowRate*value;
102 eqn.diag()[
cells[i]] +=
103 mesh().V()[
cells[i]]/set_.V()*massFlowRate;
110 void Foam::fv::massSourceBase::addSupType
113 const word& fieldName
116 addGeneralSupType(eqn, fieldName);
120 void Foam::fv::massSourceBase::addSupType
122 fvMatrix<scalar>& eqn,
123 const word& fieldName
128 if (fieldName == rhoName_)
130 const scalar massFlowRate = this->massFlowRate();
134 eqn.source()[
cells[i]] -=
135 mesh().V()[
cells[i]]/set_.V()*massFlowRate;
138 else if (fieldName == heName_ && fieldValues_.found(TName_))
140 const scalar massFlowRate = this->massFlowRate();
142 if (massFlowRate > 0)
144 if (fieldValues_.found(heName_))
147 <<
"Source " <<
name() <<
" defined for both field "
148 << heName_ <<
" and " << TName_
149 <<
". Only one of these should be present." <<
endl;
152 const basicThermo&
thermo =
153 mesh().lookupObject<basicThermo>
157 physicalProperties::typeName,
163 fieldValues_[TName_]->value<scalar>
165 mesh().time().userTimeValue()
175 eqn.source()[
cells[i]] -=
176 mesh().V()[
cells[i]]/set_.V()*massFlowRate*hs[i];
183 eqn.diag()[
cells[i]] +=
184 mesh().V()[
cells[i]]/set_.V()*massFlowRate;
190 addGeneralSupType(eqn, fieldName);
196 void Foam::fv::massSourceBase::addSupType
200 const word& fieldName
203 addSupType(eqn, fieldName);
208 void Foam::fv::massSourceBase::addSupType
213 const word& fieldName
216 addSupType(eqn, fieldName);
220 void Foam::fv::massSource::readCoeffs()
233 Foam::scalar Foam::fv::massSource::massFlowRate()
const
235 return massFlowRate_->value(mesh().time().userTimeValue());
249 fieldValues_.clear();
267 const word& modelType,
287 const word& modelType,
308 && massFlowRate() > 0
309 && !(fieldName == rhoName_)
310 && !(fieldName == heName_ && fieldValues_.found(TName_))
311 && !fieldValues_.found(fieldName)
315 <<
"No value supplied for field " << fieldName <<
" in "
329 if (fieldValues_.found(TName_))
356 set_.topoChange(map);
368 set_.distribute(map);
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Macros for easy insertion into run-time selection tables.
static autoPtr< Function1< Type > > New(const word &name, const dictionary &dict)
Selector.
word group() const
Return group (extension part of name)
static word groupName(Name name, const word &group)
A list of keyword definitions, which are a keyword followed by any number of values (e....
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
General run-time selected cell set selection class for fvMesh.
Mesh data needed to do the Finite Volume discretisation.
Finite volume model abstract base class.
const dictionary & coeffs() const
Return dictionary.
virtual bool read(const dictionary &dict)
Read source dictionary.
const fvMesh & mesh() const
Return const access to the mesh database.
virtual bool movePoints()
Update for mesh motion.
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
void readSet()
Read the set.
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
void readFieldValues()
Read the field values.
virtual bool read(const dictionary &dict)
Read source dictionary.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
massSourceBase(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
Construct from explicit source name and mesh.
virtual bool addsSupToField(const word &fieldName) const
Return true if the fvModel adds a source term to the given.
This fvModel applies a mass source to the continuity equation and to all field equations.
virtual bool read(const dictionary &dict)
Read source dictionary.
massSource(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
Construct from explicit source name and mesh.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type and name.
bool foundObject(const word &name) const
Is the named Type in registry.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Wrapper around Function1 that constructs a function for an as yet unknown primitive type....
A class for handling words, derived from string.
static const word null
An empty word.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
#define IMPLEMENT_FV_MODEL_ADD_RHO_SUP(Type, modelType)
#define IMPLEMENT_FV_MODEL_ADD_SUP(Type, modelType)
#define IMPLEMENT_FV_MODEL_ADD_ALPHA_RHO_SUP(Type, modelType)
static List< word > fieldNames
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
#define WarningInFunction
Report a warning using Foam::Warning.
addToRunTimeSelectionTable(fvConstraint, bound, dictionary)
defineTypeNameAndDebug(bound, 0)
Ostream & endl(Ostream &os)
Add newline and flush stream.
FOR_ALL_FIELD_TYPES(DefineContiguousFvWallLocationDataType)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
VolField< scalar > volScalarField
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
word name(const complex &)
Return a string representation of a complex.
UList< label > labelUList
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
fluidMulticomponentThermo & thermo