All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
BlendedInterfacialModel.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) 2014-2020 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::BlendedInterfacialModel
26 
27 Description
28 
29 SourceFiles
30  BlendedInterfacialModel.C
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #ifndef BlendedInterfacialModel_H
35 #define BlendedInterfacialModel_H
36 
37 #include "blendingMethod.H"
38 #include "phasePair.H"
39 #include "orderedPhasePair.H"
40 #include "HashPtrTable.H"
41 #include "hashedWordList.H"
42 #include "geometricZeroField.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 
49 /*---------------------------------------------------------------------------*\
50  Class BlendedInterfacialModel Declaration
51 \*---------------------------------------------------------------------------*/
52 
53 template<class ModelType>
55 :
56  public regIOobject
57 {
58  // Private Data
59 
60  //- Reference to phase 1
61  const phaseModel& phase1_;
62 
63  //- Reference to phase 2
64  const phaseModel& phase2_;
65 
66  //- Blending model
67  const blendingMethod& blending_;
68 
69  //- Model for region with no obvious dispersed phase
70  autoPtr<ModelType> model_;
71 
72  //- Model for dispersed phase 1 in continuous phase 2
73  autoPtr<ModelType> model1In2_;
74 
75  //- Model for dispersed phase 2 in continuous phase 1
76  autoPtr<ModelType> model2In1_;
77 
78  //- If true set coefficients and forces to 0 at fixed-flux BCs
79  bool correctFixedFluxBCs_;
80 
81 
82  // Private Member Functions
83 
84  //- Calculate the blending coefficients
85  template<template<class> class PatchField, class GeoMesh>
86  void calculateBlendingCoeffs
87  (
90  const bool subtract
91  ) const;
92 
93  //- Correct coeff/value on fixed flux boundary conditions
94  template<class Type, template<class> class PatchField, class GeoMesh>
95  void correctFixedFluxBCs
96  (
98  ) const;
99 
100  //- Return the blended coeff/value
101  template
102  <
103  class Type,
104  template<class> class PatchField,
105  class GeoMesh,
106  class ... Args
107  >
109  (
111  (ModelType::*method)(Args ...) const,
112  const word& name,
113  const dimensionSet& dims,
114  const bool subtract,
115  Args ... args
116  ) const;
117 
118  //- Return the blended coeff/value
119  template
120  <
121  class Type,
122  template<class> class PatchField,
123  class GeoMesh,
124  class ... Args
125  >
127  (
129  (ModelType::*method)(Args ...) const,
130  const word& name,
131  const dimensionSet& dims,
132  const bool subtract,
133  Args ... args
134  ) const;
135 
136 
137 public:
138 
139  //- Runtime type information
140  TypeName("BlendedInterfacialModel");
141 
142 
143  // Constructors
144 
145  //- Construct from two phases, blending method and three models
147  (
148  const phaseModel& phase1,
149  const phaseModel& phase2,
150  const blendingMethod& blending,
151  autoPtr<ModelType> model,
152  autoPtr<ModelType> model1In2,
153  autoPtr<ModelType> model2In1,
154  const bool correctFixedFluxBCs = true
155  );
156 
157 
158  //- Construct from the model table, dictionary and pairs
160  (
161  const phasePair::dictTable& modelTable,
162  const blendingMethod& blending,
163  const phasePair& pair,
164  const orderedPhasePair& pair1In2,
165  const orderedPhasePair& pair2In1,
166  const bool correctFixedFluxBCs = true
167  );
168 
169  //- Disallow default bitwise copy construction
171  (
173  ) = delete;
174 
175 
176  //- Destructor
178 
179 
180  // Member Functions
181 
182  //- Return the blended force coefficient
183  tmp<volScalarField> K() const;
184 
185  //- Return the blended force coefficient with a specified residual alpha
186  tmp<volScalarField> K(const scalar residualAlpha) const;
187 
188  //- Return the face blended force coefficient
189  tmp<surfaceScalarField> Kf() const;
190 
191  //- Return the blended force
192  template<class Type>
194 
195  //- Return the face blended force
196  tmp<surfaceScalarField> Ff() const;
197 
198  //- Return the blended diffusivity
199  tmp<volScalarField> D() const;
200 
201  //- Return the list of individual species that are transferred
202  bool mixture() const;
203 
204  //- Return the blended mass transfer rate
205  tmp<volScalarField> dmdtf() const;
206 
207  //- Return the list of individual species that are transferred
208  hashedWordList species() const;
209 
210  //- Return the blended mass transfer rates for individual species
212 
213  //- Dummy write for regIOobject
214  bool writeData(Ostream& os) const;
215 
216 
217  // Member Operators
218 
219  //- Disallow default bitwise assignment
220  void operator=(const BlendedInterfacialModel<ModelType>&) = delete;
221 };
222 
223 
224 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
226 #define defineBlendedInterfacialModelTypeNameAndDebug(ModelType, DebugSwitch) \
227  \
228  defineTemplateTypeNameAndDebugWithName \
229  ( \
230  BlendedInterfacialModel<ModelType>, \
231  ( \
232  word(BlendedInterfacialModel<ModelType>::typeName_()) + "<" \
233  + ModelType::typeName_() + ">" \
234  ).c_str(), \
235  DebugSwitch \
236  );
237 
238 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
239 
240 } // End namespace Foam
241 
242 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
243 
244 #ifdef NoRepository
245  #include "BlendedInterfacialModel.C"
246 #endif
247 
248 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
249 
250 #endif
251 
252 // ************************************************************************* //
tmp< volScalarField > dmdtf() const
Return the blended mass transfer rate.
const word & name() const
Return name.
Definition: IOobject.H:303
bool writeData(Ostream &os) const
Dummy write for regIOobject.
~BlendedInterfacialModel()
Destructor.
bool mixture() const
Return the list of individual species that are transferred.
tmp< volScalarField > K() const
Return the blended force coefficient.
tmp< surfaceScalarField > Ff() const
Return the face blended force.
tmp< surfaceScalarField > Kf() const
Return the face blended force coefficient.
Generic GeometricField class.
A HashTable specialisation for hashing pointers.
Definition: HashPtrTable.H:50
tmp< GeometricField< Type, fvPatchField, volMesh > > F() const
Return the blended force.
void subtract(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
scalar f1
Definition: createFields.H:15
Dimension set for the base types.
Definition: dimensionSet.H:120
A class for handling words, derived from string.
Definition: word.H:59
BlendedInterfacialModel(const phaseModel &phase1, const phaseModel &phase2, const blendingMethod &blending, autoPtr< ModelType > model, autoPtr< ModelType > model1In2, autoPtr< ModelType > model2In1, const bool correctFixedFluxBCs=true)
Construct from two phases, blending method and three models.
An STL-conforming hash table.
Definition: HashTable.H:61
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
hashedWordList species() const
Return the list of individual species that are transferred.
A wordList with hashed indices for faster lookup by name.
tmp< volScalarField > D() const
Return the blended diffusivity.
HashPtrTable< volScalarField > dmidtf() const
Return the blended mass transfer rates for individual species.
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:52
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
rDeltaTY field()
A class for managing temporary objects.
Definition: PtrList.H:53
Foam::argList args(argc, argv)
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
void operator=(const BlendedInterfacialModel< ModelType > &)=delete
Disallow default bitwise assignment.
TypeName("BlendedInterfacialModel")
Runtime type information.
Namespace for OpenFOAM.