MovingPhaseModel.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 "MovingPhaseModel.H"
27 #include "phaseSystem.H"
30 #include "slipFvPatchFields.H"
32 
33 #include "fvmDdt.H"
34 #include "fvmDiv.H"
35 #include "fvmSup.H"
36 #include "fvcDdt.H"
37 #include "fvcDiv.H"
38 #include "surfaceInterpolate.H"
39 #include "fvMatrix.H"
40 
41 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
42 
43 template<class BasePhaseModel>
46 {
47  word phiName(IOobject::groupName("phi", this->name()));
48 
49  IOobject phiHeader
50  (
51  phiName,
52  U.mesh().time().timeName(),
53  U.mesh(),
54  IOobject::NO_READ
55  );
56 
57  if (phiHeader.headerOk())
58  {
59  Info<< "Reading face flux field " << phiName << endl;
60 
61  return tmp<surfaceScalarField>
62  (
64  (
65  IOobject
66  (
67  phiName,
68  U.mesh().time().timeName(),
69  U.mesh(),
70  IOobject::MUST_READ,
71  IOobject::AUTO_WRITE
72  ),
73  U.mesh()
74  )
75  );
76  }
77  else
78  {
79  Info<< "Calculating face flux field " << phiName << endl;
80 
81  wordList phiTypes
82  (
83  U.boundaryField().size(),
84  calculatedFvPatchScalarField::typeName
85  );
86 
87  forAll(U.boundaryField(), i)
88  {
89  if
90  (
91  isA<fixedValueFvPatchVectorField>(U.boundaryField()[i])
92  || isA<slipFvPatchVectorField>(U.boundaryField()[i])
93  || isA<partialSlipFvPatchVectorField>(U.boundaryField()[i])
94  )
95  {
96  phiTypes[i] = fixedValueFvPatchScalarField::typeName;
97  }
98  }
99 
100  return tmp<surfaceScalarField>
101  (
103  (
104  IOobject
105  (
106  phiName,
107  U.mesh().time().timeName(),
108  U.mesh(),
109  IOobject::NO_READ,
110  IOobject::AUTO_WRITE
111  ),
112  fvc::interpolate(U) & U.mesh().Sf(),
113  phiTypes
114  )
115  );
116  }
117 }
118 
119 
120 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
121 
122 template<class BasePhaseModel>
124 (
125  const phaseSystem& fluid,
126  const word& phaseName,
127  const label index
128 )
129 :
130  BasePhaseModel(fluid, phaseName, index),
131  U_
132  (
133  IOobject
134  (
135  IOobject::groupName("U", this->name()),
136  fluid.mesh().time().timeName(),
137  fluid.mesh(),
138  IOobject::MUST_READ,
139  IOobject::AUTO_WRITE
140  ),
141  fluid.mesh()
142  ),
143  phi_(phi(U_)),
144  alphaPhi_
145  (
146  IOobject
147  (
148  IOobject::groupName("alphaPhi", this->name()),
149  fluid.mesh().time().timeName(),
150  fluid.mesh()
151  ),
152  fluid.mesh(),
153  dimensionedScalar("0", dimensionSet(0, 3, -1, 0, 0), 0)
154  ),
155  alphaRhoPhi_
156  (
157  IOobject
158  (
159  IOobject::groupName("alphaRhoPhi", this->name()),
160  fluid.mesh().time().timeName(),
161  fluid.mesh()
162  ),
163  fluid.mesh(),
164  dimensionedScalar("0", dimensionSet(1, 0, -1, 0, 0), 0)
165  ),
166  DUDt_
167  (
168  IOobject
169  (
170  IOobject::groupName("DUDt", this->name()),
171  fluid.mesh().time().timeName(),
172  fluid.mesh()
173  ),
174  fluid.mesh(),
176  ),
177  turbulence_
178  (
180  (
181  *this,
182  this->thermo().rho(),
183  U_,
184  alphaRhoPhi_,
185  phi_,
186  *this
187  )
188  ),
189  continuityError_
190  (
191  IOobject
192  (
193  IOobject::groupName("continuityError", this->name()),
194  fluid.mesh().time().timeName(),
195  fluid.mesh()
196  ),
197  fluid.mesh(),
199  )
200 {
201  phi_.writeOpt() = IOobject::AUTO_WRITE;
202  correctKinematics();
203 }
204 
205 
206 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
207 
208 template<class BasePhaseModel>
210 {}
211 
212 
213 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
214 
215 template<class BasePhaseModel>
217 {
219 
220  this->fluid().MRF().correctBoundaryVelocity(U_);
221 
222  volScalarField& rho = this->thermo().rho();
223 
224  continuityError_ =
225  fvc::ddt(*this, rho) + fvc::div(alphaRhoPhi_)
226  - (this->fluid().fvOptions()(*this, rho)&rho);
227 }
228 
229 
230 template<class BasePhaseModel>
232 {
233  BasePhaseModel::correctKinematics();
234 
235  DUDt_ = fvc::ddt(U_) + fvc::div(phi_, U_) - fvc::div(phi_)*U_;
236 }
237 
238 
239 template<class BasePhaseModel>
241 {
242  BasePhaseModel::correctTurbulence();
243 
244  turbulence_->correct();
245 }
246 
247 
248 template<class BasePhaseModel>
250 {
251  BasePhaseModel::correctEnergyTransport();
252  turbulence_->correctEnergyTransport();
253 }
254 
255 
256 template<class BasePhaseModel>
259 {
260  return
261  (
262  fvm::ddt(*this, this->thermo().rho(), U_)
263  + fvm::div(alphaRhoPhi_, U_)
264  - fvm::Sp(continuityError_, U_)
265  + this->fluid().MRF().DDt(*this*this->thermo().rho(), U_)
266  + turbulence_->divDevRhoReff(U_)
267  );
268 }
269 
270 
271 template<class BasePhaseModel>
274 {
275  if (DbyA_.valid())
276  {
277  return DbyA_;
278  }
279  else
280  {
281  return surfaceScalarField::null();
282  }
283 }
284 
285 
286 template<class BasePhaseModel>
288 (
289  const tmp<surfaceScalarField>& DbyA
290 )
291 {
292  DbyA_ = DbyA;
293 }
294 
295 
296 template<class BasePhaseModel>
299 {
300  return U_;
301 }
302 
303 
304 template<class BasePhaseModel>
307 {
308  return U_;
309 }
310 
311 
312 template<class BasePhaseModel>
315 {
316  return DUDt_;
317 }
318 
319 
320 template<class BasePhaseModel>
323 {
324  return continuityError_;
325 }
326 
327 
328 template<class BasePhaseModel>
331 {
332  return phi_;
333 }
334 
335 
336 template<class BasePhaseModel>
339 {
340  return phi_;
341 }
342 
343 
344 template<class BasePhaseModel>
347 {
348  return alphaPhi_;
349 }
350 
351 
352 template<class BasePhaseModel>
355 {
356  return alphaPhi_;
357 }
358 
359 
360 template<class BasePhaseModel>
363 {
364  return alphaRhoPhi_;
365 }
366 
367 
368 template<class BasePhaseModel>
371 {
372  return alphaRhoPhi_;
373 }
374 
375 
376 template<class BasePhaseModel>
379 {
380  return turbulence_;
381 }
382 
383 
384 // ************************************************************************* //
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
virtual const phaseCompressibleTurbulenceModel & turbulence() const
Return the turbulence model.
virtual tmp< volVectorField > U() const
Constant access the velocity.
virtual tmp< fvVectorMatrix > UEqn()
Return the momentum equation.
virtual tmp< surfaceScalarField > alphaPhi() const
Constant access the volumetric flux of the phase.
multiphaseSystem & fluid
Definition: createFields.H:10
Calculate the matrix for implicit and explicit sources.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
virtual ~MovingPhaseModel()
Destructor.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil)-psi)*pSat, rhoMin);#1"/home/dm2/henry/OpenFOAM/OpenFOAM-3.0.x/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H"1{alphav=max(min((rho-rholSat)/(rhovSat-rholSat), scalar(1)), scalar(0));alphal=1.0-alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:74
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
Calculate the matrix for the divergence of the given field and flux.
virtual tmp< volVectorField > DUDt() const
Return the substantive acceleration.
Calulate the matrix for the first temporal derivative.
Surface Interpolation.
virtual tmp< volScalarField > continuityError() const
Constant access the continuity error.
messageStream Info
tmp< surfaceScalarField > interpolate(const RhoType &rho)
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
dynamicFvMesh & mesh
const dimensionSet dimDensity
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
virtual void correctKinematics()
Correct the kinematics.
virtual tmp< surfaceScalarField > phi() const
Constant access the volumetric flux.
IOMRFZoneList & MRF
Calculate the first temporal derivative.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual void correctEnergyTransport()
Correct the energy transport e.g. alphat.
virtual const surfaceScalarField & DbyA() const
Return the phase diffusivity divided by the momentum coefficient.
#define forAll(list, i)
Definition: UList.H:421
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
ThermalDiffusivity< PhaseCompressibleTurbulenceModel< phaseModel > > phaseCompressibleTurbulenceModel
Typedef for phaseCompressibleTurbulenceModel.
surfaceScalarField & phi
tmp< GeometricField< Type, fvPatchField, volMesh > > DDt(const surfaceScalarField &phi, const GeometricField< Type, fvPatchField, volMesh > &psi)
Definition: fvcDDt.C:45
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
MovingPhaseModel(const phaseSystem &fluid, const word &phaseName, const label index)
static const GeometricField< scalar, fvsPatchField, surfaceMesh > & null()
Return a null geometric field.
const dimensionSet dimAcceleration
Calculate the divergence of the given field.
List< word > wordList
A List of words.
Definition: fileName.H:54
virtual tmp< surfaceScalarField > alphaRhoPhi() const
Constant access the mass flux of the phase.
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
virtual void correctTurbulence()
Correct the turbulence.
psiReactionThermo & thermo
Definition: createFields.H:32
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void correct()
Correct the phase properties other than the thermo and turbulence.
A class for managing temporary objects.
Definition: PtrList.H:118
word timeName
Definition: getTimeIndex.H:3