turbulenceIntensity.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2018-2022 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "turbulenceIntensity.H"
27 #include "momentumTransportModel.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace functionObjects
35 {
37 
39  (
43  );
44 }
45 }
46 
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
50 void Foam::functionObjects::turbulenceIntensity::writeFileHeader(const label i)
51 {
52  writeHeader(file(), "I ()");
53  writeCommented(file(), "Time");
54  writeTabbed(file(), "min");
55  writeTabbed(file(), "max");
56  writeTabbed(file(), "average");
57  file() << endl;
58 }
59 
60 
61 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
62 
64 (
65  const word& name,
66  const Time& runTime,
67  const dictionary& dict
68 )
69 :
70  fvMeshFunctionObject(name, runTime, dict),
71  logFiles(obr_, name),
72  writeLocalObjects(obr_, log)
73 {
74  read(dict);
75 }
76 
77 
78 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
79 
81 {}
82 
83 
84 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
85 
87 {
90 
91  resetName(typeName);
92  resetLocalObjectName("I");
93 
94  return true;
95 }
96 
97 
99 {
100  if (mesh_.foundType<momentumTransportModel>())
101  {
102  const momentumTransportModel& turbModel =
103  mesh_.lookupType<momentumTransportModel>();
104 
105  volScalarField uPrime(sqrt((2.0/3.0)*turbModel.k()));
106 
107  word name("I");
108 
109  return
110  store
111  (
112  name,
113  uPrime
114  /max
115  (
116  max(uPrime, mag(turbModel.U())),
118  )
119  );
120  }
121  else
122  {
124  << "Unable to find turbulence model in the "
125  << "database" << exit(FatalError);
126 
127  return false;
128  }
129 }
130 
131 
133 {
134  Log << type() << " " << name() << " write:" << nl;
135 
137 
138  logFiles::write();
139 
141  mesh_.lookupObject<volScalarField>("I");
142 
143  const scalar minTurbulenceIntensity = min(turbulenceIntensity).value();
144  const scalar maxTurbulenceIntensity = max(turbulenceIntensity).value();
145  const scalar avgTurbulenceIntensity = turbulenceIntensity.average().value();
146 
147  if (Pstream::master())
148  {
149  Log << " I : min = " << minTurbulenceIntensity
150  << ", max = " << maxTurbulenceIntensity
151  << ", average = " << avgTurbulenceIntensity << nl;
152 
153  writeTime(file());
154  file()
155  << tab << minTurbulenceIntensity
156  << tab << maxTurbulenceIntensity
157  << tab << avgTurbulenceIntensity
158  << endl;
159  }
160 
161  Log << endl;
162 
163  return true;
164 }
165 
166 
167 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base-class for Time/database functionObjects.
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
functionObject base class for creating, maintaining and writing log files e.g. integrated of averaged...
Definition: logFiles.H:60
OFstream & file()
Return access to the file (if only 1)
Definition: logFiles.C:113
virtual bool write()
Write function.
Definition: logFiles.C:173
Calculates the natural logarithm of the specified scalar field.
Definition: log.H:106
const ObjectType & lookupObject(const word &fieldName) const
Lookup object from the objectRegistry.
virtual bool read(const dictionary &)
Read optional controls.
Evaluates and writes the turbulence intensity field 'I'.
virtual bool execute()
Calculate the turbulenceIntensity field.
virtual bool write()
Write the turbulenceIntensity field.
turbulenceIntensity(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
virtual bool read(const dictionary &)
Read the turbulenceIntensity data.
void writeTabbed(Ostream &os, const string &str) const
Write a tabbed string to stream.
Definition: writeFile.C:121
void writeHeader(Ostream &os, const string &str) const
Write a commented header to stream.
Definition: writeFile.C:131
void writeCommented(Ostream &os, const string &str) const
Write a commented string to stream.
Definition: writeFile.C:110
FunctionObject base class for managing a list of objects on behalf of the inheriting function object,...
virtual bool read(const dictionary &)
Read the list of objects to be written.
virtual bool write()
Write function.
Abstract base class for turbulence models (RAS, LES and laminar).
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
const volVectorField & U() const
Access function to velocity field.
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
#define Log
Report write to Foam::Info if the local log switch is true.
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
static const char tab
Definition: Ostream.H:265
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimVelocity
error FatalError
static const char nl
Definition: Ostream.H:266
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict