twoPhaseSystem.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-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::twoPhaseSystem
26 
27 Description
28 
29 SourceFiles
30  twoPhaseSystem.C
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #ifndef twoPhaseSystem_H
35 #define twoPhaseSystem_H
36 
37 #include "IOdictionary.H"
38 #include "phaseModel.H"
39 #include "phasePair.H"
40 #include "orderedPhasePair.H"
41 #include "volFields.H"
42 #include "surfaceFields.H"
43 #include "dragModel.H"
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 namespace Foam
48 {
49 
50 class virtualMassModel;
51 class heatTransferModel;
52 class liftModel;
53 class wallLubricationModel;
54 class turbulentDispersionModel;
55 
56 class blendingMethod;
57 template<class modelType> class BlendedInterfacialModel;
58 
59 /*---------------------------------------------------------------------------*\
60  Class twoPhaseSystem Declaration
61 \*---------------------------------------------------------------------------*/
62 
63 class twoPhaseSystem
64 :
65  public IOdictionary
66 {
67  // Private Data
68 
69  //- Reference to the mesh
70  const fvMesh& mesh_;
71 
72  //- Phase model 1
73  phaseModel phase1_;
74 
75  //- Phase model 2
76  phaseModel phase2_;
77 
78  //- Total volumetric flux
79  surfaceScalarField phi_;
80 
81  //- Dilatation term
82  volScalarField dgdt_;
83 
84  //- Optional dispersion diffusivity
85  tmp<surfaceScalarField> pPrimeByA_;
86 
87  //- Unordered phase pair
88  autoPtr<phasePair> pair_;
89 
90  //- Phase pair for phase 1 dispersed in phase 2
91  autoPtr<orderedPhasePair> pair1In2_;
92 
93  //- Phase pair for phase 2 dispersed in phase 1
94  autoPtr<orderedPhasePair> pair2In1_;
95 
96  //- Blending methods
97  HashTable<autoPtr<blendingMethod>, word, word::hash> blendingMethods_;
98 
99  //- Drag model
100  autoPtr<BlendedInterfacialModel<dragModel>> drag_;
101 
102  //- Virtual mass model
103  autoPtr<BlendedInterfacialModel<virtualMassModel>> virtualMass_;
104 
105  //- Heat transfer model
106  autoPtr<BlendedInterfacialModel<heatTransferModel>> heatTransfer_;
107 
108  //- Lift model
109  autoPtr<BlendedInterfacialModel<liftModel>> lift_;
110 
111  //- Wall lubrication model
112  autoPtr<BlendedInterfacialModel<wallLubricationModel>>
113  wallLubrication_;
114 
115  //- Wall lubrication model
116  autoPtr<BlendedInterfacialModel<turbulentDispersionModel>>
117  turbulentDispersion_;
118 
119 
120  // Private Member Functions
121 
122  //- Return the mixture flux
123  tmp<surfaceScalarField> calcPhi() const;
124 
125 
126 public:
127 
128  // Constructors
129 
130  //- Construct from fvMesh
131  twoPhaseSystem(const fvMesh&, const dimensionedVector& g);
132 
133 
134  //- Destructor
135  virtual ~twoPhaseSystem();
136 
137 
138  // Member Functions
139 
140  //- Return the mixture density
141  tmp<volScalarField> rho() const;
142 
143  //- Return the mixture velocity
144  tmp<volVectorField> U() const;
145 
146  //- Return the drag coefficient
147  tmp<volScalarField> Kd() const;
148 
149  //- Return the face drag coefficient
150  tmp<surfaceScalarField> Kdf() const;
151 
152  //- Return the virtual mass coefficient
153  tmp<volScalarField> Vm() const;
154 
155  //- Return the face virtual mass coefficient
156  tmp<surfaceScalarField> Vmf() const;
157 
158  //- Return the heat transfer coefficient
159  tmp<volScalarField> Kh() const;
160 
161  //- Return the combined force (lift + wall-lubrication)
162  tmp<volVectorField> F() const;
163 
164  //- Return the combined face-force (lift + wall-lubrication)
165  tmp<surfaceScalarField> Ff() const;
166 
167  //- Return the turbulent diffusivity
168  // Multiplies the phase-fraction gradient
169  tmp<volScalarField> D() const;
170 
171  //- Solve for the two-phase-fractions
172  void solve();
173 
174  //- Correct two-phase properties other than turbulence
175  void correct();
176 
177  //- Correct two-phase turbulence
178  void correctTurbulence();
179 
180  //- Read base phaseProperties dictionary
181  bool read();
182 
183  // Access
184 
185  //- Access a sub model between a phase pair
186  template<class modelType>
187  const modelType& lookupSubModel(const phasePair& key) const;
188 
189  //- Access a sub model between two phases
190  template<class modelType>
191  const modelType& lookupSubModel
192  (
193  const phaseModel& dispersed,
194  const phaseModel& continuous
195  ) const;
196 
197  //- Return the surface tension coefficient
198  const dimensionedScalar& sigma() const;
199 
200  //- Return the mesh
201  inline const fvMesh& mesh() const;
202 
203  //- Return phase model 1
204  inline const phaseModel& phase1() const;
205 
206  //- Return non-const access to phase model 1
207  inline phaseModel& phase1();
208 
209  //- Return phase model 2
210  inline const phaseModel& phase2() const;
211 
212  //- Return non-const access to phase model 2
213  inline phaseModel& phase2();
214 
215  //- Return the phase not given as an argument
216  inline const phaseModel& otherPhase(const phaseModel& phase) const;
217 
218  //- Return the mixture flux
219  inline const surfaceScalarField& phi() const;
220 
221  //- Return non-const access to the mixture flux
222  inline surfaceScalarField& phi();
223 
224  //- Return the dilatation term
225  inline const volScalarField& dgdt() const;
226 
227  //- Return non-const access to the dilatation parameter
228  inline volScalarField& dgdt();
229 
230  //- Return non-const access to the dispersion diffusivity
231  inline tmp<surfaceScalarField>& pPrimeByA();
232 };
233 
234 
235 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
236 
237 } // End namespace Foam
238 
239 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240 
241 #include "twoPhaseSystemI.H"
242 
243 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244 
245 #ifdef NoRepository
246  #include "twoPhaseSystemTemplates.C"
247 #endif
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 #endif
252 
253 // ************************************************************************* //
Foam::surfaceFields.
phaseModel & phase1_
Phase model 1.
tmp< volVectorField > F() const
Return the combined force (lift + wall-lubrication)
const phaseModel & phase1() const
Return phase model 1.
const modelType & lookupSubModel(const phasePair &key) const
Access a sub model between a phase pair.
tmp< volScalarField > Kh() const
Return the heat transfer coefficient.
void correct()
Correct two-phase properties other than turbulence.
tmp< surfaceScalarField > & pPrimeByA()
Return non-const access to the dispersion diffusivity.
tmp< volScalarField > sigma() const
Return the surface tension coefficient.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
const fvMesh & mesh() const
Return the mesh.
tmp< volScalarField > Kd() const
Return the drag coefficient.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const phaseModel & phase2() const
Return phase model 2.
phaseModel & phase2_
Phase model 2.
const volScalarField & dgdt() const
Return the dilatation term.
twoPhaseSystem(const fvMesh &)
Construct from fvMesh.
tmp< surfaceScalarField > Kdf() const
Return the face drag coefficient.
tmp< volScalarField > D() const
Return the turbulent diffusivity.
bool read()
Read base phaseProperties dictionary.
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument.
tmp< volVectorField > U() const
Return the mixture velocity.
void correctTurbulence()
Correct two-phase turbulence.
virtual void solve()
Solve for the phase fractions.
virtual ~twoPhaseSystem()
Destructor.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< surfaceScalarField > Vmf() const
Return the face virtual mass coefficient.
tmp< surfaceScalarField > Ff() const
Return the combined face-force (lift + wall-lubrication)
const surfaceScalarField & phi() const
Return the mixture flux.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const dimensionedVector & g
tmp< volScalarField > rho() const
Return the mixture density.
Namespace for OpenFOAM.
tmp< volScalarField > Vm() const
Return the virtual mass coefficient.