basicThermo.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) 2011-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 "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 mesh.objectRegistry::lookupObjectRef<volScalarField>(name);
156 }
157 
158 
160 (
161  const fvMesh& mesh,
162  const word& phaseName
163 )
164 :
166  (
167  IOobject
168  (
169  phasePropertyName(dictName, phaseName),
170  mesh.time().constant(),
171  mesh,
174  )
175  ),
176 
177  phaseName_(phaseName),
178 
179  p_(lookupOrConstruct(mesh, "p")),
180 
181  T_
182  (
183  IOobject
184  (
185  phasePropertyName("T"),
186  mesh.time().timeName(),
187  mesh,
190  ),
191  mesh
192  ),
193 
194  alpha_
195  (
196  IOobject
197  (
198  phasePropertyName("thermo:alpha"),
199  mesh.time().timeName(),
200  mesh,
203  ),
204  mesh,
205  dimensionedScalar(dimensionSet(1, -1, -1, 0, 0), Zero)
206  ),
207 
208  dpdt_(lookupOrDefault<Switch>("dpdt", true))
209 {}
210 
211 
213 (
214  const fvMesh& mesh,
215  const dictionary& dict,
216  const word& phaseName
217 )
218 :
220  (
221  IOobject
222  (
223  phasePropertyName(dictName, phaseName),
224  mesh.time().constant(),
225  mesh,
228  ),
229  dict
230  ),
231 
232  phaseName_(phaseName),
233 
234  p_(lookupOrConstruct(mesh, "p")),
235 
236  T_
237  (
238  IOobject
239  (
240  phasePropertyName("T"),
241  mesh.time().timeName(),
242  mesh,
245  ),
246  mesh
247  ),
248 
249  alpha_
250  (
251  IOobject
252  (
253  phasePropertyName("thermo:alpha"),
254  mesh.time().timeName(),
255  mesh,
258  ),
259  mesh,
260  dimensionedScalar(dimensionSet(1, -1, -1, 0, 0), Zero)
261  )
262 {}
263 
264 
265 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
266 
268 (
269  const fvMesh& mesh,
270  const word& phaseName
271 )
272 {
273  return New<basicThermo>(mesh, phaseName);
274 }
275 
276 
277 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
278 
280 {}
281 
282 
283 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
284 
286 (
287  const fvPatchScalarField& pf
288 )
289 {
290  if (pf.db().foundObject<basicThermo>(dictName))
291  {
292  return pf.db().lookupObject<basicThermo>(dictName);
293  }
294  else
295  {
297  pf.db().lookupClass<basicThermo>();
298 
299  for
300  (
302  iter != thermos.end();
303  ++iter
304  )
305  {
306  if
307  (
308  &(iter()->he().internalField())
309  == &(pf.internalField())
310  )
311  {
312  return *iter();
313  }
314  }
315  }
316 
317  return pf.db().lookupObject<basicThermo>(dictName);
318 }
319 
320 
322 (
323  const string& app,
324  const word& a
325 ) const
326 {
327  if (!(he().name() == phasePropertyName(a)))
328  {
330  << "Supported energy type is " << phasePropertyName(a)
331  << ", thermodynamics package provides " << he().name()
332  << exit(FatalError);
333  }
334 }
335 
337 (
338  const string& app,
339  const word& a,
340  const word& b
341 ) const
342 {
343  if
344  (
345  !(
346  he().name() == phasePropertyName(a)
347  || he().name() == phasePropertyName(b)
348  )
349  )
350  {
352  << "Supported energy types are " << phasePropertyName(a)
353  << " and " << phasePropertyName(b)
354  << ", thermodynamics package provides " << he().name()
355  << exit(FatalError);
356  }
357 }
358 
360 (
361  const string& app,
362  const word& a,
363  const word& b,
364  const word& c
365 ) const
366 {
367  if
368  (
369  !(
370  he().name() == phasePropertyName(a)
371  || he().name() == phasePropertyName(b)
372  || he().name() == phasePropertyName(c)
373  )
374  )
375  {
377  << "Supported energy types are " << phasePropertyName(a)
378  << ", " << phasePropertyName(b)
379  << " and " << phasePropertyName(c)
380  << ", thermodynamics package provides " << he().name()
381  << exit(FatalError);
382  }
383 }
384 
386 (
387  const string& app,
388  const word& a,
389  const word& b,
390  const word& c,
391  const word& d
392 ) const
393 {
394  if
395  (
396  !(
397  he().name() == phasePropertyName(a)
398  || he().name() == phasePropertyName(b)
399  || he().name() == phasePropertyName(c)
400  || he().name() == phasePropertyName(d)
401  )
402  )
403  {
405  << "Supported energy types are " << phasePropertyName(a)
406  << ", " << phasePropertyName(b)
407  << ", " << phasePropertyName(c)
408  << " and " << phasePropertyName(d)
409  << ", thermodynamics package provides " << he().name()
410  << exit(FatalError);
411  }
412 }
413 
414 
416 (
417  const word& thermoName,
418  const int nCmpt
419 )
420 {
421  wordList cmpts(nCmpt);
422 
423  string::size_type beg=0, end=0, endb=0, endc=0;
424  int i = 0;
425 
426  while
427  (
428  (endb = thermoName.find('<', beg)) != string::npos
429  || (endc = thermoName.find(',', beg)) != string::npos
430  )
431  {
432  if (endb == string::npos)
433  {
434  end = endc;
435  }
436  else if ((endc = thermoName.find(',', beg)) != string::npos)
437  {
438  end = min(endb, endc);
439  }
440  else
441  {
442  end = endb;
443  }
444 
445  if (beg < end)
446  {
447  cmpts[i] = thermoName.substr(beg, end-beg);
448  cmpts[i++].replaceAll(">","");
449 
450  // If the number of number of components in the name
451  // is greater than nCmpt return an empty list
452  if (i == nCmpt)
453  {
454  return wordList();
455  }
456  }
457  beg = end + 1;
458  }
459 
460  // If the number of number of components in the name is not equal to nCmpt
461  // return an empty list
462  if (i + 1 != nCmpt)
463  {
464  return wordList();
465  }
466 
467  if (beg < thermoName.size())
468  {
469  cmpts[i] = thermoName.substr(beg, string::npos);
470  cmpts[i].replaceAll(">","");
471  }
472 
473  return cmpts;
474 }
475 
476 
478 {
479  return p_;
480 }
481 
482 
484 {
485  return p_;
486 }
487 
488 
490 {
491  return T_;
492 }
493 
494 
496 {
497  return T_;
498 }
499 
500 
502 {
503  return alpha_;
504 }
505 
506 
508 {
509  return alpha_.boundaryField()[patchi];
510 }
511 
512 
514 {
515  return regIOobject::read();
516 }
517 
518 
519 // ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:50
virtual ~basicThermo()
Destructor.
Definition: basicThermo.C:279
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
const word & name() const
Return name.
Definition: IOobject.H:295
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:112
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:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static const char *const typeName
Definition: Field.H:105
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual const volScalarField & alpha() const
Thermal diffusivity for enthalpy of mixture [kg/m/s].
Definition: basicThermo.C:501
bool foundObject(const word &name) const
Is the named Type found?
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:489
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
static autoPtr< Thermo > New(const fvMesh &, const word &phaseName=word::null)
Generic New for each of the related thermodynamics packages.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:123
virtual volScalarField & p()
Pressure [Pa].
Definition: basicThermo.C:477
An STL-conforming iterator.
Definition: HashTable.H:426
Dimension set for the base types.
Definition: dimensionSet.H:120
dynamicFvMesh & mesh
wordList heBoundaryBaseTypes()
Return the enthalpy/internal energy field boundary base types.
Definition: basicThermo.C:50
static wordList splitThermoName(const word &thermoName, const int nCmpt)
Split name of thermo package into a list of the components names.
Definition: basicThermo.C:416
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
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
static const word null
An empty word.
Definition: word.H:77
const word dictName("particleTrackDict")
static const zero Zero
Definition: zero.H:97
iterator begin()
Iterator set to the beginning of the HashTable.
Definition: HashTableI.H:411
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
volScalarField & lookupOrConstruct(const fvMesh &mesh, const char *name) const
Definition: basicThermo.C:128
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
static Table::iterator lookupThermo(const dictionary &thermoTypeDict, Table *tablePtr, const int nCmpt, const char *cmptNames[], const word &thermoTypeName)
Generic lookup for thermodynamics package thermoTypeName.
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
List< word > wordList
A List of words.
Definition: fileName.H:54
label patchi
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
basicThermo(const basicThermo &)
Construct as copy (not implemented)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:167
virtual bool read()
Read thermophysical properties dictionary.
Definition: basicThermo.C:513
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
void validate(const string &app, const word &) const
Check that the thermodynamics package is consistent.
Definition: basicThermo.C:322
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
HashTable< const Type * > lookupClass(const bool strict=false) const
Lookup and return all objects of the given Type.
Namespace for OpenFOAM.
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:332