AnisothermalPhaseModel.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 "AnisothermalPhaseModel.H"
27 #include "phaseSystem.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class BasePhaseModel>
33 (
34  const phaseSystem& fluid,
35  const word& phaseName,
36  const label index
37 )
38 :
39  BasePhaseModel(fluid, phaseName, index),
40  divU_(NULL),
41  K_
42  (
43  IOobject
44  (
45  IOobject::groupName("K", this->name()),
46  fluid.mesh().time().timeName(),
47  fluid.mesh()
48  ),
49  fluid.mesh(),
50  dimensionedScalar("K", sqr(dimVelocity), scalar(0))
51  )
52 {}
53 
54 
55 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
56 
57 template<class BasePhaseModel>
59 {}
60 
61 
62 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
63 
64 template<class BasePhaseModel>
66 {
67  BasePhaseModel::correctKinematics();
68  K_ = 0.5*magSqr(this->U());
69 }
70 
71 
72 template<class BasePhaseModel>
74 {
76 
77  this->thermo_->correct();
78 }
79 
80 
81 template<class BasePhaseModel>
84 {
85  const volScalarField& alpha = *this;
86  const surfaceScalarField& alphaPhi = this->alphaPhi();
87  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi();
88 
89  const volScalarField& contErr(this->continuityError());
90 
91  const volScalarField alphaEff(this->turbulence().alphaEff());
92 
93  volScalarField& he = this->thermo_->he();
94 
95  tmp<fvScalarMatrix> tEEqn
96  (
97  fvm::ddt(alpha, this->rho(), he) + fvm::div(alphaRhoPhi, he)
98  - fvm::Sp(contErr, he)
99 
100  + fvc::ddt(alpha, this->rho(), K_) + fvc::div(alphaRhoPhi, K_)
101  - contErr*K_
102 
104  (
105  fvc::interpolate(alpha)
107  he
108  )
109  ==
110  this->Sh()
111  );
112 
113  // Add the appropriate pressure-work term
114  if (he.name() == this->thermo_->phasePropertyName("e"))
115  {
116  tEEqn() +=
117  fvc::ddt(alpha)*this->thermo().p()
118  + fvc::div(alphaPhi, this->thermo().p());
119  }
120  else if (this->thermo_->dpdt())
121  {
122  tEEqn() -= alpha*this->fluid().dpdt();
123  }
124 
125  return tEEqn;
126 }
127 
128 
129 template<class BasePhaseModel>
131 {
132  return !this->thermo().incompressible();
133 }
134 
135 
136 template<class BasePhaseModel>
139 {
140  return divU_;
141 }
142 
143 
144 template<class BasePhaseModel>
145 void
147 (
148  const tmp<volScalarField>& divU
149 )
150 {
151  divU_ = divU;
152 }
153 
154 
155 template<class BasePhaseModel>
158 {
159  return K_;
160 }
161 
162 
163 // ************************************************************************* //
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf, const surfaceScalarField &faceFlux, Istream &schemeData)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
virtual bool compressible() const
Return true if the phase is compressible otherwise false.
scalar Sh
Definition: solveChemistry.H:2
volScalarField divU(fvc::div(fvc::absolute(phi, U)))
multiphaseSystem & fluid
Definition: createFields.H:10
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual tmp< fvScalarMatrix > heEqn()
Return the enthalpy equation.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
fluid correctThermo()
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
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
virtual void correctThermo()
Correct the thermodynamics.
autoPtr< compressible::turbulenceModel > turbulence
Definition: createFields.H:23
virtual void correctKinematics()
Correct the kinematics.
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
dynamicFvMesh & mesh
AnisothermalPhaseModel(const phaseSystem &fluid, const word &phaseName, const label index)
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
volScalarField & p
Definition: createFields.H:51
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual ~AnisothermalPhaseModel()
Destructor.
virtual const volScalarField & K() const
Return the phase kinetic energy.
surfaceScalarField alphaPhi(IOobject( "alphaPhi", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE ), phi *fvc::interpolate(alpha1))
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
psiReactionThermo & thermo
Definition: createFields.H:32
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimVelocity
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
const volScalarField & alphaEff
Definition: setAlphaEff.H:49
U
Definition: pEqn.H:82
virtual const tmp< volScalarField > & divU() const
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
word timeName
Definition: getTimeIndex.H:3