phaseModel.C
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) 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 \*---------------------------------------------------------------------------*/
25 
26 #include "phaseModel.H"
27 #include "phaseSystem.H"
28 #include "diameterModel.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(phaseModel, 0);
35  defineRunTimeSelectionTable(phaseModel, phaseSystem);
36 }
37 
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
42 (
43  const phaseSystem& fluid,
44  const word& phaseName,
45  const label index
46 )
47 :
49  (
50  IOobject
51  (
52  IOobject::groupName("alpha", phaseName),
53  fluid.mesh().time().timeName(),
54  fluid.mesh(),
55  IOobject::READ_IF_PRESENT,
56  IOobject::AUTO_WRITE
57  ),
58  fluid.mesh(),
59  dimensionedScalar("alpha", dimless, 0)
60  ),
61 
62  fluid_(fluid),
63  name_(phaseName),
64  index_(index),
65  residualAlpha_
66  (
67  "residualAlpha",
68  dimless,
69  fluid.subDict(phaseName).lookup("residualAlpha")
70  ),
71  alphaMax_(fluid.subDict(phaseName).lookupOrDefault("alphaMax", 1.0))
72 {
73  diameterModel_ = diameterModel::New(fluid.subDict(phaseName), *this);
74 }
75 
76 
78 {
79  notImplemented("phaseModel::clone() const");
80  return autoPtr<phaseModel>(NULL);
81 }
82 
83 
84 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
85 
87 {}
88 
89 
90 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
91 
93 {
94  return name_;
95 }
96 
97 
99 {
100  return name_;
101 }
102 
103 
105 {
106  return index_;
107 }
108 
109 
111 {
112  return fluid_;
113 }
114 
115 
117 {
118  return residualAlpha_;
119 }
120 
121 
122 Foam::scalar Foam::phaseModel::alphaMax() const
123 {
124  return alphaMax_;
125 }
126 
127 
129 {
130  return diameterModel_().d();
131 }
132 
133 
135 {
136  diameterModel_->correct();
137 }
138 
139 
141 {}
142 
143 
145 {}
146 
147 
149 {}
150 
151 
153 {}
154 
155 
157 {
158  return diameterModel_->read(fluid_.subDict(name_));
159 }
160 
161 
163 {
164  return false;
165 }
166 
167 
169 {
170  notImplemented("Foam::phaseModel::divU()");
171  static tmp<Foam::volScalarField> divU_(NULL);
172  return divU_;
173 }
174 
175 
176 void Foam::phaseModel::divU(const tmp<volScalarField>& divU)
177 {
178  WarningIn("phaseModel::divU(const tmp<volScalarField>& divU)")
179  << "Attempt to set the dilatation rate of an incompressible phase"
180  << endl;
181 }
182 
183 
185 {
186  notImplemented("Foam::phaseModel::K()");
187  return volScalarField::null();
188 }
189 
190 
192 {
193  return surfaceScalarField::null();
194 }
195 
196 
197 void Foam::phaseModel::DbyA(const tmp<surfaceScalarField>& DbyA)
198 {
199  WarningIn("phaseModel::DbyA(const surfaceScalarField& DbyA)")
200  << "Attempt to set the dilatation rate of an incompressible phase"
201  << endl;
202 }
203 
204 
205 // ************************************************************************* //
const word & name() const
Definition: phaseModel.H:109
static autoPtr< diameterModel > New(const dictionary &dict, const phaseModel &phase)
volScalarField divU(fvc::div(fvc::absolute(phi, U)))
const word & keyword() const
Definition: phaseModel.H:114
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:67
virtual void correctThermo()
Correct the thermodynamics.
A class for handling words, derived from string.
Definition: word.H:59
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 void correctKinematics()
Correct the kinematics.
dynamicFvMesh & mesh
const phaseSystem & fluid() const
Return the system to which this phase belongs.
Namespace for OpenFOAM.
virtual const surfaceScalarField & DbyA() const
Return the phase diffusivity divided by the momentum coefficient.
autoPtr< phaseModel > clone() const
Return clone.
void correct()
Correct the laminar viscosity.
virtual ~phaseModel()
Destructor.
Generic dimensioned Type class.
virtual const tmp< volScalarField > & divU() const
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
stressControl lookup("compactNormalStress") >> compactNormalStress
#define WarningIn(functionName)
Report a warning using Foam::Warning.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual void correctTurbulence()
Correct the turbulence.
tmp< volScalarField > d() const
virtual void correctEnergyTransport()
Correct the energy transport e.g. alphat.
const dimensionedScalar & residualAlpha() const
Return the residual phase-fraction for given phase.
static const GeometricField< scalar, fvPatchField, volMesh > & null()
Return a null geometric field.
virtual const volScalarField & K() const
Return the phase kinetic energy.
virtual bool read()
Read phase properties dictionary.
label index() const
Return the index of the phase.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
phaseModel(const word &phaseName, const volScalarField &p, const volScalarField &T)
Construct from components.
virtual bool compressible() const
Return true if the phase is compressible otherwise false.
#define notImplemented(functionName)
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:356
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
A class for managing temporary objects.
Definition: PtrList.H:118
defineTypeNameAndDebug(combustionModel, 0)
word timeName
Definition: getTimeIndex.H:3
scalar alphaMax() const
Return the maximum phase-fraction (e.g. packing limit)