multiphaseMixture.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) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 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::multiphaseMixture
29
30Description
31 Incompressible multi-phase mixture with built in solution for the
32 phase fractions with interface compression for interface-capturing.
33
34 Derived from transportModel so that it can be unused in conjunction with
35 the incompressible turbulence models.
36
37 Surface tension and contact-angle is handled for the interface
38 between each phase-pair.
39
40SourceFiles
41 multiphaseMixture.C
42
43\*---------------------------------------------------------------------------*/
44
45#ifndef multiphaseMixture_H
46#define multiphaseMixture_H
47
49#include "IOdictionary.H"
50#include "phase.H"
51#include "PtrDictionary.H"
52#include "volFields.H"
53#include "surfaceFields.H"
54
55// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
56
57namespace Foam
58{
59
60// Class forward declarations
61class alphaContactAngleFvPatchScalarField;
62
63/*---------------------------------------------------------------------------*\
64 Class multiphaseMixture Declaration
65\*---------------------------------------------------------------------------*/
68:
69 public IOdictionary,
70 public transportModel
71{
72public:
73
74 //- Symmetric pair of interface names
75 class interfacePair
76 :
77 public Pair<word>
78 {
79 public:
80
81 // Always use symmetric hashing
83
84 // Always use symmetric hashing (alias)
86
87
88 // Constructors
90 interfacePair() = default;
92 interfacePair(const word& alpha1Name, const word& alpha2Name)
93 :
94 Pair<word>(alpha1Name, alpha2Name)
95 {}
97 interfacePair(const phase& alpha1, const phase& alpha2)
98 :
100 {}
101
102
103 // Friend Operators
105 friend bool operator==
106 (
107 const interfacePair& a,
108 const interfacePair& b
109 )
110 {
111 return (0 != Pair<word>::compare(a, b));
112 }
114 friend bool operator!=
115 (
116 const interfacePair& a,
117 const interfacePair& b
118 )
119 {
120 return (!(a == b));
121 }
122 };
123
124
125private:
126
127 // Private data
128
129 //- Dictionary of phases
130 PtrDictionary<phase> phases_;
131
132 const fvMesh& mesh_;
133 const volVectorField& U_;
134 const surfaceScalarField& phi_;
135
136 surfaceScalarField rhoPhi_;
137 volScalarField alphas_;
138
139 volScalarField nu_;
140
142 sigmaTable;
143
144 sigmaTable sigmas_;
145 dimensionSet dimSigma_;
146
147 //- Stabilisation for normalisation of the interface normal
148 const dimensionedScalar deltaN_;
149
150
151 // Private member functions
152
153 void calcAlphas();
154
155 void solveAlphas(const scalar cAlpha);
156
158 (
159 const volScalarField& alpha1,
161 ) const;
162
164 (
165 const volScalarField& alpha1,
167 ) const;
168
169 void correctContactAngle
170 (
171 const phase& alpha1,
172 const phase& alpha2,
174 ) const;
175
176 void correctBoundaryContactAngle
177 (
179 label patchi,
180 const phase& alpha1,
181 const phase& alpha2,
183 ) const;
184
185 tmp<volScalarField> K(const phase& alpha1, const phase& alpha2) const;
186
187
188public:
189
190 // Constructors
191
192 //- Construct from components
194 (
195 const volVectorField& U,
197 );
198
199
200 //- Destructor
201 virtual ~multiphaseMixture() = default;
202
203
204 // Member Functions
205
206 //- Return the phases
207 const PtrDictionary<phase>& phases() const
208 {
209 return phases_;
210 }
211
212 //- Return the velocity
213 const volVectorField& U() const
214 {
215 return U_;
216 }
217
218 //- Return the volumetric flux
219 const surfaceScalarField& phi() const
220 {
221 return phi_;
222 }
224 const surfaceScalarField& rhoPhi() const
225 {
226 return rhoPhi_;
227 }
228
229 //- Return the mixture density
230 tmp<volScalarField> rho() const;
231
232 //- Return the mixture density for patch
233 tmp<scalarField> rho(const label patchi) const;
234
235 //- Return the dynamic laminar viscosity
236 tmp<volScalarField> mu() const;
237
238 //- Return the dynamic laminar viscosity for patch
239 tmp<scalarField> mu(const label patchi) const;
240
241 //- Return the face-interpolated dynamic laminar viscosity
243
244 //- Return the kinematic laminar viscosity
245 tmp<volScalarField> nu() const;
246
247 //- Return the laminar viscosity for patch
248 tmp<scalarField> nu(const label patchi) const;
249
250 //- Return the face-interpolated dynamic laminar viscosity
254
255 //- Indicator of the proximity of the interface
256 // Field values are 1 near and 0 away for the interface.
258
259 //- Solve for the mixture phase-fractions
260 void solve();
261
262 //- Correct the mixture properties
263 void correct();
264
265 //- Read base transportProperties dictionary
266 bool read();
267};
268
269
270// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
271
272} // End namespace Foam
273
274// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
275
276#endif
277
278// ************************************************************************* //
CGAL::Exact_predicates_exact_constructions_kernel K
const volScalarField & alpha1
const volScalarField & alpha2
A HashTable similar to std::unordered_map.
Definition: HashTable.H:123
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
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
Contact-angle boundary condition for multi-phase interface-capturing simulations. Used in conjunction...
const word & name() const
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
Symmetric pair of interface names.
interfacePair(const word &alpha1Name, const word &alpha2Name)
interfacePair(const phase &alpha1, const phase &alpha2)
Incompressible multi-phase mixture with built in solution for the phase fractions with interface comp...
const volVectorField & U() const
Return the velocity.
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
tmp< surfaceScalarField > surfaceTensionForce() const
tmp< scalarField > mu(const label patchi) const
Return the dynamic laminar viscosity for patch.
void correct()
Correct the mixture properties.
const surfaceScalarField & phi() const
Return the volumetric flux.
tmp< surfaceScalarField > muf() const
Return the face-interpolated dynamic laminar viscosity.
tmp< scalarField > nu(const label patchi) const
Return the laminar viscosity for patch.
virtual ~multiphaseMixture()=default
Destructor.
const surfaceScalarField & rhoPhi() const
tmp< volScalarField > nu() const
Return the kinematic laminar viscosity.
tmp< scalarField > rho(const label patchi) const
Return the mixture density for patch.
tmp< volScalarField > mu() const
Return the dynamic laminar viscosity.
tmp< volScalarField > rho() const
Return the mixture density.
const PtrDictionary< phase > & phases() const
Return the phases.
multiphaseMixture(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
void solve()
Solve for the mixture phase-fractions.
bool read()
Read base transportProperties dictionary.
tmp< surfaceScalarField > nuf() const
Return the face-interpolated dynamic laminar viscosity.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:57
A class for managing temporary objects.
Definition: tmp.H:65
Base-class for all transport models used by the incompressible turbulence models.
A class for handling words, derived from Foam::string.
Definition: word.H:68
Namespace for OpenFOAM.
volScalarField & b
Definition: createFields.H:27
Symmetric hashing functor for Pair, hashes lower value first.
Definition: Pair.H:152
Foam::surfaceFields.