phaseModel.C
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) 2011-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 \*---------------------------------------------------------------------------*/
25 
26 #include "phaseModel.H"
27 #include "twoPhaseSystem.H"
28 #include "diameterModel.H"
29 #include "fvMatrix.H"
31 #include "dragModel.H"
32 #include "heatTransferModel.H"
35 #include "slipFvPatchFields.H"
37 #include "fvcFlux.H"
38 #include "surfaceInterpolate.H"
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const twoPhaseSystem& fluid,
46  const dictionary& phaseProperties,
47  const word& phaseName
48 )
49 :
51  (
52  IOobject
53  (
54  IOobject::groupName("alpha", phaseName),
55  fluid.mesh().time().timeName(),
56  fluid.mesh(),
57  IOobject::READ_IF_PRESENT,
58  IOobject::AUTO_WRITE
59  ),
60  fluid.mesh(),
61  dimensionedScalar("alpha", dimless, 0)
62  ),
63  fluid_(fluid),
64  name_(phaseName),
65  phaseDict_
66  (
67  phaseProperties.subDict(name_)
68  ),
69  residualAlpha_
70  (
71  "residualAlpha",
72  dimless,
73  fluid.subDict(phaseName).lookup("residualAlpha")
74  ),
75  alphaMax_(phaseDict_.lookupOrDefault("alphaMax", 1.0)),
76  thermo_(rhoThermo::New(fluid.mesh(), name_)),
77  U_
78  (
79  IOobject
80  (
81  IOobject::groupName("U", name_),
82  fluid.mesh().time().timeName(),
83  fluid.mesh(),
84  IOobject::MUST_READ,
85  IOobject::AUTO_WRITE
86  ),
87  fluid.mesh()
88  ),
89  alphaPhi_
90  (
91  IOobject
92  (
93  IOobject::groupName("alphaPhi", name_),
94  fluid.mesh().time().timeName(),
95  fluid.mesh()
96  ),
97  fluid.mesh(),
98  dimensionedScalar("0", dimensionSet(0, 3, -1, 0, 0), 0)
99  ),
100  alphaRhoPhi_
101  (
102  IOobject
103  (
104  IOobject::groupName("alphaRhoPhi", name_),
105  fluid.mesh().time().timeName(),
106  fluid.mesh()
107  ),
108  fluid.mesh(),
109  dimensionedScalar("0", dimensionSet(1, 0, -1, 0, 0), 0)
110  )
111 {
112  thermo_->validate("phaseModel " + name_, "h", "e");
113 
114  const word phiName = IOobject::groupName("phi", name_);
115 
116  IOobject phiHeader
117  (
118  phiName,
119  fluid_.mesh().time().timeName(),
120  fluid_.mesh(),
122  );
123 
124  if (phiHeader.typeHeaderOk<surfaceScalarField>(true))
125  {
126  Info<< "Reading face flux field " << phiName << endl;
127 
128  phiPtr_.reset
129  (
131  (
132  IOobject
133  (
134  phiName,
135  fluid_.mesh().time().timeName(),
136  fluid_.mesh(),
139  ),
140  fluid_.mesh()
141  )
142  );
143  }
144  else
145  {
146  Info<< "Calculating face flux field " << phiName << endl;
147 
148  wordList phiTypes
149  (
150  U_.boundaryField().size(),
151  calculatedFvPatchScalarField::typeName
152  );
153 
154  forAll(U_.boundaryField(), i)
155  {
156  if
157  (
158  isA<fixedValueFvPatchVectorField>(U_.boundaryField()[i])
159  || isA<slipFvPatchVectorField>(U_.boundaryField()[i])
160  || isA<partialSlipFvPatchVectorField>(U_.boundaryField()[i])
161  )
162  {
163  phiTypes[i] = fixedValueFvsPatchScalarField::typeName;
164  }
165  }
166 
167  phiPtr_.reset
168  (
170  (
171  IOobject
172  (
173  phiName,
174  fluid_.mesh().time().timeName(),
175  fluid_.mesh(),
178  ),
179  fvc::flux(U_),
180  phiTypes
181  )
182  );
183  }
184 
185  dPtr_ = diameterModel::New
186  (
187  phaseDict_,
188  *this
189  );
190 
191  turbulence_ =
193  (
194  *this,
195  thermo_->rho(),
196  U_,
197  alphaRhoPhi_,
198  phi(),
199  *this
200  );
201 }
202 
203 
204 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
205 
207 {}
208 
209 
210 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
211 
213 {
214  return fluid_.otherPhase(*this);
215 }
216 
217 
219 {
220  return dPtr_().d();
221 }
222 
223 
226 {
227  return turbulence_();
228 }
229 
230 
233 {
234  return turbulence_();
235 }
236 
237 
239 {
240  return dPtr_->correct();
241 }
242 
243 
244 bool Foam::phaseModel::read(const dictionary& phaseProperties)
245 {
246  phaseDict_ = phaseProperties.subDict(name_);
247  return dPtr_->read(phaseDict_);
248 }
249 
250 
252 {
253  surfaceScalarField::Boundary& alphaPhiBf = alphaPhi.boundaryFieldRef();
254  const volScalarField::Boundary& alphaBf = boundaryField();
255  const surfaceScalarField::Boundary& phiBf = phi().boundaryField();
256 
257  forAll(alphaPhiBf, patchi)
258  {
259  fvsPatchScalarField& alphaPhip = alphaPhiBf[patchi];
260 
261  if (!alphaPhip.coupled())
262  {
263  alphaPhip = phiBf[patchi]*alphaBf[patchi];
264  }
265  }
266 }
267 
268 
269 // ************************************************************************* //
fvsPatchField< scalar > fvsPatchScalarField
void correctInflowOutflow(surfaceScalarField &alphaPhi) const
Ensure that the flux at inflow/outflow BCs is preserved.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const PhaseCompressibleTurbulenceModel< phaseModel > & turbulence() const
Return the turbulence model.
const surfaceScalarField & phi() const
Return the volumetric flux.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
IOobject(const word &name, const fileName &instance, const objectRegistry &registry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Definition: IOobject.C:209
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
tmp< volScalarField > d() const
phaseModel(const word &phaseName, const volScalarField &p, const volScalarField &T)
Construct from components.
static autoPtr< diameterModel > New(const dictionary &dict, const phaseModel &phase)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
dictionary()
Construct top-level dictionary null.
Definition: dictionary.C:238
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual ~phaseModel()
Destructor.
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
virtual bool read()
Read phase properties dictionary.
Calculate the face-flux of the given field.
static word groupName(Name name, const word &group)
const phaseModel & otherPhase() const
Return the other phase in this two-phase system.
word timeName
Definition: getTimeIndex.H:3
void correct()
Correct the laminar viscosity.
static autoPtr< PhaseCompressibleTurbulenceModel > New(const alphaField &alpha, const volScalarField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &trasportModel, const word &propertiesName=turbulenceModel::propertiesName)
Return a reference to the selected turbulence model.
List< word > wordList
A List of words.
Definition: fileName.H:54
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
messageStream Info
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53