All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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-2020 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  8
53 >::names[] =
54 {
55  "k",
56  "epsilon",
57  "omega",
58  "mut",
59  "muEff",
60  "alphaEff",
61  "R",
62  "devTau"
63 };
64 
65 const Foam::NamedEnum
66 <
68  8
70 
71 template<>
72 const char* Foam::NamedEnum
73 <
75  7
76 >::names[] =
77 {
78  "k",
79  "epsilon",
80  "omega",
81  "nut",
82  "nuEff",
83  "R",
84  "devSigma"
85 };
86 
87 const Foam::NamedEnum
88 <
90  7
92 
93 
94 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
95 
97 (
98  const word& name,
99  const Time& runTime,
100  const dictionary& dict
101 )
102 :
103  fvMeshFunctionObject(name, runTime, dict),
104  fieldSet_()
105 {
106  read(dict);
107 }
108 
109 
110 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
111 
113 {}
114 
115 
116 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
117 
119 {
120  return Foam::momentumTransportModel::typeName;
121 }
122 
123 
125 {
126  if (dict.found("field"))
127  {
128  fieldSet_.insert(word(dict.lookup("field")));
129  }
130  else
131  {
132  fieldSet_.insert(wordList(dict.lookup("fields")));
133  }
134 
135  if (dict.lookupOrDefault<Switch>("prefix", false))
136  {
137  prefix_ = modelName() + ':';
138  }
139  else
140  {
141  prefix_ = word::null;
142  }
143 
144  Info<< type() << " " << name() << ": ";
145  if (fieldSet_.size())
146  {
147  Info<< "storing fields:" << nl;
148  forAllConstIter(wordHashSet, fieldSet_, iter)
149  {
150  Info<< " " << prefix_ + iter.key() << nl;
151  }
152  Info<< endl;
153  }
154  else
155  {
156  Info<< "no fields requested to be stored" << nl << endl;
157  }
158 
159  return true;
160 }
161 
162 
164 {
165  if (obr_.foundObject<thermophysicalTransportModel>(modelName()))
166  {
167  const thermophysicalTransportModel& ttm =
168  obr_.lookupObject<thermophysicalTransportModel>
169  (
170  thermophysicalTransportModel::typeName
171  );
172 
174  ttm.momentumTransport();
175 
176  forAllConstIter(wordHashSet, fieldSet_, iter)
177  {
178  const word& f = iter.key();
179  switch (compressibleFieldNames_[f])
180  {
182  {
183  processField<scalar>(f, model.k());
184  break;
185  }
187  {
188  processField<scalar>(f, model.epsilon());
189  break;
190  }
191  case compressibleField::omega:
192  {
193  processField<scalar>(f, omega(model));
194  break;
195  }
196  case compressibleField::mut:
197  {
198  processField<scalar>(f, model.mut());
199  break;
200  }
201  case compressibleField::muEff:
202  {
203  processField<scalar>(f, model.muEff());
204  break;
205  }
206  case compressibleField::alphaEff:
207  {
208  processField<scalar>(f, ttm.alphaEff());
209  break;
210  }
212  {
213  processField<symmTensor>(f, model.sigma());
214  break;
215  }
216  case compressibleField::devTau:
217  {
218  processField<symmTensor>(f, model.devTau());
219  break;
220  }
221  default:
222  {
224  << "Invalid field selection" << exit(FatalError);
225  }
226  }
227  }
228  }
229  else if (obr_.foundObject<compressibleMomentumTransportModel>(modelName()))
230  {
232  obr_.lookupObject<compressibleMomentumTransportModel>(modelName());
233 
234  forAllConstIter(wordHashSet, fieldSet_, iter)
235  {
236  const word& f = iter.key();
237  switch (compressibleFieldNames_[f])
238  {
240  {
241  processField<scalar>(f, model.k());
242  break;
243  }
245  {
246  processField<scalar>(f, model.epsilon());
247  break;
248  }
249  case compressibleField::omega:
250  {
251  processField<scalar>(f, omega(model));
252  break;
253  }
254  case compressibleField::mut:
255  {
256  processField<scalar>(f, model.mut());
257  break;
258  }
259  case compressibleField::muEff:
260  {
261  processField<scalar>(f, model.muEff());
262  break;
263  }
265  {
266  processField<symmTensor>(f, model.sigma());
267  break;
268  }
269  case compressibleField::devTau:
270  {
271  processField<symmTensor>(f, model.devTau());
272  break;
273  }
274  default:
275  {
277  << "Invalid field selection" << exit(FatalError);
278  }
279  }
280  }
281  }
282  else if
283  (
284  obr_.foundObject<incompressible::momentumTransportModel>(modelName())
285  )
286  {
288  obr_.lookupObject<incompressible::momentumTransportModel>
289  (
290  modelName()
291  );
292 
293  forAllConstIter(wordHashSet, fieldSet_, iter)
294  {
295  const word& f = iter.key();
296  switch (incompressibleFieldNames_[f])
297  {
299  {
300  processField<scalar>(f, model.k());
301  break;
302  }
304  {
305  processField<scalar>(f, model.epsilon());
306  break;
307  }
308  case incompressibleField::omega:
309  {
310  processField<scalar>(f, omega(model));
311  break;
312  }
314  {
315  processField<scalar>(f, model.nut());
316  break;
317  }
318  case incompressibleField::nuEff:
319  {
320  processField<scalar>(f, model.nuEff());
321  break;
322  }
324  {
325  processField<symmTensor>(f, model.sigma());
326  break;
327  }
328  case incompressibleField::devSigma:
329  {
330  processField<symmTensor>(f, model.devSigma());
331  break;
332  }
333  default:
334  {
336  << "Invalid field selection" << exit(FatalError);
337  }
338  }
339  }
340  }
341  else
342  {
344  << "Turbulence model not found in database, deactivating"
345  << exit(FatalError);
346  }
347 
348  return true;
349 }
350 
351 
353 {
354  forAllConstIter(wordHashSet, fieldSet_, iter)
355  {
356  const word fieldName = prefix_ + iter.key();
357  writeObject(fieldName);
358  }
359 
360  return true;
361 }
362 
363 
364 // ************************************************************************* //
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:667
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual tmp< volSymmTensorField > sigma() const =0
Return the stress tensor [m^2/s^2].
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
turbulenceFields(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const compressibleMomentumTransportModel & momentumTransport() const
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
virtual tmp< volScalarField > muEff() const =0
Return the effective dynamic viscosity.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
static const NamedEnum< incompressibleField, 7 > incompressibleFieldNames_
label k
Boltzmann constant.
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.
static const NamedEnum< compressibleField, 8 > compressibleFieldNames_
bool read(const char *, int32_t &)
Definition: int32IO.C:85
A class for handling words, derived from string.
Definition: word.H:59
static const word null
An empty word.
Definition: word.H:77
virtual bool read(const dictionary &)
Read the controls.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
static const char nl
Definition: Ostream.H:260
labelList f(nPoints)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
defineTypeNameAndDebug(Qdot, 0)
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
virtual tmp< volScalarField > mut() const =0
Return the turbulence dynamic viscosity.
Abstract base class for thermophysical transport models (RAS, LES and laminar).
List< word > wordList
A List of words.
Definition: fileName.H:54
Templated abstract base class for single-phase incompressible turbulence models.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define R(A, B, C, D, E, F, K, M)
virtual tmp< volScalarField > alphaEff() const =0
Effective thermal turbulent diffusivity of mixture [kg/m/s].
virtual tmp< volSymmTensorField > devTau() const =0
Return the effective stress tensor including the laminar stress.
scalar epsilon
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual bool execute()
Calculate turbulence fields.
messageStream Info
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
virtual tmp< volSymmTensorField > devSigma() const
Return the effective stress tensor.
virtual tmp< volScalarField > nuEff() const =0
Return the effective viscosity.
Abstract base class for turbulence models (RAS, LES and laminar).
scalar nut
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812