solverTemplate.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015 OpenFOAM Foundation
9 Copyright (C) 2015-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "solverTemplate.H"
30#include "Time.H"
31#include "IOPtrList.H"
32#include "polyMesh.H"
33#include "regionProperties.H"
34
35// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36
37const Foam::Enum
38<
40>
42({
43 { solverType::stCompressible, "compressible" },
44 { solverType::stIncompressible, "incompressible" },
45 { solverType::stBuoyant, "buoyant" },
46 { solverType::stUnknown, "unknown" },
47});
48
49
50// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51
53(
54 IOobject& dictHeader,
55 const word& entryName
56) const
57{
58 if (!dictHeader.typeHeaderOk<IOdictionary>(true))
59 {
61 << "Unable to open file "
62 << dictHeader.objectPath()
63 << exit(FatalError);
64 }
65
66 IOdictionary dict(dictHeader);
67 return dict.get<word>(entryName);
68}
69
70
71Foam::dictionary Foam::solverTemplate::readFluidFieldTemplates
72(
73 const word& regionName,
74 const fileName& baseDir,
75 const dictionary& solverDict,
76 const Time& runTime
77) const
78{
79 Info<< " Reading fluid field templates";
80
81 if (!regionName.empty())
82 {
83 Info<< " for region " << regionName;
84 }
85 Info<< endl;
86
87 dictionary fieldTemplates = solverDict.subDict("fluidFields");
88
89 const fileName turbModelDir(baseDir/"models"/"turbulence");
90
91 const dictionary fieldModels(solverDict.subDict("fluidModels"));
92
93 word turbModel("laminar"); // default to laminar
94 word turbType("none");
95
96 if (fieldModels.readIfPresent("turbulenceModel", turbType))
97 {
98 if (turbType == "turbulenceModel")
99 {
100 IOdictionary turbPropDict
101 (
102 IOobject
103 (
104 // turbulenceModel::propertiesName
105 "turbulenceProperties",
108 runTime,
111 false
112 )
113 );
114
115 const word modelType(turbPropDict.get<word>("simulationType"));
116
117 if (modelType == "laminar")
118 {
119 // Leave as laminar
120 }
121 else if (modelType == "RAS")
122 {
123 // "RASModel" for v2006 and earlier
124 turbPropDict.subDict(modelType)
125 .readCompat("model", {{"RASModel", -2006}}, turbModel);
126 }
127 else if (modelType == "LES")
128 {
129 // "LESModel" for v2006 and earlier
130 turbPropDict.subDict(modelType)
131 .readCompat("model", {{"LESModel", -2006}}, turbModel);
132 }
133 else
134 {
136 << "Unhandled turbulence model option " << modelType
137 << ". Valid options are laminar, RAS, LES"
138 << exit(FatalError);
139 }
140 }
141 else
142 {
144 << "Unhandled turbulence model option " << turbType
145 << ". Valid options are turbulenceModel"
146 << exit(FatalError);
147 }
148 }
149
150 Info<< " Selecting " << turbType << ": " << turbModel << endl;
151
152 IOdictionary turbModelDict
153 (
154 IOobject
155 (
156 fileName(turbModelDir/turbModel),
157 runTime,
159 )
160 );
161
162 // Merge common fluid fields
163 fieldTemplates.merge(turbModelDict.subDict("fluidFields"));
164
165 // Merge specific compressible or incompressible fluid fields
166 switch (solverType_)
167 {
168 case stIncompressible:
169 {
170 fieldTemplates.merge(turbModelDict.subDict("incompressibleFields"));
171 break;
172 }
173 case stCompressible:
174 case stBuoyant:
175 {
176 fieldTemplates.merge(turbModelDict.subDict("compressibleFields"));
177 break;
178 }
179 default:
180 {
181 // do nothing
182 }
183 }
184
185 return fieldTemplates;
186}
187
188
189Foam::dictionary Foam::solverTemplate::readSolidFieldTemplates
190(
191 const word& regionName,
192 const dictionary& solverDict
193) const
194{
195 Info<< " Reading solid field templates for region " << regionName
196 << endl;
197
198 // nothing more to do for solids
199 return solverDict.subDict("solidFields");
200}
201
202
203void Foam::solverTemplate::setRegionProperties
204(
205 const dictionary& fieldDict,
206 const word& regionType,
207 const word& regionName,
208 const label regionI
209)
210{
211 regionTypes_[regionI] = regionType,
212 regionNames_[regionI] = regionName,
213 fieldNames_[regionI] = fieldDict.toc();
214 fieldTypes_[regionI].setSize(fieldNames_[regionI].size());
215 fieldDimensions_[regionI].setSize(fieldNames_[regionI].size());
216
217 forAll(fieldNames_[regionI], i)
218 {
219 const word& fieldName = fieldNames_[regionI][i];
220 const dictionary& dict = fieldDict.subDict(fieldName);
221
222 dict.readEntry("type", fieldTypes_[regionI][i]);
223 fieldDimensions_[regionI].set
224 (
225 i,
226 new dimensionSet(dict, "dimensions")
227 );
228 }
229}
230
231
232// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
233
235(
236 const fileName& baseDir,
237 const Time& runTime,
238 const word& solverName
239)
240:
241 solverType_(stUnknown),
242 multiRegion_(false),
243 regionTypes_(),
244 fieldNames_(),
245 fieldTypes_(),
246 fieldDimensions_()
247{
248 IOdictionary solverDict
249 (
250 IOobject
251 (
252 fileName(baseDir/"solvers"/solverName),
253 runTime,
254 IOobject::MUST_READ
255 )
256 );
257
258 Info<< "Selecting " << solverName << ": ";
259
260 solverType_ = solverTypeNames_.get("solverType", solverDict);
261 multiRegion_ = solverDict.get<bool>("multiRegion");
262
263 Info<< solverTypeNames_[solverType_];
264 if (multiRegion_)
265 {
266 Info<< ", multi-region";
267 }
268 else
269 {
270 Info<< ", single-region";
271 }
272
273 Info<< " case" << endl;
274
275 if (multiRegion_)
276 {
277 // read regionProperties
278 regionProperties rp(runTime);
279
280 const wordList fluidNames(rp["fluid"]);
281 const wordList solidNames(rp["solid"]);
282
283 const label nRegion = fluidNames.size() + solidNames.size();
284
285 regionTypes_.setSize(nRegion);
286 regionNames_.setSize(nRegion);
287 fieldNames_.setSize(nRegion);
288 fieldTypes_.setSize(nRegion);
289 fieldDimensions_.setSize(nRegion);
290
291 // read templates for solver lists of available
292 // - fields and dimensions required for the solver
293 // - models
294 label regionI = 0;
296 {
297 const dictionary fieldDict
298 (
299 readFluidFieldTemplates
300 (
301 fluidNames[i],
302 baseDir,
303 solverDict,
304 runTime
305 )
306 );
307
308 setRegionProperties(fieldDict, "fluid", fluidNames[i], regionI++);
309 }
310
312 {
313 const dictionary fieldDict
314 (
315 readSolidFieldTemplates
316 (
317 solidNames[i],
318 solverDict
319 )
320 );
321 setRegionProperties(fieldDict, "solid", solidNames[i], regionI++);
322 }
323 }
324 else
325 {
326 regionTypes_.setSize(1);
327 regionNames_.setSize(1);
328 fieldNames_.setSize(1);
329 fieldTypes_.setSize(1);
330 fieldDimensions_.setSize(1);
331
332 // read templates for solver lists of available
333 // - fields and dimensions required for the solver
334 // - models
335 const dictionary fieldDict
336 (
337 readFluidFieldTemplates
338 (
339 word::null, // use region = null for single-region cases
340 baseDir,
341 solverDict,
342 runTime
343 )
344 );
345
346 setRegionProperties(fieldDict, "fluid", word::null, 0);
347 }
348}
349
350
351// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
352
354{
355 return solverTypeNames_[solverType_];
356}
357
358
360{
361 return multiRegion_;
362}
363
364
365Foam::label Foam::solverTemplate::nRegion() const
366{
367 return regionTypes_.size();
368}
369
370
371const Foam::word& Foam::solverTemplate::regionType(const label regionI) const
372{
373 return regionTypes_[regionI];
374}
375
376
377const Foam::word& Foam::solverTemplate::regionName(const label regionI) const
378{
379 return regionNames_[regionI];
380}
381
382
384(
385 const label regionI
386) const
387{
388 return fieldNames_[regionI];
389}
390
391
393(
394 const label regionI
395) const
396{
397 return fieldTypes_[regionI];
398}
399
400
402(
403 const label regionI
404) const
405{
406 return fieldDimensions_[regionI];
407}
408
409
410// ************************************************************************* //
const wordList fluidNames(rp["fluid"])
void readFromDict()
Read old info from dict.
Definition: BFGS.C:144
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
const word & constant() const
Return constant name.
Definition: TimePathsI.H:96
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Definition: dictionary.C:460
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
virtual const wordRes & fieldNames() const noexcept
Return names of fields to probe.
Definition: probes.H:335
Class to store solver template specifications.
word type() const
Solver type name.
label nRegion() const
Return the number of regions.
const wordList & fieldTypes(const label regionI) const
Return the field types.
bool multiRegion() const
Return the multi-region flag.
const PtrList< dimensionSet > & fieldDimensions(const label regionI) const
Return the field dimensions.
static const Enum< solverType > solverTypeNames_
Solver type names.
solverType
Solver type.
A class for handling words, derived from Foam::string.
Definition: word.H:68
engineTime & runTime
regionProperties rp(runTime)
Foam::word regionName(Foam::polyMesh::defaultRegion)
const wordList solidNames(rp["solid"])
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
List< word > wordList
A List of words.
Definition: fileName.H:63
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333