triSurface.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2016-2020 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 "triSurface.H"
30 #include "Time.H"
31 #include "surfZoneList.H"
32 #include "MeshedSurface.H"
33 #include "ListOps.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(triSurface, 0);
40 }
41 
42 
43 // * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 // Helper function to print triangle info
49 static void printTriangle
50 (
51  Ostream& os,
52  const string& pre,
53  const labelledTri& f,
54  const pointField& points
55 )
56 {
57  os
58  << pre.c_str() << "vertex numbers:"
59  << f[0] << ' ' << f[1] << ' ' << f[2] << nl
60 
61  << pre.c_str() << "vertex coords :"
62  << points[f[0]] << ' ' << points[f[1]] << ' ' << points[f[2]]
63 
64  << pre.c_str() << "region :" << f.region() << nl
65  << endl;
66 }
67 
68 } // End namespace Foam
69 
70 
71 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
72 
74 {
75  fileName foamName(d.caseName() + ".ftr");
76 
77  // Search back through the time directories list to find the time
78  // closest to and lower than current time
79 
80  instantList ts = d.times();
81  label i;
82 
83  for (i=ts.size()-1; i>=0; i--)
84  {
85  if (ts[i].value() <= d.timeOutputValue())
86  {
87  break;
88  }
89  }
90 
91  // Noting that the current directory has already been searched
92  // for mesh data, start searching from the previously stored time directory
93 
94  if (i>=0)
95  {
96  for (label j=i; j>=0; j--)
97  {
98  if (isFile(d.path()/ts[j].name()/typeName/foamName))
99  {
100  if (debug)
101  {
102  Pout<< " triSurface::triSurfInstance(const Time& d)"
103  << "reading " << foamName
104  << " from " << ts[j].name()/typeName
105  << endl;
106  }
107 
108  return ts[j].name();
109  }
110  }
111  }
112 
113  if (debug)
114  {
115  Pout<< " triSurface::triSurfInstance(const Time& d)"
116  << "reading " << foamName
117  << " from constant/" << endl;
118  }
119  return d.constant();
120 }
121 
122 
123 Foam::List<Foam::labelledTri> Foam::triSurface::convertToTri
124 (
125  const faceList& faces,
126  const label defaultRegion
127 )
128 {
129  List<labelledTri> triFaces(faces.size());
130 
131  forAll(triFaces, facei)
132  {
133  const face& f = faces[facei];
134 
135  if (f.size() != 3)
136  {
138  << "Face at position " << facei
139  << " does not have three vertices:" << f
140  << abort(FatalError);
141  }
142 
143  labelledTri& tri = triFaces[facei];
144 
145  tri[0] = f[0];
146  tri[1] = f[1];
147  tri[2] = f[2];
148  tri.region() = defaultRegion;
149  }
150 
151  return triFaces;
152 }
153 
154 
155 Foam::List<Foam::labelledTri> Foam::triSurface::convertToTri
156 (
157  const triFaceList& faces,
158  const label defaultRegion
159 )
160 {
161  List<labelledTri> triFaces(faces.size());
162 
163  forAll(triFaces, facei)
164  {
165  const triFace& f = faces[facei];
166 
167  labelledTri& tri = triFaces[facei];
168 
169  tri[0] = f[0];
170  tri[1] = f[1];
171  tri[2] = f[2];
172  tri.region() = defaultRegion;
173  }
174 
175  return triFaces;
176 }
177 
178 
179 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
180 
181 // Remove non-triangles, double triangles.
182 void Foam::triSurface::checkTriangles(const bool verbose)
183 {
184  // Simple check on indices ok.
185  const label maxPointi = points().size() - 1;
186 
187  for (const auto& f : *this)
188  {
189  for (const label verti : f)
190  {
191  if (verti < 0 || verti > maxPointi)
192  {
194  << "triangle " << f
195  << " uses point indices outside point range 0.."
196  << maxPointi
197  << exit(FatalError);
198  }
199  }
200  }
201 
202  // Two phase process
203  // 1. mark invalid faces
204  // 2. pack
205  // Done to keep numbering constant in phase 1
206 
207  // List of valid triangles
208  bitSet valid(size(), true);
209 
210  forAll(*this, facei)
211  {
212  const labelledTri& f = (*this)[facei];
213 
214  if ((f[0] == f[1]) || (f[0] == f[2]) || (f[1] == f[2]))
215  {
216  // 'degenerate' triangle check
217  valid.unset(facei);
218 
219  if (verbose)
220  {
222  << "triangle " << facei
223  << " does not have three unique vertices:\n";
224  printTriangle(Warning, " ", f, points());
225  }
226  }
227  else
228  {
229  // duplicate triangle check
230  const labelList& fEdges = faceEdges()[facei];
231 
232  // Check if faceNeighbours use same points as this face.
233  // Note: discards normal information - sides of baffle are merged.
234 
235  for (const label edgei : fEdges)
236  {
237  const labelList& eFaces = edgeFaces()[edgei];
238 
239  for (const label neighbour : eFaces)
240  {
241  if (neighbour > facei)
242  {
243  // lower numbered faces already checked
244  const labelledTri& n = (*this)[neighbour];
245 
246  if
247  (
248  ((f[0] == n[0]) || (f[0] == n[1]) || (f[0] == n[2]))
249  && ((f[1] == n[0]) || (f[1] == n[1]) || (f[1] == n[2]))
250  && ((f[2] == n[0]) || (f[2] == n[1]) || (f[2] == n[2]))
251  )
252  {
253  valid.unset(facei);
254 
255  if (verbose)
256  {
258  << "triangles share the same vertices:\n"
259  << " face 1 :" << facei << endl;
260  printTriangle(Warning, " ", f, points());
261 
262  Warning
263  << endl
264  << " face 2 :"
265  << neighbour << endl;
266  printTriangle(Warning, " ", n, points());
267  }
268 
269  break;
270  }
271  }
272  }
273  }
274  }
275  }
276 
277  if (!valid.all())
278  {
279  // Compact
280  label newFacei = 0;
281  for (const label facei : valid)
282  {
283  (*this)[newFacei++] = (*this)[facei];
284  }
285 
286  if (verbose)
287  {
289  << "Removing " << size() - newFacei
290  << " illegal faces." << endl;
291  }
292  (*this).setSize(newFacei);
293 
294  // Topology can change because of renumbering
295  clearOut();
296  }
297 }
298 
299 
300 // Check/fix edges with more than two triangles
301 void Foam::triSurface::checkEdges(const bool verbose)
302 {
303  const labelListList& eFaces = edgeFaces();
304 
305  forAll(eFaces, edgei)
306  {
307  const labelList& myFaces = eFaces[edgei];
308 
309  if (myFaces.empty())
310  {
312  << "Edge " << edgei << " with vertices " << edges()[edgei]
313  << " has no edgeFaces"
314  << exit(FatalError);
315  }
316  else if (myFaces.size() > 2 && verbose)
317  {
319  << "Edge " << edgei << " with vertices " << edges()[edgei]
320  << " has more than 2 faces connected to it : " << myFaces
321  << endl;
322  }
323  }
324 }
325 
326 
327 // Returns patch info. Sets faceMap to the indexing according to patch
328 // numbers. Patch numbers start at 0.
330 Foam::triSurface::calcPatches(labelList& faceMap) const
331 {
332  // Determine the sorted order:
333  // use sortedOrder directly (the intermediate list is discarded anyhow)
334 
335  labelList regions(size());
336  forAll(regions, facei)
337  {
338  regions[facei] = operator[](facei).region();
339  }
340  sortedOrder(regions, faceMap);
341  regions.clear();
342 
343  // Extend regions
344  label maxRegion = patches_.size()-1; // for non-compacted regions
345 
346  if (faceMap.size())
347  {
348  maxRegion = max
349  (
350  maxRegion,
351  operator[](faceMap.last()).region()
352  );
353  }
354 
355  // Get new region list
356  surfacePatchList newPatches(maxRegion + 1);
357 
358  // Fill patch sizes
359  forAll(*this, facei)
360  {
361  label region = operator[](facei).region();
362 
363  newPatches[region].size()++;
364  }
365 
366 
367  // Fill rest of patch info
368 
369  label startFacei = 0;
370  forAll(newPatches, newPatchi)
371  {
372  surfacePatch& newPatch = newPatches[newPatchi];
373 
374  newPatch.index() = newPatchi;
375  newPatch.start() = startFacei;
376 
377  // Take over any information from existing patches
378  if
379  (
380  newPatchi < patches_.size()
381  && !patches_[newPatchi].name().empty()
382  )
383  {
384  newPatch.name() = patches_[newPatchi].name();
385  }
386  else
387  {
388  newPatch.name() = surfacePatch::defaultName(newPatchi);
389  }
390 
391  if
392  (
393  newPatchi < patches_.size()
394  && !patches_[newPatchi].geometricType().empty()
395  )
396  {
397  newPatch.geometricType() = patches_[newPatchi].geometricType();
398  }
399  else
400  {
401  newPatch.geometricType() = surfacePatch::emptyType;
402  }
403 
404  startFacei += newPatch.size();
405  }
406 
407  return newPatches;
408 }
409 
410 
411 void Foam::triSurface::setDefaultPatches()
412 {
414 
415  // Get names, types and sizes
416  surfacePatchList newPatches(calcPatches(faceMap));
417 
418  // Take over names and types (but not size)
419  patches_.setSize(newPatches.size());
420 
421  forAll(newPatches, patchi)
422  {
423  patches_[patchi].index() = patchi;
424  patches_[patchi].name() = newPatches[patchi].name();
425  patches_[patchi].geometricType() = newPatches[patchi].geometricType();
426  }
427 }
428 
429 
430 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
431 
433 :
435  patches_(),
436  sortedEdgeFacesPtr_(nullptr),
437  edgeOwnerPtr_(nullptr)
438 {}
439 
440 
442 :
443  MeshReference(surf, surf.points()),
444  patches_(surf.patches()),
445  sortedEdgeFacesPtr_(nullptr),
446  edgeOwnerPtr_(nullptr)
447 {}
448 
449 
451 :
452  triSurface()
453 {
454  transfer(surf);
455 }
456 
457 
459 (
460  const List<labelledTri>& triangles,
462  const pointField& pts
463 )
464 :
465  MeshReference(triangles, pts),
466  patches_(patches),
467  sortedEdgeFacesPtr_(nullptr),
468  edgeOwnerPtr_(nullptr)
469 {}
470 
471 
473 (
474  List<labelledTri>& triangles,
476  pointField& pts,
477  const bool reuse
478 )
479 :
480  MeshReference(triangles, pts, reuse),
481  patches_(patches),
482  sortedEdgeFacesPtr_(nullptr),
483  edgeOwnerPtr_(nullptr)
484 {}
485 
486 
488 (
489  const List<labelledTri>& triangles,
490  const pointField& pts
491 )
492 :
493  MeshReference(triangles, pts),
494  patches_(),
495  sortedEdgeFacesPtr_(nullptr),
496  edgeOwnerPtr_(nullptr)
497 {
498  setDefaultPatches();
499 }
500 
501 
503 (
504  const triFaceList& triangles,
505  const pointField& pts
506 )
507 :
508  MeshReference(convertToTri(triangles, 0), pts),
509  patches_(),
510  sortedEdgeFacesPtr_(nullptr),
511  edgeOwnerPtr_(nullptr)
512 {
513  setDefaultPatches();
514 }
515 
516 
518 (
519  const fileName& name,
520  const scalar scaleFactor
521 )
522 :
523  triSurface(name, name.ext(), scaleFactor)
524 {}
525 
526 
528 (
529  const fileName& name,
530  const word& fileType,
531  const scalar scaleFactor
532 )
533 :
534  triSurface()
535 {
536  read(name, fileType);
537  scalePoints(scaleFactor);
538  setDefaultPatches();
539 }
540 
541 
542 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
543 
545 {
546  clearOut();
547 }
548 
549 
550 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
551 
553 {
554  MeshReference::clearTopology();
555  sortedEdgeFacesPtr_.reset(nullptr);
556  edgeOwnerPtr_.reset(nullptr);
557 }
558 
559 
561 {
562  MeshReference::clearPatchMeshAddr();
563 }
564 
565 
567 {
568  MeshReference::clearOut();
569  clearTopology();
570  clearPatchMeshAddr();
571 }
572 
573 
575 {
576  if (this == &surf)
577  {
578  return; // Self-swap is a no-op
579  }
580 
581  clearOut();
582  surf.clearOut();
583 
584  storedFaces().swap(surf.storedFaces());
585  storedPoints().swap(surf.storedPoints());
586  patches_.swap(surf.patches());
587 }
588 
589 
591 {
592  if (!sortedEdgeFacesPtr_)
593  {
594  calcSortedEdgeFaces();
595  }
596 
597  return *sortedEdgeFacesPtr_;
598 }
599 
600 
602 {
603  if (!edgeOwnerPtr_)
604  {
605  calcEdgeOwner();
606  }
607 
608  return *edgeOwnerPtr_;
609 }
610 
611 
613 {
614  // Remove all geometry dependent data
615  sortedEdgeFacesPtr_.reset(nullptr);
616 
617  // Adapt for new point positions
618  MeshReference::movePoints(pts);
619 
620  // Copy new points
621  storedPoints() = pts;
622 }
623 
624 
626 {
627  // Remove all geometry dependent data
628  sortedEdgeFacesPtr_.reset(nullptr);
629 
630  // Adapt for new point positions
631  MeshReference::movePoints(pts);
632 
633  // Move/swap new points
634  storedPoints().swap(pts);
635 }
636 
637 
638 void Foam::triSurface::scalePoints(const scalar scaleFactor)
639 {
640  // Avoid bad or no scaling
641  if (scaleFactor > SMALL && !equal(scaleFactor, 1))
642  {
643  // Remove all geometry dependent data
644  this->clearTopology();
645 
646  // Adapt for new point positions
647  MeshReference::movePoints(pointField());
648 
649  this->storedPoints() *= scaleFactor;
650  }
651 }
652 
653 
654 // Remove non-triangles, double triangles.
655 void Foam::triSurface::cleanup(const bool verbose)
656 {
657  // Merge points (already done for STL, TRI)
658  stitchTriangles(SMALL, verbose);
659 
660  // Merging points might have changed geometric factors
661  clearOut();
662 
663  checkTriangles(verbose);
664 
665  checkEdges(verbose);
666 }
667 
668 
670 {
671  this->clearOut(); // Topology changes
672 
673  // Remove unused points while walking and renumbering faces
674  // in visit order - walk order as per localFaces()
675 
676  labelList oldToCompact(this->points().size(), -1);
677  DynamicList<label> compactPointMap(oldToCompact.size());
678 
679  for (auto& f : this->storedFaces())
680  {
681  for (label& pointi : f)
682  {
683  label compacti = oldToCompact[pointi];
684  if (compacti == -1)
685  {
686  compacti = compactPointMap.size();
687  oldToCompact[pointi] = compacti;
688  compactPointMap.append(pointi);
689  }
690  pointi = compacti;
691  }
692  }
693 
694  pointField newPoints
695  (
696  UIndirectList<point>(this->points(), compactPointMap)
697  );
698 
699  this->swapPoints(newPoints);
700 
701  if (notNull(pointMap))
702  {
703  pointMap.transfer(compactPointMap);
704  }
705 }
706 
707 
710 {
711  surfacePatchList patches(calcPatches(faceMap));
712 
713  surfZoneList zones(patches.size());
714  forAll(patches, patchi)
715  {
716  zones[patchi] = surfZone(patches[patchi]);
717  }
718 
719  return zones;
720 }
721 
722 
724 {
725  plainFaces.setSize(size());
726 
727  forAll(*this, facei)
728  {
729  plainFaces[facei] = this->operator[](facei);
730  }
731 }
732 
733 
734 // Finds area, starting at facei, delimited by borderEdge. Marks all visited
735 // faces (from face-edge-face walk) with currentZone.
737 (
738  const boolList& borderEdge,
739  const label facei,
740  const label currentZone,
742 ) const
743 {
744  // List of faces whose faceZone has been set.
745  labelList changedFaces(1, facei);
746 
747  while (true)
748  {
749  // Pick up neighbours of changedFaces
750  DynamicList<label> newChangedFaces(2*changedFaces.size());
751 
752  for (const label facei : changedFaces)
753  {
754  const labelList& fEdges = faceEdges()[facei];
755 
756  for (const label edgei : fEdges)
757  {
758  if (!borderEdge[edgei])
759  {
760  const labelList& eFaces = edgeFaces()[edgei];
761 
762  for (const label nbrFacei : eFaces)
763  {
764  if (faceZone[nbrFacei] == -1)
765  {
766  faceZone[nbrFacei] = currentZone;
767  newChangedFaces.append(nbrFacei);
768  }
769  else if (faceZone[nbrFacei] != currentZone)
770  {
772  << "Zones " << faceZone[nbrFacei]
773  << " at face " << nbrFacei
774  << " connects to zone " << currentZone
775  << " at face " << facei
776  << abort(FatalError);
777  }
778  }
779  }
780  }
781  }
782 
783  if (newChangedFaces.empty())
784  {
785  break;
786  }
787 
788  changedFaces.transfer(newChangedFaces);
789  }
790 }
791 
792 
793 // Finds areas delimited by borderEdge (or 'real' edges).
794 // Fills faceZone accordingly
795 Foam::label Foam::triSurface::markZones
796 (
797  const boolList& borderEdge,
799 ) const
800 {
801  faceZone.setSize(size());
802  faceZone = -1;
803 
804  if (borderEdge.size() != nEdges())
805  {
807  << "borderEdge boolList not same size as number of edges" << endl
808  << "borderEdge:" << borderEdge.size() << endl
809  << "nEdges :" << nEdges()
810  << exit(FatalError);
811  }
812 
813  label zoneI = 0;
814 
815  for (label startFacei = 0;; zoneI++)
816  {
817  // Find first non-coloured face
818  for (; startFacei < size(); startFacei++)
819  {
820  if (faceZone[startFacei] == -1)
821  {
822  break;
823  }
824  }
825 
826  if (startFacei >= size())
827  {
828  break;
829  }
830 
831  faceZone[startFacei] = zoneI;
832 
833  markZone(borderEdge, startFacei, zoneI, faceZone);
834  }
835 
836  return zoneI;
837 }
838 
839 
840 Foam::triSurface Foam::triSurface::subsetMeshImpl
841 (
842  const labelList& pointMap,
843  const labelList& faceMap
844 ) const
845 {
846  const pointField& locPoints = localPoints();
847  const List<labelledTri>& locFaces = localFaces();
848 
849  // Subset of points (compact)
850  pointField newPoints(UIndirectList<point>(locPoints, pointMap));
851 
852  // Inverse point mapping - same as ListOps invert() without checks
853  labelList oldToNew(locPoints.size(), -1);
854  forAll(pointMap, pointi)
855  {
856  oldToNew[pointMap[pointi]] = pointi;
857  }
858 
859  // Subset of faces
860  List<labelledTri> newFaces(UIndirectList<labelledTri>(locFaces, faceMap));
861 
862  // Renumber face node labels
863  for (auto& f : newFaces)
864  {
865  for (label& vert : f)
866  {
867  vert = oldToNew[vert];
868  }
869  }
870  oldToNew.clear();
871 
872  // Construct sub-surface
873  return triSurface(newFaces, patches(), newPoints, true);
874 }
875 
876 
879 (
880  const UList<bool>& include,
881  labelList& pointMap,
883 ) const
884 {
885  this->subsetMeshMap(include, pointMap, faceMap);
886  return this->subsetMeshImpl(pointMap, faceMap);
887 }
888 
889 
892 (
893  const bitSet& include,
894  labelList& pointMap,
896 ) const
897 {
898  this->subsetMeshMap(include, pointMap, faceMap);
899  return this->subsetMeshImpl(pointMap, faceMap);
900 }
901 
902 
905 {
906  labelList pointMap, faceMap;
907  return this->subsetMesh(include, pointMap, faceMap);
908 }
909 
910 
913 {
914  labelList pointMap, faceMap;
915  return this->subsetMesh(include, pointMap, faceMap);
916 }
917 
918 
921 (
922  const wordRes& includeNames,
923  const wordRes& excludeNames
924 ) const
925 {
926  const bitSet selectPatches
927  (
929  (
930  patches_,
931  includeNames,
932  excludeNames,
934  )
935  );
936 
937  bitSet include(this->size());
938 
939  forAll(*this, facei)
940  {
941  const label patchi = (*this)[facei].region();
942 
943  if (selectPatches.test(patchi))
944  {
945  include.set(facei);
946  }
947  }
948 
949  return this->subsetMesh(include);
950 }
951 
952 
954 {
955  clearOut(); // Topology changes
956 
957  this->storedFaces().swap(faceLst);
958 }
959 
960 
962 {
963  clearOut();
964 
965  storedFaces().transfer(surf.storedFaces());
966  storedPoints().transfer(surf.storedPoints());
967  patches_.transfer(surf.patches());
968 
969  surf.clearOut();
970 }
971 
972 
974 {
975  // Transcribe zone -> patch info
976  auto patches = ListOps::create<geometricSurfacePatch>
977  (
978  surf.surfZones(),
980  );
981 
982  // Fairly ugly, but the only way to get at the data safely
983  List<labelledTri> fcs;
984  pointField pts;
985 
986  surf.swapFaces(fcs);
987  surf.swapPoints(pts);
988 
989  clearOut();
990  surf.clear();
991 
992  triSurface s(fcs, patches, pts, true);
993  swap(s);
994 }
995 
996 
997 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
998 
1000 {
1001  clearOut();
1002 
1003  storedFaces() = surf;
1004  storedPoints() = surf.points();
1005  patches_ = surf.patches();
1006 }
1007 
1008 
1010 {
1011  transfer(surf);
1012 }
1013 
1014 
1016 {
1017  transfer(surf);
1018 }
1019 
1020 
1021 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::PrimitivePatch::points
const Field< point_type > & points() const noexcept
Return reference to global points.
Definition: PrimitivePatch.H:299
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:67
Foam::printTriangle
static void printTriangle(Ostream &os, const string &pre, const labelledTri &f, const pointField &points)
Definition: triSurface.C:50
Foam::triSurface::subsetMesh
triSurface subsetMesh(const UList< bool > &include, labelList &pointMap, labelList &faceMap) const
Return a new surface subsetted on the selected faces.
Definition: triSurface.C:879
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Foam::TimeState::timeOutputValue
scalar timeOutputValue() const
Return current time value.
Definition: TimeStateI.H:31
Foam::MeshedSurface::surfZones
const surfZoneList & surfZones() const
Const access to the surface zones.
Definition: MeshedSurface.H:429
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::faceMap
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
Definition: blockMeshMergeTopological.C:94
Foam::triSurface::clearOut
void clearOut()
Definition: triSurface.C:566
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::fileName
A class for handling file names.
Definition: fileName.H:73
s
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputSpray.H:25
Foam::bitSet::all
bool all() const
True if all bits in this bitset are set or if the set is empty.
Definition: bitSetI.H:460
Foam::bitSet
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition: bitSet.H:63
Foam::triSurface::markZone
void markZone(const boolList &borderEdge, const label facei, const label currentZone, labelList &faceZone) const
Fill faceZone with currentZone for every face reachable.
Definition: triSurface.C:737
Foam::DynamicList< label >
Foam::triSurface::triFaceFaces
void triFaceFaces(List< face > &plainFaceList) const
Create a list of faces from the triFaces.
Definition: triSurface.C:723
Foam::bitSet::unset
bitSet & unset(const bitSet &other)
Definition: bitSetI.H:612
Foam::fileName::name
static std::string name(const std::string &str)
Return basename (part beyond last /), including its extension.
Definition: fileNameI.H:199
Foam::Warning
messageStream Warning
Foam::read
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
Foam::surfacePatchList
List< surfacePatch > surfacePatchList
A List of surfacePatch.
Definition: surfacePatchList.H:47
Foam::geometricSurfacePatch::defaultName
static word defaultName(const label n=-1)
Default patch name: "patch" or "patchN".
Definition: geometricSurfacePatch.H:77
Foam::bitSet::set
void set(const bitSet &bitset)
Set specified bits from another bitset.
Definition: bitSetI.H:574
Foam::isFile
bool isFile(const fileName &name, const bool checkGzip=true, const bool followLink=true)
Does the name exist as a FILE in the file system?
Definition: MSwindows.C:658
Foam::List::append
void append(const T &val)
Append an element at the end of the list.
Definition: ListI.H:175
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::word::ext
word ext() const
Return file name extension (part after last .)
Definition: word.C:126
Foam::MeshedSurface::clear
virtual void clear()
Clear all storage.
Definition: MeshedSurface.C:598
Foam::triSurface::clearPatchMeshAddr
void clearPatchMeshAddr()
Definition: triSurface.C:560
Foam::Pout
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
triSurface.H
Foam::bitSet::test
bool test(const label pos) const
Test value at specified position, never auto-vivify entries.
Definition: bitSetI.H:520
Foam::triSurface::patches
const geometricSurfacePatchList & patches() const noexcept
Definition: triSurface.H:399
Foam::triSurface::~triSurface
virtual ~triSurface()
Destructor.
Definition: triSurface.C:544
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::triSurface::swapFaces
void swapFaces(List< labelledTri > &faceLst)
Swap the list of faces being addressed.
Definition: triSurface.C:953
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::triSurface::scalePoints
virtual void scalePoints(const scalar scaleFactor)
Scale points. A non-positive factor is ignored.
Definition: triSurface.C:638
Foam::triSurface::checkTriangles
void checkTriangles(const bool verbose)
Check/remove duplicate/degenerate triangles.
Definition: triSurface.C:182
Foam::triSurface::movePoints
virtual void movePoints(const pointField &pts)
Move points.
Definition: triSurface.C:612
Foam::Field< vector >
Foam::triSurface
Triangulated surface description with patch information.
Definition: triSurface.H:76
Foam::triSurface::swapPoints
virtual void swapPoints(pointField &pts)
Swap points. Similar to movePoints, but returns the old points.
Definition: triSurface.C:625
Foam::faceZone
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:64
Foam::OSstream::name
virtual const fileName & name() const
Get the name of the stream.
Definition: OSstream.H:107
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:222
Foam::List::transfer
void transfer(List< T > &list)
Definition: List.C:456
Foam::TimePaths::times
instantList times() const
Search the case for valid time directories.
Definition: TimePaths.C:149
Foam::triSurface::transfer
void transfer(triSurface &surf)
Alter contents by transferring (triangles, points) components.
Definition: triSurface.C:961
Foam::nameOp
Extract name (as a word) from an object, typically using its name() method.
Definition: word.H:237
Foam::IOstream::name
virtual const fileName & name() const
Return the name of the stream.
Definition: IOstream.C:40
Foam::triSurface::triSurface
triSurface()
Default construct.
Definition: triSurface.C:432
Foam::triSurface::storedPoints
pointField & storedPoints()
Non-const access to global points.
Definition: triSurface.H:189
Foam::triSurface::sortedZones
List< surfZone > sortedZones(labelList &faceMap) const
Sort faces according to zoneIds.
Definition: triSurface.C:709
Foam::stringListOps::findMatching
labelList findMatching(const StringListType &input, const wordRes &allow, const wordRes &deny=wordRes(), AccessOp aop=noOp())
Return ids for items with matching names.
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::MeshedSurface::swapPoints
void swapPoints(pointField &points)
Swap the stored points.
Definition: MeshedSurface.C:1422
Foam::notNull
bool notNull(const T *ptr)
True if ptr is not a pointer (of type T) to the nullObject.
Definition: nullObject.H:207
Foam::FatalError
error FatalError
os
OBJstream os(runTime.globalPath()/outputName)
Foam::triSurface::clearTopology
void clearTopology()
Definition: triSurface.C:552
Foam::triFaceList
List< triFace > triFaceList
list of triFaces
Definition: triFaceList.H:44
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::triSurface::edgeOwner
const labelList & edgeOwner() const
If 2 face neighbours: label of face where ordering of edge.
Definition: triSurface.C:601
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::triSurface::cleanup
void cleanup(const bool verbose)
Remove non-valid triangles.
Definition: triSurface.C:655
Time.H
Foam::triSurface::markZones
label markZones(const boolList &borderEdge, labelList &faceZone) const
(size and) fills faceZone with zone of face. Zone is area
Definition: triSurface.C:796
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::Time::caseName
const fileName & caseName() const
Return case name.
Definition: TimePathsI.H:62
Foam::geometricSurfacePatch::fromIdentifier
Helper to convert identifier types as an operation.
Definition: geometricSurfacePatch.H:91
Foam::triSurface::swap
void swap(triSurface &surf)
Definition: triSurface.C:574
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::Time::path
fileName path() const
Return path.
Definition: Time.H:358
f
labelList f(nPoints)
Foam::List< instant >
Foam::surfZone
A surface zone on a MeshedSurface.
Definition: surfZone.H:56
Foam::triSurface::storedFaces
List< labelledTri > & storedFaces()
Non-const access to the faces.
Definition: triSurface.H:195
Foam::UList
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:103
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::wordRes
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:51
Foam::geometricSurfacePatch::emptyType
static constexpr const char *const emptyType
The name for an 'empty' type.
Definition: geometricSurfacePatch.H:71
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:65
Foam::labelledTri
A triFace with additional (region) index.
Definition: labelledTri.H:57
surfZoneList.H
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::UIndirectList
A List with indirect addressing.
Definition: faMatrix.H:60
Foam::face
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:72
ListOps.H
Various functions to operate on Lists.
Foam::triSurface::sortedEdgeFaces
const labelListList & sortedEdgeFaces() const
Return edge-face addressing sorted (for edges with more than.
Definition: triSurface.C:590
Foam::sortedOrder
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::triSurface::operator=
void operator=(const triSurface &surf)
Copy assignment.
Definition: triSurface.C:999
Foam::MeshedSurface::swapFaces
void swapFaces(List< Face > &faces)
Swap the stored faces. Use with caution.
Definition: MeshedSurface.C:1409
Foam::TimePaths::constant
const word & constant() const
Return constant name.
Definition: TimePathsI.H:96
Foam::MeshedSurface
A surface geometry mesh with zone information, not to be confused with the similarly named surfaceMes...
Definition: triSurfaceTools.H:80
Foam::equal
bool equal(const T &s1, const T &s2)
Compare two values for equality.
Definition: doubleFloat.H:46
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
triFace
face triFace(3)
Foam::triSurface::compactPoints
void compactPoints(labelList &pointMap=const_cast< labelList & >(labelList::null()))
Remove unused points and renumber faces in local visit order.
Definition: triSurface.C:669
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:328
Foam::triSurface::triSurfInstance
static fileName triSurfInstance(const Time &)
Name of triSurface directory to use.
Definition: triSurface.C:73
MeshedSurface.H
Foam::PrimitivePatch<::Foam::List< labelledTri >, pointField >
Foam::triSurface::checkEdges
void checkEdges(const bool verbose)
Check triply (or more) connected edges.
Definition: triSurface.C:301