polyBoundaryMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2018-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "polyBoundaryMesh.H"
30 #include "polyMesh.H"
31 #include "primitiveMesh.H"
32 #include "processorPolyPatch.H"
33 #include "PstreamBuffers.H"
34 #include "lduSchedule.H"
35 #include "globalMeshData.H"
36 #include "stringListOps.H"
37 #include "DynamicList.H"
38 #include "PtrListOps.H"
39 #include "edgeHashes.H"
40 
41 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45  defineTypeNameAndDebug(polyBoundaryMesh, 0);
46 }
47 
48 
49 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 bool Foam::polyBoundaryMesh::hasGroupIDs() const
52 {
53  if (groupIDsPtr_)
54  {
55  // Use existing cache
56  return !groupIDsPtr_->empty();
57  }
58 
59  const polyPatchList& patches = *this;
60 
61  for (const polyPatch& p : patches)
62  {
63  if (!p.inGroups().empty())
64  {
65  return true;
66  }
67  }
68 
69  return false;
70 }
71 
72 
73 void Foam::polyBoundaryMesh::calcGroupIDs() const
74 {
75  if (groupIDsPtr_)
76  {
77  return; // Or FatalError
78  }
79 
80  groupIDsPtr_.reset(new HashTable<labelList>(16));
81  auto& groupLookup = *groupIDsPtr_;
82 
83  const polyPatchList& patches = *this;
84 
85  forAll(patches, patchi)
86  {
87  const wordList& groups = patches[patchi].inGroups();
88 
89  for (const word& groupName : groups)
90  {
91  groupLookup(groupName).append(patchi);
92  }
93  }
94 
95  // Remove groups that clash with patch names
96  forAll(patches, patchi)
97  {
98  if (groupLookup.erase(patches[patchi].name()))
99  {
101  << "Removed group '" << patches[patchi].name()
102  << "' which clashes with patch " << patchi
103  << " of the same name."
104  << endl;
105  }
106  }
107 }
108 
109 
110 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
111 
112 Foam::polyBoundaryMesh::polyBoundaryMesh
113 (
114  const IOobject& io,
115  const polyMesh& mesh
116 )
117 :
118  polyPatchList(),
119  regIOobject(io),
120  mesh_(mesh)
121 {
122  if
123  (
124  readOpt() == IOobject::MUST_READ
125  || readOpt() == IOobject::MUST_READ_IF_MODIFIED
126  )
127  {
128  // Warn for MUST_READ_IF_MODIFIED
129  warnNoRereading<polyBoundaryMesh>();
130 
131  polyPatchList& patches = *this;
132 
133  // Read polyPatchList
134  Istream& is = readStream(typeName);
135 
136  PtrList<entry> patchEntries(is);
137  patches.setSize(patchEntries.size());
138 
139  forAll(patches, patchi)
140  {
141  patches.set
142  (
143  patchi,
145  (
146  patchEntries[patchi].keyword(),
147  patchEntries[patchi].dict(),
148  patchi,
149  *this
150  )
151  );
152  }
153 
154  is.check(FUNCTION_NAME);
155 
156  close();
157  }
158 }
159 
160 
161 Foam::polyBoundaryMesh::polyBoundaryMesh
162 (
163  const IOobject& io,
164  const polyMesh& pm,
165  const label size
166 )
167 :
168  polyPatchList(size),
169  regIOobject(io),
170  mesh_(pm)
171 {}
172 
173 
174 Foam::polyBoundaryMesh::polyBoundaryMesh
175 (
176  const IOobject& io,
177  const polyMesh& pm,
178  const polyPatchList& ppl
179 )
180 :
181  polyPatchList(),
182  regIOobject(io),
183  mesh_(pm)
184 {
185  if
186  (
187  (this->readOpt() == IOobject::READ_IF_PRESENT && this->headerOk())
188  || this->readOpt() == IOobject::MUST_READ
189  || this->readOpt() == IOobject::MUST_READ_IF_MODIFIED
190  )
191  {
192  // Warn for MUST_READ_IF_MODIFIED
193  warnNoRereading<polyBoundaryMesh>();
194 
195  polyPatchList& patches = *this;
196 
197  // Read polyPatchList
198  Istream& is = readStream(typeName);
199 
200  PtrList<entry> patchEntries(is);
201  patches.setSize(patchEntries.size());
202 
203  forAll(patches, patchi)
204  {
205  patches.set
206  (
207  patchi,
209  (
210  patchEntries[patchi].keyword(),
211  patchEntries[patchi].dict(),
212  patchi,
213  *this
214  )
215  );
216  }
217 
218  is.check(FUNCTION_NAME);
219 
220  close();
221  }
222  else
223  {
224  polyPatchList& patches = *this;
225  patches.setSize(ppl.size());
226  forAll(patches, patchi)
227  {
228  patches.set(patchi, ppl[patchi].clone(*this).ptr());
229  }
230  }
231 }
232 
233 
234 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
235 
237 {
238  polyPatchList& patches = *this;
239 
240  for (polyPatch& p : patches)
241  {
242  p.clearGeom();
243  }
244 }
245 
246 
248 {
249  neighbourEdgesPtr_.clear();
250  patchIDPtr_.clear();
251  groupIDsPtr_.clear();
252 
253  polyPatchList& patches = *this;
254 
255  for (polyPatch& p : patches)
256  {
257  p.clearAddressing();
258  }
259 }
260 
261 
262 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
263 
264 void Foam::polyBoundaryMesh::calcGeometry()
265 {
267 
268  if
269  (
272  )
273  {
274  forAll(*this, patchi)
275  {
276  operator[](patchi).initGeometry(pBufs);
277  }
278 
279  pBufs.finishedSends();
280 
281  forAll(*this, patchi)
282  {
283  operator[](patchi).calcGeometry(pBufs);
284  }
285  }
287  {
288  const lduSchedule& patchSchedule = mesh().globalData().patchSchedule();
289 
290  // Dummy.
291  pBufs.finishedSends();
292 
293  for (const auto& patchEval : patchSchedule)
294  {
295  const label patchi = patchEval.patch;
296 
297  if (patchEval.init)
298  {
299  operator[](patchi).initGeometry(pBufs);
300  }
301  else
302  {
303  operator[](patchi).calcGeometry(pBufs);
304  }
305  }
306  }
307 }
308 
309 
312 {
313  const polyPatchList& patches = *this;
314 
315  UPtrList<const labelUList> list(patches.size());
316 
317  forAll(patches, patchi)
318  {
319  list.set(patchi, &patches[patchi].faceCells());
320  }
321 
322  return list;
323 }
324 
325 
328 {
329  if (Pstream::parRun())
330  {
332  << "Neighbour edge addressing not correct across parallel"
333  << " boundaries." << endl;
334  }
335 
336  if (!neighbourEdgesPtr_)
337  {
338  neighbourEdgesPtr_.reset(new List<labelPairList>(size()));
339  auto& neighbourEdges = *neighbourEdgesPtr_;
340 
341  // Initialize.
342  label nEdgePairs = 0;
343  forAll(*this, patchi)
344  {
345  const polyPatch& pp = operator[](patchi);
346 
347  neighbourEdges[patchi].setSize(pp.nEdges() - pp.nInternalEdges());
348 
349  for (labelPair& edgeInfo : neighbourEdges[patchi])
350  {
351  edgeInfo[0] = -1;
352  edgeInfo[1] = -1;
353  }
354 
355  nEdgePairs += pp.nEdges() - pp.nInternalEdges();
356  }
357 
358  // From mesh edge (expressed as a point pair so as not to construct
359  // point addressing) to patch + relative edge index.
360  EdgeMap<labelPair> pointsToEdge(nEdgePairs);
361 
362  forAll(*this, patchi)
363  {
364  const polyPatch& pp = operator[](patchi);
365 
366  const edgeList& edges = pp.edges();
367 
368  for
369  (
370  label edgei = pp.nInternalEdges();
371  edgei < edges.size();
372  edgei++
373  )
374  {
375  // Edge in patch local points
376  const edge& e = edges[edgei];
377 
378  // Edge in mesh points.
379  edge meshEdge(pp.meshPoints()[e[0]], pp.meshPoints()[e[1]]);
380 
381  auto fnd = pointsToEdge.find(meshEdge);
382 
383  if (!fnd.found())
384  {
385  // First occurrence of mesh edge. Store patch and my
386  // local index.
387  pointsToEdge.insert
388  (
389  meshEdge,
390  labelPair
391  (
392  patchi,
393  edgei - pp.nInternalEdges()
394  )
395  );
396  }
397  else
398  {
399  // Second occurrence. Store.
400  const labelPair& edgeInfo = fnd.val();
401 
402  neighbourEdges[patchi][edgei - pp.nInternalEdges()] =
403  edgeInfo;
404 
405  neighbourEdges[edgeInfo[0]][edgeInfo[1]]
406  = labelPair(patchi, edgei - pp.nInternalEdges());
407 
408  // Found all two occurrences of this edge so remove from
409  // hash to save space. Note that this will give lots of
410  // problems if the polyBoundaryMesh is multiply connected.
411  pointsToEdge.erase(meshEdge);
412  }
413  }
414  }
415 
416  if (pointsToEdge.size())
417  {
419  << "Not all boundary edges of patches match up." << nl
420  << "Is the outside of your mesh multiply connected?"
421  << abort(FatalError);
422  }
423 
424  forAll(*this, patchi)
425  {
426  const polyPatch& pp = operator[](patchi);
427 
428  const labelPairList& nbrEdges = neighbourEdges[patchi];
429 
430  forAll(nbrEdges, i)
431  {
432  const labelPair& edgeInfo = nbrEdges[i];
433 
434  if (edgeInfo[0] == -1 || edgeInfo[1] == -1)
435  {
436  const label edgei = pp.nInternalEdges() + i;
437  const edge& e = pp.edges()[edgei];
438 
440  << "Not all boundary edges of patches match up." << nl
441  << "Edge " << edgei << " on patch " << pp.name()
442  << " end points " << pp.localPoints()[e[0]] << ' '
443  << pp.localPoints()[e[1]] << " is not matched to an"
444  << " edge on any other patch." << nl
445  << "Is the outside of your mesh multiply connected?"
446  << abort(FatalError);
447  }
448  }
449  }
450  }
451 
452  return *neighbourEdgesPtr_;
453 }
454 
455 
457 {
458  if (!patchIDPtr_)
459  {
460  patchIDPtr_.reset(new labelList(mesh_.nBoundaryFaces()));
461  labelList& list = *patchIDPtr_;
462 
463  const polyPatchList& patches = *this;
464 
465  forAll(patches, patchi)
466  {
468  (
469  list,
470  patches[patchi].size(),
471  (patches[patchi].start() - mesh_.nInternalFaces())
472  ) = patchi;
473  }
474  }
475 
476  return *patchIDPtr_;
477 }
478 
479 
482 {
483  if (!groupIDsPtr_)
484  {
485  calcGroupIDs();
486  }
487 
488  return *groupIDsPtr_;
489 }
490 
491 
493 (
494  const word& groupName,
495  const labelUList& patchIDs
496 )
497 {
498  groupIDsPtr_.clear();
499 
500  polyPatchList& patches = *this;
501 
502  boolList donePatch(patches.size(), false);
503 
504  // Add to specified patches
505  for (const label patchi : patchIDs)
506  {
507  patches[patchi].inGroups().appendUniq(groupName);
508  donePatch[patchi] = true;
509  }
510 
511  // Remove from other patches
512  forAll(patches, patchi)
513  {
514  if (!donePatch[patchi])
515  {
516  wordList& groups = patches[patchi].inGroups();
517 
518  if (groups.found(groupName))
519  {
520  label newi = 0;
521  forAll(groups, i)
522  {
523  if (groups[i] != groupName)
524  {
525  groups[newi++] = groups[i];
526  }
527  }
528  groups.resize(newi);
529  }
530  }
531  }
532 }
533 
534 
536 {
537  const polyPatchList& patches = *this;
538 
539  label nonProc = 0;
540 
541  for (const polyPatch& p : patches)
542  {
543  if (isA<processorPolyPatch>(p))
544  {
545  break;
546  }
547 
548  ++nonProc;
549  }
550 
551  return nonProc;
552 }
553 
554 
556 {
557  return PtrListOps::get<word>(*this, nameOp<polyPatch>());
558 }
559 
560 
562 {
563  return PtrListOps::get<word>(*this, typeOp<polyPatch>());
564 }
565 
566 
568 {
569  return
570  PtrListOps::get<word>
571  (
572  *this,
573  [](const polyPatch& p) { return p.physicalType(); }
574  );
575 }
576 
577 
579 {
580  return
581  PtrListOps::get<label>
582  (
583  *this,
584  [](const polyPatch& p) { return p.start(); }
585  );
586 }
587 
588 
590 {
591  return
592  PtrListOps::get<label>
593  (
594  *this,
595  [](const polyPatch& p) { return p.size(); }
596  );
597 }
598 
599 
601 {
602  return
603  PtrListOps::get<labelRange>
604  (
605  *this,
606  [](const polyPatch& p) { return p.range(); }
607  );
608 }
609 
610 
611 Foam::label Foam::polyBoundaryMesh::start() const
612 {
613  return mesh_.nInternalFaces();
614 }
615 
616 
618 {
619  return mesh_.nBoundaryFaces();
620 }
621 
622 
624 {
625  return labelRange(mesh_.nInternalFaces(), mesh_.nBoundaryFaces());
626 }
627 
628 
630 {
631  if (patchi < 0)
632  {
633  return labelRange(mesh_.nInternalFaces(), 0);
634  }
635 
636  // Will fail if patchi >= size()
637  return (*this)[patchi].range();
638 }
639 
640 
642 (
643  const wordRe& matcher,
644  const bool useGroups
645 ) const
646 {
647  if (matcher.empty())
648  {
649  return labelList();
650  }
651 
652  // Only check groups if requested and they exist
653  const bool checkGroups = (useGroups && this->hasGroupIDs());
654 
655  labelHashSet ids;
656 
657  if (matcher.isPattern())
658  {
659  if (checkGroups)
660  {
661  const auto& groupLookup = groupPatchIDs();
662  forAllConstIters(groupLookup, iter)
663  {
664  if (matcher.match(iter.key()))
665  {
666  // Hash ids associated with the group
667  ids.insert(iter.val());
668  }
669  }
670  }
671 
672  if (ids.empty())
673  {
674  return PtrListOps::findMatching(*this, matcher);
675  }
676  else
677  {
678  ids.insert(PtrListOps::findMatching(*this, matcher));
679  }
680  }
681  else
682  {
683  // Literal string.
684  // Special version of above for reduced memory footprint.
685 
686  const label patchId = PtrListOps::firstMatching(*this, matcher);
687 
688  if (patchId >= 0)
689  {
690  return labelList(one{}, patchId);
691  }
692  else if (checkGroups)
693  {
694  const auto iter = groupPatchIDs().cfind(matcher);
695 
696  if (iter.found())
697  {
698  // Hash ids associated with the group
699  ids.insert(iter.val());
700  }
701  }
702  }
703 
704  return ids.sortedToc();
705 }
706 
707 
709 (
710  const wordRes& matcher,
711  const bool useGroups
712 ) const
713 {
714  if (matcher.empty())
715  {
716  return labelList();
717  }
718  else if (matcher.size() == 1)
719  {
720  return this->indices(matcher.first(), useGroups);
721  }
722 
723  labelHashSet ids;
724 
725  // Only check groups if requested and they exist
726  if (useGroups && this->hasGroupIDs())
727  {
728  ids.resize(2*this->size());
729 
730  const auto& groupLookup = groupPatchIDs();
731  forAllConstIters(groupLookup, iter)
732  {
733  if (matcher.match(iter.key()))
734  {
735  // Hash ids associated with the group
736  ids.insert(iter.val());
737  }
738  }
739  }
740 
741  if (ids.empty())
742  {
743  return PtrListOps::findMatching(*this, matcher);
744  }
745  else
746  {
747  ids.insert(PtrListOps::findMatching(*this, matcher));
748  }
749 
750  return ids.sortedToc();
751 }
752 
753 
755 {
756  if (key.empty())
757  {
758  return -1;
759  }
760  return PtrListOps::firstMatching(*this, key);
761 }
762 
763 
765 (
766  const word& patchName,
767  bool allowNotFound
768 ) const
769 {
770  if (patchName.empty())
771  {
772  return -1;
773  }
774 
775  const label patchId = PtrListOps::firstMatching(*this, patchName);
776 
777  if (patchId >= 0)
778  {
779  return patchId;
780  }
781 
782  if (!allowNotFound)
783  {
785  << "Patch '" << patchName << "' not found. "
786  << "Available patch names";
787 
788  if (polyMesh::defaultRegion != mesh_.name())
789  {
790  FatalError
791  << " in region '" << mesh_.name() << "'";
792  }
793 
794  FatalError
795  << " include: " << names() << endl
796  << exit(FatalError);
797  }
798 
799  // Patch not found
800  if (debug)
801  {
802  Pout<< "label polyBoundaryMesh::findPatchID(const word&) const"
803  << "Patch named " << patchName << " not found. "
804  << "Available patch names: " << names() << endl;
805  }
806 
807  // Not found, return -1
808  return -1;
809 }
810 
811 
812 Foam::label Foam::polyBoundaryMesh::whichPatch(const label faceIndex) const
813 {
814  // Find out which patch the current face belongs to by comparing label
815  // with patch start labels.
816  // If the face is internal, return -1;
817  // if it is off the end of the list, abort
818  if (faceIndex < mesh().nInternalFaces())
819  {
820  return -1;
821  }
822  else if (faceIndex >= mesh().nFaces())
823  {
825  << "Face " << faceIndex
826  << " out of bounds. Number of geometric faces " << mesh().nFaces()
827  << abort(FatalError);
828  }
829 
830 
831  // Patches are ordered, use binary search
832 
833  const polyPatchList& patches = *this;
834 
835  const label patchi =
836  findLower
837  (
838  patches,
839  faceIndex,
840  0,
841  // Must include the start in the comparison
842  [](const polyPatch& p, label val) { return (p.start() <= val); }
843  );
844 
845  if (patchi < 0 || !patches[patchi].range().found(faceIndex))
846  {
847  // If not in any of above, it is trouble!
849  << "Face " << faceIndex << " not found in any of the patches "
850  << flatOutput(names()) << nl
851  << "The patches appear to be inconsistent with the mesh :"
852  << " internalFaces:" << mesh().nInternalFaces()
853  << " total number of faces:" << mesh().nFaces()
854  << abort(FatalError);
855 
856  return -1;
857  }
858 
859  return patchi;
860 }
861 
862 
864 (
865  const UList<wordRe>& patchNames,
866  const bool warnNotFound,
867  const bool useGroups
868 ) const
869 {
870  const wordList allPatchNames(this->names());
871  labelHashSet ids(2*this->size());
872 
873  // Only check groups if requested and they exist
874  const bool checkGroups = (useGroups && this->hasGroupIDs());
875 
876  for (const wordRe& matcher : patchNames)
877  {
878  labelList matchIndices = findMatchingStrings(matcher, allPatchNames);
879  ids.insert(matchIndices);
880 
881  bool missed = matchIndices.empty();
882 
883  if (missed && checkGroups)
884  {
885  // Check group names
886  if (matcher.isPattern())
887  {
888  forAllConstIters(groupPatchIDs(), iter)
889  {
890  if (matcher.match(iter.key()))
891  {
892  // Hash ids associated with the group
893  ids.insert(iter.val());
894  missed = false;
895  }
896  }
897  }
898  else
899  {
900  const auto iter = groupPatchIDs().cfind(matcher);
901 
902  if (iter.found())
903  {
904  // Hash ids associated with the group
905  ids.insert(iter.val());
906  missed = false;
907  }
908  }
909  }
910 
911  if (missed && warnNotFound)
912  {
913  if (checkGroups)
914  {
916  << "Cannot find any patch or group names matching "
917  << matcher << endl;
918  }
919  else
920  {
922  << "Cannot find any patch names matching "
923  << matcher << endl;
924  }
925  }
926  }
927 
928  return ids;
929 }
930 
931 
933 (
934  const labelUList& patchIDs,
935  wordList& groups,
936  labelHashSet& nonGroupPatches
937 ) const
938 {
939  // Current matched groups
940  DynamicList<word> matchedGroups(1);
941 
942  // Current set of unmatched patches
943  nonGroupPatches = labelHashSet(patchIDs);
944 
945  const HashTable<labelList>& groupLookup = this->groupPatchIDs();
946  forAllConstIters(groupLookup, iter)
947  {
948  // Store currently unmatched patches so we can restore
949  labelHashSet oldNonGroupPatches(nonGroupPatches);
950 
951  // Match by deleting patches in group from the current set and seeing
952  // if all have been deleted.
953  labelHashSet groupPatchSet(iter.val());
954 
955  label nMatch = nonGroupPatches.erase(groupPatchSet);
956 
957  if (nMatch == groupPatchSet.size())
958  {
959  matchedGroups.append(iter.key());
960  }
961  else if (nMatch != 0)
962  {
963  // No full match. Undo.
964  nonGroupPatches.transfer(oldNonGroupPatches);
965  }
966  }
967 
968  groups.transfer(matchedGroups);
969 }
970 
971 
972 bool Foam::polyBoundaryMesh::checkParallelSync(const bool report) const
973 {
974  if (!Pstream::parRun())
975  {
976  return false;
977  }
978 
979  const polyBoundaryMesh& bm = *this;
980 
981  bool hasError = false;
982 
983  // Collect non-proc patches and check proc patches are last.
984  wordList names(bm.size());
985  wordList types(bm.size());
986 
987  label nonProci = 0;
988 
989  forAll(bm, patchi)
990  {
991  if (!isA<processorPolyPatch>(bm[patchi]))
992  {
993  if (nonProci != patchi)
994  {
995  // There is processor patch in between normal patches.
996  hasError = true;
997 
998  if (debug || report)
999  {
1000  Pout<< " ***Problem with boundary patch " << patchi
1001  << " named " << bm[patchi].name()
1002  << " of type " << bm[patchi].type()
1003  << ". The patch seems to be preceeded by processor"
1004  << " patches. This is can give problems."
1005  << endl;
1006  }
1007  }
1008  else
1009  {
1010  names[nonProci] = bm[patchi].name();
1011  types[nonProci] = bm[patchi].type();
1012  nonProci++;
1013  }
1014  }
1015  }
1016  names.setSize(nonProci);
1017  types.setSize(nonProci);
1018 
1019  List<wordList> allNames(Pstream::nProcs());
1020  allNames[Pstream::myProcNo()] = names;
1021  Pstream::gatherList(allNames);
1022  Pstream::scatterList(allNames);
1023 
1024  List<wordList> allTypes(Pstream::nProcs());
1025  allTypes[Pstream::myProcNo()] = types;
1026  Pstream::gatherList(allTypes);
1027  Pstream::scatterList(allTypes);
1028 
1029  // Have every processor check but print error on master
1030  // (in processor sequence).
1031 
1032  for (const int proci : Pstream::subProcs())
1033  {
1034  if
1035  (
1036  (allNames[proci] != allNames.first())
1037  || (allTypes[proci] != allTypes.first())
1038  )
1039  {
1040  hasError = true;
1041 
1042  if (debug || (report && Pstream::master()))
1043  {
1044  Info<< " ***Inconsistent patches across processors, "
1045  "processor 0 has patch names:"
1046  << allNames.first()
1047  << " patch types:" << allTypes.first()
1048  << " processor " << proci
1049  << " has patch names:" << allNames[proci]
1050  << " patch types:" << allTypes[proci]
1051  << endl;
1052  }
1053  }
1054  }
1055 
1056  return hasError;
1057 }
1058 
1059 
1060 bool Foam::polyBoundaryMesh::checkDefinition(const bool report) const
1061 {
1062  label nextPatchStart = mesh().nInternalFaces();
1063  const polyBoundaryMesh& bm = *this;
1064 
1065  bool hasError = false;
1066 
1067  wordHashSet patchNames(2*this->size());
1068 
1069  forAll(bm, patchi)
1070  {
1071  if (bm[patchi].start() != nextPatchStart && !hasError)
1072  {
1073  hasError = true;
1074 
1075  Info<< " ****Problem with boundary patch " << patchi
1076  << " named " << bm[patchi].name()
1077  << " of type " << bm[patchi].type()
1078  << ". The patch should start on face no " << nextPatchStart
1079  << " and the patch specifies " << bm[patchi].start()
1080  << "." << endl
1081  << "Possibly consecutive patches have this same problem."
1082  << " Suppressing future warnings." << endl;
1083  }
1084 
1085  if (!patchNames.insert(bm[patchi].name()) && !hasError)
1086  {
1087  hasError = true;
1088 
1089  Info<< " ****Duplicate boundary patch " << patchi
1090  << " named " << bm[patchi].name()
1091  << " of type " << bm[patchi].type()
1092  << "." << endl
1093  << "Suppressing future warnings." << endl;
1094  }
1095 
1096  nextPatchStart += bm[patchi].size();
1097  }
1098 
1099  reduce(hasError, orOp<bool>());
1100 
1101  if (debug || report)
1102  {
1103  if (hasError)
1104  {
1105  Pout<< " ***Boundary definition is in error." << endl;
1106  }
1107  else
1108  {
1109  Info<< " Boundary definition OK." << endl;
1110  }
1111  }
1112 
1113  return hasError;
1114 }
1115 
1116 
1118 {
1120 
1121  if
1122  (
1125  )
1126  {
1127  forAll(*this, patchi)
1128  {
1129  operator[](patchi).initMovePoints(pBufs, p);
1130  }
1131 
1132  pBufs.finishedSends();
1133 
1134  forAll(*this, patchi)
1135  {
1136  operator[](patchi).movePoints(pBufs, p);
1137  }
1138  }
1140  {
1141  const lduSchedule& patchSchedule = mesh().globalData().patchSchedule();
1142 
1143  // Dummy.
1144  pBufs.finishedSends();
1145 
1146  for (const auto& patchEval : patchSchedule)
1147  {
1148  const label patchi = patchEval.patch;
1149 
1150  if (patchEval.init)
1151  {
1152  operator[](patchi).initMovePoints(pBufs, p);
1153  }
1154  else
1155  {
1156  operator[](patchi).movePoints(pBufs, p);
1157  }
1158  }
1159  }
1160 }
1161 
1162 
1164 {
1165  neighbourEdgesPtr_.clear();
1166  patchIDPtr_.clear();
1167  groupIDsPtr_.clear();
1168 
1170 
1171  if
1172  (
1175  )
1176  {
1177  forAll(*this, patchi)
1178  {
1179  operator[](patchi).initUpdateMesh(pBufs);
1180  }
1181 
1182  pBufs.finishedSends();
1183 
1184  forAll(*this, patchi)
1185  {
1186  operator[](patchi).updateMesh(pBufs);
1187  }
1188  }
1190  {
1191  const lduSchedule& patchSchedule = mesh().globalData().patchSchedule();
1192 
1193  // Dummy.
1194  pBufs.finishedSends();
1195 
1196  for (const auto& patchEval : patchSchedule)
1197  {
1198  const label patchi = patchEval.patch;
1199 
1200  if (patchEval.init)
1201  {
1202  operator[](patchi).initUpdateMesh(pBufs);
1203  }
1204  else
1205  {
1206  operator[](patchi).updateMesh(pBufs);
1207  }
1208  }
1209  }
1210 }
1211 
1212 
1215  const labelUList& oldToNew,
1216  const bool validBoundary
1217 )
1218 {
1219  // Change order of patches
1220  polyPatchList::reorder(oldToNew);
1221 
1222  // Adapt indices
1223  polyPatchList& patches = *this;
1224 
1225  forAll(patches, patchi)
1226  {
1227  patches[patchi].index() = patchi;
1228  }
1229 
1230  if (validBoundary)
1231  {
1232  updateMesh();
1233  }
1234 }
1235 
1236 
1238 {
1239  const polyPatchList& patches = *this;
1240 
1241  os << patches.size() << nl << token::BEGIN_LIST << incrIndent << nl;
1242 
1243  for (const polyPatch& pp : patches)
1244  {
1245  os.beginBlock(pp.name());
1246  os << pp;
1247  os.endBlock();
1248  }
1249 
1251 
1253  return os.good();
1254 }
1255 
1256 
1259  IOstreamOption streamOpt,
1260  const bool valid
1261 ) const
1262 {
1264  return regIOobject::writeObject(streamOpt, valid);
1265 }
1266 
1267 
1268 // * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
1269 
1270 const Foam::polyPatch& Foam::polyBoundaryMesh::operator[]
1272  const word& patchName
1273 ) const
1274 {
1275  const label patchi = findPatchID(patchName);
1276 
1277  if (patchi < 0)
1278  {
1280  << "Patch named " << patchName << " not found." << nl
1281  << "Available patch names: " << names() << endl
1282  << abort(FatalError);
1283  }
1284 
1285  return operator[](patchi);
1286 }
1287 
1288 
1289 Foam::polyPatch& Foam::polyBoundaryMesh::operator[]
1291  const word& patchName
1292 )
1293 {
1294  const label patchi = findPatchID(patchName);
1295 
1296  if (patchi < 0)
1297  {
1299  << "Patch named " << patchName << " not found." << nl
1300  << "Available patch names: " << names() << endl
1301  << abort(FatalError);
1302  }
1303 
1304  return operator[](patchi);
1305 }
1306 
1307 
1308 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
1309 
1311 {
1312  pbm.writeData(os);
1313  return os;
1314 }
1315 
1316 
1317 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::IOstreamOption::UNCOMPRESSED
compression = false
Definition: IOstreamOption.H:79
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:67
Foam::HashTable< T, edge, Hash< edge > >::size
label size() const noexcept
The number of elements in table.
Definition: HashTableI.H:52
Foam::polyBoundaryMesh::patchSizes
labelList patchSizes() const
Return a list of patch sizes.
Definition: polyBoundaryMesh.C:589
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::wordRe::isPattern
bool isPattern() const noexcept
The wordRe is treated as a pattern, not as literal string.
Definition: wordReI.H:107
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::polyBoundaryMesh::nFaces
label nFaces() const
The number of boundary faces in the underlying mesh.
Definition: polyBoundaryMesh.C:617
Foam::polyBoundaryMesh::matchGroups
void matchGroups(const labelUList &patchIDs, wordList &groups, labelHashSet &nonGroupPatches) const
Match the patches to groups.
Definition: polyBoundaryMesh.C:933
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::lduSchedule
List< lduScheduleEntry > lduSchedule
Definition: lduSchedule.H:76
Foam::polyBoundaryMesh::physicalTypes
wordList physicalTypes() const
Return a list of physical types.
Definition: polyBoundaryMesh.C:567
Foam::List::resize
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
Foam::polyBoundaryMesh
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Definition: polyBoundaryMesh.H:63
Foam::PrimitivePatch::edges
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
Definition: PrimitivePatch.C:183
Foam::polyMesh::defaultRegion
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:318
Foam::DynamicList< word >
globalMeshData.H
Foam::Pstream::scatterList
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
Definition: gatherScatterList.C:215
Foam::SubList
A List obtained as a section of another List.
Definition: SubList.H:54
Foam::PrimitivePatch::nEdges
label nEdges() const
Number of edges in patch.
Definition: PrimitivePatch.H:322
Foam::PstreamBuffers
Buffers for inter-processor communications streams (UOPstream, UIPstream).
Definition: PstreamBuffers.H:88
Foam::polyBoundaryMesh::groupPatchIDs
const HashTable< labelList > & groupPatchIDs() const
The patch indices per patch group.
Definition: polyBoundaryMesh.C:481
Foam::polyBoundaryMesh::reorder
void reorder(const labelUList &oldToNew, const bool validBoundary)
Reorders patches. Ordering does not have to be done in.
Definition: polyBoundaryMesh.C:1214
Foam::edge
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:63
Foam::glTF::key
auto key(const Type &t) -> typename std::enable_if< std::is_enum< Type >::value, typename std::underlying_type< Type >::type >::type
Definition: foamGltfBase.H:108
Foam::HashTable< T, edge, Hash< edge > >::insert
bool insert(const edge &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
Definition: HashTableI.H:180
Foam::PtrListOps::findMatching
labelList findMatching(const UPtrList< T > &list, const UnaryMatchPredicate &matcher)
Extract list indices for all items with 'name()' that matches.
Foam::one
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
Definition: one.H:61
Foam::UPstream::master
static bool master(const label communicator=worldComm)
Am I the master process.
Definition: UPstream.H:457
PtrListOps.H
Functions to operate on Pointer Lists.
primitiveMesh.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::PtrList::set
const T * set(const label i) const
Return const pointer to element (can be nullptr),.
Definition: PtrList.H:138
Foam::UPstream::defaultCommsType
static commsTypes defaultCommsType
Default commsType.
Definition: UPstream.H:281
Foam::polyBoundaryMesh::clearGeom
void clearGeom()
Clear geometry at this level and at patches.
Definition: polyBoundaryMesh.C:236
Foam::Ostream::beginBlock
virtual Ostream & beginBlock(const keyType &kw)
Write begin block group with the given name.
Definition: Ostream.C:91
Foam::Pout
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
polyMesh.H
Foam::HashSet< label, Hash< label > >
Foam::incrIndent
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:346
Foam::polyBoundaryMesh::clearAddressing
void clearAddressing()
Clear addressing at this level and at patches.
Definition: polyBoundaryMesh.C:247
Foam::wordRe
A wordRe is a Foam::word, but can contain a regular expression for matching words or strings.
Definition: wordRe.H:80
Foam::polyBoundaryMesh::start
label start() const
The start label of the boundary faces in the polyMesh face list.
Definition: polyBoundaryMesh.C:611
Foam::polyBoundaryMesh::findIndex
label findIndex(const wordRe &key) const
Return patch index for the first match, return -1 if not found.
Definition: polyBoundaryMesh.C:754
Foam::polyBoundaryMesh::types
wordList types() const
Return a list of patch types.
Definition: polyBoundaryMesh.C:561
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::polyBoundaryMesh::names
wordList names() const
Return a list of patch names.
Definition: polyBoundaryMesh.C:555
Foam::polyBoundaryMesh::range
labelRange range() const
The face range for all boundary faces.
Definition: polyBoundaryMesh.C:623
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::IOstream::good
bool good() const noexcept
True if next operation might succeed.
Definition: IOstream.H:233
Foam::wordRes::match
bool match(const std::string &text, bool literal=false) const
Smart match as literal or regex, stopping on the first match.
Definition: wordResI.H:91
Foam::labelPair
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:54
Foam::wordList
List< word > wordList
A List of words.
Definition: fileName.H:62
Foam::operator<<
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Definition: boundaryPatch.C:83
Foam::polyBoundaryMesh::patchID
const labelList & patchID() const
Per boundary face label the patch index.
Definition: polyBoundaryMesh.C:456
Foam::polyBoundaryMesh::checkParallelSync
bool checkParallelSync(const bool report=false) const
Check whether all procs have all patches and in same order. Return.
Definition: polyBoundaryMesh.C:972
Foam::reduce
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Definition: PstreamReduceOps.H:51
Foam::UPstream::subProcs
static rangeType subProcs(const label communicator=worldComm)
Range of process indices for sub-processes.
Definition: UPstream.H:515
Foam::Field< vector >
Foam::findLower
label findLower(const ListType &input, const T &val, const label start, const ComparePredicate &comp)
Foam::Istream
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:61
Foam::HashTable< T, edge, Hash< edge > >::erase
bool erase(const iterator &iter)
Erase an entry specified by given iterator.
Definition: HashTable.C:392
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::polyBoundaryMesh::faceCells
UPtrList< const labelUList > faceCells() const
Return a list of faceCells for each patch.
Definition: polyBoundaryMesh.C:311
Foam::polyBoundaryMesh::checkDefinition
bool checkDefinition(const bool report=false) const
Check boundary definition. Return true if in error.
Definition: polyBoundaryMesh.C:1060
Foam::polyPatchList
PtrList< polyPatch > polyPatchList
container classes for polyPatch
Definition: polyPatchList.H:47
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:222
Foam::DynamicList::append
DynamicList< T, SizeMin > & append(const T &val)
Append an element to the end of this list.
Definition: DynamicListI.H:511
Foam::List::transfer
void transfer(List< T > &list)
Definition: List.C:456
Foam::IOobject::READ_IF_PRESENT
Definition: IOobject.H:187
Foam::IOstreamOption
The IOstreamOption is a simple container for options an IOstream can normally have.
Definition: IOstreamOption.H:63
Foam::PtrList::setSize
void setSize(const label newLen)
Same as resize()
Definition: PtrList.H:151
Foam::polyBoundaryMesh::patchStarts
labelList patchStarts() const
Return a list of patch start face indices.
Definition: polyBoundaryMesh.C:578
Foam::UPtrList
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
Definition: UPtrList.H:62
Foam::labelRange
A range or interval of labels defined by a start and a size.
Definition: labelRange.H:55
Foam::nameOp
Extract name (as a word) from an object, typically using its name() method.
Definition: word.H:237
Foam::polyBoundaryMesh::whichPatch
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
Definition: polyBoundaryMesh.C:812
Foam::reorder
ListType reorder(const labelUList &oldToNew, const ListType &input, const bool prune=false)
Reorder the elements of a list.
Definition: ListOpsTemplates.C:80
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:59
patchNames
wordList patchNames(nPatches)
Foam::polyBoundaryMesh::findPatchID
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
Definition: polyBoundaryMesh.C:765
Foam::PtrList::append
void append(T *ptr)
Append an element to the end of the list.
Definition: PtrListI.H:113
lduSchedule.H
Foam::PstreamBuffers::finishedSends
void finishedSends(const bool block=true)
Mark all sends as having been done.
Definition: PstreamBuffers.C:73
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::Ostream::endBlock
virtual Ostream & endBlock()
Write end block group.
Definition: Ostream.C:109
Foam::IOstream::check
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:58
Foam::FatalError
error FatalError
processorPolyPatch.H
Foam::polyPatch::New
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return pointer to a new patch created on freestore from components.
Definition: polyPatchNew.C:35
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::regIOobject::writeObject
virtual bool writeObject(IOstreamOption streamOpt, const bool valid) const
Write using stream options.
Definition: regIOobjectWrite.C:36
Foam::ProcessorTopology::patchSchedule
const lduSchedule & patchSchedule() const
Order in which the patches should be initialised/evaluated.
Definition: ProcessorTopology.H:101
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::flatOutput
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Definition: FlatOutput.H:216
Foam::decrIndent
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:353
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::HashTable< T, edge, Hash< edge > >::find
iterator find(const edge &key)
Find and return an iterator set at the hashed entry.
Definition: HashTableI.H:114
Foam::polyBoundaryMesh::movePoints
void movePoints(const pointField &p)
Correct polyBoundaryMesh after moving points.
Definition: polyBoundaryMesh.C:1117
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::PrimitivePatch::localPoints
const Field< point_type > & localPoints() const
Return pointField of points in patch.
Definition: PrimitivePatch.C:359
Foam::HashTable
A HashTable similar to std::unordered_map.
Definition: HashTable.H:105
Foam::polyBoundaryMesh::updateMesh
void updateMesh()
Correct polyBoundaryMesh after topology update.
Definition: polyBoundaryMesh.C:1163
Foam::polyBoundaryMesh::nNonProcessor
label nNonProcessor() const
The number of patches before the first processor patch.
Definition: polyBoundaryMesh.C:535
Foam::PrimitivePatch::nInternalEdges
label nInternalEdges() const
Number of internal edges.
Definition: PrimitivePatch.C:214
Foam::polyBoundaryMesh::writeData
virtual bool writeData(Ostream &os) const
writeData member function required by regIOobject
Definition: polyBoundaryMesh.C:1237
Foam::IOobject::name
const word & name() const noexcept
Return name.
Definition: IOobjectI.H:65
Foam::EdgeMap
Map from edge (expressed as its endpoints) to value. For easier forward declaration it is currently i...
Definition: EdgeMap.H:51
Foam::regIOobject
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:73
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::Pstream::gatherList
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
Definition: gatherScatterList.C:52
Foam::polyBoundaryMesh::neighbourEdges
const List< labelPairList > & neighbourEdges() const
Per patch the edges on the neighbouring patch.
Definition: polyBoundaryMesh.C:327
edgeHashes.H
Foam::UPstream::commsTypes::nonBlocking
Foam::UPstream::myProcNo
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:463
range
scalar range
Definition: LISASMDCalcMethod1.H:12
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::Pair< label >
forAllConstIters
forAllConstIters(mixture.phases(), phase)
Definition: pEqn.H:28
Foam::UPstream::commsTypes::scheduled
Foam::polyBoundaryMesh::writeObject
virtual bool writeObject(IOstreamOption streamOpt, const bool valid) const
Write using stream options.
Definition: polyBoundaryMesh.C:1258
Foam::primitiveMesh::nInternalFaces
label nInternalFaces() const noexcept
Number of internal faces.
Definition: primitiveMeshI.H:78
Foam::UPstream::parRun
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:433
Foam::polyBoundaryMesh::indices
labelList indices(const wordRe &matcher, const bool useGroups=true) const
Return (sorted) patch indices for all matches.
Definition: polyBoundaryMesh.C:642
Foam::List
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: BitOps.H:63
Foam::polyBoundaryMesh::setGroup
void setGroup(const word &groupName, const labelUList &patchIDs)
Set/add group with patches.
Definition: polyBoundaryMesh.C:493
Foam::UList< label >
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
Foam::wordRes
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:51
Foam::polyBoundaryMesh::patchSet
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool useGroups=true) const
Return the set of patch IDs corresponding to the given names.
Definition: polyBoundaryMesh.C:864
Foam::polyBoundaryMesh::patchRanges
List< labelRange > patchRanges() const
Return a list of patch ranges.
Definition: polyBoundaryMesh.C:600
PstreamBuffers.H
Foam::HashSet::insert
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
Foam::PtrListOps::firstMatching
label firstMatching(const UPtrList< T > &list, const UnaryMatchPredicate &matcher)
Find first list item with 'name()' that matches, -1 on failure.
Foam::IOobject::MUST_READ_IF_MODIFIED
Definition: IOobject.H:186
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:65
polyBoundaryMesh.H
FUNCTION_NAME
#define FUNCTION_NAME
Definition: messageStream.H:295
Foam::token::END_LIST
End list [isseparator].
Definition: token.H:156
stringListOps.H
Operations on lists of strings.
Foam::primitiveMesh::nFaces
label nFaces() const noexcept
Number of mesh faces.
Definition: primitiveMeshI.H:90
DynamicList.H
patchId
label patchId(-1)
Foam::patchIdentifier::name
const word & name() const noexcept
The patch name.
Definition: patchIdentifier.H:135
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::polyMesh::globalData
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1295
Foam::IOstreamOption::compression
compressionType compression() const noexcept
Get the stream compression.
Definition: IOstreamOption.H:312
Foam::labelHashSet
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
Definition: HashSet.H:85
Foam::PrimitivePatch::meshPoints
const labelList & meshPoints() const
Return labelList of mesh points in patch.
Definition: PrimitivePatch.C:330
Foam::token::BEGIN_LIST
Begin list [isseparator].
Definition: token.H:155
Foam::orOp
Definition: ops.H:234
Foam::PtrListOps::names
List< word > names(const UPtrList< T > &list, const UnaryMatchPredicate &matcher)
Foam::findMatchingStrings
labelList findMatchingStrings(const UnaryMatchPredicate &matcher, const UList< StringType > &input, const bool invert=false)
Extract list indices for all matches.
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::typeOp
Extract type (as a word) from an object, typically using its type() method.
Definition: word.H:248
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:328
Foam::faceCells
Smooth ATC in cells next to a set of patches supplied by type.
Definition: faceCells.H:56
Foam::wordRe::match
bool match(const std::string &text, bool literal=false) const
Smart match as regular expression or as a string.
Definition: wordReI.H:174
Foam::UPstream::nProcs
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
Definition: UPstream.H:445
Foam::UPstream::commsTypes::blocking
Foam::IOobject::MUST_READ
Definition: IOobject.H:185