SpecieMixture.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-2020 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 "SpecieMixture.H"
27 #include "fvMesh.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class MixtureType>
34 (
35  const word& psiName,
36  const dimensionSet& psiDim,
37  scalar (MixtureType::thermoType::*psiMethod)
38  (
39  const scalar,
40  const scalar
41  ) const,
42  const label speciei,
43  const volScalarField& p,
44  const volScalarField& T
45 ) const
46 {
47  const typename MixtureType::thermoType& thermo =
48  this->specieThermo(speciei);
49 
50  tmp<volScalarField> tPsi
51  (
53  (
54  IOobject::groupName(psiName, T.group()),
55  T.mesh(),
56  psiDim
57  )
58  );
59 
60  volScalarField& psi = tPsi.ref();
61 
62  forAll(p, celli)
63  {
64  psi[celli] = (thermo.*psiMethod)(p[celli], T[celli]);
65  }
66 
67  volScalarField::Boundary& psiBf = psi.boundaryFieldRef();
68 
69  forAll(psiBf, patchi)
70  {
71  const fvPatchScalarField& pp = p.boundaryField()[patchi];
72  const fvPatchScalarField& pT = T.boundaryField()[patchi];
73  fvPatchScalarField& ppsi = psiBf[patchi];
74 
75  forAll(pp, facei)
76  {
77  ppsi[facei] = (thermo.*psiMethod)(pp[facei], pT[facei]);
78  }
79  }
80 
81  return tPsi;
82 }
83 
84 
85 template<class MixtureType>
87 (
88  scalar (MixtureType::thermoType::*psiMethod)
89  (
90  const scalar,
91  const scalar
92  ) const,
93  const label speciei,
94  const scalarField& p,
95  const scalarField& T
96 ) const
97 {
98  const typename MixtureType::thermoType& thermo =
99  this->specieThermo(speciei);
100 
101  tmp<scalarField> tPsi(new scalarField(p.size()));
102 
103  scalarField& psi = tPsi.ref();
104 
105  forAll(p, facei)
106  {
107  psi[facei] = (thermo.*psiMethod)(p[facei], T[facei]);
108  }
109 
110  return tPsi;
111 }
112 
113 
114 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
115 
116 template<class MixtureType>
118 (
119  const dictionary& thermoDict,
120  const fvMesh& mesh,
121  const word& phaseName
122 )
123 :
124  MixtureType(thermoDict, mesh, phaseName)
125 {}
126 
127 
128 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
129 
130 template<class MixtureType>
131 Foam::scalar Foam::SpecieMixture<MixtureType>::Wi(const label speciei) const
132 {
133  return this->specieThermo(speciei).W();
134 }
135 
136 
137 template<class MixtureType>
138 Foam::scalar Foam::SpecieMixture<MixtureType>::Hf(const label speciei) const
139 {
140  return this->specieThermo(speciei).Hf();
141 }
142 
143 
144 template<class MixtureType>
146 (
147  const label speciei,
148  const scalar p,
149  const scalar T
150 ) const
151 {
152  return this->specieThermo(speciei).rho(p, T);
153 }
154 
155 
156 template<class MixtureType>
158 (
159  const label speciei,
160  const volScalarField& p,
161  const volScalarField& T
162 ) const
163 {
164  return volScalarFieldProperty
165  (
166  "rho",
167  dimDensity,
169  speciei,
170  p,
171  T
172  );
173 }
174 
175 
176 template<class MixtureType>
178 (
179  const label speciei,
180  const scalar p,
181  const scalar T
182 ) const
183 {
184  return this->specieThermo(speciei).Cp(p, T);
185 }
186 
187 
188 template<class MixtureType>
190 (
191  const label speciei,
192  const volScalarField& p,
193  const volScalarField& T
194 ) const
195 {
196  return volScalarFieldProperty
197  (
198  "Cp",
201  speciei,
202  p,
203  T
204  );
205 }
206 
207 
208 template<class MixtureType>
210 (
211  const label speciei,
212  const scalar p,
213  const scalar T
214 ) const
215 {
216  return this->specieThermo(speciei).HE(p, T);
217 }
218 
219 
220 template<class MixtureType>
222 (
223  const label speciei,
224  const scalarField& p,
225  const scalarField& T
226 ) const
227 {
228  return fieldProperty(&MixtureType::thermoType::HE, speciei, p, T);
229 }
230 
231 
232 template<class MixtureType>
234 (
235  const label speciei,
236  const volScalarField& p,
237  const volScalarField& T
238 ) const
239 {
240  return volScalarFieldProperty
241  (
242  "HE",
244  &MixtureType::thermoType::HE,
245  speciei,
246  p,
247  T
248  );
249 }
250 
251 
252 template<class MixtureType>
254 (
255  const label speciei,
256  const scalar p,
257  const scalar T
258 ) const
259 {
260  return this->specieThermo(speciei).Hs(p, T);
261 }
262 
263 
264 template<class MixtureType>
266 (
267  const label speciei,
268  const scalarField& p,
269  const scalarField& T
270 ) const
271 {
272  return fieldProperty(&MixtureType::thermoType::Hs, speciei, p, T);
273 }
274 
275 
276 template<class MixtureType>
278 (
279  const label speciei,
280  const volScalarField& p,
281  const volScalarField& T
282 ) const
283 {
284  return volScalarFieldProperty
285  (
286  "Hs",
289  speciei,
290  p,
291  T
292  );
293 }
294 
295 
296 template<class MixtureType>
298 (
299  const label speciei,
300  const scalar p,
301  const scalar T
302 ) const
303 {
304  return this->specieThermo(speciei).Ha(p, T);
305 }
306 
307 
308 template<class MixtureType>
310 (
311  const label speciei,
312  const scalarField& p,
313  const scalarField& T
314 ) const
315 {
316  return fieldProperty(&MixtureType::thermoType::Ha, speciei, p, T);
317 }
318 
319 
320 template<class MixtureType>
322 (
323  const label speciei,
324  const volScalarField& p,
325  const volScalarField& T
326 ) const
327 {
328  return volScalarFieldProperty
329  (
330  "Ha",
333  speciei,
334  p,
335  T
336  );
337 }
338 
339 
340 template<class MixtureType>
342 (
343  const label speciei,
344  const scalar p,
345  const scalar T
346 ) const
347 {
348  return this->specieThermo(speciei).mu(p, T);
349 }
350 
351 
352 template<class MixtureType>
354 (
355  const label speciei,
356  const volScalarField& p,
357  const volScalarField& T
358 ) const
359 {
360  return volScalarFieldProperty
361  (
362  "mu",
365  speciei,
366  p,
367  T
368  );
369 }
370 
371 
372 template<class MixtureType>
374 (
375  const label speciei,
376  const scalar p,
377  const scalar T
378 ) const
379 {
380  return this->specieThermo(speciei).kappa(p, T);
381 }
382 
383 
384 template<class MixtureType>
386 (
387  const label speciei,
388  const volScalarField& p,
389  const volScalarField& T
390 ) const
391 {
392  return volScalarFieldProperty
393  (
394  "kappa",
397  speciei,
398  p,
399  T
400  );
401 }
402 
403 
404 // ************************************************************************* //
Foam::SpecieMixture.
Definition: SpecieMixture.H:55
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual scalar Hf(const label speciei) const
Enthalpy of formation [J/kg].
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 dimensionSet dimPower
scalar Hs(const scalar p, const scalar T) const
Definition: EtoHthermo.H:11
virtual scalar HE(const label speciei, const scalar p, const scalar T) const
Enthalpy/Internal energy [J/kg].
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
virtual scalar Hs(const label speciei, const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
virtual scalar Wi(const label speciei) const
Molecular weight of the given specie [kg/kmol].
const dimensionSet dimLength
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
SpecieMixture(const dictionary &, const fvMesh &, const word &phaseName)
Construct from dictionary, mesh and phase name.
const dimensionSet dimTime
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
fvPatchField< scalar > fvPatchScalarField
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
Definition: word.H:59
const dimensionSet dimDensity
virtual scalar Ha(const label speciei, const scalar p, const scalar T) const
Absolute enthalpy [J/kg].
volScalarField scalarField(fieldObject, mesh)
const dimensionedScalar mu
Atomic mass unit.
const dimensionSet dimEnergy
const dimensionSet dimMass
Internal & ref()
Return a reference to the dimensioned internal field.
virtual scalar rho(const label speciei, const scalar p, const scalar T) const
Density [kg/m^3].
label patchi
virtual scalar mu(const label speciei, const scalar p, const scalar T) const
Dynamic viscosity [kg/m/s].
virtual tmp< volScalarField > kappa() const =0
Thermal diffusivity for temperature of mixture [W/m/K].
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
scalar Ha(const scalar p, const scalar T) const
Definition: EtoHthermo.H:20
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volScalarField > Cp() const =0
Heat capacity at constant pressure for patch [J/kg/K].
const dimensionSet dimTemperature