dynamicKEqn.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) 2011-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 
26 #include "dynamicKEqn.H"
27 #include "fvModels.H"
28 #include "fvConstraints.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace LESModels
35 {
36 
37 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
38 
39 template<class BasicMomentumTransportModel>
42 {
43  return max
44  (
45  0.5*(filter_(magSqr(this->U_)) - magSqr(filter_(this->U_))),
47  );
48 }
49 
50 
51 template<class BasicMomentumTransportModel>
53 (
54  const volSymmTensorField& D,
55  const volScalarField& KK
56 ) const
57 {
58  const volSymmTensorField LL
59  (
60  simpleFilter_(dev(filter_(sqr(this->U_)) - (sqr(filter_(this->U_)))))
61  );
62 
63  const volSymmTensorField MM
64  (
65  simpleFilter_(-2.0*this->delta()*sqrt(KK)*filter_(D))
66  );
67 
68  const volScalarField Ck
69  (
70  simpleFilter_(0.5*(LL && MM))
71  /(
72  simpleFilter_(magSqr(MM))
73  + dimensionedScalar(sqr(MM.dimensions()), vSmall)
74  )
75  );
76 
77  tmp<volScalarField> tfld = 0.5*(mag(Ck) + Ck);
78  return tfld();
79 }
80 
81 
82 template<class BasicMomentumTransportModel>
84 (
85  const volSymmTensorField& D,
86  const volScalarField& KK
87 ) const
88 {
89  const volScalarField Ce
90  (
91  simpleFilter_(this->nuEff()*(filter_(magSqr(D)) - magSqr(filter_(D))))
92  /simpleFilter_(pow(KK, 1.5)/(2.0*this->delta()))
93  );
94 
95  tmp<volScalarField> tfld = 0.5*(mag(Ce) + Ce);
96  return tfld();
97 }
98 
99 
100 template<class BasicMomentumTransportModel>
102 {
103  const volSymmTensorField D(dev(symm(fvc::grad(this->U_))));
104  return Ce(D, KK());
105 }
106 
107 
108 template<class BasicMomentumTransportModel>
110 (
111  const volSymmTensorField& D,
112  const volScalarField& KK
113 )
114 {
115  this->nut_ = Ck(D, KK)*sqrt(k_)*this->delta();
116  this->nut_.correctBoundaryConditions();
117  fvConstraints::New(this->mesh_).constrain(this->nut_);
118 }
119 
120 
121 template<class BasicMomentumTransportModel>
123 {
124  correctNut(symm(fvc::grad(this->U_)), KK());
125 }
126 
127 
128 template<class BasicMomentumTransportModel>
130 {
131  return tmp<fvScalarMatrix>
132  (
133  new fvScalarMatrix
134  (
135  k_,
136  dimVolume*this->rho_.dimensions()*k_.dimensions()
137  /dimTime
138  )
139  );
140 }
141 
142 
143 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
144 
145 template<class BasicMomentumTransportModel>
147 (
148  const alphaField& alpha,
149  const rhoField& rho,
150  const volVectorField& U,
151  const surfaceScalarField& alphaRhoPhi,
152  const surfaceScalarField& phi,
153  const transportModel& transport,
154  const word& type
155 )
156 :
158  (
159  type,
160  alpha,
161  rho,
162  U,
163  alphaRhoPhi,
164  phi,
165  transport
166  ),
167 
168  k_
169  (
170  IOobject
171  (
172  IOobject::groupName("k", this->alphaRhoPhi_.group()),
173  this->runTime_.timeName(),
174  this->mesh_,
177  ),
178  this->mesh_
179  ),
180 
181  simpleFilter_(this->mesh_),
182  filterPtr_(LESfilter::New(this->mesh_, this->coeffDict())),
183  filter_(filterPtr_())
184 {
185  bound(k_, this->kMin_);
186 
187  if (type == typeName)
188  {
189  this->printCoeffs(type);
190  }
191 }
192 
193 
194 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
195 
196 template<class BasicMomentumTransportModel>
198 {
200  {
201  filter_.read(this->coeffDict());
202 
203  return true;
204  }
205  else
206  {
207  return false;
208  }
209 }
210 
211 
212 template<class BasicMomentumTransportModel>
214 {
215  return volScalarField::New
216  (
217  IOobject::groupName("epsilon", this->alphaRhoPhi_.group()),
218  Ce()*k()*sqrt(k())/this->delta()
219  );
220 }
221 
222 
223 template<class BasicMomentumTransportModel>
225 {
226  if (!this->turbulence_)
227  {
228  return;
229  }
230 
231  // Local references
232  const alphaField& alpha = this->alpha_;
233  const rhoField& rho = this->rho_;
234  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
235  const volVectorField& U = this->U_;
236  volScalarField& nut = this->nut_;
237  const Foam::fvModels& fvModels(Foam::fvModels::New(this->mesh_));
239  (
240  Foam::fvConstraints::New(this->mesh_)
241  );
242 
244 
246 
247  tmp<volTensorField> tgradU(fvc::grad(U));
248  const volSymmTensorField D(dev(symm(tgradU())));
249  const volScalarField G(this->GName(), 2.0*nut*(tgradU() && D));
250  tgradU.clear();
251 
252  const volScalarField KK(this->KK());
253 
255  (
256  fvm::ddt(alpha, rho, k_)
257  + fvm::div(alphaRhoPhi, k_)
258  - fvm::laplacian(alpha*rho*DkEff(), k_)
259  ==
260  alpha*rho*G
261  - fvm::SuSp((2.0/3.0)*alpha*rho*divU, k_)
262  - fvm::Sp(Ce(D, KK)*alpha*rho*sqrt(k_)/this->delta(), k_)
263  + kSource()
264  + fvModels.source(alpha, rho, k_)
265  );
266 
267  kEqn.ref().relax();
268  fvConstraints.constrain(kEqn.ref());
269  solve(kEqn);
270  fvConstraints.constrain(k_);
271  bound(k_, this->kMin_);
272 
273  correctNut(D, KK);
274 }
275 
276 
277 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
278 
279 } // End namespace LESModels
280 } // End namespace Foam
281 
282 // ************************************************************************* //
Eddy viscosity LES SGS model base class.
scalar delta
virtual bool read()
Read model coefficients if they have changed.
Definition: dynamicKEqn.C:197
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:237
virtual void correct()
Correct Eddy-Viscosity and related properties.
Definition: dynamicKEqn.C:224
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
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, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
dimensionedScalar sqrt(const dimensionedScalar &ds)
label k
Boltzmann constant.
const dimensionedScalar G
Newtonian constant of gravitation.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
One equation eddy-viscosity model.
Definition: kEqn.H:71
BasicMomentumTransportModel::transportModel transportModel
Definition: dynamicKEqn.H:130
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
const dimensionSet dimTime
const dimensionSet & dimensions() const
Return dimensions.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
volScalarField Ck(const volSymmTensorField &D, const volScalarField &KK) const
Calculate Ck by filtering the velocity field U.
Definition: dynamicKEqn.C:53
Foam::fvConstraints & fvConstraints
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
Finite volume models.
Definition: fvModels.H:60
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
phi
Definition: correctPhi.H:3
volScalarField Ce() const
Definition: dynamicKEqn.C:101
dimensioned< scalar > magSqr(const dimensioned< Type > &)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
const dimensionSet dimVelocity
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:96
tmp< volScalarField > divU
Definition: alphaSuSp.H:9
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
Foam::fvModels & fvModels
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
U
Definition: pEqn.H:72
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
Finite volume constraints.
Definition: fvConstraints.H:57
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
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
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual tmp< volScalarField > epsilon() const
Return sub-grid disipation rate.
Definition: dynamicKEqn.C:213
const dimensionSet dimVolume
virtual tmp< fvScalarMatrix > kSource() const
Definition: dynamicKEqn.C:129
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > mag(const dimensioned< Type > &)
tmp< volScalarField > KK() const
Return the KK field obtained by filtering the velocity field U.
Definition: dynamicKEqn.C:41
A class for managing temporary objects.
Definition: PtrList.H:53
BasicMomentumTransportModel::alphaField alphaField
Definition: dynamicKEqn.H:128
static autoPtr< LESfilter > New(const fvMesh &, const dictionary &, const word &filterDictName="filter")
Return a reference to the selected LES filter.
Definition: LESfilter.C:41
dynamicKEqn(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Construct from components.
Definition: dynamicKEqn.C:147
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
scalar nut
Namespace for OpenFOAM.
BasicMomentumTransportModel::rhoField rhoField
Definition: dynamicKEqn.H:129