dynamicLagrangian.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-2020 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 "dynamicLagrangian.H"
27 #include "fvOptions.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace LESModels
34 {
35 
36 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
37 
38 template<class BasicMomentumTransportModel>
40 (
41  const tmp<volTensorField>& gradU
42 )
43 {
44  this->nut_ = (flm_/fmm_)*sqr(this->delta())*mag(dev(symm(gradU)));
45  this->nut_.correctBoundaryConditions();
46  fv::options::New(this->mesh_).correct(this->nut_);
47 }
48 
49 
50 template<class BasicMomentumTransportModel>
52 {
53  correctNut(fvc::grad(this->U_));
54 }
55 
56 
57 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
58 
59 template<class BasicMomentumTransportModel>
61 (
62  const alphaField& alpha,
63  const rhoField& rho,
64  const volVectorField& U,
65  const surfaceScalarField& alphaRhoPhi,
66  const surfaceScalarField& phi,
67  const transportModel& transport,
68  const word& type
69 )
70 :
72  (
73  type,
74  alpha,
75  rho,
76  U,
77  alphaRhoPhi,
78  phi,
79  transport
80  ),
81 
82  flm_
83  (
84  IOobject
85  (
86  IOobject::groupName("flm", this->alphaRhoPhi_.group()),
87  this->runTime_.timeName(),
88  this->mesh_,
91  ),
92  this->mesh_
93  ),
94  fmm_
95  (
96  IOobject
97  (
98  IOobject::groupName("fmm", this->alphaRhoPhi_.group()),
99  this->runTime_.timeName(),
100  this->mesh_,
103  ),
104  this->mesh_
105  ),
106  theta_
107  (
109  (
110  "theta",
111  this->coeffDict_,
112  1.5
113  )
114  ),
115 
116  simpleFilter_(U.mesh()),
117  filterPtr_(LESfilter::New(U.mesh(), this->coeffDict())),
118  filter_(filterPtr_()),
119 
120  flm0_("flm0", flm_.dimensions(), 0.0),
121  fmm0_("fmm0", fmm_.dimensions(), vSmall)
122 {
123  if (type == typeName)
124  {
125  this->printCoeffs(type);
126  }
127 }
128 
129 
130 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
131 
132 template<class BasicMomentumTransportModel>
134 {
136  {
137  filter_.read(this->coeffDict());
138  theta_.readIfPresent(this->coeffDict());
139 
140  return true;
141  }
142  else
143  {
144  return false;
145  }
146 }
147 
148 
149 template<class BasicMomentumTransportModel>
151 {
152  if (!this->turbulence_)
153  {
154  return;
155  }
156 
157  // Local references
158  const alphaField& alpha = this->alpha_;
159  const rhoField& rho = this->rho_;
160  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
161  const volVectorField& U = this->U_;
162  fv::options& fvOptions(fv::options::New(this->mesh_));
163 
165 
166  tmp<volTensorField> tgradU(fvc::grad(U));
167  const volTensorField& gradU = tgradU();
168 
169  volSymmTensorField S(dev(symm(gradU)));
170  volScalarField magS(mag(S));
171 
172  volVectorField Uf(filter_(U));
174  volScalarField magSf(mag(Sf));
175 
176  volSymmTensorField L(dev(filter_(sqr(U)) - (sqr(filter_(U)))));
178  (
179  2.0*sqr(this->delta())*(filter_(magS*S) - 4.0*magSf*Sf)
180  );
181 
182  volScalarField invT
183  (
184  alpha*rho*(1.0/(theta_.value()*this->delta()))*pow(flm_*fmm_, 1.0/8.0)
185  );
186 
187  volScalarField LM(L && M);
188 
189  fvScalarMatrix flmEqn
190  (
191  fvm::ddt(alpha, rho, flm_)
192  + fvm::div(alphaRhoPhi, flm_)
193  ==
194  invT*LM
195  - fvm::Sp(invT, flm_)
196  + fvOptions(alpha, rho, flm_)
197  );
198 
199  flmEqn.relax();
200  fvOptions.constrain(flmEqn);
201  flmEqn.solve();
202  fvOptions.correct(flm_);
203  bound(flm_, flm0_);
204 
205  volScalarField MM(M && M);
206 
207  fvScalarMatrix fmmEqn
208  (
209  fvm::ddt(alpha, rho, fmm_)
210  + fvm::div(alphaRhoPhi, fmm_)
211  ==
212  invT*MM
213  - fvm::Sp(invT, fmm_)
214  + fvOptions(alpha, rho, fmm_)
215  );
216 
217  fmmEqn.relax();
218  fvOptions.constrain(fmmEqn);
219  fmmEqn.solve();
220  fvOptions.correct(fmm_);
221  bound(fmm_, fmm0_);
222 
223  correctNut(gradU);
224 }
225 
226 
227 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
228 
229 } // End namespace LESModels
230 } // End namespace Foam
231 
232 // ************************************************************************* //
Eddy viscosity LES SGS model base class.
scalar delta
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
fv::options & fvOptions
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dynamicLagrangian(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.
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Finite-volume options.
Definition: fvOptions.H:52
phi
Definition: pEqn.H:104
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
BasicMomentumTransportModel::transportModel transportModel
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
autoPtr< surfaceVectorField > Uf
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
virtual bool read()
Read model coefficients if they have changed.
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:521
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
U
Definition: pEqn.H:72
BasicMomentumTransportModel::rhoField rhoField
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:101
virtual void correct()
Correct Eddy-Viscosity and related properties.
static autoPtr< LESfilter > New(const fvMesh &, const dictionary &, const word &filterDictName="filter")
Return a reference to the selected LES filter.
Definition: LESfilter.C:41
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
BasicMomentumTransportModel::alphaField alphaField
#define M(I)
Namespace for OpenFOAM.