turbulenceFields.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 | Copyright (C) 2013-2015 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"
27 #include "dictionary.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(turbulenceFields, 0);
36 
37  template<>
39  {
40  "k",
41  "epsilon",
42  "mut",
43  "muEff",
44  "alphat",
45  "alphaEff",
46  "R",
47  "devRhoReff"
48  };
51 
52  template<>
54  {
55  "k",
56  "epsilon",
57  "nut",
58  "nuEff",
59  "R",
60  "devReff"
61  };
64 
66 }
67 
68 
69 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
70 
72 {
73  if (obr_.foundObject<compressible::turbulenceModel>(modelName))
74  {
75  return true;
76  }
77  else if (obr_.foundObject<incompressible::turbulenceModel>(modelName))
78  {
79  return false;
80  }
81  else
82  {
83  WarningIn("Foam::word& Foam::turbulenceFields::compressible() const")
84  << "Turbulence model not found in database, deactivating";
85  active_ = false;
86  }
87 
88  return false;
89 }
90 
91 
92 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
93 
95 (
96  const word& name,
97  const objectRegistry& obr,
98  const dictionary& dict,
99  const bool loadFromFiles
100 )
101 :
102  name_(name),
103  obr_(obr),
104  active_(true),
105  fieldSet_()
106 {
107  // Check if the available mesh is an fvMesh otherise deactivate
108  if (isA<fvMesh>(obr_))
109  {
110  read(dict);
111  }
112  else
113  {
114  active_ = false;
115  WarningIn
116  (
117  "turbulenceFields::turbulenceFields"
118  "("
119  "const word&, "
120  "const objectRegistry&, "
121  "const dictionary&, "
122  "const bool"
123  ")"
124  ) << "No fvMesh available, deactivating " << name_
125  << endl;
126  }
127 }
128 
129 
130 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
131 
133 {}
134 
135 
136 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
137 
139 {
140  if (active_)
141  {
142  fieldSet_.insert(wordList(dict.lookup("fields")));
143 
144  Info<< type() << " " << name_ << ": ";
145  if (fieldSet_.size())
146  {
147  Info<< "storing fields:" << nl;
148  forAllConstIter(wordHashSet, fieldSet_, iter)
149  {
150  Info<< " " << modelName << ':' << iter.key() << nl;
151  }
152  Info<< endl;
153  }
154  else
155  {
156  Info<< "no fields requested to be stored" << nl << endl;
157  }
158  }
159 }
160 
161 
163 {
164  bool comp = compressible();
165 
166  if (!active_)
167  {
168  return;
169  }
170 
171  if (comp)
172  {
173  const compressible::turbulenceModel& model =
174  obr_.lookupObject<compressible::turbulenceModel>(modelName);
175 
176  forAllConstIter(wordHashSet, fieldSet_, iter)
177  {
178  const word& f = iter.key();
179  switch (compressibleFieldNames_[f])
180  {
181  case cfK:
182  {
183  processField<scalar>(f, model.k());
184  break;
185  }
186  case cfEpsilon:
187  {
188  processField<scalar>(f, model.epsilon());
189  break;
190  }
191  case cfMut:
192  {
193  processField<scalar>(f, model.mut());
194  break;
195  }
196  case cfMuEff:
197  {
198  processField<scalar>(f, model.muEff());
199  break;
200  }
201  case cfAlphat:
202  {
203  processField<scalar>(f, model.alphat());
204  break;
205  }
206  case cfAlphaEff:
207  {
208  processField<scalar>(f, model.alphaEff());
209  break;
210  }
211  case cfR:
212  {
213  processField<symmTensor>(f, model.R());
214  break;
215  }
216  case cfDevRhoReff:
217  {
218  processField<symmTensor>(f, model.devRhoReff());
219  break;
220  }
221  default:
222  {
223  FatalErrorIn("void Foam::turbulenceFields::execute()")
224  << "Invalid field selection" << abort(FatalError);
225  }
226  }
227  }
228  }
229  else
230  {
231  const incompressible::turbulenceModel& model =
232  obr_.lookupObject<incompressible::turbulenceModel>(modelName);
233 
234  forAllConstIter(wordHashSet, fieldSet_, iter)
235  {
236  const word& f = iter.key();
237  switch (incompressibleFieldNames_[f])
238  {
239  case ifK:
240  {
241  processField<scalar>(f, model.k());
242  break;
243  }
244  case ifEpsilon:
245  {
246  processField<scalar>(f, model.epsilon());
247  break;
248  }
249  case ifNut:
250  {
251  processField<scalar>(f, model.nut());
252  break;
253  }
254  case ifNuEff:
255  {
256  processField<scalar>(f, model.nuEff());
257  break;
258  }
259  case ifR:
260  {
261  processField<symmTensor>(f, model.R());
262  break;
263  }
264  case ifDevReff:
265  {
266  processField<symmTensor>(f, model.devReff());
267  break;
268  }
269  default:
270  {
271  FatalErrorIn("void Foam::turbulenceFields::execute()")
272  << "Invalid field selection" << abort(FatalError);
273  }
274  }
275  }
276  }
277 }
278 
279 
281 {
282  if (active_)
283  {
284  execute();
285  }
286 }
287 
288 
290 {}
291 
292 
294 {}
295 
296 
297 // ************************************************************************* //
virtual void end()
Execute at the final time-loop, currently does nothing.
static const NamedEnum< incompressibleField, 6 > incompressibleFieldNames_
static const NamedEnum< compressibleField, 8 > compressibleFieldNames_
labelList f(nPoints)
Templated abstract base class for single-phase incompressible turbulence models.
static const word modelName
virtual void timeSet()
Called when time was set at the end of the Time::operator++.
A class for handling words, derived from string.
Definition: word.H:59
virtual void write()
Write.
messageStream Info
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Namespace for OpenFOAM.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define WarningIn(functionName)
Report a warning using Foam::Warning.
virtual tmp< volScalarField > alphaEff() const
Return the effective turbulent thermal diffusivity for enthalpy.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:52
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual tmp< volSymmTensorField > devReff() const
Return the effective stress tensor.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
bool compressible
Definition: pEqn.H:40
List< word > wordList
A List of words.
Definition: fileName.H:54
error FatalError
virtual void read(const dictionary &)
Read the field min/max data.
Registry of regIOobjects.
A HashTable with keys but without contents.
Definition: HashSet.H:59
static const word propertiesName
Default name of the turbulence properties dictionary.
bool read(const char *, int32_t &)
Definition: int32IO.C:87
bool compressible()
Return true if compressible turbulence model is identified.
virtual void execute()
Execute, currently does nothing.
turbulenceFields(const turbulenceFields &)
Disallow default bitwise copy construct.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual ~turbulenceFields()
Destructor.
defineTypeNameAndDebug(combustionModel, 0)
virtual tmp< volScalarField > alphat() const
Return the turbulent thermal diffusivity for enthalpy [kg/m/s].