basicThermo.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) 2011-2016 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 "basicThermo.H"
31 #include "fixedJumpFvPatchFields.H"
35 
36 
37 /* * * * * * * * * * * * * * * private static data * * * * * * * * * * * * * */
38 
39 namespace Foam
40 {
41  defineTypeNameAndDebug(basicThermo, 0);
42  defineRunTimeSelectionTable(basicThermo, fvMesh);
43 }
44 
45 const Foam::word Foam::basicThermo::dictName("thermophysicalProperties");
46 
47 
48 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
49 
51 {
52  const volScalarField::Boundary& tbf =
53  this->T_.boundaryField();
54 
55  wordList hbt(tbf.size(), word::null);
56 
57  forAll(tbf, patchi)
58  {
59  if (isA<fixedJumpFvPatchScalarField>(tbf[patchi]))
60  {
61  const fixedJumpFvPatchScalarField& pf =
62  dynamic_cast<const fixedJumpFvPatchScalarField&>(tbf[patchi]);
63 
64  hbt[patchi] = pf.interfaceFieldType();
65  }
66  else if (isA<fixedJumpAMIFvPatchScalarField>(tbf[patchi]))
67  {
68  const fixedJumpAMIFvPatchScalarField& pf =
69  dynamic_cast<const fixedJumpAMIFvPatchScalarField&>
70  (
71  tbf[patchi]
72  );
73 
74  hbt[patchi] = pf.interfaceFieldType();
75  }
76  }
77 
78  return hbt;
79 }
80 
81 
83 {
84  const volScalarField::Boundary& tbf =
85  this->T_.boundaryField();
86 
87  wordList hbt = tbf.types();
88 
89  forAll(tbf, patchi)
90  {
91  if (isA<fixedValueFvPatchScalarField>(tbf[patchi]))
92  {
93  hbt[patchi] = fixedEnergyFvPatchScalarField::typeName;
94  }
95  else if
96  (
97  isA<zeroGradientFvPatchScalarField>(tbf[patchi])
98  || isA<fixedGradientFvPatchScalarField>(tbf[patchi])
99  )
100  {
101  hbt[patchi] = gradientEnergyFvPatchScalarField::typeName;
102  }
103  else if (isA<mixedFvPatchScalarField>(tbf[patchi]))
104  {
105  hbt[patchi] = mixedEnergyFvPatchScalarField::typeName;
106  }
107  else if (isA<fixedJumpFvPatchScalarField>(tbf[patchi]))
108  {
110  }
111  else if (isA<fixedJumpAMIFvPatchScalarField>(tbf[patchi]))
112  {
114  }
115  else if (tbf[patchi].type() == "energyRegionCoupledFvPatchScalarField")
116  {
117  hbt[patchi] = "energyRegionCoupledFvPatchScalarField";
118  }
119  }
120 
121  return hbt;
122 }
123 
124 
125 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
126 
128 (
129  const fvMesh& mesh,
130  const char* name
131 ) const
132 {
133  if (!mesh.objectRegistry::foundObject<volScalarField>(name))
134  {
135  volScalarField* fPtr
136  (
137  new volScalarField
138  (
139  IOobject
140  (
141  name,
142  mesh.time().timeName(),
143  mesh,
146  ),
147  mesh
148  )
149  );
150 
151  // Transfer ownership of this object to the objectRegistry
152  fPtr->store(fPtr);
153  }
154 
155  return const_cast<volScalarField&>
156  (
157  mesh.objectRegistry::lookupObject<volScalarField>(name)
158  );
159 }
160 
161 
163 (
164  const fvMesh& mesh,
165  const word& phaseName
166 )
167 :
169  (
170  IOobject
171  (
172  phasePropertyName(dictName, phaseName),
173  mesh.time().constant(),
174  mesh,
177  )
178  ),
179 
180  phaseName_(phaseName),
181 
182  p_(lookupOrConstruct(mesh, "p")),
183 
184  T_
185  (
186  IOobject
187  (
188  phasePropertyName("T"),
189  mesh.time().timeName(),
190  mesh,
193  ),
194  mesh
195  ),
196 
197  alpha_
198  (
199  IOobject
200  (
201  phasePropertyName("thermo:alpha"),
202  mesh.time().timeName(),
203  mesh,
206  ),
207  mesh,
208  dimensionSet(1, -1, -1, 0, 0)
209  ),
210 
211  dpdt_(lookupOrDefault<Switch>("dpdt", true))
212 {}
213 
214 
216 (
217  const fvMesh& mesh,
218  const dictionary& dict,
219  const word& phaseName
220 )
221 :
223  (
224  IOobject
225  (
226  phasePropertyName(dictName, phaseName),
227  mesh.time().constant(),
228  mesh,
231  ),
232  dict
233  ),
234 
235  phaseName_(phaseName),
236 
237  p_(lookupOrConstruct(mesh, "p")),
238 
239  T_
240  (
241  IOobject
242  (
243  phasePropertyName("T"),
244  mesh.time().timeName(),
245  mesh,
248  ),
249  mesh
250  ),
251 
252  alpha_
253  (
254  IOobject
255  (
256  phasePropertyName("thermo:alpha"),
257  mesh.time().timeName(),
258  mesh,
261  ),
262  mesh,
263  dimensionSet(1, -1, -1, 0, 0)
264  )
265 {}
266 
267 
268 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
269 
271 (
272  const fvMesh& mesh,
273  const word& phaseName
274 )
275 {
276  return New<basicThermo>(mesh, phaseName);
277 }
278 
279 
280 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
281 
283 {}
284 
285 
286 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
287 
289 (
290  const fvPatchScalarField& pf
291 )
292 {
293  if (pf.db().foundObject<basicThermo>(dictName))
294  {
295  return pf.db().lookupObject<basicThermo>(dictName);
296  }
297  else
298  {
300  pf.db().lookupClass<basicThermo>();
301 
302  for
303  (
305  iter != thermos.end();
306  ++iter
307  )
308  {
309  if
310  (
311  &(iter()->he().internalField())
312  == &(pf.internalField())
313  )
314  {
315  return *iter();
316  }
317  }
318  }
319 
320  return pf.db().lookupObject<basicThermo>(dictName);
321 }
322 
323 
325 (
326  const string& app,
327  const word& a
328 ) const
329 {
330  if (!(he().name() == phasePropertyName(a)))
331  {
333  << "Supported energy type is " << phasePropertyName(a)
334  << ", thermodynamics package provides " << he().name()
335  << exit(FatalError);
336  }
337 }
338 
340 (
341  const string& app,
342  const word& a,
343  const word& b
344 ) const
345 {
346  if
347  (
348  !(
349  he().name() == phasePropertyName(a)
350  || he().name() == phasePropertyName(b)
351  )
352  )
353  {
355  << "Supported energy types are " << phasePropertyName(a)
356  << " and " << phasePropertyName(b)
357  << ", thermodynamics package provides " << he().name()
358  << exit(FatalError);
359  }
360 }
361 
363 (
364  const string& app,
365  const word& a,
366  const word& b,
367  const word& c
368 ) const
369 {
370  if
371  (
372  !(
373  he().name() == phasePropertyName(a)
374  || he().name() == phasePropertyName(b)
375  || he().name() == phasePropertyName(c)
376  )
377  )
378  {
380  << "Supported energy types are " << phasePropertyName(a)
381  << ", " << phasePropertyName(b)
382  << " and " << phasePropertyName(c)
383  << ", thermodynamics package provides " << he().name()
384  << exit(FatalError);
385  }
386 }
387 
389 (
390  const string& app,
391  const word& a,
392  const word& b,
393  const word& c,
394  const word& d
395 ) const
396 {
397  if
398  (
399  !(
400  he().name() == phasePropertyName(a)
401  || he().name() == phasePropertyName(b)
402  || he().name() == phasePropertyName(c)
403  || he().name() == phasePropertyName(d)
404  )
405  )
406  {
408  << "Supported energy types are " << phasePropertyName(a)
409  << ", " << phasePropertyName(b)
410  << ", " << phasePropertyName(c)
411  << " and " << phasePropertyName(d)
412  << ", thermodynamics package provides " << he().name()
413  << exit(FatalError);
414  }
415 }
416 
417 
419 (
420  const word& thermoName,
421  const int nCmpt
422 )
423 {
424  wordList cmpts(nCmpt);
425 
426  string::size_type beg=0, end=0, endb=0, endc=0;
427  int i = 0;
428 
429  while
430  (
431  (endb = thermoName.find('<', beg)) != string::npos
432  || (endc = thermoName.find(',', beg)) != string::npos
433  )
434  {
435  if (endb == string::npos)
436  {
437  end = endc;
438  }
439  else if ((endc = thermoName.find(',', beg)) != string::npos)
440  {
441  end = min(endb, endc);
442  }
443  else
444  {
445  end = endb;
446  }
447 
448  if (beg < end)
449  {
450  cmpts[i] = thermoName.substr(beg, end-beg);
451  cmpts[i++].replaceAll(">","");
452  }
453  beg = end + 1;
454  }
455 
456  if (beg < thermoName.size())
457  {
458  cmpts[i] = thermoName.substr(beg, string::npos);
459  cmpts[i++].replaceAll(">","");
460  }
461 
462  return cmpts;
463 }
464 
465 
467 {
468  return p_;
469 }
470 
471 
473 {
474  return p_;
475 }
476 
477 
479 {
480  return T_;
481 }
482 
483 
485 {
486  return T_;
487 }
488 
489 
491 {
492  return alpha_;
493 }
494 
495 
497 {
498  return alpha_.boundaryField()[patchi];
499 }
500 
501 
503 {
504  return regIOobject::read();
505 }
506 
507 
508 // ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:56
virtual ~basicThermo()
Destructor.
Definition: basicThermo.C:282
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
virtual bool read()
Read object.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
bool foundObject(const word &name) const
Is the named Type found?
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static const char *const typeName
Definition: Field.H:94
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
static autoPtr< Thermo > New(const fvMesh &, const word &phaseName=word::null)
Generic New for each of the related thermodynamics packages.
static Table::iterator lookupThermo(const dictionary &thermoDict, Table *tablePtr)
Generic lookup for each of the related thermodynamics packages.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
virtual volScalarField & p()
Pressure [Pa].
Definition: basicThermo.C:466
An STL-conforming iterator.
Definition: HashTable.H:415
void validate(const string &app, const word &) const
Check that the thermodynamics package is consistent.
Definition: basicThermo.C:325
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Dimension set for the base types.
Definition: dimensionSet.H:118
dynamicFvMesh & mesh
wordList heBoundaryBaseTypes()
Return the enthalpy/internal energy field boundary base types.
Definition: basicThermo.C:50
const Boundary & boundaryField() const
Return const-reference to the boundary field.
static wordList splitThermoName(const word &thermoName, const int nCmpt)
Split name of thermo package into a list of the components names.
Definition: basicThermo.C:419
A class for handling words, derived from string.
Definition: word.H:59
PtrList< solidThermo > thermos(solidRegions.size())
wordList heBoundaryTypes()
Return the enthalpy/internal energy field boundary types.
Definition: basicThermo.C:82
static const word null
An empty word.
Definition: word.H:77
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
iterator begin()
Iterator set to the beginning of the HashTable.
Definition: HashTableI.H:419
An STL-conforming hash table.
Definition: HashTable.H:61
graph_traits< Graph >::vertices_size_type size_type
Definition: SloanRenumber.C:73
void store()
Transfer ownership of this object to its registry.
Definition: regIOobjectI.H:34
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:478
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:198
label patchi
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
word dictName("noiseDict")
virtual const volScalarField & alpha() const
Thermal diffusivity for enthalpy of mixture [kg/m/s].
Definition: basicThermo.C:490
basicThermo(const basicThermo &)
Construct as copy (not implemented)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual bool read()
Read thermophysical properties dictionary.
Definition: basicThermo.C:502
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
volScalarField & lookupOrConstruct(const fvMesh &mesh, const char *name) const
Definition: basicThermo.C:128
HashTable< const Type * > lookupClass(const bool strict=false) const
Lookup and return all objects of the given Type.
const word & name() const
Return name.
Definition: IOobject.H:260
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:346
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243