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-2016 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  K_
41  (
42  IOobject
43  (
44  IOobject::groupName("K", this->name()),
45  fluid.mesh().time().timeName(),
46  fluid.mesh()
47  ),
48  fluid.mesh(),
49  dimensionedScalar("K", sqr(dimVelocity), scalar(0))
50  )
51 {}
52 
53 
54 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
55 
56 template<class BasePhaseModel>
58 {}
59 
60 
61 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
62 
63 template<class BasePhaseModel>
65 {
66  BasePhaseModel::correctKinematics();
67  K_ = 0.5*magSqr(this->U());
68 }
69 
70 
71 template<class BasePhaseModel>
73 {
74  BasePhaseModel::correctThermo();
75 
76  this->thermo_->correct();
77 }
78 
79 
80 template<class BasePhaseModel>
83 {
84  const volScalarField& alpha = *this;
85  const surfaceScalarField& alphaPhi = this->alphaPhi();
86  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi();
87 
88  const volScalarField& contErr(this->continuityError());
89 
90  const volScalarField alphaEff(this->turbulence().alphaEff());
91 
92  volScalarField& he = this->thermo_->he();
93 
94  tmp<fvScalarMatrix> tEEqn
95  (
96  fvm::ddt(alpha, this->rho(), he) + fvm::div(alphaRhoPhi, he)
97  - fvm::Sp(contErr, he)
98 
99  + fvc::ddt(alpha, this->rho(), K_) + fvc::div(alphaRhoPhi, K_)
100  - contErr*K_
101 
103  (
104  fvc::interpolate(alpha)
106  he
107  )
108  ==
109  this->Sh()
110  );
111 
112  // Add the appropriate pressure-work term
113  if (he.name() == this->thermo_->phasePropertyName("e"))
114  {
115  tEEqn.ref() +=
116  fvc::ddt(alpha)*this->thermo().p()
117  + fvc::div(alphaPhi, this->thermo().p());
118  }
119  else if (this->thermo_->dpdt())
120  {
121  tEEqn.ref() -= alpha*this->fluid().dpdt();
122  }
123 
124  return tEEqn;
125 }
126 
127 
128 template<class BasePhaseModel>
130 {
131  return !this->thermo().incompressible();
132 }
133 
134 
135 template<class BasePhaseModel>
138 {
139  return K_;
140 }
141 
142 
143 // ************************************************************************* //
autoPtr< compressible::turbulenceModel > turbulence
Definition: createFields.H:23
virtual void correctThermo()
Correct the thermodynamics.
U
Definition: pEqn.H:83
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
multiphaseSystem & fluid
Definition: createFields.H:10
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
scalar Sh
Definition: solveChemistry.H:2
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
virtual const volScalarField & K() const
Return the phase kinetic energy.
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
psiReactionThermo & thermo
Definition: createFields.H:31
surfaceScalarField alphaPhi(IOobject("alphaPhi", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), phi *fvc::interpolate(alpha1))
dynamicFvMesh & mesh
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
AnisothermalPhaseModel(const phaseSystem &fluid, const word &phaseName, const label index)
virtual void correctKinematics()
Correct the kinematics.
word timeName
Definition: getTimeIndex.H:3
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual bool compressible() const
Return true if the phase is compressible otherwise false.
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
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
virtual ~AnisothermalPhaseModel()
Destructor.
volScalarField alphaEff("alphaEff", turbulence->nu()/Pr+alphat)
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:54
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual tmp< fvScalarMatrix > heEqn()
Return the enthalpy equation.
const dimensionSet dimVelocity