decompositionMethod.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-2014 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 InClass
25  decompositionMethod
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "decompositionMethod.H"
30 #include "globalIndex.H"
31 #include "syncTools.H"
32 #include "Tuple2.H"
33 #include "faceSet.H"
34 #include "regionSplit.H"
35 #include "localPointRegion.H"
36 #include "minData.H"
37 #include "FaceCellWave.H"
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43  defineTypeNameAndDebug(decompositionMethod, 0);
44  defineRunTimeSelectionTable(decompositionMethod, dictionary);
45 }
46 
47 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
48 
50 (
51  const dictionary& decompositionDict
52 )
53 {
54  word methodType(decompositionDict.lookup("method"));
55 
56  if (methodType == "scotch" && Pstream::parRun())
57  {
58  methodType = "ptscotch";
59  }
60 
61 
62  Info<< "Selecting decompositionMethod " << methodType << endl;
63 
64  dictionaryConstructorTable::iterator cstrIter =
65  dictionaryConstructorTablePtr_->find(methodType);
66 
67  if (cstrIter == dictionaryConstructorTablePtr_->end())
68  {
70  (
71  "decompositionMethod::New"
72  "(const dictionary& decompositionDict)"
73  ) << "Unknown decompositionMethod "
74  << methodType << nl << nl
75  << "Valid decompositionMethods are : " << endl
76  << dictionaryConstructorTablePtr_->sortedToc()
77  << exit(FatalError);
78  }
79 
80  return autoPtr<decompositionMethod>(cstrIter()(decompositionDict));
81 }
82 
83 
85 (
86  const polyMesh& mesh,
87  const pointField& points
88 )
89 {
90  scalarField weights(points.size(), 1.0);
91 
92  return decompose(mesh, points, weights);
93 }
94 
95 
97 (
98  const polyMesh& mesh,
99  const labelList& fineToCoarse,
100  const pointField& coarsePoints,
101  const scalarField& coarseWeights
102 )
103 {
104  CompactListList<label> coarseCellCells;
105  calcCellCells
106  (
107  mesh,
108  fineToCoarse,
109  coarsePoints.size(),
110  true, // use global cell labels
111  coarseCellCells
112  );
113 
114  // Decompose based on agglomerated points
115  labelList coarseDistribution
116  (
117  decompose
118  (
119  coarseCellCells(),
120  coarsePoints,
121  coarseWeights
122  )
123  );
124 
125  // Rework back into decomposition for original mesh_
126  labelList fineDistribution(fineToCoarse.size());
127 
128  forAll(fineDistribution, i)
129  {
130  fineDistribution[i] = coarseDistribution[fineToCoarse[i]];
131  }
132 
133  return fineDistribution;
134 }
135 
136 
138 (
139  const polyMesh& mesh,
140  const labelList& fineToCoarse,
141  const pointField& coarsePoints
142 )
143 {
144  scalarField cWeights(coarsePoints.size(), 1.0);
145 
146  return decompose
147  (
148  mesh,
149  fineToCoarse,
150  coarsePoints,
151  cWeights
152  );
153 }
154 
155 
157 (
158  const labelListList& globalCellCells,
159  const pointField& cc
160 )
161 {
162  scalarField cWeights(cc.size(), 1.0);
163 
164  return decompose(globalCellCells, cc, cWeights);
165 }
166 
167 
169 (
170  const polyMesh& mesh,
171  const labelList& agglom,
172  const label nLocalCoarse,
173  const bool parallel,
174  CompactListList<label>& cellCells
175 )
176 {
177  const labelList& faceOwner = mesh.faceOwner();
178  const labelList& faceNeighbour = mesh.faceNeighbour();
179  const polyBoundaryMesh& patches = mesh.boundaryMesh();
180 
181 
182  // Create global cell numbers
183  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
184 
185  globalIndex globalAgglom
186  (
187  nLocalCoarse,
190  parallel
191  );
192 
193 
194  // Get agglomerate owner on other side of coupled faces
195  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
196 
197  labelList globalNeighbour(mesh.nFaces()-mesh.nInternalFaces());
198 
199  forAll(patches, patchI)
200  {
201  const polyPatch& pp = patches[patchI];
202 
203  if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
204  {
205  label faceI = pp.start();
206  label bFaceI = pp.start() - mesh.nInternalFaces();
207 
208  forAll(pp, i)
209  {
210  globalNeighbour[bFaceI] = globalAgglom.toGlobal
211  (
212  agglom[faceOwner[faceI]]
213  );
214 
215  bFaceI++;
216  faceI++;
217  }
218  }
219  }
220 
221  // Get the cell on the other side of coupled patches
222  syncTools::swapBoundaryFaceList(mesh, globalNeighbour);
223 
224 
225  // Count number of faces (internal + coupled)
226  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
227 
228  // Number of faces per coarse cell
229  labelList nFacesPerCell(nLocalCoarse, 0);
230 
231  for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
232  {
233  label own = agglom[faceOwner[faceI]];
234  label nei = agglom[faceNeighbour[faceI]];
235 
236  nFacesPerCell[own]++;
237  nFacesPerCell[nei]++;
238  }
239 
240  forAll(patches, patchI)
241  {
242  const polyPatch& pp = patches[patchI];
243 
244  if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
245  {
246  label faceI = pp.start();
247  label bFaceI = pp.start()-mesh.nInternalFaces();
248 
249  forAll(pp, i)
250  {
251  label own = agglom[faceOwner[faceI]];
252 
253  label globalNei = globalNeighbour[bFaceI];
254  if
255  (
256  !globalAgglom.isLocal(globalNei)
257  || globalAgglom.toLocal(globalNei) != own
258  )
259  {
260  nFacesPerCell[own]++;
261  }
262 
263  faceI++;
264  bFaceI++;
265  }
266  }
267  }
268 
269 
270  // Fill in offset and data
271  // ~~~~~~~~~~~~~~~~~~~~~~~
272 
273  cellCells.setSize(nFacesPerCell);
274 
275  nFacesPerCell = 0;
276 
277  labelList& m = cellCells.m();
278  const labelList& offsets = cellCells.offsets();
279 
280  // For internal faces is just offsetted owner and neighbour
281  for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
282  {
283  label own = agglom[faceOwner[faceI]];
284  label nei = agglom[faceNeighbour[faceI]];
285 
286  m[offsets[own] + nFacesPerCell[own]++] = globalAgglom.toGlobal(nei);
287  m[offsets[nei] + nFacesPerCell[nei]++] = globalAgglom.toGlobal(own);
288  }
289 
290  // For boundary faces is offsetted coupled neighbour
291  forAll(patches, patchI)
292  {
293  const polyPatch& pp = patches[patchI];
294 
295  if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
296  {
297  label faceI = pp.start();
298  label bFaceI = pp.start()-mesh.nInternalFaces();
299 
300  forAll(pp, i)
301  {
302  label own = agglom[faceOwner[faceI]];
303 
304  label globalNei = globalNeighbour[bFaceI];
305 
306  if
307  (
308  !globalAgglom.isLocal(globalNei)
309  || globalAgglom.toLocal(globalNei) != own
310  )
311  {
312  m[offsets[own] + nFacesPerCell[own]++] = globalNei;
313  }
314 
315  faceI++;
316  bFaceI++;
317  }
318  }
319  }
320 
321 
322  // Check for duplicates connections between cells
323  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
324  // Done as postprocessing step since we now have cellCells.
325  label newIndex = 0;
326  labelHashSet nbrCells;
327 
328 
329  if (cellCells.size() == 0)
330  {
331  return;
332  }
333 
334  label startIndex = cellCells.offsets()[0];
335 
336  forAll(cellCells, cellI)
337  {
338  nbrCells.clear();
339  nbrCells.insert(globalAgglom.toGlobal(cellI));
340 
341  label endIndex = cellCells.offsets()[cellI+1];
342 
343  for (label i = startIndex; i < endIndex; i++)
344  {
345  if (nbrCells.insert(cellCells.m()[i]))
346  {
347  cellCells.m()[newIndex++] = cellCells.m()[i];
348  }
349  }
350  startIndex = endIndex;
351  cellCells.offsets()[cellI+1] = newIndex;
352  }
353 
354  cellCells.m().setSize(newIndex);
355 
356  //forAll(cellCells, cellI)
357  //{
358  // Pout<< "Original: Coarse cell " << cellI << endl;
359  // forAll(mesh.cellCells()[cellI], i)
360  // {
361  // Pout<< " nbr:" << mesh.cellCells()[cellI][i] << endl;
362  // }
363  // Pout<< "Compacted: Coarse cell " << cellI << endl;
364  // const labelUList cCells = cellCells[cellI];
365  // forAll(cCells, i)
366  // {
367  // Pout<< " nbr:" << cCells[i] << endl;
368  // }
369  //}
370 }
371 
372 
373 //void Foam::decompositionMethod::calcCellCells
374 //(
375 // const polyMesh& mesh,
376 // const boolList& blockedFace,
377 // const List<labelPair>& explicitConnections,
378 // const labelList& agglom,
379 // const label nLocalCoarse,
380 // const bool parallel,
381 // CompactListList<label>& cellCells
382 //)
383 //{
384 // const labelList& faceOwner = mesh.faceOwner();
385 // const labelList& faceNeighbour = mesh.faceNeighbour();
386 // const polyBoundaryMesh& patches = mesh.boundaryMesh();
387 //
388 //
389 // // Create global cell numbers
390 // // ~~~~~~~~~~~~~~~~~~~~~~~~~~
391 //
392 // globalIndex globalAgglom
393 // (
394 // nLocalCoarse,
395 // Pstream::msgType(),
396 // Pstream::worldComm,
397 // parallel
398 // );
399 //
400 //
401 // // Get agglomerate owner on other side of coupled faces
402 // // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
403 //
404 // labelList globalNeighbour(mesh.nFaces()-mesh.nInternalFaces());
405 //
406 // forAll(patches, patchI)
407 // {
408 // const polyPatch& pp = patches[patchI];
409 //
410 // if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
411 // {
412 // label faceI = pp.start();
413 // label bFaceI = pp.start() - mesh.nInternalFaces();
414 //
415 // forAll(pp, i)
416 // {
417 // globalNeighbour[bFaceI] = globalAgglom.toGlobal
418 // (
419 // agglom[faceOwner[faceI]]
420 // );
421 //
422 // bFaceI++;
423 // faceI++;
424 // }
425 // }
426 // }
427 //
428 // // Get the cell on the other side of coupled patches
429 // syncTools::swapBoundaryFaceList(mesh, globalNeighbour);
430 //
431 //
432 // // Count number of faces (internal + coupled)
433 // // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
434 //
435 // // Number of faces per coarse cell
436 // labelList nFacesPerCell(nLocalCoarse, 0);
437 //
438 // // 1. Internal faces
439 // for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
440 // {
441 // if (!blockedFace[faceI])
442 // {
443 // label own = agglom[faceOwner[faceI]];
444 // label nei = agglom[faceNeighbour[faceI]];
445 //
446 // nFacesPerCell[own]++;
447 // nFacesPerCell[nei]++;
448 // }
449 // }
450 //
451 // // 2. Coupled faces
452 // forAll(patches, patchI)
453 // {
454 // const polyPatch& pp = patches[patchI];
455 //
456 // if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
457 // {
458 // label faceI = pp.start();
459 // label bFaceI = pp.start()-mesh.nInternalFaces();
460 //
461 // forAll(pp, i)
462 // {
463 // if (!blockedFace[faceI])
464 // {
465 // label own = agglom[faceOwner[faceI]];
466 //
467 // label globalNei = globalNeighbour[bFaceI];
468 // if
469 // (
470 // !globalAgglom.isLocal(globalNei)
471 // || globalAgglom.toLocal(globalNei) != own
472 // )
473 // {
474 // nFacesPerCell[own]++;
475 // }
476 //
477 // faceI++;
478 // bFaceI++;
479 // }
480 // }
481 // }
482 // }
483 //
484 // // 3. Explicit connections between non-coupled boundary faces
485 // forAll(explicitConnections, i)
486 // {
487 // const labelPair& baffle = explicitConnections[i];
488 // label f0 = baffle.first();
489 // label f1 = baffle.second();
490 //
491 // if (!blockedFace[f0] && blockedFace[f1])
492 // {
493 // label f0Own = agglom[faceOwner[f0]];
494 // label f1Own = agglom[faceOwner[f1]];
495 //
496 // // Always count the connection between the two owner sides
497 // if (f0Own != f1Own)
498 // {
499 // nFacesPerCell[f0Own]++;
500 // nFacesPerCell[f1Own]++;
501 // }
502 //
503 // // Add any neighbour side connections
504 // if (mesh.isInternalFace(f0))
505 // {
506 // label f0Nei = agglom[faceNeighbour[f0]];
507 //
508 // if (mesh.isInternalFace(f1))
509 // {
510 // // Internal faces
511 // label f1Nei = agglom[faceNeighbour[f1]];
512 //
513 // if (f0Own != f1Nei)
514 // {
515 // nFacesPerCell[f0Own]++;
516 // nFacesPerCell[f1Nei]++;
517 // }
518 // if (f0Nei != f1Own)
519 // {
520 // nFacesPerCell[f0Nei]++;
521 // nFacesPerCell[f1Own]++;
522 // }
523 // if (f0Nei != f1Nei)
524 // {
525 // nFacesPerCell[f0Nei]++;
526 // nFacesPerCell[f1Nei]++;
527 // }
528 // }
529 // else
530 // {
531 // // f1 boundary face
532 // if (f0Nei != f1Own)
533 // {
534 // nFacesPerCell[f0Nei]++;
535 // nFacesPerCell[f1Own]++;
536 // }
537 // }
538 // }
539 // else
540 // {
541 // if (mesh.isInternalFace(f1))
542 // {
543 // label f1Nei = agglom[faceNeighbour[f1]];
544 // if (f0Own != f1Nei)
545 // {
546 // nFacesPerCell[f0Own]++;
547 // nFacesPerCell[f1Nei]++;
548 // }
549 // }
550 // }
551 // }
552 // }
553 //
554 //
555 // // Fill in offset and data
556 // // ~~~~~~~~~~~~~~~~~~~~~~~
557 //
558 // cellCells.setSize(nFacesPerCell);
559 //
560 // nFacesPerCell = 0;
561 //
562 // labelList& m = cellCells.m();
563 // const labelList& offsets = cellCells.offsets();
564 //
565 // // 1. For internal faces is just offsetted owner and neighbour
566 // for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
567 // {
568 // if (!blockedFace[faceI])
569 // {
570 // label own = agglom[faceOwner[faceI]];
571 // label nei = agglom[faceNeighbour[faceI]];
572 //
573 // m[offsets[own] + nFacesPerCell[own]++] =
574 // globalAgglom.toGlobal(nei);
575 // m[offsets[nei] + nFacesPerCell[nei]++] =
576 // globalAgglom.toGlobal(own);
577 // }
578 // }
579 //
580 // // 2. For boundary faces is offsetted coupled neighbour
581 // forAll(patches, patchI)
582 // {
583 // const polyPatch& pp = patches[patchI];
584 //
585 // if (pp.coupled() && (parallel || !isA<processorPolyPatch>(pp)))
586 // {
587 // label faceI = pp.start();
588 // label bFaceI = pp.start()-mesh.nInternalFaces();
589 //
590 // forAll(pp, i)
591 // {
592 // if (!blockedFace[faceI])
593 // {
594 // label own = agglom[faceOwner[faceI]];
595 //
596 // label globalNei = globalNeighbour[bFaceI];
597 //
598 // if
599 // (
600 // !globalAgglom.isLocal(globalNei)
601 // || globalAgglom.toLocal(globalNei) != own
602 // )
603 // {
604 // m[offsets[own] + nFacesPerCell[own]++] = globalNei;
605 // }
606 //
607 // faceI++;
608 // bFaceI++;
609 // }
610 // }
611 // }
612 // }
613 //
614 // // 3. Explicit connections between non-coupled boundary faces
615 // forAll(explicitConnections, i)
616 // {
617 // const labelPair& baffle = explicitConnections[i];
618 // label f0 = baffle.first();
619 // label f1 = baffle.second();
620 //
621 // if (!blockedFace[f0] && blockedFace[f1])
622 // {
623 // label f0Own = agglom[faceOwner[f0]];
624 // label f1Own = agglom[faceOwner[f1]];
625 //
626 // // Always count the connection between the two owner sides
627 // if (f0Own != f1Own)
628 // {
629 // m[offsets[f0Own] + nFacesPerCell[f0Own]++] =
630 // globalAgglom.toGlobal(f1Own);
631 // m[offsets[f1Own] + nFacesPerCell[f1Own]++] =
632 // globalAgglom.toGlobal(f0Own);
633 // }
634 //
635 // // Add any neighbour side connections
636 // if (mesh.isInternalFace(f0))
637 // {
638 // label f0Nei = agglom[faceNeighbour[f0]];
639 //
640 // if (mesh.isInternalFace(f1))
641 // {
642 // // Internal faces
643 // label f1Nei = agglom[faceNeighbour[f1]];
644 //
645 // if (f0Own != f1Nei)
646 // {
647 // m[offsets[f0Own] + nFacesPerCell[f0Own]++] =
648 // globalAgglom.toGlobal(f1Nei);
649 // m[offsets[f1Nei] + nFacesPerCell[f1Nei]++] =
650 // globalAgglom.toGlobal(f1Nei);
651 // }
652 // if (f0Nei != f1Own)
653 // {
654 // m[offsets[f0Nei] + nFacesPerCell[f0Nei]++] =
655 // globalAgglom.toGlobal(f1Own);
656 // m[offsets[f1Own] + nFacesPerCell[f1Own]++] =
657 // globalAgglom.toGlobal(f0Nei);
658 // }
659 // if (f0Nei != f1Nei)
660 // {
661 // m[offsets[f0Nei] + nFacesPerCell[f0Nei]++] =
662 // globalAgglom.toGlobal(f1Nei);
663 // m[offsets[f1Nei] + nFacesPerCell[f1Nei]++] =
664 // globalAgglom.toGlobal(f0Nei);
665 // }
666 // }
667 // else
668 // {
669 // // f1 boundary face
670 // if (f0Nei != f1Own)
671 // {
672 // m[offsets[f0Nei] + nFacesPerCell[f0Nei]++] =
673 // globalAgglom.toGlobal(f1Own);
674 // m[offsets[f1Own] + nFacesPerCell[f1Own]++] =
675 // globalAgglom.toGlobal(f0Nei);
676 // }
677 // }
678 // }
679 // else
680 // {
681 // if (mesh.isInternalFace(f1))
682 // {
683 // label f1Nei = agglom[faceNeighbour[f1]];
684 // if (f0Own != f1Nei)
685 // {
686 // m[offsets[f0Own] + nFacesPerCell[f0Own]++] =
687 // globalAgglom.toGlobal(f1Nei);
688 // m[offsets[f1Nei] + nFacesPerCell[f1Nei]++] =
689 // globalAgglom.toGlobal(f0Own);
690 // }
691 // }
692 // }
693 // }
694 // }
695 //
696 //
697 // // Check for duplicates connections between cells
698 // // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
699 // // Done as postprocessing step since we now have cellCells.
700 // label newIndex = 0;
701 // labelHashSet nbrCells;
702 //
703 //
704 // if (cellCells.size() == 0)
705 // {
706 // return;
707 // }
708 //
709 // label startIndex = cellCells.offsets()[0];
710 //
711 // forAll(cellCells, cellI)
712 // {
713 // nbrCells.clear();
714 // nbrCells.insert(globalAgglom.toGlobal(cellI));
715 //
716 // label endIndex = cellCells.offsets()[cellI+1];
717 //
718 // for (label i = startIndex; i < endIndex; i++)
719 // {
720 // if (nbrCells.insert(cellCells.m()[i]))
721 // {
722 // cellCells.m()[newIndex++] = cellCells.m()[i];
723 // }
724 // }
725 // startIndex = endIndex;
726 // cellCells.offsets()[cellI+1] = newIndex;
727 // }
728 //
729 // cellCells.m().setSize(newIndex);
730 //
731 // //forAll(cellCells, cellI)
732 // //{
733 // // Pout<< "Original: Coarse cell " << cellI << endl;
734 // // forAll(mesh.cellCells()[cellI], i)
735 // // {
736 // // Pout<< " nbr:" << mesh.cellCells()[cellI][i] << endl;
737 // // }
738 // // Pout<< "Compacted: Coarse cell " << cellI << endl;
739 // // const labelUList cCells = cellCells[cellI];
740 // // forAll(cCells, i)
741 // // {
742 // // Pout<< " nbr:" << cCells[i] << endl;
743 // // }
744 // //}
745 //}
746 
747 
749 (
750  const polyMesh& mesh,
751  const scalarField& cellWeights,
752 
753  //- Whether owner and neighbour should be on same processor
754  // (takes priority over explicitConnections)
755  const boolList& blockedFace,
756 
757  //- Whether whole sets of faces (and point neighbours) need to be kept
758  // on single processor
759  const PtrList<labelList>& specifiedProcessorFaces,
760  const labelList& specifiedProcessor,
761 
762  //- Additional connections between boundary faces
763  const List<labelPair>& explicitConnections
764 )
765 {
766  // Any weights specified?
767  label nWeights = returnReduce(cellWeights.size(), sumOp<label>());
768 
769  if (nWeights > 0 && cellWeights.size() != mesh.nCells())
770  {
772  (
773  "decompositionMethod::decompose\n"
774  "(\n"
775  " const polyMesh&,\n"
776  " const scalarField&,\n"
777  " const boolList&,\n"
778  " const PtrList<labelList>&,\n"
779  " const labelList&,\n"
780  " const List<labelPair>&\n"
781  ) << "Number of weights " << cellWeights.size()
782  << " differs from number of cells " << mesh.nCells()
783  << exit(FatalError);
784  }
785 
786 
787  // Any processor sets?
788  label nProcSets = 0;
789  forAll(specifiedProcessorFaces, setI)
790  {
791  nProcSets += specifiedProcessorFaces[setI].size();
792  }
793  reduce(nProcSets, sumOp<label>());
794 
795  // Any non-mesh connections?
796  label nConnections = returnReduce
797  (
798  explicitConnections.size(),
799  sumOp<label>()
800  );
801 
802  // Any faces not blocked?
803  label nUnblocked = 0;
804  forAll(blockedFace, faceI)
805  {
806  if (!blockedFace[faceI])
807  {
808  nUnblocked++;
809  }
810  }
811  reduce(nUnblocked, sumOp<label>());
812 
813 
814 
815  // Either do decomposition on cell centres or on agglomeration
816 
817  labelList finalDecomp;
818 
819 
820  if (nProcSets+nConnections+nUnblocked == 0)
821  {
822  // No constraints, possibly weights
823 
824  if (nWeights > 0)
825  {
826  finalDecomp = decompose
827  (
828  mesh,
829  mesh.cellCentres(),
830  cellWeights
831  );
832  }
833  else
834  {
835  finalDecomp = decompose(mesh, mesh.cellCentres());
836  }
837  }
838  else
839  {
840  if (debug)
841  {
842  Info<< "Constrained decomposition:" << endl
843  << " faces with same owner and neighbour processor : "
844  << nUnblocked << endl
845  << " baffle faces with same owner processor : "
846  << nConnections << endl
847  << " faces all on same processor : "
848  << nProcSets << endl << endl;
849  }
850 
851  // Determine local regions, separated by blockedFaces
852  regionSplit localRegion(mesh, blockedFace, explicitConnections, false);
853 
854 
855  if (debug)
856  {
857  Info<< "Constrained decomposition:" << endl
858  << " split into " << localRegion.nLocalRegions()
859  << " regions."
860  << endl;
861  }
862 
863  // Determine region cell centres
864  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
865 
866  // This just takes the first cell in the region. Otherwise the problem
867  // is with cyclics - if we'd average the region centre might be
868  // somewhere in the middle of the domain which might not be anywhere
869  // near any of the cells.
870 
871  pointField regionCentres(localRegion.nLocalRegions(), point::max);
872 
873  forAll(localRegion, cellI)
874  {
875  label regionI = localRegion[cellI];
876 
877  if (regionCentres[regionI] == point::max)
878  {
879  regionCentres[regionI] = mesh.cellCentres()[cellI];
880  }
881  }
882 
883  // Do decomposition on agglomeration
884  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
885 
886  scalarField regionWeights(localRegion.nLocalRegions(), 0);
887 
888  if (nWeights > 0)
889  {
890  forAll(localRegion, cellI)
891  {
892  label regionI = localRegion[cellI];
893 
894  regionWeights[regionI] += cellWeights[cellI];
895  }
896  }
897  else
898  {
899  forAll(localRegion, cellI)
900  {
901  label regionI = localRegion[cellI];
902 
903  regionWeights[regionI] += 1.0;
904  }
905  }
906 
907  finalDecomp = decompose
908  (
909  mesh,
910  localRegion,
911  regionCentres,
912  regionWeights
913  );
914 
915 
916 
917  // Implement the explicitConnections since above decompose
918  // does not know about them
919  forAll(explicitConnections, i)
920  {
921  const labelPair& baffle = explicitConnections[i];
922  label f0 = baffle.first();
923  label f1 = baffle.second();
924 
925  if (!blockedFace[f0] && !blockedFace[f1])
926  {
927  // Note: what if internal faces and owner and neighbour on
928  // different processor? So for now just push owner side
929  // proc
930 
931  const label procI = finalDecomp[mesh.faceOwner()[f0]];
932 
933  finalDecomp[mesh.faceOwner()[f1]] = procI;
934  if (mesh.isInternalFace(f1))
935  {
936  finalDecomp[mesh.faceNeighbour()[f1]] = procI;
937  }
938  }
939  else if (blockedFace[f0] != blockedFace[f1])
940  {
942  (
943  "labelList decompose\n"
944  "(\n"
945  " const polyMesh&,\n"
946  " const scalarField&,\n"
947  " const boolList&,\n"
948  " const PtrList<labelList>&,\n"
949  " const labelList&,\n"
950  " const List<labelPair>&\n"
951  ")"
952  ) << "On explicit connection between faces " << f0
953  << " and " << f1
954  << " the two blockedFace status are not equal : "
955  << blockedFace[f0] << " and " << blockedFace[f1]
956  << exit(FatalError);
957  }
958  }
959 
960 
961  // blockedFaces corresponding to processor faces need to be handled
962  // separately since not handled by local regionSplit. We need to
963  // walk now across coupled faces and make sure to move a whole
964  // global region across
965  if (Pstream::parRun())
966  {
967  // Re-do regionSplit
968 
969  // Field on cells and faces.
970  List<minData> cellData(mesh.nCells());
971  List<minData> faceData(mesh.nFaces());
972 
973  // Take over blockedFaces by seeding a negative number
974  // (so is always less than the decomposition)
975  label nUnblocked = 0;
976  forAll(blockedFace, faceI)
977  {
978  if (blockedFace[faceI])
979  {
980  faceData[faceI] = minData(-123);
981  }
982  else
983  {
984  nUnblocked++;
985  }
986  }
987 
988  // Seed unblocked faces with destination processor
989  labelList seedFaces(nUnblocked);
990  List<minData> seedData(nUnblocked);
991  nUnblocked = 0;
992 
993  forAll(blockedFace, faceI)
994  {
995  if (!blockedFace[faceI])
996  {
997  label own = mesh.faceOwner()[faceI];
998  seedFaces[nUnblocked] = faceI;
999  seedData[nUnblocked] = minData(finalDecomp[own]);
1000  nUnblocked++;
1001  }
1002  }
1003 
1004 
1005  // Propagate information inwards
1006  FaceCellWave<minData> deltaCalc
1007  (
1008  mesh,
1009  seedFaces,
1010  seedData,
1011  faceData,
1012  cellData,
1013  mesh.globalData().nTotalCells()+1
1014  );
1015 
1016  // And extract
1017  forAll(finalDecomp, cellI)
1018  {
1019  if (cellData[cellI].valid(deltaCalc.data()))
1020  {
1021  finalDecomp[cellI] = cellData[cellI].data();
1022  }
1023  }
1024  }
1025 
1026 
1027  // For specifiedProcessorFaces rework the cellToProc to enforce
1028  // all on one processor since we can't guarantee that the input
1029  // to regionSplit was a single region.
1030  // E.g. faceSet 'a' with the cells split into two regions
1031  // by a notch formed by two walls
1032  //
1033  // \ /
1034  // \ /
1035  // ---a----+-----a-----
1036  //
1037  //
1038  // Note that reworking the cellToProc might make the decomposition
1039  // unbalanced.
1040  forAll(specifiedProcessorFaces, setI)
1041  {
1042  const labelList& set = specifiedProcessorFaces[setI];
1043 
1044  label procI = specifiedProcessor[setI];
1045  if (procI == -1)
1046  {
1047  // If no processor specified use the one from the
1048  // 0th element
1049  procI = finalDecomp[mesh.faceOwner()[set[0]]];
1050  }
1051 
1052  forAll(set, fI)
1053  {
1054  const face& f = mesh.faces()[set[fI]];
1055  forAll(f, fp)
1056  {
1057  const labelList& pFaces = mesh.pointFaces()[f[fp]];
1058  forAll(pFaces, i)
1059  {
1060  label faceI = pFaces[i];
1061 
1062  finalDecomp[mesh.faceOwner()[faceI]] = procI;
1063  if (mesh.isInternalFace(faceI))
1064  {
1065  finalDecomp[mesh.faceNeighbour()[faceI]] = procI;
1066  }
1067  }
1068  }
1069  }
1070  }
1071 
1072 
1073  if (debug && Pstream::parRun())
1074  {
1075  labelList nbrDecomp;
1076  syncTools::swapBoundaryCellList(mesh, finalDecomp, nbrDecomp);
1077 
1078  const polyBoundaryMesh& patches = mesh.boundaryMesh();
1079  forAll(patches, patchI)
1080  {
1081  const polyPatch& pp = patches[patchI];
1082  if (pp.coupled())
1083  {
1084  forAll(pp, i)
1085  {
1086  label faceI = pp.start()+i;
1087  label own = mesh.faceOwner()[faceI];
1088  label bFaceI = faceI-mesh.nInternalFaces();
1089 
1090  if (!blockedFace[faceI])
1091  {
1092  label ownProc = finalDecomp[own];
1093  label nbrProc = nbrDecomp[bFaceI];
1094  if (ownProc != nbrProc)
1095  {
1096  FatalErrorIn("decompositionMethod::decompose()")
1097  << "patch:" << pp.name()
1098  << " face:" << faceI
1099  << " at:" << mesh.faceCentres()[faceI]
1100  << " ownProc:" << ownProc
1101  << " nbrProc:" << nbrProc
1102  << exit(FatalError);
1103  }
1104  }
1105  }
1106  }
1107  }
1108  }
1109  }
1110 
1111  return finalDecomp;
1112 }
1113 
1114 
1117  const polyMesh& mesh,
1118  boolList& blockedFace,
1119  PtrList<labelList>& specifiedProcessorFaces,
1120  labelList& specifiedProcessor,
1121  List<labelPair>& explicitConnections
1122 )
1123 {
1124  blockedFace.setSize(mesh.nFaces());
1125  blockedFace = true;
1126  //label nUnblocked = 0;
1127 
1128  specifiedProcessorFaces.clear();
1129  explicitConnections.clear();
1130 
1131 
1132  if (decompositionDict_.found("preservePatches"))
1133  {
1134  wordList pNames(decompositionDict_.lookup("preservePatches"));
1135 
1136  Info<< nl
1137  << "Keeping owner of faces in patches " << pNames
1138  << " on same processor. This only makes sense for cyclics." << endl;
1139 
1140  const polyBoundaryMesh& patches = mesh.boundaryMesh();
1141 
1142  forAll(pNames, i)
1143  {
1144  const label patchI = patches.findPatchID(pNames[i]);
1145 
1146  if (patchI == -1)
1147  {
1148  FatalErrorIn("decompositionMethod::decompose(const polyMesh&)")
1149  << "Unknown preservePatch " << pNames[i]
1150  << endl << "Valid patches are " << patches.names()
1151  << exit(FatalError);
1152  }
1153 
1154  const polyPatch& pp = patches[patchI];
1155 
1156  forAll(pp, i)
1157  {
1158  if (blockedFace[pp.start() + i])
1159  {
1160  blockedFace[pp.start() + i] = false;
1161  //nUnblocked++;
1162  }
1163  }
1164  }
1165  }
1166  if (decompositionDict_.found("preserveFaceZones"))
1167  {
1168  wordList zNames(decompositionDict_.lookup("preserveFaceZones"));
1169 
1170  Info<< nl
1171  << "Keeping owner and neighbour of faces in zones " << zNames
1172  << " on same processor" << endl;
1173 
1174  const faceZoneMesh& fZones = mesh.faceZones();
1175 
1176  forAll(zNames, i)
1177  {
1178  label zoneI = fZones.findZoneID(zNames[i]);
1179 
1180  if (zoneI == -1)
1181  {
1182  FatalErrorIn("decompositionMethod::decompose(const polyMesh&)")
1183  << "Unknown preserveFaceZone " << zNames[i]
1184  << endl << "Valid faceZones are " << fZones.names()
1185  << exit(FatalError);
1186  }
1187 
1188  const faceZone& fz = fZones[zoneI];
1189 
1190  forAll(fz, i)
1191  {
1192  if (blockedFace[fz[i]])
1193  {
1194  blockedFace[fz[i]] = false;
1195  //nUnblocked++;
1196  }
1197  }
1198  }
1199  }
1200 
1201  bool preserveBaffles = decompositionDict_.lookupOrDefault
1202  (
1203  "preserveBaffles",
1204  false
1205  );
1206  if (preserveBaffles)
1207  {
1208  Info<< nl
1209  << "Keeping owner of faces in baffles "
1210  << " on same processor." << endl;
1211 
1212  explicitConnections = localPointRegion::findDuplicateFacePairs(mesh);
1213  forAll(explicitConnections, i)
1214  {
1215  blockedFace[explicitConnections[i].first()] = false;
1216  blockedFace[explicitConnections[i].second()] = false;
1217  }
1218  }
1219 
1220  if
1221  (
1222  decompositionDict_.found("preservePatches")
1223  || decompositionDict_.found("preserveFaceZones")
1224  || preserveBaffles
1225  )
1226  {
1227  syncTools::syncFaceList(mesh, blockedFace, andEqOp<bool>());
1228  //reduce(nUnblocked, sumOp<label>());
1229  }
1230 
1231 
1232 
1233  // Specified processor for group of cells connected to faces
1234 
1235  label nProcSets = 0;
1236  if (decompositionDict_.found("singleProcessorFaceSets"))
1237  {
1238  List<Tuple2<word, label> > zNameAndProcs
1239  (
1240  decompositionDict_.lookup("singleProcessorFaceSets")
1241  );
1242 
1243  specifiedProcessorFaces.setSize(zNameAndProcs.size());
1244  specifiedProcessor.setSize(zNameAndProcs.size());
1245 
1246  forAll(zNameAndProcs, setI)
1247  {
1248  Info<< "Keeping all cells connected to faceSet "
1249  << zNameAndProcs[setI].first()
1250  << " on processor " << zNameAndProcs[setI].second() << endl;
1251 
1252  // Read faceSet
1253  faceSet fz(mesh, zNameAndProcs[setI].first());
1254 
1255  specifiedProcessorFaces.set(setI, new labelList(fz.sortedToc()));
1256  specifiedProcessor[setI] = zNameAndProcs[setI].second();
1257  nProcSets += fz.size();
1258  }
1259  reduce(nProcSets, sumOp<label>());
1260 
1261 
1262  // Unblock all point connected faces
1263  // 1. Mark all points on specifiedProcessorFaces
1264  boolList procFacePoint(mesh.nPoints(), false);
1265  forAll(specifiedProcessorFaces, setI)
1266  {
1267  const labelList& set = specifiedProcessorFaces[setI];
1268  forAll(set, fI)
1269  {
1270  const face& f = mesh.faces()[set[fI]];
1271  forAll(f, fp)
1272  {
1273  procFacePoint[f[fp]] = true;
1274  }
1275  }
1276  }
1277  syncTools::syncPointList(mesh, procFacePoint, orEqOp<bool>(), false);
1278 
1279  // 2. Unblock all faces on procFacePoint
1280  forAll(procFacePoint, pointI)
1281  {
1282  if (procFacePoint[pointI])
1283  {
1284  const labelList& pFaces = mesh.pointFaces()[pointI];
1285  forAll(pFaces, i)
1286  {
1287  blockedFace[pFaces[i]] = false;
1288  }
1289  }
1290  }
1291  syncTools::syncFaceList(mesh, blockedFace, andEqOp<bool>());
1292  }
1293 }
1294 
1295 
1298  const polyMesh& mesh,
1299  const scalarField& cellWeights
1300 )
1301 {
1302  boolList blockedFace;
1303  PtrList<labelList> specifiedProcessorFaces;
1304  labelList specifiedProcessor;
1305  List<labelPair> explicitConnections;
1306  setConstraints
1307  (
1308  mesh,
1309  blockedFace,
1310  specifiedProcessorFaces,
1311  specifiedProcessor,
1312  explicitConnections
1313  );
1314 
1315 
1316  // Construct decomposition method and either do decomposition on
1317  // cell centres or on agglomeration
1318 
1319  labelList finalDecomp = decompose
1320  (
1321  mesh,
1322  cellWeights, // optional weights
1323  blockedFace, // any cells to be combined
1324  specifiedProcessorFaces,// any whole cluster of cells to be kept
1325  specifiedProcessor,
1326  explicitConnections // baffles
1327  );
1328 
1329  return finalDecomp;
1330 }
1331 
1332 
1333 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:380
const faceZoneMesh & faceZones() const
Return face zone mesh.
Definition: polyMesh.H:463
bool isLocal(const label i) const
Is on local processor.
Definition: globalIndexI.H:95
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:142
label size() const
Return the primary size, i.e. the number of rows.
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
label nFaces() const
bool set(const label) const
Is element set.
Definition: PtrListI.H:107
A list of face labels.
Definition: faceSet.H:48
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
Definition: polyPatch.H:310
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
const word & name() const
Return name.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
labelList f(nPoints)
void clear()
Clear the PtrList, i.e. set size to zero deleting all the.
Definition: PtrList.C:185
virtual labelList decompose(const pointField &points, const scalarField &pointWeights)
Return for every coordinate the wanted processor number.
List< Key > sortedToc() const
Return the table of contents as a sorted list.
Definition: HashTable.C:216
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:64
scalar f1
Definition: createFields.H:28
const Type & second() const
Return second.
Definition: Pair.H:99
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
A class for handling words, derived from string.
Definition: word.H:59
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const List< label > & offsets() const
Return the offset table (= size()+1)
label nCells() const
static const Vector max
Definition: Vector.H:82
const vectorField & cellCentres() const
void setConstraints(const polyMesh &mesh, boolList &blockedFace, PtrList< labelList > &specifiedProcessorFaces, labelList &specifiedProcessor, List< labelPair > &explicitConnections)
Helper: extract constraints:
messageStream Info
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
static void syncFaceList(const polyMesh &mesh, UList< T > &l, const CombineOp &cop)
Synchronize values on all mesh faces.
Definition: syncTools.H:380
T & first()
Return the first element of the list.
Definition: UListI.H:117
Namespace for OpenFOAM.
T * data()
Return a pointer to the first data element,.
Definition: UListI.H:152
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
patches[0]
void setSize(const label nRows)
Reset size of CompactListList.
static labelPairList findDuplicateFacePairs(const polyMesh &)
Helper routine to find all baffles (two boundary faces.
wordList names() const
Return a list of patch names.
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
void clear()
Clear the list, i.e. set size to zero.
Definition: List.C:379
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Wave propagation of information through grid. Every iteration information goes through one layer of c...
Definition: FaceCellWave.H:74
void clear()
Clear all entries from table.
Definition: HashTable.C:473
Foam::polyBoundaryMesh.
#define forAll(list, i)
Definition: UList.H:421
const List< T > & m() const
Return the packed matrix of data.
label nLocalRegions() const
Return local number of regions.
Definition: regionSplit.H:212
label nTotalCells() const
Return total number of cells in decomposed mesh.
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:451
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled boundary face values.
Definition: syncTools.H:429
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
static autoPtr< decompositionMethod > New(const dictionary &decompositionDict)
Return a reference to the selected decomposition method.
const labelListList & pointFaces() const
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1073
static void syncPointList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh points.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
static void calcCellCells(const polyMesh &mesh, const labelList &agglom, const label nLocalCoarse, const bool global, CompactListList< label > &cellCells)
Helper: determine (local or global) cellCells from mesh.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
A packed storage unstructured matrix of objects of type <T> using an offset table for access...
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
label nInternalFaces() const
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
label toGlobal(const label i) const
From local to global.
Definition: globalIndexI.H:82
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
For use with FaceCellWave. Transports minimum passive data.
Definition: minData.H:52
label toLocal(const label i) const
From global to local on current processor.
Definition: globalIndexI.H:117
List< label > labelList
A List of labels.
Definition: labelList.H:56
const Type & first() const
Return first.
Definition: Pair.H:87
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: ZoneMesh.C:354
static void swapBoundaryCellList(const polyMesh &mesh, const UList< T > &cellData, List< T > &neighbourCellData)
Swap to obtain neighbour cell values for all boundary faces.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1060
wordList names() const
Return a list of zone names.
Definition: ZoneMesh.C:269
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints,-1);const cellModel &hex=*(cellModeller::lookup("hex"));labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){pointMap[i]=i;}for(label i=0;i< nPoints;i++){if(f[i] > 0.0){hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells]=cellShape(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei){if(edges[ei].mag(points)< SMALL){label start=pointMap[edges[ei].start()];while(start!=pointMap[start]){start=pointMap[start];}label end=pointMap[edges[ei].end()];while(end!=pointMap[end]){end=pointMap[end];}label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;}}cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){cellShape &cs=cellShapes[celli];forAll(cs, i){cs[i]=pointMap[cs[i]];}cs.collapse();}label bcIDs[11]={-1, 0, 2, 4,-1, 5,-1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={&wallPolyPatch::typeName,&wallPolyPatch::typeName,&wallPolyPatch::typeName,&wallPolyPatch::typeName,&symmetryPolyPatch::typeName,&wedgePolyPatch::typeName,&polyPatch::typeName,&polyPatch::typeName,&polyPatch::typeName,&polyPatch::typeName,&symmetryPolyPatch::typeName,&oldCyclicPolyPatch::typeName};enum patchTypeNames{PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={"piston","valve","liner","cylinderHead","axis","wedge","inflow","outflow","presin","presout","symmetryPlane","cyclic"};List< SLList< face > > pFaces[nBCs]
Definition: readKivaGrid.H:235
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1079
label findPatchID(const word &patchName) const
Find patch index given a name.
static label worldComm
Default communicator (all processors)
Definition: UPstream.H:261
This class separates the mesh into distinct unconnected regions, each of which is then given a label ...
Definition: regionSplit.H:116
label nPoints() const
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1200
defineTypeNameAndDebug(combustionModel, 0)
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
const vectorField & faceCentres() const
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:116