twoPhaseSystem.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-2015 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  Class which solves the volume fraction equations for two phases.
29 
30 SourceFiles
31  twoPhaseSystem.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef twoPhaseSystem_H
36 #define twoPhaseSystem_H
37 
38 #include "phaseSystem.H"
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 namespace Foam
43 {
44 
45 class dragModel;
46 class virtualMassModel;
47 
48 /*---------------------------------------------------------------------------*\
49  Class twoPhaseSystem Declaration
50 \*---------------------------------------------------------------------------*/
51 
52 class twoPhaseSystem
53 :
54  public phaseSystem
55 {
56  // Private member functions
57 
58  //- Return the drag coefficient for phase pair
59  virtual tmp<volScalarField> Kd(const phasePairKey& key) const = 0;
60 
61  //- Return the face drag coefficient for phase pair
62  virtual tmp<surfaceScalarField> Kdf(const phasePairKey& key) const = 0;
63 
64  //- Return the virtual mass coefficient for phase pair
65  virtual tmp<volScalarField> Vm(const phasePairKey& key) const = 0;
66 
67  //- Return the face virtual mass coefficient for phase pair
68  virtual tmp<surfaceScalarField> Vmf(const phasePairKey& key) const = 0;
69 
70  //- Return the combined force (lift + wall-lubrication) for phase pair
71  virtual tmp<volVectorField> F(const phasePairKey& key) const = 0;
72 
73  //- Return the combined face-force (lift + wall-lubrication)
74  // for phase pair
75  virtual tmp<surfaceScalarField> Ff(const phasePairKey& key) const = 0;
76 
77  //- Return the turbulent diffusivity for phase pair
78  // Multiplies the phase-fraction gradient
79  virtual tmp<volScalarField> D(const phasePairKey& key) const = 0;
80 
81  //- Return true if there is mass transfer for phase
82  virtual bool transfersMass(const phaseModel& phase) const = 0;
83 
84  //- Return the interfacial mass flow rate for phase pair
85  virtual tmp<volScalarField> dmdt(const phasePairKey& key) const = 0;
86 
87 
88 protected:
89 
90  // Protected data
91 
92  //- Phase model 1
94 
95  //- Phase model 2
97 
98 
99 public:
100 
101  //- Runtime type information
102  TypeName("twoPhaseSystem");
103 
104  // Declare runtime construction
105 
107  (
108  autoPtr,
110  dictionary,
111  (
112  const fvMesh& mesh
113  ),
114  (mesh)
115  );
116 
117 
118  // Constructors
119 
120  //- Construct from fvMesh
121  twoPhaseSystem(const fvMesh&);
122 
123 
124  //- Destructor
125  virtual ~twoPhaseSystem();
126 
127 
128  // Selectors
129 
131  (
132  const fvMesh& mesh
133  );
134 
135 
136  // Member Functions
137 
138  //- Constant access phase model 1
139  const phaseModel& phase1() const;
140 
141  //- Access phase model 1
142  phaseModel& phase1();
143 
144  //- Constant access phase model 2
145  const phaseModel& phase2() const;
146 
147  //- Access phase model 2
148  phaseModel& phase2();
149 
150  //- Constant access the phase not given as an argument
151  const phaseModel& otherPhase
152  (
153  const phaseModel& phase
154  ) const;
155 
156  //- Return the momentum transfer matrices
158 
159  //- Return the heat transfer matrices
160  virtual autoPtr<heatTransferTable> heatTransfer() const = 0;
161 
162  //- Return the mass transfer matrices
163  virtual autoPtr<massTransferTable> massTransfer() const = 0;
164 
165  using phaseSystem::sigma;
166 
167  //- Return the surface tension coefficient
168  tmp<volScalarField> sigma() const;
169 
170  //- Return the drag model for the given phase
171  const dragModel& drag(const phaseModel& phase) const;
172 
173  //- Return the drag coefficient
174  tmp<volScalarField> Kd() const;
175 
176  //- Return the face drag coefficient
178 
179  //- Return the virtual mass model for the given phase
180  const virtualMassModel& virtualMass(const phaseModel& phase) const;
181 
182  //- Return the virtual mass coefficient
183  tmp<volScalarField> Vm() const;
184 
185  //- Return the face virtual mass coefficient
187 
188  //- Return the combined force (lift + wall-lubrication)
189  tmp<volVectorField> F() const;
190 
191  //- Return the combined face-force (lift + wall-lubrication)
192  tmp<surfaceScalarField> Ff() const;
193 
194  //- Return the turbulent diffusivity
195  // Multiplies the phase-fraction gradient
196  tmp<volScalarField> D() const;
197 
198  //- Return true if there is mass transfer
199  bool transfersMass() const;
200 
201  //- Return the interfacial mass flow rate
202  tmp<volScalarField> dmdt() const;
203 
204  //- Solve for the phase fractions
205  virtual void solve();
206 };
207 
208 
209 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
210 
211 } // End namespace Foam
212 
213 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
214 
215 #include "twoPhaseSystemI.H"
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 #endif
220 
221 // ************************************************************************* //
phaseModel & phase1_
Phase model 1.
const fvMesh & mesh() const
Return the mesh.
tmp< volScalarField > Vm() const
Return the virtual mass coefficient.
const virtualMassModel & virtualMass(const phaseModel &phase) const
Return the virtual mass model for the given phase.
tmp< surfaceScalarField > Ff() const
Return the combined face-force (lift + wall-lubrication)
tmp< volVectorField > F() const
Return the combined force (lift + wall-lubrication)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
tmp< volScalarField > D() const
Return the turbulent diffusivity.
tmp< surfaceScalarField > Vmf() const
Return the face virtual mass coefficient.
declareRunTimeSelectionTable(autoPtr, twoPhaseSystem, dictionary,(const fvMesh &mesh),(mesh))
tmp< surfaceScalarField > Kdf() const
Return the face drag coefficient.
virtual autoPtr< heatTransferTable > heatTransfer() const =0
Return the heat transfer matrices.
Class which solves the volume fraction equations for two phases.
tmp< volScalarField > sigma() const
Return the surface tension coefficient.
TypeName("twoPhaseSystem")
Runtime type information.
virtual autoPtr< massTransferTable > massTransfer() const =0
Return the mass transfer matrices.
phaseModel & phase2_
Phase model 2.
bool transfersMass() const
Return true if there is mass transfer.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:67
twoPhaseSystem(const fvMesh &)
Construct from fvMesh.
const phaseModel & otherPhase(const phaseModel &phase) const
Constant access the phase not given as an argument.
virtual void solve()
Solve for the phase fractions.
tmp< volScalarField > Kd() const
Return the drag coefficient.
const phaseModel & phase2() const
Constant access phase model 2.
virtual ~twoPhaseSystem()
Destructor.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dragModel & drag(const phaseModel &phase) const
Return the drag model for the given phase.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
static autoPtr< twoPhaseSystem > New(const fvMesh &mesh)
tmp< volScalarField > dmdt() const
Return the interfacial mass flow rate.
const phaseModel & phase1() const
Constant access phase model 1.
A class for managing temporary objects.
Definition: PtrList.H:54
tmp< volScalarField > sigma(const phasePairKey &key) const
Return the surface tension coefficient.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
virtual autoPtr< momentumTransferTable > momentumTransfer() const =0
Return the momentum transfer matrices.
Namespace for OpenFOAM.