mergePolyMesh.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 "mergePolyMesh.H"
27 #include "Time.H"
28 #include "polyTopoChangeMap.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(mergePolyMesh, 0);
35 }
36 
37 
38 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
39 
40 Foam::label Foam::mergePolyMesh::patchIndex(const polyPatch& p)
41 {
42  // Find the patch name on the list. If the patch is already there
43  // and patch types match, return index
44  const word& pType = p.type();
45  const word& pName = p.name();
46 
47  bool nameFound = false;
48 
49  forAll(patchNames_, patchi)
50  {
51  if (patchNames_[patchi] == pName)
52  {
53  if (word(patchDicts_[patchi]["type"]) == pType)
54  {
55  // Found name and types match
56  return patchi;
57  }
58  else
59  {
60  // Found the name, but type is different
61  nameFound = true;
62  }
63  }
64  }
65 
66  // Patch not found. Append to the list
67  {
68  OStringStream os;
69  p.write(os);
70  patchDicts_.append(dictionary(IStringStream(os.str())()));
71  }
72 
73  if (nameFound)
74  {
75  // Duplicate name is not allowed. Create a composite name from the
76  // patch name and case name
77  const word& caseName = p.boundaryMesh().mesh().time().caseName();
78 
79  patchNames_.append(pName + "_" + caseName);
80 
81  Info<< "label patchIndex(const polyPatch& p) : "
82  << "Patch " << p.index() << " named "
83  << pName << " in mesh " << caseName
84  << " already exists, but patch types "
85  << " do not match.\nCreating a composite name as "
86  << patchNames_.last() << endl;
87  }
88  else
89  {
90  patchNames_.append(pName);
91  }
92 
93  return patchNames_.size() - 1;
94 }
95 
96 
97 Foam::label Foam::mergePolyMesh::zoneIndex
98 (
99  DynamicList<word>& names,
100  const word& curName
101 )
102 {
103  forAll(names, zonei)
104  {
105  if (names[zonei] == curName)
106  {
107  return zonei;
108  }
109  }
110 
111  // Not found. Add new name to the list
112  names.append(curName);
113 
114  return names.size() - 1;
115 }
116 
117 
118 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
119 
121 :
122  mesh_(mesh),
123  meshMod_(mesh_),
124  patchNames_(2*mesh_.boundaryMesh().size()),
125  patchDicts_(2*mesh_.boundaryMesh().size()),
126  pointZoneNames_(),
127  faceZoneNames_(),
128  cellZoneNames_()
129 {
130  // Insert the original patches into the list
131  wordList curPatchNames = mesh_.boundaryMesh().names();
132 
133  forAll(mesh_.boundaryMesh(), patchi)
134  {
135  patchNames_.append(mesh_.boundaryMesh()[patchi].name());
136 
137  OStringStream os;
138  mesh_.boundaryMesh()[patchi].write(os);
139  patchDicts_.append(dictionary(IStringStream(os.str())()));
140  }
141 
142  // Insert point, face and cell zones into the list
143 
144  // Point zones
145  wordList curPointZoneNames = mesh_.pointZones().toc();
146  if (curPointZoneNames.size())
147  {
148  pointZoneNames_.setCapacity(2*curPointZoneNames.size());
149  }
150 
151  forAll(curPointZoneNames, zonei)
152  {
153  pointZoneNames_.append(curPointZoneNames[zonei]);
154  }
155 
156  pointZonesAddedPoints_.setSize(pointZoneNames_.size());
157 
158  // Face zones
159  wordList curFaceZoneNames = mesh_.faceZones().toc();
160 
161  if (curFaceZoneNames.size())
162  {
163  faceZoneNames_.setCapacity(2*curFaceZoneNames.size());
164  }
165 
166  forAll(curFaceZoneNames, zonei)
167  {
168  faceZoneNames_.append(curFaceZoneNames[zonei]);
169  }
170 
171  faceZonesAddedFaces_.setSize(faceZoneNames_.size());
172 
173  // Cell zones
174  wordList curCellZoneNames = mesh_.cellZones().toc();
175 
176  if (curCellZoneNames.size())
177  {
178  cellZoneNames_.setCapacity(2*curCellZoneNames.size());
179  }
180 
181  forAll(curCellZoneNames, zonei)
182  {
183  cellZoneNames_.append(curCellZoneNames[zonei]);
184  }
185 
186  cellZonesAddedCells_.setSize(cellZoneNames_.size());
187 }
188 
189 
190 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
191 
192 void Foam::mergePolyMesh::addMesh(const polyMesh& m)
193 {
194  // Add all the points, faces and cells of the new mesh
195 
196  // Add points
197 
198  const pointField& p = m.points();
199  labelList renumberPoints(p.size());
200 
201  const pointZoneList& pz = m.pointZones();
202  labelList pointZoneIndices(pz.size());
203 
204  forAll(pz, zonei)
205  {
206  pointZoneIndices[zonei] = zoneIndex(pointZoneNames_, pz[zonei].name());
207  pointZonesAddedPoints_.setSize(pointZoneNames_.size());
208  }
209 
210  forAll(p, pointi)
211  {
212  renumberPoints[pointi] = meshMod_.addPoint
213  (
214  p[pointi], // Point to add
215  -1, // Master point (straight addition)
216  pointi < m.nPoints() // Is in cell?
217  );
218 
219  const labelList zones(pz.whichZones(pointi));
220  forAll(zones, zonei)
221  {
222  pointZonesAddedPoints_[pointZoneIndices[zonei]]
223  .insert(renumberPoints[pointi]);
224  }
225  }
226 
227  // Add cells
228 
229  const cellList& c = m.cells();
230  labelList renumberCells(c.size());
231 
232  const cellZoneList& cz = m.cellZones();
233  labelList cellZoneIndices(cz.size());
234 
235  forAll(cz, zonei)
236  {
237  cellZoneIndices[zonei] = zoneIndex(cellZoneNames_, cz[zonei].name());
238  cellZonesAddedCells_.setSize(cellZoneNames_.size());
239  }
240 
241  forAll(c, celli)
242  {
243  renumberCells[celli] = meshMod_.addCell(-1);
244 
245  const labelList zones(cz.whichZones(celli));
246  forAll(zones, zonei)
247  {
248  cellZonesAddedCells_[cellZoneIndices[zonei]]
249  .insert(renumberCells[celli]);
250  }
251  }
252 
253  // Add faces
254  const polyBoundaryMesh& bm = m.boundaryMesh();
255 
256  // Gather the patch indices
257  labelList patchIndices(bm.size());
258 
259  forAll(patchIndices, patchi)
260  {
261  patchIndices[patchi] = patchIndex(bm[patchi]);
262  }
263 
264  // Temporary: update number of allowable patches. This should be
265  // determined at the top - before adding anything.
266  meshMod_.setNumPatches(patchNames_.size());
267 
268 
269 
270  const faceZoneList& fz = m.faceZones();
271  labelList faceZoneIndices(fz.size());
272 
273  forAll(fz, zonei)
274  {
275  faceZoneIndices[zonei] = zoneIndex(faceZoneNames_, fz[zonei].name());
276  faceZonesAddedFaces_.setSize(faceZoneNames_.size());
277  }
278 
279  const faceList& f = m.faces();
280  labelList renumberFaces(f.size());
281 
282  const labelList& own = m.faceOwner();
283  const labelList& nei = m.faceNeighbour();
284 
285  label newOwn, newNei, newPatch;
286 
287  forAll(f, facei)
288  {
289  const face& curFace = f[facei];
290 
291  face newFace(curFace.size());
292 
293  forAll(curFace, pointi)
294  {
295  newFace[pointi] = renumberPoints[curFace[pointi]];
296  }
297 
298  if (debug)
299  {
300  // Check that the face is valid
301  if (min(newFace) < 0)
302  {
304  << "Error in point mapping for face " << facei
305  << ". Old face: " << curFace << " New face: " << newFace
306  << abort(FatalError);
307  }
308  }
309 
310  if (facei < m.nInternalFaces() || facei >= m.nFaces())
311  {
312  newPatch = -1;
313  }
314  else
315  {
316  newPatch = patchIndices[bm.whichPatch(facei)];
317  }
318 
319  newOwn = own[facei];
320  if (newOwn > -1) newOwn = renumberCells[newOwn];
321 
322  if (newPatch > -1)
323  {
324  newNei = -1;
325  }
326  else
327  {
328  newNei = nei[facei];
329  newNei = renumberCells[newNei];
330  }
331 
332  renumberFaces[facei] = meshMod_.addFace
333  (
334  newFace,
335  newOwn,
336  newNei,
337  -1,
338  false,
339  newPatch
340  );
341 
342  const labelList zones(fz.whichZones(facei));
343  forAll(zones, zonei)
344  {
345  const faceZone& fzi = fz[zones[zonei]];
346  const bool flip = fzi.flipMap()[fzi.localIndex(facei)];
347 
348  faceZonesAddedFaces_[faceZoneIndices[zonei]]
349  .insert(renumberFaces[facei], flip);
350  }
351  }
352 }
353 
354 
356 {
357  if (debug)
358  {
359  Info<< "patch names: " << patchNames_ << nl
360  << "patch dicts: " << patchDicts_ << nl
361  << "point zone names: " << pointZoneNames_ << nl
362  << "face zone names: " << faceZoneNames_ << nl
363  << "cell zone names: " << cellZoneNames_ << endl;
364  }
365 
366  // Add the patches if necessary
367  if (patchNames_.size() != mesh_.boundaryMesh().size())
368  {
369  if (debug)
370  {
371  Info<< "Copying old patches" << endl;
372  }
373 
374  List<polyPatch*> newPatches(patchNames_.size());
375 
376  const polyBoundaryMesh& oldPatches = mesh_.boundaryMesh();
377 
378  // Note. Re-using counter in two for loops
379  label patchi = 0;
380 
381  for (patchi = 0; patchi < oldPatches.size(); patchi++)
382  {
383  newPatches[patchi] = oldPatches[patchi].clone(oldPatches).ptr();
384  }
385 
386  if (debug)
387  {
388  Info<< "Adding new patches. " << endl;
389  }
390 
391  label endOfLastPatch =
392  patchi == 0
393  ? 0
394  : oldPatches[patchi - 1].start() + oldPatches[patchi - 1].size();
395 
396  for (; patchi < patchNames_.size(); patchi++)
397  {
398  // Add a patch
399  dictionary dict(patchDicts_[patchi]);
400  dict.set("nFaces", 0);
401  dict.set("startFace", endOfLastPatch);
402 
403  newPatches[patchi] =
404  (
406  (
407  patchNames_[patchi],
408  dict,
409  patchi,
410  oldPatches
411  ).ptr()
412  );
413  }
414 
415  mesh_.removeBoundary();
416  mesh_.addPatches(newPatches);
417  }
418 
419  // Add the zones if necessary
420  if (pointZoneNames_.size() > mesh_.pointZones().size())
421  {
422  if (debug)
423  {
424  Info<< "Adding new pointZones. " << endl;
425  }
426 
427  label nZones = mesh_.pointZones().size();
428 
429  mesh_.pointZones().setSize(pointZoneNames_.size());
430 
431  for (label zonei = nZones; zonei < pointZoneNames_.size(); zonei++)
432  {
433  mesh_.pointZones().set
434  (
435  zonei,
436  new pointZone
437  (
438  pointZoneNames_[zonei],
439  labelList(),
440  mesh_.pointZones()
441  )
442  );
443  }
444  }
445 
446  if (cellZoneNames_.size() > mesh_.cellZones().size())
447  {
448  if (debug)
449  {
450  Info<< "Adding new cellZones. " << endl;
451  }
452 
453  label nZones = mesh_.cellZones().size();
454 
455  mesh_.cellZones().setSize(cellZoneNames_.size());
456 
457  for (label zonei = nZones; zonei < cellZoneNames_.size(); zonei++)
458  {
459  mesh_.cellZones().set
460  (
461  zonei,
462  new cellZone
463  (
464  cellZoneNames_[zonei],
465  labelList(),
466  mesh_.cellZones()
467  )
468  );
469  }
470  }
471 
472  if (faceZoneNames_.size() > mesh_.faceZones().size())
473  {
474  if (debug)
475  {
476  Info<< "Adding new faceZones. " << endl;
477  }
478 
479  label nZones = mesh_.faceZones().size();
480 
481  mesh_.faceZones().setSize(faceZoneNames_.size());
482 
483  for (label zonei = nZones; zonei < faceZoneNames_.size(); zonei++)
484  {
485  mesh_.faceZones().set
486  (
487  zonei,
488  new faceZone
489  (
490  faceZoneNames_[zonei],
491  labelList(),
492  boolList(),
493  mesh_.faceZones()
494  )
495  );
496  }
497  }
498 
499  // Change mesh
500  autoPtr<polyTopoChangeMap> map(meshMod_.changeMesh(mesh_));
501 
502  // Add the new points to the pointZones in the merged mesh
503  mesh_.pointZones().insert(pointZonesAddedPoints_);
504 
505  // Add the new faces to the faceZones in the merged mesh
506  mesh_.faceZones().insert(faceZonesAddedFaces_);
507 
508  // Add the new cells to the cellZones in the merged mesh
509  mesh_.cellZones().insert(cellZonesAddedCells_);
510 
511  mesh_.topoChange(map);
512 
513  // Clear topo change for the next operation
514  meshMod_.clear();
515 }
516 
517 
518 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
T & last()
Return the last element of the list.
Definition: UListI.H:128
void set(entry *)
Assign a new entry, overwrite any existing entry.
Definition: dictionary.C:1152
void addMesh(const polyMesh &m)
Add a mesh.
mergePolyMesh(polyMesh &mesh)
Construct from polyMesh.
void merge()
Merge meshes.
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
void setNumPatches(const label nPatches)
Explicitly set the number of patches if construct-without-mesh.
label addCell(const label masterCellID)
Add cell and return new cell index.
label addPoint(const point &, const label masterPointID, const bool inCell)
Add point and return new point index.
label addFace(const face &f, const label own, const label nei, const label masterFaceID, const bool flipFaceFlux, const label patchID)
Add face to cells and return new face index.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
const dimensionedScalar c
Speed of light in a vacuum.
Namespace for OpenFOAM.
List< word > wordList
A List of words.
Definition: fileName.H:54
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
List< cell > cellList
list of cells
Definition: cellList.H:42
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
errorManip< error > abort(error &err)
Definition: errorManip.H:131
messageStream Info
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
defineTypeNameAndDebug(combustionModel, 0)
error FatalError
List< face > faceList
Definition: faceListFwd.H:41
static const char nl
Definition: Ostream.H:266
labelList f(nPoints)
dictionary dict
volScalarField & p