All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
sampledSets.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) 2011-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 "sampledSets.H"
27 #include "dictionary.H"
28 #include "Time.H"
29 #include "volFields.H"
30 #include "ListListOps.H"
31 #include "SortableList.H"
32 #include "volPointInterpolation.H"
33 #include "polyTopoChangeMap.H"
34 #include "writeFile.H"
35 #include "OFstream.H"
36 #include "OSspecific.H"
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43 namespace functionObjects
44 {
45  defineTypeNameAndDebug(sampledSets, 0);
46 
48  (
49  functionObject,
50  sampledSets,
51  dictionary
52  );
53 }
54 }
55 
56 bool Foam::functionObjects::sampledSets::verbose_ = false;
57 
58 
59 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
60 
61 void Foam::functionObjects::sampledSets::combineSampledSets()
62 {
63  masterSets_.setSize(size());
64  masterSetOrders_.setSize(size());
65 
66  forAll(*this, seti)
67  {
68  const sampledSet& s = operator[](seti);
69 
70  Tuple2<coordSet, labelList> g = s.gather();
71 
72  masterSets_.set(seti, new coordSet(g.first()));
73  masterSetOrders_[seti] = g.second();
74 
75  if (Pstream::master() && masterSets_[seti].size() == 0)
76  {
78  << "Sample set " << s.name()
79  << " has zero points." << endl;
80  }
81  }
82 }
83 
84 
85 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
86 
88 (
89  const word& name,
90  const Time& t,
91  const dictionary& dict
92 )
93 :
94  fvMeshFunctionObject(name, t, dict),
96  outputPath_(fileName::null),
97  searchEngine_(mesh_),
98  interpolationScheme_(word::null),
99  formatter_(nullptr)
100 {
101  outputPath_ =
102  mesh_.time().globalPath()/functionObjects::writeFile::outputPrefix/name;
103 
104  if (mesh_.name() != fvMesh::defaultRegion)
105  {
106  outputPath_ = outputPath_/mesh_.name();
107  }
108 
109  read(dict);
110 }
111 
112 
113 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
114 
116 {}
117 
118 
119 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
120 
122 {
123  verbose_ = verbosity;
124 }
125 
126 
128 {
129  return fields_;
130 }
131 
132 
134 {
135  return true;
136 }
137 
138 
140 {
141  if (size())
142  {
143  if (Pstream::master())
144  {
145  if (debug)
146  {
147  Pout<< "Creating directory "
148  << outputPath_/mesh_.time().timeName() << nl << endl;
149  }
150 
151  mkDir(outputPath_/mesh_.time().timeName());
152  }
153 
154  // Create a list of names of fields that are actually available
156  forAll(fields_, fieldi)
157  {
158  #define FoundFieldType(Type, nullArg) \
159  || foundObject<VolField<Type>>(fields_[fieldi])
161  {
162  fieldNames.append(fields_[fieldi]);
163  }
164  else
165  {
166  cannotFindObject(fields_[fieldi]);
167  }
168  #undef FoundFieldType
169  }
170 
171  // Create table of cached interpolations, to prevent unnecessary work
172  // when interpolating fields over multiple surfaces
173  #define DeclareInterpolations(Type, nullArg) \
174  HashPtrTable<interpolation<Type>> interpolation##Type##s;
176  #undef DeclareInterpolations
177 
178  // Sample and write the sets
179  forAll(*this, seti)
180  {
181  #define GenerateFieldTypeValues(Type, nullArg) \
182  PtrList<Field<Type>> field##Type##Values = \
183  sampleType<Type>(seti, fieldNames, interpolation##Type##s);
185  #undef GenerateFieldTypeValues
186 
187  if (Pstream::parRun())
188  {
189  if (Pstream::master() && masterSets_[seti].size())
190  {
191  formatter_->write
192  (
193  outputPath_/mesh_.time().timeName(),
194  operator[](seti).name(),
195  masterSets_[seti],
196  fieldNames
197  #define FieldTypeValuesParameter(Type, nullArg) \
198  , field##Type##Values
200  #undef FieldTypeValuesParameter
201  );
202 
203  }
204  }
205  else
206  {
207  if (operator[](seti).size())
208  {
209  formatter_->write
210  (
211  outputPath_/mesh_.time().timeName(),
212  operator[](seti).name(),
213  operator[](seti),
214  fieldNames
215  #define FieldTypeValuesParameter(Type, nullArg) \
216  , field##Type##Values
218  #undef FieldTypeValuesParameter
219  );
220  }
221  }
222  }
223  }
224 
225  return true;
226 }
227 
228 
230 {
231  dict_ = dict;
232 
233  bool setsFound = dict_.found("sets");
234  if (setsFound)
235  {
236  dict_.lookup("fields") >> fields_;
237 
238  dict.lookup("interpolationScheme") >> interpolationScheme_;
239 
240  const word writeType(dict.lookup("setFormat"));
241 
242  // Define the set formatter
243  formatter_ = setWriter::New(writeType, dict);
244 
245  PtrList<sampledSet> newList
246  (
247  dict_.lookup("sets"),
248  sampledSet::iNew(mesh_, searchEngine_)
249  );
250  transfer(newList);
251  combineSampledSets();
252 
253  if (this->size())
254  {
255  Info<< "Reading set description:" << nl;
256  forAll(*this, seti)
257  {
258  Info<< " " << operator[](seti).name() << nl;
259  }
260  Info<< endl;
261  }
262  }
263 
264  return true;
265 }
266 
267 
269 {
270  bool setsFound = dict_.found("sets");
271  if (setsFound)
272  {
273  searchEngine_.correct();
274 
275  PtrList<sampledSet> newList
276  (
277  dict_.lookup("sets"),
278  sampledSet::iNew(mesh_, searchEngine_)
279  );
280  transfer(newList);
281  combineSampledSets();
282  }
283 }
284 
285 
287 {
288  if (&mesh == &mesh_)
289  {
290  correct();
291  }
292 }
293 
294 
295 
297 (
298  const polyTopoChangeMap& map
299 )
300 {
301  if (&map.mesh() == &mesh_)
302  {
303  correct();
304  }
305 }
306 
307 
309 {
310  correct();
311 }
312 
313 
315 (
316  const polyMesh::readUpdateState state
317 )
318 {
319  if (state != polyMesh::UNCHANGED)
320  {
321  correct();
322  }
323 }
324 
325 
326 // ************************************************************************* //
static autoPtr< setWriter > New(const word &writeType, const IOstream::streamFormat writeFormat=IOstream::ASCII, const IOstream::compressionType writeCompression=IOstream::UNCOMPRESSED)
Select given write options.
Definition: setWriter.C:279
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define FieldTypeValuesParameter(Type, nullArg)
#define DeclareInterpolations(Type, nullArg)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-10/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
#define FoundFieldType(Type, nullArg)
sampledSets(const word &name, const Time &time, const dictionary &dict)
Construct from Time and dictionary.
Definition: sampledSets.C:88
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
static const fileName null
An empty fileName.
Definition: fileName.H:97
virtual bool execute()
Execute, currently does nothing.
Definition: sampledSets.C:133
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:325
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define GenerateFieldTypeValues(Type, nullArg)
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
Definition: sampledSets.C:297
fvMesh & mesh
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
Macros for easy insertion into run-time selection tables.
static List< word > fieldNames
Definition: globalFoam.H:46
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
bool read(const char *, int32_t &)
Definition: int32IO.C:85
static const word outputPrefix
Directory prefix.
Definition: writeFile.H:72
virtual bool write()
Sample and write.
Definition: sampledSets.C:139
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.
virtual void readUpdate(const polyMesh::readUpdateState state)
Update topology using the given map due to readUpdate.
Definition: sampledSets.C:315
static const word null
An empty word.
Definition: word.H:77
virtual ~sampledSets()
Destructor.
Definition: sampledSets.C:115
const polyMesh & mesh() const
Return polyMesh.
virtual bool read(const dictionary &)
Read the sampledSets.
Definition: sampledSets.C:229
Class used for the read-construction of.
Definition: sampledSet.H:135
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
Definition: sampledSets.C:286
static const char nl
Definition: Ostream.H:260
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:290
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
defineTypeNameAndDebug(Qdot, 0)
void correct()
Correct for mesh changes.
Definition: sampledSets.C:268
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
#define WarningInFunction
Report a warning using Foam::Warning.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
messageStream Info
FOR_ALL_FIELD_TYPES(DefineFvWallInfoType)
void verbose(const bool verbosity=true)
Set verbosity level.
Definition: sampledSets.C:121
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:90
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:50
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: sampledSets.C:308
virtual wordList fields() const
Return the list of fields required.
Definition: sampledSets.C:127
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