45 { volumeModeType::vmAbsolute,
"absolute" },
46 { volumeModeType::vmSpecific,
"specific" },
55 const dictionary&
dict,
56 const word& keyExplicit,
57 const word& keyImplicit
62 fieldNames_.resize_nocopy(
count);
71 driverSu_.resize(2*
count);
72 driverSp_.resize(2*
count);
76 valueExprSu_.resize(2*
count);
77 valueExprSp_.resize(2*
count);
79 fv::option::resetApplied();
82 Tuple2<Type, scalar> sourceRates;
86 for (
const entry& dEntry :
dict)
88 const word& fieldName = dEntry.keyword();
93 const dictionary& subDict = dEntry.dict();
99 (eptr = subDict.findEntry(keyExplicit, keyType::LITERAL))
108 && eptr->dict().readEntry(
"type", modelType, keyType::LITERAL)
109 && (modelType ==
"exprField")
112 const dictionary& exprDict = eptr->dict();
114 valueExprSu_.emplace_set(fieldName);
115 valueExprSu_[fieldName].readEntry(
"expression", exprDict);
135 (eptr = subDict.findEntry(keyImplicit, keyType::LITERAL))
144 && eptr->dict().readEntry(
"type", modelType, keyType::LITERAL)
145 && (modelType ==
"exprField")
148 const dictionary& exprDict = eptr->dict();
150 valueExprSp_.emplace_set(fieldName);
151 valueExprSp_[fieldName].readEntry(
"expression", exprDict);
164 Function1<scalar>::New(keyImplicit, subDict, &mesh_)
173 dEntry.readEntry(sourceRates);
180 new Function1Types::Constant<Type>
189 new Function1Types::Constant<scalar>
200 <<
"Require at least one of " 201 << keyExplicit <<
'/' << keyImplicit <<
" entries for " 202 <<
"field: " << fieldName <<
endl 206 fieldNames_[
count] = fieldName;
218 const word& modelType,
224 volumeMode_(vmAbsolute),
255 <<
">::addSup for source " << name_ <<
endl;
261 const word& fieldName = fieldNames_[fieldi];
263 const scalar tmVal = mesh_.time().timeOutputValue();
270 const auto iter1 = valueExprSu_.cfind(fieldName);
271 const auto iter2 = Su_.cfind(fieldName);
277 const auto& valueExpr = iter1.val();
285 Info<<
"Explicit expression source:" <<
nl 287 << valueExpr.c_str() <<
nl 291 auto& driver = *(driverSu_[fieldName]);
293 driver.clearVariables();
297 driver.addContextObject(
"rho", &
rho);
303 driver.parse(valueExpr);
308 const ExprResultType* ptr =
309 driver.template isResultType<ExprResultType>();
314 <<
"Expression for Su " << fieldName
315 <<
" evaluated to <" << driver.resultType()
316 <<
"> but expected <" << ExprResultType::typeName
320 else if (ptr->size() != mesh_.nCells())
323 <<
"Expression for Su " << fieldName
324 <<
" evaluated to " << ptr->size()
325 <<
" instead of " << mesh_.nCells() <<
" values" <<
endl 331 driver.removeContextObject(&
rho);
334 const Field<Type>& exprFld = ptr->primitiveField();
338 name_ + fieldName +
"Su",
340 dimensioned<Type>(SuDims,
Zero)
343 if (this->useSubMesh())
345 for (
const label celli : cells_)
347 tsu.
ref()[celli] = exprFld[celli]/VDash_;
352 tsu.
ref().field() = exprFld;
354 if (!
equal(VDash_, 1))
356 tsu.
ref().field() /= VDash_;
360 else if (iter2.good() && iter2.val()->good())
362 const dimensioned<Type> SuValue
366 iter2.val()->value(tmVal)/VDash_
369 if (
mag(SuValue.value()) <= ROOTVSMALL)
373 else if (this->useSubMesh())
377 name_ + fieldName +
"Su",
379 dimensioned<Type>(SuDims,
Zero)
381 UIndirectList<Type>(tsu.
ref(), cells_) = SuValue.value();
398 const auto iter1 = valueExprSp_.cfind(fieldName);
399 const auto iter2 = Sp_.cfind(fieldName);
401 tmp<DimensionedField<scalar, volMesh>> tsp;
405 const auto& valueExpr = iter1.val();
411 Info<<
"Implicit expression source:" <<
nl 413 << valueExpr.c_str() <<
nl 417 auto& driver = *(driverSp_[fieldName]);
419 driver.clearVariables();
423 driver.addContextObject(
"rho", &
rho);
429 driver.parse(valueExpr);
434 const ExprResultType* ptr =
435 driver.template isResultType<ExprResultType>();
440 <<
"Expression for Sp " << fieldName
441 <<
" evaluated to <" << driver.resultType()
442 <<
"> but expected <" << ExprResultType::typeName
446 else if (ptr->size() != mesh_.nCells())
449 <<
"Expression for Sp " << fieldName
450 <<
" evaluated to " << ptr->size()
451 <<
" instead of " << mesh_.nCells() <<
" values" <<
endl 457 driver.removeContextObject(&
rho);
460 const Field<scalar>& exprFld = ptr->primitiveField();
464 name_ + fieldName +
"Sp",
466 dimensioned<scalar>(SpDims,
Zero)
469 if (this->useSubMesh())
471 for (
const label celli : cells_)
473 tsp.ref()[celli] = exprFld[celli]/VDash_;
478 tsp.ref().field() = exprFld;
480 if (!
equal(VDash_, 1))
482 tsp.ref().field() /= VDash_;
486 else if (iter2.good() && iter2.val()->good())
488 const dimensioned<scalar> SpValue
492 iter2.val()->value(tmVal)/VDash_
495 if (
mag(SpValue.value()) <= ROOTVSMALL)
499 else if (this->useSubMesh())
503 name_ + fieldName +
"Sp",
505 dimensioned<scalar>(SpDims,
Zero)
507 UIndirectList<scalar>(tsp.ref(), cells_) = SpValue.value();
530 volumeMode_ = volumeModeTypeNames_.get(
"volumeMode", coeffs_);
533 if (volumeMode_ == vmAbsolute)
539 const dictionary* injectDict =
bool valid() const noexcept
Identical to good(), or bool operator.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
bool equal(const T &a, const T &b)
Compare two values for equality.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
T & ref() const
Return non-const reference to the contents of a non-null managed pointer.
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A traits class, which is primarily used for primitives.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
virtual bool read(const dictionary &dict)
Read source dictionary.
virtual bool read(const dictionary &dict)
Read source dictionary.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const dimensionSet dimVolume(pow3(dimLength))
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
const GeometricField< Type, fvPatchField, volMesh > & psi(const label i=0) const
Return psi.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
volumeExpr::parseDriver volumeExprDriver
Typedef for volumeExpr parseDriver.
static const Enum< volumeModeType > volumeModeTypeNames_
Names for volumeModeType.
A class for handling words, derived from Foam::string.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
zeroField SuSp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
int debug
Static debugging option.
static const GeometricField< scalar, fvPatchField, volMesh > & null()
Return a null geometric field.
static bool checking() noexcept
True if dimension checking is enabled (the usual default)
const dimensionSet & dimensions() const noexcept
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
volumeModeType
Options for the volume mode type.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Mesh data needed to do the Finite Volume discretisation.
virtual void addSup(fvMatrix< Type > &eqn, const label fieldi)
Add explicit contribution to incompressible equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Intermediate abstract class for handling cell-set options for the derived fvOptions.
const volScalarField & psi
Applies semi-implicit source within a specified region for Type, where <Type>=Scalar/Vector/Spherical...
A class for managing temporary objects.
SemiImplicitSource(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &dims, const Field< Type > &iField)
Return tmp field from name, mesh, dimensions, copy of internal field.
bool notNull(const T *ptr)
True if ptr is not a pointer (of type T) to the nullObject.
Calculate the finiteVolume matrix for implicit and explicit sources.
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
static constexpr const zero Zero
Global zero (0)