All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
multiphaseMixtureThermo.H
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) 2013-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 Class
25  Foam::multiphaseMixtureThermo
26 
27 Description
28 
29 SourceFiles
30  multiphaseMixtureThermo.C
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #ifndef multiphaseMixtureThermo_H
35 #define multiphaseMixtureThermo_H
36 
37 #include "phaseModel.H"
38 #include "PtrDictionary.H"
39 #include "volFields.H"
40 #include "surfaceFields.H"
41 #include "rhoThermo.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 /*---------------------------------------------------------------------------*\
49  Class multiphaseMixtureThermo Declaration
50 \*---------------------------------------------------------------------------*/
51 
53 :
55 {
56 public:
57 
58  class interfacePair
59  :
60  public Pair<word>
61  {
62  public:
63 
64  class hash
65  :
66  public Hash<interfacePair>
67  {
68  public:
69 
70  hash()
71  {}
72 
73  label operator()(const interfacePair& key) const
74  {
75  return word::hash()(key.first()) + word::hash()(key.second());
76  }
77  };
78 
79 
80  // Constructors
81 
83  {}
84 
85  interfacePair(const word& alpha1Name, const word& alpha2Name)
86  :
87  Pair<word>(alpha1Name, alpha2Name)
88  {}
89 
91  :
92  Pair<word>(alpha1.name(), alpha2.name())
93  {}
94 
95 
96  // Friend Operators
97 
98  friend bool operator==
99  (
100  const interfacePair& a,
101  const interfacePair& b
102  )
103  {
104  return
105  (
106  ((a.first() == b.first()) && (a.second() == b.second()))
107  || ((a.first() == b.second()) && (a.second() == b.first()))
108  );
109  }
110 
111  friend bool operator!=
112  (
113  const interfacePair& a,
114  const interfacePair& b
115  )
116  {
117  return (!(a == b));
118  }
119  };
120 
121 
122 private:
123 
124  // Private Data
125 
126  //- Dictionary of phases
128 
129  const fvMesh& mesh_;
130  const volVectorField& U_;
131  const surfaceScalarField& phi_;
132 
133  surfaceScalarField rhoPhi_;
134 
135  volScalarField alphas_;
136 
138  sigmaTable;
139 
140  sigmaTable sigmas_;
141  dimensionSet dimSigma_;
142 
143  //- Stabilisation for normalisation of the interface normal
144  const dimensionedScalar deltaN_;
145 
146 
147  // Private Member Functions
148 
149  void calcAlphas();
150 
151  void solveAlphas(const scalar cAlpha);
152 
154  (
155  const volScalarField& alpha1,
156  const volScalarField& alpha2
157  ) const;
158 
160  (
161  const volScalarField& alpha1,
162  const volScalarField& alpha2
163  ) const;
164 
165  void correctContactAngle
166  (
167  const phaseModel& alpha1,
168  const phaseModel& alpha2,
169  surfaceVectorField::Boundary& nHatb
170  ) const;
171 
173  (
174  const phaseModel& alpha1,
175  const phaseModel& alpha2
176  ) const;
177 
178 
179 public:
180 
181  //- Runtime type information
182  TypeName("multiphaseMixtureThermo");
183 
184 
185  // Constructors
186 
187  //- Construct from components
189  (
190  const volVectorField& U,
191  const surfaceScalarField& phi
192  );
193 
194 
195  //- Destructor
196  virtual ~multiphaseMixtureThermo()
197  {}
198 
199 
200  // Member Functions
201 
202  //- Return the phases
203  const PtrDictionary<phaseModel>& phases() const
204  {
205  return phases_;
206  }
207 
208  //- Return non-const access to the phases
210  {
211  return phases_;
212  }
213 
214  //- Return the velocity
215  const volVectorField& U() const
216  {
217  return U_;
218  }
219 
220  //- Return the volumetric flux
221  const surfaceScalarField& phi() const
222  {
223  return phi_;
224  }
226  const surfaceScalarField& rhoPhi() const
227  {
228  return rhoPhi_;
229  }
230 
231  //- Correct the thermodynamics of each phase
232  virtual void correctThermo();
233 
234  //- Update properties
235  virtual void correct();
236 
237  //- Update densities for given pressure change
238  void correctRho(const volScalarField& dp);
239 
240  //- Return the name of the thermo physics
241  virtual word thermoName() const;
242 
243  //- Return true if the equation of state is incompressible
244  // i.e. rho != f(p)
245  virtual bool incompressible() const;
246 
247  //- Return true if the equation of state is isochoric
248  // i.e. rho = const
249  virtual bool isochoric() const;
250 
251 
252  // Access to thermodynamic state variables
253 
254  //- Enthalpy/Internal energy [J/kg]
255  // Non-const access allowed for transport equations
256  virtual volScalarField& he()
257  {
259  return phases_[0].thermo().he();
260  }
261 
262  //- Enthalpy/Internal energy [J/kg]
263  virtual const volScalarField& he() const
264  {
266  return phases_[0].thermo().he();
267  }
268 
269  //- Enthalpy/Internal energy
270  // for given pressure and temperature [J/kg]
271  virtual tmp<volScalarField> he
272  (
273  const volScalarField& p,
274  const volScalarField& T
275  ) const;
276 
277  //- Enthalpy/Internal energy for cell-set [J/kg]
278  virtual tmp<scalarField> he
279  (
280  const scalarField& T,
281  const labelList& cells
282  ) const;
283 
284  //- Enthalpy/Internal energy for patch [J/kg]
285  virtual tmp<scalarField> he
286  (
287  const scalarField& T,
288  const label patchi
289  ) const;
290 
291  //- Sensible enthalpy [J/kg]
292  virtual tmp<volScalarField> hs() const;
293 
294  //- Sensible enthalpy
295  // for given pressure and temperature [J/kg]
296  virtual tmp<volScalarField> hs
297  (
298  const volScalarField& p,
299  const volScalarField& T
300  ) const;
301 
302  //- Sensible enthalpy for cell-set [J/kg]
303  virtual tmp<scalarField> hs
304  (
305  const scalarField& T,
306  const labelList& cells
307  ) const;
308 
309  //- Sensible enthalpy for patch [J/kg]
310  virtual tmp<scalarField> hs
311  (
312  const scalarField& T,
313  const label patchi
314  ) const;
315 
316  //- Absolute enthalpy [J/kg]
317  virtual tmp<volScalarField> ha() const;
318 
319  //- Absolute enthalpy
320  // for given pressure and temperature [J/kg]
321  virtual tmp<volScalarField> ha
322  (
323  const volScalarField& p,
324  const volScalarField& T
325  ) const;
326 
327  //- Absolute enthalpy for cell-set [J/kg]
328  virtual tmp<scalarField> ha
329  (
330  const scalarField& T,
331  const labelList& cells
332  ) const;
333 
334  //- Absolute enthalpy for patch [J/kg]
335  virtual tmp<scalarField> ha
336  (
337  const scalarField& T,
338  const label patchi
339  ) const;
340 
341  //- Enthalpy of formation [J/kg]
342  virtual tmp<volScalarField> hc() const;
343 
344  //- Temperature from enthalpy/internal energy
345  virtual tmp<volScalarField> THE
346  (
347  const volScalarField& h,
348  const volScalarField& p,
349  const volScalarField& T0 // starting temperature
350  ) const;
351 
352  //- Temperature from enthalpy/internal energy for cell-set
353  virtual tmp<scalarField> THE
354  (
355  const scalarField& h,
356  const scalarField& T0, // starting temperature
357  const labelList& cells
358  ) const;
359 
360  //- Temperature from enthalpy/internal energy for patch
361  virtual tmp<scalarField> THE
362  (
363  const scalarField& h,
364  const scalarField& T0, // starting temperature
365  const label patchi
366  ) const;
367 
368 
369  // Fields derived from thermodynamic state variables
370 
371  //- Heat capacity at constant pressure [J/kg/K]
372  virtual tmp<volScalarField> Cp() const;
373 
374  //- Heat capacity at constant pressure for patch [J/kg/K]
375  virtual tmp<scalarField> Cp
376  (
377  const scalarField& T,
378  const label patchi
379  ) const;
380 
381  //- Heat capacity at constant volume [J/kg/K]
382  virtual tmp<volScalarField> Cv() const;
383 
384  //- Heat capacity at constant volume for patch [J/kg/K]
385  virtual tmp<scalarField> Cv
386  (
387  const scalarField& T,
388  const label patchi
389  ) const;
390 
391  //- Gamma = Cp/Cv []
392  virtual tmp<volScalarField> gamma() const;
393 
394  //- Gamma = Cp/Cv for patch []
395  virtual tmp<scalarField> gamma
396  (
397  const scalarField& T,
398  const label patchi
399  ) const;
400 
401  //- Heat capacity at constant pressure/volume [J/kg/K]
402  virtual tmp<volScalarField> Cpv() const;
403 
404  //- Heat capacity at constant pressure/volume for patch [J/kg/K]
405  virtual tmp<scalarField> Cpv
406  (
407  const scalarField& T,
408  const label patchi
409  ) const;
410 
411  //- Molecular weight [kg/kmol]
412  virtual tmp<volScalarField> W() const;
413 
414  //- Molecular weight for patch [kg/kmol]
415  virtual tmp<scalarField> W(const label patchi) const;
416 
417 
418  // Fields derived from transport state variables
419 
420  //- Kinematic viscosity of mixture [m^2/s]
421  virtual tmp<volScalarField> nu() const;
422 
423  //- Kinematic viscosity of mixture for patch [m^2/s]
424  virtual tmp<scalarField> nu(const label patchi) const;
425 
426  //- Thermal diffusivity for temperature of mixture [W/m/K]
427  virtual tmp<volScalarField> kappa() const;
428 
429  //- Thermal diffusivity of mixture for patch [W/m/K]
430  virtual tmp<scalarField> kappa
431  (
432  const label patchi
433  ) const;
434 
435  //- Thermal diffusivity for energy of mixture [kg/m/s]
436  virtual tmp<volScalarField> alphahe() const;
437 
438  //- Thermal diffusivity for energy of mixture for patch [kg/m/s]
439  virtual tmp<scalarField> alphahe(const label patchi) const;
440 
441  //- Effective thermal diffusivity of mixture [W/m/K]
443  (
444  const volScalarField& alphat
445  ) const;
446 
447  //- Effective thermal diffusivity of mixture for patch [W/m/K]
448  virtual tmp<scalarField> kappaEff
449  (
450  const scalarField& alphat,
451  const label patchi
452  ) const;
453 
454  //- Effective thermal diffusivity of mixture [W/m/K]
456  (
457  const volScalarField& alphat
458  ) const;
459 
460  //- Effective thermal diffusivity of mixture for patch [W/m/K]
461  virtual tmp<scalarField> alphaEff
462  (
463  const scalarField& alphat,
464  const label patchi
465  ) const;
466 
467 
468  //- Return the phase-averaged reciprocal Cv
469  tmp<volScalarField> rCv() const;
470 
472 
473  //- Indicator of the proximity of the interface
474  // Field values are 1 near and 0 away for the interface.
476 
477  //- Solve for the mixture phase-fractions
478  void solve();
479 };
480 
481 
482 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
483 
484 } // End namespace Foam
485 
486 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
487 
488 #endif
489 
490 // ************************************************************************* //
Template dictionary class which manages the storage associated with it.
Definition: PtrDictionary.H:53
Foam::surfaceFields.
const volVectorField & U() const
Return the velocity.
Hashing function class, shared by all the derived classes.
Definition: string.H:92
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 surfaceScalarField & phi() const
Return the volumetric flux.
virtual tmp< volScalarField > gamma() const
Gamma = Cp/Cv [].
virtual tmp< volScalarField > hc() const
Enthalpy of formation [J/kg].
virtual tmp< volScalarField > ha() const
Absolute enthalpy [J/kg].
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
volScalarField & alpha1(mixture.alpha1())
virtual ~multiphaseMixtureThermo()
Destructor.
const dimensionedScalar h
Planck constant.
CGAL::Exact_predicates_exact_constructions_kernel K
virtual tmp< volScalarField > Cp() const
Heat capacity at constant pressure [J/kg/K].
alpha2
Definition: alphaEqn.H:115
void correctRho(const volScalarField &dp)
Update densities for given pressure change.
void solve()
Solve for the mixture phase-fractions.
Dimension set for the base types.
Definition: dimensionSet.H:120
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
const fileName & name() const
Return the dictionary name.
Definition: dictionary.H:109
const cellShapeList & cells
virtual tmp< volScalarField > alphahe() const
Thermal diffusivity for energy of mixture [kg/m/s].
tmp< volScalarField > rCv() const
Return the phase-averaged reciprocal Cv.
const PtrDictionary< phaseModel > & phases() const
Return the phases.
A class for handling words, derived from string.
Definition: word.H:59
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:449
virtual tmp< volScalarField > nu() const
Kinematic viscosity of mixture [m^2/s].
label operator()(const interfacePair &key) const
const Type & second() const
Return second.
Definition: Pair.H:110
TypeName("multiphaseMixtureThermo")
Runtime type information.
virtual tmp< volScalarField > kappaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [W/m/K].
multiphaseMixtureThermo(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
label patchi
virtual tmp< volScalarField > W() const
Molecular weight [kg/kmol].
const surfaceScalarField & rhoPhi() const
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Hash function class for primitives. All non-primitives used to hash entries on hash tables likely nee...
Definition: Hash.H:52
virtual tmp< volScalarField > hs() const
Sensible enthalpy [J/kg].
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
virtual tmp< volScalarField > Cv() const
Heat capacity at constant volume [J/kg/K].
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< surfaceScalarField > surfaceTensionForce() const
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [W/m/K].
virtual tmp< volScalarField > Cpv() const
Heat capacity at constant pressure/volume [J/kg/K].
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:370
virtual tmp< volScalarField > kappa() const
Thermal diffusivity for temperature of mixture [W/m/K].
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
const Type & first() const
Return first.
Definition: Pair.H:98
virtual bool isochoric() const
Return true if the equation of state is isochoric.
virtual word thermoName() const
Return the name of the thermo physics.
virtual volScalarField & p()
Pressure [Pa].
Definition: fluidThermo.C:97
virtual void correctThermo()
Correct the thermodynamics of each phase.
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
Namespace for OpenFOAM.
virtual tmp< volScalarField > THE(const volScalarField &h, const volScalarField &p, const volScalarField &T0) const
Temperature from enthalpy/internal energy.
scalar T0
Definition: createFields.H:22
virtual bool incompressible() const
Return true if the equation of state is incompressible.
virtual void correct()
Update properties.