heheuPsiThermo.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-2022 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& CpCells = this->Cp_.primitiveFieldRef();
42  scalarField& CvCells = this->Cv_.primitiveFieldRef();
43  scalarField& psiCells = this->psi_.primitiveFieldRef();
44  scalarField& muCells = this->mu_.primitiveFieldRef();
45  scalarField& kappaCells = this->kappa_.primitiveFieldRef();
46 
47  forAll(TCells, celli)
48  {
49  const typename MixtureType::thermoMixtureType& thermoMixture =
50  this->cellThermoMixture(celli);
51 
52  const typename MixtureType::transportMixtureType& transportMixture =
53  this->cellTransportMixture(celli, thermoMixture);
54 
55  TCells[celli] = thermoMixture.THE
56  (
57  hCells[celli],
58  pCells[celli],
59  TCells[celli]
60  );
61 
62  CpCells[celli] = thermoMixture.Cp(pCells[celli], TCells[celli]);
63  CvCells[celli] = thermoMixture.Cv(pCells[celli], TCells[celli]);
64  psiCells[celli] = thermoMixture.psi(pCells[celli], TCells[celli]);
65 
66  muCells[celli] = transportMixture.mu(pCells[celli], TCells[celli]);
67  kappaCells[celli] =
68  transportMixture.kappa(pCells[celli], TCells[celli]);
69 
70  TuCells[celli] = this->cellReactants(celli).THE
71  (
72  heuCells[celli],
73  pCells[celli],
74  TuCells[celli]
75  );
76  }
77 
78  volScalarField::Boundary& pBf =
79  this->p_.boundaryFieldRef();
80 
81  volScalarField::Boundary& TBf =
82  this->T_.boundaryFieldRef();
83 
84  volScalarField::Boundary& TuBf =
85  this->Tu_.boundaryFieldRef();
86 
87  volScalarField::Boundary& CpBf =
88  this->Cp_.boundaryFieldRef();
89 
90  volScalarField::Boundary& CvBf =
91  this->Cv_.boundaryFieldRef();
92 
93  volScalarField::Boundary& psiBf =
94  this->psi_.boundaryFieldRef();
95 
96  volScalarField::Boundary& heBf =
97  this->he().boundaryFieldRef();
98 
99  volScalarField::Boundary& heuBf =
100  this->heu().boundaryFieldRef();
101 
102  volScalarField::Boundary& muBf =
103  this->mu_.boundaryFieldRef();
104 
105  volScalarField::Boundary& kappaBf =
106  this->kappa_.boundaryFieldRef();
107 
108  forAll(this->T_.boundaryField(), patchi)
109  {
110  fvPatchScalarField& pp = pBf[patchi];
111  fvPatchScalarField& pT = TBf[patchi];
112  fvPatchScalarField& pTu = TuBf[patchi];
113  fvPatchScalarField& pCp = CpBf[patchi];
114  fvPatchScalarField& pCv = CvBf[patchi];
115  fvPatchScalarField& ppsi = psiBf[patchi];
116  fvPatchScalarField& phe = heBf[patchi];
117  fvPatchScalarField& pheu = heuBf[patchi];
118  fvPatchScalarField& pmu = muBf[patchi];
119  fvPatchScalarField& pkappa = kappaBf[patchi];
120 
121  if (pT.fixesValue())
122  {
123  forAll(pT, facei)
124  {
125  const typename MixtureType::thermoMixtureType&
126  thermoMixture = this->patchFaceThermoMixture(patchi, facei);
127 
128  const typename MixtureType::transportMixtureType&
129  transportMixture =
130  this->patchFaceTransportMixture
131  (patchi, facei, thermoMixture);
132 
133  phe[facei] = thermoMixture.HE(pp[facei], pT[facei]);
134 
135  pCp[facei] = thermoMixture.Cp(pp[facei], pT[facei]);
136  pCv[facei] = thermoMixture.Cv(pp[facei], pT[facei]);
137  ppsi[facei] = thermoMixture.psi(pp[facei], pT[facei]);
138  pmu[facei] = transportMixture.mu(pp[facei], pT[facei]);
139  pkappa[facei] = transportMixture.kappa(pp[facei], pT[facei]);
140  }
141  }
142  else
143  {
144  forAll(pT, facei)
145  {
146  const typename MixtureType::thermoMixtureType&
147  thermoMixture = this->patchFaceThermoMixture(patchi, facei);
148 
149  const typename MixtureType::transportMixtureType&
150  transportMixture =
151  this->patchFaceTransportMixture
152  (patchi, facei, thermoMixture);
153 
154  pT[facei] = thermoMixture.THE(phe[facei], pp[facei], pT[facei]);
155 
156  pCp[facei] = thermoMixture.Cp(pp[facei], pT[facei]);
157  pCv[facei] = thermoMixture.Cv(pp[facei], pT[facei]);
158  ppsi[facei] = thermoMixture.psi(pp[facei], pT[facei]);
159  pmu[facei] = transportMixture.mu(pp[facei], pT[facei]);
160  pkappa[facei] = transportMixture.kappa(pp[facei], pT[facei]);
161 
162  pTu[facei] =
163  this->patchFaceReactants(patchi, facei)
164  .THE(pheu[facei], pp[facei], pTu[facei]);
165  }
166  }
167  }
168 }
169 
170 
171 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
172 
173 template<class BasicPsiThermo, class MixtureType>
175 (
176  const fvMesh& mesh,
177  const word& phaseName
178 )
179 :
181  Tu_
182  (
183  IOobject
184  (
185  "Tu",
186  mesh.time().timeName(),
187  mesh,
188  IOobject::MUST_READ,
189  IOobject::AUTO_WRITE
190  ),
191  mesh
192  ),
193  heu_
194  (
195  IOobject
196  (
197  MixtureType::thermoType::heName() + 'u',
198  mesh.time().timeName(),
199  mesh,
200  IOobject::NO_READ,
201  IOobject::NO_WRITE
202  ),
203  this->volScalarFieldProperty
204  (
205  MixtureType::thermoType::heName() + 'u',
207  &MixtureType::cellReactants,
208  &MixtureType::patchFaceReactants,
209  &MixtureType::thermoMixtureType::HE,
210  this->p_,
211  this->Tu_
212  ),
213  this->heuBoundaryTypes()
214  )
215 {
216  this->heuBoundaryCorrection(this->heu_);
217 
218  calculate();
219 
220  this->psi_.oldTime(); // Switch on saving old time
221 }
222 
223 
224 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
225 
226 template<class BasicPsiThermo, class MixtureType>
228 {}
229 
230 
231 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
232 
233 template<class BasicPsiThermo, class MixtureType>
235 {
236  if (debug)
237  {
238  InfoInFunction << endl;
239  }
240 
241  // force the saving of the old-time values
242  this->psi_.oldTime();
243 
244  calculate();
245 
246  if (debug)
247  {
248  Info<< " Finished" << endl;
249  }
250 }
251 
252 
253 template<class BasicPsiThermo, class MixtureType>
256 (
257  const scalarField& Tu,
258  const labelList& cells
259 ) const
260 {
261  return this->cellSetProperty
262  (
263  &MixtureType::cellReactants,
264  &MixtureType::thermoMixtureType::HE,
265  cells,
266  UIndirectList<scalar>(this->p_, cells),
267  Tu
268  );
269 }
270 
271 
272 template<class BasicPsiThermo, class MixtureType>
275 (
276  const scalarField& Tu,
277  const label patchi
278 ) const
279 {
280  return this->patchFieldProperty
281  (
282  &MixtureType::patchFaceReactants,
283  &MixtureType::thermoMixtureType::HE,
284  patchi,
285  this->p_.boundaryField()[patchi],
286  Tu
287  );
288 }
289 
290 
291 template<class BasicPsiThermo, class MixtureType>
294 {
295  return this->volScalarFieldProperty
296  (
297  "Tb",
299  &MixtureType::cellProducts,
300  &MixtureType::patchFaceProducts,
301  &MixtureType::thermoMixtureType::THE,
302  this->he_,
303  this->p_,
304  this->T_
305  );
306 }
307 
308 
309 template<class BasicPsiThermo, class MixtureType>
312 {
313  return this->volScalarFieldProperty
314  (
315  "psiu",
316  this->psi_.dimensions(),
317  &MixtureType::cellReactants,
318  &MixtureType::patchFaceReactants,
320  this->p_,
321  this->Tu_
322  );
323 }
324 
325 
326 template<class BasicPsiThermo, class MixtureType>
329 {
330  const volScalarField Tb(this->Tb());
331 
332  return this->volScalarFieldProperty
333  (
334  "psib",
335  this->psi_.dimensions(),
336  &MixtureType::cellProducts,
337  &MixtureType::patchFaceProducts,
339  this->p_,
340  Tb
341  );
342 }
343 
344 
345 template<class BasicPsiThermo, class MixtureType>
348 {
349  return this->volScalarFieldProperty
350  (
351  "muu",
353  &MixtureType::cellReactants,
354  &MixtureType::patchFaceReactants,
356  this->p_,
357  this->Tu_
358  );
359 }
360 
361 
362 template<class BasicPsiThermo, class MixtureType>
365 {
366  const volScalarField Tb(this->Tb());
367 
368  return this->volScalarFieldProperty
369  (
370  "mub",
372  &MixtureType::cellProducts,
373  &MixtureType::patchFaceProducts,
375  this->p_,
376  Tb
377  );
378 }
379 
380 
381 // ************************************************************************* //
virtual tmp< volScalarField > mub() const
Dynamic viscosity of burnt gas [kg/m/s].
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual tmp< volScalarField > psiu() const
Unburnt gas compressibility [s^2/m^2].
virtual tmp< volScalarField > Tb() const
Burnt gas temperature [K].
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
virtual tmp< volScalarField > psib() const
Burnt gas compressibility [s^2/m^2].
fvMesh & mesh
virtual tmp< volScalarField > muu() const
Dynamic viscosity of unburnt gas [kg/m/s].
label calculate(const fvMesh &mesh, const labelHashSet &patchIDs, const scalar minFaceFraction, GeometricField< scalar, PatchField, GeoMesh > &distance)
Calculate distance data from patches.
fvPatchField< scalar > fvPatchScalarField
A class for handling words, derived from string.
Definition: word.H:59
const volScalarField & psi
const dimensionSet dimDynamicViscosity
virtual void correct()
Update properties.
volScalarField scalarField(fieldObject, mesh)
virtual ~heheuPsiThermo()
Destructor.
const dimensionedScalar mu
Atomic mass unit.
thermo he()
const dimensionSet dimEnergy
const dimensionSet dimMass
label patchi
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
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:98
virtual volScalarField & heu()
Unburnt gas enthalpy [J/kg].
heheuPsiThermo(const fvMesh &, const word &phaseName)
Construct from mesh and phase name.
const dimensionSet dimTemperature
#define InfoInFunction
Report an information message using Foam::Info.