tetrahedronI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2018 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "triangle.H"
30#include "IOstreams.H"
31#include "tetPoints.H"
32#include "plane.H"
33
34// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35
36template<class Point, class PointRef>
38(
39 const Point& a,
40 const Point& b,
41 const Point& c,
42 const Point& d
43)
44:
45 a_(a),
46 b_(b),
47 c_(c),
48 d_(d)
49{}
50
51
52template<class Point, class PointRef>
54(
55 const UList<Point>& points,
56 const FixedList<label, 4>& indices
57)
58:
59 a_(points[indices[0]]),
60 b_(points[indices[1]]),
61 c_(points[indices[2]]),
62 d_(points[indices[3]])
63{}
64
65
66template<class Point, class PointRef>
68{
69 is >> *this;
70}
71
72
73// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
74
75template<class Point, class PointRef>
77{
78 return a_;
79}
80
81
82template<class Point, class PointRef>
84{
85 return b_;
86}
87
88
89template<class Point, class PointRef>
91{
92 return c_;
93}
94
95
96template<class Point, class PointRef>
98{
99 return d_;
100}
101
102
103template<class Point, class PointRef>
105(
106 const label facei
107) const
108{
109 // Warning. Ordering of faces needs to be the same for a tetrahedron
110 // class, a tetrahedron cell shape model and a tetCell
111
112 if (facei == 0)
113 {
114 return triPointRef(b_, c_, d_);
115 }
116 else if (facei == 1)
117 {
118 return triPointRef(a_, d_, c_);
119 }
120 else if (facei == 2)
121 {
122 return triPointRef(a_, b_, d_);
123 }
124 else if (facei == 3)
125 {
126 return triPointRef(a_, c_, b_);
127 }
128
130 << "Face index (" << facei << ") out of range 0..3\n"
131 << abort(FatalError);
132 return triPointRef(b_, c_, d_);
133}
134
135
136template<class Point, class PointRef>
138{
139 return triangle<Point, PointRef>(b_, c_, d_).areaNormal();
140}
141
142
143template<class Point, class PointRef>
145{
146 return triangle<Point, PointRef>(a_, d_, c_).areaNormal();
147}
148
149
150template<class Point, class PointRef>
152{
153 return triangle<Point, PointRef>(a_, b_, d_).areaNormal();
154}
155
156
157template<class Point, class PointRef>
159{
160 return triangle<Point, PointRef>(a_, c_, b_).areaNormal();
161}
162
163
164template<class Point, class PointRef>
166{
167 return 0.25*(a_ + b_ + c_ + d_);
168}
169
170
171template<class Point, class PointRef>
173{
174 return (1.0/6.0)*(((b_ - a_) ^ (c_ - a_)) & (d_ - a_));
175}
176
177
178template<class Point, class PointRef>
180{
181 vector a = b_ - a_;
182 vector b = c_ - a_;
183 vector c = d_ - a_;
184
185 scalar lambda = magSqr(c) - (a & c);
186 scalar mu = magSqr(b) - (a & b);
187
188 vector ba = b ^ a;
189 vector ca = c ^ a;
190
191 vector num = lambda*ba - mu*ca;
192 scalar denom = (c & ba);
193
194 if (Foam::mag(denom) < ROOTVSMALL)
195 {
196 // Degenerate tetrahedron, returning centre instead of circumCentre.
197
198 return centre();
199 }
200
201 return a_ + 0.5*(a + num/denom);
202}
203
204
205template<class Point, class PointRef>
207{
208 vector a = b_ - a_;
209 vector b = c_ - a_;
210 vector c = d_ - a_;
211
212 scalar lambda = magSqr(c) - (a & c);
213 scalar mu = magSqr(b) - (a & b);
214
215 vector ba = b ^ a;
216 vector ca = c ^ a;
217
218 vector num = lambda*ba - mu*ca;
219 scalar denom = (c & ba);
220
221 if (Foam::mag(denom) < ROOTVSMALL)
222 {
223 // Degenerate tetrahedron, returning GREAT for circumRadius.
224 return GREAT;
225 }
226
227 return Foam::mag(0.5*(a + num/denom));
228}
229
230
231template<class Point, class PointRef>
233{
234 return
235 mag()
236 /(
237 8.0/(9.0*sqrt(3.0))
238 *pow3(min(circumRadius(), GREAT))
239 + ROOTVSMALL
240 );
241}
242
243
244template<class Point, class PointRef>
246(
248) const
249{
250 return barycentricToPoint(barycentric01(rndGen));
251}
252
253
254template<class Point, class PointRef>
256(
257 const barycentric& bary
258) const
259{
260 return bary[0]*a_ + bary[1]*b_ + bary[2]*c_ + bary[3]*d_;
261}
262
263
264template<class Point, class PointRef>
266(
267 const point& pt
268) const
269{
270 barycentric bary;
271 pointToBarycentric(pt, bary);
272 return bary;
273}
274
275
276template<class Point, class PointRef>
278(
279 const point& pt,
280 barycentric& bary
281) const
282{
283 // Reference:
284 // http://en.wikipedia.org/wiki/Barycentric_coordinate_system_(mathematics)
285
286 vector e0(a_ - d_);
287 vector e1(b_ - d_);
288 vector e2(c_ - d_);
289
290 tensor t
291 (
292 e0.x(), e1.x(), e2.x(),
293 e0.y(), e1.y(), e2.y(),
294 e0.z(), e1.z(), e2.z()
295 );
296
297 scalar detT = det(t);
298
299 if (Foam::mag(detT) < SMALL)
300 {
301 // Degenerate tetrahedron, returning 1/4 barycentric coordinates
302
303 bary = barycentric(0.25, 0.25, 0.25, 0.25);
304
305 return detT;
306 }
307
308 vector res = inv(t, detT) & (pt - d_);
309
310 bary[0] = res.x();
311 bary[1] = res.y();
312 bary[2] = res.z();
313 bary[3] = 1 - cmptSum(res);
314
315 return detT;
316}
317
318
319template<class Point, class PointRef>
321(
322 const point& p
323) const
324{
325 // Adapted from:
326 // Real-time collision detection, Christer Ericson, 2005, p142-144
327
328 // Assuming initially that the point is inside all of the
329 // halfspaces, so closest to itself.
330
331 point closestPt = p;
332
333 scalar minOutsideDistance = VGREAT;
334
335 bool inside = true;
336
337 if (((p - b_) & Sa()) >= 0)
338 {
339 // p is outside halfspace plane of tri
341
342 inside = false;
343
344 if (info.distance() < minOutsideDistance)
345 {
346 closestPt = info.rawPoint();
347
348 minOutsideDistance = info.distance();
349 }
350 }
351
352 if (((p - a_) & Sb()) >= 0)
353 {
354 // p is outside halfspace plane of tri
356
357 inside = false;
358
359 if (info.distance() < minOutsideDistance)
360 {
361 closestPt = info.rawPoint();
362
363 minOutsideDistance = info.distance();
364 }
365 }
366
367 if (((p - a_) & Sc()) >= 0)
368 {
369 // p is outside halfspace plane of tri
371
372 inside = false;
373
374 if (info.distance() < minOutsideDistance)
375 {
376 closestPt = info.rawPoint();
377
378 minOutsideDistance = info.distance();
379 }
380 }
381
382 if (((p - a_) & Sd()) >= 0)
383 {
384 // p is outside halfspace plane of tri
386
387 inside = false;
388
389 if (info.distance() < minOutsideDistance)
390 {
391 closestPt = info.rawPoint();
392
393 minOutsideDistance = info.distance();
394 }
395 }
396
397 // If the point is inside, then the distance to the closest point
398 // is zero
399 if (inside)
400 {
401 minOutsideDistance = 0;
402 }
403
404 return pointHit
405 (
406 inside,
407 closestPt,
408 minOutsideDistance,
409 !inside
410 );
411}
412
413
414template<class Point, class PointRef>
416{
417 // For robustness, assuming that the point is in the tet unless
418 // "definitively" shown otherwise by obtaining a positive dot
419 // product greater than a tolerance of SMALL.
420
421 // The tet is defined: tet(Cc, tetBasePt, pA, pB) where the area normal
422 // vectors and base points for the half-space planes are:
423 // area[0] = Sa();
424 // area[1] = Sb();
425 // area[2] = Sc();
426 // area[3] = Sd();
427 // planeBase[0] = tetBasePt = b_
428 // planeBase[1] = ptA = c_
429 // planeBase[2] = tetBasePt = b_
430 // planeBase[3] = tetBasePt = b_
431
432 vector n = Zero;
433
434 {
435 // 0, a
436 const point& basePt = b_;
437
438 n = Sa();
439 n /= (Foam::mag(n) + VSMALL);
440
441 if (((pt - basePt) & n) > SMALL)
442 {
443 return false;
444 }
445 }
446
447 {
448 // 1, b
449 const point& basePt = c_;
450
451 n = Sb();
452 n /= (Foam::mag(n) + VSMALL);
453
454 if (((pt - basePt) & n) > SMALL)
455 {
456 return false;
457 }
458 }
459
460 {
461 // 2, c
462 const point& basePt = b_;
463
464 n = Sc();
465 n /= (Foam::mag(n) + VSMALL);
466
467 if (((pt - basePt) & n) > SMALL)
468 {
469 return false;
470 }
471 }
472
473 {
474 // 3, d
475 const point& basePt = b_;
476
477 n = Sd();
478 n /= (Foam::mag(n) + VSMALL);
479
480 if (((pt - basePt) & n) > SMALL)
481 {
482 return false;
483 }
484 }
485
486 return true;
487}
488
489
490template<class Point, class PointRef>
492(
493 const tetPoints&
494)
495{}
496
497
498template<class Point, class PointRef>
500:
501 vol_(0.0)
502{}
503
504
505template<class Point, class PointRef>
507(
508 const tetPoints& tet
509)
510{
511 vol_ += tet.tet().mag();
512}
513
514
515template<class Point, class PointRef>
517(
519 label& nTets
520)
521:
522 tets_(tets),
523 nTets_(nTets)
524{}
525
526
527template<class Point, class PointRef>
529(
530 const tetPoints& tet
531)
532{
533 tets_[nTets_++] = tet;
534}
535
536
537template<class Point, class PointRef>
539(
540 const FixedList<scalar, 4>& d,
541 const tetPoints& t,
542 const label negI,
543 const label posI
544)
545{
546 return
547 (d[posI]*t[negI] - d[negI]*t[posI])
548 / (-d[negI]+d[posI]);
549}
550
551
552template<class Point, class PointRef>
553template<class TetOp>
555(
557 TetOp& op
558)
559{
560 op(tetPoints(points[1], points[3], points[2], points[0]));
561 op(tetPoints(points[1], points[2], points[3], points[4]));
562 op(tetPoints(points[4], points[2], points[3], points[5]));
563}
564
565
566template<class Point, class PointRef>
567template<class AboveTetOp, class BelowTetOp>
569(
570 const plane& pln,
571 const tetPoints& tet,
572 AboveTetOp& aboveOp,
573 BelowTetOp& belowOp
574)
575{
576 // Distance to plane
577 FixedList<scalar, 4> d;
578 label nPos = 0;
579 forAll(tet, i)
580 {
581 d[i] = pln.signedDistance(tet[i]);
582 if (d[i] > 0)
583 {
584 ++nPos;
585 }
586 }
587
588 if (nPos == 4)
589 {
590 aboveOp(tet);
591 }
592 else if (nPos == 3)
593 {
594 // Sliced into below tet and above prism.
595 // Prism gets split into two tets.
596
597 // Find the below tet
598 label i0 = -1;
599 forAll(d, i)
600 {
601 if (d[i] <= 0)
602 {
603 i0 = i;
604 break;
605 }
606 }
607
608 label i1 = d.fcIndex(i0);
609 label i2 = d.fcIndex(i1);
610 label i3 = d.fcIndex(i2);
611
612 point p01(planeIntersection(d, tet, i0, i1));
613 point p02(planeIntersection(d, tet, i0, i2));
614 point p03(planeIntersection(d, tet, i0, i3));
615
616 // i0 = tetCell vertex 0: p01,p02,p03 outwards pointing triad
617 // ,, 1 : ,, inwards pointing triad
618 // ,, 2 : ,, outwards pointing triad
619 // ,, 3 : ,, inwards pointing triad
620
621 //Pout<< "Split 3pos tet " << tet << " d:" << d << " into" << nl;
622
623 if (i0 == 0 || i0 == 2)
624 {
625 tetPoints t(tet[i0], p01, p02, p03);
626 //Pout<< " belowtet:" << t << " around i0:" << i0 << endl;
627 //checkTet(t, "nPos 3, belowTet i0==0 or 2");
628 belowOp(t);
629
630 // Prism
631 FixedList<point, 6> p
632 (
633 {
634 tet[i1],
635 tet[i3],
636 tet[i2],
637 p01,
638 p03,
639 p02
640 }
641 );
642
643 //Pout<< " aboveprism:" << p << endl;
644 decomposePrism(p, aboveOp);
645 }
646 else
647 {
648 tetPoints t(p01, p02, p03, tet[i0]);
649 //Pout<< " belowtet:" << t << " around i0:" << i0 << endl;
650 //checkTet(t, "nPos 3, belowTet i0==1 or 3");
651 belowOp(t);
652
653 // Prism
654 FixedList<point, 6> p
655 (
656 {
657 tet[i3],
658 tet[i1],
659 tet[i2],
660 p03,
661 p01,
662 p02
663 }
664 );
665 //Pout<< " aboveprism:" << p << endl;
666 decomposePrism(p, aboveOp);
667 }
668 }
669 else if (nPos == 2)
670 {
671 // Tet cut into two prisms. Determine the positive one.
672 label pos0 = -1;
673 label pos1 = -1;
674 forAll(d, i)
675 {
676 if (d[i] > 0)
677 {
678 if (pos0 == -1)
679 {
680 pos0 = i;
681 }
682 else
683 {
684 pos1 = i;
685 }
686 }
687 }
688
689 //Pout<< "Split 2pos tet " << tet << " d:" << d
690 // << " around pos0:" << pos0 << " pos1:" << pos1
691 // << " neg0:" << neg0 << " neg1:" << neg1 << " into" << nl;
692
693 const edge posEdge(pos0, pos1);
694
695 if (posEdge == edge(0, 1))
696 {
697 point p02(planeIntersection(d, tet, 0, 2));
698 point p03(planeIntersection(d, tet, 0, 3));
699 point p12(planeIntersection(d, tet, 1, 2));
700 point p13(planeIntersection(d, tet, 1, 3));
701 // Split the resulting prism
702 {
703 FixedList<point, 6> p
704 (
705 {
706 tet[0],
707 p02,
708 p03,
709 tet[1],
710 p12,
711 p13
712 }
713 );
714
715 //Pout<< " 01 aboveprism:" << p << endl;
716 decomposePrism(p, aboveOp);
717 }
718 {
719 FixedList<point, 6> p
720 (
721 {
722 tet[2],
723 p02,
724 p12,
725 tet[3],
726 p03,
727 p13
728 }
729 );
730
731 //Pout<< " 01 belowprism:" << p << endl;
732 decomposePrism(p, belowOp);
733 }
734 }
735 else if (posEdge == edge(1, 2))
736 {
737 point p01(planeIntersection(d, tet, 0, 1));
738 point p13(planeIntersection(d, tet, 1, 3));
739 point p02(planeIntersection(d, tet, 0, 2));
740 point p23(planeIntersection(d, tet, 2, 3));
741 // Split the resulting prism
742 {
743 FixedList<point, 6> p
744 (
745 {
746 tet[1],
747 p01,
748 p13,
749 tet[2],
750 p02,
751 p23
752 }
753 );
754
755 //Pout<< " 12 aboveprism:" << p << endl;
756 decomposePrism(p, aboveOp);
757 }
758 {
759 FixedList<point, 6> p
760 (
761 {
762 tet[3],
763 p23,
764 p13,
765 tet[0],
766 p02,
767 p01
768 }
769 );
770
771 //Pout<< " 12 belowprism:" << p << endl;
772 decomposePrism(p, belowOp);
773 }
774 }
775 else if (posEdge == edge(2, 0))
776 {
777 point p01(planeIntersection(d, tet, 0, 1));
778 point p03(planeIntersection(d, tet, 0, 3));
779 point p12(planeIntersection(d, tet, 1, 2));
780 point p23(planeIntersection(d, tet, 2, 3));
781 // Split the resulting prism
782 {
783 FixedList<point, 6> p
784 (
785 {
786 tet[2],
787 p12,
788 p23,
789 tet[0],
790 p01,
791 p03
792 }
793 );
794
795 //Pout<< " 20 aboveprism:" << p << endl;
796 decomposePrism(p, aboveOp);
797 }
798 {
799 FixedList<point, 6> p
800 (
801 {
802 tet[1],
803 p12,
804 p01,
805 tet[3],
806 p23,
807 p03
808 }
809 );
810
811 //Pout<< " 20 belowprism:" << p << endl;
812 decomposePrism(p, belowOp);
813 }
814 }
815 else if (posEdge == edge(0, 3))
816 {
817 point p01(planeIntersection(d, tet, 0, 1));
818 point p02(planeIntersection(d, tet, 0, 2));
819 point p13(planeIntersection(d, tet, 1, 3));
820 point p23(planeIntersection(d, tet, 2, 3));
821 // Split the resulting prism
822 {
823 FixedList<point, 6> p
824 (
825 {
826 tet[3],
827 p23,
828 p13,
829 tet[0],
830 p02,
831 p01
832 }
833 );
834
835 //Pout<< " 03 aboveprism:" << p << endl;
836 decomposePrism(p, aboveOp);
837 }
838 {
839 FixedList<point, 6> p
840 (
841 {
842 tet[2],
843 p23,
844 p02,
845 tet[1],
846 p13,
847 p01
848 }
849 );
850
851 //Pout<< " 03 belowprism:" << p << endl;
852 decomposePrism(p, belowOp);
853 }
854 }
855 else if (posEdge == edge(1, 3))
856 {
857 point p01(planeIntersection(d, tet, 0, 1));
858 point p12(planeIntersection(d, tet, 1, 2));
859 point p03(planeIntersection(d, tet, 0, 3));
860 point p23(planeIntersection(d, tet, 2, 3));
861 // Split the resulting prism
862 {
863 FixedList<point, 6> p
864 (
865 {
866 tet[1],
867 p12,
868 p01,
869 tet[3],
870 p23,
871 p03
872 }
873 );
874
875 //Pout<< " 13 aboveprism:" << p << endl;
876 decomposePrism(p, aboveOp);
877 }
878 {
879 FixedList<point, 6> p
880 (
881 {
882 tet[2],
883 p12,
884 p23,
885 tet[0],
886 p01,
887 p03
888 }
889 );
890
891 //Pout<< " 13 belowprism:" << p << endl;
892 decomposePrism(p, belowOp);
893 }
894 }
895 else if (posEdge == edge(2, 3))
896 {
897 point p02(planeIntersection(d, tet, 0, 2));
898 point p12(planeIntersection(d, tet, 1, 2));
899 point p03(planeIntersection(d, tet, 0, 3));
900 point p13(planeIntersection(d, tet, 1, 3));
901 // Split the resulting prism
902 {
903 FixedList<point, 6> p
904 (
905 {
906 tet[2],
907 p02,
908 p12,
909 tet[3],
910 p03,
911 p13
912 }
913 );
914
915 //Pout<< " 23 aboveprism:" << p << endl;
916 decomposePrism(p, aboveOp);
917 }
918 {
919 FixedList<point, 6> p
920 (
921 {
922 tet[0],
923 p02,
924 p03,
925 tet[1],
926 p12,
927 p13
928 }
929 );
930
931 //Pout<< " 23 belowprism:" << p << endl;
932 decomposePrism(p, belowOp);
933 }
934 }
935 else
936 {
938 << "Missed edge:" << posEdge
939 << abort(FatalError);
940 }
941 }
942 else if (nPos == 1)
943 {
944 // Find the positive tet
945 label i0 = -1;
946 forAll(d, i)
947 {
948 if (d[i] > 0)
949 {
950 i0 = i;
951 break;
952 }
953 }
954
955 label i1 = d.fcIndex(i0);
956 label i2 = d.fcIndex(i1);
957 label i3 = d.fcIndex(i2);
958
959 point p01(planeIntersection(d, tet, i0, i1));
960 point p02(planeIntersection(d, tet, i0, i2));
961 point p03(planeIntersection(d, tet, i0, i3));
962
963 //Pout<< "Split 1pos tet " << tet << " d:" << d << " into" << nl;
964
965 if (i0 == 0 || i0 == 2)
966 {
967 tetPoints t(tet[i0], p01, p02, p03);
968 //Pout<< " abovetet:" << t << " around i0:" << i0 << endl;
969 //checkTet(t, "nPos 1, aboveTets i0==0 or 2");
970 aboveOp(t);
971
972 // Prism
973 FixedList<point, 6> p
974 (
975 {
976 tet[i1],
977 tet[i3],
978 tet[i2],
979 p01,
980 p03,
981 p02
982 }
983 );
984
985 //Pout<< " belowprism:" << p << endl;
986 decomposePrism(p, belowOp);
987 }
988 else
989 {
990 tetPoints t(p01, p02, p03, tet[i0]);
991 //Pout<< " abovetet:" << t << " around i0:" << i0 << endl;
992 //checkTet(t, "nPos 1, aboveTets i0==1 or 3");
993 aboveOp(t);
994
995 // Prism
996 FixedList<point, 6> p
997 (
998 {
999 tet[i3],
1000 tet[i1],
1001 tet[i2],
1002 p03,
1003 p01,
1004 p02
1005 }
1006 );
1007
1008 //Pout<< " belowprism:" << p << endl;
1009 decomposePrism(p, belowOp);
1010 }
1011 }
1012 else // nPos == 0
1013 {
1014 belowOp(tet);
1015 }
1016}
1017
1018
1019template<class Point, class PointRef>
1020template<class AboveTetOp, class BelowTetOp>
1022(
1023 const plane& pl,
1024 AboveTetOp& aboveOp,
1025 BelowTetOp& belowOp
1026) const
1027{
1028 tetSliceWithPlane(pl, tetPoints(a_, b_, c_, d_), aboveOp, belowOp);
1029}
1030
1031
1032// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
1033
1034template<class Point, class PointRef>
1035inline Foam::Istream& Foam::operator>>
1036(
1037 Istream& is,
1039)
1040{
1041 is.readBegin("tetrahedron");
1042 is >> t.a_ >> t.b_ >> t.c_ >> t.d_;
1043 is.readEnd("tetrahedron");
1044
1045 is.check(FUNCTION_NAME);
1046
1047 return is;
1048}
1049
1050
1051template<class Point, class PointRef>
1052inline Foam::Ostream& Foam::operator<<
1053(
1054 Ostream& os,
1055 const tetrahedron<Point, PointRef>& t
1056)
1057{
1058 os << nl
1060 << t.a_ << token::SPACE
1061 << t.b_ << token::SPACE
1062 << t.c_ << token::SPACE
1063 << t.d_
1064 << token::END_LIST;
1065
1066 return os;
1067}
1068
1069
1070// ************************************************************************* //
CGAL::Point_3< K > Point
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
label n
A 1D vector of objects of type <T> with a fixed length <N>.
Definition: FixedList.H:81
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:64
bool readBegin(const char *funcName)
Begin read of data chunk, starts with '('.
Definition: Istream.C:111
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Describes the interaction of a face and a point. It carries the info of a successful hit and (if succ...
Definition: PointHit.H:54
const point_type & rawPoint() const noexcept
The point, no checks.
Definition: PointHit.H:172
scalar distance() const noexcept
Return distance to hit.
Definition: PointHit.H:139
Random number generator.
Definition: Random.H:60
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: UList.H:94
const Cmpt & z() const
Access to the vector z component.
Definition: VectorI.H:85
const Cmpt & y() const
Access to the vector y component.
Definition: VectorI.H:79
const Cmpt & x() const
Access to the vector x component.
Definition: VectorI.H:73
friend Ostream & operator(Ostream &, const faMatrix< Type > &)
Geometric class that creates a 3D plane and can return the intersection point between a line and the ...
Definition: plane.H:95
Tet storage. Null constructable (unfortunately tetrahedron<point, point> is not)
Definition: tetPoints.H:56
Store resulting tets.
Definition: tetrahedron.H:120
A tetrahedron primitive.
Definition: tetrahedron.H:87
Point circumCentre() const
Return circum-centre.
Definition: tetrahedronI.H:179
barycentric pointToBarycentric(const point &pt) const
Calculate the barycentric coordinates from the given point.
Definition: tetrahedronI.H:266
vector Sd() const
Face area normal for side d.
Definition: tetrahedronI.H:158
triPointRef tri(const label facei) const
Return i-th face.
Definition: tetrahedronI.H:105
scalar circumRadius() const
Return circum-radius.
Definition: tetrahedronI.H:206
const Point & a() const
Return vertices.
Definition: tetrahedronI.H:76
Point centre() const
Return centre (centroid)
Definition: tetrahedronI.H:165
vector Sc() const
Face area normal for side c.
Definition: tetrahedronI.H:151
bool inside(const point &pt) const
Return true if point is inside tetrahedron.
Definition: tetrahedronI.H:415
vector Sb() const
Face area normal for side b.
Definition: tetrahedronI.H:144
const Point & c() const
Definition: tetrahedronI.H:90
scalar mag() const
Return volume.
Definition: tetrahedronI.H:172
pointHit nearestPoint(const point &p) const
Return nearest point to p on tetrahedron. Is p itself.
Definition: tetrahedronI.H:321
Point randomPoint(Random &rndGen) const
Return a random point in the tetrahedron from a.
Definition: tetrahedronI.H:246
void sliceWithPlane(const plane &pl, AboveTetOp &aboveOp, BelowTetOp &belowOp) const
Decompose tet into tets above and below plane.
Point barycentricToPoint(const barycentric &bary) const
Calculate the point from the given barycentric coordinates.
Definition: tetrahedronI.H:256
const Point & b() const
Definition: tetrahedronI.H:83
scalar quality() const
Return quality: Ratio of tetrahedron and circum-sphere.
Definition: tetrahedronI.H:232
const Point & d() const
Definition: tetrahedronI.H:97
vector Sa() const
Face area normal for side a.
Definition: tetrahedronI.H:137
@ BEGIN_LIST
Begin list [isseparator].
Definition: token.H:155
@ END_LIST
End list [isseparator].
Definition: token.H:156
@ SPACE
Space [isspace].
Definition: token.H:125
A triangle primitive used to calculate face normals and swept volumes.
Definition: triangle.H:80
vector areaNormal() const
The area normal - with magnitude equal to area of triangle.
Definition: triangleI.H:112
pointHit nearestPoint(const point &p) const
Return nearest point to p on triangle.
Definition: triangleI.H:671
volScalarField & p
const volScalarField & mu
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
const pointField & points
#define FUNCTION_NAME
dimensionedScalar det(const dimensionedSphericalTensor &dt)
barycentric barycentric01(Random &rndGen)
Generate a random barycentric coordinate within the unit tetrahedron.
Definition: barycentric.C:73
Barycentric< scalar > barycentric
A scalar version of the templated Barycentric.
Definition: barycentric.H:51
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar pow3(const dimensionedScalar &ds)
vector point
Point is a vector.
Definition: point.H:43
Cmpt cmptSum(const SphericalTensor< Cmpt > &st)
Return the sum of components of a SphericalTensor.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
triangle< point, const point & > triPointRef
A triangle using referred points.
Definition: triangle.H:71
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
errorManip< error > abort(error &err)
Definition: errorManip.H:144
PointHit< point > pointHit
A PointIndexHit for 3D points.
Definition: pointHit.H:44
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
volScalarField & b
Definition: createFields.H:27
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Random rndGen
Definition: createFields.H:23