LESModel.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) 2013-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 "LESModel.H"
28 
29 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
30 
31 template<class BasicMomentumTransportModel>
33 {
34  if (printCoeffs_)
35  {
36  Info<< coeffDict_.dictName() << coeffDict_ << endl;
37  }
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
43 template<class BasicMomentumTransportModel>
45 (
46  const word& type,
47  const alphaField& alpha,
48  const rhoField& rho,
49  const volVectorField& U,
50  const surfaceScalarField& alphaRhoPhi,
51  const surfaceScalarField& phi,
52  const viscosity& viscosity
53 )
54 :
55  BasicMomentumTransportModel
56  (
57  type,
58  alpha,
59  rho,
60  U,
61  alphaRhoPhi,
62  phi,
63  viscosity
64  ),
65 
66  LESDict_(this->subOrEmptyDict("LES")),
67  turbulence_(LESDict_.lookup("turbulence")),
68  printCoeffs_(LESDict_.lookupOrDefault<Switch>("printCoeffs", false)),
69  coeffDict_(LESDict_.optionalSubDict(type + "Coeffs")),
70 
71  kMin_
72  (
74  (
75  "kMin",
76  LESDict_,
78  small
79  )
80  ),
81 
82  epsilonMin_
83  (
85  (
86  "epsilonMin",
87  LESDict_,
88  kMin_.dimensions()/dimTime,
89  small
90  )
91  ),
92 
93  omegaMin_
94  (
96  (
97  "omegaMin",
98  LESDict_,
100  small
101  )
102  ),
103 
104  viscosityModel_
105  (
106  coeffDict_.found("viscosityModel")
108  (
109  coeffDict_,
110  viscosity,
111  U
112  )
113  : autoPtr<laminarModels::generalisedNewtonianViscosityModel>
114  (
115  new laminarModels::generalisedNewtonianViscosityModels::Newtonian
116  (
117  coeffDict_,
118  viscosity,
119  U
120  )
121  )
122  ),
123 
124  delta_
125  (
126  LESdelta::New
127  (
128  IOobject::groupName("delta", alphaRhoPhi.group()),
129  *this,
130  LESDict_
131  )
132  )
133 {
134  // Force the construction of the mesh deltaCoeffs which may be needed
135  // for the construction of the derived models and BCs
136  this->mesh_.deltaCoeffs();
137 }
138 
139 
140 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
141 
142 template<class BasicMomentumTransportModel>
145 (
146  const alphaField& alpha,
147  const rhoField& rho,
148  const volVectorField& U,
149  const surfaceScalarField& alphaRhoPhi,
150  const surfaceScalarField& phi,
151  const viscosity& viscosity
152 )
153 {
154  const IOdictionary modelDict
155  (
156  momentumTransportModel::readModelDict
157  (
158  U.db(),
159  alphaRhoPhi.group()
160  )
161  );
162 
163  const word modelType =
164  modelDict.subDict("LES").lookupBackwardsCompatible<word>
165  (
166  {"model", "LESModel"}
167  );
168 
169  Info<< "Selecting LES turbulence model " << modelType << endl;
170 
171  typename dictionaryConstructorTable::iterator cstrIter =
172  dictionaryConstructorTablePtr_->find(modelType);
173 
174  if (cstrIter == dictionaryConstructorTablePtr_->end())
175  {
177  << "Unknown LESModel type "
178  << modelType << nl << nl
179  << "Valid LESModel types:" << endl
180  << dictionaryConstructorTablePtr_->sortedToc()
181  << exit(FatalError);
182  }
183 
184  return autoPtr<LESModel>
185  (
186  cstrIter()(alpha, rho, U, alphaRhoPhi, phi, viscosity)
187  );
188 }
189 
190 
191 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
192 
193 template<class BasicMomentumTransportModel>
195 {
197  {
198  LESDict_ <<= this->subDict("LES");
199  LESDict_.lookup("turbulence") >> turbulence_;
200 
201  coeffDict_ <<= LESDict_.optionalSubDict(type() + "Coeffs");
202 
203  delta_().read(LESDict_);
204 
205  kMin_.readIfPresent(LESDict_);
206 
207  return true;
208  }
209  else
210  {
211  return false;
212  }
213 }
214 
215 
216 template<class BasicMomentumTransportModel>
218 {
219  viscosityModel_->correct();
220  delta_().correct();
222 }
223 
224 
225 // ************************************************************************* //
const char *const group
Group name for atomic constants.
static word group(const word &name)
Return group (extension part of name)
Definition: IOobject.C:134
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
U
Definition: pEqn.H:72
error FatalError
LESModel(const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity)
Construct from components.
Definition: LESModel.C:45
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-10/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Abstract base class for LES deltas.
Definition: LESdelta.H:50
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
virtual void printCoeffs(const word &type)
Print model coefficients.
Definition: LESModel.C:32
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
const dimensionSet dimless
Generic dimensioned Type class.
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:1002
const dimensionSet dimTime
A class for handling words, derived from string.
Definition: word.H:59
Info<< "Reading field p\"<< endl;volScalarField p(IOobject("p", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);pressureReference pressureReference(p, simple.dict());mesh.schemes().setFluxRequired(p.name());Info<< "Reading field pa\"<< endl;volScalarField pa(IOobject("pa", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field Ua\"<< endl;volVectorField Ua(IOobject("Ua", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);# 65 "/home/ubuntu/OpenFOAM-10/applications/solvers/incompressible/adjointShapeOptimisationFoam/createFields.H" 2label paRefCell=0;scalar paRefValue=0.0;setRefCell(pa, simple.dict(), paRefCell, paRefValue);mesh.schemes().setFluxRequired(pa.name());autoPtr< viscosityModel > viscosity(viscosityModel::New(mesh))
phi
Definition: correctPhi.H:3
static autoPtr< LESModel > New(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity)
Return a reference to the selected LES model.
Definition: LESModel.C:145
static const char nl
Definition: Ostream.H:260
const dimensionSet dimVelocity
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: LESModel.C:217
BasicMomentumTransportModel::rhoField rhoField
Definition: LESModel.H:107
virtual bool read()
Read model coefficients if they have changed.
Definition: LESModel.C:194
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
messageStream Info
BasicMomentumTransportModel::alphaField alphaField
Definition: LESModel.H:106
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
ITstream & lookupBackwardsCompatible(const wordList &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream, trying a list of keywords.
Definition: dictionary.C:875