TTK
Loading...
Searching...
No Matches
AbstractTriangulation.h
Go to the documentation of this file.
1
11
12#pragma once
13
14// base code includes
15#include <Cache.h>
16#include <Geometry.h>
17#include <Wrapper.h>
18
19#include <array>
20#include <ostream>
21#include <unordered_map>
22
23#ifdef TTK_ENABLE_KAMIKAZE
24#define TTK_TRIANGULATION_INTERNAL(NAME) NAME
25#else
26#define TTK_TRIANGULATION_INTERNAL(NAME) NAME##Internal
27#endif // TTK_ENABLE_KAMIKAZE
28
29#define ttkTemplateMacroCase(triangulationType, triangulationClass, call) \
30 case triangulationType: { \
31 using TTK_TT = triangulationClass; \
32 call; \
33 }; break
34
35#define ttkTemplateMacro(triangulationType, call) \
36 switch(triangulationType) { \
37 ttkTemplateMacroCase( \
38 ttk::Triangulation::Type::EXPLICIT, ttk::ExplicitTriangulation, call); \
39 ttkTemplateMacroCase( \
40 ttk::Triangulation::Type::IMPLICIT, ttk::ImplicitNoPreconditions, call); \
41 ttkTemplateMacroCase(ttk::Triangulation::Type::HYBRID_IMPLICIT, \
42 ttk::ImplicitWithPreconditions, call); \
43 ttkTemplateMacroCase( \
44 ttk::Triangulation::Type::PERIODIC, ttk::PeriodicNoPreconditions, call); \
45 ttkTemplateMacroCase(ttk::Triangulation::Type::HYBRID_PERIODIC, \
46 ttk::PeriodicWithPreconditions, call); \
47 ttkTemplateMacroCase( \
48 ttk::Triangulation::Type::COMPACT, ttk::CompactTriangulation, call); \
49 }
50
51namespace ttk {
52
53 // forward declaration of the ttk::dcg::DiscreteGradient class to
54 // give it access to the gradient cache (with a `friend`
55 // declaration)
56 namespace dcg {
57 class DiscreteGradient;
58 }
59
61
62 public:
64
66
71
74 int clear();
75
78 size_t footprint(size_t size = 0) const;
79
101 virtual inline int getCellEdge(const SimplexId &cellId,
102 const int &localEdgeId,
103 SimplexId &edgeId) const {
104
105#ifndef TTK_ENABLE_KAMIKAZE
106 // initialize output variable before early return
107 edgeId = -1;
108
110 return -1;
111#endif
112
113 if(getDimensionality() == 1)
114 return getCellNeighbor(cellId, localEdgeId, edgeId);
115
116 else if(getDimensionality() == 2)
117 return getTriangleEdgeInternal(cellId, localEdgeId, edgeId);
118
119 return getCellEdgeInternal(cellId, localEdgeId, edgeId);
120 }
121
138 virtual inline SimplexId getCellEdgeNumber(const SimplexId &cellId) const {
139#ifndef TTK_ENABLE_KAMIKAZE
141 return -1;
142#endif
143 if(getDimensionality() == 1)
144 return getCellNeighborNumber(cellId);
145
146 else if(getDimensionality() == 2)
147 return getTriangleEdgeNumber(cellId);
148
149 return getCellEdgeNumberInternal(cellId);
150 }
151
180 virtual inline const std::vector<std::vector<SimplexId>> *getCellEdges() {
181#ifndef TTK_ENABLE_KAMIKAZE
183 return nullptr;
184#endif
185 if(getDimensionality() == 1)
186 return getCellNeighbors();
187
188 else if(getDimensionality() == 2)
190
191 return getCellEdgesInternal();
192 }
193
213 virtual inline int getCellNeighbor(const SimplexId &cellId,
214 const int &localNeighborId,
215 SimplexId &neighborId) const {
216
217#ifndef TTK_ENABLE_KAMIKAZE
218 // initialize output variable before early return
219 neighborId = -1;
220
222 return -1;
223#endif
224 return getCellNeighborInternal(cellId, localNeighborId, neighborId);
225 }
226
241 virtual inline SimplexId
242 getCellNeighborNumber(const SimplexId &cellId) const {
243
244#ifndef TTK_ENABLE_KAMIKAZE
246 return -1;
247#endif
248 return getCellNeighborNumberInternal(cellId);
249 }
250
276 virtual inline const std::vector<std::vector<SimplexId>> *
278#ifndef TTK_ENABLE_KAMIKAZE
280 return nullptr;
281#endif
283 }
284
307 virtual inline int getCellTriangle(const SimplexId &cellId,
308 const int &localTriangleId,
309 SimplexId &triangleId) const {
310#ifndef TTK_ENABLE_KAMIKAZE
311 // initialize output variable before early return
312 triangleId = -1;
313
314 if(getDimensionality() == 1)
315 return -1;
316
318 return -2;
319#endif
320 if(getDimensionality() == 2)
321 return getCellNeighbor(cellId, localTriangleId, triangleId);
322
323 return getCellTriangleInternal(cellId, localTriangleId, triangleId);
324 }
325
343 virtual inline SimplexId
344 getCellTriangleNumber(const SimplexId &cellId) const {
345#ifndef TTK_ENABLE_KAMIKAZE
346 if(getDimensionality() == 1)
347 return -1;
348
350 return -2;
351#endif
352 if(getDimensionality() == 2)
353 return getCellNeighborNumber(cellId);
354
355 return getCellTriangleNumberInternal(cellId);
356 }
357
389 virtual inline const std::vector<std::vector<SimplexId>> *
391#ifndef TTK_ENABLE_KAMIKAZE
392 if(getDimensionality() == 1)
393 return nullptr;
394
396 return nullptr;
397#endif
398 if(getDimensionality() == 2)
399 return getCellNeighbors();
400
402 }
403
415 virtual inline int getCellVertex(const SimplexId &cellId,
416 const int &localVertexId,
417 SimplexId &vertexId) const {
418
419 return getCellVertexInternal(cellId, localVertexId, vertexId);
420 }
421
428 virtual inline SimplexId
429 getCellVertexNumber(const SimplexId &cellId) const {
430 return getCellVertexNumberInternal(cellId);
431 }
432
436 virtual inline int getDimensionality() const {
438 }
439
464 virtual inline const std::vector<std::array<SimplexId, 2>> *getEdges() {
465#ifndef TTK_ENABLE_KAMIKAZE
466 if(getDimensionality() == 1)
467 return nullptr;
468
470 return nullptr;
471#endif
472 return getEdgesInternal();
473 }
474
494 virtual inline int getEdgeLink(const SimplexId &edgeId,
495 const int &localLinkId,
496 SimplexId &linkId) const {
497#ifndef TTK_ENABLE_KAMIKAZE
498 // initialize output variable before early return
499 linkId = -1;
500
501 if(getDimensionality() == 1)
502 return -1;
503
505 return -2;
506#endif
507 return getEdgeLinkInternal(edgeId, localLinkId, linkId);
508 }
509
524 virtual inline SimplexId getEdgeLinkNumber(const SimplexId &edgeId) const {
525#ifndef TTK_ENABLE_KAMIKAZE
526 if(getDimensionality() == 1)
527 return -1;
528
530 return -2;
531#endif
532 return getEdgeLinkNumberInternal(edgeId);
533 }
534
558 virtual inline const std::vector<std::vector<SimplexId>> *getEdgeLinks() {
559#ifndef TTK_ENABLE_KAMIKAZE
560 if(getDimensionality() == 1)
561 return nullptr;
562
564 return nullptr;
565#endif
566 return getEdgeLinksInternal();
567 }
568
591 virtual inline int getEdgeStar(const SimplexId &edgeId,
592 const int &localStarId,
593 SimplexId &starId) const {
594
595#ifndef TTK_ENABLE_KAMIKAZE
596 // initialize output variable before early return
597 starId = -1;
598
599 if(getDimensionality() == 1)
600 return -1;
601
603 return -2;
604#endif
605 return getEdgeStarInternal(edgeId, localStarId, starId);
606 }
607
625 virtual inline SimplexId getEdgeStarNumber(const SimplexId &edgeId) const {
626#ifndef TTK_ENABLE_KAMIKAZE
627 if(getDimensionality() == 1)
628 return -1;
629
631 return -2;
632#endif
633 return getEdgeStarNumberInternal(edgeId);
634 }
635
664 virtual inline const std::vector<std::vector<SimplexId>> *getEdgeStars() {
665#ifndef TTK_ENABLE_KAMIKAZE
666 if(getDimensionality() == 1)
667 return nullptr;
668
670 return nullptr;
671#endif
672 return getEdgeStarsInternal();
673 }
674
694 virtual inline int getEdgeTriangle(const SimplexId &edgeId,
695 const int &localTriangleId,
696 SimplexId &triangleId) const {
697#ifndef TTK_ENABLE_KAMIKAZE
698 // initialize output variable before early return
699 triangleId = -1;
700
701 if(getDimensionality() == 1)
702 return -1;
703
705 return -2;
706#endif
707 if(getDimensionality() == 2)
708 return getEdgeStar(edgeId, localTriangleId, triangleId);
709
710 return getEdgeTriangleInternal(edgeId, localTriangleId, triangleId);
711 }
712
727 virtual inline SimplexId
728 getEdgeTriangleNumber(const SimplexId &edgeId) const {
729#ifndef TTK_ENABLE_KAMIKAZE
730 if(getDimensionality() == 1)
731 return -1;
732
734 return -2;
735#endif
736
737 if(getDimensionality() == 2)
738 return getEdgeStarNumber(edgeId);
739
740 return getEdgeTriangleNumberInternal(edgeId);
741 }
742
768 virtual inline const std::vector<std::vector<SimplexId>> *
770#ifndef TTK_ENABLE_KAMIKAZE
771 if(getDimensionality() == 1)
772 return nullptr;
773
775 return nullptr;
776#endif
777
778 if(getDimensionality() == 2)
779 return getEdgeStars();
780
782 }
783
801 virtual inline int getEdgeVertex(const SimplexId &edgeId,
802 const int &localVertexId,
803 SimplexId &vertexId) const {
804#ifndef TTK_ENABLE_KAMIKAZE
805 // initialize output variable before early return
806 vertexId = -1;
807
809 return -2;
810#endif
811 if(getDimensionality() == 1)
812 return getCellVertex(edgeId, localVertexId, vertexId);
813
814 return getEdgeVertexInternal(edgeId, localVertexId, vertexId);
815 }
816
819 inline int getEdgeVertexNumber(const SimplexId ttkNotUsed(edgeId)) const {
820 return 2;
821 }
822
830 virtual inline const std::array<SimplexId, 3> &getGridDimensions() const {
831 return this->gridDimensions_;
832 }
833
839 virtual inline SimplexId getNumberOfCells() const {
841 }
842
856 virtual inline SimplexId getNumberOfEdges() const {
857#ifndef TTK_ENABLE_KAMIKAZE
859 return -2;
860#endif
861 if(getDimensionality() == 1)
862 return getNumberOfCells();
863
865 }
866
880 virtual inline SimplexId getNumberOfTriangles() const {
881#ifndef TTK_ENABLE_KAMIKAZE
882 if(getDimensionality() == 1)
883 return -1;
884
886 return -2;
887#endif
888 if(getDimensionality() == 2)
889 return getNumberOfCells();
890
892 }
893
896 virtual inline SimplexId getNumberOfVertices() const {
898 }
899
924 virtual inline const std::vector<std::array<SimplexId, 3>> *getTriangles() {
925#ifndef TTK_ENABLE_KAMIKAZE
927 return nullptr;
928#endif
929 return getTrianglesInternal();
930 }
931
950 virtual inline int getTriangleEdge(const SimplexId &triangleId,
951 const int &localEdgeId,
952 SimplexId &edgeId) const {
953#ifndef TTK_ENABLE_KAMIKAZE
954 // initialize output variable before early return
955 edgeId = -1;
956
957 if(getDimensionality() == 1)
958 return -1;
959
961 return -2;
962#endif
963
964 return getTriangleEdgeInternal(triangleId, localEdgeId, edgeId);
965 }
966
981 virtual inline SimplexId
982 getTriangleEdgeNumber(const SimplexId &triangleId) const {
983#ifndef TTK_ENABLE_KAMIKAZE
984 if(getDimensionality() == 1)
985 return -1;
986
988 return -2;
989#endif
990
991 return getTriangleEdgeNumberInternal(triangleId);
992 }
993
1019 virtual inline const std::vector<std::vector<SimplexId>> *
1021#ifndef TTK_ENABLE_KAMIKAZE
1022 if(getDimensionality() == 1)
1023 return nullptr;
1024
1026 return nullptr;
1027#endif
1028
1029 return getTriangleEdgesInternal();
1030 }
1031
1051 virtual inline int getTriangleLink(const SimplexId &triangleId,
1052 const int &localLinkId,
1053 SimplexId &linkId) const {
1054#ifndef TTK_ENABLE_KAMIKAZE
1055 // initialize output variable before early return
1056 linkId = -1;
1057
1058 if(getDimensionality() != 3)
1059 return -1;
1060
1062 return -2;
1063#endif
1064 return getTriangleLinkInternal(triangleId, localLinkId, linkId);
1065 }
1066
1082 virtual inline SimplexId
1083 getTriangleLinkNumber(const SimplexId &triangleId) const {
1084#ifndef TTK_ENABLE_KAMIKAZE
1085 if(getDimensionality() != 3)
1086 return -1;
1087
1089 return -2;
1090#endif
1091 return getTriangleLinkNumberInternal(triangleId);
1092 }
1093
1118 virtual inline const std::vector<std::vector<SimplexId>> *
1120#ifndef TTK_ENABLE_KAMIKAZE
1121 if(getDimensionality() != 3)
1122 return nullptr;
1123
1125 return nullptr;
1126#endif
1127 return getTriangleLinksInternal();
1128 }
1129
1150 virtual inline int getTriangleStar(const SimplexId &triangleId,
1151 const int &localStarId,
1152 SimplexId &starId) const {
1153#ifndef TTK_ENABLE_KAMIKAZE
1154 // initialize output variable before early return
1155 starId = -1;
1156
1157 if(getDimensionality() != 3)
1158 return -1;
1159
1161 return -2;
1162#endif
1163 return getTriangleStarInternal(triangleId, localStarId, starId);
1164 }
1165
1180 virtual inline SimplexId
1181 getTriangleStarNumber(const SimplexId &triangleId) const {
1182#ifndef TTK_ENABLE_KAMIKAZE
1183 if(getDimensionality() != 3)
1184 return -1;
1185
1187 return -2;
1188#endif
1189 return getTriangleStarNumberInternal(triangleId);
1190 }
1191
1217 virtual inline const std::vector<std::vector<SimplexId>> *
1219#ifndef TTK_ENABLE_KAMIKAZE
1220 if(getDimensionality() != 3)
1221 return nullptr;
1222
1224 return nullptr;
1225#endif
1226 return getTriangleStarsInternal();
1227 }
1228
1246 virtual inline int getTriangleVertex(const SimplexId &triangleId,
1247 const int &localVertexId,
1248 SimplexId &vertexId) const {
1249#ifndef TTK_ENABLE_KAMIKAZE
1250 // initialize output variable before early return
1251 vertexId = -1;
1252
1253 if(getDimensionality() == 1)
1254 return -1;
1255
1257 return -2;
1258#endif
1259 if(getDimensionality() == 2)
1260 return getCellVertex(triangleId, localVertexId, vertexId);
1261
1262 return getTriangleVertexInternal(triangleId, localVertexId, vertexId);
1263 }
1264
1267 inline int
1269 return 3;
1270 }
1271
1292 virtual inline int getVertexEdge(const SimplexId &vertexId,
1293 const int &localEdgeId,
1294 SimplexId &edgeId) const {
1295
1296#ifndef TTK_ENABLE_KAMIKAZE
1297 // initialize output variable before early return
1298 edgeId = -1;
1299
1301 return -1;
1302#endif
1303 if(getDimensionality() == 1)
1304 return getVertexStar(vertexId, localEdgeId, edgeId);
1305
1306 return getVertexEdgeInternal(vertexId, localEdgeId, edgeId);
1307 }
1308
1323 virtual inline SimplexId
1324 getVertexEdgeNumber(const SimplexId &vertexId) const {
1325#ifndef TTK_ENABLE_KAMIKAZE
1327 return -1;
1328#endif
1329 if(getDimensionality() == 1)
1330 return getVertexStarNumber(vertexId);
1331
1332 return getVertexEdgeNumberInternal(vertexId);
1333 }
1334
1361 virtual inline const std::vector<std::vector<SimplexId>> *getVertexEdges() {
1362#ifndef TTK_ENABLE_KAMIKAZE
1364 return nullptr;
1365#endif
1366 if(getDimensionality() == 1)
1367 return getVertexStars();
1368
1369 return getVertexEdgesInternal();
1370 }
1371
1391 virtual inline int getVertexLink(const SimplexId &vertexId,
1392 const int &localLinkId,
1393 SimplexId &linkId) const {
1394#ifndef TTK_ENABLE_KAMIKAZE
1395 // initialize output variable before early return
1396 linkId = -1;
1397
1399 return -1;
1400#endif
1401 return getVertexLinkInternal(vertexId, localLinkId, linkId);
1402 }
1403
1418 virtual inline SimplexId
1419 getVertexLinkNumber(const SimplexId &vertexId) const {
1420#ifndef TTK_ENABLE_KAMIKAZE
1422 return -1;
1423#endif
1424 return getVertexLinkNumberInternal(vertexId);
1425 }
1426
1451 virtual inline const std::vector<std::vector<SimplexId>> *getVertexLinks() {
1452#ifndef TTK_ENABLE_KAMIKAZE
1454 return nullptr;
1455#endif
1456 return getVertexLinksInternal();
1457 }
1458
1475 virtual inline int getVertexNeighbor(const SimplexId &vertexId,
1476 const int &localNeighborId,
1477 SimplexId &neighborId) const {
1478#ifndef TTK_ENABLE_KAMIKAZE
1479 // initialize output variable before early return
1480 neighborId = -1;
1481
1483 return -1;
1484#endif
1485 return getVertexNeighborInternal(vertexId, localNeighborId, neighborId);
1486 }
1487
1499 virtual inline SimplexId
1500 getVertexNeighborNumber(const SimplexId &vertexId) const {
1501#ifndef TTK_ENABLE_KAMIKAZE
1503 return -1;
1504#endif
1505 return getVertexNeighborNumberInternal(vertexId);
1506 }
1507
1531 virtual inline const std::vector<std::vector<SimplexId>> *
1533#ifndef TTK_ENABLE_KAMIKAZE
1535 return nullptr;
1536#endif
1538 }
1539
1546 virtual inline int getVertexPoint(const SimplexId &vertexId,
1547 float &x,
1548 float &y,
1549 float &z) const {
1550 return getVertexPointInternal(vertexId, x, y, z);
1551 }
1552
1573 virtual inline int getVertexStar(const SimplexId &vertexId,
1574 const int &localStarId,
1575 SimplexId &starId) const {
1576#ifndef TTK_ENABLE_KAMIKAZE
1577 // initialize output variable before early return
1578 starId = -1;
1579
1581 return -1;
1582#endif
1583 return getVertexStarInternal(vertexId, localStarId, starId);
1584 }
1585
1600 virtual inline SimplexId
1601 getVertexStarNumber(const SimplexId &vertexId) const {
1602#ifndef TTK_ENABLE_KAMIKAZE
1604 return -1;
1605#endif
1606 return getVertexStarNumberInternal(vertexId);
1607 }
1608
1634 virtual inline const std::vector<std::vector<SimplexId>> *getVertexStars() {
1635#ifndef TTK_ENABLE_KAMIKAZE
1637 return nullptr;
1638#endif
1639 return getVertexStarsInternal();
1640 }
1641
1662 virtual inline int getVertexTriangle(const SimplexId &vertexId,
1663 const int &localTriangleId,
1664 SimplexId &triangleId) const {
1665#ifndef TTK_ENABLE_KAMIKAZE
1666 // initialize output variable before early return
1667 triangleId = -1;
1668
1669 if(getDimensionality() == 1)
1670 return -1;
1671
1673 return -2;
1674#endif
1675 if(getDimensionality() == 2)
1676 return getVertexStar(vertexId, localTriangleId, triangleId);
1677
1678 return getVertexTriangleInternal(vertexId, localTriangleId, triangleId);
1679 }
1680
1695 virtual inline SimplexId
1696 getVertexTriangleNumber(const SimplexId &vertexId) const {
1697#ifndef TTK_ENABLE_KAMIKAZE
1698 if(getDimensionality() == 1)
1699 return -1;
1700
1702 return -2;
1703#endif
1704
1705 if(getDimensionality() == 2)
1706 return getVertexStarNumber(vertexId);
1707
1708 return getVertexTriangleNumberInternal(vertexId);
1709 }
1710
1736 virtual inline const std::vector<std::vector<SimplexId>> *
1738#ifndef TTK_ENABLE_KAMIKAZE
1739 if(getDimensionality() == 1)
1740 return nullptr;
1741
1743 return nullptr;
1744#endif
1745 if(getDimensionality() == 2)
1746 return getVertexStars();
1747
1749 }
1750
1752 inline bool hasPeriodicBoundaries() const {
1754 }
1755
1775 virtual inline bool isEdgeOnBoundary(const SimplexId &edgeId) const {
1776#ifndef TTK_ENABLE_KAMIKAZE
1778 return false;
1779#endif
1780 return isEdgeOnBoundaryInternal(edgeId);
1781 }
1782
1785 virtual inline bool isEmpty() const {
1786 return true;
1787 }
1788
1792 virtual inline bool isManifold() const {
1793#ifndef TTK_ENABLE_KAMIKAZE
1794 if(!this->hasPreconditionedManifold())
1795 return false;
1796#endif
1797 return this->isManifold_;
1798 }
1799
1820 virtual inline int preconditionManifold() {
1821
1822 if(!this->hasPreconditionedManifold_) {
1824 this->hasPreconditionedManifold_ = true;
1825 }
1826 return 0;
1827 }
1828
1849 virtual inline bool
1850 isTriangleOnBoundary(const SimplexId &triangleId) const {
1851#ifndef TTK_ENABLE_KAMIKAZE
1853 return false;
1854#endif
1855 return isTriangleOnBoundaryInternal(triangleId);
1856 }
1857
1875 virtual inline bool isVertexOnBoundary(const SimplexId &vertexId) const {
1876#ifndef TTK_ENABLE_KAMIKAZE
1878 return false;
1879#endif
1880 return isVertexOnBoundaryInternal(vertexId);
1881 }
1882
1896 virtual inline int preconditionBoundaryEdges() {
1897
1902 }
1903 return 0;
1904 }
1905
1919 virtual inline int preconditionBoundaryTriangles() {
1920
1925 }
1926 return 0;
1927 }
1928
1942 virtual inline int preconditionBoundaryVertices() {
1943
1947 }
1948 return 0;
1949 }
1950
1966 virtual inline int preconditionCellEdges() {
1967
1970 if(getDimensionality() == 1)
1972
1975 }
1976 return 0;
1977 }
1978
1996 virtual inline int preconditionCellNeighbors() {
1997
2001 }
2002 return 0;
2003 }
2004
2022 virtual inline int preconditionCellTriangles() {
2023
2026
2027#ifndef TTK_ENABLE_KAMIKAZE
2028 if(getDimensionality() == 1)
2029 return -1;
2030#endif
2031 if(getDimensionality() == 2)
2033
2036 }
2037 return 0;
2038 }
2039
2057 virtual inline int preconditionEdges() {
2058
2062 }
2063 return 0;
2064 }
2065
2083 virtual inline int preconditionEdgeLinks() {
2084
2086#ifndef TTK_ENABLE_KAMIKAZE
2087 if(getDimensionality() == 1)
2088 return -1;
2089#endif
2093 }
2094 return 0;
2095 }
2096
2114 virtual inline int preconditionEdgeStars() {
2115
2117#ifndef TTK_ENABLE_KAMIKAZE
2118 if(getDimensionality() == 1)
2119 return -1;
2120#endif
2124 }
2125 return 0;
2126 }
2127
2145 virtual inline int preconditionEdgeTriangles() {
2146
2149
2150#ifndef TTK_ENABLE_KAMIKAZE
2151 if(getDimensionality() == 1)
2152 return -1;
2153#endif
2154
2155 if(getDimensionality() == 2)
2156 return preconditionEdgeStars();
2157
2161 }
2162
2163 return 0;
2164 }
2165
2183 virtual inline int preconditionTriangles() {
2184
2187
2188#ifndef TTK_ENABLE_KAMIKAZE
2189 if(getDimensionality() == 1)
2190 return -1;
2191#endif
2192 if(getDimensionality() == 2)
2193 return 0;
2194
2196 }
2197
2198 return 0;
2199 }
2200
2218 virtual inline int preconditionTriangleEdges() {
2219
2222
2223#ifndef TTK_ENABLE_KAMIKAZE
2224 if(getDimensionality() == 1)
2225 return -1;
2226#endif
2227 if(getDimensionality() == 2)
2228 return preconditionCellEdges();
2229
2233 }
2234 return 0;
2235 }
2236
2254 virtual inline int preconditionTriangleLinks() {
2255
2257#ifndef TTK_ENABLE_KAMIKAZE
2258 if(getDimensionality() != 3)
2259 return -2;
2260#endif
2264 }
2265 return 0;
2266 }
2267
2285 virtual inline int preconditionTriangleStars() {
2286
2288#ifndef TTK_ENABLE_KAMIKAZE
2289 if(getDimensionality() != 3)
2290 return -1;
2291#endif
2295 }
2296 return 0;
2297 }
2298
2316 virtual inline int preconditionVertexEdges() {
2317
2320
2321 if(getDimensionality() == 1)
2322 return preconditionVertexStars();
2323
2326 }
2327 return 0;
2328 }
2329
2347 virtual inline int preconditionVertexLinks() {
2348
2352 }
2353 return 0;
2354 }
2355
2373 virtual inline int preconditionVertexNeighbors() {
2374
2378 }
2379 return 0;
2380 }
2381
2399 virtual inline int preconditionVertexStars() {
2400
2404 }
2405 return 0;
2406 }
2407
2425 virtual inline int preconditionVertexTriangles() {
2426
2429
2430#ifndef TTK_ENABLE_KAMIKAZE
2431 if(getDimensionality() == 1)
2432 return -1;
2433#endif
2434
2435 if(getDimensionality() == 2)
2436 return preconditionVertexStars();
2437
2440 }
2441 return 0;
2442 }
2443
2447 inline int getEdgeIncenter(const SimplexId edgeId,
2448 float incenter[3]) const {
2449 std::array<SimplexId, 2> vertexId{};
2450 getEdgeVertex(edgeId, 0, vertexId[0]);
2451 getEdgeVertex(edgeId, 1, vertexId[1]);
2452
2453 std::array<float, 6> p{};
2454 getVertexPoint(vertexId[0], p[0], p[1], p[2]);
2455 getVertexPoint(vertexId[1], p[3], p[4], p[5]);
2456
2457 incenter[0] = 0.5 * p[0] + 0.5 * p[3];
2458 incenter[1] = 0.5 * p[1] + 0.5 * p[4];
2459 incenter[2] = 0.5 * p[2] + 0.5 * p[5];
2460
2461 return 0;
2462 }
2463
2467 inline int getTriangleIncenter(const SimplexId triangleId,
2468 float incenter[3]) const {
2469
2470 std::array<SimplexId, 3> vertexId{};
2471 if(getDimensionality() == 2) {
2472 getCellVertex(triangleId, 0, vertexId[0]);
2473 getCellVertex(triangleId, 1, vertexId[1]);
2474 getCellVertex(triangleId, 2, vertexId[2]);
2475 } else if(getDimensionality() == 3) {
2476 getTriangleVertex(triangleId, 0, vertexId[0]);
2477 getTriangleVertex(triangleId, 1, vertexId[1]);
2478 getTriangleVertex(triangleId, 2, vertexId[2]);
2479 }
2480
2481 std::array<float, 9> p{};
2482 getVertexPoint(vertexId[0], p[0], p[1], p[2]);
2483 getVertexPoint(vertexId[1], p[3], p[4], p[5]);
2484 getVertexPoint(vertexId[2], p[6], p[7], p[8]);
2485
2486 std::array<float, 3> d{};
2487 d[0] = Geometry::distance(&p[3], &p[6]);
2488 d[1] = Geometry::distance(&p[0], &p[6]);
2489 d[2] = Geometry::distance(&p[0], &p[3]);
2490 const float sum = d[0] + d[1] + d[2];
2491
2492 d[0] = d[0] / sum;
2493 d[1] = d[1] / sum;
2494 d[2] = d[2] / sum;
2495
2496 incenter[0] = d[0] * p[0] + d[1] * p[3] + d[2] * p[6];
2497 incenter[1] = d[0] * p[1] + d[1] * p[4] + d[2] * p[7];
2498 incenter[2] = d[0] * p[2] + d[1] * p[5] + d[2] * p[8];
2499
2500 return 0;
2501 }
2502
2507 inline int getTetraIncenter(const SimplexId tetraId,
2508 float incenter[3]) const {
2509 incenter[0] = 0.0f;
2510 incenter[1] = 0.0f;
2511 incenter[2] = 0.0f;
2512
2513 std::array<float, 3> p{};
2514 for(int i = 0; i < 4; ++i) {
2515 SimplexId triangleId;
2516 getCellTriangle(tetraId, i, triangleId);
2517 getTriangleIncenter(triangleId, p.data());
2518 incenter[0] += p[0];
2519 incenter[1] += p[1];
2520 incenter[2] += p[2];
2521 }
2522
2523 incenter[0] /= 4.0f;
2524 incenter[1] /= 4.0f;
2525 incenter[2] /= 4.0f;
2526
2527 return 0;
2528 }
2529
2533 inline int getCellIncenter(const SimplexId cellid,
2534 const int dim,
2535 float incenter[3]) const {
2536 switch(dim) {
2537 case 0:
2538 getVertexPoint(cellid, incenter[0], incenter[1], incenter[2]);
2539 break;
2540 case 1:
2541 getEdgeIncenter(cellid, incenter);
2542 break;
2543 case 2:
2544 getTriangleIncenter(cellid, incenter);
2545 break;
2546 case 3:
2547 getTetraIncenter(cellid, incenter);
2548 break;
2549 }
2550 return 0;
2551 }
2552
2553 virtual inline int getCellVTKID(const int &ttkId, int &vtkId) const {
2554
2555#ifndef TTK_ENABLE_KAMIKAZE
2556 // initialize output variable before early return
2557 vtkId = -1;
2558#endif
2559 return getCellVTKIDInternal(ttkId, vtkId);
2560 }
2561
2562#ifdef TTK_ENABLE_MPI
2563
2564 // "vtkGhostType" on points & cells
2565
2566 inline void setVertexGhostArray(const unsigned char *const data) {
2567 this->vertexGhost_ = data;
2568 }
2569
2570 inline void setCellGhostArray(const unsigned char *const data) {
2571 this->cellGhost_ = data;
2572 }
2573
2584 virtual int preconditionGlobalBoundary() {
2585
2586 if(!hasPreconditionedGlobalBoundary_) {
2587 preconditionGlobalBoundaryInternal();
2588 hasPreconditionedGlobalBoundary_ = true;
2589 }
2590 return 0;
2591 }
2592
2593 virtual inline int preconditionGlobalBoundaryInternal() {
2594 return 0;
2595 }
2596
2597 // public preconditions methods (used in Triangulation/ttkAlgorithm)
2598 virtual int preconditionDistributedVertices() {
2599 return 0;
2600 }
2601 virtual int preconditionDistributedCells() {
2602 return 0;
2603 }
2604
2605 // This method should be called to initialize and populate the
2606 // ghostCellsPerOwner_ and the remoteGhostCells_ attributes.
2607
2608 virtual int preconditionExchangeGhostCells() {
2609 return 0;
2610 }
2611
2612 // This method should be called to initialize and populate the
2613 // ghostVerticesPerOwner_ and the remoteGhostVertices_ attributes.
2614
2615 virtual int preconditionExchangeGhostVertices() {
2616 return 0;
2617 }
2618
2619 virtual int setVertexRankArray(const int *ttkNotUsed(rankArray)) {
2620 return 0;
2621 }
2622
2623 virtual int setCellRankArray(const int *ttkNotUsed(rankArray)) {
2624 return 0;
2625 }
2626
2627 virtual int preconditionEdgeRankArray() {
2628 return 0;
2629 }
2630
2631 virtual int preconditionTriangleRankArray() {
2632 return 0;
2633 }
2634
2635 // global <-> local id mappings
2636
2637 virtual inline SimplexId getVertexGlobalId(const SimplexId lvid) const {
2638 if(!ttk::isRunningWithMPI()) {
2639 return lvid;
2640 }
2641#ifndef TTK_ENABLE_KAMIKAZE
2642 if(this->getDimensionality() != 1 && this->getDimensionality() != 2
2643 && this->getDimensionality() != 3) {
2644 this->printErr("Only 1D, 2D and 3D datasets are supported");
2645 return -1;
2646 }
2647 if(!this->hasPreconditionedDistributedVertices_) {
2648 this->printErr("VertexGlobalId query without pre-process!");
2649 this->printErr(
2650 "Please call preconditionDistributedVertices() in a pre-process.");
2651 return -1;
2652 }
2653 if(lvid < 0 || lvid >= this->getNumberOfVertices()) {
2654 return -1;
2655 }
2656#endif // TTK_ENABLE_KAMIKAZE
2657 return this->getVertexGlobalIdInternal(lvid);
2658 }
2659 virtual inline SimplexId getVertexLocalId(const SimplexId gvid) const {
2660
2661 if(!ttk::isRunningWithMPI()) {
2662 return gvid;
2663 }
2664#ifndef TTK_ENABLE_KAMIKAZE
2665 if(this->getDimensionality() != 1 && this->getDimensionality() != 2
2666 && this->getDimensionality() != 3) {
2667 this->printErr("Only 1D, 2D and 3D datasets are supported");
2668 return -1;
2669 }
2670 if(!this->hasPreconditionedDistributedVertices_) {
2671 this->printErr("VertexLocalId query without pre-process!");
2672 this->printErr(
2673 "Please call preconditionDistributedVertices() in a pre-process.");
2674 return -1;
2675 }
2676#endif // TTK_ENABLE_KAMIKAZE
2677 return this->getVertexLocalIdInternal(gvid);
2678 }
2679
2680 virtual inline SimplexId getCellGlobalId(const SimplexId lcid) const {
2681#ifndef TTK_ENABLE_KAMIKAZE
2682 if(this->getDimensionality() != 1 && this->getDimensionality() != 2
2683 && this->getDimensionality() != 3) {
2684 this->printErr("Only 1D, 2D and 3D datasets are supported");
2685 return -1;
2686 }
2687 if(!this->hasPreconditionedDistributedCells_) {
2688 this->printErr("CellGlobalId query without pre-process!");
2689 this->printErr(
2690 "Please call preconditionDistributedCells() in a pre-process.");
2691 return -1;
2692 }
2693 if(lcid < 0 || lcid >= this->getNumberOfCells()) {
2694 return -1;
2695 }
2696#endif // TTK_ENABLE_KAMIKAZE
2697 if(!ttk::isRunningWithMPI()) {
2698 return lcid;
2699 }
2700 return this->getCellGlobalIdInternal(lcid);
2701 }
2702 virtual inline SimplexId getCellLocalId(const SimplexId gcid) const {
2703#ifndef TTK_ENABLE_KAMIKAZE
2704 if(this->getDimensionality() != 1 && this->getDimensionality() != 2
2705 && this->getDimensionality() != 3) {
2706 this->printErr("Only 1D, 2D and 3D datasets are supported");
2707 return -1;
2708 }
2709 if(!this->hasPreconditionedDistributedCells_) {
2710 this->printErr("CellLocalId query without pre-process!");
2711 this->printErr(
2712 "Please call preconditionDistributedCells() in a pre-process.");
2713 return -1;
2714 }
2715#endif // TTK_ENABLE_KAMIKAZE
2716 if(!ttk::isRunningWithMPI()) {
2717 return gcid;
2718 }
2719 return this->getCellLocalIdInternal(gcid);
2720 }
2721
2722 virtual inline SimplexId getEdgeGlobalId(const SimplexId leid) const {
2723 const auto dim{this->getDimensionality()};
2724#ifndef TTK_ENABLE_KAMIKAZE
2725 if(dim != 1 && dim != 2 && dim != 3) {
2726 this->printErr("Only 1D, 2D and 3D datasets are supported");
2727 return -1;
2728 }
2729 if(!this->hasPreconditionedDistributedEdges_) {
2730 this->printErr("EdgeGlobalId query without pre-process!");
2731 this->printErr(
2732 "Please call preconditionDistributedEdges() in a pre-process.");
2733 return -1;
2734 }
2735 if(leid < 0 || leid >= this->getNumberOfEdges()) {
2736 return -1;
2737 }
2738#endif // TTK_ENABLE_KAMIKAZE
2739 if(!ttk::isRunningWithMPI()) {
2740 return leid;
2741 }
2742 if(dim == 2 || dim == 3) {
2743 return this->getEdgeGlobalIdInternal(leid);
2744 } else if(dim == 1) {
2745 return this->getCellGlobalIdInternal(leid);
2746 }
2747 return -1;
2748 }
2749 virtual inline SimplexId getEdgeLocalId(const SimplexId geid) const {
2750 const auto dim{this->getDimensionality()};
2751#ifndef TTK_ENABLE_KAMIKAZE
2752 if(dim != 1 && dim != 2 && dim != 3) {
2753 this->printErr("Only 1D, 2D and 3D datasets are supported");
2754 return -1;
2755 }
2756 if(!this->hasPreconditionedDistributedEdges_) {
2757 this->printErr("EdgeLocalId query without pre-process!");
2758 this->printErr(
2759 "Please call preconditionDistributedEdges() in a pre-process.");
2760 return -1;
2761 }
2762#endif // TTK_ENABLE_KAMIKAZE
2763 if(!ttk::isRunningWithMPI()) {
2764 return geid;
2765 }
2766 if(dim == 2 || dim == 3) {
2767 return this->getEdgeLocalIdInternal(geid);
2768 } else if(dim == 1) {
2769 return this->getCellLocalIdInternal(geid);
2770 }
2771 return -1;
2772 }
2773
2774 virtual inline SimplexId getTriangleGlobalId(const SimplexId ltid) const {
2775 const auto dim{this->getDimensionality()};
2776#ifndef TTK_ENABLE_KAMIKAZE
2777 if(dim != 3 && dim != 2) {
2778 this->printErr("Only 2D and 3D datasets are supported");
2779 return -1;
2780 }
2781 if(!this->hasPreconditionedDistributedEdges_) {
2782 this->printErr("TriangleGlobalId query without pre-process!");
2783 this->printErr(
2784 "Please call preconditionDistributedTriangles() in a pre-process.");
2785 return -1;
2786 }
2787 if(ltid < 0 || ltid >= this->getNumberOfTriangles()) {
2788 return -1;
2789 }
2790#endif // TTK_ENABLE_KAMIKAZE
2791 if(!ttk::isRunningWithMPI()) {
2792 return ltid;
2793 }
2794 if(dim == 3) {
2795 return this->getTriangleGlobalIdInternal(ltid);
2796 } else if(dim == 2) {
2797 return this->getCellGlobalIdInternal(ltid);
2798 }
2799 return -1;
2800 }
2801 virtual inline SimplexId getTriangleLocalId(const SimplexId gtid) const {
2802 const auto dim{this->getDimensionality()};
2803#ifndef TTK_ENABLE_KAMIKAZE
2804 if(dim != 3 && dim != 2) {
2805 this->printErr("Only 2D and 3D datasets are supported");
2806 return -1;
2807 }
2808 if(!this->hasPreconditionedDistributedEdges_) {
2809 this->printErr("TriangleLocalId query without pre-process!");
2810 this->printErr(
2811 "Please call preconditionDistributedTriangles() in a pre-process.");
2812 return -1;
2813 }
2814#endif // TTK_ENABLE_KAMIKAZE
2815 if(!ttk::isRunningWithMPI()) {
2816 return gtid;
2817 }
2818 if(dim == 3) {
2819 return this->getTriangleLocalIdInternal(gtid);
2820 } else if(dim == 2) {
2821 return this->getCellLocalIdInternal(gtid);
2822 }
2823 return -1;
2824 }
2825
2826 virtual inline int getVertexRank(const SimplexId lvid) const {
2827
2828 if(!ttk::isRunningWithMPI()) {
2829 return 0;
2830 }
2831#ifndef TTK_ENABLE_KAMIKAZE
2832 if(this->getDimensionality() != 1 && this->getDimensionality() != 2
2833 && this->getDimensionality() != 3) {
2834 this->printErr("Only 1D, 2D and 3D datasets are supported");
2835 return -1;
2836 }
2837 if(!this->hasPreconditionedDistributedVertices_) {
2838 this->printErr("VertexRankId query without pre-process!");
2839 this->printErr(
2840 "Please call preconditionDistributedVertices() in a pre-process.");
2841 return -1;
2842 }
2843 if(lvid < 0 || lvid >= this->getNumberOfVertices()) {
2844 return -1;
2845 }
2846#endif // TTK_ENABLE_KAMIKAZE
2847 return this->getVertexRankInternal(lvid);
2848 }
2849
2850 virtual inline int getCellRank(const SimplexId lcid) const {
2851#ifndef TTK_ENABLE_KAMIKAZE
2852 if(this->getDimensionality() != 1 && this->getDimensionality() != 2
2853 && this->getDimensionality() != 3) {
2854 this->printErr("Only 1D, 2D and 3D datasets are supported");
2855 return -1;
2856 }
2857 if(!this->hasPreconditionedDistributedCells_) {
2858 this->printErr("CellRankId query without pre-process!");
2859 this->printErr(
2860 "Please call preconditionDistributedCells() in a pre-process.");
2861 return -1;
2862 }
2863 if(lcid < 0 || lcid >= this->getNumberOfCells()) {
2864 return -1;
2865 }
2866#endif // TTK_ENABLE_KAMIKAZE
2867 if(!ttk::isRunningWithMPI()) {
2868 return lcid;
2869 }
2870 return this->getCellRankInternal(lcid);
2871 }
2872
2873 virtual inline const std::vector<int> &getNeighborRanks() const {
2874 return this->neighborRanks_;
2875 }
2876 virtual inline std::vector<int> &getNeighborRanks() {
2877 return this->neighborRanks_;
2878 }
2879
2880 virtual inline const std::vector<std::array<ttk::SimplexId, 6>> &
2881 getNeighborVertexBBoxes() const {
2882 return this->neighborVertexBBoxes_;
2883 }
2884
2885 virtual inline const std::vector<std::vector<SimplexId>> &
2886 getGhostCellsPerOwner() const {
2887 if(!hasPreconditionedExchangeGhostCells_) {
2888 printErr("The ghostCellsOwner_ attribute has not been populated!");
2889 printErr(
2890 "Please call preconditionExchangeGhostCells in a pre-process.");
2891 }
2892 return this->ghostCellsPerOwner_;
2893 }
2894
2895 virtual inline const std::vector<std::vector<SimplexId>> &
2896 getRemoteGhostCells() const {
2897 if(!hasPreconditionedExchangeGhostCells_) {
2898 printErr("The remoteGhostCells_ attribute has not been populated!");
2899 printErr(
2900 "Please call preconditionExchangeGhostCells in a pre-process.");
2901 }
2902 return this->remoteGhostCells_;
2903 }
2904
2905 virtual inline const std::vector<std::vector<SimplexId>> &
2906 getGhostVerticesPerOwner() const {
2907 if(!hasPreconditionedExchangeGhostVertices_) {
2908 printErr(
2909 "The ghostVerticesPerOwner_ attribute has not been populated!");
2910 printErr(
2911 "Please call preconditionExchangeGhostVertices in a pre-process.");
2912 }
2913 return this->ghostVerticesPerOwner_;
2914 }
2915
2916 virtual inline const std::vector<std::vector<SimplexId>> &
2917 getRemoteGhostVertices() const {
2918 if(!hasPreconditionedExchangeGhostVertices_) {
2919 printErr("The remoteGhostVertices_ attribute has not been populated!");
2920 printErr(
2921 "Please call preconditionExchangeGhostVertices in a pre-process.");
2922 }
2923 return this->remoteGhostVertices_;
2924 }
2925
2926 virtual inline void setHasPreconditionedDistributedVertices(bool flag) {
2927 this->hasPreconditionedDistributedVertices_ = flag;
2928 }
2929
2930 virtual inline bool hasPreconditionedDistributedVertices() const {
2931 return this->hasPreconditionedDistributedVertices_;
2932 }
2933 virtual inline bool hasPreconditionedDistributedCells() const {
2934 return this->hasPreconditionedDistributedCells_;
2935 }
2936
2937 inline int getDistributedGlobalCellId(const SimplexId &localCellId,
2938 const int &cellDim,
2939 SimplexId &globalCellId) const {
2940 if(ttk::hasInitializedMPI()) {
2941 switch(cellDim) {
2942 case 0:
2943 globalCellId = this->getVertexGlobalIdInternal(localCellId);
2944 break;
2945 case 1:
2946 globalCellId = this->getEdgeGlobalIdInternal(localCellId);
2947 break;
2948 case 2:
2949 if(getDimensionality() == 2) {
2950 globalCellId = this->getCellGlobalIdInternal(localCellId);
2951 break;
2952 } else {
2953 globalCellId = this->getTriangleGlobalIdInternal(localCellId);
2954 break;
2955 }
2956 case 3: {
2957 globalCellId = this->getCellGlobalIdInternal(localCellId);
2958 break;
2959 }
2960 default:
2961 globalCellId = -1;
2962 break;
2963 }
2964 } else {
2965 globalCellId = localCellId;
2966 }
2967 return 0;
2968 }
2969
2970 protected:
2971 virtual inline SimplexId
2972 getVertexGlobalIdInternal(const SimplexId ttkNotUsed(lvid)) const {
2973 return -1;
2974 }
2975
2976 virtual inline SimplexId
2977 getVertexLocalIdInternal(const SimplexId ttkNotUsed(gvid)) const {
2978 return -1;
2979 }
2980
2981 // overriden in ImplicitTriangulation &
2982 // PeriodicImplicitTriangulation (where cellGid_ refers to
2983 // squares/cubes and not triangles/tetrahedron)
2984 virtual inline SimplexId
2985 getCellGlobalIdInternal(const SimplexId ttkNotUsed(lcid)) const {
2986 return -1;
2987 }
2988
2989 virtual inline SimplexId
2990 getCellLocalIdInternal(const SimplexId ttkNotUsed(gcid)) const {
2991 return -1;
2992 }
2993
2994 virtual inline SimplexId
2995 getEdgeGlobalIdInternal(const SimplexId ttkNotUsed(leid)) const {
2996 return -1;
2997 }
2998
2999 virtual inline SimplexId
3000 getEdgeLocalIdInternal(const SimplexId ttkNotUsed(geid)) const {
3001 return -1;
3002 }
3003
3004 virtual inline SimplexId
3005 getTriangleGlobalIdInternal(const SimplexId ttkNotUsed(ltid)) const {
3006 return -1;
3007 }
3008
3009 virtual inline SimplexId
3010 getTriangleLocalIdInternal(const SimplexId ttkNotUsed(gtid)) const {
3011 return -1;
3012 }
3013
3014 virtual inline int
3015 getVertexRankInternal(const SimplexId ttkNotUsed(lvid)) const {
3016 return -1;
3017 }
3018
3019 virtual inline int
3020 getCellRankInternal(const SimplexId ttkNotUsed(lcid)) const {
3021 return -1;
3022 }
3023
3024 // these protected methods should not be exposed since they are
3025 // called by their non-MPI-aware conterparts
3026 virtual inline bool
3027 isVertexOnGlobalBoundaryInternal(const SimplexId ttkNotUsed(lvid)) const {
3028 return false;
3029 }
3030 virtual inline bool
3031 isEdgeOnGlobalBoundaryInternal(const SimplexId ttkNotUsed(leid)) const {
3032 return false;
3033 }
3034 virtual inline bool isTriangleOnGlobalBoundaryInternal(
3035 const SimplexId ttkNotUsed(ltid)) const {
3036 return false;
3037 }
3038
3039#endif // TTK_ENABLE_MPI
3040
3041 protected:
3042 virtual int getCellEdgeInternal(const SimplexId &ttkNotUsed(cellId),
3043 const int &ttkNotUsed(localEdgeId),
3044 SimplexId &ttkNotUsed(edgeId)) const {
3045 return 0;
3046 }
3047
3048 virtual inline SimplexId
3050 return 0;
3051 }
3052
3053 virtual inline const std::vector<std::vector<SimplexId>> *
3055 return nullptr;
3056 }
3057
3058 virtual inline int
3060 const int &ttkNotUsed(localNeighborId),
3061 SimplexId &ttkNotUsed(neighborId)) const {
3062 return 0;
3063 }
3064
3065 virtual inline SimplexId
3067 return 0;
3068 }
3069
3070 virtual inline const std::vector<std::vector<SimplexId>> *
3072 return nullptr;
3073 }
3074
3075 virtual inline int
3077 const int &ttkNotUsed(localTriangleId),
3078 SimplexId &ttkNotUsed(triangleId)) const {
3079 return 0;
3080 }
3081
3082 virtual inline SimplexId
3084 return 0;
3085 }
3086
3087 virtual inline const std::vector<std::vector<SimplexId>> *
3089 return nullptr;
3090 }
3091
3092 virtual inline int
3094 const int &ttkNotUsed(localVertexId),
3095 SimplexId &ttkNotUsed(vertexId)) const {
3096 return 0;
3097 }
3098
3099 virtual inline SimplexId
3101 return 0;
3102 }
3103
3104 virtual inline int getDimensionalityInternal() const {
3105 return 0;
3106 }
3107
3108 virtual inline const std::vector<std::array<SimplexId, 2>> *
3110 return nullptr;
3111 }
3112
3113 virtual inline int
3115 const int &ttkNotUsed(localLinkId),
3116 SimplexId &ttkNotUsed(linkId)) const {
3117 return 0;
3118 }
3119
3120 virtual inline SimplexId
3122 return 0;
3123 }
3124
3125 virtual inline const std::vector<std::vector<SimplexId>> *
3127 return nullptr;
3128 }
3129
3130 virtual inline int
3132 const int &ttkNotUsed(localStarId),
3133 SimplexId &ttkNotUsed(starId)) const {
3134 return 0;
3135 }
3136
3137 virtual inline SimplexId
3139 return 0;
3140 }
3141
3142 virtual inline const std::vector<std::vector<SimplexId>> *
3144 return nullptr;
3145 }
3146
3147 virtual inline int
3149 const int &ttkNotUsed(localTriangleId),
3150 SimplexId &ttkNotUsed(triangleId)) const {
3151 return 0;
3152 }
3153
3154 virtual inline SimplexId
3156 return 0;
3157 }
3158
3159 virtual inline const std::vector<std::vector<SimplexId>> *
3161 return nullptr;
3162 }
3163
3164 virtual inline int
3166 const int &ttkNotUsed(localVertexId),
3167 SimplexId &ttkNotUsed(vertexId)) const {
3168 return 0;
3169 }
3170
3171 virtual inline SimplexId getNumberOfCellsInternal() const {
3172 return 0;
3173 }
3174
3175 virtual inline SimplexId getNumberOfEdgesInternal() const {
3176 return 0;
3177 }
3178
3180 return 0;
3181 }
3182
3184 return 0;
3185 }
3186
3187 virtual inline const std::vector<std::array<SimplexId, 3>> *
3189 return nullptr;
3190 }
3191
3192 virtual inline int
3194 const int &ttkNotUsed(localEdgeId),
3195 SimplexId &ttkNotUsed(edgeId)) const {
3196 return 0;
3197 }
3198
3200 const SimplexId &ttkNotUsed(triangleId)) const {
3201 return 0;
3202 }
3203
3204 virtual inline const std::vector<std::vector<SimplexId>> *
3206 return nullptr;
3207 }
3208
3209 virtual inline int
3211 const int &ttkNotUsed(localLinkId),
3212 SimplexId &ttkNotUsed(linkId)) const {
3213 return 0;
3214 }
3215
3217 const SimplexId &ttkNotUsed(triangleId)) const {
3218 return 0;
3219 }
3220
3221 virtual inline const std::vector<std::vector<SimplexId>> *
3223 return nullptr;
3224 }
3225
3226 virtual inline int
3228 const int &ttkNotUsed(localStarId),
3229 SimplexId &ttkNotUsed(starId)) const {
3230 return 0;
3231 }
3232
3234 const SimplexId &ttkNotUsed(triangleId)) const {
3235 return 0;
3236 }
3237
3238 virtual inline const std::vector<std::vector<SimplexId>> *
3240 return nullptr;
3241 }
3242
3243 virtual inline int
3245 const int &ttkNotUsed(localVertexId),
3246 SimplexId &ttkNotUsed(vertexId)) const {
3247 return 0;
3248 }
3249
3250 virtual inline int
3252 const int &ttkNotUsed(localEdgeId),
3253 SimplexId &ttkNotUsed(edgeId)) const {
3254 return 0;
3255 }
3256
3257 virtual inline SimplexId
3259 return 0;
3260 }
3261
3262 virtual inline const std::vector<std::vector<SimplexId>> *
3264 return nullptr;
3265 }
3266
3267 virtual inline int
3269 const int &ttkNotUsed(localLinkId),
3270 SimplexId &ttkNotUsed(linkId)) const {
3271 return 0;
3272 }
3273
3274 virtual inline SimplexId
3276 return 0;
3277 }
3278
3279 virtual inline const std::vector<std::vector<SimplexId>> *
3281 return nullptr;
3282 }
3283
3284 virtual inline int
3286 const int &ttkNotUsed(localNeighborId),
3287 SimplexId &ttkNotUsed(neighborId)) const {
3288 return 0;
3289 }
3290
3292 const SimplexId &ttkNotUsed(vertexId)) const {
3293 return 0;
3294 }
3295
3296 virtual inline const std::vector<std::vector<SimplexId>> *
3298 return nullptr;
3299 }
3300
3301 virtual inline int
3303 float &ttkNotUsed(x),
3304 float &ttkNotUsed(y),
3305 float &ttkNotUsed(z)) const {
3306 return 0;
3307 }
3308
3309 virtual inline int
3311 const int &ttkNotUsed(localStarId),
3312 SimplexId &ttkNotUsed(starId)) const {
3313 return 0;
3314 }
3315
3316 virtual inline SimplexId
3318 return 0;
3319 }
3320
3321 virtual inline const std::vector<std::vector<SimplexId>> *
3323 return nullptr;
3324 }
3325
3326 virtual inline int
3328 const int &ttkNotUsed(localTriangleId),
3329 SimplexId &ttkNotUsed(triangleId)) const {
3330 return 0;
3331 }
3332
3334 const SimplexId &ttkNotUsed(vertexId)) const {
3335 return 0;
3336 }
3337
3338 virtual inline const std::vector<std::vector<SimplexId>> *
3340 return nullptr;
3341 }
3342
3343 virtual inline bool
3345 return false;
3346 }
3347
3348 virtual inline bool isTriangleOnBoundaryInternal(
3349 const SimplexId &ttkNotUsed(triangleId)) const {
3350 return false;
3351 }
3352
3353 virtual inline bool
3355 return false;
3356 }
3357
3359
3360#ifndef TTK_ENABLE_KAMIKAZE
3362
3363 printMsg("BoundaryEdge query without pre-process!",
3365 printMsg("Please call preconditionBoundaryEdges() in a pre-process.",
3367 }
3368#endif
3369
3371 }
3372
3374#ifndef TTK_ENABLE_KAMIKAZE
3376
3377 printMsg("BoundaryTriangle query without pre-process!",
3379 printMsg(
3380 "Please call preconditionBoundaryTriangles() in a pre-process.",
3382 }
3383#endif
3385 }
3386
3388#ifndef TTK_ENABLE_KAMIKAZE
3390
3391 printMsg("BoundaryVertex query without pre-process!",
3393 printMsg("Please call preconditionBoundaryVertices() in a pre-process.",
3395 }
3396#endif
3398 }
3399
3400 inline bool hasPreconditionedCellEdges() const {
3401
3402#ifndef TTK_ENABLE_KAMIKAZE
3405
3406 printMsg("CellEdge query without pre-process!", debug::Priority::ERROR,
3407 debug::LineMode::NEW, std::cerr);
3408 printMsg("Please call preconditionCellEdges() in a pre-process.",
3410 return false;
3411 }
3412#endif
3413 return true;
3414 }
3415
3417#ifndef TTK_ENABLE_KAMIKAZE
3419
3420 printMsg("CellNeighbor query without pre-process!",
3422 printMsg("Please call preconditionCellNeighbors() in a pre-process.",
3424 }
3425#endif
3427 }
3428
3430#ifndef TTK_ENABLE_KAMIKAZE
3432 || ((getDimensionality() == 3)
3434
3435 printMsg("CellTriangle query without pre-process!",
3437 printMsg("Please call preconditionCellTriangles() in a pre-process.",
3439
3440 return false;
3441 }
3442#endif
3443 return true;
3444 }
3445
3446 inline bool hasPreconditionedEdgeLinks() const {
3447#ifndef TTK_ENABLE_KAMIKAZE
3449
3450 printMsg("EdgeLink query without pre-process!", debug::Priority::ERROR,
3451 debug::LineMode::NEW, std::cerr);
3452 printMsg("Please call preconditionEdgeLinks() in a pre-process.",
3454 }
3455#endif
3457 }
3458
3459 inline bool hasPreconditionedEdgeStars() const {
3460#ifndef TTK_ENABLE_KAMIKAZE
3462
3463 printMsg("EdgeStar query without pre-process!", debug::Priority::ERROR,
3464 debug::LineMode::NEW, std::cerr);
3465 printMsg("Please call preconditionEdgeStars() in a pre-process.",
3467 }
3468#endif
3470 }
3471
3473#ifndef TTK_ENABLE_KAMIKAZE
3475 || ((getDimensionality() == 3)
3477
3478 printMsg("EdgeTriangle query without pre-process!",
3480 printMsg("Please call preconditionEdgeTriangles() in a pre-process.",
3482
3483 return false;
3484 }
3485#endif
3486 return true;
3487 }
3488
3489 inline bool hasPreconditionedEdges() const {
3490#ifndef TTK_ENABLE_KAMIKAZE
3491 if((getDimensionality() != 1) && (!hasPreconditionedEdges_)) {
3492
3493 printMsg("Edge query without pre-process!", debug::Priority::ERROR,
3494 debug::LineMode::NEW, std::cerr);
3495 printMsg("Please call preconditionEdges() in a pre-process.",
3497
3498 return false;
3499 }
3500#endif
3501 return true;
3502 }
3503
3504 inline bool hasPreconditionedTriangles() const {
3505#ifndef TTK_ENABLE_KAMIKAZE
3507
3508 printMsg("Triangle query without pre-process!", debug::Priority::ERROR,
3509 debug::LineMode::NEW, std::cerr);
3510 printMsg("Please call preconditionTriangles() in a pre-process.",
3512
3513 return false;
3514 }
3515#endif
3516 return true;
3517 }
3518
3520#ifndef TTK_ENABLE_KAMIKAZE
3522 || ((getDimensionality() == 3)
3524
3525 printMsg("TriangleEdge query without pre-process!",
3527 printMsg("Please call preconditionTriangleEdges() in a pre-process.",
3529
3530 return false;
3531 }
3532#endif
3533 return true;
3534 }
3535
3537#ifndef TTK_ENABLE_KAMIKAZE
3539
3540 printMsg("TriangleLink query without pre-process!",
3542 printMsg("Please call preconditionTriangleLinks() in a pre-process.",
3544 }
3545#endif
3547 }
3548
3550#ifndef TTK_ENABLE_KAMIKAZE
3552
3553 printMsg("TriangleStar query without pre-process!",
3555 printMsg("Please call preconditionTriangleStars() in a pre-process.",
3557 }
3558#endif
3560 }
3561
3562 inline bool hasPreconditionedVertexEdges() const {
3563#ifndef TTK_ENABLE_KAMIKAZE
3566
3567 printMsg("VertexEdge query without pre-process!",
3569 printMsg("Please call preconditionVertexEdges() in a pre-process.",
3571
3572 return false;
3573 }
3574#endif
3575 return true;
3576 }
3577
3578 inline bool hasPreconditionedVertexLinks() const {
3579#ifndef TTK_ENABLE_KAMIKAZE
3581
3582 printMsg("VertexLink query without pre-process!",
3584 printMsg("Please call preconditionVertexLinks() in a pre-process.",
3586 }
3587#endif
3589 }
3590
3592#ifndef TTK_ENABLE_KAMIKAZE
3594
3595 printMsg("VertexNeighbor query without pre-process!",
3597 printMsg("Please call preconditionVertexNeighbors() in a pre-process.",
3599 }
3600#endif
3602 }
3603
3604 inline bool hasPreconditionedVertexStars() const {
3605#ifndef TTK_ENABLE_KAMIKAZE
3607
3608 printMsg("VertexStar query without pre-process!",
3610 printMsg("Please call preconditionVertexStars() in a pre-process.",
3612 }
3613#endif
3615 }
3616
3618#ifndef TTK_ENABLE_KAMIKAZE
3620 || ((getDimensionality() == 3)
3622
3623 printMsg("VertexTriangle query without pre-process!",
3625 printMsg("Please call preconditionVertexTriangles() in a pre-process.",
3627
3628 return false;
3629 }
3630#endif
3631 return true;
3632 }
3633
3634 inline bool hasPreconditionedManifold() const {
3635#ifndef TTK_ENABLE_KAMIKAZE
3636 if(!this->hasPreconditionedManifold_) {
3637 this->printMsg("isManifold query without pre-process!",
3639 this->printMsg("Please call preconditionManifold() in a pre-process.",
3641 }
3642#endif // TTK_ENABLE_KAMIKAZE
3643 return this->hasPreconditionedManifold_;
3644 }
3645
3646 // empty wrapping to VTK for now
3647 bool needsToAbort() override {
3648 return false;
3649 }
3650
3652 return 0;
3653 }
3654
3656 return 0;
3657 }
3658
3660 return 0;
3661 }
3662
3663 virtual inline int preconditionCellEdgesInternal() {
3664 return 0;
3665 }
3666
3668 return 0;
3669 }
3670
3672 return 0;
3673 }
3674
3675 virtual inline int preconditionEdgesInternal() {
3676 return 0;
3677 }
3678
3679 virtual inline int preconditionEdgeLinksInternal() {
3680 return 0;
3681 }
3682
3683 virtual inline int preconditionEdgeStarsInternal() {
3684 return 0;
3685 }
3686
3688 return 0;
3689 }
3690
3691 virtual inline int preconditionTrianglesInternal() {
3692 return 0;
3693 }
3694
3696 return 0;
3697 }
3698
3700 return 0;
3701 }
3702
3704 return 0;
3705 }
3706
3707 virtual inline int preconditionVertexEdgesInternal() {
3708 return 0;
3709 }
3710
3711 virtual inline int preconditionVertexLinksInternal() {
3712 return 0;
3713 }
3714
3716 return 0;
3717 }
3718
3719 virtual inline int preconditionVertexStarsInternal() {
3720 return 0;
3721 }
3722
3724 return 0;
3725 }
3726
3727 virtual inline int preconditionManifoldInternal() {
3728 return 0;
3729 }
3730
3731 virtual inline int getCellVTKIDInternal(const int &ttkId,
3732 int &vtkId) const {
3733#ifndef TTK_ENABLE_KAMIKAZE
3734 if(ttkId < 0) {
3735 return -1;
3736 }
3737#endif
3738 vtkId = ttkId;
3739 return 0;
3740 }
3741
3742 template <class itemType>
3743 size_t tableFootprint(const std::vector<itemType> &table,
3744 const std::string &tableName = "",
3745 std::ostream &stream = std::cout) const {
3746
3747 std::stringstream msg;
3748 if((table.size()) && (tableName.length())) {
3749 msg << tableName << ": " << table.size() * sizeof(itemType) << " bytes";
3750 printMsg(
3751 msg.str(), debug::Priority::INFO, debug::LineMode::NEW, stream);
3752 }
3753
3754 return table.size() * sizeof(itemType);
3755 }
3756
3757 template <class itemType>
3758 size_t tableTableFootprint(const std::vector<std::vector<itemType>> &table,
3759 const std::string &tableName = "",
3760 std::ostream &stream = std::cout) const;
3761
3762 int updateProgress(const float &ttkNotUsed(progress)) override {
3763 return 0;
3764 }
3765
3777
3778 bool isManifold_{true};
3779
3780 std::array<SimplexId, 3> gridDimensions_;
3781
3783 std::vector<std::array<SimplexId, 6>> tetraEdgeList_;
3784 std::vector<std::vector<SimplexId>> cellNeighborList_;
3785 std::vector<std::array<SimplexId, 4>> tetraTriangleList_;
3786 std::vector<std::vector<SimplexId>> edgeLinkList_;
3787 std::vector<std::array<SimplexId, 2>> edgeList_;
3788 std::vector<std::vector<SimplexId>> edgeStarList_;
3789 std::vector<std::vector<SimplexId>> edgeTriangleList_;
3790 std::vector<std::array<SimplexId, 3>> triangleList_;
3791 std::vector<std::array<SimplexId, 3>> triangleEdgeList_;
3792 std::vector<std::vector<SimplexId>> triangleLinkList_;
3793 std::vector<std::vector<SimplexId>> triangleStarList_;
3794 std::vector<std::vector<SimplexId>> vertexEdgeList_;
3795 std::vector<std::vector<SimplexId>> vertexLinkList_;
3796 std::vector<std::vector<SimplexId>> vertexNeighborList_;
3797 std::vector<std::vector<SimplexId>> vertexStarList_;
3798 std::vector<std::vector<SimplexId>> vertexTriangleList_;
3799
3800 // keep compatibility between getCellEdges(), getCellTriangles(),
3801 // getCellNeighbors() and getTriangleEdges()
3802 std::vector<std::vector<SimplexId>> cellEdgeVector_{};
3803 std::vector<std::vector<SimplexId>> cellTriangleVector_{};
3804 std::vector<std::vector<SimplexId>> triangleEdgeVector_{};
3805
3806#ifdef TTK_ENABLE_MPI
3807
3808 // precondition methods for distributed meshes
3809
3810 virtual int preconditionDistributedEdges() {
3811 return 0;
3812 }
3813 virtual int preconditionDistributedTriangles() {
3814 return 0;
3815 }
3816
3817 // "vtkGhostType" PointData array
3818 const unsigned char *vertexGhost_{};
3819 // "vtkGhostType" CellData array
3820 const unsigned char *cellGhost_{};
3821
3822 // list of neighboring ranks (sharing ghost cells to current rank)
3823 std::vector<int> neighborRanks_{};
3824 // global ids of (local) ghost cells per each MPI (neighboring) rank
3825 std::vector<std::vector<SimplexId>> ghostCellsPerOwner_{};
3826 // global ids of local (owned) cells that are ghost cells of other
3827 // (neighboring) ranks (per MPI rank)
3828 std::vector<std::vector<SimplexId>> remoteGhostCells_{};
3829 // global ids of (local) ghost vertices per each MPI (neighboring) rank
3830 std::vector<std::vector<SimplexId>> ghostVerticesPerOwner_{};
3831 // global ids of local (owned) vertices that are ghost cells of other
3832 // (neighboring) ranks (per MPI rank)
3833 std::vector<std::vector<SimplexId>> remoteGhostVertices_{};
3834 // hold the neighboring ranks vertex bounding boxes (metaGrid_ coordinates)
3835 std::vector<std::array<SimplexId, 6>> neighborVertexBBoxes_{};
3836 // hold the neighboring ranks cells bounding boxes (metaGrid_ coordinates)
3837 std::vector<std::array<SimplexId, 6>> neighborCellBBoxes_{};
3838
3839 bool hasPreconditionedDistributedCells_{false};
3840 bool hasPreconditionedExchangeGhostCells_{false};
3841 bool hasPreconditionedDistributedEdges_{false};
3842 bool hasPreconditionedDistributedTriangles_{false};
3843 bool hasPreconditionedDistributedVertices_{false};
3844 bool hasPreconditionedExchangeGhostVertices_{false};
3845 bool hasPreconditionedGlobalBoundary_{false};
3846
3847#endif // TTK_ENABLE_MPI
3848
3849 // only ttk::dcg::DiscreteGradient should use what's defined below.
3851
3852#ifdef TTK_ENABLE_DCG_OPTIMIZE_MEMORY
3853 // might no be sufficient for Rips complexes, where a vertex can
3854 // have more than 128 neighbors
3855 using gradIdType = char;
3856#else
3858#endif
3859
3873 using gradientType = std::array<std::vector<gradIdType>, 6>;
3882 using gradientKeyType = std::pair<const void *, size_t>;
3883 /*
3884 * @brief Type for caching Discrete Gradient internal data structure.
3885 *
3886 * Uses the ttk::LRUCache with \ref gradientKeytype as key and
3887 * \ref gradientType as value types.
3888 */
3890 /*
3891 * @brief Access to the gradientCache_ mutable member variable
3892 *
3893 * \warning This should only be used by the
3894 * ttk::dcg::DiscreteGradient class.
3895 */
3897 return &this->gradientCache_;
3898 }
3899
3900 // store, for each triangulation object and per offset field, a
3901 // reference to the discrete gradient internal structure
3903 };
3904} // namespace ttk
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
Definition: BaseClass.h:47
AbstractTriangulation is an interface class that defines an interface for efficient traversal methods...
virtual const std::vector< std::vector< SimplexId > > * getCellNeighborsInternal()
virtual const std::vector< std::vector< SimplexId > > * getVertexNeighbors()
virtual int getVertexPoint(const SimplexId &vertexId, float &x, float &y, float &z) const
std::vector< std::vector< SimplexId > > vertexStarList_
virtual SimplexId getVertexEdgeNumberInternal(const SimplexId &ttkNotUsed(vertexId)) const
std::pair< const void *, size_t > gradientKeyType
Key type for gradientCacheType.
virtual int getEdgeLink(const SimplexId &edgeId, const int &localLinkId, SimplexId &linkId) const
virtual int getEdgeVertexInternal(const SimplexId &ttkNotUsed(edgeId), const int &ttkNotUsed(localVertexId), SimplexId &ttkNotUsed(vertexId)) const
virtual SimplexId getVertexTriangleNumber(const SimplexId &vertexId) const
virtual const std::vector< std::vector< SimplexId > > * getCellNeighbors()
virtual int getVertexTriangle(const SimplexId &vertexId, const int &localTriangleId, SimplexId &triangleId) const
virtual const std::vector< std::vector< SimplexId > > * getCellEdgesInternal()
virtual int getVertexEdge(const SimplexId &vertexId, const int &localEdgeId, SimplexId &edgeId) const
int getCellIncenter(const SimplexId cellid, const int dim, float incenter[3]) const
virtual SimplexId getEdgeTriangleNumber(const SimplexId &edgeId) const
virtual const std::array< SimplexId, 3 > & getGridDimensions() const
virtual bool isTriangleOnBoundary(const SimplexId &triangleId) const
virtual SimplexId getVertexTriangleNumberInternal(const SimplexId &ttkNotUsed(vertexId)) const
virtual const std::vector< std::vector< SimplexId > > * getVertexStars()
int getEdgeVertexNumber(const SimplexId ttkNotUsed(edgeId)) const
virtual const std::vector< std::vector< SimplexId > > * getEdgeStars()
virtual int getTriangleLinkInternal(const SimplexId &ttkNotUsed(triangleId), const int &ttkNotUsed(localLinkId), SimplexId &ttkNotUsed(linkId)) const
virtual SimplexId getNumberOfEdgesInternal() const
virtual const std::vector< std::vector< SimplexId > > * getVertexLinks()
virtual const std::vector< std::array< SimplexId, 2 > > * getEdgesInternal()
virtual const std::vector< std::vector< SimplexId > > * getVertexTriangles()
virtual int getTriangleStar(const SimplexId &triangleId, const int &localStarId, SimplexId &starId) const
virtual SimplexId getVertexNeighborNumber(const SimplexId &vertexId) const
bool hasPeriodicBoundaries() const
Returns true if the grid uses period boundary conditions.
virtual int getEdgeTriangle(const SimplexId &edgeId, const int &localTriangleId, SimplexId &triangleId) const
virtual const std::vector< std::vector< SimplexId > > * getTriangleEdgesInternal()
std::vector< std::array< SimplexId, 3 > > triangleList_
virtual int getCellEdgeInternal(const SimplexId &ttkNotUsed(cellId), const int &ttkNotUsed(localEdgeId), SimplexId &ttkNotUsed(edgeId)) const
virtual int getEdgeVertex(const SimplexId &edgeId, const int &localVertexId, SimplexId &vertexId) const
virtual const std::vector< std::vector< SimplexId > > * getVertexStarsInternal()
virtual int preconditionVertexTrianglesInternal()
virtual int getEdgeLinkInternal(const SimplexId &ttkNotUsed(edgeId), const int &ttkNotUsed(localLinkId), SimplexId &ttkNotUsed(linkId)) const
AbstractTriangulation & operator=(AbstractTriangulation &&)=default
virtual SimplexId getNumberOfVerticesInternal() const
virtual const std::vector< std::array< SimplexId, 2 > > * getEdges()
virtual bool isVertexOnBoundaryInternal(const SimplexId &ttkNotUsed(vertexId)) const
virtual int getTriangleEdge(const SimplexId &triangleId, const int &localEdgeId, SimplexId &edgeId) const
virtual SimplexId getTriangleLinkNumberInternal(const SimplexId &ttkNotUsed(triangleId)) const
virtual SimplexId getNumberOfCellsInternal() const
std::vector< bool > boundaryTriangles_
std::vector< std::array< SimplexId, 4 > > tetraTriangleList_
int updateProgress(const float &ttkNotUsed(progress)) override
virtual SimplexId getTriangleLinkNumber(const SimplexId &triangleId) const
AbstractTriangulation & operator=(const AbstractTriangulation &)=default
std::vector< std::vector< SimplexId > > triangleLinkList_
virtual int getTriangleStarInternal(const SimplexId &ttkNotUsed(triangleId), const int &ttkNotUsed(localStarId), SimplexId &ttkNotUsed(starId)) const
virtual bool isVertexOnBoundary(const SimplexId &vertexId) const
virtual const std::vector< std::vector< SimplexId > > * getVertexEdgesInternal()
virtual const std::vector< std::vector< SimplexId > > * getTriangleLinks()
virtual SimplexId getCellEdgeNumber(const SimplexId &cellId) const
std::vector< std::vector< SimplexId > > vertexNeighborList_
std::vector< std::vector< SimplexId > > vertexLinkList_
virtual int getVertexTriangleInternal(const SimplexId &ttkNotUsed(vertexId), const int &ttkNotUsed(localTriangleId), SimplexId &ttkNotUsed(triangleId)) const
gradientCacheType * getGradientCacheHandler() const
virtual SimplexId getEdgeStarNumberInternal(const SimplexId &ttkNotUsed(edgeId)) const
virtual const std::vector< std::vector< SimplexId > > * getEdgeTrianglesInternal()
virtual bool isEdgeOnBoundary(const SimplexId &edgeId) const
virtual int getVertexNeighborInternal(const SimplexId &ttkNotUsed(vertexId), const int &ttkNotUsed(localNeighborId), SimplexId &ttkNotUsed(neighborId)) const
virtual int getCellNeighborInternal(const SimplexId &ttkNotUsed(cellId), const int &ttkNotUsed(localNeighborId), SimplexId &ttkNotUsed(neighborId)) const
AbstractTriangulation(const AbstractTriangulation &)=default
virtual const std::vector< std::vector< SimplexId > > * getEdgeTriangles()
virtual SimplexId getCellVertexNumberInternal(const SimplexId &ttkNotUsed(cellId)) const
virtual int getCellVTKIDInternal(const int &ttkId, int &vtkId) const
virtual const std::vector< std::array< SimplexId, 3 > > * getTriangles()
virtual int getTriangleLink(const SimplexId &triangleId, const int &localLinkId, SimplexId &linkId) const
virtual const std::vector< std::array< SimplexId, 3 > > * getTrianglesInternal()
virtual SimplexId getEdgeLinkNumberInternal(const SimplexId &ttkNotUsed(edgeId)) const
virtual int getVertexPointInternal(const SimplexId &ttkNotUsed(vertexId), float &ttkNotUsed(x), float &ttkNotUsed(y), float &ttkNotUsed(z)) const
virtual const std::vector< std::vector< SimplexId > > * getEdgeLinks()
virtual int getTriangleVertexInternal(const SimplexId &ttkNotUsed(triangleId), const int &ttkNotUsed(localVertexId), SimplexId &ttkNotUsed(vertexId)) const
std::vector< std::vector< SimplexId > > edgeLinkList_
virtual int getVertexLinkInternal(const SimplexId &ttkNotUsed(vertexId), const int &ttkNotUsed(localLinkId), SimplexId &ttkNotUsed(linkId)) const
std::vector< bool > boundaryVertices_
std::vector< std::vector< SimplexId > > vertexEdgeList_
std::vector< std::vector< SimplexId > > triangleEdgeVector_
std::vector< std::vector< SimplexId > > triangleStarList_
virtual SimplexId getNumberOfCells() const
virtual SimplexId getCellTriangleNumber(const SimplexId &cellId) const
virtual SimplexId getNumberOfVertices() const
virtual SimplexId getVertexNeighborNumberInternal(const SimplexId &ttkNotUsed(vertexId)) const
virtual int preconditionVertexNeighborsInternal()
virtual SimplexId getNumberOfTrianglesInternal() const
virtual int getEdgeStarInternal(const SimplexId &ttkNotUsed(edgeId), const int &ttkNotUsed(localStarId), SimplexId &ttkNotUsed(starId)) const
virtual SimplexId getCellNeighborNumber(const SimplexId &cellId) const
virtual int getDimensionality() const
virtual int getCellEdge(const SimplexId &cellId, const int &localEdgeId, SimplexId &edgeId) const
virtual int getEdgeStar(const SimplexId &edgeId, const int &localStarId, SimplexId &starId) const
virtual int getVertexNeighbor(const SimplexId &vertexId, const int &localNeighborId, SimplexId &neighborId) const
virtual SimplexId getCellEdgeNumberInternal(const SimplexId &ttkNotUsed(cellId)) const
std::vector< std::vector< SimplexId > > cellNeighborList_
virtual SimplexId getEdgeLinkNumber(const SimplexId &edgeId) const
virtual const std::vector< std::vector< SimplexId > > * getVertexEdges()
virtual const std::vector< std::vector< SimplexId > > * getCellTrianglesInternal()
virtual bool isTriangleOnBoundaryInternal(const SimplexId &ttkNotUsed(triangleId)) const
virtual int getCellTriangleInternal(const SimplexId &ttkNotUsed(cellId), const int &ttkNotUsed(localTriangleId), SimplexId &ttkNotUsed(triangleId)) const
virtual int getVertexStar(const SimplexId &vertexId, const int &localStarId, SimplexId &starId) const
virtual SimplexId getTriangleStarNumber(const SimplexId &triangleId) const
std::vector< std::array< SimplexId, 2 > > edgeList_
virtual int getCellVertex(const SimplexId &cellId, const int &localVertexId, SimplexId &vertexId) const
virtual SimplexId getCellVertexNumber(const SimplexId &cellId) const
int getTriangleIncenter(const SimplexId triangleId, float incenter[3]) const
virtual const std::vector< std::vector< SimplexId > > * getTriangleStars()
virtual const std::vector< std::vector< SimplexId > > * getTriangleLinksInternal()
virtual SimplexId getVertexLinkNumberInternal(const SimplexId &ttkNotUsed(vertexId)) const
std::vector< std::array< SimplexId, 6 > > tetraEdgeList_
std::array< SimplexId, 3 > gridDimensions_
virtual SimplexId getVertexStarNumberInternal(const SimplexId &ttkNotUsed(vertexId)) const
std::vector< std::vector< SimplexId > > edgeStarList_
virtual const std::vector< std::vector< SimplexId > > * getVertexTrianglesInternal()
virtual SimplexId getNumberOfTriangles() const
size_t tableFootprint(const std::vector< itemType > &table, const std::string &tableName="", std::ostream &stream=std::cout) const
virtual int getTriangleVertex(const SimplexId &triangleId, const int &localVertexId, SimplexId &vertexId) const
virtual int getTriangleEdgeInternal(const SimplexId &ttkNotUsed(triangleId), const int &ttkNotUsed(localEdgeId), SimplexId &ttkNotUsed(edgeId)) const
virtual SimplexId getVertexStarNumber(const SimplexId &vertexId) const
virtual int getCellVTKID(const int &ttkId, int &vtkId) const
virtual SimplexId getCellTriangleNumberInternal(const SimplexId &ttkNotUsed(cellId)) const
virtual SimplexId getCellNeighborNumberInternal(const SimplexId &ttkNotUsed(cellId)) const
virtual SimplexId getEdgeTriangleNumberInternal(const SimplexId &ttkNotUsed(edgeId)) const
virtual int preconditionBoundaryTrianglesInternal()
virtual bool isEdgeOnBoundaryInternal(const SimplexId &ttkNotUsed(edgeId)) const
int getTriangleVertexNumber(const SimplexId ttkNotUsed(triangleId)) const
virtual int getEdgeTriangleInternal(const SimplexId &ttkNotUsed(edgeId), const int &ttkNotUsed(localTriangleId), SimplexId &ttkNotUsed(triangleId)) const
int getEdgeIncenter(const SimplexId edgeId, float incenter[3]) const
virtual int getCellTriangle(const SimplexId &cellId, const int &localTriangleId, SimplexId &triangleId) const
virtual const std::vector< std::vector< SimplexId > > * getTriangleEdges()
virtual const std::vector< std::vector< SimplexId > > * getEdgeStarsInternal()
virtual int getVertexEdgeInternal(const SimplexId &ttkNotUsed(vertexId), const int &ttkNotUsed(localEdgeId), SimplexId &ttkNotUsed(edgeId)) const
virtual SimplexId getVertexLinkNumber(const SimplexId &vertexId) const
virtual const std::vector< std::vector< SimplexId > > * getVertexLinksInternal()
virtual SimplexId getVertexEdgeNumber(const SimplexId &vertexId) const
virtual SimplexId getTriangleEdgeNumberInternal(const SimplexId &ttkNotUsed(triangleId)) const
std::array< std::vector< gradIdType >, 6 > gradientType
Discrete gradient internal struct.
std::vector< std::vector< SimplexId > > vertexTriangleList_
virtual SimplexId getNumberOfEdges() const
std::vector< std::vector< SimplexId > > cellEdgeVector_
size_t footprint(size_t size=0) const
std::vector< std::vector< SimplexId > > cellTriangleVector_
size_t tableTableFootprint(const std::vector< std::vector< itemType > > &table, const std::string &tableName="", std::ostream &stream=std::cout) const
virtual const std::vector< std::vector< SimplexId > > * getCellTriangles()
int getTetraIncenter(const SimplexId tetraId, float incenter[3]) const
virtual const std::vector< std::vector< SimplexId > > * getTriangleStarsInternal()
virtual SimplexId getTriangleEdgeNumber(const SimplexId &triangleId) const
virtual int preconditionBoundaryVerticesInternal()
virtual int getCellNeighbor(const SimplexId &cellId, const int &localNeighborId, SimplexId &neighborId) const
virtual const std::vector< std::vector< SimplexId > > * getVertexNeighborsInternal()
virtual int getVertexLink(const SimplexId &vertexId, const int &localLinkId, SimplexId &linkId) const
std::vector< std::array< SimplexId, 3 > > triangleEdgeList_
virtual int getVertexStarInternal(const SimplexId &ttkNotUsed(vertexId), const int &ttkNotUsed(localStarId), SimplexId &ttkNotUsed(starId)) const
virtual SimplexId getEdgeStarNumber(const SimplexId &edgeId) const
virtual int getCellVertexInternal(const SimplexId &ttkNotUsed(cellId), const int &ttkNotUsed(localVertexId), SimplexId &ttkNotUsed(vertexId)) const
std::vector< std::vector< SimplexId > > edgeTriangleList_
virtual const std::vector< std::vector< SimplexId > > * getEdgeLinksInternal()
virtual SimplexId getTriangleStarNumberInternal(const SimplexId &ttkNotUsed(triangleId)) const
virtual const std::vector< std::vector< SimplexId > > * getCellEdges()
AbstractTriangulation(AbstractTriangulation &&)=default
virtual int getDimensionalityInternal() const
int printMsg(const std::string &msg, const debug::Priority &priority=debug::Priority::INFO, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cout) const
Definition: Debug.h:118
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition: Debug.h:149
Wrapper class to wrap ttk code.
Definition: Wrapper.h:14
TTK discreteGradient processing package.
T distance(const T *p0, const T *p1, const int &dimension=3)
Definition: Geometry.cpp:344
The Topology ToolKit.
int SimplexId
Identifier type for simplices of any dimension.
Definition: DataTypes.h:22