nonUnityLewisEddyDiffusivity.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-2021 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 "fvcDiv.H"
28 #include "fvcLaplacian.H"
29 #include "fvcSnGrad.H"
30 #include "fvmSup.H"
31 #include "surfaceInterpolate.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace turbulenceThermophysicalTransportModels
38 {
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
42 template<class TurbulenceThermophysicalTransportModel>
45 (
46  const momentumTransportModel& momentumTransport,
47  const thermoModel& thermo
48 )
49 :
51  (
52  typeName,
53  momentumTransport,
54  thermo,
55  false
56  ),
57 
58  Sct_("Sct", dimless, this->coeffDict_)
59 {}
60 
61 
62 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
63 
64 template<class TurbulenceThermophysicalTransportModel>
65 bool
67 {
68  if
69  (
71  ::read()
72  )
73  {
74  Sct_.read(this->coeffDict());
75 
76  return true;
77  }
78  else
79  {
80  return false;
81  }
82 }
83 
84 
85 template<class TurbulenceThermophysicalTransportModel>
88 {
90  (
92  (
94  (
95  "q",
96  this->momentumTransport().alphaRhoPhi().group()
97  ),
98  -fvc::interpolate(this->alpha()*this->kappaEff())
99  *fvc::snGrad(this->thermo().T())
100  )
101  );
102 
103  const basicSpecieMixture& composition = this->thermo().composition();
104  const PtrList<volScalarField>& Y = composition.Y();
105 
106  if (Y.size())
107  {
108  surfaceScalarField hGradY
109  (
111  (
112  "hGradY",
113  Y[0].mesh(),
115  )
116  );
117 
118  forAll(Y, i)
119  {
120  const volScalarField hi
121  (
122  composition.Hs(i, this->thermo().p(), this->thermo().T())
123  );
124 
125  hGradY += fvc::interpolate(hi)*fvc::snGrad(Y[i]);
126  }
127 
128  tmpq.ref() -=
130  (
131  this->alpha()
132  *this->thermo().alphaEff((this->Prt_/Sct_)*this->alphat())
133  )*hGradY;
134  }
135 
136  return tmpq;
137 }
138 
139 
140 template<class TurbulenceThermophysicalTransportModel>
143 (
145 ) const
146 {
147  tmp<fvScalarMatrix> tmpDivq
148  (
149  fvm::Su
150  (
151  -fvc::laplacian(this->alpha()*this->kappaEff(), this->thermo().T()),
152  he
153  )
154  );
155 
156  const basicSpecieMixture& composition = this->thermo().composition();
157  const PtrList<volScalarField>& Y = composition.Y();
158 
159  tmpDivq.ref() -=
160  correction(fvm::laplacian(this->alpha()*this->alphaEff(), he));
161 
162  surfaceScalarField hGradY
163  (
165  (
166  "hGradY",
167  he.mesh(),
169  )
170  );
171 
172  forAll(Y, i)
173  {
174  const volScalarField hi
175  (
176  composition.Hs(i, this->thermo().p(), this->thermo().T())
177  );
178 
179  hGradY += fvc::interpolate(hi)*fvc::snGrad(Y[i]);
180  }
181 
182  tmpDivq.ref() -=
183  fvc::div
184  (
186  (
187  this->alpha()
188  *this->thermo().alphaEff((this->Prt_/Sct_)*this->alphat())
189  )*hGradY*he.mesh().magSf()
190  );
191 
192  return tmpDivq;
193 }
194 
195 
196 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
197 
198 } // End namespace turbulenceThermophysicalTransportModels
199 } // End namespace Foam
200 
201 // ************************************************************************* //
const char *const group
Group name for atomic constants.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
fluidReactionThermo & thermo
Definition: createFields.H:28
virtual tmp< fvScalarMatrix > divq(volScalarField &he) const
Return the source term for the energy equation.
basicSpecieMixture & composition
tmp< fvMatrix< Type > > Su(const DimensionedField< Type, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
static tmp< GeometricField< scalar, fvsPatchField, surfaceMesh > > New(const word &name, const Internal &, const PtrList< fvsPatchField< scalar >> &)
Return a temporary field constructed from name,.
Calculate the snGrad of the given volField.
Eddy-diffusivity based energy gradient heat flux model for RAS or LES of turbulent flow...
const dimensionSet dimless
Specialisation of basicMixture for a mixture consisting of a number for molecular species...
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
const dimensionSet dimLength
const dimensionSet & dimensions() const
Return dimensions.
TurbulenceThermophysicalTransportModel::momentumTransportModel momentumTransportModel
bool read(const char *, int32_t &)
Definition: int32IO.C:85
dynamicFvMesh & mesh
nonUnityLewisEddyDiffusivity(const momentumTransportModel &momentumTransport, const thermoModel &thermo)
Construct from a momentum transport model and a thermo model.
virtual tmp< surfaceScalarField > q() const
Return the heat flux [W/m^2].
static word groupName(Name name, const word &group)
Calculate the laplacian of the given field.
Calculate the divergence of the given field.
thermo he()
const Mesh & mesh() const
Return mesh.
const dimensionSet dimEnergy
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const dimensionSet dimMass
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.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
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
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
PtrList< volScalarField > & Y
virtual scalar Hs(const label speciei, const scalar p, const scalar T) const =0
Sensible enthalpy [J/kg].
A class for managing temporary objects.
Definition: PtrList.H:53
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.