triangleI.H
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-2018 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 "IOstreams.H"
27 #include "pointHit.H"
28 #include "mathematicalConstants.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
32 template<class Point, class PointRef>
34 (
35  const Point& a,
36  const Point& b,
37  const Point& c
38 )
39 :
40  a_(a),
41  b_(b),
42  c_(c)
43 {}
44 
45 
46 template<class Point, class PointRef>
48 (
49  const UList<Point>& points,
50  const FixedList<label, 3>& indices
51 )
52 :
53  a_(points[indices[0]]),
54  b_(points[indices[1]]),
55  c_(points[indices[2]])
56 {}
57 
58 
59 
60 template<class Point, class PointRef>
62 {
63  is >> *this;
64 }
65 
66 
67 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
68 
69 template<class Point, class PointRef>
70 inline const Point& Foam::triangle<Point, PointRef>::a() const
71 {
72  return a_;
73 }
74 
75 template<class Point, class PointRef>
76 inline const Point& Foam::triangle<Point, PointRef>::b() const
77 {
78  return b_;
79 }
80 
81 template<class Point, class PointRef>
82 inline const Point& Foam::triangle<Point, PointRef>::c() const
83 {
84  return c_;
85 }
86 
87 
88 template<class Point, class PointRef>
90 {
91  return (1.0/3.0)*(a_ + b_ + c_);
92 }
93 
94 
95 template<class Point, class PointRef>
97 {
98  return 0.5*((b_ - a_)^(c_ - a_));
99 }
100 
101 
102 template<class Point, class PointRef>
103 inline Foam::scalar Foam::triangle<Point, PointRef>::mag() const
104 {
105  return Foam::mag(area());
106 }
107 
108 
109 template<class Point, class PointRef>
111 {
112  const vector a = area();
113  const scalar maga = Foam::mag(a);
114  return maga > 0 ? a/maga : Zero;
115 }
116 
117 
118 template<class Point, class PointRef>
120 {
121  scalar d1 = (c_ - a_) & (b_ - a_);
122  scalar d2 = -(c_ - b_) & (b_ - a_);
123  scalar d3 = (c_ - a_) & (c_ - b_);
124 
125  scalar c1 = d2*d3;
126  scalar c2 = d3*d1;
127  scalar c3 = d1*d2;
128 
129  scalar c = c1 + c2 + c3;
130 
131  if (Foam::mag(c) < rootVSmall)
132  {
133  // Degenerate triangle, returning centre instead of circumCentre.
134 
135  return centre();
136  }
137 
138  return
139  (
140  ((c2 + c3)*a_ + (c3 + c1)*b_ + (c1 + c2)*c_)/(2*c)
141  );
142 }
143 
144 
145 template<class Point, class PointRef>
147 {
148  scalar d1 = (c_ - a_) & (b_ - a_);
149  scalar d2 = -(c_ - b_) & (b_ - a_);
150  scalar d3 = (c_ - a_) & (c_ - b_);
151 
152  scalar denom = d2*d3 + d3*d1 + d1*d2;
153 
154  if (Foam::mag(denom) < vSmall)
155  {
156  // Degenerate triangle, returning great for circumRadius.
157 
158  return great;
159  }
160  else
161  {
162  scalar a = (d1 + d2)*(d2 + d3)*(d3 + d1) / denom;
163 
164  return 0.5*Foam::sqrt(min(great, max(0, a)));
165  }
166 }
167 
168 
169 template<class Point, class PointRef>
170 inline Foam::scalar Foam::triangle<Point, PointRef>::quality() const
171 {
172  scalar c = circumRadius();
173 
174  if (c < rootVSmall)
175  {
176  // zero circumRadius, something has gone wrong.
177  return small;
178  }
179 
180  return mag()/(Foam::sqr(c)*3.0*sqrt(3.0)/4.0);
181 }
182 
183 
184 template<class Point, class PointRef>
186 (
187  const triangle& t
188 ) const
189 {
190  return (1.0/12.0)*
191  (
192  ((t.a_ - a_) & ((b_ - a_)^(c_ - a_)))
193  + ((t.b_ - b_) & ((c_ - b_)^(t.a_ - b_)))
194  + ((c_ - t.c_) & ((t.b_ - t.c_)^(t.a_ - t.c_)))
195 
196  + ((t.a_ - a_) & ((b_ - a_)^(c_ - a_)))
197  + ((b_ - t.b_) & ((t.a_ - t.b_)^(t.c_ - t.b_)))
198  + ((c_ - t.c_) & ((b_ - t.c_)^(t.a_ - t.c_)))
199  );
200 }
201 
202 
203 template<class Point, class PointRef>
205 (
206  PointRef refPt,
207  scalar density
208 ) const
209 {
210  Point aRel = a_ - refPt;
211  Point bRel = b_ - refPt;
212  Point cRel = c_ - refPt;
213 
214  tensor V
215  (
216  aRel.x(), aRel.y(), aRel.z(),
217  bRel.x(), bRel.y(), bRel.z(),
218  cRel.x(), cRel.y(), cRel.z()
219  );
220 
221  scalar a = Foam::mag((b_ - a_)^(c_ - a_));
222 
223  tensor S = 1/24.0*(tensor::one + I);
224 
225  return
226  (
227  a*I/24.0*
228  (
229  (aRel & aRel)
230  + (bRel & bRel)
231  + (cRel & cRel)
232  + ((aRel + bRel + cRel) & (aRel + bRel + cRel))
233  )
234  - a*(V.T() & S & V)
235  )
236  *density;
237 }
238 
239 
240 template<class Point, class PointRef>
242 {
243  return barycentricToPoint(barycentric2D01(rndGen));
244 }
245 
246 
247 template<class Point, class PointRef>
249 (
250  const barycentric2D& bary
251 ) const
252 {
253  return bary[0]*a_ + bary[1]*b_ + bary[2]*c_;
254 }
255 
256 
257 template<class Point, class PointRef>
259 (
260  const point& pt
261 ) const
262 {
263  barycentric2D bary;
264  pointToBarycentric(pt, bary);
265  return bary;
266 }
267 
268 
269 template<class Point, class PointRef>
271 (
272  const point& pt,
273  barycentric2D& bary
274 ) const
275 {
276  // Reference:
277  // Real-time collision detection, Christer Ericson, 2005, p47-48
278 
279  vector v0 = b_ - a_;
280  vector v1 = c_ - a_;
281  vector v2 = pt - a_;
282 
283  scalar d00 = v0 & v0;
284  scalar d01 = v0 & v1;
285  scalar d11 = v1 & v1;
286  scalar d20 = v2 & v0;
287  scalar d21 = v2 & v1;
288 
289  scalar denom = d00*d11 - d01*d01;
290 
291  if (Foam::mag(denom) < small)
292  {
293  // Degenerate triangle, returning 1/3 barycentric coordinates.
294 
295  bary = barycentric2D(1.0/3.0, 1.0/3.0, 1.0/3.0);
296 
297  return denom;
298  }
299 
300  bary[1] = (d11*d20 - d01*d21)/denom;
301  bary[2] = (d00*d21 - d01*d20)/denom;
302  bary[0] = 1.0 - bary[1] - bary[2];
303 
304  return denom;
305 }
306 
307 
308 template<class Point, class PointRef>
310 (
311  const point& p,
312  const vector& q,
313  const intersection::algorithm alg,
314  const intersection::direction dir
315 ) const
316 {
317  // Express triangle in terms of baseVertex (point a_) and
318  // two edge vectors
319  const vector E0 = b_ - a_;
320  const vector E1 = c_ - a_;
321 
322  // Initialize intersection to miss.
323  pointHit inter(p);
324 
325  vector n(0.5*(E0 ^ E1));
326  const scalar magArea = Foam::mag(n);
327 
328  if (magArea < vSmall)
329  {
330  // Ineligible miss.
331  inter.setMiss(false);
332 
333  // The miss point is the nearest point on the triangle. Take any one.
334  inter.setPoint(a_);
335 
336  // The distance to the miss is the distance between the
337  // original point and plane of intersection. No normal so use
338  // distance from p to a_
339  inter.setDistance(Foam::mag(a_ - p));
340 
341  return inter;
342  }
343 
344  const vector q1 = q/Foam::mag(q);
345 
346  if (dir == intersection::CONTACT_SPHERE)
347  {
348  n /= magArea;
349 
350  return ray(p, q1 - n, alg, intersection::VECTOR);
351  }
352 
353  // Intersection point with triangle plane
354  point pInter;
355 
356  // Is intersection point inside triangle
357  bool hit;
358  {
359  // Reuse the fast ray intersection routine below in FULL_RAY
360  // mode since the original intersection routine has rounding problems.
361  pointHit fastInter = intersection(p, q1, intersection::FULL_RAY);
362  hit = fastInter.hit();
363 
364  if (hit)
365  {
366  pInter = fastInter.rawPoint();
367  }
368  else
369  {
370  // Calculate intersection of ray with triangle plane
371  vector v = a_ - p;
372  pInter = p + (q1&v)*q1;
373  }
374  }
375 
376  // Distance to intersection point
377  const scalar dist = q1 & (pInter - p);
378 
379  const scalar planarPointTol =
380  Foam::min
381  (
382  Foam::min
383  (
384  Foam::mag(E0),
385  Foam::mag(E1)
386  ),
387  Foam::mag(c_ - b_)
388  )*intersection::planarTol();
389 
390  bool eligible =
391  alg == intersection::FULL_RAY
392  || (alg == intersection::HALF_RAY && dist > -planarPointTol)
393  || (
394  alg == intersection::VISIBLE
395  && ((q1 & area()) < -vSmall)
396  );
397 
398  if (hit && eligible)
399  {
400  // Hit. Set distance to intersection.
401  inter.setHit();
402  inter.setPoint(pInter);
403  inter.setDistance(dist);
404  }
405  else
406  {
407  // Miss or ineligible hit.
408  inter.setMiss(eligible);
409 
410  // The miss point is the nearest point on the triangle
411  inter.setPoint(nearestPoint(p).rawPoint());
412 
413  // The distance to the miss is the distance between the
414  // original point and plane of intersection
415  inter.setDistance(Foam::mag(pInter - p));
416  }
417 
418  return inter;
419 }
420 
421 
422 // From "Fast, Minimum Storage Ray/Triangle Intersection"
423 // Moeller/Trumbore.
424 template<class Point, class PointRef>
426 (
427  const point& orig,
428  const vector& dir,
429  const intersection::algorithm alg,
430  const scalar tol
431 ) const
432 {
433  const vector edge1 = b_ - a_;
434  const vector edge2 = c_ - a_;
435 
436  // begin calculating determinant - also used to calculate U parameter
437  const vector pVec = dir ^ edge2;
438 
439  // if determinant is near zero, ray lies in plane of triangle
440  const scalar det = edge1 & pVec;
441 
442  // Initialise to miss
443  pointHit intersection(false, Zero, great, false);
444 
445  if (alg == intersection::VISIBLE)
446  {
447  // Culling branch
448  if (det < rootVSmall)
449  {
450  // Ray on wrong side of triangle. Return miss
451  return intersection;
452  }
453  }
454  else if (alg == intersection::HALF_RAY || alg == intersection::FULL_RAY)
455  {
456  // Non-culling branch
457  if (det > -rootVSmall && det < rootVSmall)
458  {
459  // Ray parallel to triangle. Return miss
460  return intersection;
461  }
462  }
463 
464  const scalar inv_det = 1.0 / det;
465 
466  /* calculate distance from a_ to ray origin */
467  const vector tVec = orig-a_;
468 
469  /* calculate U parameter and test bounds */
470  const scalar u = (tVec & pVec)*inv_det;
471 
472  if (u < -tol || u > 1.0+tol)
473  {
474  // return miss
475  return intersection;
476  }
477 
478  /* prepare to test V parameter */
479  const vector qVec = tVec ^ edge1;
480 
481  /* calculate V parameter and test bounds */
482  const scalar v = (dir & qVec) * inv_det;
483 
484  if (v < -tol || u + v > 1.0+tol)
485  {
486  // return miss
487  return intersection;
488  }
489 
490  /* calculate t, scale parameters, ray intersects triangle */
491  const scalar t = (edge2 & qVec) * inv_det;
492 
493  if (alg == intersection::HALF_RAY && t < -tol)
494  {
495  // Wrong side of orig. Return miss
496  return intersection;
497  }
498 
499  intersection.setHit();
500  intersection.setPoint(a_ + u*edge1 + v*edge2);
501  intersection.setDistance(t);
502 
503  return intersection;
504 }
505 
506 
507 template<class Point, class PointRef>
509 (
510  const point& p,
511  label& nearType,
512  label& nearLabel
513 ) const
514 {
515  // Adapted from:
516  // Real-time collision detection, Christer Ericson, 2005, p136-142
517 
518  // Check if P in vertex region outside A
519  vector ab = b_ - a_;
520  vector ac = c_ - a_;
521  vector ap = p - a_;
522 
523  scalar d1 = ab & ap;
524  scalar d2 = ac & ap;
525 
526  if (d1 <= 0.0 && d2 <= 0.0)
527  {
528  // barycentric coordinates (1, 0, 0)
529 
530  nearType = POINT;
531  nearLabel = 0;
532  return pointHit(false, a_, Foam::mag(a_ - p), true);
533  }
534 
535  // Check if P in vertex region outside B
536  vector bp = p - b_;
537  scalar d3 = ab & bp;
538  scalar d4 = ac & bp;
539 
540  if (d3 >= 0.0 && d4 <= d3)
541  {
542  // barycentric coordinates (0, 1, 0)
543 
544  nearType = POINT;
545  nearLabel = 1;
546  return pointHit(false, b_, Foam::mag(b_ - p), true);
547  }
548 
549  // Check if P in edge region of AB, if so return projection of P onto AB
550  scalar vc = d1*d4 - d3*d2;
551 
552  if (vc <= 0.0 && d1 >= 0.0 && d3 <= 0.0)
553  {
554  if ((d1 - d3) < rootVSmall)
555  {
556  // Degenerate triangle, for d1 = d3, a_ and b_ are likely coincident
557  nearType = POINT;
558  nearLabel = 0;
559  return pointHit(false, a_, Foam::mag(a_ - p), true);
560  }
561 
562  // barycentric coordinates (1-v, v, 0)
563  scalar v = d1/(d1 - d3);
564 
565  point nearPt = a_ + v*ab;
566  nearType = EDGE;
567  nearLabel = 0;
568  return pointHit(false, nearPt, Foam::mag(nearPt - p), true);
569  }
570 
571  // Check if P in vertex region outside C
572  vector cp = p - c_;
573  scalar d5 = ab & cp;
574  scalar d6 = ac & cp;
575 
576  if (d6 >= 0.0 && d5 <= d6)
577  {
578  // barycentric coordinates (0, 0, 1)
579 
580  nearType = POINT;
581  nearLabel = 2;
582  return pointHit(false, c_, Foam::mag(c_ - p), true);
583  }
584 
585  // Check if P in edge region of AC, if so return projection of P onto AC
586  scalar vb = d5*d2 - d1*d6;
587 
588  if (vb <= 0.0 && d2 >= 0.0 && d6 <= 0.0)
589  {
590  if ((d2 - d6) < rootVSmall)
591  {
592  // Degenerate triangle, for d2 = d6, a_ and c_ are likely coincident
593  nearType = POINT;
594  nearLabel = 0;
595  return pointHit(false, a_, Foam::mag(a_ - p), true);
596  }
597 
598  // barycentric coordinates (1-w, 0, w)
599  scalar w = d2/(d2 - d6);
600 
601  point nearPt = a_ + w*ac;
602  nearType = EDGE;
603  nearLabel = 2;
604  return pointHit(false, nearPt, Foam::mag(nearPt - p), true);
605  }
606 
607  // Check if P in edge region of BC, if so return projection of P onto BC
608  scalar va = d3*d6 - d5*d4;
609 
610  if (va <= 0.0 && (d4 - d3) >= 0.0 && (d5 - d6) >= 0.0)
611  {
612  if (((d4 - d3) + (d5 - d6)) < rootVSmall)
613  {
614  // Degenerate triangle, for (d4 - d3) = (d6 - d5), b_ and c_ are
615  // likely coincident
616  nearType = POINT;
617  nearLabel = 1;
618  return pointHit(false, b_, Foam::mag(b_ - p), true);
619  }
620 
621  // barycentric coordinates (0, 1-w, w)
622  scalar w = (d4 - d3)/((d4 - d3) + (d5 - d6));
623 
624  point nearPt = b_ + w*(c_ - b_);
625  nearType = EDGE;
626  nearLabel = 1;
627  return pointHit(false, nearPt, Foam::mag(nearPt - p), true);
628  }
629 
630  // P inside face region. Compute Q through its barycentric
631  // coordinates (u, v, w)
632 
633  if ((va + vb + vc) < rootVSmall)
634  {
635  // Degenerate triangle, return the centre because no edge or points are
636  // closest
637  point nearPt = centre();
638  nearType = NONE,
639  nearLabel = -1;
640  return pointHit(true, nearPt, Foam::mag(nearPt - p), false);
641  }
642 
643  scalar denom = 1.0/(va + vb + vc);
644  scalar v = vb * denom;
645  scalar w = vc * denom;
646 
647  // = u*a + v*b + w*c, u = va*denom = 1.0 - v - w
648 
649  point nearPt = a_ + ab*v + ac*w;
650  nearType = NONE,
651  nearLabel = -1;
652  return pointHit(true, nearPt, Foam::mag(nearPt - p), false);
653 }
654 
655 
656 template<class Point, class PointRef>
658 (
659  const point& p
660 ) const
661 {
662  // Dummy labels
663  label nearType = -1;
664  label nearLabel = -1;
665 
666  return nearestPointClassify(p, nearType, nearLabel);
667 }
668 
669 
670 template<class Point, class PointRef>
672 (
673  const point& p,
674  label& nearType,
675  label& nearLabel
676 ) const
677 {
678  return nearestPointClassify(p, nearType, nearLabel).hit();
679 }
680 
681 
682 template<class Point, class PointRef>
684 (
685  const linePointRef& ln,
686  pointHit& lnInfo
687 ) const
688 {
689  vector q = ln.vec();
690  pointHit triInfo
691  (
693  (
694  ln.start(),
695  q,
696  intersection::FULL_RAY
697  )
698  );
699 
700  if (triInfo.hit())
701  {
702  // Line hits triangle. Find point on line.
703  if (triInfo.distance() > 1)
704  {
705  // Hit beyond endpoint
706  lnInfo.setMiss(true);
707  lnInfo.setPoint(ln.end());
708  scalar dist = Foam::mag(triInfo.hitPoint()-lnInfo.missPoint());
709  lnInfo.setDistance(dist);
710  triInfo.setMiss(true);
711  triInfo.setDistance(dist);
712  }
713  else if (triInfo.distance() < 0)
714  {
715  // Hit beyond startpoint
716  lnInfo.setMiss(true);
717  lnInfo.setPoint(ln.start());
718  scalar dist = Foam::mag(triInfo.hitPoint()-lnInfo.missPoint());
719  lnInfo.setDistance(dist);
720  triInfo.setMiss(true);
721  triInfo.setDistance(dist);
722  }
723  else
724  {
725  // Hit on line
726  lnInfo.setHit();
727  lnInfo.setPoint(triInfo.hitPoint());
728  lnInfo.setDistance(0.0);
729  triInfo.setDistance(0.0);
730  }
731  }
732  else
733  {
734  // Line skips triangle. See which triangle edge it gets closest to
735 
736  point nearestEdgePoint;
737  point nearestLinePoint;
738  // label minEdgeIndex = 0;
739  scalar minDist = ln.nearestDist
740  (
741  linePointRef(a_, b_),
742  nearestLinePoint,
743  nearestEdgePoint
744  );
745 
746  {
747  point linePoint;
748  point triEdgePoint;
749  scalar dist = ln.nearestDist
750  (
751  linePointRef(b_, c_),
752  linePoint,
753  triEdgePoint
754  );
755  if (dist < minDist)
756  {
757  minDist = dist;
758  nearestEdgePoint = triEdgePoint;
759  nearestLinePoint = linePoint;
760  // minEdgeIndex = 1;
761  }
762  }
763 
764  {
765  point linePoint;
766  point triEdgePoint;
767  scalar dist = ln.nearestDist
768  (
769  linePointRef(c_, a_),
770  linePoint,
771  triEdgePoint
772  );
773  if (dist < minDist)
774  {
775  minDist = dist;
776  nearestEdgePoint = triEdgePoint;
777  nearestLinePoint = linePoint;
778  // minEdgeIndex = 2;
779  }
780  }
781 
782  lnInfo.setDistance(minDist);
783  triInfo.setDistance(minDist);
784  triInfo.setMiss(false);
785  triInfo.setPoint(nearestEdgePoint);
786 
787  // Convert point on line to pointHit
788  if (Foam::mag(nearestLinePoint-ln.start()) < small)
789  {
790  lnInfo.setMiss(true);
791  lnInfo.setPoint(ln.start());
792  }
793  else if (Foam::mag(nearestLinePoint-ln.end()) < small)
794  {
795  lnInfo.setMiss(true);
796  lnInfo.setPoint(ln.end());
797  }
798  else
799  {
800  lnInfo.setHit();
801  lnInfo.setPoint(nearestLinePoint);
802  }
803  }
804  return triInfo;
805 }
806 
807 
808 // * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
809 
810 template<class Point, class PointRef>
811 inline Foam::Istream& Foam::operator>>
812 (
813  Istream& is,
815 )
816 {
817  is.readBegin("triangle");
818  is >> t.a_ >> t.b_ >> t.c_;
819  is.readEnd("triangle");
820 
821  is.check("Istream& operator>>(Istream&, triangle&)");
822  return is;
823 }
824 
825 
826 template<class Point, class PointRef>
827 inline Foam::Ostream& Foam::operator<<
828 (
829  Ostream& os,
831 )
832 {
833  os << nl
834  << token::BEGIN_LIST
835  << t.a_ << token::SPACE
836  << t.b_ << token::SPACE
837  << t.c_
838  << token::END_LIST;
839 
840  return os;
841 }
842 
843 
844 // ************************************************************************* //
barycentric2D barycentric2D01(Random &rndGen)
Generate a random barycentric coordinate within the unit triangle.
Definition: barycentric2D.C:31
scalar sweptVol(const triangle &t) const
Return swept-volume.
Definition: triangleI.H:186
Istream & readBegin(const char *funcName)
Definition: Istream.C:86
A triangle primitive used to calculate face areas and swept volumes.
Definition: triangle.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
A line primitive.
Definition: line.H:56
const Point & c() const
Return third vertex.
Definition: triangleI.H:82
Barycentric2D< scalar > barycentric2D
A scalar version of the templated Barycentric2D.
Definition: barycentric2D.H:45
barycentric2D pointToBarycentric(const point &pt) const
Calculate the barycentric coordinates from the given point.
Definition: triangleI.H:259
scalar mag() const
Return scalar magnitude.
Definition: triangleI.H:103
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
scalar quality() const
Return quality: Ratio of triangle and circum-circle.
Definition: triangleI.H:170
dimensionedSymmTensor sqr(const dimensionedVector &dv)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
const dimensionedScalar c2
Second radiation constant: default SI units: [m.K].
dimensionedScalar sqrt(const dimensionedScalar &ds)
pointHit nearestPointClassify(const point &p, label &nearType, label &nearLabel) const
Find the nearest point to p on the triangle and classify it:
Definition: triangleI.H:509
dimensionedScalar det(const dimensionedSphericalTensor &dt)
scalar minDist(const List< pointIndexHit > &hitList)
const Point & missPoint() const
Return miss point.
Definition: PointHit.H:145
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
Templated 2D Barycentric derived from VectorSpace. Has 3 components, one of which is redundant...
Definition: Barycentric2D.H:50
vector area() const
Return vector area.
Definition: triangleI.H:96
Random rndGen(label(0))
bool classify(const point &p, label &nearType, label &nearLabel) const
Classify nearest point to p in triangle plane.
Definition: triangleI.H:672
pointHit nearestPoint(const point &p) const
Return nearest point to p on triangle.
Definition: triangleI.H:658
void setHit()
Definition: PointHit.H:169
Foam::intersection.
Definition: intersection.H:49
void setDistance(const scalar d)
Definition: PointHit.H:186
line< point, const point & > linePointRef
Line using referred points.
Definition: linePointRef.H:45
pointHit ray(const point &p, const vector &q, const intersection::algorithm=intersection::FULL_RAY, const intersection::direction dir=intersection::VECTOR) const
Return point intersection with a ray.
Definition: triangleI.H:310
static const Identity< scalar > I
Definition: Identity.H:93
const Point & hitPoint() const
Return hit point.
Definition: PointHit.H:126
PointRef start() const
Return first vertex.
Definition: lineI.H:60
Point centre() const
Return centre (centroid)
Definition: triangleI.H:89
bool hit() const
Is there a hit.
Definition: PointHit.H:120
static const zero Zero
Definition: zero.H:97
scalar circumRadius() const
Return circum-radius.
Definition: triangleI.H:146
tensor inertia(PointRef refPt=vector::zero, scalar density=1.0) const
Return the inertia tensor, with optional reference.
Definition: triangleI.H:205
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:61
Random number generator.
Definition: Random.H:57
vector normal() const
Return unit normal.
Definition: triangleI.H:110
void setMiss(const bool eligible)
Definition: PointHit.H:175
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
triangle(const Point &a, const Point &b, const Point &c)
Construct from three points.
Definition: triangleI.H:34
const Point & rawPoint() const
Return point with no checking.
Definition: PointHit.H:158
static const char nl
Definition: Ostream.H:265
const volScalarField & cp
Point barycentricToPoint(const barycentric2D &bary) const
Calculate the point from the given barycentric coordinates.
Definition: triangleI.H:249
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const Point & a() const
Return first vertex.
Definition: triangleI.H:70
Point vec() const
Return start-end vector.
Definition: lineI.H:87
Point circumCentre() const
Return circum-centre.
Definition: triangleI.H:119
const Point & b() const
Return second vertex.
Definition: triangleI.H:76
PointRef end() const
Return second vertex.
Definition: lineI.H:66
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
Point randomPoint(Random &rndGen) const
Return a random point on the triangle from a uniform.
Definition: triangleI.H:241
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
PointHit< point > pointHit
Definition: pointHit.H:41
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
volScalarField & p
pointHit intersection(const point &p, const vector &q, const intersection::algorithm alg, const scalar tol=0.0) const
Fast intersection with a ray.
Definition: triangleI.H:426
void setPoint(const Point &p)
Definition: PointHit.H:181
PointHit< Point > nearestDist(const Point &p) const
Return nearest distance to line from a given point.
Definition: lineI.H:95