LESThermophysicalTransportModel.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) 2020-2022 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 
28 
29 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
30 
31 template<class BasicThermophysicalTransportModel>
33 <
34  BasicThermophysicalTransportModel
35 >::printCoeffs
36 (
37  const word& type)
38 {
39  if (printCoeffs_)
40  {
41  Info<< coeffDict_.dictName() << coeffDict_ << endl;
42  }
43 }
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 template<class BasicThermophysicalTransportModel>
50 <
51  BasicThermophysicalTransportModel
53 (
54  const word& type,
55  const momentumTransportModel& momentumTransport,
56  const thermoModel& thermo
57 )
58 :
59  BasicThermophysicalTransportModel(momentumTransport, thermo),
60  LESDict_(this->subOrEmptyDict("LES")),
61  printCoeffs_(LESDict_.lookupOrDefault<Switch>("printCoeffs", false)),
62  coeffDict_(LESDict_.optionalSubDict(type + "Coeffs"))
63 {}
64 
65 
66 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
67 
68 template<class BasicThermophysicalTransportModel>
70 <
72  <
73  BasicThermophysicalTransportModel
74  >
75 >
77 <
78  BasicThermophysicalTransportModel
79 >::New
80 (
81  const momentumTransportModel& momentumTransport,
82  const thermoModel& thermo
83 )
84 {
86  (
88  (
89  thermophysicalTransportModel::typeName,
90  momentumTransport.alphaRhoPhi().group()
91  ),
92  momentumTransport.time().constant(),
93  momentumTransport.mesh(),
96  false
97  );
98 
99  if (header.headerOk())
100  {
101  IOdictionary modelDict(header);
102 
103  const word modelType(modelDict.subDict("LES").lookup( "model"));
104 
105  Info<< "Selecting LES thermophysical transport model "
106  << modelType << endl;
107 
108  typename dictionaryConstructorTable::iterator cstrIter =
109  dictionaryConstructorTablePtr_->find(modelType);
110 
111  if (cstrIter == dictionaryConstructorTablePtr_->end())
112  {
114  << "Unknown LES thermophysical transport model "
115  << modelType << nl << nl
116  << "Available models:" << endl
117  << dictionaryConstructorTablePtr_->sortedToc()
118  << exit(FatalError);
119  }
120 
121  return autoPtr<LESThermophysicalTransportModel>
122  (
123  cstrIter()(momentumTransport, thermo)
124  );
125  }
126  else
127  {
128  typedef
129  turbulenceThermophysicalTransportModels::unityLewisEddyDiffusivity
130  <
131  LESThermophysicalTransportModel
132  <
133  BasicThermophysicalTransportModel
134  >
135  > LESunityLewisEddyDiffusivity;
136 
137  Info<< "Selecting default LES thermophysical transport model "
138  << LESunityLewisEddyDiffusivity::typeName << endl;
139 
140  return autoPtr<LESThermophysicalTransportModel>
141  (
142  new LESunityLewisEddyDiffusivity
143  (
144  LESunityLewisEddyDiffusivity::typeName,
145  momentumTransport,
146  thermo,
147  true
148  )
149  );
150  }
151 }
152 
153 
154 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
155 
156 template<class BasicThermophysicalTransportModel>
158 <
159  BasicThermophysicalTransportModel
161 {
163  {
164  LESDict_ <<= this->subDict("LES");
165 
166  coeffDict_ <<= LESDict_.optionalSubDict(type() + "Coeffs");
167 
168  return true;
169  }
170  else
171  {
172  return false;
173  }
174 }
175 
176 
177 template<class BasicThermophysicalTransportModel>
179 predict()
180 {
181  BasicThermophysicalTransportModel::predict();
182 }
183 
184 
185 template<class BasicThermophysicalTransportModel>
187 correct()
188 {
190 }
191 
192 
193 // ************************************************************************* //
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
static word group(const word &name)
Return group (extension part of name)
Definition: IOobject.C:134
static word groupName(Name name, const word &group)
Templated abstract base class for LES thermophysical transport models.
virtual void correct()
Solve the thermophysical transport model equations.
virtual void predict()
Predict the LES transport coefficients if possible.
BasicThermophysicalTransportModel::thermoModel thermoModel
BasicThermophysicalTransportModel::momentumTransportModel momentumTransportModel
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
static const word & constant()
Return constant name.
Definition: TimePaths.H:122
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Abstract base class for turbulence models (RAS, LES and laminar).
const surfaceScalarField & alphaRhoPhi() const
Access function to phase flux field.
Templated form of IOobject providing type information for file reading and header type checking.
Definition: IOobject.H:531
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
bool read(const char *, int32_t &)
Definition: int32IO.C:85
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
error FatalError
static const char nl
Definition: Ostream.H:266
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
fluidMulticomponentThermo & thermo
Definition: createFields.H:31