multiphaseMixture.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) 2011-2019 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::multiphaseMixture
26 
27 Description
28  Incompressible multi-phase mixture with built in solution for the
29  phase fractions with interface compression for interface-capturing.
30 
31  Derived from transportModel so that it can be used in conjunction with
32  the incompressible turbulence models.
33 
34  Surface tension and contact-angle is handled for the interface
35  between each phase-pair.
36 
37 SourceFiles
38  multiphaseMixture.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef multiphaseMixture_H
43 #define multiphaseMixture_H
44 
46 #include "IOdictionary.H"
47 #include "phase.H"
48 #include "PtrDictionary.H"
49 #include "volFields.H"
50 #include "surfaceFields.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 /*---------------------------------------------------------------------------*\
58  Class multiphaseMixture Declaration
59 \*---------------------------------------------------------------------------*/
60 
62 :
63  public IOdictionary,
64  public transportModel
65 {
66 public:
67 
68  class interfacePair
69  :
70  public Pair<word>
71  {
72  public:
73 
74  class hash
75  :
76  public Hash<interfacePair>
77  {
78  public:
79 
80  hash()
81  {}
82 
83  label operator()(const interfacePair& key) const
84  {
85  return word::hash()(key.first()) + word::hash()(key.second());
86  }
87  };
88 
89 
90  // Constructors
91 
93  {}
94 
95  interfacePair(const word& alpha1Name, const word& alpha2Name)
96  :
97  Pair<word>(alpha1Name, alpha2Name)
98  {}
99 
100  interfacePair(const phase& alpha1, const phase& alpha2)
101  :
102  Pair<word>(alpha1.name(), alpha2.name())
103  {}
104 
105 
106  // Friend Operators
107 
108  friend bool operator==
109  (
110  const interfacePair& a,
111  const interfacePair& b
112  )
113  {
114  return
115  (
116  ((a.first() == b.first()) && (a.second() == b.second()))
117  || ((a.first() == b.second()) && (a.second() == b.first()))
118  );
119  }
120 
121  friend bool operator!=
122  (
123  const interfacePair& a,
124  const interfacePair& b
125  )
126  {
127  return (!(a == b));
128  }
129  };
130 
131 
132 private:
133 
134  // Private Data
135 
136  //- Dictionary of phases
137  PtrDictionary<phase> phases_;
138 
139  const fvMesh& mesh_;
140  const volVectorField& U_;
141  const surfaceScalarField& phi_;
142 
143  surfaceScalarField rhoPhi_;
144  volScalarField alphas_;
145 
146  volScalarField nu_;
147 
149  sigmaTable;
150 
151  sigmaTable sigmas_;
152  dimensionSet dimSigma_;
153 
154  //- Stabilisation for normalisation of the interface normal
155  const dimensionedScalar deltaN_;
156 
157  //- Conversion factor for degrees into radians
158  static const scalar convertToRad;
159 
160 
161  // Private Member Functions
162 
163  void calcAlphas();
164 
165  void solveAlphas(const scalar cAlpha);
166 
168  (
169  const volScalarField& alpha1,
170  const volScalarField& alpha2
171  ) const;
172 
174  (
175  const volScalarField& alpha1,
176  const volScalarField& alpha2
177  ) const;
178 
179  void correctContactAngle
180  (
181  const phase& alpha1,
182  const phase& alpha2,
183  surfaceVectorField::Boundary& nHatb
184  ) const;
185 
186  tmp<volScalarField> K(const phase& alpha1, const phase& alpha2) const;
187 
188 
189 public:
190 
191  // Constructors
192 
193  //- Construct from components
195  (
196  const volVectorField& U,
197  const surfaceScalarField& phi
198  );
199 
200 
201  //- Destructor
202  virtual ~multiphaseMixture()
203  {}
204 
205 
206  // Member Functions
207 
208  //- Return the phases
209  const PtrDictionary<phase>& phases() const
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  //- Return the mixture density
232  tmp<volScalarField> rho() const;
233 
234  //- Return the mixture density for patch
235  tmp<scalarField> rho(const label patchi) const;
236 
237  //- Return the dynamic laminar viscosity
238  tmp<volScalarField> mu() const;
239 
240  //- Return the dynamic laminar viscosity for patch
241  tmp<scalarField> mu(const label patchi) const;
242 
243  //- Return the face-interpolated dynamic laminar viscosity
245 
246  //- Return the kinematic laminar viscosity
247  tmp<volScalarField> nu() const;
248 
249  //- Return the laminar viscosity for patch
250  tmp<scalarField> nu(const label patchi) const;
251 
252  //- Return the face-interpolated dynamic laminar viscosity
254 
256 
257  //- Indicator of the proximity of the interface
258  // Field values are 1 near and 0 away for the interface.
260 
261  //- Solve for the mixture phase-fractions
262  void solve();
263 
264  //- Correct the mixture properties
265  void correct();
266 
267  //- Read base transportProperties dictionary
268  bool read();
269 };
270 
271 
272 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
273 
274 } // End namespace Foam
275 
276 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
277 
278 #endif
279 
280 // ************************************************************************* //
Template dictionary class which manages the storage associated with it.
Definition: PtrDictionary.H:53
Foam::surfaceFields.
Hashing function class, shared by all the derived classes.
Definition: string.H:90
tmp< surfaceScalarField > muf() const
Return the face-interpolated dynamic laminar viscosity.
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.
Incompressible multi-phase mixture with built in solution for the phase fractions with interface comp...
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
const PtrDictionary< phase > & phases() const
Return the phases.
virtual ~multiphaseMixture()
Destructor.
CGAL::Exact_predicates_exact_constructions_kernel K
multiphaseMixture(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
tmp< volScalarField > rho() const
Return the mixture density.
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
bool read()
Read base transportProperties dictionary.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
A class for handling words, derived from string.
Definition: word.H:59
void correct()
Correct the mixture properties.
const volScalarField & alpha1
const word & name() const
Name function is needed to disambiguate those inherited.
alpha2
Definition: alphaEqn.H:115
tmp< volScalarField > nu() const
Return the kinematic laminar viscosity.
tmp< volScalarField > mu() const
Return the dynamic laminar viscosity.
void solve()
Solve for the mixture phase-fractions.
const Type & second() const
Return second.
Definition: Pair.H:99
tmp< surfaceScalarField > nuf() const
Return the face-interpolated dynamic laminar viscosity.
label patchi
Base-class for all transport models used by the incompressible turbulence models. ...
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
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:52
const volVectorField & U() const
Return the velocity.
label operator()(const interfacePair &key) const
A class for managing temporary objects.
Definition: PtrList.H:53
const Type & first() const
Return first.
Definition: Pair.H:87
const surfaceScalarField & rhoPhi() const
Namespace for OpenFOAM.
tmp< surfaceScalarField > surfaceTensionForce() const