turbulenceFields.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) 2013-2018 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 "turbulenceFields.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace functionObjects
36 {
37  defineTypeNameAndDebug(turbulenceFields, 0);
38 
40  (
41  functionObject,
42  turbulenceFields,
43  dictionary
44  );
45 }
46 }
47 
48 template<>
49 const char* Foam::NamedEnum
50 <
52  9
53 >::names[] =
54 {
55  "k",
56  "epsilon",
57  "omega",
58  "mut",
59  "muEff",
60  "alphat",
61  "alphaEff",
62  "R",
63  "devRhoReff"
64 };
65 
66 const Foam::NamedEnum
67 <
69  9
71 
72 template<>
73 const char* Foam::NamedEnum
74 <
76  7
77 >::names[] =
78 {
79  "k",
80  "epsilon",
81  "omega",
82  "nut",
83  "nuEff",
84  "R",
85  "devReff"
86 };
87 
88 const Foam::NamedEnum
89 <
91  7
93 
95 (
97 );
98 
99 
100 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
101 
103 {
104  if (obr_.foundObject<compressible::turbulenceModel>(modelName))
105  {
106  return true;
107  }
108  else if (obr_.foundObject<incompressible::turbulenceModel>(modelName))
109  {
110  return false;
111  }
112  else
113  {
115  << "Turbulence model not found in database, deactivating"
116  << exit(FatalError);
117  }
118 
119  return false;
120 }
121 
122 
123 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
124 
125 Foam::functionObjects::turbulenceFields::turbulenceFields
126 (
127  const word& name,
128  const Time& runTime,
129  const dictionary& dict
130 )
131 :
132  fvMeshFunctionObject(name, runTime, dict),
133  fieldSet_()
134 {
135  read(dict);
136 }
137 
138 
139 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
140 
142 {}
143 
144 
145 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
146 
148 {
149  if (dict.found("field"))
150  {
151  fieldSet_.insert(word(dict.lookup("field")));
152  }
153  else
154  {
155  fieldSet_.insert(wordList(dict.lookup("fields")));
156  }
157 
158  Info<< type() << " " << name() << ": ";
159  if (fieldSet_.size())
160  {
161  Info<< "storing fields:" << nl;
162  forAllConstIter(wordHashSet, fieldSet_, iter)
163  {
164  Info<< " " << modelName << ':' << iter.key() << nl;
165  }
166  Info<< endl;
167  }
168  else
169  {
170  Info<< "no fields requested to be stored" << nl << endl;
171  }
172 
173  return true;
174 }
175 
176 
178 {
179  bool comp = compressible();
180 
181  if (comp)
182  {
183  const compressible::turbulenceModel& model =
184  obr_.lookupObject<compressible::turbulenceModel>(modelName);
185 
186  forAllConstIter(wordHashSet, fieldSet_, iter)
187  {
188  const word& f = iter.key();
189  switch (compressibleFieldNames_[f])
190  {
191  case cfK:
192  {
193  processField<scalar>(f, model.k());
194  break;
195  }
196  case cfEpsilon:
197  {
198  processField<scalar>(f, model.epsilon());
199  break;
200  }
201  case cfOmega:
202  {
203  processField<scalar>(f, omega(model));
204  break;
205  }
206  case cfMut:
207  {
208  processField<scalar>(f, model.mut());
209  break;
210  }
211  case cfMuEff:
212  {
213  processField<scalar>(f, model.muEff());
214  break;
215  }
216  case cfAlphat:
217  {
218  processField<scalar>(f, model.alphat());
219  break;
220  }
221  case cfAlphaEff:
222  {
223  processField<scalar>(f, model.alphaEff());
224  break;
225  }
226  case cfR:
227  {
228  processField<symmTensor>(f, model.R());
229  break;
230  }
231  case cfDevRhoReff:
232  {
233  processField<symmTensor>(f, model.devRhoReff());
234  break;
235  }
236  default:
237  {
239  << "Invalid field selection" << abort(FatalError);
240  }
241  }
242  }
243  }
244  else
245  {
246  const incompressible::turbulenceModel& model =
247  obr_.lookupObject<incompressible::turbulenceModel>(modelName);
248 
249  forAllConstIter(wordHashSet, fieldSet_, iter)
250  {
251  const word& f = iter.key();
252  switch (incompressibleFieldNames_[f])
253  {
254  case ifK:
255  {
256  processField<scalar>(f, model.k());
257  break;
258  }
259  case ifEpsilon:
260  {
261  processField<scalar>(f, model.epsilon());
262  break;
263  }
264  case ifOmega:
265  {
266  processField<scalar>(f, omega(model));
267  break;
268  }
269  case ifNut:
270  {
271  processField<scalar>(f, model.nut());
272  break;
273  }
274  case ifNuEff:
275  {
276  processField<scalar>(f, model.nuEff());
277  break;
278  }
279  case ifR:
280  {
281  processField<symmTensor>(f, model.R());
282  break;
283  }
284  case ifDevReff:
285  {
286  processField<symmTensor>(f, model.devReff());
287  break;
288  }
289  default:
290  {
292  << "Invalid field selection" << abort(FatalError);
293  }
294  }
295  }
296  }
297 
298  return true;
299 }
300 
301 
303 {
304  forAllConstIter(wordHashSet, fieldSet_, iter)
305  {
306  const word fieldName = modelName + ':' + iter.key();
307  writeObject(fieldName);
308  }
309 
310  return true;
311 }
312 
313 
314 // ************************************************************************* //
A HashTable with keys but without contents.
Definition: HashSet.H:59
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:431
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
virtual tmp< volSymmTensorField > R() const =0
Return the Reynolds stress tensor.
static const NamedEnum< incompressibleField, 7 > incompressibleFieldNames_
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Macros for easy insertion into run-time selection tables.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
bool compressible()
Return true if compressible turbulence model is identified.
virtual tmp< volScalarField > alphaEff() const
Effective thermal turbulent diffusivity of mixture [kg/m/s].
Templated abstract base class for single-phase incompressible turbulence models.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
virtual bool read(const dictionary &)
Read the controls.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
virtual tmp< volSymmTensorField > devReff() const
Return the effective stress tensor.
static const char nl
Definition: Ostream.H:265
virtual tmp< volScalarField > nuEff() const =0
Return the effective viscosity.
labelList f(nPoints)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
List< word > wordList
A List of words.
Definition: fileName.H:54
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
virtual bool execute()
Calculate turbulence fields.
messageStream Info
virtual tmp< volScalarField > alphat() const
Turbulent thermal diffusivity for enthalpy [kg/m/s].
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
static const NamedEnum< compressibleField, 9 > compressibleFieldNames_
addToRunTimeSelectionTable(functionObject, add, dictionary)
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576