heThermo.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 "heThermo.H"
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
32 template<class BasicThermo, class MixtureType>
35 {
36  volScalarField::Boundary& hBf = h.boundaryFieldRef();
37 
38  forAll(hBf, patchi)
39  {
40  if (isA<gradientEnergyFvPatchScalarField>(hBf[patchi]))
41  {
42  refCast<gradientEnergyFvPatchScalarField>(hBf[patchi]).gradient()
43  = hBf[patchi].fvPatchField::snGrad();
44  }
45  else if (isA<mixedEnergyFvPatchScalarField>(hBf[patchi]))
46  {
47  refCast<mixedEnergyFvPatchScalarField>(hBf[patchi]).refGrad()
48  = hBf[patchi].fvPatchField::snGrad();
49  }
50  }
51 }
52 
53 
54 template<class BasicThermo, class MixtureType>
56 {
57  scalarField& heCells = he_.primitiveFieldRef();
58  const scalarField& pCells = this->p_;
59  const scalarField& TCells = this->T_;
60 
61  forAll(heCells, celli)
62  {
63  heCells[celli] =
64  this->cellMixture(celli).HE(pCells[celli], TCells[celli]);
65  }
66 
67  volScalarField::Boundary& heBf = he_.boundaryFieldRef();
68 
69  forAll(heBf, patchi)
70  {
71  heBf[patchi] == he
72  (
73  this->p_.boundaryField()[patchi],
74  this->T_.boundaryField()[patchi],
75  patchi
76  );
77  }
78 
79  this->heBoundaryCorrection(he_);
80 }
81 
82 
83 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
84 
85 
86 template<class BasicThermo, class MixtureType>
88 (
89  const fvMesh& mesh,
90  const word& phaseName
91 )
92 :
93  BasicThermo(mesh, phaseName),
94  MixtureType(*this, mesh, phaseName),
95 
96  he_
97  (
98  IOobject
99  (
100  BasicThermo::phasePropertyName
101  (
102  MixtureType::thermoType::heName()
103  ),
104  mesh.time().timeName(),
105  mesh,
106  IOobject::NO_READ,
107  IOobject::NO_WRITE
108  ),
109  mesh,
111  this->heBoundaryTypes(),
112  this->heBoundaryBaseTypes()
113  )
114 {
115  init();
116 }
117 
118 
119 template<class BasicThermo, class MixtureType>
121 (
122  const fvMesh& mesh,
123  const dictionary& dict,
124  const word& phaseName
125 )
126 :
127  BasicThermo(mesh, dict, phaseName),
128  MixtureType(*this, mesh, phaseName),
129 
130  he_
131  (
132  IOobject
133  (
134  BasicThermo::phasePropertyName
135  (
136  MixtureType::thermoType::heName()
137  ),
138  mesh.time().timeName(),
139  mesh,
140  IOobject::NO_READ,
141  IOobject::NO_WRITE
142  ),
143  mesh,
145  this->heBoundaryTypes(),
146  this->heBoundaryBaseTypes()
147  )
148 {
149  init();
150 }
151 
152 
153 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
154 
155 template<class BasicThermo, class MixtureType>
157 {}
158 
159 
160 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
161 
162 template<class BasicThermo, class MixtureType>
164 (
165  const volScalarField& p,
166  const volScalarField& T
167 ) const
168 {
169  const fvMesh& mesh = this->T_.mesh();
170 
172  (
174  (
175  "he",
176  mesh,
177  he_.dimensions()
178  )
179  );
180 
181  volScalarField& he = the.ref();
182  scalarField& heCells = he.primitiveFieldRef();
183  const scalarField& pCells = p;
184  const scalarField& TCells = T;
185 
186  forAll(heCells, celli)
187  {
188  heCells[celli] =
189  this->cellMixture(celli).HE(pCells[celli], TCells[celli]);
190  }
191 
192  volScalarField::Boundary& heBf = he.boundaryFieldRef();
193 
194  forAll(heBf, patchi)
195  {
196  scalarField& hep = heBf[patchi];
197  const scalarField& pp = p.boundaryField()[patchi];
198  const scalarField& Tp = T.boundaryField()[patchi];
199 
200  forAll(hep, facei)
201  {
202  hep[facei] =
203  this->patchFaceMixture(patchi, facei).HE(pp[facei], Tp[facei]);
204  }
205  }
206 
207  return the;
208 }
209 
210 
211 template<class BasicThermo, class MixtureType>
213 (
214  const scalarField& p,
215  const scalarField& T,
216  const labelList& cells
217 ) const
218 {
219  tmp<scalarField> the(new scalarField(T.size()));
220  scalarField& he = the.ref();
221 
222  forAll(T, celli)
223  {
224  he[celli] = this->cellMixture(cells[celli]).HE(p[celli], T[celli]);
225  }
226 
227  return the;
228 }
229 
230 
231 template<class BasicThermo, class MixtureType>
233 (
234  const scalarField& p,
235  const scalarField& T,
236  const label patchi
237 ) const
238 {
239  tmp<scalarField> the(new scalarField(T.size()));
240  scalarField& he = the.ref();
241 
242  forAll(T, facei)
243  {
244  he[facei] =
245  this->patchFaceMixture(patchi, facei).HE(p[facei], T[facei]);
246  }
247 
248  return the;
249 }
250 
251 
252 template<class BasicThermo, class MixtureType>
255 {
256  const fvMesh& mesh = this->T_.mesh();
257 
259  (
261  (
262  "hc",
263  mesh,
264  he_.dimensions()
265  )
266  );
267 
268  volScalarField& hcf = thc.ref();
269  scalarField& hcCells = hcf.primitiveFieldRef();
270 
271  forAll(hcCells, celli)
272  {
273  hcCells[celli] = this->cellMixture(celli).Hc();
274  }
275 
276  volScalarField::Boundary& hcfBf = hcf.boundaryFieldRef();
277 
278  forAll(hcfBf, patchi)
279  {
280  scalarField& hcp = hcfBf[patchi];
281 
282  forAll(hcp, facei)
283  {
284  hcp[facei] = this->patchFaceMixture(patchi, facei).Hc();
285  }
286  }
287 
288  return thc;
289 }
290 
291 
292 template<class BasicThermo, class MixtureType>
294 (
295  const scalarField& p,
296  const scalarField& T,
297  const label patchi
298 ) const
299 {
300  tmp<scalarField> tCp(new scalarField(T.size()));
301  scalarField& cp = tCp.ref();
302 
303  forAll(T, facei)
304  {
305  cp[facei] =
306  this->patchFaceMixture(patchi, facei).Cp(p[facei], T[facei]);
307  }
308 
309  return tCp;
310 }
311 
312 
313 template<class BasicThermo, class MixtureType>
316 {
317  const fvMesh& mesh = this->T_.mesh();
318 
320  (
322  (
323  "Cp",
324  mesh,
326  )
327  );
328 
329  volScalarField& cp = tCp.ref();
330 
331  forAll(this->T_, celli)
332  {
333  cp[celli] =
334  this->cellMixture(celli).Cp(this->p_[celli], this->T_[celli]);
335  }
336 
337  volScalarField::Boundary& cpBf = cp.boundaryFieldRef();
338 
339  forAll(cpBf, patchi)
340  {
341  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
342  const fvPatchScalarField& pT = this->T_.boundaryField()[patchi];
343  fvPatchScalarField& pCp = cpBf[patchi];
344 
345  forAll(pT, facei)
346  {
347  pCp[facei] =
348  this->patchFaceMixture(patchi, facei).Cp(pp[facei], pT[facei]);
349  }
350  }
351 
352  return tCp;
353 }
354 
355 
356 template<class BasicThermo, class MixtureType>
359 (
360  const scalarField& p,
361  const scalarField& T,
362  const label patchi
363 ) const
364 {
365  tmp<scalarField> tCv(new scalarField(T.size()));
366  scalarField& cv = tCv.ref();
367 
368  forAll(T, facei)
369  {
370  cv[facei] =
371  this->patchFaceMixture(patchi, facei).Cv(p[facei], T[facei]);
372  }
373 
374  return tCv;
375 }
376 
377 
378 template<class BasicThermo, class MixtureType>
381 {
382  const fvMesh& mesh = this->T_.mesh();
383 
385  (
387  (
388  "Cv",
389  mesh,
391  )
392  );
393 
394  volScalarField& cv = tCv.ref();
395 
396  forAll(this->T_, celli)
397  {
398  cv[celli] =
399  this->cellMixture(celli).Cv(this->p_[celli], this->T_[celli]);
400  }
401 
402  volScalarField::Boundary& cvBf = cv.boundaryFieldRef();
403 
404  forAll(cvBf, patchi)
405  {
406  cvBf[patchi] = Cv
407  (
408  this->p_.boundaryField()[patchi],
409  this->T_.boundaryField()[patchi],
410  patchi
411  );
412  }
413 
414  return tCv;
415 }
416 
417 
418 template<class BasicThermo, class MixtureType>
420 (
421  const scalarField& p,
422  const scalarField& T,
423  const label patchi
424 ) const
425 {
426  tmp<scalarField> tgamma(new scalarField(T.size()));
427  scalarField& gamma = tgamma.ref();
428 
429  forAll(T, facei)
430  {
431  gamma[facei] =
432  this->patchFaceMixture(patchi, facei).gamma(p[facei], T[facei]);
433  }
434 
435  return tgamma;
436 }
437 
438 
439 template<class BasicThermo, class MixtureType>
442 {
443  const fvMesh& mesh = this->T_.mesh();
444 
445  tmp<volScalarField> tgamma
446  (
448  (
449  "gamma",
450  mesh,
451  dimless
452  )
453  );
454 
455  volScalarField& gamma = tgamma.ref();
456 
457  forAll(this->T_, celli)
458  {
459  gamma[celli] =
460  this->cellMixture(celli).gamma(this->p_[celli], this->T_[celli]);
461  }
462 
463  volScalarField::Boundary& gammaBf = gamma.boundaryFieldRef();
464 
465  forAll(gammaBf, patchi)
466  {
467  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
468  const fvPatchScalarField& pT = this->T_.boundaryField()[patchi];
469  fvPatchScalarField& pgamma = gammaBf[patchi];
470 
471  forAll(pT, facei)
472  {
473  pgamma[facei] = this->patchFaceMixture(patchi, facei).gamma
474  (
475  pp[facei],
476  pT[facei]
477  );
478  }
479  }
480 
481  return tgamma;
482 }
483 
484 
485 template<class BasicThermo, class MixtureType>
487 (
488  const scalarField& p,
489  const scalarField& T,
490  const label patchi
491 ) const
492 {
493  tmp<scalarField> tCpv(new scalarField(T.size()));
494  scalarField& Cpv = tCpv.ref();
495 
496  forAll(T, facei)
497  {
498  Cpv[facei] =
499  this->patchFaceMixture(patchi, facei).Cpv(p[facei], T[facei]);
500  }
501 
502  return tCpv;
503 }
504 
505 
506 template<class BasicThermo, class MixtureType>
509 {
510  const fvMesh& mesh = this->T_.mesh();
511 
513  (
515  (
516  "Cpv",
517  mesh,
519  )
520  );
521 
522  volScalarField& Cpv = tCpv.ref();
523 
524  forAll(this->T_, celli)
525  {
526  Cpv[celli] =
527  this->cellMixture(celli).Cpv(this->p_[celli], this->T_[celli]);
528  }
529 
530  volScalarField::Boundary& CpvBf = Cpv.boundaryFieldRef();
531 
532  forAll(CpvBf, patchi)
533  {
534  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
535  const fvPatchScalarField& pT = this->T_.boundaryField()[patchi];
536  fvPatchScalarField& pCpv = CpvBf[patchi];
537 
538  forAll(pT, facei)
539  {
540  pCpv[facei] =
541  this->patchFaceMixture(patchi, facei).Cpv(pp[facei], pT[facei]);
542  }
543  }
544 
545  return tCpv;
546 }
547 
548 
549 template<class BasicThermo, class MixtureType>
551 (
552  const scalarField& p,
553  const scalarField& T,
554  const label patchi
555 ) const
556 {
557  tmp<scalarField> tCpByCpv(new scalarField(T.size()));
558  scalarField& CpByCpv = tCpByCpv.ref();
559 
560  forAll(T, facei)
561  {
562  CpByCpv[facei] =
563  this->patchFaceMixture(patchi, facei).CpByCpv(p[facei], T[facei]);
564  }
565 
566  return tCpByCpv;
567 }
568 
569 
570 template<class BasicThermo, class MixtureType>
573 {
574  const fvMesh& mesh = this->T_.mesh();
575 
576  tmp<volScalarField> tCpByCpv
577  (
579  (
580  "CpByCpv",
581  mesh,
582  dimless
583  )
584  );
585 
586  volScalarField& CpByCpv = tCpByCpv.ref();
587 
588  forAll(this->T_, celli)
589  {
590  CpByCpv[celli] = this->cellMixture(celli).CpByCpv
591  (
592  this->p_[celli],
593  this->T_[celli]
594  );
595  }
596 
597  volScalarField::Boundary& CpByCpvBf =
598  CpByCpv.boundaryFieldRef();
599 
600  forAll(CpByCpvBf, patchi)
601  {
602  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
603  const fvPatchScalarField& pT = this->T_.boundaryField()[patchi];
604  fvPatchScalarField& pCpByCpv = CpByCpvBf[patchi];
605 
606  forAll(pT, facei)
607  {
608  pCpByCpv[facei] = this->patchFaceMixture(patchi, facei).CpByCpv
609  (
610  pp[facei],
611  pT[facei]
612  );
613  }
614  }
615 
616  return tCpByCpv;
617 }
618 
619 
620 template<class BasicThermo, class MixtureType>
622 (
623  const scalarField& h,
624  const scalarField& p,
625  const scalarField& T0,
626  const labelList& cells
627 ) const
628 {
629  tmp<scalarField> tT(new scalarField(h.size()));
630  scalarField& T = tT.ref();
631 
632  forAll(h, celli)
633  {
634  T[celli] =
635  this->cellMixture(cells[celli]).THE(h[celli], p[celli], T0[celli]);
636  }
637 
638  return tT;
639 }
640 
641 
642 template<class BasicThermo, class MixtureType>
644 (
645  const scalarField& h,
646  const scalarField& p,
647  const scalarField& T0,
648  const label patchi
649 ) const
650 {
651 
652  tmp<scalarField> tT(new scalarField(h.size()));
653  scalarField& T = tT.ref();
654  forAll(h, facei)
655  {
656  T[facei] = this->patchFaceMixture
657  (
658  patchi,
659  facei
660  ).THE(h[facei], p[facei], T0[facei]);
661  }
662 
663  return tT;
664 }
665 
666 
667 template<class BasicThermo, class MixtureType>
669 (
670 ) const
671 {
672  const fvMesh& mesh = this->T_.mesh();
673 
675  (
677  (
678  "W",
679  mesh,
681  )
682  );
683 
684  volScalarField& W = tW.ref();
685  scalarField& WCells = W.primitiveFieldRef();
686 
687  forAll(WCells, celli)
688  {
689  WCells[celli] = this->cellMixture(celli).W();
690  }
691 
692  volScalarField::Boundary& WBf = W.boundaryFieldRef();
693 
694  forAll(WBf, patchi)
695  {
696  scalarField& Wp = WBf[patchi];
697  forAll(Wp, facei)
698  {
699  Wp[facei] = this->patchFaceMixture(patchi, facei).W();
700  }
701  }
702 
703  return tW;
704 }
705 
706 
707 template<class BasicThermo, class MixtureType>
709 (
710  const label patchi
711 ) const
712 {
713  const fvMesh& mesh = this->T_.mesh();
714 
716  scalarField& W = tW.ref();
717  forAll(W, facei)
718  {
719  W[facei] = this->patchFaceMixture(patchi, facei).W();
720  }
721 
722  return tW;
723 }
724 
725 
726 template<class BasicThermo, class MixtureType>
729 {
730  tmp<Foam::volScalarField> kappa(Cp()*this->alpha_);
731  kappa.ref().rename("kappa");
732  return kappa;
733 }
734 
735 
736 template<class BasicThermo, class MixtureType>
738 (
739  const label patchi
740 ) const
741 {
742  return
743  Cp
744  (
745  this->p_.boundaryField()[patchi],
746  this->T_.boundaryField()[patchi],
747  patchi
748  )*this->alpha_.boundaryField()[patchi];
749 }
750 
751 
752 template<class BasicThermo, class MixtureType>
755 {
756  tmp<Foam::volScalarField> alphaEff(this->CpByCpv()*this->alpha_);
757  alphaEff.ref().rename("alphahe");
758  return alphaEff;
759 }
760 
761 
762 template<class BasicThermo, class MixtureType>
765 {
766  return
767  this->CpByCpv
768  (
769  this->p_.boundaryField()[patchi],
770  this->T_.boundaryField()[patchi],
771  patchi
772  )
773  *this->alpha_.boundaryField()[patchi];
774 }
775 
776 
777 template<class BasicThermo, class MixtureType>
780 (
781  const volScalarField& alphat
782 ) const
783 {
784  tmp<Foam::volScalarField> kappaEff(Cp()*(this->alpha_ + alphat));
785  kappaEff.ref().rename("kappaEff");
786  return kappaEff;
787 }
788 
789 
790 template<class BasicThermo, class MixtureType>
793 (
794  const scalarField& alphat,
795  const label patchi
796 ) const
797 {
798  return
799  Cp
800  (
801  this->p_.boundaryField()[patchi],
802  this->T_.boundaryField()[patchi],
803  patchi
804  )
805  *(
806  this->alpha_.boundaryField()[patchi]
807  + alphat
808  );
809 }
810 
811 
812 template<class BasicThermo, class MixtureType>
815 (
816  const volScalarField& alphat
817 ) const
818 {
819  tmp<Foam::volScalarField> alphaEff(this->CpByCpv()*(this->alpha_ + alphat));
820  alphaEff.ref().rename("alphaEff");
821  return alphaEff;
822 }
823 
824 
825 template<class BasicThermo, class MixtureType>
828 (
829  const scalarField& alphat,
830  const label patchi
831 ) const
832 {
833  return
834  this->CpByCpv
835  (
836  this->p_.boundaryField()[patchi],
837  this->T_.boundaryField()[patchi],
838  patchi
839  )
840  *(
841  this->alpha_.boundaryField()[patchi]
842  + alphat
843  );
844 }
845 
846 
847 template<class BasicThermo, class MixtureType>
849 {
850  if (BasicThermo::read())
851  {
852  MixtureType::read(*this);
853  return true;
854  }
855  else
856  {
857  return false;
858  }
859 }
860 
861 
862 // ************************************************************************* //
virtual tmp< volScalarField > kappa() const
Thermal diffusivity for temperature of mixture [W/m/K].
Definition: heThermo.C:728
volScalarField & he
Definition: YEEqn.H:50
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
virtual tmp< scalarField > THE(const scalarField &he, const scalarField &p, const scalarField &T0, const labelList &cells) const
Temperature from enthalpy/internal energy for cell-set.
Definition: heThermo.C:622
scalar Cv(const scalar p, const scalar T) const
Definition: HtoEthermo.H:2
#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
virtual tmp< volScalarField > gamma() const
Gamma = Cp/Cv [].
Definition: heThermo.C:441
virtual tmp< volScalarField > CpByCpv() const
Heat capacity ratio [].
Definition: heThermo.C:572
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
virtual ~heThermo()
Destructor.
Definition: heThermo.C:156
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
virtual tmp< volScalarField > Cpv() const
Heat capacity at constant pressure/volume [J/kg/K].
Definition: heThermo.C:508
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
Definition: heThermo.H:147
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
virtual tmp< volScalarField > hc() const
Chemical enthalpy [J/kg].
Definition: heThermo.C:254
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
virtual tmp< volScalarField > Cp() const
Heat capacity at constant pressure [J/kg/K].
Definition: heThermo.C:315
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
virtual tmp< volScalarField > alphahe() const
Thermal diffusivity for energy of mixture [kg/m/s].
Definition: heThermo.C:754
volScalarField scalarField(fieldObject, mesh)
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Effective thermal turbulent diffusivity of mixture [kg/m/s].
Definition: heThermo.C:815
const volScalarField & cp
const dimensionSet dimMoles(0, 0, 0, 0, 1, 0, 0)
Definition: dimensionSets.H:53
virtual void rename(const word &newName)
Rename.
Definition: regIOobject.C:410
const volScalarField & T
const dimensionSet dimEnergy
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
virtual bool read()
Read thermophysical properties dictionary.
Definition: heThermo.C:848
label patchi
virtual tmp< volScalarField > kappaEff(const volScalarField &) const
Effective thermal turbulent diffusivity for temperature.
Definition: heThermo.C:780
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
virtual tmp< volScalarField > Cv() const
Heat capacity at constant volume [J/kg/K].
Definition: heThermo.C:380
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
Enthalpy/Internal energy for a mixture.
Definition: heThermo.H:49
const scalarList W(::W(thermo))
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
scalar Cp(const scalar p, const scalar T) const
Definition: EtoHthermo.H:2
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void heBoundaryCorrection(volScalarField &he)
Correct the enthalpy/internal energy field boundaries.
Definition: heThermo.C:34
virtual tmp< volScalarField > W() const
Molecular weight [kg/kmol].
Definition: heThermo.C:669