edgeIntersections.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 "edgeIntersections.H"
27 #include "triSurfaceSearch.H"
28 #include "labelPairLookup.H"
29 #include "OFstream.H"
30 #include "HashSet.H"
31 #include "triSurface.H"
32 #include "pointIndexHit.H"
33 #include "treeDataTriSurface.H"
34 #include "indexedOctree.H"
35 #include "meshTools.H"
36 #include "plane.H"
37 #include "Random.H"
38 #include "unitConversion.H"
39 #include "treeBoundBox.H"
40 
41 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 defineTypeNameAndDebug(edgeIntersections, 0);
46 
47 scalar edgeIntersections::alignedCos_ = cos(degToRad(89.0));
48 }
49 
50 
51 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
52 
53 void Foam::edgeIntersections::checkEdges(const triSurface& surf)
54 {
55  const pointField& localPoints = surf.localPoints();
56  const edgeList& edges = surf.edges();
57  const labelListList& edgeFaces = surf.edgeFaces();
58 
59  treeBoundBox bb(localPoints);
60 
61  scalar minSize = SMALL * bb.minDim();
62 
63  forAll(edges, edgeI)
64  {
65  const edge& e = edges[edgeI];
66 
67  scalar eMag = e.mag(localPoints);
68 
69  if (eMag < minSize)
70  {
71  WarningIn
72  (
73  "Foam::edgeIntersections::checkEdges(const triSurface& surf)"
74  ) << "Edge " << edgeI << " vertices " << e
75  << " coords:" << localPoints[e[0]] << ' '
76  << localPoints[e[1]] << " is very small compared to bounding"
77  << " box dimensions " << bb << endl
78  << "This might lead to problems in intersection"
79  << endl;
80  }
81 
82  if (edgeFaces[edgeI].size() == 1)
83  {
84  WarningIn
85  (
86  "Foam::edgeIntersections::checkEdges(const triSurface& surf)"
87  ) << "Edge " << edgeI << " vertices " << e
88  << " coords:" << localPoints[e[0]] << ' '
89  << localPoints[e[1]] << " has only one face connected to it:"
90  << edgeFaces[edgeI] << endl
91  << "This might lead to problems in intersection"
92  << endl;
93  }
94  }
95 }
96 
97 
98 // Update intersections for selected edges.
99 void Foam::edgeIntersections::intersectEdges
100 (
101  const triSurface& surf1,
102  const pointField& points1, // surf1 meshPoints (not localPoints!)
103  const triSurfaceSearch& querySurf2,
104  const scalarField& surf1PointTol, // surf1 tolerance per point
105  const labelList& edgeLabels
106 )
107 {
108  const triSurface& surf2 = querySurf2.surface();
109  const vectorField& normals2 = surf2.faceNormals();
110 
111  const labelList& meshPoints = surf1.meshPoints();
112 
113  if (debug)
114  {
115  Pout<< "Calculating intersection of " << edgeLabels.size() << " edges"
116  << " out of " << surf1.nEdges() << " with " << surf2.size()
117  << " triangles ..." << endl;
118  }
119 
120  pointField start(edgeLabels.size());
121  pointField end(edgeLabels.size());
122  vectorField edgeDirs(edgeLabels.size());
123 
124  // Go through all edges, calculate intersections
125  forAll(edgeLabels, i)
126  {
127  label edgeI = edgeLabels[i];
128 
129  if (debug)// && (i % 1000 == 0))
130  {
131  Pout<< "Intersecting edge " << edgeI << " with surface" << endl;
132  }
133 
134  const edge& e = surf1.edges()[edgeI];
135 
136  const point& pStart = points1[meshPoints[e.start()]];
137  const point& pEnd = points1[meshPoints[e.end()]];
138 
139  const vector eVec(pEnd - pStart);
140  const vector n(eVec/(mag(eVec) + VSMALL));
141 
142  // Start tracking somewhat before pStart and up to somewhat after p1.
143  // Note that tolerances here are smaller than those used to classify
144  // hit below.
145  // This will cause this hit to be marked as degenerate and resolved
146  // later on.
147  start[i] = pStart - 0.5*surf1PointTol[e[0]]*n;
148  end[i] = pEnd + 0.5*surf1PointTol[e[1]]*n;
149 
150  edgeDirs[i] = n;
151  }
152 
153  List<List<pointIndexHit> > edgeIntersections;
154  querySurf2.findLineAll
155  (
156  start,
157  end,
158  edgeIntersections
159  );
160 
161  label nHits = 0;
162 
163  // Classify the hits
164  forAll(edgeLabels, i)
165  {
166  const label edgeI = edgeLabels[i];
167 
168  labelList& intersectionTypes = classification_[edgeI];
169  intersectionTypes.setSize(edgeIntersections[i].size(), -1);
170 
171  this->operator[](edgeI).transfer(edgeIntersections[i]);
172 
173  forAll(intersectionTypes, hitI)
174  {
175  const pointIndexHit& pHit = this->operator[](edgeI)[hitI];
176  label& hitType = intersectionTypes[hitI];
177 
178  if (!pHit.hit())
179  {
180  continue;
181  }
182 
183  const edge& e = surf1.edges()[edgeI];
184 
185  // Classify point on surface1 edge.
186  if (mag(pHit.hitPoint() - start[i]) < surf1PointTol[e[0]])
187  {
188  // Intersection is close to edge start
189  hitType = 0;
190  }
191  else if (mag(pHit.hitPoint() - end[i]) < surf1PointTol[e[1]])
192  {
193  // Intersection is close to edge end
194  hitType = 1;
195  }
196  else if (mag(edgeDirs[i] & normals2[pHit.index()]) < alignedCos_)
197  {
198  // Edge is almost coplanar with the face
199 
200  Pout<< "Flat angle edge:" << edgeI
201  << " face:" << pHit.index()
202  << " cos:" << mag(edgeDirs[i] & normals2[pHit.index()])
203  << endl;
204 
205  hitType = 2;
206  }
207 
208  if (debug)
209  {
210  Info<< " hit " << pHit << " classify = " << hitType << endl;
211  }
212 
213  nHits++;
214  }
215  }
216 
217  if (debug)
218  {
219  Pout<< "Found " << nHits << " intersections of edges with surface ..."
220  << endl;
221  }
222 }
223 
224 
225 // If edgeI intersections are close to endpoint of edge shift endpoints
226 // slightly along edge
227 // Updates
228 // - points1 with new endpoint position
229 // - affectedEdges with all edges affected by moving the point
230 // Returns true if changed anything.
231 bool Foam::edgeIntersections::inlinePerturb
232 (
233  const triSurface& surf1,
234  const scalarField& surf1PointTol, // surf1 tolerance per point
235  const label edgeI,
236  Random& rndGen,
237  pointField& points1,
238  boolList& affectedEdges
239 ) const
240 {
241  bool hasPerturbed = false;
242 
243  // Check if edge close to endpoint. Note that we only have to check
244  // the intersection closest to the edge endpoints (i.e. first and last in
245  // edgeEdgs)
246 
247  const labelList& edgeEnds = classification_[edgeI];
248 
249  if (edgeEnds.size())
250  {
251  bool perturbStart = false;
252  bool perturbEnd = false;
253 
254  // Check first intersection.
255  if (edgeEnds.first() == 0)
256  {
257  perturbStart = true;
258  }
259 
260  if (edgeEnds.last() == 1)
261  {
262  perturbEnd = true;
263  }
264 
265 
266  if (perturbStart || perturbEnd)
267  {
268  const edge& e = surf1.edges()[edgeI];
269 
270  label v0 = surf1.meshPoints()[e[0]];
271  label v1 = surf1.meshPoints()[e[1]];
272 
273  vector eVec(points1[v1] - points1[v0]);
274  vector n = eVec/mag(eVec);
275 
276  if (perturbStart)
277  {
278  // Perturb with something (hopefully) larger than tolerance.
279  scalar t = 4.0*(rndGen.scalar01() - 0.5);
280  points1[v0] += t*surf1PointTol[e[0]]*n;
281 
282  const labelList& pEdges = surf1.pointEdges()[e[0]];
283 
284  forAll(pEdges, i)
285  {
286  affectedEdges[pEdges[i]] = true;
287  }
288  }
289  if (perturbEnd)
290  {
291  // Perturb with something larger than tolerance.
292  scalar t = 4.0*(rndGen.scalar01() - 0.5);
293  points1[v1] += t*surf1PointTol[e[1]]*n;
294 
295  const labelList& pEdges = surf1.pointEdges()[e[1]];
296 
297  forAll(pEdges, i)
298  {
299  affectedEdges[pEdges[i]] = true;
300  }
301  }
302  hasPerturbed = true;
303  }
304  }
305 
306  return hasPerturbed;
307 }
308 
309 
310 // Perturb single edge endpoint when perpendicular to face
311 bool Foam::edgeIntersections::rotatePerturb
312 (
313  const triSurface& surf1,
314  const scalarField& surf1PointTol, // surf1 tolerance per point
315  const label edgeI,
316 
317  Random& rndGen,
318  pointField& points1,
319  boolList& affectedEdges
320 ) const
321 {
322  const labelList& meshPoints = surf1.meshPoints();
323 
324  const labelList& edgeEnds = classification_[edgeI];
325 
326  bool hasPerturbed = false;
327 
328  forAll(edgeEnds, i)
329  {
330  if (edgeEnds[i] == 2)
331  {
332  const edge& e = surf1.edges()[edgeI];
333 
334  // Endpoint to modify. Choose either start or end.
335  label pointI = e[rndGen.bit()];
336  //label pointI = e[0];
337 
338  // Generate random vector slightly larger than tolerance.
339  vector rndVec = rndGen.vector01() - vector(0.5, 0.5, 0.5);
340 
341  // Make sure rndVec only perp to edge
342  vector n(points1[meshPoints[e[1]]] - points1[meshPoints[e[0]]]);
343  scalar magN = mag(n) + VSMALL;
344  n /= magN;
345 
346  rndVec -= n*(n & rndVec);
347 
348  // Normalize
349  rndVec /= mag(rndVec) + VSMALL;
350 
351  // Scale to be moved by tolerance.
352  rndVec *= 0.01*magN;
353 
354  Pout<< "rotating: shifting endpoint " << meshPoints[pointI]
355  << " of edge:" << edgeI << " verts:"
356  << points1[meshPoints[e[0]]] << ' '
357  << points1[meshPoints[e[1]]]
358  << " by " << rndVec
359  << " tol:" << surf1PointTol[pointI] << endl;
360 
361  points1[meshPoints[pointI]] += rndVec;
362 
363  // Mark edges affected by change to point
364  const labelList& pEdges = surf1.pointEdges()[pointI];
365 
366  forAll(pEdges, i)
367  {
368  affectedEdges[pEdges[i]] = true;
369  }
370 
371  hasPerturbed = true;
372 
373  // Enough done for current edge; no need to test other intersections
374  // of this edge.
375  break;
376  }
377  }
378 
379  return hasPerturbed;
380 }
381 
382 
383 // Perturb edge when close to face
384 bool Foam::edgeIntersections::offsetPerturb
385 (
386  const triSurface& surf1,
387  const triSurface& surf2,
388  const label edgeI,
389 
390  Random& rndGen,
391  pointField& points1,
392  boolList& affectedEdges
393 ) const
394 {
395  const labelList& meshPoints = surf1.meshPoints();
396 
397  const edge& e = surf1.edges()[edgeI];
398 
399  const List<pointIndexHit>& hits = operator[](edgeI);
400 
401 
402  bool hasPerturbed = false;
403 
404  // For all hits on edge
405  forAll(hits, i)
406  {
407  const pointIndexHit& pHit = hits[i];
408 
409  // Classify point on face of surface2
410  label surf2FaceI = pHit.index();
411 
412  const triSurface::FaceType& f2 = surf2.localFaces()[surf2FaceI];
413  const pointField& surf2Pts = surf2.localPoints();
414 
415  const point ctr = f2.centre(surf2Pts);
416 
417  label nearType, nearLabel;
418 
419  f2.nearestPointClassify(pHit.hitPoint(), surf2Pts, nearType, nearLabel);
420 
421  if (nearType == triPointRef::POINT || nearType == triPointRef::EDGE)
422  {
423  // Shift edge towards tri centre
424  vector offset = 0.01*rndGen.scalar01()*(ctr - pHit.hitPoint());
425 
426  // shift e[0]
427  points1[meshPoints[e[0]]] += offset;
428 
429  // Mark edges affected by change to e0
430  const labelList& pEdges0 = surf1.pointEdges()[e[0]];
431 
432  forAll(pEdges0, i)
433  {
434  affectedEdges[pEdges0[i]] = true;
435  }
436 
437  // shift e[1]
438  points1[meshPoints[e[1]]] += offset;
439 
440  // Mark edges affected by change to e1
441  const labelList& pEdges1 = surf1.pointEdges()[e[1]];
442 
443  forAll(pEdges1, i)
444  {
445  affectedEdges[pEdges1[i]] = true;
446  }
447 
448  hasPerturbed = true;
449 
450  // No need to test any other hits on this edge
451  break;
452  }
453  }
454 
455  return hasPerturbed;
456 }
457 
458 
459 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
460 
461 // Construct null
463 :
464  List<List<pointIndexHit> >(),
465  classification_()
466 {}
467 
468 
469 // Construct from surface and tolerance
471 (
472  const triSurface& surf1,
473  const triSurfaceSearch& query2,
474  const scalarField& surf1PointTol
475 )
476 :
477  List<List<pointIndexHit> >(surf1.nEdges()),
478  classification_(surf1.nEdges())
479 {
480  checkEdges(surf1);
481  checkEdges(query2.surface());
482 
483  // Current set of edges to test
484  labelList edgesToTest(surf1.nEdges());
485 
486  // Start off with all edges
487  forAll(edgesToTest, i)
488  {
489  edgesToTest[i] = i;
490  }
491 
492  // Determine intersections for edgesToTest
493  intersectEdges
494  (
495  surf1,
496  surf1.points(), // surf1 meshPoints (not localPoints!)
497  query2,
498  surf1PointTol,
499  edgesToTest
500  );
501 }
502 
503 
504 // Construct from components
506 (
507  const List<List<pointIndexHit> >& intersections,
508  const labelListList& classification
509 )
510 :
511  List<List<pointIndexHit> >(intersections),
512  classification_(classification)
513 {}
514 
515 
516 // * * * * * * * * * * * * * * * Static Functions * * * * * * * * * * * * * //
517 
519 {
520  const pointField& localPoints = surf.localPoints();
521  const labelListList& pointEdges = surf.pointEdges();
522  const edgeList& edges = surf.edges();
523 
524  scalarField minLen(localPoints.size());
525 
526  forAll(minLen, pointI)
527  {
528  const labelList& pEdges = pointEdges[pointI];
529 
530  scalar minDist = GREAT;
531 
532  forAll(pEdges, i)
533  {
534  minDist = min(minDist, edges[pEdges[i]].mag(localPoints));
535  }
536 
537  minLen[pointI] = minDist;
538  }
539  return minLen;
540 }
541 
542 
543 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
544 
546 (
547  const label nIters,
548  const triSurface& surf1,
549  const triSurfaceSearch& query2,
550  const scalarField& surf1PointTol,
551  pointField& points1
552 )
553 {
554  const triSurface& surf2 = query2.surface();
555 
556  Random rndGen(356574);
557 
558  // Current set of edges to (re)test
559  labelList edgesToTest(surf1.nEdges());
560 
561  // Start off with all edges
562  forAll(edgesToTest, i)
563  {
564  edgesToTest[i] = i;
565  }
566 
567 
568  label iter = 0;
569 
570  for (; iter < nIters; iter++)
571  {
572  // Go through all edges to (re)test and perturb points if they are
573  // degenerate hits. Mark off edges that need to be recalculated.
574 
575  boolList affectedEdges(surf1.nEdges(), false);
576  label nShifted = 0;
577  label nRotated = 0;
578  label nOffset = 0;
579 
580  forAll(edgesToTest, i)
581  {
582  label edgeI = edgesToTest[i];
583 
584  // If edge not already marked for retesting
585  if (!affectedEdges[edgeI])
586  {
587  // 1. Check edges close to endpoint and perturb if necessary.
588 
589  bool shiftedEdgeEndPoints =
590  inlinePerturb
591  (
592  surf1,
593  surf1PointTol,
594  edgeI,
595  rndGen,
596  points1,
597  affectedEdges
598  );
599 
600  nShifted += (shiftedEdgeEndPoints ? 1 : 0);
601 
602  if (!shiftedEdgeEndPoints)
603  {
604  bool rotatedEdge =
605  rotatePerturb
606  (
607  surf1,
608  surf1PointTol,
609  edgeI,
610  rndGen,
611  points1,
612  affectedEdges
613  );
614 
615  nRotated += (rotatedEdge ? 1 : 0);
616 
617  if (!rotatedEdge)
618  {
619  // 2. we're sure now that the edge actually pierces the
620  // face. Now check the face for intersections close its
621  // points/edges
622 
623  bool offsetEdgePoints =
624  offsetPerturb
625  (
626  surf1,
627  surf2,
628  edgeI,
629  rndGen,
630  points1,
631  affectedEdges
632  );
633 
634  nOffset += (offsetEdgePoints ? 1 : 0);
635  }
636  }
637  }
638  }
639 
640  if (debug)
641  {
642  Pout<< "Edges to test : " << nl
643  << " total:" << edgesToTest.size() << nl
644  << " resolved by:" << nl
645  << " shifting : " << nShifted << nl
646  << " rotating : " << nRotated << nl
647  << " offsetting : " << nOffset << nl
648  << endl;
649  }
650 
651 
652  if (nShifted == 0 && nRotated == 0 && nOffset == 0)
653  {
654  // Nothing changed in current iteration. Current hit pattern is
655  // without any degenerates.
656  break;
657  }
658 
659  // Repack affected edges
660  labelList newEdgesToTest(surf1.nEdges());
661  label newEdgeI = 0;
662 
663  forAll(affectedEdges, edgeI)
664  {
665  if (affectedEdges[edgeI])
666  {
667  newEdgesToTest[newEdgeI++] = edgeI;
668  }
669  }
670  newEdgesToTest.setSize(newEdgeI);
671 
672  if (debug)
673  {
674  Pout<< "Edges to test:" << nl
675  << " was : " << edgesToTest.size() << nl
676  << " is : " << newEdgesToTest.size() << nl
677  << endl;
678  }
679 
680  // Transfer and test.
681  edgesToTest.transfer(newEdgesToTest);
682 
683  if (edgesToTest.empty())
684  {
685  FatalErrorIn("perturb") << "oops" << abort(FatalError);
686  }
687 
688  // Re intersect moved edges.
689  intersectEdges
690  (
691  surf1,
692  points1, // surf1 meshPoints (not localPoints!)
693  query2,
694  surf1PointTol,
695  edgesToTest
696  );
697  }
698 
699  return iter;
700 }
701 
702 
704 (
705  const edgeIntersections& subInfo,
706  const labelList& edgeMap,
707  const labelList& faceMap,
708  const bool merge
709 )
710 {
711  forAll(subInfo, subI)
712  {
713  const List<pointIndexHit>& subHits = subInfo[subI];
714  const labelList& subClass = subInfo.classification()[subI];
715 
716  label edgeI = edgeMap[subI];
717  List<pointIndexHit>& intersections = operator[](edgeI);
718  labelList& intersectionTypes = classification_[edgeI];
719 
720  // Count unique hits. Assume edge can hit face only once
721  label sz = 0;
722  if (merge)
723  {
724  sz = intersections.size();
725  }
726 
727  label nNew = 0;
728  forAll(subHits, i)
729  {
730  const pointIndexHit& subHit = subHits[i];
731 
732  bool foundFace = false;
733  for (label interI = 0; interI < sz; interI++)
734  {
735  if (intersections[interI].index() == faceMap[subHit.index()])
736  {
737  foundFace = true;
738  break;
739  }
740  }
741  if (!foundFace)
742  {
743  nNew++;
744  }
745  }
746 
747 
748  intersections.setSize(sz+nNew);
749  intersectionTypes.setSize(sz+nNew);
750  nNew = sz;
751 
752  forAll(subHits, i)
753  {
754  const pointIndexHit& subHit = subHits[i];
755 
756  bool foundFace = false;
757  for (label interI = 0; interI < sz; interI++)
758  {
759  if (intersections[interI].index() == faceMap[subHit.index()])
760  {
761  foundFace = true;
762  break;
763  }
764  }
765 
766  if (!foundFace)
767  {
768  intersections[nNew] = pointIndexHit
769  (
770  subHit.hit(),
771  subHit.rawPoint(),
772  faceMap[subHit.index()]
773  );
774  intersectionTypes[nNew] = subClass[i];
775  nNew++;
776  }
777  }
778  }
779 }
780 
781 
782 // ************************************************************************* //
const labelListList & classification() const
For every intersection the classification status.
Simple random number generator.
Definition: Random.H:49
cachedRandom rndGen(label(0),-1)
static scalarField minEdgeLength(const triSurface &surf)
Calculate min edge length for every surface point.
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
dimensioned< scalar > mag(const dimensioned< Type > &)
const Point & rawPoint() const
Return point with no checking.
Triangulated surface description with patch information.
Definition: triSurface.H:57
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
messageStream Info
bool hit() const
Is there a hit.
label index() const
Return index.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Namespace for OpenFOAM.
void merge(const edgeIntersections &, const labelList &edgeMap, const labelList &faceMap, const bool merge=true)
Merge (or override) edge intersection for a subset.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
label n
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
T & operator[](const label)
Return element of UList.
Definition: UListI.H:163
const Field< PointType > & localPoints() const
Return pointField of points in patch.
static const char nl
Definition: Ostream.H:260
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
#define WarningIn(functionName)
Report a warning using Foam::Warning.
Field< vector > vectorField
Specialisation of Field<T> for vector.
#define forAll(list, i)
Definition: UList.H:421
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
dimensionedScalar cos(const dimensionedScalar &ds)
Unit conversion functions.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const Vector< Cmpt > & centre(const Foam::List< Vector< Cmpt > > &) const
Return *this (used for point which is a typedef to Vector<scalar>.
Definition: VectorI.H:106
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
List< edge > edgeList
Definition: edgeList.H:38
error FatalError
const labelListList & pointEdges() const
Return point-edge addressing.
List< label > labelList
A List of labels.
Definition: labelList.H:56
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
edgeIntersections()
Construct null.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
label nEdges() const
Return number of edges in patch.
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
label removeDegenerates(const label nIters, const triSurface &surf1, const triSurfaceSearch &query2, const scalarField &surf1PointTol, pointField &points1)
Resolve ties. Shuffles points so all edge - face intersections.
const triSurface & surface() const
Return reference to the surface.
const Field< PointType > & points() const
Return reference to global points.
defineTypeNameAndDebug(combustionModel, 0)
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
Helper class to search on triSurface.