phaseModel.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-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Class
27 Foam::phaseModel
28
29SourceFiles
30 phaseModel.C
31
32\*---------------------------------------------------------------------------*/
33
34#ifndef phaseModel_H
35#define phaseModel_H
36
37#include "dictionary.H"
38#include "dimensionedScalar.H"
39#include "volFields.H"
40#include "surfaceFields.H"
41#include "transportModel.H"
42#include "rhoThermo.H"
43
44// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45
46namespace Foam
47{
48
49// Forward declarations
50class twoPhaseSystem;
51class diameterModel;
52
53template<class Phase>
54class PhaseCompressibleTurbulenceModel;
55
56
57/*---------------------------------------------------------------------------*\
58 Class phaseModel Declaration
59\*---------------------------------------------------------------------------*/
60
61class phaseModel
62:
63 public volScalarField,
64 public transportModel
65{
66 // Private data
67
68 //- Reference to the twoPhaseSystem to which this phase belongs
69 const twoPhaseSystem& fluid_;
70
71 //- Name of phase
72 word name_;
73
74 dictionary phaseDict_;
75
76 //- Return the residual phase-fraction for given phase
77 // Used to stabilize the phase momentum as the phase-fraction -> 0
78 dimensionedScalar residualAlpha_;
79
80 //- Optional maximum phase-fraction (e.g. packing limit)
81 scalar alphaMax_;
82
83 //- Thermophysical properties
84 autoPtr<rhoThermo> thermo_;
85
86 //- Velocity
88
89 //- Volumetric flux of the phase
90 surfaceScalarField alphaPhi_;
91
92 //- Mass flux of the phase
93 surfaceScalarField alphaRhoPhi_;
94
95 //- Volumetric flux of the phase
96 autoPtr<surfaceScalarField> phiPtr_;
97
98 //- Diameter model
99 autoPtr<diameterModel> dPtr_;
100
101 //- Turbulence model
102 autoPtr<PhaseCompressibleTurbulenceModel<phaseModel>> turbulence_;
103
104
105public:
106
107 // Constructors
108
110 (
111 const twoPhaseSystem& fluid,
112 const dictionary& phaseProperties,
113 const word& phaseName
114 );
115
116
117 //- Destructor
118 virtual ~phaseModel();
119
120
121 // Member Functions
122
123 //- Return the name of this phase
124 const word& name() const
125 {
126 return name_;
127 }
128
129 //- Return the twoPhaseSystem to which this phase belongs
130 const twoPhaseSystem& fluid() const
131 {
132 return fluid_;
133 }
134
135 //- Return the other phase in this two-phase system
136 const phaseModel& otherPhase() const;
137
138 //- Return the residual phase-fraction for given phase
139 // Used to stabilize the phase momentum as the phase-fraction -> 0
140 const dimensionedScalar& residualAlpha() const
141 {
142 return residualAlpha_;
143 }
144
145 //- Optional maximum phase-fraction (e.g. packing limit)
146 // Defaults to 1
147 scalar alphaMax() const
148 {
149 return alphaMax_;
150 }
151
152 //- Return the Sauter-mean diameter
153 tmp<volScalarField> d() const;
154
155 //- Return the turbulence model
157 turbulence() const;
158
159 //- Return non-const access to the turbulence model
160 // for correction
162 turbulence();
163
164 //- Return the thermophysical model
165 const rhoThermo& thermo() const
166 {
167 return *thermo_;
168 }
169
170 //- Return non-const access to the thermophysical model
171 // for correction
173 {
174 return *thermo_;
175 }
176
177 //- Return the laminar viscosity
178 tmp<volScalarField> nu() const
179 {
180 return thermo_->nu();
181 }
182
183 //- Return the laminar viscosity for patch
184 tmp<scalarField> nu(const label patchi) const
185 {
186 return thermo_->nu(patchi);
187 }
188
189 //- Return the laminar dynamic viscosity
190 tmp<volScalarField> mu() const
191 {
192 return thermo_->mu();
193 }
194
195 //- Return the laminar dynamic viscosity for patch
196 tmp<scalarField> mu(const label patchi) const
197 {
198 return thermo_->mu(patchi);
199 }
200
201 //- Return the thermal conductivity on a patch
202 tmp<scalarField> kappa(const label patchi) const
203 {
204 return thermo_->kappa(patchi);
205 }
206
207 //- Return the thermal conductivity
209 {
210 return thermo_->kappa();
211 }
212
213 //- Thermal diffusivity for energy of mixture [kg/m/s]
215 {
216 return thermo_->alphahe();
217 }
218
219 //- Thermal diffusivity for energy of mixture for patch [kg/m/s]
220 tmp<scalarField> alphahe(const label patchi) const
221 {
222 return thermo_->alphahe(patchi);
223 }
224
225 //- Return the laminar thermal conductivity
227 (
228 const volScalarField& alphat
229 ) const
230 {
231 return thermo_->kappaEff(alphat);
232 }
233
234 //- Return the laminar thermal conductivity on a patch
236 (
237 const scalarField& alphat,
238 const label patchi
239 ) const
240 {
241 return thermo_->kappaEff(alphat, patchi);
242 }
243
244 //- Return the laminar thermal diffusivity for enthalpy
246 {
247 return thermo_->alpha();
248 }
249
250 //- Return the laminar thermal diffusivity for enthalpy on a patch
251 tmp<scalarField> alpha(const label patchi) const
252 {
253 return thermo_->alpha(patchi);
254 }
255
256 //- Return the effective thermal diffusivity for enthalpy
258 (
259 const volScalarField& alphat
260 ) const
261 {
262 return thermo_->alphaEff(alphat);
263 }
264
265 //- Return the effective thermal diffusivity for enthalpy on a patch
267 (
268 const scalarField& alphat,
269 const label patchi
270 ) const
271 {
272 return thermo_->alphaEff(alphat, patchi);
273 }
274
275 //- Return the specific heat capacity
276 tmp<volScalarField> Cp() const
277 {
278 return thermo_->Cp();
279 }
280
281 //- Return the density
282 const volScalarField& rho() const
283 {
284 return thermo_->rho();
285 }
286
287 //- Return the velocity
288 const volVectorField& U() const
289 {
290 return U_;
291 }
292
293 //- Return non-const access to the velocity
294 // Used in the momentum equation
296 {
297 return U_;
298 }
299
300 //- Return the volumetric flux
301 const surfaceScalarField& phi() const
302 {
303 return *phiPtr_;
304 }
305
306 //- Return non-const access to the volumetric flux
308 {
309 return *phiPtr_;
310 }
311
312 //- Return the volumetric flux of the phase
313 const surfaceScalarField& alphaPhi() const
314 {
315 return alphaPhi_;
316 }
317
318 //- Return non-const access to the volumetric flux of the phase
320 {
321 return alphaPhi_;
322 }
323
324 //- Return the mass flux of the phase
325 const surfaceScalarField& alphaRhoPhi() const
326 {
327 return alphaRhoPhi_;
328 }
329
330 //- Return non-const access to the mass flux of the phase
332 {
333 return alphaRhoPhi_;
334 }
335
336 //- Ensure that the flux at inflow/outflow BCs is preserved
338
339 //- Correct the phase properties
340 // other than the thermodynamics and turbulence
341 // which have special treatment
342 void correct();
343
344 //- Read phaseProperties dictionary
345 virtual bool read(const dictionary& phaseProperties);
346
347 //- Dummy Read for transportModel
348 virtual bool read()
349 {
350 return true;
351 }
352};
353
354
355// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
356
357} // End namespace Foam
358
359// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
360
361#endif
362
363// ************************************************************************* //
Templated abstract base class for multiphase compressible turbulence models.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:61
const rhoThermo & thermo() const
Return the thermophysical model.
Definition: phaseModel.H:164
const volVectorField & U() const
Return the velocity.
Definition: phaseModel.H:287
tmp< scalarField > alpha(const label patchi) const
Return the laminar thermal diffusivity for enthalpy on a patch.
Definition: phaseModel.H:250
virtual bool read(const dictionary &phaseProperties)
Read phaseProperties dictionary.
const phaseSystem & fluid() const
Return the system to which this phase belongs.
Definition: phaseModel.C:148
const PhaseCompressibleTurbulenceModel< phaseModel > & turbulence() const
Return the turbulence model.
Definition: phaseModel.C:238
const volScalarField & rho() const
Return the density.
Definition: phaseModel.H:281
surfaceScalarField & alphaRhoPhi()
Return non-const access to the mass flux of the phase.
Definition: phaseModel.H:330
tmp< scalarField > kappaEff(const scalarField &alphat, const label patchi) const
Return the laminar thermal conductivity on a patch.
Definition: phaseModel.H:235
tmp< volScalarField > Cp() const
Return the specific heat capacity.
Definition: phaseModel.H:275
const surfaceScalarField & alphaRhoPhi() const
Return the mass flux of the phase.
Definition: phaseModel.H:324
const dimensionedScalar & residualAlpha() const
Return the residual phase-fraction for given phase.
Definition: phaseModel.C:154
tmp< volScalarField > kappa() const
Return the thermal conductivity.
Definition: phaseModel.H:207
tmp< scalarField > mu(const label patchi) const
Return the laminar dynamic viscosity for patch.
Definition: phaseModel.H:195
tmp< scalarField > alphahe(const label patchi) const
Thermal diffusivity for energy of mixture for patch [kg/m/s].
Definition: phaseModel.H:219
tmp< volScalarField > d() const
Return the Sauter-mean diameter.
virtual tmp< surfaceScalarField > alphaRhoPhi() const =0
Return the mass flux of the phase.
void correct()
Correct the phase properties.
Definition: phaseModel.C:216
virtual ~phaseModel()
Destructor.
Definition: phaseModel.C:202
phaseModel(const word &phaseName, const dictionary &phaseDict, const fvMesh &mesh)
Definition: phaseModel.C:40
const surfaceScalarField & phi() const
Return the volumetric flux.
Definition: phaseModel.H:300
virtual tmp< volScalarField > alphahe() const =0
Thermal diffusivity for energy of mixture [kg/m/s].
const word & name() const
Return the name of this phase.
Definition: phaseModel.H:123
tmp< volScalarField > kappaEff(const volScalarField &alphat) const
Return the laminar thermal conductivity.
Definition: phaseModel.H:226
surfaceScalarField & phi()
Return non-const access to the volumetric flux.
Definition: phaseModel.H:306
scalar alphaMax() const
Return the maximum phase-fraction (e.g. packing limit)
Definition: phaseModel.C:160
const dimensionedScalar & rho() const
Definition: phaseModel.H:171
const phaseModel & otherPhase() const
Return the other phase in this two-phase system.
Definition: phaseModel.C:218
tmp< scalarField > kappa(const label patchi) const
Return the thermal conductivity on a patch.
Definition: phaseModel.H:201
virtual tmp< volScalarField > alpha() const =0
Thermal diffusivity for enthalpy of mixture [kg/m/s].
tmp< scalarField > nu(const label patchi) const
Return the laminar viscosity for patch.
Definition: phaseModel.H:183
const dimensionedScalar & nu() const
Return the laminar viscosity.
Definition: phaseModel.H:156
tmp< volScalarField > d() const
Definition: phaseModel.C:258
tmp< volScalarField > nu() const
Return the laminar viscosity.
Definition: phaseModel.H:177
tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Return the effective thermal diffusivity for enthalpy.
Definition: phaseModel.H:257
const surfaceScalarField & alphaPhi() const
Return the volumetric flux of the phase.
Definition: phaseModel.H:312
tmp< volScalarField > mu() const
Return the laminar dynamic viscosity.
Definition: phaseModel.H:189
const dimensionedScalar & Cp() const
Definition: phaseModel.H:166
tmp< scalarField > alphaEff(const scalarField &alphat, const label patchi) const
Return the effective thermal diffusivity for enthalpy on a patch.
Definition: phaseModel.H:266
surfaceScalarField & alphaPhi()
Return non-const access to the volumetric flux of the phase.
Definition: phaseModel.H:318
virtual tmp< volScalarField > kappaEff() const =0
Effective thermal turbulent diffusivity for temperature.
virtual bool read()
Dummy Read for transportModel.
Definition: phaseModel.H:347
tmp< volScalarField > alpha() const
Return the laminar thermal diffusivity for enthalpy.
Definition: phaseModel.H:244
void correctInflowOutflow(surfaceScalarField &alphaPhi) const
Ensure that the flux at inflow/outflow BCs is preserved.
tmp< volScalarField > alphahe() const
Thermal diffusivity for energy of mixture [kg/m/s].
Definition: phaseModel.H:213
virtual ~phaseModel()
Destructor.
const word & name() const
Definition: phaseModel.H:144
virtual const rhoThermo & thermo() const =0
Return the thermophysical model.
virtual tmp< volScalarField > alphaEff() const =0
Effective thermal turbulent diffusivity of mixture [kg/m/s].
virtual tmp< volScalarField > mu() const =0
Return the laminar dynamic viscosity.
const twoPhaseSystem & fluid() const
Return the twoPhaseSystem to which this phase belongs.
Definition: phaseModel.H:129
virtual bool read()
Read phase properties dictionary.
Definition: phaseModel.C:200
const dimensionedScalar & kappa() const
Definition: phaseModel.H:161
Helper class to manage multi-specie phase properties.
Basic thermodynamic properties based on density.
Definition: rhoThermo.H:58
A class for managing temporary objects.
Definition: tmp.H:65
Class which solves the volume fraction equations for two phases.
A class for handling words, derived from Foam::string.
Definition: word.H:68
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Foam::surfaceFields.