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 | Copyright (C) 2013-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 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 #include "psiThermo.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 
49 /*---------------------------------------------------------------------------*\
50  Class multiphaseMixtureThermo Declaration
51 \*---------------------------------------------------------------------------*/
52 
54 :
55  public psiThermo
56 {
57 public:
58 
59  class interfacePair
60  :
61  public Pair<word>
62  {
63  public:
64 
65  class hash
66  :
67  public Hash<interfacePair>
68  {
69  public:
70 
71  hash()
72  {}
73 
74  label operator()(const interfacePair& key) const
75  {
76  return word::hash()(key.first()) + word::hash()(key.second());
77  }
78  };
79 
80 
81  // Constructors
82 
84  {}
85 
86  interfacePair(const word& alpha1Name, const word& alpha2Name)
87  :
88  Pair<word>(alpha1Name, alpha2Name)
89  {}
90 
92  :
93  Pair<word>(alpha1.name(), alpha2.name())
94  {}
95 
96 
97  // Friend Operators
98 
99  friend bool operator==
100  (
101  const interfacePair& a,
102  const interfacePair& b
103  )
104  {
105  return
106  (
107  ((a.first() == b.first()) && (a.second() == b.second()))
108  || ((a.first() == b.second()) && (a.second() == b.first()))
109  );
110  }
111 
112  friend bool operator!=
113  (
114  const interfacePair& a,
115  const interfacePair& b
116  )
117  {
118  return (!(a == b));
119  }
120  };
121 
122 
123 private:
124 
125  // Private data
126 
127  //- Dictionary of phases
129 
130  const fvMesh& mesh_;
131  const volVectorField& U_;
132  const surfaceScalarField& phi_;
133 
134  surfaceScalarField rhoPhi_;
135 
136  volScalarField alphas_;
137 
139  sigmaTable;
140 
141  sigmaTable sigmas_;
142  dimensionSet dimSigma_;
143 
144  //- Stabilisation for normalisation of the interface normal
145  const dimensionedScalar deltaN_;
146 
147  //- Conversion factor for degrees into radians
148  static const scalar convertToRad;
149 
150 
151  // Private member functions
152 
153  void calcAlphas();
154 
155  void solveAlphas(const scalar cAlpha);
156 
158  (
159  const volScalarField& alpha1,
160  const volScalarField& alpha2
161  ) const;
162 
164  (
165  const volScalarField& alpha1,
166  const volScalarField& alpha2
167  ) const;
168 
169  void correctContactAngle
170  (
171  const phaseModel& alpha1,
172  const phaseModel& alpha2,
173  surfaceVectorField::Boundary& nHatb
174  ) const;
175 
177  (
178  const phaseModel& alpha1,
179  const phaseModel& alpha2
180  ) const;
181 
182 
183 public:
184 
185  //- Runtime type information
186  TypeName("multiphaseMixtureThermo");
187 
188 
189  // Constructors
190 
191  //- Construct from components
193  (
194  const volVectorField& U,
195  const surfaceScalarField& phi
196  );
197 
198 
199  //- Destructor
200  virtual ~multiphaseMixtureThermo()
201  {}
202 
203 
204  // Member Functions
205 
206  //- Return the phases
207  const PtrDictionary<phaseModel>& phases() const
208  {
209  return phases_;
210  }
211 
212  //- Return non-const access to the phases
214  {
215  return phases_;
216  }
217 
218  //- Return the velocity
219  const volVectorField& U() const
220  {
221  return U_;
222  }
223 
224  //- Return the volumetric flux
225  const surfaceScalarField& phi() const
226  {
227  return phi_;
228  }
230  const surfaceScalarField& rhoPhi() const
231  {
232  return rhoPhi_;
233  }
234 
235  //- Update properties
236  virtual void correct();
237 
238  //- Update densities for given pressure change
239  void correctRho(const volScalarField& dp);
240 
241  //- Return true if the equation of state is incompressible
242  // i.e. rho != f(p)
243  virtual bool incompressible() const;
244 
245  //- Return true if the equation of state is isochoric
246  // i.e. rho = const
247  virtual bool isochoric() const;
248 
249 
250  // Access to thermodynamic state variables
251 
252  //- Enthalpy/Internal energy [J/kg]
253  // Non-const access allowed for transport equations
254  virtual volScalarField& he()
255  {
257  return phases_[0].thermo().he();
258  }
259 
260  //- Enthalpy/Internal energy [J/kg]
261  virtual const volScalarField& he() const
262  {
264  return phases_[0].thermo().he();
265  }
266 
267  //- Enthalpy/Internal energy
268  // for given pressure and temperature [J/kg]
269  virtual tmp<volScalarField> he
270  (
271  const volScalarField& p,
272  const volScalarField& T
273  ) const;
274 
275  //- Enthalpy/Internal energy for cell-set [J/kg]
276  virtual tmp<scalarField> he
277  (
278  const scalarField& p,
279  const scalarField& T,
280  const labelList& cells
281  ) const;
282 
283  //- Enthalpy/Internal energy for patch [J/kg]
284  virtual tmp<scalarField> he
285  (
286  const scalarField& p,
287  const scalarField& T,
288  const label patchi
289  ) const;
290 
291  //- Chemical enthalpy [J/kg]
292  virtual tmp<volScalarField> hc() const;
293 
294  //- Temperature from enthalpy/internal energy for cell-set
295  virtual tmp<scalarField> THE
296  (
297  const scalarField& h,
298  const scalarField& p,
299  const scalarField& T0, // starting temperature
300  const labelList& cells
301  ) const;
302 
303  //- Temperature from enthalpy/internal energy for patch
304  virtual tmp<scalarField> THE
305  (
306  const scalarField& h,
307  const scalarField& p,
308  const scalarField& T0, // starting temperature
309  const label patchi
310  ) const;
311 
312 
313  // Fields derived from thermodynamic state variables
314 
315  //- Density [kg/m^3]
316  virtual tmp<volScalarField> rho() const;
317 
318  //- Density for patch [kg/m^3]
319  virtual tmp<scalarField> rho(const label patchi) const;
320 
321  //- Heat capacity at constant pressure [J/kg/K]
322  virtual tmp<volScalarField> Cp() const;
323 
324  //- Heat capacity at constant pressure for patch [J/kg/K]
325  virtual tmp<scalarField> Cp
326  (
327  const scalarField& p,
328  const scalarField& T,
329  const label patchi
330  ) const;
331 
332  //- Heat capacity at constant volume [J/kg/K]
333  virtual tmp<volScalarField> Cv() const;
334 
335  //- Heat capacity at constant volume for patch [J/kg/K]
336  virtual tmp<scalarField> Cv
337  (
338  const scalarField& p,
339  const scalarField& T,
340  const label patchi
341  ) const;
342 
343  //- Gamma = Cp/Cv []
344  virtual tmp<volScalarField> gamma() const;
345 
346  //- Gamma = Cp/Cv for patch []
347  virtual tmp<scalarField> gamma
348  (
349  const scalarField& p,
350  const scalarField& T,
351  const label patchi
352  ) const;
353 
354  //- Heat capacity at constant pressure/volume [J/kg/K]
355  virtual tmp<volScalarField> Cpv() const;
356 
357  //- Heat capacity at constant pressure/volume for patch [J/kg/K]
358  virtual tmp<scalarField> Cpv
359  (
360  const scalarField& p,
361  const scalarField& T,
362  const label patchi
363  ) const;
364 
365  //- Heat capacity ratio []
366  virtual tmp<volScalarField> CpByCpv() const;
367 
368  //- Heat capacity ratio for patch []
369  virtual tmp<scalarField> CpByCpv
370  (
371  const scalarField& p,
372  const scalarField& T,
373  const label patchi
374  ) const;
375 
376 
377  // Fields derived from transport state variables
378 
379  //- Kinematic viscosity of mixture [m^2/s]
380  virtual tmp<volScalarField> nu() const;
381 
382  //- Kinematic viscosity of mixture for patch [m^2/s]
383  virtual tmp<scalarField> nu(const label patchi) const;
384 
385  //- Thermal diffusivity for temperature of mixture [J/m/s/K]
386  virtual tmp<volScalarField> kappa() const;
387 
388  //- Thermal diffusivity of mixture for patch [J/m/s/K]
389  virtual tmp<scalarField> kappa
390  (
391  const label patchi
392  ) const;
393 
394  //- Effective thermal diffusivity of mixture [J/m/s/K]
396  (
397  const volScalarField& alphat
398  ) const;
399 
400  //- Effective thermal diffusivity of mixture for patch [J/m/s/K]
401  virtual tmp<scalarField> kappaEff
402  (
403  const scalarField& alphat,
404  const label patchi
405  ) const;
406 
407  //- Effective thermal diffusivity of mixture [J/m/s/K]
409  (
410  const volScalarField& alphat
411  ) const;
412 
413  //- Effective thermal diffusivity of mixture for patch [J/m/s/K]
414  virtual tmp<scalarField> alphaEff
415  (
416  const scalarField& alphat,
417  const label patchi
418  ) const;
419 
420 
421  //- Return the phase-averaged reciprocal Cv
422  tmp<volScalarField> rCv() const;
423 
425 
426  //- Indicator of the proximity of the interface
427  // Field values are 1 near and 0 away for the interface.
429 
430  //- Solve for the mixture phase-fractions
431  void solve();
432 };
433 
434 
435 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
436 
437 } // End namespace Foam
438 
439 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
440 
441 #endif
442 
443 // ************************************************************************* //
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:90
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
Chemical enthalpy [J/kg].
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:489
virtual ~multiphaseMixtureThermo()
Destructor.
CGAL::Exact_predicates_exact_constructions_kernel K
virtual tmp< volScalarField > Cp() const
Heat capacity at constant pressure [J/kg/K].
virtual tmp< volScalarField > CpByCpv() const
Heat capacity ratio [].
virtual volScalarField & p()
Pressure [Pa].
Definition: basicThermo.C:477
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 cellShapeList & cells
virtual tmp< volScalarField > rho() const
Density [kg/m^3].
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
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
const word & name() const
Name function is needed to disambiguate those inherited.
alpha2
Definition: alphaEqn.H:112
virtual tmp< volScalarField > nu() const
Kinematic viscosity of mixture [m^2/s].
label operator()(const interfacePair &key) const
volScalarField & alpha1
const Type & second() const
Return second.
Definition: Pair.H:99
TypeName("multiphaseMixtureThermo")
Runtime type information.
virtual tmp< volScalarField > kappaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
multiphaseMixtureThermo(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
label patchi
volScalarField & h
Planck constant.
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:54
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
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.
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 [J/m/s/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:366
virtual tmp< volScalarField > kappa() const
Thermal diffusivity for temperature of mixture [J/m/s/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:87
Basic thermodynamic properties based on compressibility.
Definition: psiThermo.H:49
virtual bool isochoric() const
Return true if the equation of state is isochoric.
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
Namespace for OpenFOAM.
virtual bool incompressible() const
Return true if the equation of state is incompressible.
virtual void correct()
Update properties.