RASModel.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-2018 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 "RASModel.H"
27 
28 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
29 
30 template<class BasicTurbulenceModel>
32 {
33  if (printCoeffs_)
34  {
35  Info<< coeffDict_.dictName() << coeffDict_ << endl;
36  }
37 }
38 
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
42 template<class BasicTurbulenceModel>
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  const word& propertiesName
53 )
54 :
55  BasicTurbulenceModel
56  (
57  type,
58  alpha,
59  rho,
60  U,
61  alphaRhoPhi,
62  phi,
63  transport,
64  propertiesName
65  ),
66 
67  RASDict_(this->subOrEmptyDict("RAS")),
68  turbulence_(RASDict_.lookup("turbulence")),
69  printCoeffs_(RASDict_.lookupOrDefault<Switch>("printCoeffs", false)),
70  coeffDict_(RASDict_.optionalSubDict(type + "Coeffs")),
71 
72  kMin_
73  (
75  (
76  "kMin",
77  RASDict_,
79  small
80  )
81  ),
82 
83  epsilonMin_
84  (
86  (
87  "epsilonMin",
88  RASDict_,
89  kMin_.dimensions()/dimTime,
90  small
91  )
92  ),
93 
94  omegaMin_
95  (
97  (
98  "omegaMin",
99  RASDict_,
101  small
102  )
103  )
104 {
105  // Force the construction of the mesh deltaCoeffs which may be needed
106  // for the construction of the derived models and BCs
107  this->mesh_.deltaCoeffs();
108 }
109 
110 
111 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
112 
113 template<class BasicTurbulenceModel>
116 (
117  const alphaField& alpha,
118  const rhoField& rho,
119  const volVectorField& U,
120  const surfaceScalarField& alphaRhoPhi,
121  const surfaceScalarField& phi,
122  const transportModel& transport,
123  const word& propertiesName
124 )
125 {
126  // get model name, but do not register the dictionary
127  // otherwise it is registered in the database twice
128  const word modelType
129  (
131  (
132  IOobject
133  (
134  IOobject::groupName(propertiesName, alphaRhoPhi.group()),
135  U.time().constant(),
136  U.db(),
137  IOobject::MUST_READ_IF_MODIFIED,
138  IOobject::NO_WRITE,
139  false
140  )
141  ).subDict("RAS").lookup("RASModel")
142  );
143 
144  Info<< "Selecting RAS turbulence model " << modelType << endl;
145 
146  typename dictionaryConstructorTable::iterator cstrIter =
147  dictionaryConstructorTablePtr_->find(modelType);
148 
149  if (cstrIter == dictionaryConstructorTablePtr_->end())
150  {
152  << "Unknown RASModel type "
153  << modelType << nl << nl
154  << "Valid RASModel types:" << endl
155  << dictionaryConstructorTablePtr_->sortedToc()
156  << exit(FatalError);
157  }
158 
159  return autoPtr<RASModel>
160  (
161  cstrIter()(alpha, rho, U, alphaRhoPhi, phi, transport, propertiesName)
162  );
163 }
164 
165 
166 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
167 
168 template<class BasicTurbulenceModel>
170 {
172  {
173  RASDict_ <<= this->subDict("RAS");
174  RASDict_.lookup("turbulence") >> turbulence_;
175 
176  coeffDict_ <<= RASDict_.optionalSubDict(type() + "Coeffs");
177 
178  kMin_.readIfPresent(RASDict_);
179  epsilonMin_.readIfPresent(RASDict_);
180  omegaMin_.readIfPresent(RASDict_);
181 
182  return true;
183  }
184  else
185  {
186  return false;
187  }
188 }
189 
190 
191 template<class BasicTurbulenceModel>
193 {
195 }
196 
197 
198 // ************************************************************************* //
static word group(const word &name)
Return group (extension part of name)
Definition: IOobject.C:176
surfaceScalarField & phi
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:319
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
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
Generic dimensioned Type class.
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
virtual bool read()
Read model coefficients if they have changed.
Definition: RASModel.C:169
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:699
RASModel(const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName)
Construct from components.
Definition: RASModel.C:44
A class for handling words, derived from string.
Definition: word.H:59
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
virtual void printCoeffs(const word &type)
Print model coefficients.
Definition: RASModel.C:31
static const char nl
Definition: Ostream.H:265
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-7/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
U
Definition: pEqn.H:72
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
static autoPtr< RASModel > New(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName)
Return a reference to the selected RAS model.
Definition: RASModel.C:116
const Time & time() const
Return time.
Definition: IOobject.C:360
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:354
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: RASModel.C:192
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583
const dimensionSet dimVelocity