heheuPsiThermo.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 "heheuPsiThermo.H"
27 #include "fvMesh.H"
29 
30 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
31 
32 template<class BasicPsiThermo, class MixtureType>
34 {
35  const scalarField& hCells = this->he_;
36  const scalarField& heuCells = this->heu_;
37  const scalarField& pCells = this->p_;
38 
39  scalarField& TCells = this->T_.primitiveFieldRef();
40  scalarField& TuCells = this->Tu_.primitiveFieldRef();
41  scalarField& psiCells = this->psi_.primitiveFieldRef();
42  scalarField& muCells = this->mu_.primitiveFieldRef();
43  scalarField& alphaCells = this->alpha_.primitiveFieldRef();
44 
45  forAll(TCells, celli)
46  {
47  const typename MixtureType::thermoType& mixture_ =
48  this->cellMixture(celli);
49 
50  TCells[celli] = mixture_.THE
51  (
52  hCells[celli],
53  pCells[celli],
54  TCells[celli]
55  );
56 
57  psiCells[celli] = mixture_.psi(pCells[celli], TCells[celli]);
58 
59  muCells[celli] = mixture_.mu(pCells[celli], TCells[celli]);
60  alphaCells[celli] = mixture_.alphah(pCells[celli], TCells[celli]);
61 
62  TuCells[celli] = this->cellReactants(celli).THE
63  (
64  heuCells[celli],
65  pCells[celli],
66  TuCells[celli]
67  );
68  }
69 
70  volScalarField::Boundary& pBf =
71  this->p_.boundaryFieldRef();
72 
73  volScalarField::Boundary& TBf =
74  this->T_.boundaryFieldRef();
75 
76  volScalarField::Boundary& TuBf =
77  this->Tu_.boundaryFieldRef();
78 
79  volScalarField::Boundary& psiBf =
80  this->psi_.boundaryFieldRef();
81 
82  volScalarField::Boundary& heBf =
83  this->he().boundaryFieldRef();
84 
85  volScalarField::Boundary& heuBf =
86  this->heu().boundaryFieldRef();
87 
88  volScalarField::Boundary& muBf =
89  this->mu_.boundaryFieldRef();
90 
91  volScalarField::Boundary& alphaBf =
92  this->alpha_.boundaryFieldRef();
93 
94  forAll(this->T_.boundaryField(), patchi)
95  {
96  fvPatchScalarField& pp = pBf[patchi];
97  fvPatchScalarField& pT = TBf[patchi];
98  fvPatchScalarField& pTu = TuBf[patchi];
99  fvPatchScalarField& ppsi = psiBf[patchi];
100  fvPatchScalarField& phe = heBf[patchi];
101  fvPatchScalarField& pheu = heuBf[patchi];
102  fvPatchScalarField& pmu = muBf[patchi];
103  fvPatchScalarField& palpha = alphaBf[patchi];
104 
105  if (pT.fixesValue())
106  {
107  forAll(pT, facei)
108  {
109  const typename MixtureType::thermoType& mixture_ =
110  this->patchFaceMixture(patchi, facei);
111 
112  phe[facei] = mixture_.HE(pp[facei], pT[facei]);
113 
114  ppsi[facei] = mixture_.psi(pp[facei], pT[facei]);
115  pmu[facei] = mixture_.mu(pp[facei], pT[facei]);
116  palpha[facei] = mixture_.alphah(pp[facei], pT[facei]);
117  }
118  }
119  else
120  {
121  forAll(pT, facei)
122  {
123  const typename MixtureType::thermoType& mixture_ =
124  this->patchFaceMixture(patchi, facei);
125 
126  pT[facei] = mixture_.THE(phe[facei], pp[facei], pT[facei]);
127 
128  ppsi[facei] = mixture_.psi(pp[facei], pT[facei]);
129  pmu[facei] = mixture_.mu(pp[facei], pT[facei]);
130  palpha[facei] = mixture_.alphah(pp[facei], pT[facei]);
131 
132  pTu[facei] =
133  this->patchFaceReactants(patchi, facei)
134  .THE(pheu[facei], pp[facei], pTu[facei]);
135  }
136  }
137  }
138 }
139 
140 
141 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
142 
143 template<class BasicPsiThermo, class MixtureType>
145 (
146  const fvMesh& mesh,
147  const word& phaseName
148 )
149 :
151  Tu_
152  (
153  IOobject
154  (
155  "Tu",
156  mesh.time().timeName(),
157  mesh,
158  IOobject::MUST_READ,
159  IOobject::AUTO_WRITE
160  ),
161  mesh
162  ),
163 
164  heu_
165  (
166  IOobject
167  (
168  MixtureType::thermoType::heName() + 'u',
169  mesh.time().timeName(),
170  mesh,
171  IOobject::NO_READ,
172  IOobject::NO_WRITE
173  ),
174  mesh,
175  dimensionSet(0, 2, -2, 0, 0),
176  this->heuBoundaryTypes()
177  )
178 {
179  scalarField& heuCells = this->heu_.primitiveFieldRef();
180  const scalarField& pCells = this->p_;
181  const scalarField& TuCells = this->Tu_;
182 
183  forAll(heuCells, celli)
184  {
185  heuCells[celli] = this->cellReactants(celli).HE
186  (
187  pCells[celli],
188  TuCells[celli]
189  );
190  }
191 
192  volScalarField::Boundary& heuBf = heu_.boundaryFieldRef();
193 
194  forAll(heuBf, patchi)
195  {
196  fvPatchScalarField& pheu = heuBf[patchi];
197  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
198  const fvPatchScalarField& pTu = this->Tu_.boundaryField()[patchi];
199 
200  forAll(pheu, facei)
201  {
202  pheu[facei] = this->patchFaceReactants(patchi, facei).HE
203  (
204  pp[facei],
205  pTu[facei]
206  );
207  }
208  }
209 
210  this->heuBoundaryCorrection(this->heu_);
211 
212  calculate();
213  this->psi_.oldTime(); // Switch on saving old time
214 }
215 
216 
217 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
218 
219 template<class BasicPsiThermo, class MixtureType>
221 {}
222 
223 
224 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
225 
226 template<class BasicPsiThermo, class MixtureType>
228 {
229  if (debug)
230  {
231  InfoInFunction << endl;
232  }
233 
234  // force the saving of the old-time values
235  this->psi_.oldTime();
236 
237  calculate();
238 
239  if (debug)
240  {
241  Info<< " Finished" << endl;
242  }
243 }
244 
245 
246 template<class BasicPsiThermo, class MixtureType>
249 (
250  const scalarField& p,
251  const scalarField& Tu,
252  const labelList& cells
253 ) const
254 {
255  tmp<scalarField> theu(new scalarField(Tu.size()));
256  scalarField& heu = theu.ref();
257 
258  forAll(Tu, celli)
259  {
260  heu[celli] = this->cellReactants(cells[celli]).HE(p[celli], Tu[celli]);
261  }
262 
263  return theu;
264 }
265 
266 
267 template<class BasicPsiThermo, class MixtureType>
270 (
271  const scalarField& p,
272  const scalarField& Tu,
273  const label patchi
274 ) const
275 {
276  tmp<scalarField> theu(new scalarField(Tu.size()));
277  scalarField& heu = theu.ref();
278 
279  forAll(Tu, facei)
280  {
281  heu[facei] =
282  this->patchFaceReactants(patchi, facei).HE(p[facei], Tu[facei]);
283  }
284 
285  return theu;
286 }
287 
288 
289 template<class BasicPsiThermo, class MixtureType>
292 {
294  (
295  new volScalarField
296  (
297  IOobject
298  (
299  "Tb",
300  this->T_.time().timeName(),
301  this->T_.db(),
302  IOobject::NO_READ,
303  IOobject::NO_WRITE,
304  false
305  ),
306  this->T_
307  )
308  );
309 
310  volScalarField& Tb_ = tTb.ref();
311  scalarField& TbCells = Tb_.primitiveFieldRef();
312  const scalarField& pCells = this->p_;
313  const scalarField& TCells = this->T_;
314  const scalarField& hCells = this->he_;
315 
316  forAll(TbCells, celli)
317  {
318  TbCells[celli] = this->cellProducts(celli).THE
319  (
320  hCells[celli],
321  pCells[celli],
322  TCells[celli]
323  );
324  }
325 
326  volScalarField::Boundary& TbBf = Tb_.boundaryFieldRef();
327 
328  forAll(TbBf, patchi)
329  {
330  fvPatchScalarField& pTb = TbBf[patchi];
331 
332  const fvPatchScalarField& ph = this->he_.boundaryField()[patchi];
333  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
334  const fvPatchScalarField& pT = this->T_.boundaryField()[patchi];
335 
336  forAll(pTb, facei)
337  {
338  pTb[facei] =
339  this->patchFaceProducts(patchi, facei)
340  .THE(ph[facei], pp[facei], pT[facei]);
341  }
342  }
343 
344  return tTb;
345 }
346 
347 
348 template<class BasicPsiThermo, class MixtureType>
351 {
352  tmp<volScalarField> tpsiu
353  (
354  new volScalarField
355  (
356  IOobject
357  (
358  "psiu",
359  this->psi_.time().timeName(),
360  this->psi_.db(),
361  IOobject::NO_READ,
362  IOobject::NO_WRITE,
363  false
364  ),
365  this->psi_.mesh(),
366  this->psi_.dimensions()
367  )
368  );
369 
370  volScalarField& psiu = tpsiu.ref();
371  scalarField& psiuCells = psiu.primitiveFieldRef();
372  const scalarField& TuCells = this->Tu_;
373  const scalarField& pCells = this->p_;
374 
375  forAll(psiuCells, celli)
376  {
377  psiuCells[celli] =
378  this->cellReactants(celli).psi(pCells[celli], TuCells[celli]);
379  }
380 
381  volScalarField::Boundary& psiuBf = psiu.boundaryFieldRef();
382 
383  forAll(psiuBf, patchi)
384  {
385  fvPatchScalarField& ppsiu = psiuBf[patchi];
386 
387  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
388  const fvPatchScalarField& pTu = this->Tu_.boundaryField()[patchi];
389 
390  forAll(ppsiu, facei)
391  {
392  ppsiu[facei] =
393  this->
394  patchFaceReactants(patchi, facei).psi(pp[facei], pTu[facei]);
395  }
396  }
397 
398  return tpsiu;
399 }
400 
401 
402 template<class BasicPsiThermo, class MixtureType>
405 {
406  tmp<volScalarField> tpsib
407  (
408  new volScalarField
409  (
410  IOobject
411  (
412  "psib",
413  this->psi_.time().timeName(),
414  this->psi_.db(),
415  IOobject::NO_READ,
416  IOobject::NO_WRITE,
417  false
418  ),
419  this->psi_.mesh(),
420  this->psi_.dimensions()
421  )
422  );
423 
424  volScalarField& psib = tpsib.ref();
425  scalarField& psibCells = psib.primitiveFieldRef();
426  const volScalarField Tb_(Tb());
427  const scalarField& TbCells = Tb_;
428  const scalarField& pCells = this->p_;
429 
430  forAll(psibCells, celli)
431  {
432  psibCells[celli] =
433  this->cellReactants(celli).psi(pCells[celli], TbCells[celli]);
434  }
435 
436  volScalarField::Boundary& psibBf = psib.boundaryFieldRef();
437 
438  forAll(psibBf, patchi)
439  {
440  fvPatchScalarField& ppsib = psibBf[patchi];
441 
442  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
443  const fvPatchScalarField& pTb = Tb_.boundaryField()[patchi];
444 
445  forAll(ppsib, facei)
446  {
447  ppsib[facei] =
448  this->patchFaceReactants
449  (patchi, facei).psi(pp[facei], pTb[facei]);
450  }
451  }
452 
453  return tpsib;
454 }
455 
456 
457 template<class BasicPsiThermo, class MixtureType>
460 {
462  (
463  new volScalarField
464  (
465  IOobject
466  (
467  "muu",
468  this->T_.time().timeName(),
469  this->T_.db(),
470  IOobject::NO_READ,
471  IOobject::NO_WRITE,
472  false
473  ),
474  this->T_.mesh(),
475  dimensionSet(1, -1, -1, 0, 0)
476  )
477  );
478 
479  volScalarField& muu_ = tmuu.ref();
480  scalarField& muuCells = muu_.primitiveFieldRef();
481  const scalarField& pCells = this->p_;
482  const scalarField& TuCells = this->Tu_;
483 
484  forAll(muuCells, celli)
485  {
486  muuCells[celli] = this->cellReactants(celli).mu
487  (
488  pCells[celli],
489  TuCells[celli]
490  );
491  }
492 
493  volScalarField::Boundary& muuBf = muu_.boundaryFieldRef();
494 
495  forAll(muuBf, patchi)
496  {
497  fvPatchScalarField& pMuu = muuBf[patchi];
498  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
499  const fvPatchScalarField& pTu = this->Tu_.boundaryField()[patchi];
500 
501  forAll(pMuu, facei)
502  {
503  pMuu[facei] = this->patchFaceReactants(patchi, facei).mu
504  (
505  pp[facei],
506  pTu[facei]
507  );
508  }
509  }
510 
511  return tmuu;
512 }
513 
514 
515 template<class BasicPsiThermo, class MixtureType>
518 {
520  (
521  new volScalarField
522  (
523  IOobject
524  (
525  "mub",
526  this->T_.time().timeName(),
527  this->T_.db(),
528  IOobject::NO_READ,
529  IOobject::NO_WRITE,
530  false
531  ),
532  this->T_.mesh(),
533  dimensionSet(1, -1, -1, 0, 0)
534  )
535  );
536 
537  volScalarField& mub_ = tmub.ref();
538  scalarField& mubCells = mub_.primitiveFieldRef();
539  const volScalarField Tb_(Tb());
540  const scalarField& pCells = this->p_;
541  const scalarField& TbCells = Tb_;
542 
543  forAll(mubCells, celli)
544  {
545  mubCells[celli] = this->cellProducts(celli).mu
546  (
547  pCells[celli],
548  TbCells[celli]
549  );
550  }
551 
552  volScalarField::Boundary& mubBf = mub_.boundaryFieldRef();
553 
554  forAll(mubBf, patchi)
555  {
556  fvPatchScalarField& pMub = mubBf[patchi];
557  const fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
558  const fvPatchScalarField& pTb = Tb_.boundaryField()[patchi];
559 
560  forAll(pMub, facei)
561  {
562  pMub[facei] = this->patchFaceProducts(patchi, facei).mu
563  (
564  pp[facei],
565  pTb[facei]
566  );
567  }
568  }
569 
570  return tmub;
571 }
572 
573 
574 // ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:51
virtual tmp< volScalarField > mub() const
Dynamic viscosity of burnt gas [kg/ms].
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
virtual tmp< volScalarField > psiu() const
Unburnt gas compressibility [s^2/m^2].
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 > Tb() const
Burnt gas temperature [K].
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
virtual tmp< volScalarField > psib() const
Burnt gas compressibility [s^2/m^2].
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual tmp< volScalarField > muu() const
Dynamic viscosity of unburnt gas [kg/ms].
Dimension set for the base types.
Definition: dimensionSet.H:120
dynamicFvMesh & mesh
fvPatchField< scalar > fvPatchScalarField
A class for handling words, derived from string.
Definition: word.H:59
virtual void correct()
Update properties.
volScalarField scalarField(fieldObject, mesh)
virtual ~heheuPsiThermo()
Destructor.
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
Internal & ref()
Return a reference to the dimensioned internal field.
label patchi
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
Enthalpy/Internal energy for a mixture.
Definition: heThermo.H:49
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
messageStream Info
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
virtual volScalarField & heu()
Unburnt gas enthalpy [J/kg].
#define InfoInFunction
Report an information message using Foam::Info.