Go to the documentation of this file.
41 namespace functionObjects
62 const bool isNew = !result1;
152 resultName1_(
"ObukhovLength"),
153 resultName2_(
"Ustar"),
161 dict.getCheckOrDefault<scalar>
165 [&](
const scalar
x){ return x > SMALL; }
186 UName_ =
dict.getOrDefault<
word>(
"U",
"U");
187 resultName1_ =
dict.getOrDefault<
word>(
"ObukhovLength",
"ObukhovLength");
188 resultName2_ =
dict.getOrDefault<
word>(
"Ustar",
"Ustar");
190 if (UName_ !=
"U" && resultName1_ ==
"ObukhovLength")
192 resultName1_ +=
'(' + UName_ +
')';
195 if (UName_ !=
"U" && resultName1_ ==
"Ustar")
197 resultName2_ +=
'(' + UName_ +
')';
200 rhoRef_ =
dict.getOrDefault<scalar>(
"rhoRef", 1.0);
201 kappa_ =
dict.getOrDefault<scalar>(
"kappa", 0.4);
202 beta_.value() =
dict.getOrDefault<scalar>(
"beta", 3
e-3);
216 if (isNew)
Log <<
" (new)" <<
nl <<
endl;
224 const auto* ioptr1 = mesh_.cfindObject<
regIOobject>(resultName1_);
225 const auto* ioptr2 = mesh_.cfindObject<
regIOobject>(resultName2_);
230 <<
" writing field " << ioptr1->name() <<
nl
231 <<
" writing field " << ioptr2->
name() <<
endl;
243 mesh_.thisDb().checkOut(resultName1_);
244 mesh_.thisDb().checkOut(resultName2_);
250 if (&mpm.
mesh() == &mesh_)
252 removeObukhovLength();
261 removeObukhovLength();
word resultName1_
Name of the output field for ObukhovLength.
Type * getObjectPtr(const word &name, const bool recursive=false) const
bool calcOL()
Hard-coded Obukhov length field and friction velocity.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
word resultName2_
Name of the output field for Ustar.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
const dimensionSet dimVelocity
virtual bool read(const dictionary &dict)
Read the data.
bool read(const char *buf, int32_t &val)
Same as readInt32.
static word timeName(const scalar t, const int precision=precision_)
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
scalar rhoRef_
Reference density (to convert from kinematic to static pressure)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar sign(const dimensionedScalar &ds)
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
Abstract base-class for Time/database function objects.
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
word UName_
Name of velocity field.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
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
void removeObukhovLength()
Remove (checkOut) the output fields from the object registry.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
ObukhovLength(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
virtual bool write(const bool valid=true) const
Write using setting from DB.
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual bool write()
Write the output fields.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar kappa_
von Kármán constant [-]
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const Type & lookupObject(const word &name, const bool recursive=false) const
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual void movePoints(const polyMesh &m)
Update for mesh point-motion.
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.
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
const word & name() const noexcept
Return name.
dimensionedScalar beta_
Thermal expansion coefficient [1/K].
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const Type * findObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
Computes the Obukhov length field and associated friction velocity field.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
virtual bool execute()
Calculate the output fields.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
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.
const dimensionedScalar e
Elementary charge.
defineTypeNameAndDebug(ObukhovLength, 0)
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const Time & time() const
Return the top-level database.
const dimensionSet & dimensions() const
Return const reference to dimensions.
const polyMesh & mesh() const
Return polyMesh.
static const gravity & New(const Time &runTime)
Construct on Time.
const dimensionSet dimless
Dimensionless.
virtual void updateMesh(const mapPolyMesh &mpm)
Update for changes of mesh.
const dimensionedVector g_
Gravitational acceleration vector [m/s2].