phaseModel.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) 2015-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::phaseModel
26 
27 SourceFiles
28  phaseModel.C
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #ifndef phaseModel_H
33 #define phaseModel_H
34 
35 #include "dictionary.H"
36 #include "dimensionedScalar.H"
37 #include "volFields.H"
38 #include "surfaceFields.H"
39 #include "fvMatricesFwd.H"
40 #include "rhoThermo.H"
41 #include "phaseCompressibleTurbulenceModelFwd.H"
42 #include "runTimeSelectionTables.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 
49 class phaseSystem;
50 class diameterModel;
51 
52 /*---------------------------------------------------------------------------*\
53  Class phaseModel Declaration
54 \*---------------------------------------------------------------------------*/
55 
56 class phaseModel
57 :
58  public volScalarField
59 {
60  // Private data
61 
62  //- Reference to the phaseSystem to which this phase belongs
63  const phaseSystem& fluid_;
64 
65  //- Name of phase
66  word name_;
67 
68  //- Index of phase
69  label index_;
70 
71  //- Return the residual phase-fraction for given phase
72  // Used to stabilize the phase momentum as the phase-fraction -> 0
73  dimensionedScalar residualAlpha_;
74 
75  //- Optional maximum phase-fraction (e.g. packing limit)
76  scalar alphaMax_;
77 
78  //- Diameter model
79  autoPtr<diameterModel> diameterModel_;
80 
81 
82 public:
83 
84  //- Runtime type information
85  ClassName("phaseModel");
86 
87 
88  // Declare runtime construction
89 
91  (
92  autoPtr,
93  phaseModel,
94  phaseSystem,
95  (
96  const phaseSystem& fluid,
97  const word& phaseName,
98  const label index
99  ),
100  (fluid, phaseName, index)
101  );
102 
103 
104  // Constructors
105 
106  phaseModel
107  (
108  const phaseSystem& fluid,
109  const word& phaseName,
110  const label index
111  );
112 
113  //- Return clone
114  autoPtr<phaseModel> clone() const;
115 
116 
117  // Selectors
118 
119  static autoPtr<phaseModel> New
120  (
121  const phaseSystem& fluid,
122  const word& phaseName,
123  const label index
124  );
125 
126  //- Return a pointer to a new phase created on freestore
127  // from Istream
128  class iNew
129  {
130  const phaseSystem& fluid_;
131  mutable label indexCounter_;
132 
133  public:
134 
136  (
137  const phaseSystem& fluid
138  )
139  :
140  fluid_(fluid),
141  indexCounter_(-1)
142  {}
145  {
146  indexCounter_++;
147  return autoPtr<phaseModel>
148  (
149  phaseModel::New(fluid_, word(is), indexCounter_)
150  );
151  }
152  };
153 
154 
155  //- Destructor
156  virtual ~phaseModel();
157 
158 
159  // Member Functions
160 
161  //- Return the name of this phase
162  const word& name() const;
163 
164  //- Return the name of the phase for use as the keyword in PtrDictionary
165  const word& keyword() const;
166 
167  //- Return the index of the phase
168  label index() const;
169 
170  //- Return the system to which this phase belongs
171  const phaseSystem& fluid() const;
172 
173  //- Return the residual phase-fraction for given phase
174  // Used to stabilize the phase momentum as the phase-fraction -> 0
175  const dimensionedScalar& residualAlpha() const;
176 
177  //- Return the maximum phase-fraction (e.g. packing limit)
178  scalar alphaMax() const;
179 
180  //- Return the Sauter-mean diameter
181  tmp<volScalarField> d() const;
182 
183  //- Return const-reference to diameterModel of the phase
184  const autoPtr<diameterModel>& dPtr() const;
185 
186  //- Correct the phase properties
187  virtual void correct();
188 
189  //- Correct the kinematics
190  virtual void correctKinematics();
191 
192  //- Correct the thermodynamics
193  virtual void correctThermo();
194 
195  //- Correct the turbulence
196  virtual void correctTurbulence();
197 
198  //- Correct the energy transport
199  virtual void correctEnergyTransport();
200 
201  //- Ensure that the flux at inflow/outflow BCs is preserved
203 
204  //- Read phase properties dictionary
205  virtual bool read();
206 
207 
208  // Compressibility (variable density)
209 
210  //- Return true if the phase is compressible otherwise false
211  virtual bool compressible() const = 0;
212 
213  //- Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
214  virtual tmp<volScalarField> divU() const = 0;
215 
216  //- Set the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
217  virtual void divU(tmp<volScalarField> divU) = 0;
218 
219 
220  // Thermo
221 
222  //- Return whether the phase is isothermal
223  virtual bool isothermal() const = 0;
224 
225  //- Return the enthalpy equation
226  virtual tmp<fvScalarMatrix> heEqn() = 0;
227 
228  //- Return the thermophysical model
229  virtual const rhoThermo& thermo() const = 0;
230 
231  //- Access the thermophysical model
232  virtual rhoThermo& thermoRef() = 0;
233 
234  //- Return the density field
235  virtual tmp<volScalarField> rho() const = 0;
236 
237 
238  // Species
239 
240  //- Return whether the phase is pure (i.e., not multi-component)
241  virtual bool pure() const = 0;
242 
243  //- Return the species fraction equation
244  virtual tmp<fvScalarMatrix> YiEqn(volScalarField& Yi) = 0;
245 
246  //- Return the species mass fractions
247  virtual const PtrList<volScalarField>& Y() const = 0;
248 
249  //- Return a species mass fraction by name
250  virtual const volScalarField& Y(const word& name) const = 0;
251 
252  //- Access the species mass fractions
253  virtual PtrList<volScalarField>& YRef() = 0;
254 
255  //- Return the active species mass fractions
256  virtual const UPtrList<volScalarField>& YActive() const = 0;
257 
258  //- Access the active species mass fractions
259  virtual UPtrList<volScalarField>& YActiveRef() = 0;
260 
261 
262  // Momentum
263 
264  //- Return whether the phase is stationary
265  virtual bool stationary() const = 0;
266 
267  //- Return the momentum equation
268  virtual tmp<fvVectorMatrix> UEqn() = 0;
269 
270  //- Return the momentum equation for the face-based algorithm
271  virtual tmp<fvVectorMatrix> UfEqn() = 0;
272 
273  //- Return the velocity
274  virtual tmp<volVectorField> U() const = 0;
275 
276  //- Access the velocity
277  virtual volVectorField& URef() = 0;
278 
279  //- Return the volumetric flux
280  virtual tmp<surfaceScalarField> phi() const = 0;
281 
282  //- Access the volumetric flux
283  virtual surfaceScalarField& phiRef() = 0;
284 
285  //- Return the volumetric flux of the phase
286  virtual tmp<surfaceScalarField> alphaPhi() const = 0;
287 
288  //- Access the volumetric flux of the phase
289  virtual surfaceScalarField& alphaPhiRef() = 0;
290 
291  //- Return the mass flux of the phase
292  virtual tmp<surfaceScalarField> alphaRhoPhi() const = 0;
293 
294  //- Access the mass flux of the phase
295  virtual surfaceScalarField& alphaRhoPhiRef() = 0;
296 
297  //- Return the substantive acceleration
298  virtual tmp<volVectorField> DUDt() const = 0;
299 
300  //- Return the substantive acceleration on the faces
301  virtual tmp<surfaceScalarField> DUDtf() const = 0;
302 
303  //- Return the continuity error
304  virtual tmp<volScalarField> continuityError() const = 0;
305 
306  //- Return the continuity error due to the flow field
307  virtual tmp<volScalarField> continuityErrorFlow() const = 0;
308 
309  //- Return the continuity error due to any sources
310  virtual tmp<volScalarField> continuityErrorSources() const = 0;
311 
312  //- Return the phase kinetic energy
313  virtual tmp<volScalarField> K() const = 0;
314 
315 
316  // Transport
317 
318  //- Return the laminar dynamic viscosity
319  virtual tmp<volScalarField> mu() const = 0;
320 
321  //- Return the laminar dynamic viscosity on a patch
322  virtual tmp<scalarField> mu(const label patchi) const = 0;
323 
324  //- Return the laminar kinematic viscosity
325  virtual tmp<volScalarField> nu() const = 0;
326 
327  //- Return the laminar kinematic viscosity on a patch
328  virtual tmp<scalarField> nu(const label patchi) const = 0;
329 
330  //- Thermal diffusivity for enthalpy of mixture [kg/m/s]
331  virtual tmp<volScalarField> alpha() const = 0;
332 
333  //- Thermal diffusivity for enthalpy of mixture for patch [kg/m/s]
334  virtual tmp<scalarField> alpha(const label patchi) const = 0;
335 
336  //- Thermal diffusivity for temperature of mixture [J/m/s/K]
337  virtual tmp<volScalarField> kappa() const = 0;
338 
339  //- Thermal diffusivity for temperature of mixture
340  // for patch [J/m/s/K]
341  virtual tmp<scalarField> kappa(const label patchi) const = 0;
342 
343  //- Thermal diffusivity for energy of mixture [kg/m/s]
344  virtual tmp<volScalarField> alphahe() const = 0;
345 
346  //- Thermal diffusivity for energy of mixture for patch [kg/m/s]
347  virtual tmp<scalarField> alphahe(const label patchi) const = 0;
348 
349  //- Effective thermal turbulent diffusivity for temperature
350  // of mixture [J/m/s/K]
352  (
353  const volScalarField& alphat
354  ) const = 0;
355 
356  //- Effective thermal turbulent diffusivity for temperature
357  // of mixture for patch [J/m/s/K]
358  virtual tmp<scalarField> kappaEff
359  (
360  const scalarField& alphat,
361  const label patchi
362  ) const = 0;
363 
364  //- Effective thermal turbulent diffusivity of mixture [kg/m/s]
366  (
367  const volScalarField& alphat
368  ) const = 0;
369 
370  //- Effective thermal turbulent diffusivity of mixture
371  // for patch [kg/m/s]
372  virtual tmp<scalarField> alphaEff
373  (
374  const scalarField& alphat,
375  const label patchi
376  ) const = 0;
377 
378 
379  // Turbulence
380 
381  //- Return the turbulent dynamic viscosity
382  virtual tmp<volScalarField> mut() const = 0;
383 
384  //- Return the effective dynamic viscosity
385  virtual tmp<volScalarField> muEff() const = 0;
386 
387  //- Return the turbulent kinematic viscosity
388  virtual tmp<volScalarField> nut() const = 0;
389 
390  //- Return the effective kinematic viscosity
391  virtual tmp<volScalarField> nuEff() const = 0;
392 
393  //- Effective thermal turbulent diffusivity for temperature
394  // of mixture [J/m/s/K]
395  virtual tmp<volScalarField> kappaEff() const = 0;
396 
397  //- Effective thermal turbulent diffusivity for temperature
398  // of mixture for patch [J/m/s/K]
399  virtual tmp<scalarField> kappaEff(const label patchi) const = 0;
400 
401  //- Effective thermal turbulent diffusivity of mixture [kg/m/s]
402  virtual tmp<volScalarField> alphaEff() const = 0;
403 
404  //- Effective thermal turbulent diffusivity of mixture
405  // for patch [kg/m/s]
406  virtual tmp<scalarField> alphaEff(const label patchi) const = 0;
407 
408  //- Return the turbulent kinetic energy
409  virtual tmp<volScalarField> k() const = 0;
410 
411  //- Return the phase-pressure'
412  // (derivative of phase-pressure w.r.t. phase-fraction)
413  virtual tmp<volScalarField> pPrime() const = 0;
414 };
415 
416 
417 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
418 
419 } // End namespace Foam
420 
421 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
422 
423 #endif
424 
425 // ************************************************************************* //
Foam::surfaceFields.
virtual tmp< volScalarField > mut() const =0
Return the turbulent dynamic viscosity.
virtual surfaceScalarField & phiRef()=0
Access the volumetric flux.
virtual UPtrList< volScalarField > & YActiveRef()=0
Access the active species mass fractions.
void correctInflowOutflow(surfaceScalarField &alphaPhi) const
Ensure that the flux at inflow/outflow BCs is preserved.
label index() const
Return the index of the phase.
static autoPtr< phaseModel > New(const phaseSystem &fluid, const word &phaseName, const label index)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
virtual tmp< volScalarField > continuityError() const =0
Return the continuity error.
virtual surfaceScalarField & alphaRhoPhiRef()=0
Access the mass flux of the phase.
virtual tmp< volScalarField > divU() const =0
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
virtual bool compressible() const =0
Return true if the phase is compressible otherwise false.
virtual tmp< volScalarField > alpha() const =0
Thermal diffusivity for enthalpy of mixture [kg/m/s].
virtual tmp< volScalarField > alphaEff() const =0
Effective thermal turbulent diffusivity of mixture [kg/m/s].
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
autoPtr< phaseModel > operator()(Istream &is) const
Definition: phaseModel.H:100
virtual void correctKinematics()
Correct the kinematics.
autoPtr< phaseModel > clone() const
Return clone.
const word & name() const
Definition: phaseModel.H:109
virtual tmp< volVectorField > DUDt() const =0
Return the substantive acceleration.
virtual tmp< volScalarField > muEff() const =0
Return the effective dynamic viscosity.
tmp< volScalarField > d() const
phaseModel(const word &phaseName, const volScalarField &p, const volScalarField &T)
Construct from components.
virtual bool isothermal() const =0
Return whether the phase is isothermal.
const surfaceScalarField & phi() const
Definition: phaseModel.H:190
virtual tmp< volScalarField > continuityErrorSources() const =0
Return the continuity error due to any sources.
iNew(const volScalarField &p, const volScalarField &T)
Definition: phaseModel.H:91
virtual tmp< fvScalarMatrix > YiEqn(volScalarField &Yi)=0
Return the species fraction equation.
const surfaceScalarField & alphaPhi() const
Definition: phaseModel.H:200
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual tmp< fvScalarMatrix > heEqn()=0
Return the enthalpy equation.
const volVectorField & U() const
Definition: phaseModel.H:170
virtual ~phaseModel()
Destructor.
virtual tmp< volScalarField > mu() const =0
Return the laminar dynamic viscosity.
virtual PtrList< volScalarField > & YRef()=0
Access the species mass fractions.
virtual bool read()
Read phase properties dictionary.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:67
virtual tmp< volScalarField > k() const =0
Return the turbulent kinetic energy.
A class for handling words, derived from string.
Definition: word.H:59
virtual const UPtrList< volScalarField > & YActive() const =0
Return the active species mass fractions.
const rhoThermo & thermo() const
Return const-access to phase rhoThermo.
Definition: phaseModel.H:120
virtual bool pure() const =0
Return whether the phase is pure (i.e., not multi-component)
virtual tmp< volScalarField > continuityErrorFlow() const =0
Return the continuity error due to the flow field.
const word & keyword() const
Definition: phaseModel.H:114
virtual tmp< volScalarField > K() const =0
Return the phase kinetic energy.
virtual tmp< volScalarField > pPrime() const =0
Return the phase-pressure&#39;.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: UPtrList.H:54
const dimensionedScalar & residualAlpha() const
Return the residual phase-fraction for given phase.
void correct()
Correct the laminar viscosity.
virtual const PtrList< volScalarField > & Y() const =0
Return the species mass fractions.
const dimensionedScalar & kappa() const
Definition: phaseModel.H:155
scalar alphaMax() const
Return the maximum phase-fraction (e.g. packing limit)
const dimensionedScalar & nu() const
Return the laminar viscosity.
Definition: phaseModel.H:150
virtual tmp< fvVectorMatrix > UEqn()=0
Return the momentum equation.
virtual tmp< surfaceScalarField > alphaRhoPhi() const =0
Return the mass flux of the phase.
virtual tmp< surfaceScalarField > DUDtf() const =0
Return the substantive acceleration on the faces.
Forward declarations of fvMatrix specializations.
label patchi
virtual surfaceScalarField & alphaPhiRef()=0
Access the volumetric flux of the phase.
ClassName("phaseModel")
Runtime type information.
const dimensionedScalar & rho() const
Definition: phaseModel.H:165
virtual tmp< volScalarField > nuEff() const =0
Return the effective kinematic viscosity.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
const phaseSystem & fluid() const
Return the system to which this phase belongs.
virtual tmp< volScalarField > alphahe() const =0
Thermal diffusivity for energy of mixture [kg/m/s].
const autoPtr< diameterModel > & dPtr() const
Return const-reference to diameterModel of the phase.
virtual tmp< volScalarField > kappaEff() const =0
Effective thermal turbulent diffusivity for temperature.
virtual bool stationary() const =0
Return whether the phase is stationary.
Basic thermodynamic properties based on density.
Definition: rhoThermo.H:49
virtual void correctThermo()
Correct the thermodynamics.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Macros to ease declaration of run-time selection tables.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volScalarField > nut() const =0
Return the turbulent kinematic viscosity.
virtual void correctEnergyTransport()
Correct the energy transport.
virtual void correctTurbulence()
Correct the turbulence.
declareRunTimeSelectionTable(autoPtr, phaseModel, phaseSystem,(const phaseSystem &fluid, const word &phaseName, const label index),(fluid, phaseName, index))
virtual tmp< fvVectorMatrix > UfEqn()=0
Return the momentum equation for the face-based algorithm.
virtual volVectorField & URef()=0
Access the velocity.
Namespace for OpenFOAM.
virtual rhoThermo & thermoRef()=0
Access the thermophysical model.