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