polyBoundaryMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 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 "polyBoundaryMesh.H"
27 #include "polyMesh.H"
28 #include "primitiveMesh.H"
29 #include "processorPolyPatch.H"
30 #include "stringListOps.H"
31 #include "PstreamBuffers.H"
32 #include "lduSchedule.H"
33 #include "globalMeshData.H"
34 #include "stringListOps.H"
35 
36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40  defineTypeNameAndDebug(polyBoundaryMesh, 0);
41 }
42 
43 
44 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45 
46 Foam::polyBoundaryMesh::polyBoundaryMesh
47 (
48  const IOobject& io,
49  const polyMesh& mesh
50 )
51 :
52  polyPatchList(),
53  regIOobject(io),
54  mesh_(mesh)
55 {
56  if
57  (
58  readOpt() == IOobject::MUST_READ
59  || readOpt() == IOobject::MUST_READ_IF_MODIFIED
60  )
61  {
62  if (readOpt() == IOobject::MUST_READ_IF_MODIFIED)
63  {
65  << "Specified IOobject::MUST_READ_IF_MODIFIED but class"
66  << " does not support automatic rereading."
67  << endl;
68  }
69 
70 
71  polyPatchList& patches = *this;
72 
73  // Read polyPatchList
74  Istream& is = readStream(typeName);
75 
76  PtrList<entry> patchEntries(is);
77  patches.setSize(patchEntries.size());
78 
79  forAll(patches, patchi)
80  {
81  patches.set
82  (
83  patchi,
85  (
86  patchEntries[patchi].keyword(),
87  patchEntries[patchi].dict(),
88  patchi,
89  *this
90  )
91  );
92  }
93 
94  // Check state of IOstream
95  is.check
96  (
97  "polyBoundaryMesh::polyBoundaryMesh"
98  "(const IOobject&, const polyMesh&)"
99  );
100 
101  close();
102  }
103 }
104 
105 
106 Foam::polyBoundaryMesh::polyBoundaryMesh
107 (
108  const IOobject& io,
109  const polyMesh& pm,
110  const label size
111 )
112 :
113  polyPatchList(size),
114  regIOobject(io),
115  mesh_(pm)
116 {}
117 
118 
119 Foam::polyBoundaryMesh::polyBoundaryMesh
120 (
121  const IOobject& io,
122  const polyMesh& pm,
123  const polyPatchList& ppl
124 )
125 :
126  polyPatchList(),
127  regIOobject(io),
128  mesh_(pm)
129 {
130  if
131  (
132  (this->readOpt() == IOobject::READ_IF_PRESENT && this->headerOk())
133  || this->readOpt() == IOobject::MUST_READ
134  || this->readOpt() == IOobject::MUST_READ_IF_MODIFIED
135  )
136  {
137 
138  if (readOpt() == IOobject::MUST_READ_IF_MODIFIED)
139  {
141  << "Specified IOobject::MUST_READ_IF_MODIFIED but class"
142  << " does not support automatic rereading."
143  << endl;
144  }
145 
146  polyPatchList& patches = *this;
147 
148  // Read polyPatchList
149  Istream& is = readStream(typeName);
150 
151  PtrList<entry> patchEntries(is);
152  patches.setSize(patchEntries.size());
153 
154  forAll(patches, patchi)
155  {
156  patches.set
157  (
158  patchi,
160  (
161  patchEntries[patchi].keyword(),
162  patchEntries[patchi].dict(),
163  patchi,
164  *this
165  )
166  );
167  }
168 
169  // Check state of IOstream
170  is.check
171  (
172  "polyBoundaryMesh::polyBoundaryMesh"
173  "(const IOobject&, const polyMesh&, const polyPatchList&)"
174  );
175 
176  close();
177  }
178  else
179  {
180  polyPatchList& patches = *this;
181  patches.setSize(ppl.size());
182  forAll(patches, patchi)
183  {
184  patches.set(patchi, ppl[patchi].clone(*this).ptr());
185  }
186  }
187 }
188 
189 
190 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
191 
193 {}
194 
195 
197 {
198  forAll(*this, patchi)
199  {
200  if (this->set(patchi))
201  {
202  operator[](patchi).clearGeom();
203  }
204  }
205 }
206 
207 
209 {
210  nbrEdgesPtr_.clear();
211  patchIDPtr_.clear();
212  patchFaceIDPtr_.clear();
213  groupPatchIDsPtr_.clear();
214 
215  forAll(*this, patchi)
216  {
217  if (this->set(patchi))
218  {
219  operator[](patchi).clearAddressing();
220  }
221  }
222 }
223 
224 
225 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
226 
227 void Foam::polyBoundaryMesh::calcGeometry()
228 {
230 
231  if
232  (
235  )
236  {
237  forAll(*this, patchi)
238  {
239  operator[](patchi).initCalcGeometry(pBufs);
240  }
241 
242  pBufs.finishedSends();
243 
244  forAll(*this, patchi)
245  {
246  operator[](patchi).calcGeometry(pBufs);
247  }
248  }
250  {
251  const lduSchedule& patchSchedule = mesh().globalData().patchSchedule();
252 
253  // Dummy.
254  pBufs.finishedSends();
255 
256  forAll(patchSchedule, patchEvali)
257  {
258  const label patchi = patchSchedule[patchEvali].patch;
259 
260  if (patchSchedule[patchEvali].init)
261  {
262  operator[](patchi).initCalcGeometry(pBufs);
263  }
264  else
265  {
266  operator[](patchi).calcGeometry(pBufs);
267  }
268  }
269  }
270 }
271 
272 
275 {
276  if (Pstream::parRun())
277  {
279  << "Neighbour edge addressing not correct across parallel"
280  << " boundaries." << endl;
281  }
282 
283  if (!nbrEdgesPtr_.valid())
284  {
285  nbrEdgesPtr_.reset(new List<labelPairList>(size()));
286  List<labelPairList>& nbrEdges = nbrEdgesPtr_();
287 
288  // Initialise.
289  label nEdgePairs = 0;
290  forAll(*this, patchi)
291  {
292  const polyPatch& pp = operator[](patchi);
293 
294  nbrEdges[patchi].setSize(pp.nEdges() - pp.nInternalEdges());
295 
296  forAll(nbrEdges[patchi], i)
297  {
298  labelPair& edgeInfo = nbrEdges[patchi][i];
299 
300  edgeInfo[0] = -1;
301  edgeInfo[1] = -1;
302  }
303 
304  nEdgePairs += pp.nEdges() - pp.nInternalEdges();
305  }
306 
307  // From mesh edge (expressed as a point pair so as not to construct
308  // point addressing) to patch + relative edge index.
309  HashTable<labelPair, edge, Hash<edge>> pointsToEdge(nEdgePairs);
310 
311  forAll(*this, patchi)
312  {
313  const polyPatch& pp = operator[](patchi);
314 
315  const edgeList& edges = pp.edges();
316 
317  for
318  (
319  label edgei = pp.nInternalEdges();
320  edgei < edges.size();
321  edgei++
322  )
323  {
324  // Edge in patch local points
325  const edge& e = edges[edgei];
326 
327  // Edge in mesh points.
328  edge meshEdge(pp.meshPoints()[e[0]], pp.meshPoints()[e[1]]);
329 
331  pointsToEdge.find(meshEdge);
332 
333  if (fnd == pointsToEdge.end())
334  {
335  // First occurrence of mesh edge. Store patch and my
336  // local index.
337  pointsToEdge.insert
338  (
339  meshEdge,
340  labelPair
341  (
342  patchi,
343  edgei - pp.nInternalEdges()
344  )
345  );
346  }
347  else
348  {
349  // Second occurrence. Store.
350  const labelPair& edgeInfo = fnd();
351 
352  nbrEdges[patchi][edgei - pp.nInternalEdges()] =
353  edgeInfo;
354 
355  nbrEdges[edgeInfo[0]][edgeInfo[1]]
356  = labelPair(patchi, edgei - pp.nInternalEdges());
357 
358  // Found all two occurrences of this edge so remove from
359  // hash to save space. Note that this will give lots of
360  // problems if the polyBoundaryMesh is multiply connected.
361  pointsToEdge.erase(meshEdge);
362  }
363  }
364  }
365 
366  if (pointsToEdge.size())
367  {
369  << "Not all boundary edges of patches match up." << nl
370  << "Is the outside of your mesh multiply connected?"
371  << abort(FatalError);
372  }
373 
374  forAll(*this, patchi)
375  {
376  const polyPatch& pp = operator[](patchi);
377 
378  const labelPairList& nbrEdgesp = nbrEdges[patchi];
379 
380  forAll(nbrEdgesp, i)
381  {
382  const labelPair& edgeInfo = nbrEdgesp[i];
383 
384  if (edgeInfo[0] == -1 || edgeInfo[1] == -1)
385  {
386  label edgeI = pp.nInternalEdges() + i;
387  const edge& e = pp.edges()[edgeI];
388 
390  << "Not all boundary edges of patches match up." << nl
391  << "Edge " << edgeI << " on patch " << pp.name()
392  << " end points " << pp.localPoints()[e[0]] << ' '
393  << pp.localPoints()[e[1]] << " is not matched to an"
394  << " edge on any other patch." << nl
395  << "Is the outside of your mesh multiply connected?"
396  << abort(FatalError);
397  }
398  }
399  }
400  }
401 
402  return nbrEdgesPtr_();
403 }
404 
405 
407 {
408  if (!patchIDPtr_.valid())
409  {
410  patchIDPtr_.reset
411  (
412  new labelList
413  (
414  mesh_.nFaces()
415  - mesh_.nInternalFaces()
416  )
417  );
418  labelList& patchID = patchIDPtr_();
419 
420  const polyBoundaryMesh& bm = *this;
421 
422  forAll(bm, patchi)
423  {
424  label bFacei = bm[patchi].start() - mesh_.nInternalFaces();
425  forAll(bm[patchi], i)
426  {
427  patchID[bFacei++] = patchi;
428  }
429  }
430  }
431  return patchIDPtr_();
432 }
433 
434 
436 {
437  if (!patchFaceIDPtr_.valid())
438  {
439  patchFaceIDPtr_.reset
440  (
441  new labelList
442  (
443  mesh_.nFaces()
444  - mesh_.nInternalFaces()
445  )
446  );
447  labelList& patchFaceID = patchFaceIDPtr_();
448 
449  const polyBoundaryMesh& bm = *this;
450 
451  forAll(bm, patchi)
452  {
453  label bFacei = bm[patchi].start() - mesh_.nInternalFaces();
454  forAll(bm[patchi], i)
455  {
456  patchFaceID[bFacei++] = i;
457  }
458  }
459  }
460  return patchFaceIDPtr_();
461 }
462 
463 
466 {
467  if (!groupPatchIDsPtr_.valid())
468  {
469  groupPatchIDsPtr_.reset(new HashTable<labelList, word>(10));
470  HashTable<labelList, word>& groupPatchIDs = groupPatchIDsPtr_();
471 
472  const polyBoundaryMesh& bm = *this;
473 
474  forAll(bm, patchi)
475  {
476  const wordList& groups = bm[patchi].inGroups();
477 
478  forAll(groups, i)
479  {
480  const word& name = groups[i];
481 
482  HashTable<labelList, word>::iterator iter = groupPatchIDs.find
483  (
484  name
485  );
486 
487  if (iter != groupPatchIDs.end())
488  {
489  iter().append(patchi);
490  }
491  else
492  {
493  groupPatchIDs.insert(name, labelList(1, patchi));
494  }
495  }
496  }
497 
498  // Remove patch names from patchGroups
499  forAll(bm, patchi)
500  {
501  if (groupPatchIDs.erase(bm[patchi].name()))
502  {
504  << "Removing patchGroup '" << bm[patchi].name()
505  << "' which clashes with patch " << patchi
506  << " of the same name."
507  << endl;
508  }
509  }
510  }
511 
512  return groupPatchIDsPtr_();
513 }
514 
515 
517 (
518  const word& groupName,
519  const labelList& patchIDs
520 )
521 {
522  groupPatchIDsPtr_.clear();
523 
524  polyPatchList& patches = *this;
525 
526  boolList donePatch(patches.size(), false);
527 
528  // Add to specified patches
529  forAll(patchIDs, i)
530  {
531  label patchi = patchIDs[i];
532  polyPatch& pp = patches[patchi];
533 
534  if (!pp.inGroup(groupName))
535  {
536  pp.inGroups().append(groupName);
537  }
538  donePatch[patchi] = true;
539  }
540 
541  // Remove from other patches
542  forAll(patches, patchi)
543  {
544  if (!donePatch[patchi])
545  {
546  polyPatch& pp = patches[patchi];
547 
548  label newI = 0;
549  if (pp.inGroup(groupName))
550  {
551  wordList& groups = pp.inGroups();
552 
553  forAll(groups, i)
554  {
555  if (groups[i] != groupName)
556  {
557  groups[newI++] = groups[i];
558  }
559  }
560  groups.setSize(newI);
561  }
562  }
563  }
564 }
565 
566 
568 {
569  const polyPatchList& patches = *this;
570 
571  wordList t(patches.size());
572 
573  forAll(patches, patchi)
574  {
575  t[patchi] = patches[patchi].name();
576  }
577 
578  return t;
579 }
580 
581 
583 {
584  const polyPatchList& patches = *this;
585 
586  wordList t(patches.size());
587 
588  forAll(patches, patchi)
589  {
590  t[patchi] = patches[patchi].type();
591  }
592 
593  return t;
594 }
595 
596 
598 {
599  const polyPatchList& patches = *this;
600 
601  wordList t(patches.size());
602 
603  forAll(patches, patchi)
604  {
605  t[patchi] = patches[patchi].physicalType();
606  }
607 
608  return t;
609 }
610 
611 
613 (
614  const wordRe& key,
615  const bool usePatchGroups
616 ) const
617 {
618  DynamicList<label> indices;
619 
620  if (!key.empty())
621  {
622  if (key.isPattern())
623  {
624  indices = findStrings(key, this->names());
625 
626  if (usePatchGroups && groupPatchIDs().size())
627  {
628  labelHashSet indexSet(indices);
629 
630  const wordList allGroupNames = groupPatchIDs().toc();
631  labelList groupIDs = findStrings(key, allGroupNames);
632  forAll(groupIDs, i)
633  {
634  const word& grpName = allGroupNames[groupIDs[i]];
635  const labelList& patchIDs = groupPatchIDs()[grpName];
636  forAll(patchIDs, j)
637  {
638  if (indexSet.insert(patchIDs[j]))
639  {
640  indices.append(patchIDs[j]);
641  }
642  }
643  }
644  }
645  }
646  else
647  {
648  // Literal string. Special version of above to avoid
649  // unnecessary memory allocations
650 
651  indices.setCapacity(1);
652  forAll(*this, i)
653  {
654  if (key == operator[](i).name())
655  {
656  indices.append(i);
657  break;
658  }
659  }
660 
661  if (usePatchGroups && groupPatchIDs().size())
662  {
664  groupPatchIDs().find(key);
665 
666  if (iter != groupPatchIDs().end())
667  {
668  labelHashSet indexSet(indices);
669 
670  const labelList& patchIDs = iter();
671  forAll(patchIDs, j)
672  {
673  if (indexSet.insert(patchIDs[j]))
674  {
675  indices.append(patchIDs[j]);
676  }
677  }
678  }
679  }
680  }
681  }
682 
683  return move(indices);
684 }
685 
686 
688 {
689  if (!key.empty())
690  {
691  if (key.isPattern())
692  {
693  labelList indices = this->findIndices(key);
694 
695  // return first element
696  if (!indices.empty())
697  {
698  return indices[0];
699  }
700  }
701  else
702  {
703  forAll(*this, i)
704  {
705  if (key == operator[](i).name())
706  {
707  return i;
708  }
709  }
710  }
711  }
712 
713  // not found
714  return -1;
715 }
716 
717 
719 {
720  const polyPatchList& patches = *this;
721 
722  forAll(patches, patchi)
723  {
724  if (patches[patchi].name() == patchName)
725  {
726  return patchi;
727  }
728  }
729 
730  // Patch not found
731  if (debug)
732  {
733  Pout<< "label polyBoundaryMesh::findPatchID(const word&) const"
734  << "Patch named " << patchName << " not found. "
735  << "List of available patch names: " << names() << endl;
736  }
737 
738  // Not found, return -1
739  return -1;
740 }
741 
742 
744 {
745  // Find out which patch the current face belongs to by comparing label
746  // with patch start labels.
747  // If the face is internal, return -1;
748  // if it is off the end of the list, abort
749  if (faceIndex < mesh().nInternalFaces())
750  {
751  return -1;
752  }
753  else if (faceIndex >= mesh().nFaces())
754  {
756  << "given label " << faceIndex
757  << " greater than the number of geometric faces " << mesh().nFaces()
758  << abort(FatalError);
759  }
760 
761 
762  forAll(*this, patchi)
763  {
764  const polyPatch& bp = operator[](patchi);
765 
766  if
767  (
768  faceIndex >= bp.start()
769  && faceIndex < bp.start() + bp.size()
770  )
771  {
772  return patchi;
773  }
774  }
775 
776  // If not in any of above, it is trouble!
778  << "Cannot find face " << faceIndex << " in any of the patches "
779  << names() << nl
780  << "It seems your patches are not consistent with the mesh :"
781  << " internalFaces:" << mesh().nInternalFaces()
782  << " total number of faces:" << mesh().nFaces()
783  << abort(FatalError);
784 
785  return -1;
786 }
787 
788 
790 (
791  const UList<wordRe>& patchNames,
792  const bool warnNotFound,
793  const bool usePatchGroups
794 ) const
795 {
796  const wordList allPatchNames(this->names());
797  labelHashSet ids(size());
798 
799  forAll(patchNames, i)
800  {
801  const wordRe& patchName = patchNames[i];
802 
803  // Treat the given patch names as wild-cards and search the set
804  // of all patch names for matches
805  labelList patchIDs = findStrings(patchName, allPatchNames);
806 
807  forAll(patchIDs, j)
808  {
809  ids.insert(patchIDs[j]);
810  }
811 
812  if (patchIDs.empty())
813  {
814  if (usePatchGroups)
815  {
816  const wordList allGroupNames = groupPatchIDs().toc();
817 
818  // Regard as group name
819  labelList groupIDs = findStrings(patchName, allGroupNames);
820 
821  forAll(groupIDs, i)
822  {
823  const word& name = allGroupNames[groupIDs[i]];
824  const labelList& extraPatchIDs = groupPatchIDs()[name];
825 
826  forAll(extraPatchIDs, extraI)
827  {
828  ids.insert(extraPatchIDs[extraI]);
829  }
830  }
831 
832  if (groupIDs.empty() && warnNotFound)
833  {
835  << "Cannot find any patch or group names matching "
836  << patchName
837  << endl;
838  }
839  }
840  else if (warnNotFound)
841  {
843  << "Cannot find any patch names matching " << patchName
844  << endl;
845  }
846  }
847  }
848 
849  return ids;
850 }
851 
852 
854 (
855  const labelUList& patchIDs,
856  wordList& groups,
857  labelHashSet& nonGroupPatches
858 ) const
859 {
860  // Current matched groups
861  DynamicList<word> matchedGroups(1);
862 
863  // Current set of unmatched patches
864  nonGroupPatches = labelHashSet(patchIDs);
865 
866  const HashTable<labelList, word>& groupPatchIDs = this->groupPatchIDs();
867  for
868  (
870  groupPatchIDs.begin();
871  iter != groupPatchIDs.end();
872  ++iter
873  )
874  {
875  // Store currently unmatched patches so we can restore
876  labelHashSet oldNonGroupPatches(nonGroupPatches);
877 
878  // Match by deleting patches in group from the current set and seeing
879  // if all have been deleted.
880  labelHashSet groupPatchSet(iter());
881 
882  label nMatch = nonGroupPatches.erase(groupPatchSet);
883 
884  if (nMatch == groupPatchSet.size())
885  {
886  matchedGroups.append(iter.key());
887  }
888  else if (nMatch != 0)
889  {
890  // No full match. Undo.
891  nonGroupPatches.transfer(oldNonGroupPatches);
892  }
893  }
894 
895  groups.transfer(matchedGroups);
896 }
897 
898 
899 bool Foam::polyBoundaryMesh::checkParallelSync(const bool report) const
900 {
901  if (!Pstream::parRun())
902  {
903  return false;
904  }
905 
906 
907  const polyBoundaryMesh& bm = *this;
908 
909  bool hasError = false;
910 
911  // Collect non-proc patches and check proc patches are last.
912  wordList names(bm.size());
913  wordList types(bm.size());
914 
915  label nonProci = 0;
916 
917  forAll(bm, patchi)
918  {
919  if (!isA<processorPolyPatch>(bm[patchi]))
920  {
921  if (nonProci != patchi)
922  {
923  // There is processor patch in between normal patches.
924  hasError = true;
925 
926  if (debug || report)
927  {
928  Pout<< " ***Problem with boundary patch " << patchi
929  << " named " << bm[patchi].name()
930  << " of type " << bm[patchi].type()
931  << ". The patch seems to be preceded by processor"
932  << " patches. This is can give problems."
933  << endl;
934  }
935  }
936  else
937  {
938  names[nonProci] = bm[patchi].name();
939  types[nonProci] = bm[patchi].type();
940  nonProci++;
941  }
942  }
943  }
944  names.setSize(nonProci);
945  types.setSize(nonProci);
946 
947  List<wordList> allNames(Pstream::nProcs());
948  allNames[Pstream::myProcNo()] = names;
949  Pstream::gatherList(allNames);
950  Pstream::scatterList(allNames);
951 
952  List<wordList> allTypes(Pstream::nProcs());
953  allTypes[Pstream::myProcNo()] = types;
954  Pstream::gatherList(allTypes);
955  Pstream::scatterList(allTypes);
956 
957  // Have every processor check but only master print error.
958 
959  for (label proci = 1; proci < allNames.size(); ++proci)
960  {
961  if
962  (
963  (allNames[proci] != allNames[0])
964  || (allTypes[proci] != allTypes[0])
965  )
966  {
967  hasError = true;
968 
969  if (debug || (report && Pstream::master()))
970  {
971  Info<< " ***Inconsistent patches across processors, "
972  "processor 0 has patch names:" << allNames[0]
973  << " patch types:" << allTypes[0]
974  << " processor " << proci << " has patch names:"
975  << allNames[proci]
976  << " patch types:" << allTypes[proci]
977  << endl;
978  }
979  }
980  }
981 
982  return hasError;
983 }
984 
985 
986 bool Foam::polyBoundaryMesh::checkDefinition(const bool report) const
987 {
988  label nextPatchStart = mesh().nInternalFaces();
989  const polyBoundaryMesh& bm = *this;
990 
991  bool hasError = false;
992 
993  HashSet<word> patchNames(2*size());
994 
995  forAll(bm, patchi)
996  {
997  if (bm[patchi].start() != nextPatchStart && !hasError)
998  {
999  hasError = true;
1000 
1001  Info<< " ****Problem with boundary patch " << patchi
1002  << " named " << bm[patchi].name()
1003  << " of type " << bm[patchi].type()
1004  << ". The patch should start on face no " << nextPatchStart
1005  << " and the patch specifies " << bm[patchi].start()
1006  << "." << endl
1007  << "Possibly consecutive patches have this same problem."
1008  << " Suppressing future warnings." << endl;
1009  }
1010 
1011  if (!patchNames.insert(bm[patchi].name()) && !hasError)
1012  {
1013  hasError = true;
1014 
1015  Info<< " ****Duplicate boundary patch " << patchi
1016  << " named " << bm[patchi].name()
1017  << " of type " << bm[patchi].type()
1018  << "." << endl
1019  << "Suppressing future warnings." << endl;
1020  }
1021 
1022  nextPatchStart += bm[patchi].size();
1023  }
1024 
1025  reduce(hasError, orOp<bool>());
1026 
1027  if (debug || report)
1028  {
1029  if (hasError)
1030  {
1031  Pout<< " ***Boundary definition is in error." << endl;
1032  }
1033  else
1034  {
1035  Info<< " Boundary definition OK." << endl;
1036  }
1037  }
1038 
1039  return hasError;
1040 }
1041 
1042 
1044 {
1046 
1047  if
1048  (
1051  )
1052  {
1053  forAll(*this, patchi)
1054  {
1055  operator[](patchi).initMovePoints(pBufs, p);
1056  }
1057 
1058  pBufs.finishedSends();
1059 
1060  forAll(*this, patchi)
1061  {
1062  operator[](patchi).movePoints(pBufs, p);
1063  }
1064  }
1066  {
1067  const lduSchedule& patchSchedule = mesh().globalData().patchSchedule();
1068 
1069  // Dummy.
1070  pBufs.finishedSends();
1071 
1072  forAll(patchSchedule, patchEvali)
1073  {
1074  const label patchi = patchSchedule[patchEvali].patch;
1075 
1076  if (patchSchedule[patchEvali].init)
1077  {
1078  operator[](patchi).initMovePoints(pBufs, p);
1079  }
1080  else
1081  {
1082  operator[](patchi).movePoints(pBufs, p);
1083  }
1084  }
1085  }
1086 }
1087 
1088 
1090 {
1091  nbrEdgesPtr_.clear();
1092  patchIDPtr_.clear();
1093  patchFaceIDPtr_.clear();
1094  groupPatchIDsPtr_.clear();
1095 
1097 
1098  if
1099  (
1102  )
1103  {
1104  forAll(*this, patchi)
1105  {
1106  operator[](patchi).initTopoChange(pBufs);
1107  }
1108 
1109  pBufs.finishedSends();
1110 
1111  forAll(*this, patchi)
1112  {
1113  operator[](patchi).topoChange(pBufs);
1114  }
1115  }
1117  {
1118  const lduSchedule& patchSchedule = mesh().globalData().patchSchedule();
1119 
1120  // Dummy.
1121  pBufs.finishedSends();
1122 
1123  forAll(patchSchedule, patchEvali)
1124  {
1125  const label patchi = patchSchedule[patchEvali].patch;
1126 
1127  if (patchSchedule[patchEvali].init)
1128  {
1129  operator[](patchi).initTopoChange(pBufs);
1130  }
1131  else
1132  {
1133  operator[](patchi).topoChange(pBufs);
1134  }
1135  }
1136  }
1137 }
1138 
1139 
1142  const wordUList& newNames,
1143  const bool validBoundary
1144 )
1145 {
1146  polyPatchList& patches = *this;
1147  forAll(patches, patchi)
1148  {
1149  if (patches.set(patchi))
1150  {
1151  patches[patchi].rename(newNames);
1152  }
1153  }
1154 
1155  if (validBoundary)
1156  {
1157  topoChange();
1158  }
1159 }
1160 
1161 
1164  const labelUList& newToOld,
1165  const bool validBoundary
1166 )
1167 {
1168  // Change order of patches
1169  polyPatchList::shuffle(newToOld);
1170 
1171  // Adapt indices
1172  polyPatchList& patches = *this;
1173  forAll(patches, patchi)
1174  {
1175  if (patches.set(patchi))
1176  {
1177  patches[patchi].reorder(newToOld);
1178  }
1179  }
1180 
1181  if (validBoundary)
1182  {
1183  topoChange();
1184  }
1185 }
1186 
1187 
1189 {
1190  const polyPatchList& patches = *this;
1191 
1192  os << patches.size() << nl << token::BEGIN_LIST << incrIndent << nl;
1193 
1194  forAll(patches, patchi)
1195  {
1196  os << indent << patches[patchi].name() << nl
1197  << indent << token::BEGIN_BLOCK << nl
1198  << incrIndent << patches[patchi] << decrIndent
1199  << indent << token::END_BLOCK << endl;
1200  }
1201 
1202  os << decrIndent << token::END_LIST;
1203 
1204  // Check state of IOstream
1205  os.check("polyBoundaryMesh::writeData(Ostream& os) const");
1206 
1207  return os.good();
1208 }
1209 
1210 
1216  const bool write
1217 ) const
1218 {
1219  return regIOobject::writeObject(fmt, ver, IOstream::UNCOMPRESSED, write);
1220 }
1221 
1222 
1223 // * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
1224 
1227  const word& patchName
1228 ) const
1229 {
1230  const label patchi = findPatchID(patchName);
1231 
1232  if (patchi < 0)
1233  {
1235  << "Patch named " << patchName << " not found." << nl
1236  << "Available patch names: " << names() << endl
1237  << abort(FatalError);
1238  }
1239 
1240  return operator[](patchi);
1241 }
1242 
1243 
1246  const word& patchName
1247 )
1248 {
1249  const label patchi = findPatchID(patchName);
1250 
1251  if (patchi < 0)
1252  {
1254  << "Patch named " << patchName << " not found." << nl
1255  << "Available patch names: " << names() << endl
1256  << abort(FatalError);
1257  }
1258 
1259  return operator[](patchi);
1260 }
1261 
1262 
1263 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
1264 
1266 {
1267  pbm.writeData(os);
1268  return os;
1269 }
1270 
1271 
1272 // ************************************************************************* //
const fvPatchList & patches
PtrList< polyPatch > polyPatchList
container classes for polyPatch
Definition: polyPatchList.H:45
void reorder(const labelUList &oldToNew)
Reorders elements. Ordering does not have to be done in.
Definition: PtrList.C:197
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:325
bool writeData(Ostream &) const
writeData member function required by regIOobject
const word & name() const
Return name.
const word & name() const
Return name.
Definition: IOobject.H:315
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
void finishedSends(const bool block=true)
Mark all sends as having been done. This will start receives.
An STL-conforming const_iterator.
Definition: HashTable.H:481
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:112
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
label nInternalFaces() const
bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool write=true) const
Write using given format, version and form uncompression.
void renamePatches(const wordUList &newNames, const bool validBoundary)
Rename the patches. If validBoundary is set this calls topoChange()
label nFaces() const
bool inGroup(const word &) const
Test if in group.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
const lduSchedule & patchSchedule() const
Order in which the patches should be initialised/evaluated.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const labelList & patchID() const
Per boundary face label the patch index.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
void reorderPatches(const labelUList &newToOld, const bool validBoundary)
Reorders the patches. Ordering does not have to be done in.
label findPatchID(const word &patchName) const
Find patch index given a name.
label nInternalEdges() const
Number of internal edges.
void clearGeom()
Clear geometry at this level and at patches.
const wordList & inGroups() const
Return the optional groups patch belongs to.
Operations on lists of strings.
void transfer(HashTable< T, Key, Hash > &)
Transfer the contents of the argument table into this table.
Definition: HashTable.C:513
const Field< PointType > & localPoints() const
Return pointField of points in patch.
wordList types() const
Return a list of patch types.
labelList findIndices(const ListType &, typename ListType::const_reference, const label start=0)
Find all occurrences of given element. Linear search.
fvMesh & mesh
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:111
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:330
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
bool erase(const iterator &)
Erase a hashedEntry specified by given iterator.
Definition: HashTable.C:371
bool findStrings(const wordReListMatcher &matcher, const std::string &str)
Return true if string matches one of the regular expressions.
Definition: stringListOps.H:52
An STL iterator.
Definition: UPtrList.H:196
friend Ostream & operator(Ostream &, const UPtrList< T > &)
Write UPtrList to Ostream.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:211
An STL-conforming iterator.
Definition: HashTable.H:426
void matchGroups(const labelUList &patchIDs, wordList &groups, labelHashSet &nonGroupPatches) const
Match the patches to groups. Returns all the (fully matched) groups.
T clone(const T &t)
Definition: List.H:54
void topoChange()
Correct polyBoundaryMesh after topology update.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
~polyBoundaryMesh()
Destructor.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
void movePoints(const pointField &)
Correct polyBoundaryMesh after moving points.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
A class for handling words, derived from string.
Definition: word.H:59
static bool isPattern(const string &)
Test string for regular expression meta characters.
Definition: wordReI.H:34
bool checkDefinition(const bool report=false) const
Check boundary definition. Return true if in error.
wordList names() const
Return a list of patch names.
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
void setCapacity(const label)
Alter the size of the underlying storage.
Definition: DynamicListI.H:130
wordList patchNames(nPatches)
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:86
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
Pair< label > labelPair
Label pair.
Definition: labelPair.H:48
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1435
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
List< label > labelList
A List of labels.
Definition: labelList.H:56
wordList physicalTypes() const
Return a list of physical types.
iterator begin()
Iterator set to the beginning of the HashTable.
Definition: HashTableI.H:411
An STL-conforming hash table.
Definition: HashTable.H:61
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const List< labelPairList > & nbrEdges() const
Per patch the edges on the neighbouring patch. Is for every external.
A wordRe is a word, but can also have a regular expression for matching words.
Definition: wordRe.H:74
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
Foam::polyBoundaryMesh.
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:193
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
bool checkParallelSync(const bool report=false) const
Check whether all procs have all patches and in same order. Return.
label nEdges() const
Return number of edges in patch.
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:235
Buffers for inter-processor communications streams (UOPstream, UIPstream).
const labelList & patchFaceID() const
Per boundary face label the patch face index.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
edge meshEdge(const PrimitivePatch< FaceList, PointField > &p, const label edgei)
static commsTypes defaultCommsType
Default commsType.
Definition: UPstream.H:272
void setSize(const label)
Reset size of List.
Definition: List.C:281
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
void clearAddressing()
Clear addressing at this level and at patches.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:306
Ostream & operator<<(Ostream &, const ensightPart &)
Version number type.
Definition: IOstream.H:96
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:52
messageStream Info
virtual bool writeObject(IOstream::streamFormat, IOstream::versionNumber, IOstream::compressionType, const bool write) const
Write using given format, version and compression.
label findIndex(const wordRe &) const
Return patch index for the first match, return -1 if not found.
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
void shuffle(const labelUList &newToOld)
Reorders elements. Ordering does not have to be done in.
Definition: PtrList.C:248
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
volScalarField & p
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:228
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
labelList findIndices(const wordRe &, const bool usePatchGroups=true) const
Return patch indices for all matches. Optionally matches patchGroups.
void setGroup(const word &groupName, const labelList &patchIDs)
Set/add group with patches.
Namespace for OpenFOAM.
const HashTable< labelList, word > & groupPatchIDs() const
Per patch group the patch indices.