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-2023 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 "polyMeshMap.H"
35 #include "polyDistributionMap.H"
36 #include "writeFile.H"
37 #include "OFstream.H"
38 #include "OSspecific.H"
40 
41 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 namespace functionObjects
46 {
48 
50  (
54  );
55 }
56 }
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 void Foam::functionObjects::sampledSets::correct()
86 {
87  bool setsFound = dict_.found("sets");
88  if (setsFound)
89  {
90  searchEngine_.correct();
91 
92  PtrList<sampledSet> newList
93  (
94  dict_.lookup("sets"),
95  sampledSet::iNew(mesh_, searchEngine_)
96  );
97  transfer(newList);
98  combineSampledSets();
99  }
100 }
101 
102 
103 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
104 
106 (
107  const word& name,
108  const Time& t,
109  const dictionary& dict
110 )
111 :
113  PtrList<sampledSet>(),
114  outputPath_
115  (
116  mesh_.time().globalPath()
117  /writeFile::outputPrefix
118  /(mesh_.name() != polyMesh::defaultRegion ? mesh_.name() : word())
119  /name
120  ),
121  searchEngine_(mesh_),
122  interpolationScheme_(word::null),
123  formatter_(nullptr)
124 {
125  read(dict);
126 }
127 
128 
129 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
130 
132 {}
133 
134 
135 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
136 
138 {
139  dict_ = dict;
140 
141  bool setsFound = dict_.found("sets");
142  if (setsFound)
143  {
144  dict_.lookup("fields") >> fields_;
145 
146  dict.lookup("interpolationScheme") >> interpolationScheme_;
147 
148  const word writeType(dict.lookup("setFormat"));
149 
150  // Define the set formatter
151  formatter_ = setWriter::New(writeType, dict);
152 
153  PtrList<sampledSet> newList
154  (
155  dict_.lookup("sets"),
156  sampledSet::iNew(mesh_, searchEngine_)
157  );
158  transfer(newList);
159  combineSampledSets();
160 
161  if (this->size())
162  {
163  Info<< "Reading set description:" << nl;
164  forAll(*this, seti)
165  {
166  Info<< " " << operator[](seti).name() << nl;
167  }
168  Info<< endl;
169  }
170  }
171 
172  return true;
173 }
174 
175 
177 {
178  return fields_;
179 }
180 
181 
183 {
184  return true;
185 }
186 
187 
189 {
190  if (size())
191  {
192  if (Pstream::master())
193  {
194  if (debug)
195  {
196  Pout<< "Creating directory "
197  << outputPath_/mesh_.time().name() << nl << endl;
198  }
199 
200  mkDir(outputPath_/mesh_.time().name());
201  }
202 
203  // Create a list of names of fields that are actually available
205  forAll(fields_, fieldi)
206  {
207  #define FoundFieldType(Type, nullArg) \
208  || foundObject<VolField<Type>>(fields_[fieldi])
210  {
211  fieldNames.append(fields_[fieldi]);
212  }
213  else
214  {
215  cannotFindObject(fields_[fieldi]);
216  }
217  #undef FoundFieldType
218  }
219 
220  // Create table of cached interpolations, to prevent unnecessary work
221  // when interpolating fields over multiple surfaces
222  #define DeclareInterpolations(Type, nullArg) \
223  HashPtrTable<interpolation<Type>> interpolation##Type##s;
225  #undef DeclareInterpolations
226 
227  // Sample and write the sets
228  forAll(*this, seti)
229  {
230  #define GenerateFieldTypeValues(Type, nullArg) \
231  PtrList<Field<Type>> field##Type##Values = \
232  sampleType<Type>(seti, fieldNames, interpolation##Type##s);
234  #undef GenerateFieldTypeValues
235 
236  if (Pstream::parRun())
237  {
238  if (Pstream::master() && masterSets_[seti].size())
239  {
240  formatter_->write
241  (
242  outputPath_/mesh_.time().name(),
243  operator[](seti).name(),
244  masterSets_[seti],
245  fieldNames
246  #define FieldTypeValuesParameter(Type, nullArg) \
247  , field##Type##Values
250  );
251 
252  }
253  }
254  else
255  {
256  if (operator[](seti).size())
257  {
258  formatter_->write
259  (
260  outputPath_/mesh_.time().name(),
261  operator[](seti).name(),
262  operator[](seti),
263  fieldNames
264  #define FieldTypeValuesParameter(Type, nullArg) \
265  , field##Type##Values
268  );
269  }
270  }
271  }
272  }
273 
274  return true;
275 }
276 
277 
279 {
280  if (&mesh == &mesh_)
281  {
282  correct();
283  }
284 }
285 
286 
287 
289 (
290  const polyTopoChangeMap& map
291 )
292 {
293  if (&map.mesh() == &mesh_)
294  {
295  correct();
296  }
297 }
298 
299 
301 {
302  if (&map.mesh() == &mesh_)
303  {
304  correct();
305  }
306 }
307 
308 
310 (
311  const polyDistributionMap& map
312 )
313 {
314  if (&map.mesh() == &mesh_)
315  {
316  correct();
317  }
318 }
319 
320 
321 // ************************************************************************* //
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
void setSize(const label)
Reset size of List.
Definition: List.C:281
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:85
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
const T & operator[](const label) const
Return element const reference.
Definition: UPtrListI.H:105
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base-class for Time/database functionObjects.
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
sampledSets(const word &name, const Time &time, const dictionary &dict)
Construct from Time and dictionary.
Definition: sampledSets.C:106
virtual wordList fields() const
Return the list of fields required.
Definition: sampledSets.C:176
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
Definition: sampledSets.C:289
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
Definition: sampledSets.C:310
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: sampledSets.C:300
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
Definition: sampledSets.C:278
virtual ~sampledSets()
Destructor.
Definition: sampledSets.C:131
virtual bool execute()
Execute, currently does nothing.
Definition: sampledSets.C:182
virtual bool write()
Sample and write.
Definition: sampledSets.C:188
virtual bool read(const dictionary &)
Read the sampledSets.
Definition: sampledSets.C:137
Writes run time, CPU time and clock time and optionally the CPU and clock times per time step.
functionObject base class for writing single files
Definition: writeFile.H:56
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
const polyMesh & mesh() const
Return polyMesh.
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
const polyMesh & mesh() const
Return polyMesh.
Definition: polyMeshMap.H:75
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const polyMesh & mesh() const
Return polyMesh.
Class used for the read-construction of.
Definition: sampledSet.H:136
Holds list of sampling points which is filled at construction time. Various implementations of this b...
Definition: sampledSet.H:64
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:286
A class for handling words, derived from string.
Definition: word.H:62
static List< word > fieldNames
Definition: globalFoam.H:46
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
#define WarningInFunction
Report a warning using Foam::Warning.
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
Namespace for OpenFOAM.
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:290
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
messageStream Info
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
FOR_ALL_FIELD_TYPES(makeFieldSourceTypedef)
static const char nl
Definition: Ostream.H:266
#define FoundFieldType(Type, nullArg)
#define FieldTypeValuesParameter(Type, nullArg)
#define DeclareInterpolations(Type, nullArg)
#define GenerateFieldTypeValues(Type, nullArg)
dictionary dict