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 | Copyright (C) 2011-2013 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 "triSurface.H"
27 #include "demandDrivenData.H"
28 #include "IFstream.H"
29 #include "OFstream.H"
30 #include "Time.H"
31 #include "boundBox.H"
32 #include "SortableList.H"
33 #include "PackedBoolList.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 defineTypeNameAndDebug(triSurface, 0);
40 }
41 
42 
44 {
45  fileName foamName(d.caseName() + ".ftr");
46 
47  // Search back through the time directories list to find the time
48  // closest to and lower than current time
49 
50  instantList ts = d.times();
51  label i;
52 
53  for (i=ts.size()-1; i>=0; i--)
54  {
55  if (ts[i].value() <= d.timeOutputValue())
56  {
57  break;
58  }
59  }
60 
61  // Noting that the current directory has already been searched
62  // for mesh data, start searching from the previously stored time directory
63 
64  if (i>=0)
65  {
66  for (label j=i; j>=0; j--)
67  {
68  if (isFile(d.path()/ts[j].name()/typeName/foamName))
69  {
70  if (debug)
71  {
72  Pout<< " triSurface::triSurfInstance(const Time& d)"
73  << "reading " << foamName
74  << " from " << ts[j].name()/typeName
75  << endl;
76  }
77 
78  return ts[j].name();
79  }
80  }
81  }
82 
83  if (debug)
84  {
85  Pout<< " triSurface::triSurfInstance(const Time& d)"
86  << "reading " << foamName
87  << " from constant/" << endl;
88  }
89  return d.constant();
90 }
91 
92 
93 Foam::List<Foam::labelledTri> Foam::triSurface::convertToTri
94 (
95  const faceList& faces,
96  const label defaultRegion
97 )
98 {
99  List<labelledTri> triFaces(faces.size());
100 
101  forAll(triFaces, faceI)
102  {
103  const face& f = faces[faceI];
104 
105  if (f.size() != 3)
106  {
108  (
109  "triSurface::convertToTri"
110  "(const faceList&, const label)"
111  ) << "Face at position " << faceI
112  << " does not have three vertices:" << f
113  << abort(FatalError);
114  }
115 
116  labelledTri& tri = triFaces[faceI];
117 
118  tri[0] = f[0];
119  tri[1] = f[1];
120  tri[2] = f[2];
121  tri.region() = defaultRegion;
122  }
123 
124  return triFaces;
125 }
126 
127 
128 Foam::List<Foam::labelledTri> Foam::triSurface::convertToTri
129 (
130  const triFaceList& faces,
131  const label defaultRegion
132 )
133 {
134  List<labelledTri> triFaces(faces.size());
135 
136  forAll(triFaces, faceI)
137  {
138  const triFace& f = faces[faceI];
139 
140  labelledTri& tri = triFaces[faceI];
141 
142  tri[0] = f[0];
143  tri[1] = f[1];
144  tri[2] = f[2];
145  tri.region() = defaultRegion;
146  }
147 
148  return triFaces;
149 }
150 
151 
152 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
153 
154 void Foam::triSurface::printTriangle
155 (
156  Ostream& os,
157  const string& pre,
158  const labelledTri& f,
159  const pointField& points
160 )
161 {
162  os
163  << pre.c_str() << "vertex numbers:"
164  << f[0] << ' ' << f[1] << ' ' << f[2] << endl
165  << pre.c_str() << "vertex coords :"
166  << points[f[0]] << ' ' << points[f[1]] << ' ' << points[f[2]]
167  << pre.c_str() << "region :" << f.region() << endl
168  << endl;
169 }
170 
171 
172 Foam::string Foam::triSurface::getLineNoComment(IFstream& is)
173 {
174  string line;
175  do
176  {
177  is.getLine(line);
178  }
179  while ((line.empty() || line[0] == '#') && is.good());
180 
181  return line;
182 }
183 
184 
185 // Remove non-triangles, double triangles.
186 void Foam::triSurface::checkTriangles(const bool verbose)
187 {
188  // Simple check on indices ok.
189  const label maxPointI = points().size() - 1;
190 
191  forAll(*this, faceI)
192  {
193  const triSurface::FaceType& f = (*this)[faceI];
194 
195  forAll(f, fp)
196  {
197  if (f[fp] < 0 || f[fp] > maxPointI)
198  {
199  FatalErrorIn("triSurface::checkTriangles(bool)")
200  << "triangle " << f
201  << " uses point indices outside point range 0.."
202  << maxPointI
203  << exit(FatalError);
204  }
205  }
206  }
207 
208  // Two phase process
209  // 1. mark invalid faces
210  // 2. pack
211  // Done to keep numbering constant in phase 1
212 
213  // List of valid triangles
214  boolList valid(size(), true);
215  bool hasInvalid = false;
216 
217  forAll(*this, faceI)
218  {
219  const labelledTri& f = (*this)[faceI];
220 
221  if ((f[0] == f[1]) || (f[0] == f[2]) || (f[1] == f[2]))
222  {
223  // 'degenerate' triangle check
224  valid[faceI] = false;
225  hasInvalid = true;
226 
227  if (verbose)
228  {
229  WarningIn
230  (
231  "triSurface::checkTriangles(bool verbose)"
232  ) << "triangle " << faceI
233  << " does not have three unique vertices:\n";
234  printTriangle(Warning, " ", f, points());
235  }
236  }
237  else
238  {
239  // duplicate triangle check
240  const labelList& fEdges = faceEdges()[faceI];
241 
242  // Check if faceNeighbours use same points as this face.
243  // Note: discards normal information - sides of baffle are merged.
244 
245  forAll(fEdges, fp)
246  {
247  const labelList& eFaces = edgeFaces()[fEdges[fp]];
248 
249  forAll(eFaces, i)
250  {
251  label neighbour = eFaces[i];
252 
253  if (neighbour > faceI)
254  {
255  // lower numbered faces already checked
256  const labelledTri& n = (*this)[neighbour];
257 
258  if
259  (
260  ((f[0] == n[0]) || (f[0] == n[1]) || (f[0] == n[2]))
261  && ((f[1] == n[0]) || (f[1] == n[1]) || (f[1] == n[2]))
262  && ((f[2] == n[0]) || (f[2] == n[1]) || (f[2] == n[2]))
263  )
264  {
265  valid[faceI] = false;
266  hasInvalid = true;
267 
268  if (verbose)
269  {
270  WarningIn
271  (
272  "triSurface::checkTriangles(bool verbose)"
273  ) << "triangles share the same vertices:\n"
274  << " face 1 :" << faceI << endl;
275  printTriangle(Warning, " ", f, points());
276 
277  Warning
278  << endl
279  << " face 2 :"
280  << neighbour << endl;
281  printTriangle(Warning, " ", n, points());
282  }
283 
284  break;
285  }
286  }
287  }
288  }
289  }
290  }
291 
292  if (hasInvalid)
293  {
294  // Pack
295  label newFaceI = 0;
296  forAll(*this, faceI)
297  {
298  if (valid[faceI])
299  {
300  const labelledTri& f = (*this)[faceI];
301  (*this)[newFaceI++] = f;
302  }
303  }
304 
305  if (verbose)
306  {
307  WarningIn
308  (
309  "triSurface::checkTriangles(bool verbose)"
310  ) << "Removing " << size() - newFaceI
311  << " illegal faces." << endl;
312  }
313  (*this).setSize(newFaceI);
314 
315  // Topology can change because of renumbering
316  clearOut();
317  }
318 }
319 
320 
321 // Check/fix edges with more than two triangles
322 void Foam::triSurface::checkEdges(const bool verbose)
323 {
324  const labelListList& eFaces = edgeFaces();
325 
326  forAll(eFaces, edgeI)
327  {
328  const labelList& myFaces = eFaces[edgeI];
329 
330  if (myFaces.empty())
331  {
332  FatalErrorIn("triSurface::checkEdges(bool verbose)")
333  << "Edge " << edgeI << " with vertices " << edges()[edgeI]
334  << " has no edgeFaces"
335  << exit(FatalError);
336  }
337  else if (myFaces.size() > 2 && verbose)
338  {
339  WarningIn
340  (
341  "triSurface::checkEdges(bool verbose)"
342  ) << "Edge " << edgeI << " with vertices " << edges()[edgeI]
343  << " has more than 2 faces connected to it : " << myFaces
344  << endl;
345  }
346  }
347 }
348 
349 
350 // Read triangles, points from Istream
351 bool Foam::triSurface::read(Istream& is)
352 {
353  is >> patches_ >> storedPoints() >> storedFaces();
354 
355  return true;
356 }
357 
358 
359 // Read from file in given format
360 bool Foam::triSurface::read
361 (
362  const fileName& name,
363  const word& ext,
364  const bool check
365 )
366 {
367  if (check && !exists(name))
368  {
370  (
371  "triSurface::read(const fileName&, const word&, const bool)"
372  ) << "Cannnot read " << name << exit(FatalError);
373  }
374 
375  if (ext == "gz")
376  {
377  fileName unzipName = name.lessExt();
378 
379  // Do not check for existence. Let IFstream do the unzipping.
380  return read(unzipName, unzipName.ext(), false);
381  }
382  else if (ext == "ftr")
383  {
384  return read(IFstream(name)());
385  }
386  else if (ext == "stl")
387  {
388  return readSTL(name);
389  }
390  else if (ext == "stlb")
391  {
392  return readSTLBINARY(name);
393  }
394  else if (ext == "gts")
395  {
396  return readGTS(name);
397  }
398  else if (ext == "obj")
399  {
400  return readOBJ(name);
401  }
402  else if (ext == "off")
403  {
404  return readOFF(name);
405  }
406  else if (ext == "tri")
407  {
408  return readTRI(name);
409  }
410  else if (ext == "ac")
411  {
412  return readAC(name);
413  }
414  else if (ext == "nas")
415  {
416  return readNAS(name);
417  }
418  else if (ext == "vtk")
419  {
420  return readVTK(name);
421  }
422  else
423  {
425  (
426  "triSurface::read(const fileName&, const word&)"
427  ) << "unknown file extension " << ext
428  << ". Supported extensions are '.ftr', '.stl', '.stlb', '.gts'"
429  << ", '.obj', '.ac', '.off', '.nas', '.tri' and '.vtk'"
430  << exit(FatalError);
431 
432  return false;
433  }
434 }
435 
436 
437 // Write to file in given format
438 void Foam::triSurface::write
439 (
440  const fileName& name,
441  const word& ext,
442  const bool sort
443 ) const
444 {
445  if (ext == "ftr")
446  {
447  return write(OFstream(name)());
448  }
449  else if (ext == "stl")
450  {
451  return writeSTLASCII(sort, OFstream(name)());
452  }
453  else if (ext == "stlb")
454  {
455  ofstream outFile(name.c_str(), std::ios::binary);
456 
457  writeSTLBINARY(outFile);
458  }
459  else if (ext == "gts")
460  {
461  return writeGTS(sort, OFstream(name)());
462  }
463  else if (ext == "obj")
464  {
465  writeOBJ(sort, OFstream(name)());
466  }
467  else if (ext == "off")
468  {
469  writeOFF(sort, OFstream(name)());
470  }
471  else if (ext == "vtk")
472  {
473  writeVTK(sort, OFstream(name)());
474  }
475  else if (ext == "tri")
476  {
477  writeTRI(sort, OFstream(name)());
478  }
479  else if (ext == "dx")
480  {
481  writeDX(sort, OFstream(name)());
482  }
483  else if (ext == "ac")
484  {
485  writeAC(OFstream(name)());
486  }
487  else if (ext == "smesh")
488  {
489  writeSMESH(sort, OFstream(name)());
490  }
491  else
492  {
494  (
495  "triSurface::write(const fileName&, const word&, const bool)"
496  ) << "unknown file extension " << ext
497  << " for file " << name
498  << ". Supported extensions are '.ftr', '.stl', '.stlb', "
499  << "'.gts', '.obj', '.vtk'"
500  << ", '.off', '.dx', '.smesh', '.ac' and '.tri'"
501  << exit(FatalError);
502  }
503 }
504 
505 
506 // Returns patch info. Sets faceMap to the indexing according to patch
507 // numbers. Patch numbers start at 0.
508 Foam::surfacePatchList Foam::triSurface::calcPatches(labelList& faceMap) const
509 {
510  // Sort according to region numbers of labelledTri
511  SortableList<label> sortedRegion(size());
512 
513  forAll(sortedRegion, faceI)
514  {
515  sortedRegion[faceI] = operator[](faceI).region();
516  }
517  sortedRegion.sort();
518 
519  faceMap = sortedRegion.indices();
520 
521  // Extend regions
522  label maxRegion = patches_.size()-1; // for non-compacted regions
523 
524  if (faceMap.size())
525  {
526  maxRegion = max
527  (
528  maxRegion,
529  operator[](faceMap.last()).region()
530  );
531  }
532 
533  // Get new region list
534  surfacePatchList newPatches(maxRegion + 1);
535 
536  // Fill patch sizes
537  forAll(*this, faceI)
538  {
539  label region = operator[](faceI).region();
540 
541  newPatches[region].size()++;
542  }
543 
544 
545  // Fill rest of patch info
546 
547  label startFaceI = 0;
548  forAll(newPatches, newPatchI)
549  {
550  surfacePatch& newPatch = newPatches[newPatchI];
551 
552  newPatch.index() = newPatchI;
553 
554  label oldPatchI = newPatchI;
555 
556  // start of patch
557  newPatch.start() = startFaceI;
558 
559 
560  // Take over any information from existing patches
561  if ((oldPatchI < patches_.size()) && (patches_[oldPatchI].name() != ""))
562  {
563  newPatch.name() = patches_[oldPatchI].name();
564  }
565  else
566  {
567  newPatch.name() = word("patch") + name(newPatchI);
568  }
569 
570  if
571  (
572  (oldPatchI < patches_.size())
573  && (patches_[oldPatchI].geometricType() != "")
574  )
575  {
576  newPatch.geometricType() = patches_[oldPatchI].geometricType();
577  }
578  else
579  {
580  newPatch.geometricType() = "empty";
581  }
582 
583  startFaceI += newPatch.size();
584  }
585 
586  return newPatches;
587 }
588 
589 
590 void Foam::triSurface::setDefaultPatches()
591 {
593 
594  // Get names, types and sizes
595  surfacePatchList newPatches(calcPatches(faceMap));
596 
597  // Take over names and types (but not size)
598  patches_.setSize(newPatches.size());
599 
600  forAll(newPatches, patchI)
601  {
602  patches_[patchI].index() = patchI;
603  patches_[patchI].name() = newPatches[patchI].name();
604  patches_[patchI].geometricType() = newPatches[patchI].geometricType();
605  }
606 }
607 
608 
609 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
610 
612 :
614  patches_(0),
615  sortedEdgeFacesPtr_(NULL),
616  edgeOwnerPtr_(NULL)
617 {}
618 
619 
620 
622 (
623  const List<labelledTri>& triangles,
624  const geometricSurfacePatchList& patches,
625  const pointField& points
626 )
627 :
628  ParentType(triangles, points),
629  patches_(patches),
630  sortedEdgeFacesPtr_(NULL),
631  edgeOwnerPtr_(NULL)
632 {}
633 
634 
636 (
637  List<labelledTri>& triangles,
638  const geometricSurfacePatchList& patches,
639  pointField& points,
640  const bool reUse
641 )
642 :
643  ParentType(triangles, points, reUse),
644  patches_(patches),
645  sortedEdgeFacesPtr_(NULL),
646  edgeOwnerPtr_(NULL)
647 {}
648 
649 
651 (
652  const Xfer<List<labelledTri> >& triangles,
653  const geometricSurfacePatchList& patches,
654  const Xfer<List<point> >& points
655 )
656 :
657  ParentType(triangles, points),
658  patches_(patches),
659  sortedEdgeFacesPtr_(NULL),
660  edgeOwnerPtr_(NULL)
661 {}
662 
663 
665 (
666  const List<labelledTri>& triangles,
667  const pointField& points
668 )
669 :
670  ParentType(triangles, points),
671  patches_(),
672  sortedEdgeFacesPtr_(NULL),
673  edgeOwnerPtr_(NULL)
674 {
675  setDefaultPatches();
676 }
677 
678 
680 (
681  const triFaceList& triangles,
682  const pointField& points
683 )
684 :
685  ParentType(convertToTri(triangles, 0), points),
686  patches_(),
687  sortedEdgeFacesPtr_(NULL),
688  edgeOwnerPtr_(NULL)
689 {
690  setDefaultPatches();
691 }
692 
693 
695 :
697  patches_(),
698  sortedEdgeFacesPtr_(NULL),
699  edgeOwnerPtr_(NULL)
700 {
701  word ext = name.ext();
702 
703  read(name, ext);
704 
705  setDefaultPatches();
706 }
707 
708 
710 :
712  patches_(),
713  sortedEdgeFacesPtr_(NULL),
714  edgeOwnerPtr_(NULL)
715 {
716  read(is);
717 
718  setDefaultPatches();
719 }
720 
721 
723 :
725  patches_(),
726  sortedEdgeFacesPtr_(NULL),
727  edgeOwnerPtr_(NULL)
728 {
729  fileName foamFile(d.caseName() + ".ftr");
730 
731  fileName foamPath(d.path()/triSurfInstance(d)/typeName/foamFile);
732 
733  IFstream foamStream(foamPath);
734 
735  read(foamStream);
736 
737  setDefaultPatches();
738 }
739 
740 
742 :
743  ParentType(ts, ts.points()),
744  patches_(ts.patches()),
745  sortedEdgeFacesPtr_(NULL),
746  edgeOwnerPtr_(NULL)
747 {}
748 
749 
750 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
751 
753 {
754  clearOut();
755 }
756 
757 
758 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
759 
761 {
763  deleteDemandDrivenData(sortedEdgeFacesPtr_);
764  deleteDemandDrivenData(edgeOwnerPtr_);
765 }
766 
767 
769 {
771 }
772 
773 
775 {
777 
778  clearTopology();
780 }
781 
782 
784 {
785  if (!sortedEdgeFacesPtr_)
786  {
787  calcSortedEdgeFaces();
788  }
789 
790  return *sortedEdgeFacesPtr_;
791 }
792 
793 
795 {
796  if (!edgeOwnerPtr_)
797  {
798  calcEdgeOwner();
799  }
800 
801  return *edgeOwnerPtr_;
802 }
803 
804 
806 {
807  // Remove all geometry dependent data
808  deleteDemandDrivenData(sortedEdgeFacesPtr_);
809 
810  // Adapt for new point position
811  ParentType::movePoints(newPoints);
812 
813  // Copy new points
814  storedPoints() = newPoints;
815 }
816 
817 
818 void Foam::triSurface::scalePoints(const scalar scaleFactor)
819 {
820  // avoid bad scaling
821  if (scaleFactor > 0 && scaleFactor != 1.0)
822  {
823  // Remove all geometry dependent data
824  clearTopology();
825 
826  // Adapt for new point position
828 
829  storedPoints() *= scaleFactor;
830  }
831 }
832 
833 
834 // Remove non-triangles, double triangles.
835 void Foam::triSurface::cleanup(const bool verbose)
836 {
837  // Merge points (already done for STL, TRI)
838  stitchTriangles(SMALL, verbose);
839 
840  // Merging points might have changed geometric factors
841  clearOut();
842 
843  checkTriangles(verbose);
844 
845  checkEdges(verbose);
846 }
847 
848 
849 // Finds area, starting at faceI, delimited by borderEdge. Marks all visited
850 // faces (from face-edge-face walk) with currentZone.
852 (
853  const boolList& borderEdge,
854  const label faceI,
855  const label currentZone,
856  labelList& faceZone
857 ) const
858 {
859  // List of faces whose faceZone has been set.
860  labelList changedFaces(1, faceI);
861 
862  while (true)
863  {
864  // Pick up neighbours of changedFaces
865  DynamicList<label> newChangedFaces(2*changedFaces.size());
866 
867  forAll(changedFaces, i)
868  {
869  label faceI = changedFaces[i];
870 
871  const labelList& fEdges = faceEdges()[faceI];
872 
873  forAll(fEdges, i)
874  {
875  label edgeI = fEdges[i];
876 
877  if (!borderEdge[edgeI])
878  {
879  const labelList& eFaces = edgeFaces()[edgeI];
880 
881  forAll(eFaces, j)
882  {
883  label nbrFaceI = eFaces[j];
884 
885  if (faceZone[nbrFaceI] == -1)
886  {
887  faceZone[nbrFaceI] = currentZone;
888  newChangedFaces.append(nbrFaceI);
889  }
890  else if (faceZone[nbrFaceI] != currentZone)
891  {
893  (
894  "triSurface::markZone(const boolList&,"
895  "const label, const label, labelList&) const"
896  )
897  << "Zones " << faceZone[nbrFaceI]
898  << " at face " << nbrFaceI
899  << " connects to zone " << currentZone
900  << " at face " << faceI
901  << abort(FatalError);
902  }
903  }
904  }
905  }
906  }
907 
908  if (newChangedFaces.empty())
909  {
910  break;
911  }
912 
913  changedFaces.transfer(newChangedFaces);
914  }
915 }
916 
917 
918 // Finds areas delimited by borderEdge (or 'real' edges).
919 // Fills faceZone accordingly
921 (
922  const boolList& borderEdge,
923  labelList& faceZone
924 ) const
925 {
926  faceZone.setSize(size());
927  faceZone = -1;
928 
929  if (borderEdge.size() != nEdges())
930  {
932  (
933  "triSurface::markZones"
934  "(const boolList&, labelList&)"
935  )
936  << "borderEdge boolList not same size as number of edges" << endl
937  << "borderEdge:" << borderEdge.size() << endl
938  << "nEdges :" << nEdges()
939  << exit(FatalError);
940  }
941 
942  label zoneI = 0;
943 
944  for (label startFaceI = 0;; zoneI++)
945  {
946  // Find first non-coloured face
947  for (; startFaceI < size(); startFaceI++)
948  {
949  if (faceZone[startFaceI] == -1)
950  {
951  break;
952  }
953  }
954 
955  if (startFaceI >= size())
956  {
957  break;
958  }
959 
960  faceZone[startFaceI] = zoneI;
961 
962  markZone(borderEdge, startFaceI, zoneI, faceZone);
963  }
964 
965  return zoneI;
966 }
967 
968 
970 (
971  const boolList& include,
972  labelList& pointMap,
973  labelList& faceMap
974 ) const
975 {
976  const List<labelledTri>& locFaces = localFaces();
977 
978  label faceI = 0;
979  label pointI = 0;
980 
981  faceMap.setSize(locFaces.size());
982 
983  pointMap.setSize(nPoints());
984 
985  boolList pointHad(nPoints(), false);
986 
987  forAll(include, oldFaceI)
988  {
989  if (include[oldFaceI])
990  {
991  // Store new faces compact
992  faceMap[faceI++] = oldFaceI;
993 
994  // Renumber labels for face
995  const triSurface::FaceType& f = locFaces[oldFaceI];
996 
997  forAll(f, fp)
998  {
999  label labI = f[fp];
1000  if (!pointHad[labI])
1001  {
1002  pointHad[labI] = true;
1003  pointMap[pointI++] = labI;
1004  }
1005  }
1006  }
1007  }
1008 
1009  // Trim
1010  faceMap.setSize(faceI);
1011  pointMap.setSize(pointI);
1012 }
1013 
1014 
1017  const boolList& include,
1018  labelList& pointMap,
1019  labelList& faceMap
1020 ) const
1021 {
1022  const pointField& locPoints = localPoints();
1023  const List<labelledTri>& locFaces = localFaces();
1024 
1025  // Fill pointMap, faceMap
1026  subsetMeshMap(include, pointMap, faceMap);
1027 
1028 
1029  // Create compact coordinate list and forward mapping array
1030  pointField newPoints(pointMap.size());
1031  labelList oldToNew(locPoints.size());
1032  forAll(pointMap, pointi)
1033  {
1034  newPoints[pointi] = locPoints[pointMap[pointi]];
1035  oldToNew[pointMap[pointi]] = pointi;
1036  }
1037 
1038  // Renumber triangle node labels and compact
1039  List<labelledTri> newTriangles(faceMap.size());
1040 
1041  forAll(faceMap, facei)
1042  {
1043  // Get old vertex labels
1044  const labelledTri& tri = locFaces[faceMap[facei]];
1045 
1046  newTriangles[facei][0] = oldToNew[tri[0]];
1047  newTriangles[facei][1] = oldToNew[tri[1]];
1048  newTriangles[facei][2] = oldToNew[tri[2]];
1049  newTriangles[facei].region() = tri.region();
1050  }
1051 
1052  // Construct subsurface
1053  return triSurface(newTriangles, patches(), newPoints, true);
1054 }
1055 
1056 
1057 void Foam::triSurface::write
1059  const fileName& name,
1060  const bool sortByRegion
1061 ) const
1062 {
1063  write(name, name.ext(), sortByRegion);
1064 }
1065 
1066 
1067 void Foam::triSurface::write(Ostream& os) const
1068 {
1069  os << patches() << endl;
1070 
1071  //Note: Write with global point numbering
1072  os << points() << nl
1073  << static_cast<const List<labelledTri>&>(*this) << endl;
1074 
1075  // Check state of Ostream
1076  os.check("triSurface::write(Ostream&)");
1077 }
1078 
1079 
1080 void Foam::triSurface::write(const Time& d) const
1081 {
1082  fileName foamFile(d.caseName() + ".ftr");
1083 
1084  fileName foamPath(d.path()/triSurfInstance(d)/typeName/foamFile);
1085 
1086  OFstream foamStream(foamPath);
1087 
1088  write(foamStream);
1089 }
1090 
1091 
1093 {
1094  // Unfortunately nPoints constructs meshPoints() so do compact version
1095  // ourselves.
1096  PackedBoolList pointIsUsed(points().size());
1097 
1098  label nPoints = 0;
1100 
1101  forAll(*this, faceI)
1102  {
1103  const triSurface::FaceType& f = operator[](faceI);
1104 
1105  forAll(f, fp)
1106  {
1107  label pointI = f[fp];
1108  if (pointIsUsed.set(pointI, 1))
1109  {
1110  bb.min() = ::Foam::min(bb.min(), points()[pointI]);
1111  bb.max() = ::Foam::max(bb.max(), points()[pointI]);
1112  nPoints++;
1113  }
1114  }
1115  }
1116 
1117  os << "Triangles : " << size() << endl
1118  << "Vertices : " << nPoints << endl
1119  << "Bounding Box : " << bb << endl;
1120 }
1121 
1122 
1123 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1124 
1126 {
1128  clearOut();
1129  storedPoints() = ts.points();
1130  patches_ = ts.patches();
1131 }
1132 
1133 
1134 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
1135 
1137 {
1138  sm.write(os);
1139  return os;
1140 }
1141 
1142 
1143 // ************************************************************************* //
Output to file stream.
Definition: OFstream.H:81
const pointField & points
word ext() const
Return file name extension (part after last .)
Definition: fileName.C:323
word name() const
Return file name (part beyond last /)
Definition: fileName.C:206
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
const labelList & indices() const
Return the list of sorted indices. Updated every sort.
Definition: SortableList.H:90
static const boundBox invertedBox
A very large inverted boundBox: min/max == +/- VGREAT.
Definition: boundBox.H:79
void sort(UList< T > &)
Definition: UList.C:107
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
fileName path() const
Return path.
Definition: Time.H:281
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:390
const labelList & edgeOwner() const
If 2 face neighbours: label of face where ordering of edge.
Definition: triSurface.C:794
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
labelList f(nPoints)
const labelListList & edgeFaces() const
Return edge-face addressing.
void writeVTK(OFstream &os, const Type &value)
bool empty() const
Return true if the UList is empty (ie, size() is zero).
Definition: UListI.H:313
void clearPatchMeshAddr()
Definition: triSurface.C:768
void clearTopology()
Definition: triSurface.C:760
T & last()
Return the last element of the list.
Definition: UListI.H:131
void set(const PackedList< 1 > &)
Set specified bits.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
A simple container for copying or transferring objects of type <T>.
Definition: Xfer.H:85
Triangulated surface description with patch information.
Definition: triSurface.H:57
void deleteDemandDrivenData(DataPtr &dataPtr)
A bit-packed bool list.
label markZones(const boolList &borderEdge, labelList &faceZone) const
(size and) fills faceZone with zone of face. Zone is area
Definition: triSurface.C:921
A class for handling words, derived from string.
Definition: word.H:59
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A class for handling character strings derived from std::string.
Definition: string.H:74
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
triSurface()
Construct null.
Definition: triSurface.C:611
bool exists(const fileName &, const bool checkGzip=true)
Does the name exist (as DIRECTORY or FILE) in the file system?
Definition: POSIX.C:609
const word & geometricType() const
Return the type of the patch.
Ostream & operator<<(Ostream &, const edgeMesh &)
Definition: edgeMeshIO.C:133
Input from file stream.
Definition: IFstream.H:81
scalar timeOutputValue() const
Return current time value.
Definition: TimeState.C:67
void sort()
(stable) sort the list (if changed after construction time)
Definition: SortableList.C:95
label nPoints() const
Return number of points supporting patch faces.
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Namespace for OpenFOAM.
A line primitive.
Definition: line.H:56
runTime write()
messageStream Warning
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
label n
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
void checkEdges(const bool verbose)
Check triply (or more) connected edges.
Definition: triSurface.C:322
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
T & operator[](const label)
Return element of UList.
Definition: UListI.H:163
const Field< point > & localPoints() const
Return pointField of points in patch.
Triangle with additional region number.
Definition: labelledTri.H:49
static const char nl
Definition: Ostream.H:260
label region() const
Return region label.
Definition: labelledTriI.H:68
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const geometricSurfacePatchList & patches() const
Definition: triSurface.H:306
#define WarningIn(functionName)
Report a warning using Foam::Warning.
A list that is sorted upon construction or when explicitly requested with the sort() method...
Definition: List.H:65
void subsetMeshMap(const boolList &include, labelList &pointMap, labelList &faceMap) const
&#39;Create&#39; sub mesh, including only faces for which
Definition: triSurface.C:970
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define forAll(list, i)
Definition: UList.H:421
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
triSurface subsetMesh(const boolList &include, labelList &pointMap, labelList &faceMap) const
Return new surface. Returns pointMap, faceMap from.
Definition: triSurface.C:1016
label start() const
Return start label of this patch in the polyMesh face list.
Definition: surfacePatch.H:102
Template functions to aid in the implementation of demand driven data.
void operator=(const triSurface &)
Definition: triSurface.C:1125
const point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
errorManip< error > abort(error &err)
Definition: errorManip.H:131
ISstream & getLine(string &)
Raw, low-level getline into a string function.
Definition: ISstreamI.H:77
pointField & storedPoints()
Non-const access to global points.
Definition: triSurface.H:212
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:88
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:55
static fileName triSurfInstance(const Time &)
Name of triSurface directory to use.
Definition: triSurface.C:43
const word & name() const
Return name.
virtual const fileName & name() const
Return the name of the stream.
Definition: IOstream.H:297
error FatalError
const labelListList & sortedEdgeFaces() const
Return edge-face addressing sorted (for edges with more than.
Definition: triSurface.C:783
A class for handling file names.
Definition: fileName.H:69
const labelListList & faceEdges() const
Return face-edge addressing.
label size() const
Return the number of elements in the UList.
fileName lessExt() const
Return file name without extension (part before last .)
Definition: fileName.C:307
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:68
virtual void movePoints(const pointField &)
Move points.
Definition: triSurface.C:805
void checkTriangles(const bool verbose)
Check/remove duplicate/degenerate triangles.
Definition: triSurface.C:186
bool read(const char *, int32_t &)
Definition: int32IO.C:87
const List< labelledTri > & localFaces() const
Return patch faces addressing into local point list.
label size() const
Return size of this patch in the polyMesh face list.
Definition: surfacePatch.H:114
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:97
instantList times() const
Search the case for valid time directories.
Definition: Time.C:758
void cleanup(const bool verbose)
Remove non-valid triangles.
Definition: triSurface.C:835
~triSurface()
Destructor.
Definition: triSurface.C:752
bool isFile(const fileName &, const bool checkGzip=true)
Does the name exist as a FILE in the file system?
Definition: POSIX.C:623
const fileName & caseName() const
Return case name.
Definition: Time.H:275
virtual void movePoints(const Field< point > &)
Correct patch after moving points.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
label nEdges() const
Return number of edges in patch.
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:209
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:852
void writeStats(Ostream &) const
Write some statistics.
Definition: triSurface.C:1092
label index() const
Return the index of this patch in the boundaryMesh.
const Field< point > & points() const
Return reference to global points.
virtual void scalePoints(const scalar)
Scale points. A non-positive factor is ignored.
Definition: triSurface.C:818
defineTypeNameAndDebug(combustionModel, 0)
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
&#39;Patch&#39; on surface as subset of triSurface.
Definition: surfacePatch.H:51
void operator=(const UList< T > &)
Assignment from UList operator. Takes linear time.
Definition: List.C:429