40 namespace functionObjects
50 void Foam::functionObjects::momentum::purgeFields()
58 template<
class GeoField>
60 Foam::functionObjects::momentum::newField
86 void Foam::functionObjects::momentum::calc()
99 autoPtr<volVectorField> tmomentum, tAngularMom, tAngularVel;
103 const auto&
U = lookupObject<volVectorField>(UName_);
104 const auto*
rhoPtr = findObject<volScalarField>(rhoName_);
122 auto* pmomentum = getObjectPtr<volVectorField>(scopedName(
"momentum"));
127 pmomentum = tmomentum.get();
129 auto& momentum = *pmomentum;
133 momentum.ref() = (
U * mesh_.V() * (*rhoPtr));
137 momentum.ref() = (
U * mesh_.V() * rhoRef);
139 momentum.correctBoundaryConditions();
146 getObjectPtr<volVectorField>(scopedName(
"angularMomentum"));
148 if (hasCsys_ && !pAngularMom)
152 pAngularMom = tAngularMom.get();
154 else if (!pAngularMom)
158 pAngularMom = pmomentum;
160 auto& angularMom = *pAngularMom;
167 getObjectPtr<volVectorField>(scopedName(
"angularVelocity"));
174 newField<volVectorField>(
"angularVelocity",
dimVelocity);
175 pAngularVel = tAngularVel.get();
177 auto& angularVel = *pAngularVel;
182 angularVel.primitiveFieldRef() =
183 csys_.invTransform(mesh_.cellCentres(),
U.internalField());
185 angularVel.correctBoundaryConditions();
189 angularMom.ref() = (angularVel * mesh_.V() * (*rhoPtr));
193 angularMom.ref() = (angularVel * mesh_.V() * rhoRef);
196 angularMom.correctBoundaryConditions();
203 sumAngularMom_ =
Zero;
207 for (label celli=0; celli < mesh_.nCells(); ++celli)
209 sumMomentum_ += momentum[celli];
210 sumAngularMom_ += angularMom[celli];
215 for (
const label celli :
cellIDs())
217 sumMomentum_ += momentum[celli];
218 sumAngularMom_ += angularMom[celli];
231 if (!writeToFile() || writtenHeader_)
239 writeHeaderValue(
os,
"origin", csys_.origin());
240 writeHeaderValue(
os,
"axis", csys_.e3());
252 "Selection " + regionTypeNames_[regionType_]
253 +
" = " + regionName_
258 writeCommented(
os,
"Time");
259 writeTabbed(
os,
"(momentum_x momentum_y momentum_z)");
263 writeTabbed(
os,
"(momentum_r momentum_rtheta momentum_axis)");
266 writeTabbed(
os,
"volume");
269 writtenHeader_ =
true;
280 if (!foundObject<volVectorField>(UName_))
283 <<
"Could not find U: " << UName_ <<
" in database" 288 const auto* pPtr = cfindObject<volScalarField>(pName_);
294 if (!foundObject<volScalarField>(rhoName_))
297 <<
"Could not find rho:" << rhoName_
312 Info<<
" Sum of Momentum";
316 Info<<
' ' << regionTypeNames_[regionType_]
317 <<
' ' << regionName_;
321 <<
" linear : " << sumMomentum_ <<
nl;
325 Info<<
" angular : " << sumAngularMom_ <<
nl;
333 writeCurrentTime(
os);
335 os <<
tab << sumMomentum_;
339 os <<
tab << sumAngularMom_;
358 volRegion(fvMeshFunctionObject::mesh_,
dict),
359 writeFile(mesh_,
name, typeName,
dict),
361 sumAngularMom_(
Zero),
368 writeMomentum_(false),
369 writeVelocity_(false),
370 writePosition_(false),
389 fvMeshFunctionObject(
name, obr,
dict),
390 volRegion(fvMeshFunctionObject::mesh_,
dict),
391 writeFile(mesh_,
name, typeName,
dict),
393 sumAngularMom_(
Zero),
400 writeMomentum_(false),
401 writeVelocity_(false),
402 writePosition_(false),
426 UName_ =
dict.getOrDefault<
word>(
"U",
"U");
427 pName_ =
dict.getOrDefault<
word>(
"p",
"p");
428 rhoName_ =
dict.getOrDefault<
word>(
"rho",
"rho");
429 rhoRef_ =
dict.getOrDefault<scalar>(
"rhoRef", 1.0);
430 hasCsys_ =
dict.getOrDefault(
"cylindrical",
false);
437 writeMomentum_ =
dict.getOrDefault(
"writeMomentum",
false);
438 writeVelocity_ =
dict.getOrDefault(
"writeVelocity",
false);
439 writePosition_ =
dict.getOrDefault(
"writePosition",
false);
441 Info<<
"Integrating for selection: " 442 << regionTypeNames_[regionType_]
443 <<
" (" << regionName_ <<
")" <<
nl;
447 Info<<
" Momentum fields will be written" <<
endl;
467 Info<<
" Angular velocity will be written" <<
endl;
471 newField<volVectorField>(
"angularVelocity",
dimVelocity)
477 Info<<
" Angular position will be written" <<
endl;
491 writeFileHeader(file());
499 setResult(
"momentum_x", sumMomentum_[0]);
500 setResult(
"momentum_y", sumMomentum_[1]);
501 setResult(
"momentum_z", sumMomentum_[2]);
503 setResult(
"momentum_r", sumAngularMom_[0]);
504 setResult(
"momentum_rtheta", sumAngularMom_[1]);
505 setResult(
"momentum_axis", sumAngularMom_[2]);
513 if (writeMomentum_ || (hasCsys_ && (writeVelocity_ || writePosition_)))
515 Log <<
"Writing fields" <<
nl;
519 fieldPtr = findObject<volVectorField>(scopedName(
"momentum"));
520 if (fieldPtr) fieldPtr->
write();
522 fieldPtr = findObject<volVectorField>(scopedName(
"angularMomentum"));
523 if (fieldPtr) fieldPtr->
write();
525 fieldPtr = findObject<volVectorField>(scopedName(
"angularVelocity"));
526 if (fieldPtr) fieldPtr->
write();
528 if (hasCsys_ && writePosition_)
533 auto cyl_r = newField<volScalarField>(
"cyl_r",
dimLength);
534 auto cyl_t = newField<volScalarField>(
"cyl_theta",
dimless);
535 auto cyl_z = newField<volScalarField>(
"cyl_z",
dimLength);
539 const auto&
pts = mesh_.cellCentres();
540 const label len =
pts.size();
546 for (label i=0; i < len; ++i)
557 const polyBoundaryMesh&
pbm = mesh_.boundaryMesh();
561 if (isA<emptyPolyPatch>(
pbm[patchi]))
566 const auto&
pts =
pbm[patchi].faceCentres();
567 const label len =
pts.size();
569 UList<scalar>& r = cyl_r->boundaryFieldRef(
false)[patchi];
570 UList<scalar>& t = cyl_t->boundaryFieldRef(
false)[patchi];
571 UList<scalar>& z = cyl_z->boundaryFieldRef(
false)[patchi];
573 for (label i=0; i < len; ++i)
const polyBoundaryMesh & pbm
defineTypeNameAndDebug(ObukhovLength, 0)
static bool initialised_(false)
static void writeHeader(Ostream &os, const word &fieldName)
dimensionedScalar log(const dimensionedScalar &ds)
A cylindrical coordinate system (r-theta-z). The coordinate system angle theta is always in radians...
errorManipArg< error, int > exit(error &err, const int errNo=1)
momentum(const word &name, const Time &runTime, const dictionary &dict, const bool readFields=true)
Construct from Time and dictionary.
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...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual bool read(const dictionary &dict)
Read from dictionary.
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool store()
Register object with its registry and transfer ownership to the registry.
virtual bool execute()
Calculate and report the integral momentum.
constexpr char tab
The tab '\t' character(0x09)
virtual void updateMesh(const mapPolyMesh &)
Update for changes of mesh.
bool checkOut(regIOobject *io) const
Remove a regIOobject from registry and free memory if the object is ownedByRegistry. A nullptr is ignored.
Ignore writing from objectRegistry::writeObject()
const dimensionSet dimless
Dimensionless.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
bool useAllCells() const noexcept
Use all cells, not the cellIDs.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const NameMatchPredicate &selectedFields, DynamicList< regIOobject *> &storedObjects)
Read the selected GeometricFields of the templated type and store on the objectRegistry.
A class for handling words, derived from Foam::string.
void writeValues(Ostream &os)
Write momentum data.
virtual bool write(const bool writeOnProc=true) const
Write using setting from DB.
const dimensionSet dimPressure
void initialise()
Initialise the fields.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
An Ostream is an abstract base class for all output systems (streams, files, token lists...
scalar V() const
Return total volume of the selected region.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
OBJstream os(runTime.globalPath()/outputName)
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
virtual bool read(const dictionary &dict)
Read.
const dimensionSet dimDensity
vector point
Point is a vector.
Info<< "Reading mechanical properties\"<< endl;IOdictionary mechanicalProperties(IOobject("mechanicalProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));const dictionary &rhoDict(mechanicalProperties.subDict("rho"));word rhoType(rhoDict.get< word >"type"));autoPtr< volScalarField > rhoPtr
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
bool update()
Update the cached values as required.
static bool master(const label communicator=worldComm)
True if process corresponds to the master rank in the communicator.
const objectRegistry & obr_
Reference to the region objectRegistry.
virtual bool write()
Write the momentum, possibly angular momentum and velocity.
messageStream Info
Information stream (stdout output on master, null elsewhere)
virtual bool read(const dictionary &dict)
Read optional controls.
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
virtual void updateMesh(const mapPolyMesh &)
Update for changes of mesh.
Mesh consisting of general polyhedral cells.
virtual bool read(const dictionary &)
Read the momentum data.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
virtual void writeFileHeader(Ostream &os)
Output file header information.
Registry of regIOobjects.
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
Defines the attributes of an object for which implicit objectRegistry management is supported...
word scopedName(const word &name) const
Return a scoped (prefixed) name.
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVelocity