refinementFeatures.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-2014 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 "refinementFeatures.H"
27 #include "Time.H"
28 #include "Tuple2.H"
29 #include "DynamicField.H"
30 
31 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
32 
33 void Foam::refinementFeatures::read
34 (
35  const objectRegistry& io,
36  const PtrList<dictionary>& featDicts
37 )
38 {
39  forAll(featDicts, featI)
40  {
41  const dictionary& dict = featDicts[featI];
42 
43  fileName featFileName(dict.lookup("file"));
44 
45 
46  // Try reading extendedEdgeMesh first
47 
48  IOobject extFeatObj
49  (
50  featFileName, // name
51  io.time().constant(), // instance
52  "extendedFeatureEdgeMesh", // local
53  io.time(), // registry
56  false
57  );
58 
59  const fileName fName(extFeatObj.filePath());
60 
61  if (!fName.empty() && extendedEdgeMesh::canRead(fName))
62  {
63  autoPtr<extendedEdgeMesh> eMeshPtr = extendedEdgeMesh::New
64  (
65  fName
66  );
67 
68  Info<< "Read extendedFeatureEdgeMesh " << extFeatObj.name()
69  << nl << incrIndent;
70  eMeshPtr().writeStats(Info);
71  Info<< decrIndent << endl;
72 
73  set(featI, new extendedFeatureEdgeMesh(extFeatObj, eMeshPtr()));
74  }
75  else
76  {
77  // Try reading edgeMesh
78 
79  IOobject featObj
80  (
81  featFileName, // name
82  io.time().constant(), // instance
83  "triSurface", // local
84  io.time(), // registry
87  false
88  );
89 
90  const fileName fName(featObj.filePath());
91 
92  if (fName.empty())
93  {
95  (
96  "refinementFeatures::read"
97  "(const objectRegistry&"
98  ", const PtrList<dictionary>&)",
99  dict
100  ) << "Could not open " << featObj.objectPath()
101  << exit(FatalIOError);
102  }
103 
104 
105  // Read as edgeMesh
106  autoPtr<edgeMesh> eMeshPtr = edgeMesh::New(fName);
107  const edgeMesh& eMesh = eMeshPtr();
108 
109  Info<< "Read edgeMesh " << featObj.name() << nl
110  << incrIndent;
111  eMesh.writeStats(Info);
112  Info<< decrIndent << endl;
113 
114 
115  // Analyse for feature points. These are all classified as mixed
116  // points for lack of anything better
117  const labelListList& pointEdges = eMesh.pointEdges();
118 
119  labelList oldToNew(eMesh.points().size(), -1);
120  DynamicField<point> newPoints(eMesh.points().size());
121  forAll(pointEdges, pointI)
122  {
123  if (pointEdges[pointI].size() > 2)
124  {
125  oldToNew[pointI] = newPoints.size();
126  newPoints.append(eMesh.points()[pointI]);
127  }
128  //else if (pointEdges[pointI].size() == 2)
129  //MEJ: do something based on a feature angle?
130  }
131  label nFeatures = newPoints.size();
132  forAll(oldToNew, pointI)
133  {
134  if (oldToNew[pointI] == -1)
135  {
136  oldToNew[pointI] = newPoints.size();
137  newPoints.append(eMesh.points()[pointI]);
138  }
139  }
140 
141 
142  const edgeList& edges = eMesh.edges();
143  edgeList newEdges(edges.size());
144  forAll(edges, edgeI)
145  {
146  const edge& e = edges[edgeI];
147  newEdges[edgeI] = edge
148  (
149  oldToNew[e[0]],
150  oldToNew[e[1]]
151  );
152  }
153 
154  // Construct an extendedEdgeMesh with
155  // - all points on more than 2 edges : mixed feature points
156  // - all edges : external edges
157 
158  extendedEdgeMesh eeMesh
159  (
160  newPoints, // pts
161  newEdges, // eds
162  0, // (point) concaveStart
163  0, // (point) mixedStart
164  nFeatures, // (point) nonFeatureStart
165  edges.size(), // (edge) internalStart
166  edges.size(), // (edge) flatStart
167  edges.size(), // (edge) openStart
168  edges.size(), // (edge) multipleStart
169  vectorField(0), // normals
170  List<extendedEdgeMesh::sideVolumeType>(0),// normalVolumeTypes
171  vectorField(0), // edgeDirections
172  labelListList(0), // normalDirections
173  labelListList(0), // edgeNormals
174  labelListList(0), // featurePointNormals
175  labelListList(0), // featurePointEdges
176  identity(newEdges.size()) // regionEdges
177  );
178 
179  //Info<< "Constructed extendedFeatureEdgeMesh " << featObj.name()
180  // << nl << incrIndent;
181  //eeMesh.writeStats(Info);
182  //Info<< decrIndent << endl;
183 
184  set(featI, new extendedFeatureEdgeMesh(featObj, eeMesh));
185  }
186 
187  const extendedEdgeMesh& eMesh = operator[](featI);
188 
189  if (dict.found("levels"))
190  {
191  List<Tuple2<scalar, label> > distLevels(dict["levels"]);
192 
193  if (dict.size() < 1)
194  {
196  (
197  "refinementFeatures::read"
198  "(const objectRegistry&"
199  ", const PtrList<dictionary>&)"
200  ) << " : levels should be at least size 1" << endl
201  << "levels : " << dict["levels"]
202  << exit(FatalError);
203  }
204 
205  distances_[featI].setSize(distLevels.size());
206  levels_[featI].setSize(distLevels.size());
207 
208  forAll(distLevels, j)
209  {
210  distances_[featI][j] = distLevels[j].first();
211  levels_[featI][j] = distLevels[j].second();
212 
213  // Check in incremental order
214  if (j > 0)
215  {
216  if
217  (
218  (distances_[featI][j] <= distances_[featI][j-1])
219  || (levels_[featI][j] > levels_[featI][j-1])
220  )
221  {
223  (
224  "refinementFeatures::read"
225  "(const objectRegistry&"
226  ", const PtrList<dictionary>&)"
227  ) << " : Refinement should be specified in order"
228  << " of increasing distance"
229  << " (and decreasing refinement level)." << endl
230  << "Distance:" << distances_[featI][j]
231  << " refinementLevel:" << levels_[featI][j]
232  << exit(FatalError);
233  }
234  }
235  }
236  }
237  else
238  {
239  // Look up 'level' for single level
240  levels_[featI] = labelList(1, readLabel(dict.lookup("level")));
241  distances_[featI] = scalarField(1, 0.0);
242  }
243 
244  Info<< "Refinement level according to distance to "
245  << featFileName << " (" << eMesh.points().size() << " points, "
246  << eMesh.edges().size() << " edges)." << endl;
247  forAll(levels_[featI], j)
248  {
249  Info<< " level " << levels_[featI][j]
250  << " for all cells within " << distances_[featI][j]
251  << " metre." << endl;
252  }
253  }
254 }
255 
256 
257 void Foam::refinementFeatures::buildTrees(const label featI)
258 {
259  const extendedEdgeMesh& eMesh = operator[](featI);
260  const pointField& points = eMesh.points();
261  const edgeList& edges = eMesh.edges();
262 
263  // Calculate bb of all points
264  treeBoundBox bb(points);
265 
266  // Random number generator. Bit dodgy since not exactly random ;-)
267  Random rndGen(65431);
268 
269  // Slightly extended bb. Slightly off-centred just so on symmetric
270  // geometry there are less face/edge aligned items.
271  bb = bb.extend(rndGen, 1e-4);
272  bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
273  bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
274 
275  edgeTrees_.set
276  (
277  featI,
278  new indexedOctree<treeDataEdge>
279  (
280  treeDataEdge
281  (
282  false, // do not cache bb
283  edges,
284  points,
285  identity(edges.size())
286  ),
287  bb, // overall search domain
288  8, // maxLevel
289  10, // leafsize
290  3.0 // duplicity
291  )
292  );
293 
294 
295  labelList featurePoints(identity(eMesh.nonFeatureStart()));
296 
297  pointTrees_.set
298  (
299  featI,
300  new indexedOctree<treeDataPoint>
301  (
302  treeDataPoint(points, featurePoints),
303  bb, // overall search domain
304  8, // maxLevel
305  10, // leafsize
306  3.0 // duplicity
307  )
308  );
309 }
310 
311 
314 {
315  if (!regionEdgeTreesPtr_.valid())
316  {
317  regionEdgeTreesPtr_.reset
318  (
320  );
321  PtrList<indexedOctree<treeDataEdge> >& trees = regionEdgeTreesPtr_();
322 
323  forAll(*this, featI)
324  {
325  const extendedEdgeMesh& eMesh = operator[](featI);
326  const pointField& points = eMesh.points();
327  const edgeList& edges = eMesh.edges();
328 
329  // Calculate bb of all points
330  treeBoundBox bb(points);
331 
332  // Random number generator. Bit dodgy since not exactly random ;-)
333  Random rndGen(65431);
334 
335  // Slightly extended bb. Slightly off-centred just so on symmetric
336  // geometry there are less face/edge aligned items.
337  bb = bb.extend(rndGen, 1e-4);
338  bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
339  bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
340 
341  trees.set
342  (
343  featI,
345  (
347  (
348  false, // do not cache bb
349  edges,
350  points,
351  eMesh.regionEdges()
352  ),
353  bb, // overall search domain
354  8, // maxLevel
355  10, // leafsize
356  3.0 // duplicity
357  )
358  );
359  }
360  }
361  return regionEdgeTreesPtr_();
362 }
363 
364 
365 // Find maximum level of a shell.
366 void Foam::refinementFeatures::findHigherLevel
367 (
368  const pointField& pt,
369  const label featI,
370  labelList& maxLevel
371 ) const
372 {
373  const labelList& levels = levels_[featI];
374 
375  const scalarField& distances = distances_[featI];
376 
377  // Collect all those points that have a current maxLevel less than
378  // (any of) the shell. Also collect the furthest distance allowable
379  // to any shell with a higher level.
380 
381  pointField candidates(pt.size());
382  labelList candidateMap(pt.size());
383  scalarField candidateDistSqr(pt.size());
384  label candidateI = 0;
385 
386  forAll(maxLevel, pointI)
387  {
388  forAllReverse(levels, levelI)
389  {
390  if (levels[levelI] > maxLevel[pointI])
391  {
392  candidates[candidateI] = pt[pointI];
393  candidateMap[candidateI] = pointI;
394  candidateDistSqr[candidateI] = sqr(distances[levelI]);
395  candidateI++;
396  break;
397  }
398  }
399  }
400  candidates.setSize(candidateI);
401  candidateMap.setSize(candidateI);
402  candidateDistSqr.setSize(candidateI);
403 
404  // Do the expensive nearest test only for the candidate points.
405  const indexedOctree<treeDataEdge>& tree = edgeTrees_[featI];
406 
407  List<pointIndexHit> nearInfo(candidates.size());
408  forAll(candidates, candidateI)
409  {
410  nearInfo[candidateI] = tree.findNearest
411  (
412  candidates[candidateI],
413  candidateDistSqr[candidateI]
414  );
415  }
416 
417  // Update maxLevel
418  forAll(nearInfo, candidateI)
419  {
420  if (nearInfo[candidateI].hit())
421  {
422  // Check which level it actually is in.
423  label minDistI = findLower
424  (
425  distances,
426  mag(nearInfo[candidateI].hitPoint()-candidates[candidateI])
427  );
428 
429  label pointI = candidateMap[candidateI];
430 
431  // pt is inbetween shell[minDistI] and shell[minDistI+1]
432  maxLevel[pointI] = levels[minDistI+1];
433  }
434  }
435 }
436 
437 
438 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
439 
441 (
442  const objectRegistry& io,
443  const PtrList<dictionary>& featDicts
444 )
445 :
447  distances_(featDicts.size()),
448  levels_(featDicts.size()),
449  edgeTrees_(featDicts.size()),
450  pointTrees_(featDicts.size())
451 {
452  // Read features
453  read(io, featDicts);
454 
455  // Search engines
456  forAll(*this, i)
457  {
458  buildTrees(i);
459  }
460 }
461 
462 
463 //Foam::refinementFeatures::refinementFeatures
464 //(
465 // const objectRegistry& io,
466 // const PtrList<dictionary>& featDicts,
467 // const scalar minCos
468 //)
469 //:
470 // PtrList<extendedFeatureEdgeMesh>(featDicts.size()),
471 // distances_(featDicts.size()),
472 // levels_(featDicts.size()),
473 // edgeTrees_(featDicts.size()),
474 // pointTrees_(featDicts.size())
475 //{
476 // // Read features
477 // read(io, featDicts);
478 //
479 // // Search engines
480 // forAll(*this, i)
481 // {
482 // const edgeMesh& eMesh = operator[](i);
483 // const pointField& points = eMesh.points();
484 // const edgeList& edges = eMesh.edges();
485 // const labelListList& pointEdges = eMesh.pointEdges();
486 //
487 // DynamicList<label> featurePoints;
488 // forAll(pointEdges, pointI)
489 // {
490 // const labelList& pEdges = pointEdges[pointI];
491 // if (pEdges.size() > 2)
492 // {
493 // featurePoints.append(pointI);
494 // }
495 // else if (pEdges.size() == 2)
496 // {
497 // // Check the angle
498 // const edge& e0 = edges[pEdges[0]];
499 // const edge& e1 = edges[pEdges[1]];
500 //
501 // const point& p = points[pointI];
502 // const point& p0 = points[e0.otherVertex(pointI)];
503 // const point& p1 = points[e1.otherVertex(pointI)];
504 //
505 // vector v0 = p-p0;
506 // scalar v0Mag = mag(v0);
507 //
508 // vector v1 = p1-p;
509 // scalar v1Mag = mag(v1);
510 //
511 // if
512 // (
513 // v0Mag > SMALL
514 // && v1Mag > SMALL
515 // && ((v0/v0Mag & v1/v1Mag) < minCos)
516 // )
517 // {
518 // featurePoints.append(pointI);
519 // }
520 // }
521 // }
522 //
523 // Info<< "Detected " << featurePoints.size()
524 // << " featurePoints out of " << points.size()
525 // << " points on feature " << i //eMesh.name()
526 // << " when using feature cos " << minCos << endl;
527 //
528 // buildTrees(i, featurePoints);
529 // }
530 //}
531 
532 
533 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
534 
536 (
537  const pointField& samples,
538  const scalarField& nearestDistSqr,
539  labelList& nearFeature,
540  List<pointIndexHit>& nearInfo,
541  vectorField& nearNormal
542 ) const
543 {
544  nearFeature.setSize(samples.size());
545  nearFeature = -1;
546  nearInfo.setSize(samples.size());
547  nearInfo = pointIndexHit();
548  nearNormal.setSize(samples.size());
549  nearNormal = vector::zero;
550 
551  forAll(edgeTrees_, featI)
552  {
553  const indexedOctree<treeDataEdge>& tree = edgeTrees_[featI];
554 
555  if (tree.shapes().size() > 0)
556  {
557  forAll(samples, sampleI)
558  {
559  const point& sample = samples[sampleI];
560 
561  scalar distSqr;
562  if (nearInfo[sampleI].hit())
563  {
564  distSqr = magSqr(nearInfo[sampleI].hitPoint()-sample);
565  }
566  else
567  {
568  distSqr = nearestDistSqr[sampleI];
569  }
570 
571  pointIndexHit info = tree.findNearest(sample, distSqr);
572 
573  if (info.hit())
574  {
575  nearFeature[sampleI] = featI;
576  nearInfo[sampleI] = pointIndexHit
577  (
578  info.hit(),
579  info.hitPoint(),
580  tree.shapes().edgeLabels()[info.index()]
581  );
582 
583  const treeDataEdge& td = tree.shapes();
584  const edge& e = td.edges()[nearInfo[sampleI].index()];
585  nearNormal[sampleI] = e.vec(td.points());
586  nearNormal[sampleI] /= mag(nearNormal[sampleI])+VSMALL;
587  }
588  }
589  }
590  }
591 }
592 
593 
595 (
596  const pointField& samples,
597  const scalarField& nearestDistSqr,
598  labelList& nearFeature,
599  List<pointIndexHit>& nearInfo,
600  vectorField& nearNormal
601 ) const
602 {
603  nearFeature.setSize(samples.size());
604  nearFeature = -1;
605  nearInfo.setSize(samples.size());
606  nearInfo = pointIndexHit();
607  nearNormal.setSize(samples.size());
608  nearNormal = vector::zero;
609 
610 
611  const PtrList<indexedOctree<treeDataEdge> >& regionTrees =
612  regionEdgeTrees();
613 
614  forAll(regionTrees, featI)
615  {
616  const indexedOctree<treeDataEdge>& regionTree = regionTrees[featI];
617 
618  forAll(samples, sampleI)
619  {
620  const point& sample = samples[sampleI];
621 
622  scalar distSqr;
623  if (nearInfo[sampleI].hit())
624  {
625  distSqr = magSqr(nearInfo[sampleI].hitPoint()-sample);
626  }
627  else
628  {
629  distSqr = nearestDistSqr[sampleI];
630  }
631 
632  // Find anything closer than current best
633  pointIndexHit info = regionTree.findNearest(sample, distSqr);
634 
635  if (info.hit())
636  {
637  const treeDataEdge& td = regionTree.shapes();
638 
639  nearFeature[sampleI] = featI;
640  nearInfo[sampleI] = pointIndexHit
641  (
642  info.hit(),
643  info.hitPoint(),
644  regionTree.shapes().edgeLabels()[info.index()]
645  );
646 
647  const edge& e = td.edges()[nearInfo[sampleI].index()];
648  nearNormal[sampleI] = e.vec(td.points());
649  nearNormal[sampleI] /= mag(nearNormal[sampleI])+VSMALL;
650  }
651  }
652  }
653 }
654 
655 
656 //void Foam::refinementFeatures::findNearestPoint
657 //(
658 // const pointField& samples,
659 // const scalarField& nearestDistSqr,
660 // labelList& nearFeature,
661 // labelList& nearIndex
662 //) const
663 //{
664 // nearFeature.setSize(samples.size());
665 // nearFeature = -1;
666 // nearIndex.setSize(samples.size());
667 // nearIndex = -1;
668 //
669 // forAll(pointTrees_, featI)
670 // {
671 // const indexedOctree<treeDataPoint>& tree = pointTrees_[featI];
672 //
673 // if (tree.shapes().pointLabels().size() > 0)
674 // {
675 // forAll(samples, sampleI)
676 // {
677 // const point& sample = samples[sampleI];
678 //
679 // scalar distSqr;
680 // if (nearFeature[sampleI] != -1)
681 // {
682 // label nearFeatI = nearFeature[sampleI];
683 // const indexedOctree<treeDataPoint>& nearTree =
684 // pointTrees_[nearFeatI];
685 // label featPointI =
686 // nearTree.shapes().pointLabels()[nearIndex[sampleI]];
687 // const point& featPt =
688 // operator[](nearFeatI).points()[featPointI];
689 // distSqr = magSqr(featPt-sample);
690 // }
691 // else
692 // {
693 // distSqr = nearestDistSqr[sampleI];
694 // }
695 //
696 // pointIndexHit info = tree.findNearest(sample, distSqr);
697 //
698 // if (info.hit())
699 // {
700 // nearFeature[sampleI] = featI;
701 // nearIndex[sampleI] = info.index();
702 // }
703 // }
704 // }
705 // }
706 //}
707 
708 
710 (
711  const pointField& samples,
712  const scalarField& nearestDistSqr,
713  labelList& nearFeature,
714  List<pointIndexHit>& nearInfo
715 ) const
716 {
717  nearFeature.setSize(samples.size());
718  nearFeature = -1;
719  nearInfo.setSize(samples.size());
720  nearInfo = pointIndexHit();
721 
722  forAll(pointTrees_, featI)
723  {
724  const indexedOctree<treeDataPoint>& tree = pointTrees_[featI];
725 
726  if (tree.shapes().pointLabels().size() > 0)
727  {
728  forAll(samples, sampleI)
729  {
730  const point& sample = samples[sampleI];
731 
732  scalar distSqr;
733  if (nearFeature[sampleI] != -1)
734  {
735  distSqr = magSqr(nearInfo[sampleI].hitPoint()-sample);
736  }
737  else
738  {
739  distSqr = nearestDistSqr[sampleI];
740  }
741 
742  pointIndexHit info = tree.findNearest(sample, distSqr);
743 
744  if (info.hit())
745  {
746  nearFeature[sampleI] = featI;
747  nearInfo[sampleI] = pointIndexHit
748  (
749  info.hit(),
750  info.hitPoint(),
751  tree.shapes().pointLabels()[info.index()]
752  );
753  }
754  }
755  }
756  }
757 }
758 
759 
760 void Foam::refinementFeatures::findHigherLevel
761 (
762  const pointField& pt,
763  const labelList& ptLevel,
764  labelList& maxLevel
765 ) const
766 {
767  // Maximum level of any shell. Start off with level of point.
768  maxLevel = ptLevel;
769 
770  forAll(*this, featI)
771  {
772  findHigherLevel(pt, featI, maxLevel);
773  }
774 }
775 
776 
778 {
779  scalar overallMax = -GREAT;
780  forAll(distances_, featI)
781  {
782  overallMax = max(overallMax, max(distances_[featI]));
783  }
784  return overallMax;
785 }
786 
787 
788 // ************************************************************************* //
const PtrList< indexedOctree< treeDataEdge > > & regionEdgeTrees() const
Simple random number generator.
Definition: Random.H:49
cachedRandom rndGen(label(0),-1)
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
void findNearestPoint(const pointField &samples, const scalarField &nearestDistSqr, labelList &nearFeature, List< pointIndexHit > &nearInfo) const
Find nearest feature point. Sets.
label size() const
Return the number of elements in the PtrList.
bool set(const label) const
Is element set.
Definition: PtrListI.H:107
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
vector point
Point is a vector.
Definition: point.H:41
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
Tuple2< scalar, label > nearInfo
Private class for finding nearest.
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const Type & shapes() const
Reference to shape.
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
const pointField & points() const
Return points.
Definition: edgeMeshI.H:39
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:75
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
const edgeList & edges() const
Definition: treeDataEdge.H:169
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
messageStream Info
static autoPtr< extendedEdgeMesh > New(const fileName &, const word &ext)
Select constructed from filename (explicit extension)
bool hit() const
Is there a hit.
const labelList & regionEdges() const
Return the feature edges which are on the boundary between.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
T & first()
Return the first element of the list.
Definition: UListI.H:117
label index() const
Return index.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
label readLabel(Istream &is)
Definition: label.H:64
const Point & hitPoint() const
Return hit point.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static bool canRead(const fileName &, const bool verbose=false)
Can we read this file format?
static autoPtr< edgeMesh > New(const fileName &, const word &ext)
Select constructed from filename (explicit extension)
Definition: edgeMeshNew.C:31
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:228
static const char nl
Definition: Ostream.H:260
const double e
Elementary charge.
Definition: doubleFloat.H:78
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
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:235
IOerror FatalIOError
treeBoundBox extend(Random &, const scalar s) const
Return slightly wider bounding box.
const pointField & points() const
Definition: treeDataEdge.H:174
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
refinementFeatures(const objectRegistry &io, const PtrList< dictionary > &featDicts)
Construct from description.
#define forAll(list, i)
Definition: UList.H:421
const List< scalarField > & distances() const
Per featureEdgeMesh the list of ranges.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
const point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
label findLower(const ListType &, typename ListType::const_reference, const label stary, const BinaryOp &bop)
Find last element < given value in sorted list and return index,.
vector vec(const pointField &) const
Return the vector (end - start)
Definition: edgeI.H:157
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
List< edge > edgeList
Definition: edgeList.H:38
Holds data for octree to work on an edges subset.
Definition: treeDataEdge.H:53
void findNearestRegionEdge(const pointField &samples, const scalarField &nearestDistSqr, labelList &nearFeature, List< pointIndexHit > &nearInfo, vectorField &nearNormal) const
Find nearest point on nearest region edge. Sets.
#define forAllReverse(list, i)
Definition: UList.H:424
error FatalError
Registry of regIOobjects.
List< label > labelList
A List of labels.
Definition: labelList.H:56
scalar maxDistance() const
Highest distance of all features.
static const Vector zero
Definition: Vector.H:80
const labelListList & levels() const
Per featureEdgeMesh the list of level.
Description of feature edges and points.
void findNearestEdge(const pointField &samples, const scalarField &nearestDistSqr, labelList &nearFeature, List< pointIndexHit > &nearInfo, vectorField &nearNormal) const
Find nearest point on nearest feature edge. Sets.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
const extendedFeatureEdgeMesh & operator[](const label) const
Return element const reference.
const edgeList & edges() const
Return edges.
Definition: edgeMeshI.H:45