fvModels.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) 2021-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 
26 #include "fvModels.H"
27 #include "fvConstraint.H"
28 #include "fvMesh.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
35 }
36 
37 
38 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
39 
40 Foam::IOobject Foam::fvModels::createIOobject
41 (
42  const fvMesh& mesh
43 ) const
44 {
45  typeIOobject<IOdictionary> io
46  (
47  typeName,
48  mesh.time().constant(),
49  mesh,
52  );
53 
54  if (io.headerOk())
55  {
56  Info<< "Creating fvModels from "
57  << io.instance()/io.name() << nl
58  << endl;
59 
60  io.readOpt() = IOobject::MUST_READ_IF_MODIFIED;
61  return io;
62  }
63  else
64  {
65  // For backward-compatibility
66  // check if the fvOptions file is in constant
67  io.rename("fvOptions");
68 
69  if (io.headerOk())
70  {
71  Warning
72  << "Creating fvModels from "
73  << io.instance()/io.name() << nl
74  << endl;
75 
76  io.readOpt() = IOobject::MUST_READ_IF_MODIFIED;
77  return io;
78  }
79  else
80  {
81  // For backward-compatibility
82  // check if the fvOptions file is in system
83  io.instance() = mesh.time().system();
84 
85  if (io.headerOk())
86  {
87  Warning
88  << "Creating fvModels from "
89  << io.instance()/io.name()
90  << " rather than constant/fvModels"
91  << endl;
92 
93  io.readOpt() = IOobject::MUST_READ_IF_MODIFIED;
94  return io;
95  }
96  else
97  {
98  io.rename(typeName);
99  io.readOpt() = IOobject::NO_READ;
100  return io;
101  }
102  }
103  }
104 }
105 
106 
107 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
108 
109 void Foam::fvModels::checkApplied() const
110 {
111  if (mesh().time().timeIndex() > checkTimeIndex_)
112  {
113  const PtrListDictionary<fvModel>& modelList(*this);
114 
115  forAll(modelList, i)
116  {
117  const fvModel& model = modelList[i];
118 
119  wordHashSet notAddSupFields(model.addSupFields());
120  notAddSupFields -= addSupFields_[i];
121 
122  forAllConstIter(wordHashSet, notAddSupFields, iter)
123  {
125  << "Model " << model.name()
126  << " defined for field " << iter.key()
127  << " but never used" << endl;
128  }
129  }
130 
131  checkTimeIndex_ = mesh().time().timeIndex();
132  }
133 }
134 
135 
136 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
137 
139 (
140  const fvMesh& mesh
141 )
142 :
144  (
145  mesh,
146  createIOobject(mesh)
147  ),
149  checkTimeIndex_(mesh.time().timeIndex() + 1),
150  addSupFields_()
151 {
152  readHeaderOk(IOstream::ASCII, typeName);
153 
154  const bool readFromFvModels(IOobject::name() == typeName);
155 
156  const dictionary& dict(*this);
157 
158  // Count number of active fvModels
159  label count = 0;
160  forAllConstIter(dictionary, dict, iter)
161  {
162  if (iter().isDict())
163  {
164  count++;
165  }
166  }
167 
169 
170  addSupFields_.setSize(count);
171 
172  label i = 0;
173  forAllConstIter(dictionary, dict, iter)
174  {
175  if (iter().isDict())
176  {
177  const word& name = iter().keyword();
178  const dictionary& modelDict = iter().dict();
179 
180  const word modelType(modelDict.lookup("type"));
181 
182  if
183  (
184  readFromFvModels
185  || !fvConstraint::dictionaryConstructorTablePtr_->found
186  (
187  modelType
188  )
189  )
190  {
192  (
193  i,
194  name,
195  fvModel::New(name, modelDict, mesh).ptr()
196  );
197 
198  addSupFields_.set(i, new wordHashSet());
199 
200  i++;
201  }
202  }
203  }
204 
206  addSupFields_.setSize(i);
207 
208  if (readFromFvModels)
209  {
210  // Add file watch on the fvModels dictionary for
211  // MUST_READ_IF_MODIFIED
212  addWatch();
213  }
214  else
215  {
216  // If the fvOptions file was read re-name to fvModels
217  rename(typeName);
218  }
219 }
220 
221 
222 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
223 
224 bool Foam::fvModels::addsSupToField(const word& fieldName) const
225 {
226  const PtrListDictionary<fvModel>& modelList(*this);
227 
228  forAll(modelList, i)
229  {
230  if (modelList[i].addsSupToField(fieldName))
231  {
232  return true;
233  }
234  }
235 
236  return false;
237 }
238 
239 
240 Foam::scalar Foam::fvModels::maxDeltaT() const
241 {
242  const PtrListDictionary<fvModel>& modelList(*this);
243 
244  scalar maxDeltaT = great;
245 
246  forAll(modelList, i)
247  {
248  maxDeltaT = min(maxDeltaT, modelList[i].maxDeltaT());
249  }
250 
251  return maxDeltaT;
252 }
253 
254 
256 {
257  PtrListDictionary<fvModel>& modelList(*this);
258 
259  forAll(modelList, i)
260  {
261  modelList[i].preUpdateMesh();
262  }
263 }
264 
265 
267 {
268  bool allOk = true;
269 
270  PtrListDictionary<fvModel>& modelList(*this);
271 
272  forAll(modelList, i)
273  {
274  allOk = allOk && modelList[i].movePoints();
275  }
276 
277  return allOk;
278 }
279 
280 
282 {
283  PtrListDictionary<fvModel>& modelList(*this);
284 
285  forAll(modelList, i)
286  {
287  modelList[i].topoChange(map);
288  }
289 }
290 
291 
293 {
294  PtrListDictionary<fvModel>& modelList(*this);
295 
296  forAll(modelList, i)
297  {
298  modelList[i].mapMesh(map);
299  }
300 }
301 
302 
304 {
305  PtrListDictionary<fvModel>& modelList(*this);
306 
307  forAll(modelList, i)
308  {
309  modelList[i].distribute(map);
310  }
311 }
312 
313 
315 {
316  is >> *this;
317  return !is.bad();
318 }
319 
320 
322 {
323  dictionary::write(os, false);
324  return os.good();
325 }
326 
327 
329 {
331  {
332  checkTimeIndex_ = mesh().time().timeIndex() + 1;
333 
334  bool allOk = true;
335 
336  PtrListDictionary<fvModel>& modelList(*this);
337 
338  forAll(modelList, i)
339  {
340  fvModel& model = modelList[i];
341  bool ok = model.read(subDict(model.name()));
342  allOk = (allOk && ok);
343  }
344  return allOk;
345  }
346  else
347  {
348  return false;
349  }
350 }
351 
352 
354 {
355  PtrListDictionary<fvModel>& modelList(*this);
356 
357  forAll(modelList, i)
358  {
359  modelList[i].correct();
360  }
361 }
362 
363 
364 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
365 
366 Foam::Ostream& Foam::operator<<
367 (
368  Ostream& os,
369  const fvModels& models
370 )
371 {
372  models.writeData(os);
373  return os;
374 }
375 
376 
377 // ************************************************************************* //
void write(Ostream &, const bool subDict=true) const
Write dictionary, normally with sub-dictionary formatting.
Definition: dictionaryIO.C:207
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
const word & name() const
Return name.
Definition: IOobject.H:315
virtual scalar maxDeltaT() const
Return the maximum time-step for stable operation.
Definition: fvModels.C:240
bool bad() const
Return true if stream is corrupted.
Definition: IOstream.H:348
const word & name() const
Return const access to the source name.
Definition: fvModelI.H:28
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
autoPtr< T > set(const label, const word &key, T *)
Set element to pointer provided and return old element.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: fvModel.C:164
virtual void correct()
Correct the fvModels.
Definition: fvModels.C:353
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
scalar maxDeltaT
label count(const ListType &l, typename ListType::const_reference x)
Count the number of occurrences of a value in a list.
virtual void preUpdateMesh()
Prepare for mesh update.
Definition: fvModels.C:255
Finite volume model abstract base class.
Definition: fvModel.H:57
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
fvMesh & mesh
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:330
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Definition: MeshObject.H:87
A class for handling words, derived from string.
Definition: word.H:59
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
static autoPtr< fvModel > New(const word &name, const dictionary &dict, const fvMesh &mesh)
Return a reference to the selected fvModel.
Definition: fvModel.C:94
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
Definition: fvModels.C:303
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
Definition: fvModels.C:281
Finite volume models.
Definition: fvModels.H:60
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
HashSet wordHashSet
A HashSet with word keys.
Definition: HashSet.H:208
virtual bool readData(Istream &)
ReadData function which reads the fvModels dictionary.
Definition: fvModels.C:314
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
messageStream Warning
virtual bool writeData(Ostream &os) const
Write data to Ostream.
Definition: fvModels.C:321
Foam::fvModels & fvModels
virtual bool movePoints()
Update for mesh motion.
Definition: fvModels.C:266
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:28
#define WarningInFunction
Report a warning using Foam::Warning.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
label timeIndex
Definition: getTimeIndex.H:4
messageStream Info
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: fvModels.C:292
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:50
virtual bool addsSupToField(const word &fieldName) const
Return true if an fvModel adds a source term to the given.
Definition: fvModels.C:224
fvModels(const fvMesh &mesh)
Construct from components with list of field names.
Definition: fvModels.C:139
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
bool found
virtual bool read()
Read the fvModels dictionary if it has changed.
Definition: fvModels.C:328
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864