pointZoneSet.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-2024 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 "pointZoneSet.H"
27 #include "polyTopoChangeMap.H"
28 #include "polyMesh.H"
29 #include "processorPolyPatch.H"
30 #include "cyclicPolyPatch.H"
31 
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
42 
46 
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
51 {
53  sortedOrder(addressing_, order);
54  inplaceReorder(order, addressing_);
55 
57  pointSet::resize(2*addressing_.size());
58  forAll(addressing_, i)
59  {
60  pointSet::insert(addressing_[i]);
61  }
62 }
63 
64 
65 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
66 
68 (
69  const polyMesh& mesh,
70  const word& name,
71  readOption r,
72  writeOption w
73 )
74 :
75  pointSet(mesh, name, 1000), // do not read pointSet
76  mesh_(mesh),
77  addressing_(0)
78 {
79  const pointZoneList& pointZones = mesh.pointZones();
80  label zoneID = pointZones.findIndex(name);
81 
82  if
83  (
86  || (r == IOobject::READ_IF_PRESENT && zoneID != -1)
87  )
88  {
89  const pointZone& fz = pointZones[zoneID];
90  addressing_ = fz;
91  }
92 
93  updateSet();
94 
95  check(mesh.nPoints());
96 }
97 
98 
100 (
101  const polyMesh& mesh,
102  const word& name,
103  const label size,
104  writeOption w
105 )
106 :
107  pointSet(mesh, name, size, w),
108  mesh_(mesh),
109  addressing_(0)
110 {
111  updateSet();
112 }
113 
114 
116 (
117  const polyMesh& mesh,
118  const word& name,
119  const topoSet& set,
120  writeOption w
121 )
122 :
123  pointSet(mesh, name, set.size(), w),
124  mesh_(mesh),
125  addressing_(refCast<const pointZoneSet>(set).addressing())
126 {
127  updateSet();
128 }
129 
130 
131 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
132 
134 {}
135 
136 
137 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
138 
139 void pointZoneSet::invert(const label maxLen)
140 {
141  // Count
142  label n = 0;
143 
144  for (label pointi = 0; pointi < maxLen; pointi++)
145  {
146  if (!found(pointi))
147  {
148  n++;
149  }
150  }
151 
152  // Fill
153  addressing_.setSize(n);
154  n = 0;
155 
156  for (label pointi = 0; pointi < maxLen; pointi++)
157  {
158  if (!found(pointi))
159  {
160  addressing_[n] = pointi;
161  n++;
162  }
163  }
164  updateSet();
165 }
166 
167 
169 {
170  DynamicList<label> newAddressing(addressing_.size());
171 
172  const pointZoneSet& fSet = refCast<const pointZoneSet>(set);
173 
174  forAll(fSet.addressing(), i)
175  {
176  label pointi = fSet.addressing()[i];
177 
178  if (found(pointi))
179  {
180  newAddressing.append(pointi);
181  }
182  }
183 
184  addressing_.transfer(newAddressing);
185  updateSet();
186 }
187 
188 
190 {
191  DynamicList<label> newAddressing(addressing_);
192 
193  const pointZoneSet& fSet = refCast<const pointZoneSet>(set);
194 
195  forAll(fSet.addressing(), i)
196  {
197  label pointi = fSet.addressing()[i];
198 
199  if (!found(pointi))
200  {
201  newAddressing.append(pointi);
202  }
203  }
204 
205  addressing_.transfer(newAddressing);
206  updateSet();
207 }
208 
209 
211 {
212  DynamicList<label> newAddressing(addressing_.size());
213 
214  const pointZoneSet& fSet = refCast<const pointZoneSet>(set);
215 
216  forAll(addressing_, i)
217  {
218  label pointi = addressing_[i];
219 
220  if (!fSet.found(pointi))
221  {
222  // Not found in fSet so add
223  newAddressing.append(pointi);
224  }
225  }
226 
227  addressing_.transfer(newAddressing);
228  updateSet();
229 }
230 
231 
232 void pointZoneSet::sync(const polyMesh& mesh)
233 {
234  pointSet::sync(mesh);
235 
236  // Take over contents of pointSet into addressing.
237  addressing_ = sortedToc();
238  updateSet();
239 }
240 
241 
243 {
244  return mesh.nPoints();
245 }
246 
247 
249 (
253  const bool write
254 ) const
255 {
256  // Write shadow pointSet
257  word oldTypeName = typeName;
258  const_cast<word&>(type()) = pointSet::typeName;
259  bool ok = pointSet::writeObject(s, v, c, write);
260  const_cast<word&>(type()) = oldTypeName;
261 
262  // Modify pointZone
263  pointZoneList& pointZones = const_cast<polyMesh&>(mesh_).pointZones();
264  label zoneID = pointZones.findIndex(name());
265 
266  if (zoneID == -1)
267  {
268  zoneID = pointZones.size();
269 
270  pointZones.setSize(zoneID+1);
271  pointZones.set
272  (
273  zoneID,
274  new pointZone
275  (
276  name(),
277  addressing_,
278  pointZones
279  )
280  );
281  }
282  else
283  {
284  pointZones[zoneID] = addressing_;
285  }
286 
287  return ok && pointZones.write(write);
288 }
289 
290 
292 {
293  // pointZone
294  labelList newAddressing(addressing_.size());
295 
296  label n = 0;
297  forAll(addressing_, i)
298  {
299  label pointi = addressing_[i];
300  label newPointi = map.reversePointMap()[pointi];
301  if (newPointi >= 0)
302  {
303  newAddressing[n] = newPointi;
304  n++;
305  }
306  }
307  newAddressing.setSize(n);
308 
309  addressing_.transfer(newAddressing);
310 
311  updateSet();
312 }
313 
314 
316 (
317  Ostream& os,
318  const primitiveMesh& mesh,
319  const label maxLen
320 ) const
321 {
322  pointSet::writeDebug(os, mesh, maxLen);
323 }
324 
325 
326 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
327 
328 } // End namespace Foam
329 
330 // ************************************************************************* //
label n
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
bool insert(const label &key)
Insert a new entry.
Definition: HashSet.H:109
bool set(const label &key)
Same as insert (cannot overwrite nil content)
Definition: HashSet.H:123
List< Key > sortedToc() const
Return the table of contents as a sorted list.
Definition: HashTable.C:242
void clearStorage()
Clear the table entries and the table itself.
Definition: HashTable.C:566
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:138
void resize(const label newSize)
Resize the hash table for efficiency.
Definition: HashTable.C:506
readOption
Enumeration defining the read options.
Definition: IOobject.H:117
@ MUST_READ_IF_MODIFIED
Definition: IOobject.H:119
const word & name() const
Return name.
Definition: IOobject.H:310
writeOption
Enumeration defining the write options.
Definition: IOobject.H:126
Version number type.
Definition: IOstream.H:97
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:87
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:194
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void setSize(const label)
Reset size of List.
Definition: List.C:281
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
autoPtr< T > set(const label, const word &key, T *)
Set element to pointer provided and return old element.
label findIndex(const word &key) const
Return the index of the given the key or -1 if not found.
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
A set of point labels.
Definition: pointSet.H:51
virtual void sync(const polyMesh &mesh)
Sync set across coupled patches. Adds coupled points to set.
Definition: pointSet.C:113
virtual void writeDebug(Ostream &os, const primitiveMesh &, const label maxLen) const
Update any stored data for new labels.
Definition: pointSet.C:160
Like pointSet but -reads data from pointZone -updates pointZone when writing.
Definition: pointZoneSet.H:53
virtual ~pointZoneSet()
Destructor.
Definition: pointZoneSet.C:133
virtual label maxSize(const polyMesh &mesh) const
Return max index+1.
Definition: pointZoneSet.C:242
pointZoneSet(const polyMesh &mesh, const word &name, readOption r=MUST_READ, writeOption w=NO_WRITE)
Construct from objectRegistry and name.
Definition: pointZoneSet.C:68
virtual void invert(const label maxLen)
Invert contents. (insert all members 0..maxLen-1 which were not in.
Definition: pointZoneSet.C:139
virtual bool writeObject(IOstream::streamFormat, IOstream::versionNumber, IOstream::compressionType, const bool write=true) const
Write pointZone.
Definition: pointZoneSet.C:249
virtual void deleteSet(const topoSet &set)
Delete elements present in set.
Definition: pointZoneSet.C:210
virtual void sync(const polyMesh &mesh)
Sync pointZoneSet across coupled patches.
Definition: pointZoneSet.C:232
virtual void addSet(const topoSet &set)
Add elements present in set.
Definition: pointZoneSet.C:189
const labelList & addressing() const
Definition: pointZoneSet.H:104
virtual void writeDebug(Ostream &os, const primitiveMesh &, const label maxLen) const
Write maxLen items with label and coordinates.
Definition: pointZoneSet.C:316
virtual void subset(const topoSet &set)
Subset contents. Only elements present in both sets remain.
Definition: pointZoneSet.C:168
virtual void topoChange(const polyTopoChangeMap &map)
Update any stored data for new labels.
Definition: pointZoneSet.C:291
void updateSet()
Sort addressing and make pointSet part consistent with addressing.
Definition: pointZoneSet.C:50
A subset of mesh points. The labels of points in the zone can be obtained from the addressing() list.
Definition: pointZone.H:59
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
const pointZoneList & pointZones() const
Return point zones.
Definition: polyMesh.H:440
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const labelList & reversePointMap() const
Reverse point map.
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:75
label nPoints() const
virtual bool writeObject(IOstream::streamFormat, IOstream::versionNumber, IOstream::compressionType, const bool write) const
Write using given format, version and compression.
virtual bool write(const bool write=true) const
Write using setting from DB.
General set of labels of mesh quantity (points, cells, faces).
Definition: topoSet.H:65
void check(const label maxLabel)
Check validity of contents.
A class for handling words, derived from string.
Definition: word.H:62
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))
const dimensionedScalar c
Speed of light in a vacuum.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
Namespace for OpenFOAM.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:129
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
int order(const scalar s)
defineTypeNameAndDebug(combustionModel, 0)
void sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
void inplaceReorder(const labelUList &oldToNew, ListType &)
Inplace reorder the elements of a list.
label newPointi
Definition: readKivaGrid.H:495