proudmanAcousticPower.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2019-2021 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "proudmanAcousticPower.H"
29 #include "volFields.H"
30 #include "basicThermo.H"
31 #include "turbulenceModel.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace functionObjects
39 {
40  defineTypeNameAndDebug(proudmanAcousticPower, 0);
42  (
43  functionObject,
44  proudmanAcousticPower,
45  dictionary
46  );
47 }
48 }
49 
50 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
53 Foam::functionObjects::proudmanAcousticPower::rhoScale
54 (
55  const tmp<volScalarField>& fld
56 ) const
57 {
58  const auto* thermoPtr = getObjectPtr<basicThermo>(basicThermo::dictName);
59 
60  if (thermoPtr)
61  {
62  return fld*thermoPtr->rho();
63  }
64 
65  if (rhoInf_.value() < 0)
66  {
68  << type() << " " << name() << ": "
69  << "Incompressible calculation assumed, but no reference density "
70  << "set. Please set the entry 'rhoInf' to an appropriate value"
71  << nl
72  << exit(FatalError);
73  }
74 
75  return rhoInf_*fld;
76 }
77 
78 
80 Foam::functionObjects::proudmanAcousticPower::a() const
81 {
82  const auto* thermoPtr = getObjectPtr<basicThermo>(basicThermo::dictName);
83 
84  if (thermoPtr)
85  {
86  const basicThermo& thermo = *thermoPtr;
87  return sqrt(thermo.gamma()*thermo.p()/thermo.rho());
88  }
89 
90  return
92  (
93  IOobject
94  (
95  scopedName("a"),
96  mesh_.time().timeName(),
97  mesh_,
99  ),
100  mesh_,
101  aRef_
102  );
103 }
104 
105 
107 Foam::functionObjects::proudmanAcousticPower::k() const
108 {
109  if (kName_ != "none")
110  {
111  return lookupObject<volScalarField>(kName_);
112  }
113 
114  const auto& turb =
115  lookupObject<turbulenceModel>(turbulenceModel::propertiesName);
116 
117  return turb.k();
118 }
119 
120 
122 Foam::functionObjects::proudmanAcousticPower::epsilon() const
123 {
124  if (epsilonName_ != "none")
125  {
126  return lookupObject<volScalarField>(epsilonName_);
127  }
128 
129  if (omegaName_ != "none")
130  {
131  // Construct epsilon on-the-fly
132  const auto& omega = lookupObject<volScalarField>(omegaName_);
133  const scalar betaStar = 0.09;
134  return betaStar*k()*omega;
135  }
136 
137  const auto& turb =
138  lookupObject<turbulenceModel>(turbulenceModel::propertiesName);
139 
140  return turb.epsilon();
141 }
142 
143 
144 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
145 
147 (
148  const word& name,
149  const Time& runTime,
150  const dictionary& dict
151 )
152 :
154  alphaEps_(0.1),
155  rhoInf_("0", dimDensity, -1),
156  aRef_(dimVelocity, Zero),
157  kName_("none"),
158  epsilonName_("none"),
159  omegaName_("none")
160 {
161  read(dict);
162 
163  volScalarField* PAPtr
164  (
165  new volScalarField
166  (
167  IOobject
168  (
169  scopedName("P_A"),
170  mesh_.time().timeName(),
171  mesh_,
174  ),
175  mesh_,
177  )
178  );
179 
180  PAPtr->store();
181 
182  volScalarField* LPPtr
183  (
184  new volScalarField
185  (
186  IOobject
187  (
188  scopedName("L_P"),
189  mesh_.time().timeName(),
190  mesh_,
193  ),
194  mesh_,
196  )
197  );
198 
199  LPPtr->store();
200 }
201 
202 
203 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
204 
206 {
208  {
209  Info<< type() << " " << name() << nl;
210 
211  dict.readIfPresent("alphaEps", alphaEps_);
212  rhoInf_.readIfPresent("rhoInf", dict);
213  aRef_.readIfPresent("aRef", dict);
214 
215  if (dict.readIfPresent("k", kName_))
216  {
217  Info<< " k field: " << kName_ << endl;
218  }
219  else
220  {
221  Info<< " k field from turbulence model" << endl;
222  }
223 
224  if (dict.readIfPresent("epsilon", epsilonName_))
225  {
226  Info<< " epsilon field: " << epsilonName_ << endl;
227  }
228  else
229  {
230  Info<< " epsilon field from turbulence model (if needed)"
231  << endl;
232  }
233 
234  if (dict.readIfPresent("omega", omegaName_))
235  {
236  Info<< " omega field: " << omegaName_ << endl;
237  }
238  else
239  {
240  Info<< " omega field from turbulence model (if needed)" << endl;
241  }
242 
243  if (epsilonName_ != "none" && omegaName_ != "none")
244  {
246  << "either epsilon or omega field names can be set but not both"
247  << exit(FatalIOError);
248  }
249 
250  Info<< endl;
251 
252  return true;
253  }
254 
255  return false;
256 }
257 
258 
260 {
261  const volScalarField Mt(sqrt(2*k())/a());
262 
263  auto& P_A = mesh_.lookupObjectRef<volScalarField>(scopedName("P_A"));
264 
265  P_A = rhoScale(alphaEps_*epsilon()*pow5(Mt));
266 
267  auto& L_P = mesh_.lookupObjectRef<volScalarField>(scopedName("L_P"));
268 
269  L_P = 10*log10(P_A/dimensionedScalar("PRef", dimPower/dimVolume, 1e-12));
270 
271  return true;
272 }
273 
274 
276 {
277  Log << type() << " " << name() << " write:" << nl;
278 
279  const auto& P_A = mesh_.lookupObject<volScalarField>(scopedName("P_A"));
280 
281  Log << " writing field " << P_A.name() << nl;
282 
283  P_A.write();
284 
285  const auto& L_P = mesh_.lookupObject<volScalarField>(scopedName("L_P"));
286 
287  Log << " writing field " << L_P.name() << nl;
288 
289  L_P.write();
290 
291  Log << endl;
292 
293  return true;
294 }
295 
296 
297 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
volFields.H
runTime
engineTime & runTime
Definition: createEngineTime.H:13
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::functionObjects::proudmanAcousticPower::execute
virtual bool execute()
Calculate the Proudman acoustic power.
Definition: proudmanAcousticPower.C:259
basicThermo.H
Log
#define Log
Definition: PDRblock.C:35
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::dimVelocity
const dimensionSet dimVelocity
Foam::dimDensity
const dimensionSet dimDensity
Foam::read
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
Foam::turbulenceModel::propertiesName
static const word propertiesName
Default name of the turbulence properties dictionary.
Definition: turbulenceModel.H:100
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Foam::FatalIOError
IOerror FatalIOError
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::functionObjects::fvMeshFunctionObject
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Definition: fvMeshFunctionObject.H:64
Foam::functionObjects::proudmanAcousticPower::write
virtual bool write()
Write the Proudman acoustic power.
Definition: proudmanAcousticPower.C:275
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::log10
dimensionedScalar log10(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:263
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
fld
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))
Definition: gmvOutputLagrangian.H:23
Foam::dimPower
const dimensionSet dimPower
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::pow5
dimensionedScalar pow5(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:111
Foam::functionObjects::regionFunctionObject::read
virtual bool read(const dictionary &dict)
Read optional controls.
Definition: regionFunctionObject.C:173
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::functionObject::name
const word & name() const noexcept
Return the name of this functionObject.
Definition: functionObject.C:143
Foam::functionObjects::proudmanAcousticPower::proudmanAcousticPower
proudmanAcousticPower(const word &name, const Time &runTime, const dictionary &)
Construct from Time and dictionary.
Definition: proudmanAcousticPower.C:147
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::functionObjects::addToRunTimeSelectionTable
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
Foam::functionObject::type
virtual const word & type() const =0
Runtime type information.
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:590
k
label k
Boltzmann constant.
Definition: LISASMDCalcMethod2.H:41
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
proudmanAcousticPower.H
Foam::tmp::New
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
Foam::functionObjects::defineTypeNameAndDebug
defineTypeNameAndDebug(ObukhovLength, 0)
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
FatalIOErrorInFunction
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
epsilon
scalar epsilon
Definition: evaluateNearWall.H:7
Foam::dimVolume
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
Foam::GeometricField< scalar, fvPatchField, volMesh >
turb
compressible::turbulenceModel & turb
Definition: setRegionFluidFields.H:10
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::functionObjects::proudmanAcousticPower::read
virtual bool read(const dictionary &)
Read the Proudman acoustic power data.
Definition: proudmanAcousticPower.C:205
turbulenceModel.H
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:189
Foam::dictionary::dictName
word dictName() const
The local dictionary name (final part of scoped name)
Definition: dictionaryI.H:60