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  phaseName_(dict.lookupOrDefault<word>("phase", word::null))
106 {
107  read(dict);
108 }
109 
110 
111 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
112 
114 {}
115 
116 
117 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118 
120 {
121  if (dict.found("field"))
122  {
123  fieldSet_.insert(word(dict.lookup("field")));
124  }
125  else
126  {
127  fieldSet_.insert(wordList(dict.lookup("fields")));
128  }
129 
130  if (dict.lookupOrDefault<Switch>("prefix", false))
131  {
132  prefix_ = momentumTransportModel::typeName + ':';
133  }
134  else
135  {
136  prefix_ = word::null;
137  }
138 
139  Info<< type() << " " << name() << ": ";
140  if (fieldSet_.size())
141  {
142  Info<< "storing fields:" << nl;
143  forAllConstIter(wordHashSet, fieldSet_, iter)
144  {
145  Info<< " "
146  << IOobject::groupName(prefix_ + iter.key(), phaseName_) << nl;
147  }
148  Info<< endl;
149  }
150  else
151  {
152  Info<< "no fields requested to be stored" << nl << endl;
153  }
154 
155  return true;
156 }
157 
158 
160 {
161  const word modelName
162  (
163  IOobject::groupName(momentumTransportModel::typeName, phaseName_)
164  );
165 
166  if (obr_.foundObject<thermophysicalTransportModel>(modelName))
167  {
168  const thermophysicalTransportModel& ttm =
169  obr_.lookupObject<thermophysicalTransportModel>
170  (
171  thermophysicalTransportModel::typeName
172  );
173 
175  ttm.momentumTransport();
176 
177  forAllConstIter(wordHashSet, fieldSet_, iter)
178  {
179  const word& f = iter.key();
180  switch (compressibleFieldNames_[f])
181  {
183  {
184  processField<scalar>(f, model.k());
185  break;
186  }
188  {
189  processField<scalar>(f, model.epsilon());
190  break;
191  }
192  case compressibleField::omega:
193  {
194  processField<scalar>(f, omega(model));
195  break;
196  }
197  case compressibleField::mut:
198  {
199  processField<scalar>(f, model.mut());
200  break;
201  }
202  case compressibleField::muEff:
203  {
204  processField<scalar>(f, model.muEff());
205  break;
206  }
207  case compressibleField::alphaEff:
208  {
209  processField<scalar>(f, ttm.alphaEff());
210  break;
211  }
213  {
214  processField<symmTensor>(f, model.sigma());
215  break;
216  }
217  case compressibleField::devTau:
218  {
219  processField<symmTensor>(f, model.devTau());
220  break;
221  }
222  default:
223  {
225  << "Invalid field selection" << exit(FatalError);
226  }
227  }
228  }
229  }
230  else if (obr_.foundObject<compressibleMomentumTransportModel>(modelName))
231  {
233  obr_.lookupObject<compressibleMomentumTransportModel>(modelName);
234 
235  forAllConstIter(wordHashSet, fieldSet_, iter)
236  {
237  const word& f = iter.key();
238  switch (compressibleFieldNames_[f])
239  {
241  {
242  processField<scalar>(f, model.k());
243  break;
244  }
246  {
247  processField<scalar>(f, model.epsilon());
248  break;
249  }
250  case compressibleField::omega:
251  {
252  processField<scalar>(f, omega(model));
253  break;
254  }
255  case compressibleField::mut:
256  {
257  processField<scalar>(f, model.mut());
258  break;
259  }
260  case compressibleField::muEff:
261  {
262  processField<scalar>(f, model.muEff());
263  break;
264  }
266  {
267  processField<symmTensor>(f, model.sigma());
268  break;
269  }
270  case compressibleField::devTau:
271  {
272  processField<symmTensor>(f, model.devTau());
273  break;
274  }
275  default:
276  {
278  << "Invalid field selection" << exit(FatalError);
279  }
280  }
281  }
282  }
283  else if
284  (
285  obr_.foundObject<incompressible::momentumTransportModel>(modelName)
286  )
287  {
289  obr_.lookupObject<incompressible::momentumTransportModel>
290  (
291  modelName
292  );
293 
294  forAllConstIter(wordHashSet, fieldSet_, iter)
295  {
296  const word& f = iter.key();
297  switch (incompressibleFieldNames_[f])
298  {
300  {
301  processField<scalar>(f, model.k());
302  break;
303  }
305  {
306  processField<scalar>(f, model.epsilon());
307  break;
308  }
309  case incompressibleField::omega:
310  {
311  processField<scalar>(f, omega(model));
312  break;
313  }
315  {
316  processField<scalar>(f, model.nut());
317  break;
318  }
319  case incompressibleField::nuEff:
320  {
321  processField<scalar>(f, model.nuEff());
322  break;
323  }
325  {
326  processField<symmTensor>(f, model.sigma());
327  break;
328  }
329  case incompressibleField::devSigma:
330  {
331  processField<symmTensor>(f, model.devSigma());
332  break;
333  }
334  default:
335  {
337  << "Invalid field selection" << exit(FatalError);
338  }
339  }
340  }
341  }
342  else
343  {
345  << "Turbulence model not found in database, deactivating"
346  << exit(FatalError);
347  }
348 
349  return true;
350 }
351 
352 
354 {
355  forAllConstIter(wordHashSet, fieldSet_, iter)
356  {
357  const word fieldName
358  (
359  IOobject::groupName(prefix_ + iter.key(), phaseName_)
360  );
361  writeObject(fieldName);
362  }
363 
364  return true;
365 }
366 
367 
368 // ************************************************************************* //
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:643
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:156
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:323
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 word groupName(Name name, const word &group)
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
Specialisation 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:844