Go to the documentation of this file.
39 namespace functionObjects
49 void Foam::functionObjects::momentum::purgeFields()
57 template<
class GeoField>
59 Foam::functionObjects::momentum::newField
84 void Foam::functionObjects::momentum::calc()
101 const auto&
U = lookupObject<volVectorField>(UName_);
102 const auto*
rhoPtr = findObject<volScalarField>(rhoName_);
120 auto* pmomentum = getObjectPtr<volVectorField>(scopedName(
"momentum"));
125 pmomentum = tmomentum.
get();
127 auto& momentum = *pmomentum;
131 momentum.ref() = (
U * mesh_.V() * (*rhoPtr));
135 momentum.ref() = (
U * mesh_.V() * rhoRef);
137 momentum.correctBoundaryConditions();
144 getObjectPtr<volVectorField>(scopedName(
"angularMomentum"));
146 if (hasCsys_ && !pAngularMom)
150 pAngularMom = tAngularMom.
get();
152 else if (!pAngularMom)
156 pAngularMom = pmomentum;
158 auto& angularMom = *pAngularMom;
165 getObjectPtr<volVectorField>(scopedName(
"angularVelocity"));
172 newField<volVectorField>(
"angularVelocity",
dimVelocity);
173 pAngularVel = tAngularVel.
get();
175 auto& angularVel = *pAngularVel;
180 angularVel.primitiveFieldRef() =
181 csys_.invTransform(mesh_.cellCentres(),
U.internalField());
183 angularVel.correctBoundaryConditions();
187 angularMom.ref() = (angularVel * mesh_.V() * (*rhoPtr));
191 angularMom.ref() = (angularVel * mesh_.V() * rhoRef);
194 angularMom.correctBoundaryConditions();
201 sumAngularMom_ =
Zero;
208 for (
const label celli : cellIDs())
210 sumMomentum_ += momentum[celli];
211 sumAngularMom_ += angularMom[celli];
217 for (label celli=0; celli < mesh_.nCells(); ++celli)
219 sumMomentum_ += momentum[celli];
220 sumAngularMom_ += angularMom[celli];
226 reduce(sumMomentum_, sumOp<vector>());
227 reduce(sumAngularMom_, sumOp<vector>());
235 if (!writeToFile() || writtenHeader_)
243 writeHeaderValue(os,
"origin", csys_.origin());
244 writeHeaderValue(os,
"axis", csys_.e3());
251 if (regionType_ != vrtAll)
256 "Selection " + regionTypeNames_[regionType_]
257 +
" = " + regionName_
262 writeCommented(os,
"Time");
263 writeTabbed(os,
"(momentum_x momentum_y momentum_z)");
267 writeTabbed(os,
"(momentum_r momentum_rtheta momentum_axis)");
270 writeTabbed(os,
"volume");
273 writtenHeader_ =
true;
284 if (!foundObject<volVectorField>(UName_))
287 <<
"Could not find U: " << UName_ <<
" in database"
292 const auto* pPtr = cfindObject<volScalarField>(pName_);
298 if (!foundObject<volScalarField>(rhoName_))
301 <<
"Could not find rho:" << rhoName_
316 Info<<
" Sum of Momentum";
318 if (regionType_ != vrtAll)
320 Info<<
' ' << regionTypeNames_[regionType_]
321 <<
' ' << regionName_;
325 <<
" linear : " << sumMomentum_ <<
nl;
329 Info<<
" angular : " << sumAngularMom_ <<
nl;
337 writeCurrentTime(os);
339 os <<
tab << sumMomentum_;
343 os <<
tab << sumAngularMom_;
365 sumAngularMom_(
Zero),
372 writeMomentum_(
false),
373 writeVelocity_(
false),
374 writePosition_(
false),
397 sumAngularMom_(
Zero),
404 writeMomentum_(
false),
405 writeVelocity_(
false),
406 writePosition_(
false),
425 initialised_ =
false;
430 UName_ =
dict.getOrDefault<
word>(
"U",
"U");
431 pName_ =
dict.getOrDefault<
word>(
"p",
"p");
432 rhoName_ =
dict.getOrDefault<
word>(
"rho",
"rho");
433 rhoRef_ =
dict.getOrDefault<scalar>(
"rhoRef", 1.0);
434 hasCsys_ =
dict.getOrDefault(
"cylindrical",
false);
441 writeMomentum_ =
dict.getOrDefault(
"writeMomentum",
false);
442 writeVelocity_ =
dict.getOrDefault(
"writeVelocity",
false);
443 writePosition_ =
dict.getOrDefault(
"writePosition",
false);
445 Info<<
"Integrating for selection: "
446 << regionTypeNames_[regionType_]
447 <<
" (" << regionName_ <<
")" <<
nl;
451 Info<<
" Momentum fields will be written" <<
endl;
453 mesh_.objectRegistry::store
460 mesh_.objectRegistry::store
471 Info<<
" Angular velocity will be written" <<
endl;
473 mesh_.objectRegistry::store
475 newField<volVectorField>(
"angularVelocity",
dimVelocity)
481 Info<<
" Angular position will be written" <<
endl;
495 writeFileHeader(file());
503 setResult(
"momentum_x", sumMomentum_[0]);
504 setResult(
"momentum_y", sumMomentum_[1]);
505 setResult(
"momentum_z", sumMomentum_[2]);
507 setResult(
"momentum_r", sumAngularMom_[0]);
508 setResult(
"momentum_rtheta", sumAngularMom_[1]);
509 setResult(
"momentum_axis", sumAngularMom_[2]);
517 if (writeMomentum_ || (hasCsys_ && (writeVelocity_ || writePosition_)))
519 Log <<
"Writing fields" <<
nl;
523 fieldPtr = findObject<volVectorField>(scopedName(
"momentum"));
524 if (fieldPtr) fieldPtr->write();
526 fieldPtr = findObject<volVectorField>(scopedName(
"angularMomentum"));
527 if (fieldPtr) fieldPtr->write();
529 fieldPtr = findObject<volVectorField>(scopedName(
"angularVelocity"));
530 if (fieldPtr) fieldPtr->write();
532 if (hasCsys_ && writePosition_)
537 auto cyl_r = newField<volScalarField>(
"cyl_r",
dimLength);
538 auto cyl_t = newField<volScalarField>(
"cyl_theta",
dimless);
539 auto cyl_z = newField<volScalarField>(
"cyl_z",
dimLength);
543 const auto& pts = mesh_.cellCentres();
544 const label len = pts.size();
550 for (label i=0; i < len; ++i)
552 point p(csys_.localPosition(pts[i]));
565 const auto& pts = pbm[patchi].faceCentres();
566 const label len = pts.size();
572 for (label i=0; i < len; ++i)
574 point p(csys_.localPosition(pts[i]));
virtual void updateMesh(const mapPolyMesh &)
Update for changes of mesh.
static autoPtr< T > New(Args &&... args)
Construct autoPtr of T with forwarding arguments.
virtual bool read(const dictionary &)
Read the momentum data.
const dimensionSet dimPressure
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const objectRegistry & obr_
Reference to the region objectRegistry.
virtual bool write()
Write the momentum, possibly angular momentum and velocity.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
bool update()
Update the cached values as required.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Volume (cell) region selection class.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVelocity
const dimensionSet dimDensity
bool read(const char *buf, int32_t &val)
Same as readInt32.
word scopedName(const word &base) const
Return a scoped name.
Ostream & endl(Ostream &os)
Add newline and flush stream.
momentum(const word &name, const Time &runTime, const dictionary &dict, const bool readFields=true)
Construct from Time and dictionary.
static void writeHeader(Ostream &os, const word &fieldName)
Dimension set for the base types.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual void updateMesh(const mapPolyMesh &)
Update for changes of mesh.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
Info<< "Reading mechanical properties\n"<< 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
virtual bool read(const dictionary &dict)
Read from dictionary.
Registry of regIOobjects.
virtual bool read(const dictionary &dict)
Read.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual bool execute()
Calculate and report the integral momentum.
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
T * get() noexcept
Return pointer to managed object without nullptr checking.
A cylindrical coordinate system (r-theta-z). The coordinate system angle theta is always in radians.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool read(const dictionary &dict)
Read optional controls.
Macros for easy insertion into run-time selection tables.
Generic dimensioned Type class.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
virtual void writeFileHeader(Ostream &os)
Output file header information.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static bool master(const label communicator=0)
Am I the master process.
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
void initialise()
Initialise the fields.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
void writeValues(Ostream &os)
Write momentum data.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const fvMesh & mesh_
Reference to the fvMesh.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
defineTypeNameAndDebug(ObukhovLength, 0)
Reads fields from the time directories and adds them to the mesh database for further post-processing...
bool checkOut(regIOobject *io) const
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Computes the natural logarithm of an input volScalarField.
Base class for writing single files from the function objects.
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.