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"
27 
28 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
29 
30 template<class BasicMomentumTransportModel>
32 {
33  if (printCoeffs_)
34  {
35  Info<< coeffDict_.dictName() << coeffDict_ << endl;
36  }
37 }
38 
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
42 template<class BasicMomentumTransportModel>
44 (
45  const word& type,
46  const alphaField& alpha,
47  const rhoField& rho,
48  const volVectorField& U,
49  const surfaceScalarField& alphaRhoPhi,
50  const surfaceScalarField& phi,
51  const transportModel& transport
52 )
53 :
54  BasicMomentumTransportModel
55  (
56  type,
57  alpha,
58  rho,
59  U,
60  alphaRhoPhi,
61  phi,
62  transport
63  ),
64 
65  LESDict_(this->subOrEmptyDict("LES")),
66  turbulence_(LESDict_.lookup("turbulence")),
67  printCoeffs_(LESDict_.lookupOrDefault<Switch>("printCoeffs", false)),
68  coeffDict_(LESDict_.optionalSubDict(type + "Coeffs")),
69 
70  kMin_
71  (
73  (
74  "kMin",
75  LESDict_,
77  small
78  )
79  ),
80 
81  epsilonMin_
82  (
84  (
85  "epsilonMin",
86  LESDict_,
87  kMin_.dimensions()/dimTime,
88  small
89  )
90  ),
91 
92  omegaMin_
93  (
95  (
96  "omegaMin",
97  LESDict_,
99  small
100  )
101  ),
102 
103  delta_
104  (
106  (
107  IOobject::groupName("delta", alphaRhoPhi.group()),
108  *this,
109  LESDict_
110  )
111  )
112 {
113  // Force the construction of the mesh deltaCoeffs which may be needed
114  // for the construction of the derived models and BCs
115  this->mesh_.deltaCoeffs();
116 }
117 
118 
119 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
120 
121 template<class BasicMomentumTransportModel>
124 (
125  const alphaField& alpha,
126  const rhoField& rho,
127  const volVectorField& U,
128  const surfaceScalarField& alphaRhoPhi,
129  const surfaceScalarField& phi,
130  const transportModel& transport
131 )
132 {
133  const IOdictionary modelDict
134  (
135  momentumTransportModel::readModelDict
136  (
137  U.db(),
138  alphaRhoPhi.group()
139  )
140  );
141 
142  const word modelType =
143  modelDict.subDict("LES").lookupBackwardsCompatible<word>
144  (
145  {"model", "LESModel"}
146  );
147 
148  Info<< "Selecting LES turbulence model " << modelType << endl;
149 
150  typename dictionaryConstructorTable::iterator cstrIter =
151  dictionaryConstructorTablePtr_->find(modelType);
152 
153  if (cstrIter == dictionaryConstructorTablePtr_->end())
154  {
156  << "Unknown LESModel type "
157  << modelType << nl << nl
158  << "Valid LESModel types:" << endl
159  << dictionaryConstructorTablePtr_->sortedToc()
160  << exit(FatalError);
161  }
162 
163  return autoPtr<LESModel>
164  (
165  cstrIter()(alpha, rho, U, alphaRhoPhi, phi, transport)
166  );
167 }
168 
169 
170 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
171 
172 template<class BasicMomentumTransportModel>
174 {
176  {
177  LESDict_ <<= this->subDict("LES");
178  LESDict_.lookup("turbulence") >> turbulence_;
179 
180  coeffDict_ <<= LESDict_.optionalSubDict(type() + "Coeffs");
181 
182  delta_().read(LESDict_);
183 
184  kMin_.readIfPresent(LESDict_);
185 
186  return true;
187  }
188  else
189  {
190  return false;
191  }
192 }
193 
194 
195 template<class BasicMomentumTransportModel>
197 {
198  delta_().correct();
200 }
201 
202 
203 // ************************************************************************* //
static word group(const word &name)
Return group (extension part of name)
Definition: IOobject.C:134
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
dimensionedSymmTensor sqr(const dimensionedVector &dv)
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:31
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:982
const dimensionSet dimTime
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
Definition: word.H:59
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-9/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
phi
Definition: correctPhi.H:3
LESModel(const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport)
Construct from components.
Definition: LESModel.C:44
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:196
BasicMomentumTransportModel::rhoField rhoField
Definition: LESModel.H:99
BasicMomentumTransportModel::transportModel transportModel
Definition: LESModel.H:100
U
Definition: pEqn.H:72
virtual bool read()
Read model coefficients if they have changed.
Definition: LESModel.C:173
static autoPtr< LESModel > New(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport)
Return a reference to the selected LES model.
Definition: LESModel.C:124
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:98
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
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:855