unityLewisEddyDiffusivity.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) 2020-2022 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 
27 #include "fvcSnGrad.H"
28 #include "fvmLaplacian.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace turbulenceThermophysicalTransportModels
35 {
36 
37 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
38 
39 template<class TurbulenceThermophysicalTransportModel>
42 {
43  alphat_ =
44  this->momentumTransport().rho()
45  *this->momentumTransport().nut()/Prt_;
46  alphat_.correctBoundaryConditions();
47 }
48 
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
52 template<class TurbulenceThermophysicalTransportModel>
55 (
56  const momentumTransportModel& momentumTransport,
57  const thermoModel& thermo,
58  const bool allowDefaultPrt
59 )
60 :
62  (
63  typeName,
64  momentumTransport,
65  thermo,
66  allowDefaultPrt
67  )
68 {
69  this->printCoeffs(typeName);
70 }
71 
72 
73 template<class TurbulenceThermophysicalTransportModel>
76 (
77  const word& type,
78  const momentumTransportModel& momentumTransport,
79  const thermoModel& thermo,
80  const bool allowDefaultPrt
81 )
82 :
83  TurbulenceThermophysicalTransportModel
84  (
85  type,
86  momentumTransport,
87  thermo
88  ),
89 
90  Prt_
91  (
92  allowDefaultPrt
93  ? dimensioned<scalar>::lookupOrAddToDict
94  (
95  "Prt",
96  this->coeffDict_,
97  1
98  )
99  : dimensioned<scalar>
100  (
101  "Prt",
102  dimless,
103  this->coeffDict_
104  )
105  ),
106 
107  alphat_
108  (
109  IOobject
110  (
111  IOobject::groupName
112  (
113  "alphat",
114  this->momentumTransport().alphaRhoPhi().group()
115  ),
116  momentumTransport.time().name(),
117  momentumTransport.mesh(),
118  IOobject::MUST_READ,
119  IOobject::AUTO_WRITE
120  ),
121  momentumTransport.mesh()
122  )
123 {
124  if (type == typeName)
125  {
126  this->printCoeffs(type);
127  }
128 }
129 
130 
131 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
132 
133 template<class TurbulenceThermophysicalTransportModel>
135 {
137  {
138  Prt_.readIfPresent(this->coeffDict());
139 
140  return true;
141  }
142  else
143  {
144  return false;
145  }
146 }
147 
148 
149 template<class TurbulenceThermophysicalTransportModel>
150 tmp<surfaceScalarField>
152 {
154  (
156  (
157  "q",
158  this->momentumTransport().alphaRhoPhi().group()
159  ),
160  -fvc::interpolate(this->alphaEff()*this->alpha())
161  *fvc::snGrad(this->thermo().he())
162  );
163 }
164 
165 
166 template<class TurbulenceThermophysicalTransportModel>
169 (
171 ) const
172 {
173  return -fvm::laplacian(this->alpha()*this->alphaEff(), he);
174 }
175 
176 
177 template<class TurbulenceThermophysicalTransportModel>
180 (
181  const volScalarField& Yi
182 ) const
183 {
185  (
187  (
188  "j(" + Yi.name() + ')',
189  this->momentumTransport().alphaRhoPhi().group()
190  ),
191  -fvc::interpolate(this->DEff(Yi)*this->alpha())*fvc::snGrad(Yi)
192  );
193 }
194 
195 
196 template<class TurbulenceThermophysicalTransportModel>
199 (
200  volScalarField& Yi
201 ) const
202 {
203  return -fvm::laplacian(this->alpha()*this->DEff(Yi), Yi);
204 }
205 
206 
207 template<class TurbulenceThermophysicalTransportModel>
209 predict()
210 {
211  TurbulenceThermophysicalTransportModel::predict();
212  correctAlphat();
213 }
214 
215 
216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
217 
218 } // End namespace turbulenceThermophysicalTransportModels
219 } // End namespace Foam
220 
221 // ************************************************************************* //
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const word & name() const
Return name.
Definition: IOobject.H:310
static word groupName(Name name, const word &group)
Generic dimensioned Type class.
A class for managing temporary objects.
Definition: tmp.H:55
Eddy-diffusivity based energy gradient heat flux model for RAS or LES of turbulent flow....
unityLewisEddyDiffusivity(const momentumTransportModel &momentumTransport, const thermoModel &thermo, const bool allowDefaultPrt=false)
Construct from a momentum transport model and a thermo model.
virtual tmp< fvScalarMatrix > divq(volScalarField &he) const
Return the source term for the energy equation.
virtual tmp< fvScalarMatrix > divj(volScalarField &Yi) const
Return the source term for the given specie mass-fraction equation.
virtual void predict()
Correct the unityLewisEddyDiffusivity viscosity.
virtual tmp< surfaceScalarField > j(const volScalarField &Yi) const
Return the specie flux for the given specie mass-fraction [kg/m^2/s].
TurbulenceThermophysicalTransportModel::momentumTransportModel momentumTransportModel
virtual tmp< surfaceScalarField > q() const
Return the heat flux [W/m^2].
A class for handling words, derived from string.
Definition: word.H:62
Calculate the snGrad of the given volField.
Calculate the matrix for the laplacian of the field.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
const char *const group
Group name for atomic constants.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
Namespace for OpenFOAM.
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
thermo he()
fluidMulticomponentThermo & thermo
Definition: createFields.H:31