multiphaseMixtureThermo.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2013-2017 OpenFOAM Foundation
9 Copyright (C) 2019-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Class
28 Foam::multiphaseMixtureThermo
29
30Description
31
32SourceFiles
33 multiphaseMixtureThermo.C
34
35\*---------------------------------------------------------------------------*/
36
37#ifndef multiphaseMixtureThermo_H
38#define multiphaseMixtureThermo_H
39
40#include "phaseModel.H"
41#include "PtrDictionary.H"
42#include "volFields.H"
43#include "surfaceFields.H"
44#include "rhoThermo.H"
45#include "psiThermo.H"
46
47// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48
49namespace Foam
50{
51
52/*---------------------------------------------------------------------------*\
53 Class multiphaseMixtureThermo Declaration
54\*---------------------------------------------------------------------------*/
57:
58 public psiThermo
59{
60public:
61
62 //- Symmetric pair of interface names
63 class interfacePair
64 :
65 public Pair<word>
66 {
67 public:
68
69 // Always use symmetric hashing
71
72 // Always use symmetric hashing (alias)
74
75
76 // Constructors
78 interfacePair() = default;
80 interfacePair(const word& alpha1Name, const word& alpha2Name)
81 :
82 Pair<word>(alpha1Name, alpha2Name)
83 {}
86 :
88 {}
89
90
91 // Friend Operators
93 friend bool operator==
94 (
95 const interfacePair& a,
96 const interfacePair& b
97 )
98 {
99 return (0 != Pair<word>::compare(a, b));
100 }
102 friend bool operator!=
103 (
104 const interfacePair& a,
105 const interfacePair& b
106 )
107 {
108 return (!(a == b));
109 }
110 };
111
112
113private:
114
115 // Private Data
116
117 //- Dictionary of phases
119
120 const fvMesh& mesh_;
121 const volVectorField& U_;
122 const surfaceScalarField& phi_;
123
124 surfaceScalarField rhoPhi_;
125
126 volScalarField alphas_;
127
129 sigmaTable;
130
131 sigmaTable sigmas_;
132 dimensionSet dimSigma_;
133
134 //- Stabilisation for normalisation of the interface normal
135 const dimensionedScalar deltaN_;
136
137
138 // Private Member Functions
139
140 void calcAlphas();
141
142 void solveAlphas(const scalar cAlpha);
143
145 (
146 const volScalarField& alpha1,
148 ) const;
149
151 (
152 const volScalarField& alpha1,
154 ) const;
155
156 void correctContactAngle
157 (
158 const phaseModel& alpha1,
159 const phaseModel& alpha2,
161 ) const;
162
164 (
165 const phaseModel& alpha1,
166 const phaseModel& alpha2
167 ) const;
168
169
170public:
171
172 //- Runtime type information
173 TypeName("multiphaseMixtureThermo");
174
175
176 // Constructors
177
178 //- Construct from components
180 (
181 const volVectorField& U,
183 );
184
185
186 //- Destructor
187 virtual ~multiphaseMixtureThermo() = default;
188
189
190 // Member Functions
191
192 //- Return the phases
193 const PtrDictionary<phaseModel>& phases() const
194 {
195 return phases_;
196 }
197
198 //- Return non-const access to the phases
200 {
201 return phases_;
202 }
203
204 //- Return the velocity
205 const volVectorField& U() const
206 {
207 return U_;
208 }
209
210 //- Return the volumetric flux
211 const surfaceScalarField& phi() const
212 {
213 return phi_;
214 }
216 const surfaceScalarField& rhoPhi() const
217 {
218 return rhoPhi_;
219 }
220
221 //- Update properties
222 virtual void correct();
223
224 //- Update densities for given pressure change
225 void correctRho(const volScalarField& dp);
226
227 //- Return the name of the thermo physics
228 virtual word thermoName() const;
229
230 //- Return true if the equation of state is incompressible
231 // i.e. rho != f(p)
232 virtual bool incompressible() const;
233
234 //- Return true if the equation of state is isochoric
235 // i.e. rho = const
236 virtual bool isochoric() const;
237
238
239 // Access to thermodynamic state variables
240
241 //- Enthalpy/Internal energy [J/kg]
242 // Non-const access allowed for transport equations
243 virtual volScalarField& he()
244 {
246 return phases_[0].thermo().he();
247 }
248
249 //- Enthalpy/Internal energy [J/kg]
250 virtual const volScalarField& he() const
251 {
253 return phases_[0].thermo().he();
254 }
255
256 //- Enthalpy/Internal energy
257 // for given pressure and temperature [J/kg]
258 virtual tmp<volScalarField> he
259 (
260 const volScalarField& p,
261 const volScalarField& T
262 ) const;
263
264 //- Enthalpy/Internal energy for cell-set [J/kg]
265 virtual tmp<scalarField> he
266 (
267 const scalarField& p,
268 const scalarField& T,
269 const labelList& cells
270 ) const;
271
272 //- Enthalpy/Internal energy for patch [J/kg]
273 virtual tmp<scalarField> he
274 (
275 const scalarField& p,
276 const scalarField& T,
277 const label patchi
278 ) const;
279
280 //- Chemical enthalpy [J/kg]
281 virtual tmp<volScalarField> hc() const;
282
283 //- Temperature from enthalpy/internal energy for cell-set
284 virtual tmp<scalarField> THE
285 (
286 const scalarField& h,
287 const scalarField& p,
288 const scalarField& T0, // starting temperature
289 const labelList& cells
290 ) const;
291
292 //- Temperature from enthalpy/internal energy for patch
293 virtual tmp<scalarField> THE
294 (
295 const scalarField& h,
296 const scalarField& p,
297 const scalarField& T0, // starting temperature
298 const label patchi
299 ) const;
300
301
302 // Fields derived from thermodynamic state variables
303
304 //- Density [kg/m^3]
305 virtual tmp<volScalarField> rho() const;
306
307 //- Density for patch [kg/m^3]
308 virtual tmp<scalarField> rho(const label patchi) const;
309
310 //- Heat capacity at constant pressure [J/kg/K]
311 virtual tmp<volScalarField> Cp() const;
312
313 //- Heat capacity at constant pressure for patch [J/kg/K]
314 virtual tmp<scalarField> Cp
315 (
316 const scalarField& p,
317 const scalarField& T,
318 const label patchi
319 ) const;
320
321 //- Heat capacity using pressure and temperature
322 virtual tmp<scalarField> Cp
323 (
324 const scalarField& p,
325 const scalarField& T,
326 const labelList& cells
327 ) const
328 {
330 return tmp<scalarField>::New(p);
331 }
332
333 //- Heat capacity at constant volume [J/kg/K]
334 virtual tmp<volScalarField> Cv() const;
335
336 //- Heat capacity at constant volume for patch [J/kg/K]
337 virtual tmp<scalarField> Cv
338 (
339 const scalarField& p,
340 const scalarField& T,
341 const label patchi
342 ) const;
343
344 //- Density from pressure and temperature
346 (
347 const scalarField& p,
348 const scalarField& T,
349 const labelList& cells
350 ) const
351 {
353 return tmp<scalarField>::New(p);
354 }
355
356 //- Gamma = Cp/Cv []
357 virtual tmp<volScalarField> gamma() const;
358
359 //- Gamma = Cp/Cv for patch []
360 virtual tmp<scalarField> gamma
361 (
362 const scalarField& p,
363 const scalarField& T,
364 const label patchi
365 ) const;
366
367 //- Heat capacity at constant pressure/volume [J/kg/K]
368 virtual tmp<volScalarField> Cpv() const;
369
370 //- Heat capacity at constant pressure/volume for patch [J/kg/K]
371 virtual tmp<scalarField> Cpv
372 (
373 const scalarField& p,
374 const scalarField& T,
375 const label patchi
376 ) const;
377
378 //- Heat capacity ratio []
379 virtual tmp<volScalarField> CpByCpv() const;
380
381 //- Heat capacity ratio for patch []
383 (
384 const scalarField& p,
385 const scalarField& T,
386 const label patchi
387 ) const;
388
389 //- Molecular weight [kg/kmol]
390 virtual tmp<volScalarField> W() const;
391
392
393 // Fields derived from transport state variables
394
395 //- Kinematic viscosity of mixture [m^2/s]
396 virtual tmp<volScalarField> nu() const;
397
398 //- Kinematic viscosity of mixture for patch [m^2/s]
399 virtual tmp<scalarField> nu(const label patchi) const;
400
401 //- Thermal diffusivity for temperature of mixture [J/m/s/K]
402 virtual tmp<volScalarField> kappa() const;
403
404 //- Thermal diffusivity of mixture for patch [J/m/s/K]
405 virtual tmp<scalarField> kappa
406 (
407 const label patchi
408 ) const;
409
410
411 //- Thermal diffusivity for energy of mixture [kg/m/s]
412 virtual tmp<volScalarField> alphahe() const;
413
414 //- Thermal diffusivity for energy of mixture for patch [kg/m/s]
415 virtual tmp<scalarField> alphahe(const label patchi) const;
416
417 //- Effective thermal diffusivity of mixture [J/m/s/K]
419 (
420 const volScalarField& alphat
421 ) const;
422
423 //- Effective thermal diffusivity of mixture for patch [J/m/s/K]
425 (
426 const scalarField& alphat,
427 const label patchi
428 ) const;
429
430 //- Effective thermal diffusivity of mixture [J/m/s/K]
432 (
433 const volScalarField& alphat
434 ) const;
435
436 //- Effective thermal diffusivity of mixture for patch [J/m/s/K]
438 (
439 const scalarField& alphat,
440 const label patchi
441 ) const;
442
443
444 //- Return the phase-averaged reciprocal Cv
445 tmp<volScalarField> rCv() const;
448
449 //- Indicator of the proximity of the interface
450 // Field values are 1 near and 0 away for the interface.
452
453 //- Solve for the mixture phase-fractions
454 void solve();
455};
456
457
458// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
459
460} // End namespace Foam
461
462// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
463
464#endif
465
466// ************************************************************************* //
CGAL::Exact_predicates_exact_constructions_kernel K
const volScalarField & alpha1
const volScalarField & alpha2
A HashTable similar to std::unordered_map.
Definition: HashTable.H:123
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: Pair.H:69
Template dictionary class which manages the storage associated with it.
Definition: PtrDictionary.H:55
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
const word & name() const
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:622
virtual volScalarField & p()
Pressure [Pa].
Definition: basicThermo.C:610
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Definition: dimensionSet.H:109
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
interfacePair(const word &alpha1Name, const word &alpha2Name)
interfacePair(const phaseModel &alpha1, const phaseModel &alpha2)
const volVectorField & U() const
Return the velocity.
virtual ~multiphaseMixtureThermo()=default
Destructor.
virtual tmp< scalarField > he(const scalarField &p, const scalarField &T, const label patchi) const
Enthalpy/Internal energy for patch [J/kg].
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
virtual word thermoName() const
Return the name of the thermo physics.
tmp< surfaceScalarField > surfaceTensionForce() const
virtual tmp< volScalarField > Cv() const
Heat capacity at constant volume [J/kg/K].
virtual tmp< volScalarField > rho() const
Density [kg/m^3].
virtual tmp< volScalarField > W() const
Molecular weight [kg/kmol].
virtual tmp< volScalarField > gamma() const
Gamma = Cp/Cv [].
virtual tmp< scalarField > kappa(const label patchi) const
Thermal diffusivity of mixture for patch [J/m/s/K].
virtual bool incompressible() const
Return true if the equation of state is incompressible.
virtual tmp< volScalarField > Cpv() const
Heat capacity at constant pressure/volume [J/kg/K].
PtrDictionary< phaseModel > & phases()
Return non-const access to the phases.
virtual tmp< scalarField > CpByCpv(const scalarField &p, const scalarField &T, const label patchi) const
Heat capacity ratio for patch [].
virtual tmp< scalarField > alphahe(const label patchi) const
Thermal diffusivity for energy of mixture for patch [kg/m/s].
virtual tmp< volScalarField > kappaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
virtual tmp< scalarField > rho(const label patchi) const
Density for patch [kg/m^3].
virtual tmp< scalarField > THE(const scalarField &h, const scalarField &p, const scalarField &T0, const labelList &cells) const
Temperature from enthalpy/internal energy for cell-set.
const surfaceScalarField & phi() const
Return the volumetric flux.
virtual tmp< scalarField > Cv(const scalarField &p, const scalarField &T, const label patchi) const
Heat capacity at constant volume for patch [J/kg/K].
TypeName("multiphaseMixtureThermo")
Runtime type information.
virtual tmp< volScalarField > kappa() const
Thermal diffusivity for temperature of mixture [J/m/s/K].
virtual tmp< scalarField > he(const scalarField &p, const scalarField &T, const labelList &cells) const
Enthalpy/Internal energy for cell-set [J/kg].
virtual tmp< scalarField > Cp(const scalarField &p, const scalarField &T, const labelList &cells) const
Heat capacity using pressure and temperature.
virtual tmp< scalarField > Cpv(const scalarField &p, const scalarField &T, const label patchi) const
Heat capacity at constant pressure/volume for patch [J/kg/K].
const surfaceScalarField & rhoPhi() const
tmp< volScalarField > rCv() const
Return the phase-averaged reciprocal Cv.
virtual tmp< scalarField > alphaEff(const scalarField &alphat, const label patchi) const
Effective thermal diffusivity of mixture for patch [J/m/s/K].
const PtrDictionary< phaseModel > & phases() const
Return the phases.
virtual tmp< scalarField > rhoEoS(const scalarField &p, const scalarField &T, const labelList &cells) const
Density from pressure and temperature.
virtual tmp< volScalarField > Cp() const
Heat capacity at constant pressure [J/kg/K].
virtual tmp< scalarField > kappaEff(const scalarField &alphat, const label patchi) const
Effective thermal diffusivity of mixture for patch [J/m/s/K].
virtual void correct()
Update properties.
multiphaseMixtureThermo(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
virtual tmp< scalarField > Cp(const scalarField &p, const scalarField &T, const label patchi) const
Heat capacity at constant pressure for patch [J/kg/K].
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
virtual tmp< volScalarField > hc() const
Chemical enthalpy [J/kg].
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
virtual tmp< volScalarField > nu() const
Kinematic viscosity of mixture [m^2/s].
virtual tmp< scalarField > gamma(const scalarField &p, const scalarField &T, const label patchi) const
Gamma = Cp/Cv for patch [].
virtual tmp< scalarField > THE(const scalarField &h, const scalarField &p, const scalarField &T0, const label patchi) const
Temperature from enthalpy/internal energy for patch.
virtual const volScalarField & he() const
Enthalpy/Internal energy [J/kg].
void correctRho(const volScalarField &dp)
Update densities for given pressure change.
virtual tmp< scalarField > nu(const label patchi) const
Kinematic viscosity of mixture for patch [m^2/s].
virtual tmp< volScalarField > alphahe() const
Thermal diffusivity for energy of mixture [kg/m/s].
virtual tmp< volScalarField > CpByCpv() const
Heat capacity ratio [].
void solve()
Solve for the mixture phase-fractions.
virtual bool isochoric() const
Return true if the equation of state is isochoric.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:61
Basic thermodynamic properties based on compressibility.
Definition: psiThermo.H:58
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:517
const cellShapeList & cells
kappaEff
Definition: TEqn.H:10
Namespace for OpenFOAM.
volScalarField & h
volScalarField & b
Definition: createFields.H:27
scalar T0
Definition: createFields.H:22
Symmetric hashing functor for Pair, hashes lower value first.
Definition: Pair.H:152
Foam::surfaceFields.
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73