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-2018 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 
41 #include "geometricZeroField.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 /*---------------------------------------------------------------------------*\
49  Class BlendedInterfacialModel Declaration
50 \*---------------------------------------------------------------------------*/
51 
52 template<class modelType>
53 class BlendedInterfacialModel
54 {
55  // Private data
56 
57  //- Unordered phase pair
58  const phasePair& pair_;
59 
60  //- Ordered phase pair for dispersed phase 1 in continuous phase 2
61  const orderedPhasePair& pair1In2_;
62 
63  //- Ordered phase pair for dispersed phase 2 in continuous phase 1
64  const orderedPhasePair& pair2In1_;
65 
66  //- Model for region with no obvious dispersed phase
67  autoPtr<modelType> model_;
68 
69  //- Model for dispersed phase 1 in continuous phase 2
70  autoPtr<modelType> model1In2_;
71 
72  //- Model for dispersed phase 2 in continuous phase 1
73  autoPtr<modelType> model2In1_;
74 
75  //- Blending model
76  const blendingMethod& blending_;
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  //- Disallow default bitwise copy construct
85  BlendedInterfacialModel(const BlendedInterfacialModel<modelType>&);
86 
87  //- Disallow default bitwise assignment
88  void operator=(const BlendedInterfacialModel<modelType>&);
89 
90  //- Correct coeff/value on fixed flux boundary conditions
91  template<class GeometricField>
92  void correctFixedFluxBCs(GeometricField& field) const;
93 
94 
95 public:
96 
97  // Constructors
98 
99  //- Construct from the model table, dictionary and pairs
100  BlendedInterfacialModel
101  (
102  const phasePair::dictTable& modelTable,
103  const blendingMethod& blending,
104  const phasePair& pair,
105  const orderedPhasePair& pair1In2,
106  const orderedPhasePair& pair2In1,
107  const bool correctFixedFluxBCs = true
108  );
109 
110 
111  //- Destructor
113 
114 
115  // Member Functions
116 
117  //- Return true if a model is specified for the supplied phase
118  bool hasModel(const phaseModel& phase) const;
119 
120  //- Return the model for the supplied phase
121  const modelType& phaseModel(const phaseModel& phase) const;
122 
123  //- Return the blended force coefficient
124  tmp<volScalarField> K() const;
125 
126  //- Return the face blended force coefficient
127  tmp<surfaceScalarField> Kf() const;
128 
129  //- Return the blended force
130  template<class Type>
131  tmp<GeometricField<Type, fvPatchField, volMesh>> F() const;
132 
133  //- Return the face blended force
134  tmp<surfaceScalarField> Ff() const;
135 
136  //- Return the blended diffusivity
137  tmp<volScalarField> D() const;
138 };
139 
140 
141 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
142 
143 } // End namespace Foam
144 
145 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
146 
147 #ifdef NoRepository
148  #include "BlendedInterfacialModel.C"
149 #endif
150 
151 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
152 
153 #endif
154 
155 // ************************************************************************* //
HashTable< dictionary, phasePairKey, phasePairKey::hash > dictTable
Dictionary hash table.
Definition: phasePair.H:59
bool hasModel(const phaseModel &phase) const
Return true if a model is specified for the supplied phase.
const modelType & phaseModel(const phaseModel &phase) const
Return the model for the supplied phase.
~BlendedInterfacialModel()
Destructor.
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.
tmp< GeometricField< Type, fvPatchField, volMesh > > F() const
Return the blended force.
tmp< volScalarField > D() const
Return the blended diffusivity.
Namespace for OpenFOAM.