Go to the documentation of this file.
38 namespace functionObjects
44 proudmanAcousticPower,
53 Foam::functionObjects::proudmanAcousticPower::rhoScale
55 const tmp<volScalarField>&
fld
62 return fld*thermoPtr->rho();
65 if (rhoInf_.
value() < 0)
69 <<
"Incompressible calculation assumed, but no reference density "
70 <<
"set. Please set the entry 'rhoInf' to an appropriate value"
80 Foam::functionObjects::proudmanAcousticPower::a()
const
86 const basicThermo&
thermo = *thermoPtr;
96 mesh_.time().timeName(),
107 Foam::functionObjects::proudmanAcousticPower::k()
const
109 if (kName_ !=
"none")
111 return lookupObject<volScalarField>(kName_);
122 Foam::functionObjects::proudmanAcousticPower::epsilon()
const
124 if (epsilonName_ !=
"none")
126 return lookupObject<volScalarField>(epsilonName_);
129 if (omegaName_ !=
"none")
132 const auto& omega = lookupObject<volScalarField>(omegaName_);
133 const scalar betaStar = 0.09;
134 return betaStar*
k()*omega;
140 return turb.epsilon();
158 epsilonName_(
"none"),
170 mesh_.time().timeName(),
189 mesh_.time().timeName(),
211 dict.readIfPresent(
"alphaEps", alphaEps_);
212 rhoInf_.readIfPresent(
"rhoInf",
dict);
213 aRef_.readIfPresent(
"aRef",
dict);
215 if (
dict.readIfPresent(
"k", kName_))
217 Info<<
" k field: " << kName_ <<
endl;
221 Info<<
" k field from turbulence model" <<
endl;
224 if (
dict.readIfPresent(
"epsilon", epsilonName_))
226 Info<<
" epsilon field: " << epsilonName_ <<
endl;
230 Info<<
" epsilon field from turbulence model (if needed)"
234 if (
dict.readIfPresent(
"omega", omegaName_))
236 Info<<
" omega field: " << omegaName_ <<
endl;
240 Info<<
" omega field from turbulence model (if needed)" <<
endl;
243 if (epsilonName_ !=
"none" && omegaName_ !=
"none")
246 <<
"either epsilon or omega field names can be set but not both"
263 auto& P_A = mesh_.lookupObjectRef<
volScalarField>(scopedName(
"P_A"));
267 auto& L_P = mesh_.lookupObjectRef<
volScalarField>(scopedName(
"L_P"));
279 const auto& P_A = mesh_.lookupObject<
volScalarField>(scopedName(
"P_A"));
281 Log <<
" writing field " << P_A.name() <<
nl;
285 const auto& L_P = mesh_.lookupObject<
volScalarField>(scopedName(
"L_P"));
287 Log <<
" writing field " << L_P.name() <<
nl;
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual bool execute()
Calculate the Proudman acoustic power.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
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.
static const word propertiesName
Default name of the turbulence properties dictionary.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual bool write()
Write the Proudman acoustic power.
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensionedScalar log10(const dimensionedScalar &ds)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const dimensionSet dimPower
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
dimensionedScalar pow5(const dimensionedScalar &ds)
virtual bool read(const dictionary &dict)
Read optional controls.
Macros for easy insertion into run-time selection tables.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const word & name() const noexcept
Return the name of this functionObject.
proudmanAcousticPower(const word &name, const Time &runTime, const dictionary &)
Construct from Time and dictionary.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
virtual const word & type() const =0
Runtime type information.
dimensionedScalar sqrt(const dimensionedScalar &ds)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
label k
Boltzmann constant.
const dimensionedScalar e
Elementary charge.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
defineTypeNameAndDebug(ObukhovLength, 0)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
const dimensionSet dimVolume(pow3(dimLength))
compressible::turbulenceModel & turb
virtual bool read(const dictionary &)
Read the Proudman acoustic power data.
const dimensionSet dimless
Dimensionless.
word dictName() const
The local dictionary name (final part of scoped name)