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