extendedEdgeMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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 "extendedEdgeMesh.H"
27 #include "surfaceFeatures.H"
28 #include "triSurface.H"
29 #include "Random.H"
30 #include "Time.H"
31 #include "OBJstream.H"
32 #include "DynamicField.H"
33 #include "edgeMeshFormatsCore.H"
34 #include "IOmanip.H"
35 
36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40  defineTypeNameAndDebug(extendedEdgeMesh, 0);
41 
42  template<>
43  const char* Foam::NamedEnum
44  <
46  4
47  >::names[] =
48  {
49  "convex",
50  "concave",
51  "mixed",
52  "nonFeature"
53  };
54 
55  template<>
56  const char* Foam::NamedEnum
57  <
59  6
60  >::names[] =
61  {
62  "external",
63  "internal",
64  "flat",
65  "open",
66  "multiple",
67  "none"
68  };
69 
70  template<>
71  const char* Foam::NamedEnum
72  <
74  4
75  >::names[] =
76  {
77  "inside",
78  "outside",
79  "both",
80  "neither"
81  };
82 }
83 
86 
89 
92 
94  Foam::cos(degToRad(0.1));
95 
96 
97 Foam::label Foam::extendedEdgeMesh::convexStart_ = 0;
98 
99 
100 Foam::label Foam::extendedEdgeMesh::externalStart_ = 0;
101 
102 
103 Foam::label Foam::extendedEdgeMesh::nPointTypes = 4;
104 
105 
106 Foam::label Foam::extendedEdgeMesh::nEdgeTypes = 5;
107 
108 
110 {
111  return wordHashSet(*fileExtensionConstructorTablePtr_);
112 }
113 
114 
116 {
117  return wordHashSet(*writefileExtensionMemberFunctionTablePtr_);
118 }
119 
120 
121 
122 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
123 
125 (
126  const word& ext,
127  const bool verbose
128 )
129 {
130  return edgeMeshFormatsCore::checkSupport
131  (
132  readTypes(),
133  ext,
134  verbose,
135  "reading"
136  );
137 }
138 
139 
141 (
142  const word& ext,
143  const bool verbose
144 )
145 {
146  return edgeMeshFormatsCore::checkSupport
147  (
148  writeTypes(),
149  ext,
150  verbose,
151  "writing"
152  );
153 }
154 
155 
157 (
158  const fileName& name,
159  const bool verbose
160 )
161 {
162  word ext = name.ext();
163  if (ext == "gz")
164  {
165  ext = name.lessExt().ext();
166  }
167  return canReadType(ext, verbose);
168 }
169 
170 
171 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
172 
175 (
176  label ptI
177 ) const
178 {
179  const labelList& ptEds(pointEdges()[ptI]);
180 
181  label nPtEds = ptEds.size();
182  label nExternal = 0;
183  label nInternal = 0;
184 
185  if (nPtEds == 0)
186  {
187  // There are no edges attached to the point, this is a problem
188  return NONFEATURE;
189  }
190 
191  forAll(ptEds, i)
192  {
193  edgeStatus edStat = getEdgeStatus(ptEds[i]);
194 
195  if (edStat == EXTERNAL)
196  {
197  nExternal++;
198  }
199  else if (edStat == INTERNAL)
200  {
201  nInternal++;
202  }
203  }
204 
205  if (nExternal == nPtEds)
206  {
207  return CONVEX;
208  }
209  else if (nInternal == nPtEds)
210  {
211  return CONCAVE;
212  }
213  else
214  {
215  return MIXED;
216  }
217 }
218 
219 
222 (
223  const List<vector>& norms,
224  const labelList& edNorms,
225  const vector& fC0tofC1
226 )
227 {
228  label nEdNorms = edNorms.size();
229 
230  if (nEdNorms == 1)
231  {
232  return OPEN;
233  }
234  else if (nEdNorms == 2)
235  {
236  const vector& n0(norms[edNorms[0]]);
237  const vector& n1(norms[edNorms[1]]);
238 
239  if ((n0 & n1) > cosNormalAngleTol_)
240  {
241  return FLAT;
242  }
243  else if ((fC0tofC1 & n0) > 0.0)
244  {
245  return INTERNAL;
246  }
247  else
248  {
249  return EXTERNAL;
250  }
251  }
252  else if (nEdNorms > 2)
253  {
254  return MULTIPLE;
255  }
256  else
257  {
258  // There is a problem - the edge has no normals
259  return NONE;
260  }
261 }
262 
263 
264 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
265 
267 :
268  edgeMesh(pointField(0), edgeList(0)),
269  concaveStart_(0),
270  mixedStart_(0),
271  nonFeatureStart_(0),
272  internalStart_(0),
273  flatStart_(0),
274  openStart_(0),
275  multipleStart_(0),
276  normals_(0),
277  normalVolumeTypes_(0),
278  edgeDirections_(0),
279  normalDirections_(0),
280  edgeNormals_(0),
281  featurePointNormals_(0),
282  featurePointEdges_(0),
283  regionEdges_(0),
284  pointTree_(),
285  edgeTree_(),
286  edgeTreesByType_()
287 {}
288 
289 
291 :
292  edgeMesh(fem),
293  concaveStart_(fem.concaveStart()),
294  mixedStart_(fem.mixedStart()),
295  nonFeatureStart_(fem.nonFeatureStart()),
296  internalStart_(fem.internalStart()),
297  flatStart_(fem.flatStart()),
298  openStart_(fem.openStart()),
299  multipleStart_(fem.multipleStart()),
300  normals_(fem.normals()),
301  normalVolumeTypes_(fem.normalVolumeTypes()),
302  edgeDirections_(fem.edgeDirections()),
303  normalDirections_(fem.normalDirections()),
304  edgeNormals_(fem.edgeNormals()),
305  featurePointNormals_(fem.featurePointNormals()),
306  featurePointEdges_(fem.featurePointEdges()),
307  regionEdges_(fem.regionEdges()),
308  pointTree_(),
309  edgeTree_(),
311 {}
312 
313 
315 :
316  edgeMesh(move(fem)),
317  concaveStart_(fem.concaveStart()),
318  mixedStart_(fem.mixedStart()),
319  nonFeatureStart_(fem.nonFeatureStart()),
320  internalStart_(fem.internalStart()),
321  flatStart_(fem.flatStart()),
322  openStart_(fem.openStart()),
323  multipleStart_(fem.multipleStart()),
324  normals_(move(fem.normals())),
325  normalVolumeTypes_(move(fem.normalVolumeTypes())),
326  edgeDirections_(move(fem.edgeDirections())),
327  normalDirections_(move(fem.normalDirections())),
328  edgeNormals_(move(fem.edgeNormals())),
329  featurePointNormals_(move(fem.featurePointNormals())),
330  featurePointEdges_(move(fem.featurePointEdges())),
331  regionEdges_(move(fem.regionEdges())),
332  pointTree_(),
333  edgeTree_(),
335 {}
336 
337 
339 {
340  is >> *this;
341 }
342 
343 
345 (
346  pointField&& pointLst,
347  edgeList&& edgeLst
348 )
349 :
350  edgeMesh(move(pointLst), move(edgeLst)),
351  concaveStart_(0),
352  mixedStart_(0),
353  nonFeatureStart_(0),
354  internalStart_(0),
355  flatStart_(0),
356  openStart_(0),
357  multipleStart_(0),
358  normals_(0),
360  edgeDirections_(0),
362  edgeNormals_(0),
365  regionEdges_(0),
366  pointTree_(),
367  edgeTree_(),
369 {}
370 
371 
373 (
374  const surfaceFeatures& sFeat,
375  const boolList& surfBaffleRegions
376 )
377 :
378  edgeMesh(pointField(0), edgeList(0)),
379  concaveStart_(-1),
380  mixedStart_(-1),
381  nonFeatureStart_(-1),
382  internalStart_(-1),
383  flatStart_(-1),
384  openStart_(-1),
385  multipleStart_(-1),
386  normals_(0),
388  edgeDirections_(0),
390  edgeNormals_(0),
393  regionEdges_(0),
394  pointTree_(),
395  edgeTree_(),
397 {
398  // Extract and reorder the data from surfaceFeatures
399  const triSurface& surf = sFeat.surface();
400  const labelList& featureEdges = sFeat.featureEdges();
401  const labelList& featurePoints = sFeat.featurePoints();
402 
403  // Get a labelList of all the featureEdges that are region edges
404  const labelList regionFeatureEdges(identity(sFeat.nRegionEdges()));
405 
407  (
408  surf,
409  featureEdges,
410  regionFeatureEdges,
411  featurePoints
412  );
413 
414  const labelListList& edgeFaces = surf.edgeFaces();
415 
417 
418  // Noting when the normal of a face has been used so not to duplicate
419  labelList faceMap(surf.size(), -1);
420 
421  label nAdded = 0;
422 
423  forAll(featureEdges, i)
424  {
425  label sFEI = featureEdges[i];
426 
427  // Pick up the faces adjacent to the feature edge
428  const labelList& eFaces = edgeFaces[sFEI];
429 
430  forAll(eFaces, j)
431  {
432  label eFI = eFaces[j];
433 
434  // Check to see if the points have been already used
435  if (faceMap[eFI] == -1)
436  {
437  normalVolumeTypes_[nAdded++] =
438  (
439  surfBaffleRegions[surf[eFI].region()]
440  ? BOTH
441  : INSIDE
442  );
443 
444  faceMap[eFI] = nAdded - 1;
445  }
446  }
447  }
448 }
449 
450 
452 (
454  const labelList& featureEdges,
455  const labelList& regionFeatureEdges,
456  const labelList& featurePoints
457 )
458 :
459  edgeMesh(pointField(0), edgeList(0)),
460  concaveStart_(-1),
461  mixedStart_(-1),
462  nonFeatureStart_(-1),
463  internalStart_(-1),
464  flatStart_(-1),
465  openStart_(-1),
466  multipleStart_(-1),
467  normals_(0),
469  edgeDirections_(0),
471  edgeNormals_(0),
474  regionEdges_(0),
475  pointTree_(),
476  edgeTree_(),
478 {
480  (
481  surf,
482  featureEdges,
483  regionFeatureEdges,
484  featurePoints
485  );
486 }
487 
488 
490 (
491  const pointField& pts,
492  const edgeList& eds,
493  label concaveStart,
494  label mixedStart,
495  label nonFeatureStart,
496  label internalStart,
497  label flatStart,
498  label openStart,
499  label multipleStart,
500  const vectorField& normals,
501  const List<sideVolumeType>& normalVolumeTypes,
502  const vectorField& edgeDirections,
503  const labelListList& normalDirections,
504  const labelListList& edgeNormals,
505  const labelListList& featurePointNormals,
506  const labelListList& featurePointEdges,
507  const labelList& regionEdges
508 )
509 :
510  edgeMesh(pts, eds),
511  concaveStart_(concaveStart),
512  mixedStart_(mixedStart),
513  nonFeatureStart_(nonFeatureStart),
514  internalStart_(internalStart),
515  flatStart_(flatStart),
516  openStart_(openStart),
517  multipleStart_(multipleStart),
518  normals_(normals),
519  normalVolumeTypes_(normalVolumeTypes),
520  edgeDirections_(edgeDirections),
521  normalDirections_(normalDirections),
522  edgeNormals_(edgeNormals),
523  featurePointNormals_(featurePointNormals),
524  featurePointEdges_(featurePointEdges),
525  regionEdges_(regionEdges),
526  pointTree_(),
527  edgeTree_(),
529 {}
530 
531 
533 (
534  const fileName& name,
535  const word& ext
536 )
537 :
538  edgeMesh(pointField(0), edgeList(0)),
539  concaveStart_(0),
540  mixedStart_(0),
541  nonFeatureStart_(0),
542  internalStart_(0),
543  flatStart_(0),
544  openStart_(0),
545  multipleStart_(0),
546  normals_(0),
548  edgeDirections_(0),
550  edgeNormals_(0),
553  regionEdges_(0),
554  pointTree_(),
555  edgeTree_(),
557 {
558  read(name, ext);
559 }
560 
561 
563 :
564  edgeMesh(pointField(0), edgeList(0)),
565  concaveStart_(0),
566  mixedStart_(0),
567  nonFeatureStart_(0),
568  internalStart_(0),
569  flatStart_(0),
570  openStart_(0),
571  multipleStart_(0),
572  normals_(0),
574  edgeDirections_(0),
576  edgeNormals_(0),
579  regionEdges_(0),
580  pointTree_(),
581  edgeTree_(),
583 {
584  read(name);
585 }
586 
587 
588 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
589 
591 {}
592 
593 
594 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
595 
597 {
598  word ext = name.ext();
599  if (ext == "gz")
600  {
601  fileName unzipName = name.lessExt();
602  return read(unzipName, unzipName.ext());
603  }
604  else
605  {
606  return read(name, ext);
607  }
608 }
609 
610 
611 // Read from file in given format
613 (
614  const fileName& name,
615  const word& ext
616 )
617 {
618  // read via selector mechanism
619  transfer(New(name, ext)());
620  return true;
621 }
622 
623 
625 (
626  const point& sample,
627  scalar searchDistSqr,
628  pointIndexHit& info
629 ) const
630 {
631  info = pointTree().findNearest
632  (
633  sample,
634  searchDistSqr
635  );
636 }
637 
638 
640 (
641  const point& sample,
642  scalar searchDistSqr,
643  pointIndexHit& info
644 ) const
645 {
646  info = edgeTree().findNearest
647  (
648  sample,
649  searchDistSqr
650  );
651 }
652 
653 
655 (
656  const pointField& samples,
657  const scalarField& searchDistSqr,
658  pointIndexHitList& info
659 ) const
660 {
661  info.setSize(samples.size());
662 
663  forAll(samples, i)
664  {
666  (
667  samples[i],
668  searchDistSqr[i],
669  info[i]
670  );
671  }
672 }
673 
674 
676 (
677  const point& sample,
678  const scalarField& searchDistSqr,
679  pointIndexHitList& info
680 ) const
681 {
683 
684  info.setSize(edgeTrees.size());
685 
686  labelList sliceStarts(edgeTrees.size());
687 
688  sliceStarts[0] = externalStart_;
689  sliceStarts[1] = internalStart_;
690  sliceStarts[2] = flatStart_;
691  sliceStarts[3] = openStart_;
692  sliceStarts[4] = multipleStart_;
693 
694  forAll(edgeTrees, i)
695  {
696  info[i] = edgeTrees[i].findNearest
697  (
698  sample,
699  searchDistSqr[i]
700  );
701 
702  // The index returned by the indexedOctree is local to the slice of
703  // edges it was supplied with, return the index to the value in the
704  // complete edge list
705 
706  info[i].setIndex(info[i].index() + sliceStarts[i]);
707  }
708 }
709 
710 
712 (
713  const point& sample,
714  scalar searchRadiusSqr,
715  pointIndexHitList& info
716 ) const
717 {
718  // Pick up all the feature points that intersect the search sphere
719  labelList elems = pointTree().findSphere
720  (
721  sample,
722  searchRadiusSqr
723  );
724 
725  DynamicList<pointIndexHit> dynPointHit(elems.size());
726 
727  forAll(elems, elemI)
728  {
729  label index = elems[elemI];
730  label ptI = pointTree().shapes().pointLabels()[index];
731  const point& pt = points()[ptI];
732 
733  pointIndexHit nearHit(true, pt, index);
734 
735  dynPointHit.append(nearHit);
736  }
737 
738  info.transfer(dynPointHit);
739 }
740 
741 
743 (
744  const point& sample,
745  const scalar searchRadiusSqr,
746  pointIndexHitList& info
747 ) const
748 {
750 
751  info.setSize(edgeTrees.size());
752 
753  labelList sliceStarts(edgeTrees.size());
754 
755  sliceStarts[0] = externalStart_;
756  sliceStarts[1] = internalStart_;
757  sliceStarts[2] = flatStart_;
758  sliceStarts[3] = openStart_;
759  sliceStarts[4] = multipleStart_;
760 
761  DynamicList<pointIndexHit> dynEdgeHit(edgeTrees.size()*3);
762 
763  // Loop over all the feature edge types
764  forAll(edgeTrees, i)
765  {
766  // Pick up all the edges that intersect the search sphere
767  labelList elems = edgeTrees[i].findSphere
768  (
769  sample,
770  searchRadiusSqr
771  );
772 
773  forAll(elems, elemI)
774  {
775  label index = elems[elemI];
776  label edgeI = edgeTrees[i].shapes().edgeLabels()[index];
777  const edge& e = edges()[edgeI];
778 
779  pointHit hitPoint = e.line(points()).nearestDist(sample);
780 
781  label hitIndex = index + sliceStarts[i];
782 
783  pointIndexHit nearHit
784  (
785  hitPoint.hit(),
786  hitPoint.rawPoint(),
787  hitIndex
788  );
789 
790  dynEdgeHit.append(nearHit);
791  }
792  }
793 
794  info.transfer(dynEdgeHit);
795 }
796 
797 
799 (
800  const pointIndexHitList& hitList
801 ) const
802 {
803  scalar minDist = GREAT;
804 
805  for
806  (
807  label hi1 = 0;
808  hi1 < hitList.size() - 1;
809  ++hi1
810  )
811  {
812  const pointIndexHit& pHit1 = hitList[hi1];
813 
814  if (pHit1.hit())
815  {
816  const edge& e1 = edges()[pHit1.index()];
817 
818  for
819  (
820  label hi2 = hi1 + 1;
821  hi2 < hitList.size();
822  ++hi2
823  )
824  {
825  const pointIndexHit& pHit2 = hitList[hi2];
826 
827  if (pHit2.hit())
828  {
829  const edge& e2 = edges()[pHit2.index()];
830 
831  // Don't refine if the edges are connected to each other
832  if (!e1.connected(e2))
833  {
834  scalar curDist =
835  mag(pHit1.hitPoint() - pHit2.hitPoint());
836 
837  minDist = min(curDist, minDist);
838  }
839  }
840  }
841  }
842  }
843 
844  return minDist;
845 }
846 
847 
850 {
851  if (pointTree_.empty())
852  {
853  // Slightly extended bb. Slightly off-centred just so on symmetric
854  // geometry there are less face/edge aligned items.
855  treeBoundBox bb(treeBoundBox(points()).extend(1e-4));
856 
857  const labelList featurePointLabels = identity(nonFeatureStart_);
858 
859  pointTree_.reset
860  (
862  (
864  (
865  points(),
866  featurePointLabels
867  ),
868  bb, // bb
869  8, // maxLevel
870  10, // leafsize
871  3.0 // duplicity
872  )
873  );
874  }
875 
876  return pointTree_();
877 }
878 
879 
882 {
883  if (edgeTree_.empty())
884  {
885  // Slightly extended bb. Slightly off-centred just so on symmetric
886  // geometry there are less face/edge aligned items.
887  treeBoundBox bb(treeBoundBox(points()).extend(1e-4));
888 
889  labelList allEdges(identity(edges().size()));
890 
891  edgeTree_.reset
892  (
894  (
896  (
897  false, // cachebb
898  edges(), // edges
899  points(), // points
900  allEdges // selected edges
901  ),
902  bb, // bb
903  8, // maxLevel
904  10, // leafsize
905  3.0 // duplicity
906  )
907  );
908  }
909 
910  return edgeTree_();
911 }
912 
913 
916 {
917  if (edgeTreesByType_.size() == 0)
918  {
919  edgeTreesByType_.setSize(nEdgeTypes);
920 
921  // Slightly extended bb. Slightly off-centred just so on symmetric
922  // geometry there are less face/edge aligned items.
923  treeBoundBox bb(treeBoundBox(points()).extend(1e-4));
924 
925  labelListList sliceEdges(nEdgeTypes);
926 
927  // External edges
928  sliceEdges[0] =
930 
931  // Internal edges
932  sliceEdges[1] = identity(flatStart_ - internalStart_) + internalStart_;
933 
934  // Flat edges
935  sliceEdges[2] = identity(openStart_ - flatStart_) + flatStart_;
936 
937  // Open edges
938  sliceEdges[3] = identity(multipleStart_ - openStart_) + openStart_;
939 
940  // Multiple edges
941  sliceEdges[4] =
943 
945  {
946  edgeTreesByType_.set
947  (
948  i,
950  (
952  (
953  false, // cachebb
954  edges(), // edges
955  points(), // points
956  sliceEdges[i] // selected edges
957  ),
958  bb, // bb
959  8, // maxLevel
960  10, // leafsize
961  3.0 // duplicity
962  )
963  );
964  }
965  }
966 
967  return edgeTreesByType_;
968 }
969 
970 
972 {
973  edgeMesh::transfer(mesh);
974 
976  mixedStart_ = mesh.mixedStart_;
979  flatStart_ = mesh.flatStart_;
980  openStart_ = mesh.openStart_;
982  normals_.transfer(mesh.normals_);
990  pointTree_ = mesh.pointTree_;
991  edgeTree_ = mesh.edgeTree_;
992  edgeTreesByType_.transfer(mesh.edgeTreesByType_);
993 }
994 
995 
997 {
998  edgeMesh::clear();
999  concaveStart_ = 0;
1000  mixedStart_ = 0;
1001  nonFeatureStart_ = 0;
1002  internalStart_ = 0;
1003  flatStart_ = 0;
1004  openStart_ = 0;
1005  multipleStart_ = 0;
1006  normals_.clear();
1010  edgeNormals_.clear();
1013  regionEdges_.clear();
1014  pointTree_.clear();
1015  edgeTree_.clear();
1016  edgeTreesByType_.clear();
1017 }
1018 
1019 
1021 {
1022  // Points
1023  // ~~~~~~
1024 
1025  // From current points into combined points
1026  labelList reversePointMap(points().size());
1027  labelList reverseFemPointMap(fem.points().size());
1028 
1029  label newPointi = 0;
1030  for (label i = 0; i < concaveStart(); i++)
1031  {
1032  reversePointMap[i] = newPointi++;
1033  }
1034  for (label i = 0; i < fem.concaveStart(); i++)
1035  {
1036  reverseFemPointMap[i] = newPointi++;
1037  }
1038 
1039  // Concave
1040  label newConcaveStart = newPointi;
1041  for (label i = concaveStart(); i < mixedStart(); i++)
1042  {
1043  reversePointMap[i] = newPointi++;
1044  }
1045  for (label i = fem.concaveStart(); i < fem.mixedStart(); i++)
1046  {
1047  reverseFemPointMap[i] = newPointi++;
1048  }
1049 
1050  // Mixed
1051  label newMixedStart = newPointi;
1052  for (label i = mixedStart(); i < nonFeatureStart(); i++)
1053  {
1054  reversePointMap[i] = newPointi++;
1055  }
1056  for (label i = fem.mixedStart(); i < fem.nonFeatureStart(); i++)
1057  {
1058  reverseFemPointMap[i] = newPointi++;
1059  }
1060 
1061  // Non-feature
1062  label newNonFeatureStart = newPointi;
1063  for (label i = nonFeatureStart(); i < points().size(); i++)
1064  {
1065  reversePointMap[i] = newPointi++;
1066  }
1067  for (label i = fem.nonFeatureStart(); i < fem.points().size(); i++)
1068  {
1069  reverseFemPointMap[i] = newPointi++;
1070  }
1071 
1072  pointField newPoints(newPointi);
1073  newPoints.rmap(points(), reversePointMap);
1074  newPoints.rmap(fem.points(), reverseFemPointMap);
1075 
1076 
1077  // Edges
1078  // ~~~~~
1079 
1080  // From current edges into combined edges
1081  labelList reverseEdgeMap(edges().size());
1082  labelList reverseFemEdgeMap(fem.edges().size());
1083 
1084  // External
1085  label newEdgeI = 0;
1086  for (label i = 0; i < internalStart(); i++)
1087  {
1088  reverseEdgeMap[i] = newEdgeI++;
1089  }
1090  for (label i = 0; i < fem.internalStart(); i++)
1091  {
1092  reverseFemEdgeMap[i] = newEdgeI++;
1093  }
1094 
1095  // Internal
1096  label newInternalStart = newEdgeI;
1097  for (label i = internalStart(); i < flatStart(); i++)
1098  {
1099  reverseEdgeMap[i] = newEdgeI++;
1100  }
1101  for (label i = fem.internalStart(); i < fem.flatStart(); i++)
1102  {
1103  reverseFemEdgeMap[i] = newEdgeI++;
1104  }
1105 
1106  // Flat
1107  label newFlatStart = newEdgeI;
1108  for (label i = flatStart(); i < openStart(); i++)
1109  {
1110  reverseEdgeMap[i] = newEdgeI++;
1111  }
1112  for (label i = fem.flatStart(); i < fem.openStart(); i++)
1113  {
1114  reverseFemEdgeMap[i] = newEdgeI++;
1115  }
1116 
1117  // Open
1118  label newOpenStart = newEdgeI;
1119  for (label i = openStart(); i < multipleStart(); i++)
1120  {
1121  reverseEdgeMap[i] = newEdgeI++;
1122  }
1123  for (label i = fem.openStart(); i < fem.multipleStart(); i++)
1124  {
1125  reverseFemEdgeMap[i] = newEdgeI++;
1126  }
1127 
1128  // Multiple
1129  label newMultipleStart = newEdgeI;
1130  for (label i = multipleStart(); i < edges().size(); i++)
1131  {
1132  reverseEdgeMap[i] = newEdgeI++;
1133  }
1134  for (label i = fem.multipleStart(); i < fem.edges().size(); i++)
1135  {
1136  reverseFemEdgeMap[i] = newEdgeI++;
1137  }
1138 
1139  edgeList newEdges(newEdgeI);
1140  forAll(edges(), i)
1141  {
1142  const edge& e = edges()[i];
1143  newEdges[reverseEdgeMap[i]] = edge
1144  (
1145  reversePointMap[e[0]],
1146  reversePointMap[e[1]]
1147  );
1148  }
1149  forAll(fem.edges(), i)
1150  {
1151  const edge& e = fem.edges()[i];
1152  newEdges[reverseFemEdgeMap[i]] = edge
1153  (
1154  reverseFemPointMap[e[0]],
1155  reverseFemPointMap[e[1]]
1156  );
1157  }
1158 
1159  pointField newEdgeDirections(newEdgeI);
1160  newEdgeDirections.rmap(edgeDirections(), reverseEdgeMap);
1161  newEdgeDirections.rmap(fem.edgeDirections(), reverseFemEdgeMap);
1162 
1163 
1164 
1165 
1166  // Normals
1167  // ~~~~~~~
1168 
1169  // Combine normals
1170  DynamicField<point> newNormals(normals().size()+fem.normals().size());
1171  newNormals.append(normals());
1172  newNormals.append(fem.normals());
1173 
1174 
1175  // Combine and re-index into newNormals
1176  labelListList newEdgeNormals(edgeNormals().size()+fem.edgeNormals().size());
1177  UIndirectList<labelList>(newEdgeNormals, reverseEdgeMap) =
1178  edgeNormals();
1179  UIndirectList<labelList>(newEdgeNormals, reverseFemEdgeMap) =
1180  fem.edgeNormals();
1181  forAll(reverseFemEdgeMap, i)
1182  {
1183  label mapI = reverseFemEdgeMap[i];
1184  labelList& en = newEdgeNormals[mapI];
1185  forAll(en, j)
1186  {
1187  en[j] += normals().size();
1188  }
1189  }
1190 
1191 
1192  // Combine and re-index into newFeaturePointNormals
1193  labelListList newFeaturePointNormals
1194  (
1195  featurePointNormals().size()
1196  + fem.featurePointNormals().size()
1197  );
1198 
1199  // Note: featurePointNormals only go up to nonFeatureStart
1200  UIndirectList<labelList>
1201  (
1202  newFeaturePointNormals,
1203  SubList<label>(reversePointMap, featurePointNormals().size())
1204  ) = featurePointNormals();
1205  UIndirectList<labelList>
1206  (
1207  newFeaturePointNormals,
1208  SubList<label>(reverseFemPointMap, fem.featurePointNormals().size())
1209  ) = fem.featurePointNormals();
1210  forAll(fem.featurePointNormals(), i)
1211  {
1212  label mapI = reverseFemPointMap[i];
1213  labelList& fn = newFeaturePointNormals[mapI];
1214  forAll(fn, j)
1215  {
1216  fn[j] += normals().size();
1217  }
1218  }
1219 
1220 
1221  // Combine regionEdges
1222  DynamicList<label> newRegionEdges
1223  (
1224  regionEdges().size()
1225  + fem.regionEdges().size()
1226  );
1227  forAll(regionEdges(), i)
1228  {
1229  newRegionEdges.append(reverseEdgeMap[regionEdges()[i]]);
1230  }
1231  forAll(fem.regionEdges(), i)
1232  {
1233  newRegionEdges.append(reverseFemEdgeMap[fem.regionEdges()[i]]);
1234  }
1235 
1236 
1237  // Assign
1238  // ~~~~~~
1239 
1240  // Transfer
1241  concaveStart_ = newConcaveStart;
1242  mixedStart_ = newMixedStart;
1243  nonFeatureStart_ = newNonFeatureStart;
1244 
1245  // Reset points and edges
1246  reset(move(newPoints), move(newEdges));
1247 
1248  // Transfer
1249  internalStart_ = newInternalStart;
1250  flatStart_ = newFlatStart;
1251  openStart_ = newOpenStart;
1252  multipleStart_ = newMultipleStart;
1253 
1254  edgeDirections_.transfer(newEdgeDirections);
1255 
1256  normals_.transfer(newNormals);
1257  edgeNormals_.transfer(newEdgeNormals);
1258  featurePointNormals_.transfer(newFeaturePointNormals);
1259 
1260  regionEdges_.transfer(newRegionEdges);
1261 
1262  pointTree_.clear();
1263  edgeTree_.clear();
1264  edgeTreesByType_.clear();
1265 }
1266 
1267 
1269 {
1270  // Points
1271  // ~~~~~~
1272 
1273  // From current points into new points
1274  labelList reversePointMap(identity(points().size()));
1275 
1276  // Flip convex and concave points
1277 
1278  label newPointi = 0;
1279  // Concave points become convex
1280  for (label i = concaveStart(); i < mixedStart(); i++)
1281  {
1282  reversePointMap[i] = newPointi++;
1283  }
1284  // Convex points become concave
1285  label newConcaveStart = newPointi;
1286  for (label i = 0; i < concaveStart(); i++)
1287  {
1288  reversePointMap[i] = newPointi++;
1289  }
1290 
1291 
1292  // Edges
1293  // ~~~~~~
1294 
1295  // From current edges into new edges
1296  labelList reverseEdgeMap(identity(edges().size()));
1297 
1298  // Flip external and internal edges
1299 
1300  label newEdgeI = 0;
1301  // Internal become external
1302  for (label i = internalStart(); i < flatStart(); i++)
1303  {
1304  reverseEdgeMap[i] = newEdgeI++;
1305  }
1306  // External become internal
1307  label newInternalStart = newEdgeI;
1308  for (label i = 0; i < internalStart(); i++)
1309  {
1310  reverseEdgeMap[i] = newEdgeI++;
1311  }
1312 
1313 
1314  pointField newPoints(points().size());
1315  newPoints.rmap(points(), reversePointMap);
1316 
1317  edgeList newEdges(edges().size());
1318  forAll(edges(), i)
1319  {
1320  const edge& e = edges()[i];
1321  newEdges[reverseEdgeMap[i]] = edge
1322  (
1323  reversePointMap[e[0]],
1324  reversePointMap[e[1]]
1325  );
1326  }
1327 
1328 
1329  // Normals are flipped
1330  // ~~~~~~~~~~~~~~~~~~~
1331 
1332  pointField newEdgeDirections(edges().size());
1333  newEdgeDirections.rmap(-1.0*edgeDirections(), reverseEdgeMap);
1334 
1335  pointField newNormals(-1.0*normals());
1336 
1337  labelListList newEdgeNormals(edgeNormals().size());
1338  UIndirectList<labelList>(newEdgeNormals, reverseEdgeMap) = edgeNormals();
1339 
1340  labelListList newFeaturePointNormals(featurePointNormals().size());
1341 
1342  // Note: featurePointNormals only go up to nonFeatureStart
1344  (
1345  newFeaturePointNormals,
1346  SubList<label>(reversePointMap, featurePointNormals().size())
1347  ) = featurePointNormals();
1348 
1349  labelList newRegionEdges(regionEdges().size());
1350  forAll(regionEdges(), i)
1351  {
1352  newRegionEdges[i] = reverseEdgeMap[regionEdges()[i]];
1353  }
1354 
1355  // Transfer
1356  concaveStart_ = newConcaveStart;
1357 
1358  // Reset points and edges
1359  reset(move(newPoints), move(newEdges));
1360 
1361  // Transfer
1362  internalStart_ = newInternalStart;
1363 
1364  edgeDirections_.transfer(newEdgeDirections);
1365  normals_.transfer(newNormals);
1366  edgeNormals_.transfer(newEdgeNormals);
1367  featurePointNormals_.transfer(newFeaturePointNormals);
1368  regionEdges_.transfer(newRegionEdges);
1369 
1370  pointTree_.clear();
1371  edgeTree_.clear();
1372  edgeTreesByType_.clear();
1373 }
1374 
1375 
1378  const fileName& prefix,
1379  const bool verbose
1380 ) const
1381 {
1382  Info<< nl << "Writing extendedEdgeMesh components to " << prefix
1383  << endl;
1384 
1385  edgeMesh::write(prefix + "_edgeMesh.obj");
1386 
1387  if (!verbose) return;
1388 
1389  OBJstream convexFtPtStr(prefix + "_convexFeaturePts.obj");
1390  Info<< "Writing convex feature points to " << convexFtPtStr.name() << endl;
1391 
1392  for(label i = 0; i < concaveStart_; i++)
1393  {
1394  convexFtPtStr.write(points()[i]);
1395  }
1396 
1397  OBJstream concaveFtPtStr(prefix + "_concaveFeaturePts.obj");
1398  Info<< "Writing concave feature points to "
1399  << concaveFtPtStr.name() << endl;
1400 
1401  for(label i = concaveStart_; i < mixedStart_; i++)
1402  {
1403  convexFtPtStr.write(points()[i]);
1404  }
1405 
1406  OBJstream mixedFtPtStr(prefix + "_mixedFeaturePts.obj");
1407  Info<< "Writing mixed feature points to " << mixedFtPtStr.name() << endl;
1408 
1409  for(label i = mixedStart_; i < nonFeatureStart_; i++)
1410  {
1411  mixedFtPtStr.write(points()[i]);
1412  }
1413 
1414  OBJstream mixedFtPtStructureStr(prefix + "_mixedFeaturePtsStructure.obj");
1415  Info<< "Writing mixed feature point structure to "
1416  << mixedFtPtStructureStr.name() << endl;
1417 
1418  for(label i = mixedStart_; i < nonFeatureStart_; i++)
1419  {
1420  const labelList& ptEds = pointEdges()[i];
1421 
1422  forAll(ptEds, j)
1423  {
1424  const edge& e = edges()[ptEds[j]];
1425  mixedFtPtStructureStr.write
1426  (
1427  linePointRef(points()[e[0]],
1428  points()[e[1]])
1429  );
1430  }
1431  }
1432 
1433  OBJstream externalStr(prefix + "_externalEdges.obj");
1434  Info<< "Writing external edges to " << externalStr.name() << endl;
1435 
1436  for (label i = externalStart_; i < internalStart_; i++)
1437  {
1438  const edge& e = edges()[i];
1439  externalStr.write(linePointRef(points()[e[0]], points()[e[1]]));
1440  }
1441 
1442  OBJstream internalStr(prefix + "_internalEdges.obj");
1443  Info<< "Writing internal edges to " << internalStr.name() << endl;
1444 
1445  for (label i = internalStart_; i < flatStart_; i++)
1446  {
1447  const edge& e = edges()[i];
1448  internalStr.write(linePointRef(points()[e[0]], points()[e[1]]));
1449  }
1450 
1451  OBJstream flatStr(prefix + "_flatEdges.obj");
1452  Info<< "Writing flat edges to " << flatStr.name() << endl;
1453 
1454  for (label i = flatStart_; i < openStart_; i++)
1455  {
1456  const edge& e = edges()[i];
1457  flatStr.write(linePointRef(points()[e[0]], points()[e[1]]));
1458  }
1459 
1460  OBJstream openStr(prefix + "_openEdges.obj");
1461  Info<< "Writing open edges to " << openStr.name() << endl;
1462 
1463  for (label i = openStart_; i < multipleStart_; i++)
1464  {
1465  const edge& e = edges()[i];
1466  openStr.write(linePointRef(points()[e[0]], points()[e[1]]));
1467  }
1468 
1469  OBJstream multipleStr(prefix + "_multipleEdges.obj");
1470  Info<< "Writing multiple edges to " << multipleStr.name() << endl;
1471 
1472  for (label i = multipleStart_; i < edges().size(); i++)
1473  {
1474  const edge& e = edges()[i];
1475  multipleStr.write(linePointRef(points()[e[0]], points()[e[1]]));
1476  }
1477 
1478  OBJstream regionStr(prefix + "_regionEdges.obj");
1479  Info<< "Writing region edges to " << regionStr.name() << endl;
1480 
1481  forAll(regionEdges_, i)
1482  {
1483  const edge& e = edges()[regionEdges_[i]];
1484  regionStr.write(linePointRef(points()[e[0]], points()[e[1]]));
1485  }
1486 
1487  OBJstream edgeDirsStr(prefix + "_edgeDirections.obj");
1488  Info<< "Writing edge directions to " << edgeDirsStr.name() << endl;
1489 
1491  {
1492  const vector& eVec = edgeDirections_[i];
1493  const edge& e = edges()[i];
1494 
1495  edgeDirsStr.write
1496  (
1497  linePointRef(points()[e.start()], eVec + points()[e.start()])
1498  );
1499  }
1500 }
1501 
1502 
1504 {
1506 
1507  os << indent << "point classification :" << nl;
1508  os << incrIndent;
1509  os << indent << "convex feature points : "
1511  << nl;
1512  os << indent << "concave feature points : "
1513  << setw(8) << mixedStart_-concaveStart_
1514  << nl;
1515  os << indent << "mixed feature points : "
1517  << nl;
1518  os << indent << "other (non-feature) points : "
1519  << setw(8) << points().size()-nonFeatureStart_
1520  << nl;
1521  os << decrIndent;
1522 
1523  os << indent << "edge classification :" << nl;
1524  os << incrIndent;
1525  os << indent << "external (convex angle) edges : "
1527  << nl;
1528  os << indent << "internal (concave angle) edges : "
1529  << setw(8) << flatStart_-internalStart_
1530  << nl;
1531  os << indent << "flat region edges : "
1532  << setw(8) << openStart_-flatStart_
1533  << nl;
1534  os << indent << "open edges : "
1535  << setw(8) << multipleStart_-openStart_
1536  << nl;
1537  os << indent << "multiply connected edges : "
1538  << setw(8) << edges().size()-multipleStart_
1539  << nl;
1540  os << decrIndent;
1541 }
1542 
1543 
1544 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
1545 
1546 Foam::Istream& Foam::operator>>
1548  Istream& is,
1550 )
1551 {
1552  label type;
1553  is >> type;
1554 
1555  vt = static_cast<Foam::extendedEdgeMesh::sideVolumeType>(type);
1556 
1557  // Check state of Istream
1558  is.check("operator>>(Istream&, sideVolumeType&)");
1559 
1560  return is;
1561 }
1562 
1563 
1564 Foam::Ostream& Foam::operator<<
1566  Ostream& os,
1568 )
1569 {
1570  os << static_cast<label>(vt);
1571 
1572  return os;
1573 }
1574 
1575 
1577 {
1578  // fileFormats::extendedEdgeMeshFormat::write(os, em.points_, em.edges_);
1579  os << "// points" << nl
1580  << em.points() << nl
1581  << "// edges" << nl
1582  << em.edges() << nl
1583  << "// concaveStart mixedStart nonFeatureStart" << nl
1584  << em.concaveStart_ << token::SPACE
1585  << em.mixedStart_ << token::SPACE
1586  << em.nonFeatureStart_ << nl
1587  << "// internalStart flatStart openStart multipleStart" << nl
1588  << em.internalStart_ << token::SPACE
1589  << em.flatStart_ << token::SPACE
1590  << em.openStart_ << token::SPACE
1591  << em.multipleStart_ << nl
1592  << "// normals" << nl
1593  << em.normals_ << nl
1594  << "// normal volume types" << nl
1595  << em.normalVolumeTypes_ << nl
1596  << "// normalDirections" << nl
1597  << em.normalDirections_ << nl
1598  << "// edgeNormals" << nl
1599  << em.edgeNormals_ << nl
1600  << "// featurePointNormals" << nl
1601  << em.featurePointNormals_ << nl
1602  << "// featurePointEdges" << nl
1603  << em.featurePointEdges_ << nl
1604  << "// regionEdges" << nl
1605  << em.regionEdges_
1606  << endl;
1607 
1608  // Check state of Ostream
1609  os.check("Ostream& operator<<(Ostream&, const extendedEdgeMesh&)");
1610 
1611  return os;
1612 }
1613 
1614 
1616 {
1617  // fileFormats::extendedEdgeMeshFormat::read(is, em.points_, em.edges_);
1618  is >> static_cast<edgeMesh&>(em)
1619  >> em.concaveStart_
1620  >> em.mixedStart_
1621  >> em.nonFeatureStart_
1622  >> em.internalStart_
1623  >> em.flatStart_
1624  >> em.openStart_
1625  >> em.multipleStart_
1626  >> em.normals_
1627  >> em.normalVolumeTypes_
1628  >> em.normalDirections_
1629  >> em.edgeNormals_
1630  >> em.featurePointNormals_
1631  >> em.featurePointEdges_
1632  >> em.regionEdges_;
1633 
1634  // Check state of Istream
1635  is.check("Istream& operator>>(Istream&, extendedEdgeMesh&)");
1636 
1637  return is;
1638 }
1639 
1640 
1641 // ************************************************************************* //
extendedEdgeMesh()
Construct null.
label nonFeatureStart() const
Return the index of the start of the non-feature points.
static autoPtr< extendedEdgeMesh > New(const fileName &, const word &ext)
Select constructed from filename (explicit extension)
void transfer(edgeMesh &)
Transfer the contents of the argument and annul the argument.
Definition: edgeMesh.C:197
const labelListList & featurePointNormals() const
Return the indices of the normals that are adjacent to the.
A HashTable with keys but without contents.
Definition: HashSet.H:59
static const scalar GREAT
Definition: scalar.H:111
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A class for handling file names.
Definition: fileName.H:79
Holds (reference to) pointField. Encapsulation of data needed for octree searches. Used for searching for nearest point. No bounding boxes around points. Only overlaps and calcNearest are implemented, rest makes little sense.
Definition: treeDataPoint.H:59
static wordHashSet readTypes()
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
vectorField normals_
Normals of the features, to be referred to by index by both feature.
label nRegionEdges() const
Return number of region edges.
const labelListList & edgeNormals() const
Return the indices of the normals that are adjacent to the.
void nearestFeaturePoint(const point &sample, scalar searchDistSqr, pointIndexHit &info) const
Find nearest surface edge for the sample point.
static bool canWriteType(const word &ext, const bool verbose=false)
Can we write this file format type?
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
void sortPointsAndEdges(const Patch &, const labelList &featureEdges, const labelList &regionFeatureEdges, const labelList &feaurePoints)
PtrList< indexedOctree< treeDataEdge > > edgeTreesByType_
Individual search trees for each type of edge.
label mixedStart_
Index of the start of the mixed type feature points.
virtual void writeStats(Ostream &os) const
Dump some information.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Description of feature edges and points.
label multipleStart() const
Return the index of the start of the multiply-connected feature.
label openStart_
Index of the start of the open feature edges.
static const Foam::NamedEnum< sideVolumeType, 4 > sideVolumeTypeNames_
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Holds data for octree to work on an edges subset.
Definition: treeDataEdge.H:53
virtual void clear()
Clear all storage.
static const Foam::NamedEnum< edgeStatus, 6 > edgeStatusNames_
const vectorField & normals() const
Return the normals of the surfaces adjacent to the feature edges.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
label concaveStart_
Index of the start of the concave feature points.
label openStart() const
Return the index of the start of the open feature edges.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
List< sideVolumeType > normalVolumeTypes_
Type per normal: which side of normal to mesh.
void add(const extendedEdgeMesh &)
Add extendedEdgeMesh. No filtering of duplicates.
const labelList & featureEdges() const
Return feature edge list.
const labelList & regionEdges() const
Return the feature edges which are on the boundary between.
virtual Ostream & write(const char)
Write character.
Definition: OBJstream.C:81
fvMesh & mesh
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:120
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
scalar minDist(const List< pointIndexHit > &hitList)
word ext() const
Return file name extension (part after last .)
Definition: fileName.C:299
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
edgeMesh()
Construct null.
Definition: edgeMesh.C:123
label internalStart_
Index of the start of the internal feature edges.
void nearestFeatureEdgeByType(const point &sample, const scalarField &searchDistSqr, pointIndexHitList &info) const
Find the nearest point on each type of feature edge.
labelListList normalDirections_
Starting directions for the edges.
A list of faces which address into the list of points.
const Point & hitPoint() const
Return hit point.
A List obtained as a section of another List.
Definition: SubList.H:53
linePointRef line(const pointField &) const
Return edge line.
Definition: edgeI.H:187
const labelListList & pointEdges() const
Return edges.
Definition: edgeMeshI.H:74
void allNearestFeatureEdges(const point &sample, const scalar searchRadiusSqr, pointIndexHitList &info) const
Find all the feature edges within searchDistSqr of sample.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
void flipNormals()
Flip normals. All concave become convex, all internal external.
label flatStart() const
Return the index of the start of the flat feature edges.
dimensionedScalar cos(const dimensionedScalar &ds)
line< point, const point & > linePointRef
Line using referred points.
Definition: linePointRef.H:45
labelList regionEdges_
Feature edges which are on the boundary between regions.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
List< edge > edgeList
Definition: edgeList.H:38
static label convexStart_
Index of the start of the convex feature points - static as 0.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
label internalStart() const
Return the index of the start of the internal feature edges.
static label nEdgeTypes
Number of possible feature edge types (i.e. number of slices)
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
A class for handling words, derived from string.
Definition: word.H:59
label flatStart_
Index of the start of the flat feature edges.
Istream & operator>>(Istream &, directionInfo &)
const edgeList & edges() const
Return edges.
Definition: edgeMeshI.H:68
const triSurface & surface() const
bool read(const fileName &, const word &ext)
Read from file. Chooses reader based on explicit extension.
~extendedEdgeMesh()
Destructor.
bool hit() const
Is there a hit.
sideVolumeType
Normals point to the outside.
static label externalStart_
Index of the start of the external feature edges - static as 0.
void transfer(extendedEdgeMesh &)
Transfer the contents of the argument and annul the argument.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
const labelListList & edgeFaces() const
Return edge-face addressing.
bool hit() const
Is there a hit.
Definition: PointHit.H:120
Dynamically sized Field.
Definition: DynamicField.H:49
HashSet wordHashSet
A HashSet with word keys.
Definition: HashSet.H:208
autoPtr< indexedOctree< treeDataEdge > > edgeTree_
Search tree for all edges.
const indexedOctree< treeDataEdge > & edgeTree() const
Demand driven construction of octree for boundary edges.
static scalar cosNormalAngleTol_
Angular closeness tolerance for treating normals as the same.
label mixedStart() const
Return the index of the start of the mixed type feature points.
bool connected(const edge &a) const
Return true if connected to given edge.
Definition: edgeI.H:103
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
const Point & rawPoint() const
Return point with no checking.
Definition: PointHit.H:158
Istream and Ostream manipulators taking arguments.
OFstream which keeps track of vertices.
Definition: OBJstream.H:53
static const char nl
Definition: Ostream.H:260
const labelList & featurePoints() const
Return feature point list.
const pointField & points() const
Return points.
Definition: edgeMeshI.H:62
defineTypeNameAndDebug(combustionModel, 0)
Points connected by edges.
Definition: edgeMesh.H:69
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:235
virtual void writeStats(Ostream &) const
Definition: edgeMeshIO.C:117
const PtrList< indexedOctree< treeDataEdge > > & edgeTreesByType() const
Demand driven construction of octree for boundary edges by type.
virtual void reset(pointField &&points, edgeList &&edges)
Reset primitive data (points, edges)
Definition: edgeMesh.C:175
static edgeStatus classifyEdge(const List< vector > &norms, const labelList &edNorms, const vector &fC0tofC1)
Classify the type of feature edge. Requires face centre 0 to face.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
static bool canReadType(const word &ext, const bool verbose=false)
Can we read this file format?
virtual void clear()
Clear all storage.
Definition: edgeMesh.C:166
void setSize(const label)
Reset size of List.
Definition: List.C:281
const vectorField & edgeDirections() const
Return the edgeDirection vectors.
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:362
label multipleStart_
Index of the start of the multiply-connected feature edges.
static const Foam::NamedEnum< pointStatus, 4 > pointStatusNames_
label newPointi
Definition: readKivaGrid.H:501
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
A List with indirect addressing.
Definition: fvMatrix.H:106
Ostream & operator<<(Ostream &, const ensightPart &)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
static bool canRead(const fileName &, const bool verbose=false)
Can we read this file format?
static void write(const fileName &, const edgeMesh &)
Write to file.
Definition: edgeMeshIO.C:87
fileName lessExt() const
Return file name without extension (part before last .)
Definition: fileName.C:284
labelListList featurePointEdges_
Indices of feature edges attached to feature points. The edges are.
void writeObj(const fileName &prefix, const bool verbose=true) const
Write all components of the extendedEdgeMesh as obj files.
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
messageStream Info
const indexedOctree< treeDataPoint > & pointTree() const
Demand driven construction of octree for feature points.
dimensioned< scalar > mag(const dimensioned< Type > &)
static wordHashSet writeTypes()
scalar minDisconnectedDist(const pointIndexHitList &hitList) const
Return the minimum distance between disconnected edges.
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
labelListList featurePointNormals_
Indices of the normals that are adjacent to the feature points.
labelListList edgeNormals_
Indices of the normals that are adjacent to the feature edges.
Triangulated surface description with patch information.
Definition: triSurface.H:66
label nonFeatureStart_
Index of the start of the non-feature points.
void nearestFeatureEdge(const point &sample, scalar searchDistSqr, pointIndexHit &info) const
Find nearest surface edge for the sample point.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:228
vectorField edgeDirections_
Flat and open edges require the direction of the edge.
autoPtr< indexedOctree< treeDataPoint > > pointTree_
Search tree for all feature points.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
label index() const
Return index.
PointHit< Point > nearestDist(const Point &p) const
Return nearest distance to line from a given point.
Definition: lineI.H:95
DynamicField< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
label start() const
Return start vertex label.
Definition: edgeI.H:81
Holds feature edges/points of surface.
void allNearestFeaturePoints(const point &sample, scalar searchRadiusSqr, pointIndexHitList &info) const
Find all the feature points within searchDistSqr of sample.
Namespace for OpenFOAM.
label concaveStart() const
Return the index of the start of the concave feature points.
pointStatus classifyFeaturePoint(label ptI) const
Classify the type of feature point. Requires valid stored member.