createShellMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2019-2022 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "createShellMesh.H"
30 #include "polyTopoChange.H"
31 #include "meshTools.H"
32 #include "mapPolyMesh.H"
33 #include "polyAddPoint.H"
34 #include "polyAddFace.H"
35 #include "polyModifyFace.H"
36 #include "polyAddCell.H"
37 #include "labelPair.H"
38 #include "indirectPrimitivePatch.H"
39 #include "mapDistribute.H"
40 #include "globalMeshData.H"
41 #include "PatchTools.H"
42 #include "globalIndex.H"
43 
44 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
49 
50 template<>
51 class minEqOp<labelPair>
52 {
53 public:
54  void operator()(labelPair& x, const labelPair& y) const
55  {
56  x[0] = min(x[0], y[0]);
57  x[1] = min(x[1], y[1]);
58  }
59 };
60 }
61 
62 
63 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
64 
65 // Synchronise edges
66 void Foam::createShellMesh::syncEdges
67 (
68  const globalMeshData& globalData,
69 
70  const labelList& patchEdges,
71  const labelList& coupledEdges,
72  const bitSet& sameEdgeOrientation,
73  const bool syncNonCollocated,
74 
75  bitSet& isChangedEdge,
76  DynamicList<label>& changedEdges,
77  labelPairList& allEdgeData
78 )
79 {
80  const mapDistribute& map = globalData.globalEdgeSlavesMap();
81  const bitSet& cppOrientation = globalData.globalEdgeOrientation();
82 
83  // Convert patch-edge data into cpp-edge data
84  labelPairList cppEdgeData
85  (
86  map.constructSize(),
88  );
89 
90  forAll(patchEdges, i)
91  {
92  label patchEdgeI = patchEdges[i];
93  label coupledEdgeI = coupledEdges[i];
94 
95  if (isChangedEdge[patchEdgeI])
96  {
97  const labelPair& data = allEdgeData[patchEdgeI];
98 
99  // Patch-edge data needs to be converted into coupled-edge data
100  // (optionally flipped) and consistent in orientation with
101  // other coupled edge (optionally flipped)
102  if (sameEdgeOrientation[i] == cppOrientation[coupledEdgeI])
103  {
104  cppEdgeData[coupledEdgeI] = data;
105  }
106  else
107  {
108  cppEdgeData[coupledEdgeI] = labelPair(data[1], data[0]);
109  }
110  }
111  }
112 
113  // Synchronise
114  globalData.syncData
115  (
116  cppEdgeData,
117  globalData.globalEdgeSlaves(),
118  (
119  syncNonCollocated
120  ? globalData.globalEdgeTransformedSlaves() // transformed elems
121  : labelListList(globalData.globalEdgeSlaves().size()) //no transformed
122  ),
123  map,
124  minEqOp<labelPair>()
125  );
126 
127  // Back from cpp-edge to patch-edge data
128  forAll(patchEdges, i)
129  {
130  label patchEdgeI = patchEdges[i];
131  label coupledEdgeI = coupledEdges[i];
132 
133  if (cppEdgeData[coupledEdgeI] != labelPair(labelMax, labelMax))
134  {
135  const labelPair& data = cppEdgeData[coupledEdgeI];
136 
137  if (sameEdgeOrientation[i] == cppOrientation[coupledEdgeI])
138  {
139  allEdgeData[patchEdgeI] = data;
140  }
141  else
142  {
143  allEdgeData[patchEdgeI] = labelPair(data[1], data[0]);
144  }
145 
146  if (!isChangedEdge[patchEdgeI])
147  {
148  changedEdges.append(patchEdgeI);
149  isChangedEdge.set(patchEdgeI);
150  }
151  }
152  }
153 }
154 
155 
157 (
158  const globalMeshData& globalData,
159  const primitiveFacePatch& patch,
160  const bitSet& nonManifoldEdge,
161  const bool syncNonCollocated,
162 
163  faceList& pointGlobalRegions,
164  faceList& pointLocalRegions,
165  labelList& localToGlobalRegion
166 )
167 {
168  const indirectPrimitivePatch& cpp = globalData.coupledPatch();
169 
170  // Calculate correspondence between patch and globalData.coupledPatch.
171  labelList patchEdges;
172  labelList coupledEdges;
173  bitSet sameEdgeOrientation;
175  (
176  cpp,
177  patch,
178 
179  coupledEdges,
180  patchEdges,
181  sameEdgeOrientation
182  );
183 
184 
185  // Initial unique regions
186  // ~~~~~~~~~~~~~~~~~~~~~~
187  // These get merged later on across connected edges.
188 
189  // 1. Count
190  label nMaxRegions = 0;
191  forAll(patch.localFaces(), facei)
192  {
193  const face& f = patch.localFaces()[facei];
194  nMaxRegions += f.size();
195  }
196 
197  const globalIndex globalRegions(nMaxRegions);
198 
199  // 2. Assign unique regions
200  label nRegions = 0;
201 
202  pointGlobalRegions.setSize(patch.size());
203  forAll(pointGlobalRegions, facei)
204  {
205  const face& f = patch.localFaces()[facei];
206  labelList& pRegions = pointGlobalRegions[facei];
207  pRegions.setSize(f.size());
208  forAll(pRegions, fp)
209  {
210  pRegions[fp] = globalRegions.toGlobal(nRegions++);
211  }
212  }
213 
214 
215  DynamicList<label> changedEdges(patch.nEdges());
216  labelPairList allEdgeData(patch.nEdges(), labelPair(labelMax, labelMax));
217  bitSet isChangedEdge(patch.nEdges());
218 
219 
220  // Fill initial seed
221  // ~~~~~~~~~~~~~~~~~
222 
223  forAll(patch.edgeFaces(), edgeI)
224  {
225  if (!nonManifoldEdge[edgeI])
226  {
227  // Take over value from one face only.
228  const edge& e = patch.edges()[edgeI];
229  label facei = patch.edgeFaces()[edgeI][0];
230  const face& f = patch.localFaces()[facei];
231 
232  label fp0 = f.find(e[0]);
233  label fp1 = f.find(e[1]);
234  allEdgeData[edgeI] = labelPair
235  (
236  pointGlobalRegions[facei][fp0],
237  pointGlobalRegions[facei][fp1]
238  );
239  if (!isChangedEdge[edgeI])
240  {
241  changedEdges.append(edgeI);
242  isChangedEdge.set(edgeI);
243  }
244  }
245  }
246 
247 
248  syncEdges
249  (
250  globalData,
251 
252  patchEdges,
253  coupledEdges,
254  sameEdgeOrientation,
255  syncNonCollocated,
256 
257  isChangedEdge,
258  changedEdges,
259  allEdgeData
260  );
261 
262 
263  // Edge-Face-Edge walk across patch
264  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
265  // Across edge minimum regions win
266 
267  while (true)
268  {
269  // From edge to face
270  // ~~~~~~~~~~~~~~~~~
271 
272  DynamicList<label> changedFaces(patch.size());
273  bitSet isChangedFace(patch.size());
274 
275  forAll(changedEdges, changedI)
276  {
277  label edgeI = changedEdges[changedI];
278  const labelPair& edgeData = allEdgeData[edgeI];
279 
280  const edge& e = patch.edges()[edgeI];
281  const labelList& eFaces = patch.edgeFaces()[edgeI];
282 
283  forAll(eFaces, i)
284  {
285  label facei = eFaces[i];
286  const face& f = patch.localFaces()[facei];
287 
288  // Combine edgeData with face data
289  label fp0 = f.find(e[0]);
290  if (pointGlobalRegions[facei][fp0] > edgeData[0])
291  {
292  pointGlobalRegions[facei][fp0] = edgeData[0];
293  if (!isChangedFace[facei])
294  {
295  isChangedFace.set(facei);
296  changedFaces.append(facei);
297  }
298  }
299 
300  label fp1 = f.find(e[1]);
301  if (pointGlobalRegions[facei][fp1] > edgeData[1])
302  {
303  pointGlobalRegions[facei][fp1] = edgeData[1];
304  if (!isChangedFace[facei])
305  {
306  isChangedFace.set(facei);
307  changedFaces.append(facei);
308  }
309  }
310  }
311  }
312 
313 
314  if (returnReduceAnd(changedFaces.empty()))
315  {
316  break;
317  }
318 
319 
320  // From face to edge
321  // ~~~~~~~~~~~~~~~~~
322 
323  isChangedEdge = false;
324  changedEdges.clear();
325 
326  forAll(changedFaces, i)
327  {
328  label facei = changedFaces[i];
329  const face& f = patch.localFaces()[facei];
330  const labelList& fEdges = patch.faceEdges()[facei];
331 
332  forAll(fEdges, fp)
333  {
334  label edgeI = fEdges[fp];
335 
336  if (!nonManifoldEdge[edgeI])
337  {
338  const edge& e = patch.edges()[edgeI];
339  label fp0 = f.find(e[0]);
340  label region0 = pointGlobalRegions[facei][fp0];
341  label fp1 = f.find(e[1]);
342  label region1 = pointGlobalRegions[facei][fp1];
343 
344  if
345  (
346  (allEdgeData[edgeI][0] > region0)
347  || (allEdgeData[edgeI][1] > region1)
348  )
349  {
350  allEdgeData[edgeI] = labelPair(region0, region1);
351  if (!isChangedEdge[edgeI])
352  {
353  changedEdges.append(edgeI);
354  isChangedEdge.set(edgeI);
355  }
356  }
357  }
358  }
359  }
360 
361  syncEdges
362  (
363  globalData,
364 
365  patchEdges,
366  coupledEdges,
367  sameEdgeOrientation,
368  syncNonCollocated,
369 
370  isChangedEdge,
371  changedEdges,
372  allEdgeData
373  );
374 
375 
376  if (returnReduceAnd(changedEdges.empty()))
377  {
378  break;
379  }
380  }
381 
382 
383 
384  // Assign local regions
385  // ~~~~~~~~~~~~~~~~~~~~
386 
387  // Calculate addressing from global region back to local region
388  pointLocalRegions.setSize(patch.size());
389  Map<label> globalToLocalRegion(globalRegions.localSize()/4);
390  DynamicList<label> dynLocalToGlobalRegion(globalToLocalRegion.size());
391  forAll(patch.localFaces(), facei)
392  {
393  const face& f = patch.localFaces()[facei];
394  face& pRegions = pointLocalRegions[facei];
395  pRegions.setSize(f.size());
396  forAll(f, fp)
397  {
398  label globalRegionI = pointGlobalRegions[facei][fp];
399 
400  const auto fnd = globalToLocalRegion.cfind(globalRegionI);
401 
402  if (fnd.good())
403  {
404  // Already encountered this global region. Assign same local one
405  pRegions[fp] = fnd();
406  }
407  else
408  {
409  // Region not yet seen. Create new one
410  label localRegionI = globalToLocalRegion.size();
411  pRegions[fp] = localRegionI;
412  globalToLocalRegion.insert(globalRegionI, localRegionI);
413  dynLocalToGlobalRegion.append(globalRegionI);
414  }
415  }
416  }
417  localToGlobalRegion.transfer(dynLocalToGlobalRegion);
418 }
419 
420 
421 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
422 
423 Foam::createShellMesh::createShellMesh
424 (
425  const primitiveFacePatch& patch,
426  const faceList& pointRegions,
427  const labelList& regionPoints
428 )
429 :
430  patch_(patch),
431  pointRegions_(pointRegions),
432  regionPoints_(regionPoints)
433 {
434  if (pointRegions_.size() != patch_.size())
435  {
437  << "nFaces:" << patch_.size()
438  << " pointRegions:" << pointRegions.size()
439  << exit(FatalError);
440  }
441 }
442 
443 
444 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
445 
447 (
448  const pointField& firstLayerDisp,
449  const scalar expansionRatio,
450  const label nLayers,
451  const labelList& topPatchID,
452  const labelList& bottomPatchID,
453  const labelListList& extrudeEdgePatches,
454  polyTopoChange& meshMod
455 )
456 {
457  if (firstLayerDisp.size() != regionPoints_.size())
458  {
460  << "nRegions:" << regionPoints_.size()
461  << " firstLayerDisp:" << firstLayerDisp.size()
462  << exit(FatalError);
463  }
464 
465  if
466  (
467  topPatchID.size() != patch_.size()
468  && bottomPatchID.size() != patch_.size()
469  )
470  {
472  << "nFaces:" << patch_.size()
473  << " topPatchID:" << topPatchID.size()
474  << " bottomPatchID:" << bottomPatchID.size()
475  << exit(FatalError);
476  }
477 
478  if (extrudeEdgePatches.size() != patch_.nEdges())
479  {
481  << "nEdges:" << patch_.nEdges()
482  << " extrudeEdgePatches:" << extrudeEdgePatches.size()
483  << exit(FatalError);
484  }
485 
486 
487 
488  // From cell to patch (trivial)
489  DynamicList<label> cellToFaceMap(nLayers*patch_.size());
490  // From face to patch+turning index
491  DynamicList<label> faceToFaceMap
492  (
493  (nLayers+1)*(patch_.size()+patch_.nEdges())
494  );
495  // From face to patch edge index
496  DynamicList<label> faceToEdgeMap(nLayers*(patch_.nEdges()+patch_.nEdges()));
497  // From point to patch point index
498  DynamicList<label> pointToPointMap((nLayers+1)*patch_.nPoints());
499 
500 
501  // Introduce new cell for every face
502  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
503 
504  labelList addedCells(nLayers*patch_.size());
505  forAll(patch_, facei)
506  {
507  for (label layerI = 0; layerI < nLayers; layerI++)
508  {
509  addedCells[nLayers*facei+layerI] = meshMod.addCell
510  (
511  -1, // masterPointID
512  -1, // masterEdgeID
513  -1, // masterFaceID
514  cellToFaceMap.size(), // masterCellID
515  -1 // zoneID
516  );
517  cellToFaceMap.append(facei);
518  }
519  }
520 
521 
522  // Introduce original points
523  // ~~~~~~~~~~~~~~~~~~~~~~~~~
524 
525  // Original point numbers in local point ordering so no need to store.
526  forAll(patch_.localPoints(), pointi)
527  {
528  //label addedPointi =
529  meshMod.addPoint
530  (
531  patch_.localPoints()[pointi], // point
532  pointToPointMap.size(), // masterPointID
533  -1, // zoneID
534  true // inCell
535  );
536  pointToPointMap.append(pointi);
537 
538  //Pout<< "Added bottom point " << addedPointi
539  // << " at " << patch_.localPoints()[pointi]
540  // << " from point " << pointi
541  // << endl;
542  }
543 
544 
545  // Introduce new points (one for every region)
546  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
547 
548  labelList addedPoints(nLayers*regionPoints_.size());
549  forAll(regionPoints_, regionI)
550  {
551  label pointi = regionPoints_[regionI];
552 
553  point pt = patch_.localPoints()[pointi];
554  point disp = firstLayerDisp[regionI];
555  for (label layerI = 0; layerI < nLayers; layerI++)
556  {
557  pt += disp;
558 
559  addedPoints[nLayers*regionI+layerI] = meshMod.addPoint
560  (
561  pt, // point
562  pointToPointMap.size(), // masterPointID - used only addressing
563  -1, // zoneID
564  true // inCell
565  );
566  pointToPointMap.append(pointi);
567 
568  disp *= expansionRatio;
569  }
570  }
571 
572 
573  // Add face on bottom side
574  forAll(patch_.localFaces(), facei)
575  {
576  meshMod.addFace
577  (
578  patch_.localFaces()[facei].reverseFace(),// vertices
579  addedCells[nLayers*facei], // own
580  -1, // nei
581  -1, // masterPointID
582  -1, // masterEdgeID
583  faceToFaceMap.size(), // masterFaceID : current facei
584  true, // flipFaceFlux
585  bottomPatchID[facei], // patchID
586  -1, // zoneID
587  false // zoneFlip
588  );
589  faceToFaceMap.append(-facei-1); // points to flipped original face
590  faceToEdgeMap.append(-1);
591 
592  //const face newF(patch_.localFaces()[facei].reverseFace());
593  //Pout<< "Added bottom face "
594  // << newF
595  // << " coords:" << UIndirectList<point>(meshMod.points(), newF)
596  // << " own " << addedCells[facei]
597  // << " patch:" << bottomPatchID[facei]
598  // << " at " << patch_.faceCentres()[facei]
599  // << endl;
600  }
601 
602  // Add inbetween faces and face on top
603  forAll(patch_.localFaces(), facei)
604  {
605  // Get face in original ordering
606  const face& f = patch_.localFaces()[facei];
607 
608  face newF(f.size());
609 
610  for (label layerI = 0; layerI < nLayers; layerI++)
611  {
612  // Pick up point based on region and layer
613  forAll(f, fp)
614  {
615  label region = pointRegions_[facei][fp];
616  newF[fp] = addedPoints[region*nLayers+layerI];
617  }
618 
619  label own = addedCells[facei*nLayers+layerI];
620  label nei;
621  label patchi;
622  if (layerI == nLayers-1)
623  {
624  nei = -1;
625  patchi = topPatchID[facei];
626  }
627  else
628  {
629  nei = addedCells[facei*nLayers+layerI+1];
630  patchi = -1;
631  }
632 
633  meshMod.addFace
634  (
635  newF, // vertices
636  own, // own
637  nei, // nei
638  -1, // masterPointID
639  -1, // masterEdgeID
640  faceToFaceMap.size(), // masterFaceID : current facei
641  false, // flipFaceFlux
642  patchi, // patchID
643  -1, // zoneID
644  false // zoneFlip
645  );
646  faceToFaceMap.append(facei+1); // unflipped
647  faceToEdgeMap.append(-1);
648 
649  //Pout<< "Added inbetween face " << newF
650  // << " coords:" << UIndirectList<point>(meshMod.points(), newF)
651  // << " at layer " << layerI
652  // << " own " << own
653  // << " nei " << nei
654  // << " at " << patch_.faceCentres()[facei]
655  // << endl;
656  }
657  }
658 
659 
660  // Add side faces
661  // ~~~~~~~~~~~~~~
662  // Note that we loop over edges multiple times so for edges with
663  // two cyclic faces they get added in two passes (for correct ordering)
664 
665  // Pass1. Internal edges and first face of other edges
666  forAll(extrudeEdgePatches, edgeI)
667  {
668  const labelList& eFaces = patch_.edgeFaces()[edgeI];
669  const labelList& ePatches = extrudeEdgePatches[edgeI];
670 
671  if (ePatches.size() == 0)
672  {
673  // Internal face
674  if (eFaces.size() != 2)
675  {
677  << "edge:" << edgeI
678  << " not internal but does not have side-patches defined."
679  << exit(FatalError);
680  }
681  }
682  else
683  {
684  if (eFaces.size() != ePatches.size())
685  {
687  << "external/feature edge:" << edgeI
688  << " has " << eFaces.size() << " connected extruded faces "
689  << " but only " << ePatches.size()
690  << " boundary faces defined." << exit(FatalError);
691  }
692  }
693 
694 
695 
696  // Make face pointing in to eFaces[0] so out of new master face
697  const face& f = patch_.localFaces()[eFaces[0]];
698  const edge& e = patch_.edges()[edgeI];
699 
700  label fp0 = f.find(e[0]);
701  label fp1 = f.fcIndex(fp0);
702 
703  if (f[fp1] != e[1])
704  {
705  fp1 = fp0;
706  fp0 = f.rcIndex(fp1);
707  }
708 
709  face newF(4);
710 
711  for (label layerI = 0; layerI < nLayers; layerI++)
712  {
713  label region0 = pointRegions_[eFaces[0]][fp0];
714  label region1 = pointRegions_[eFaces[0]][fp1];
715 
716  // Pick up points with correct normal
717  if (layerI == 0)
718  {
719  newF[0] = f[fp0];
720  newF[1] = f[fp1];
721  newF[2] = addedPoints[nLayers*region1+layerI];
722  newF[3] = addedPoints[nLayers*region0+layerI];
723  }
724  else
725  {
726  newF[0] = addedPoints[nLayers*region0+layerI-1];
727  newF[1] = addedPoints[nLayers*region1+layerI-1];
728  newF[2] = addedPoints[nLayers*region1+layerI];
729  newF[3] = addedPoints[nLayers*region0+layerI];
730  }
731 
732  // Optionally rotate so e[0] is always 0th vertex. Note that
733  // this normally is automatically done by coupled face ordering
734  // but with NOORDERING we have to do it ourselves.
735  if (f[fp0] != e[0])
736  {
737  // rotate one back to get newF[1] (originating from e[0])
738  // into newF[0]
739  label v0 = newF[0];
740  for (label i = 0; i < newF.size()-1; i++)
741  {
742  newF[i] = newF[newF.fcIndex(i)];
743  }
744  newF.last() = v0;
745  }
746 
747 
748  label minCelli = addedCells[nLayers*eFaces[0]+layerI];
749  label maxCelli;
750  label patchi;
751  if (ePatches.size() == 0)
752  {
753  maxCelli = addedCells[nLayers*eFaces[1]+layerI];
754  if (minCelli > maxCelli)
755  {
756  // Swap
757  std::swap(minCelli, maxCelli);
758  newF.flip();
759  }
760  patchi = -1;
761  }
762  else
763  {
764  maxCelli = -1;
765  patchi = ePatches[0];
766  }
767 
768  //{
769  // Pout<< "Adding from face:" << patch_.faceCentres()[eFaces[0]]
770  // << " from edge:"
771  // << patch_.localPoints()[f[fp0]]
772  // << patch_.localPoints()[f[fp1]]
773  // << " at layer:" << layerI
774  // << " with new points:" << newF
775  // << " locations:"
776  // << UIndirectList<point>(meshMod.points(), newF)
777  // << " own:" << minCelli
778  // << " nei:" << maxCelli
779  // << endl;
780  //}
781 
782 
783  // newF already outwards pointing.
784  meshMod.addFace
785  (
786  newF, // vertices
787  minCelli, // own
788  maxCelli, // nei
789  -1, // masterPointID
790  -1, // masterEdgeID
791  faceToFaceMap.size(), // masterFaceID
792  false, // flipFaceFlux
793  patchi, // patchID
794  -1, // zoneID
795  false // zoneFlip
796  );
797  faceToFaceMap.append(0);
798  faceToEdgeMap.append(edgeI);
799  }
800  }
801 
802  // Pass2. Other faces of boundary edges
803  forAll(extrudeEdgePatches, edgeI)
804  {
805  const labelList& eFaces = patch_.edgeFaces()[edgeI];
806  const labelList& ePatches = extrudeEdgePatches[edgeI];
807 
808  if (ePatches.size() >= 2)
809  {
810  for (label i = 1; i < ePatches.size(); i++)
811  {
812  // Extrude eFaces[i]
813  label minFacei = eFaces[i];
814 
815  // Make face pointing in to eFaces[0] so out of new master face
816  const face& f = patch_.localFaces()[minFacei];
817 
818  const edge& e = patch_.edges()[edgeI];
819  label fp0 = f.find(e[0]);
820  label fp1 = f.fcIndex(fp0);
821 
822  if (f[fp1] != e[1])
823  {
824  fp1 = fp0;
825  fp0 = f.rcIndex(fp1);
826  }
827 
828  face newF(4);
829  for (label layerI = 0; layerI < nLayers; layerI++)
830  {
831  label region0 = pointRegions_[minFacei][fp0];
832  label region1 = pointRegions_[minFacei][fp1];
833 
834  if (layerI == 0)
835  {
836  newF[0] = f[fp0];
837  newF[1] = f[fp1];
838  newF[2] = addedPoints[nLayers*region1+layerI];
839  newF[3] = addedPoints[nLayers*region0+layerI];
840  }
841  else
842  {
843  newF[0] = addedPoints[nLayers*region0+layerI-1];
844  newF[1] = addedPoints[nLayers*region1+layerI-1];
845  newF[2] = addedPoints[nLayers*region1+layerI];
846  newF[3] = addedPoints[nLayers*region0+layerI];
847  }
848 
849 
850  // Optionally rotate so e[0] is always 0th vertex. Note that
851  // this normally is automatically done by coupled face
852  // ordering but with NOORDERING we have to do it ourselves.
853  if (f[fp0] != e[0])
854  {
855  // rotate one back to get newF[1] (originating
856  // from e[0]) into newF[0].
857  label v0 = newF[0];
858  for (label i = 0; i < newF.size()-1; i++)
859  {
860  newF[i] = newF[newF.fcIndex(i)];
861  }
862  newF.last() = v0;
863  }
865  //{
866  // Pout<< "Adding from MULTI face:"
867  // << patch_.faceCentres()[minFacei]
868  // << " from edge:"
869  // << patch_.localPoints()[f[fp0]]
870  // << patch_.localPoints()[f[fp1]]
871  // << " at layer:" << layerI
872  // << " with new points:" << newF
873  // << " locations:"
874  // << UIndirectList<point>(meshMod.points(), newF)
875  // << endl;
876  //}
877 
878  // newF already outwards pointing.
879  meshMod.addFace
880  (
881  newF, // vertices
882  addedCells[nLayers*minFacei+layerI], // own
883  -1, // nei
884  -1, // masterPointID
885  -1, // masterEdgeID
886  faceToFaceMap.size(), // masterFaceID
887  false, // flipFaceFlux
888  ePatches[i], // patchID
889  -1, // zoneID
890  false // zoneFlip
891  );
892  faceToFaceMap.append(0);
893  faceToEdgeMap.append(edgeI);
894  }
895  }
896  }
897  }
898 
899 
900  cellToFaceMap_.transfer(cellToFaceMap);
901  faceToFaceMap_.transfer(faceToFaceMap);
902  faceToEdgeMap_.transfer(faceToEdgeMap);
903  pointToPointMap_.transfer(pointToPointMap);
904 }
905 
906 
908 {
909  inplaceReorder(map.reverseCellMap(), cellToFaceMap_);
910  inplaceReorder(map.reverseFaceMap(), faceToFaceMap_);
911  inplaceReorder(map.reverseFaceMap(), faceToEdgeMap_);
912  inplaceReorder(map.reversePointMap(), pointToPointMap_);
913 }
914 
915 
916 // ************************************************************************* //
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:116
void set(const bitSet &bitset)
Set specified bits from another bitset.
Definition: bitSetI.H:498
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:125
error FatalError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL ERROR&#39; header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:608
const mapDistribute & globalEdgeSlavesMap() const
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
void inplaceReorder(const labelUList &oldToNew, ListType &input, const bool prune=false)
Inplace reorder the elements of a list.
static void matchEdges(const PrimitivePatch< FaceList1, PointField1 > &p1, const PrimitivePatch< FaceList2, PointField2 > &p2, labelList &p1EdgeLabels, labelList &p2EdgeLabels, bitSet &sameOrientation)
Find corresponding edges on patches sharing the same points.
bool empty() const noexcept
True if List is empty (ie, size() is zero)
Definition: UList.H:675
const bitSet & globalEdgeOrientation() const
Is my edge same orientation as master edge.
label constructSize() const noexcept
Constructed data size.
List< labelPair > labelPairList
List of labelPair.
Definition: labelPair.H:33
const labelListList & globalEdgeSlaves() const
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field...
List< labelList > labelListList
List of labelList.
Definition: labelList.H:38
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:421
static void syncData(List< Type > &elems, const labelListList &slaves, const labelListList &transformedSlaves, const mapDistribute &slavesMap, const globalIndexAndTransform &, const CombineOp &cop, const TransformOp &top)
Helper: synchronise data with transforms.
const labelList & reverseFaceMap() const noexcept
Reverse face map.
Definition: mapPolyMesh.H:620
label fcIndex(const label i) const noexcept
The forward circular index. The next index in the list which returns to the first at the end of the l...
Definition: UListI.H:90
List< face > faceList
List of faces.
Definition: faceListFwd.H:39
scalar y
A list of faces which address into the list of points.
bool returnReduceAnd(const bool value, const label comm=UPstream::worldComm)
Perform logical (and) MPI Allreduce on a copy. Uses UPstream::reduceAnd.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:38
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
void setSize(const label n)
Alias for resize()
Definition: List.H:320
const labelList & reverseCellMap() const noexcept
Reverse cell map.
Definition: mapPolyMesh.H:656
Creates mesh by extruding a patch.
const labelList & reversePointMap() const noexcept
Reverse point map.
Definition: mapPolyMesh.H:582
void append(const T &val)
Copy append an element to the end of this list.
Definition: DynamicList.H:584
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:26
label find(const T &val) const
Find index of the first occurrence of the value.
Definition: UList.C:173
static void calcPointRegions(const globalMeshData &globalData, const primitiveFacePatch &patch, const bitSet &nonManifoldEdge, const bool syncNonCollocated, faceList &pointGlobalRegions, faceList &pointLocalRegions, labelList &localToGlobalRegion)
Helper: calculate point regions. The point region is the.
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:51
void updateMesh(const mapPolyMesh &)
Update any locally stored mesh information.
defineTypeNameAndDebug(combustionModel, 0)
void clear() noexcept
Clear the addressed list, i.e. set the size to zero.
Definition: DynamicListI.H:405
labelList f(nPoints)
PrimitivePatch< List< face >, const pointField & > primitiveFacePatch
A PrimitivePatch with List storage for the faces, const reference for the point field.
T & last()
Access last element of the list, position [size()-1].
Definition: UList.H:876
void operator()(T &x, const T &y) const
Definition: ops.H:76
Class containing processor-to-processor mapping information.
vector point
Point is a vector.
Definition: point.H:37
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition: bitSet.H:59
const std::string patch
OpenFOAM patch number as a std::string.
constexpr label labelMax
Definition: label.H:55
void setRefinement(const pointField &firstLayerThickness, const scalar expansionRatio, const label nLayers, const labelList &topPatchID, const labelList &bottomPatchID, const labelListList &extrudeEdgePatches, polyTopoChange &meshMod)
Play commands into polyTopoChange to create layer mesh.
List< label > labelList
A List of labels.
Definition: List.H:62
label rcIndex(const label i) const noexcept
The reverse circular index. The previous index in the list which returns to the last at the beginning...
Definition: UListI.H:97
Namespace for OpenFOAM.
const labelListList & globalEdgeTransformedSlaves() const