TTK
Loading...
Searching...
No Matches
FTMTree_MT.cpp
Go to the documentation of this file.
1
13
14#include "FTMTree_MT.h"
15
16#include <stack>
17
18#define PRIOR(x)
19// #define PRIOR(x) priority(x)
20
21#ifdef __INTEL_COMPILER
22#define HIGHER
23#endif
24
25#ifndef TTK_ENABLE_OPENMP
26#define HIGHER
27#endif
28
29#ifdef TTK_ENABLE_OMP_PRIORITY
30#define OPTIONAL_PRIORITY(value) priority(value)
31#else
32#define OPTIONAL_PRIORITY(value)
33#endif
34
35using namespace std;
36using namespace ttk;
37using namespace ftm;
38
39FTMTree_MT::FTMTree_MT(const std::shared_ptr<Params> &params,
40 const std::shared_ptr<Scalars> &scalars,
41 TreeType type)
42 : params_(params), scalars_(scalars) {
43
44 this->setDebugMsgPrefix("FTMtree_MT");
45
46 mt_data_.treeType = type;
47}
48
50 this->clear();
51}
52
54
55 // remove UF data structures
56 if(!mt_data_.ufs.empty()) {
57 sort(mt_data_.ufs.begin(), mt_data_.ufs.end());
58 auto it = unique(mt_data_.ufs.begin(), mt_data_.ufs.end());
59 mt_data_.ufs.resize(std::distance(mt_data_.ufs.begin(), it));
60 }
61
62 // if (mt_data_.propagation) {
63 // Already cleaned by ufs
64 // sort(mt_data_.propagation->begin(), mt_data_.propagation->end());
65 // auto it = unique(mt_data_.propagation->begin(),
66 // mt_data_.propagation->end());
67 // mt_data_.propagation->resize(std::distance(mt_data_.propagation->begin(),
68 // it)); for (auto* addr : *mt_data_.propagation) if(addr) delete addr;
69 // }
70
71 // remove containers
73 mt_data_.superArcs.reset();
74 }
75 if(mt_data_.nodes) {
76 mt_data_.nodes.reset();
77 }
78 if(mt_data_.roots) {
79 mt_data_.roots.reset();
80 }
81 mt_data_.leaves.clear();
82 mt_data_.vert2tree.clear();
83 mt_data_.trunkSegments.clear();
84 mt_data_.visitOrder.clear();
85 mt_data_.ufs.clear();
86
87 if(mt_data_.states) {
88 mt_data_.states.reset();
89 }
90 mt_data_.propagation.clear();
91 mt_data_.valences.clear();
92 mt_data_.openedNodes.clear();
93
94#ifdef TTK_ENABLE_FTM_TREE_STATS_TIME
95 mt_data_.activeTasksStats.clear();
96#endif
97
98 this->params_.reset();
99 this->scalars_.reset();
100}
101
103
104 const idSuperArc nbArcs = mt_data_.superArcs->size();
105
106 // Make reserve
107
108 // SuperArc i correspond to segment i,
109 // one arc correspond to one segment
110 vector<SimplexId> sizes(nbArcs);
111
112 // get the size of each segment
113 const idSuperArc arcChunkSize = getChunkSize(nbArcs);
114 const idSuperArc arcChunkNb = getChunkCount(nbArcs);
115 for(idSuperArc arcChunkId = 0; arcChunkId < arcChunkNb; ++arcChunkId) {
116 // WHY shared(sizes) is needed ??
117#ifdef TTK_ENABLE_OPENMP
118#pragma omp task firstprivate(arcChunkId) shared(sizes) \
119 OPTIONAL_PRIORITY(isPrior())
120#endif
121 {
122 const idSuperArc lowerBound = arcChunkId * arcChunkSize;
123 const idSuperArc upperBound
124 = min(nbArcs, (arcChunkId + 1) * arcChunkSize);
125 for(idSuperArc a = lowerBound; a < upperBound; ++a) {
126 sizes[a]
127 = max(SimplexId{0}, (*mt_data_.superArcs)[a].getNbVertSeen() - 1);
128 }
129 }
130 }
131#ifdef TTK_ENABLE_OPENMP
132#pragma omp taskwait
133#endif
134
135 // change segments size using the created vector
137
138 Timer segmentsSet;
139
140 // Fill segments using vert2tree
141
142 // current status of the segmentation of this arc
143 vector<SimplexId> posSegm(nbArcs, 0);
144
145 // Segments are connex region of geometry forming
146 // the segmentation (sorted in ascending order)
147 const SimplexId nbVert = scalars_->size;
148 const SimplexId chunkSize = getChunkSize();
149 const SimplexId chunkNb = getChunkCount();
150 for(SimplexId chunkId = 0; chunkId < chunkNb; ++chunkId) {
151#ifdef TTK_ENABLE_OPENMP
152#pragma omp task firstprivate(chunkId) shared(posSegm) \
153 OPTIONAL_PRIORITY(isPrior())
154#endif
155 {
156 const SimplexId lowerBound = chunkId * chunkSize;
157 const SimplexId upperBound = min(nbVert, (chunkId + 1) * chunkSize);
158 for(SimplexId i = lowerBound; i < upperBound; ++i) {
159 const auto vert = scalars_->sortedVertices[i];
160 if(isCorrespondingArc(vert)) {
161 idSuperArc const sa = getCorrespondingSuperArcId(vert);
162 SimplexId vertToAdd;
163 if(mt_data_.visitOrder[vert] != nullVertex) {
164 // Opposite order for Split Tree
165 vertToAdd = mt_data_.visitOrder[vert];
166 if(isST())
167 vertToAdd = getSuperArc(sa)->getNbVertSeen() - vertToAdd - 2;
168 mt_data_.segments_[sa][vertToAdd] = vert;
169 } else if(mt_data_.trunkSegments.empty()) {
170 // MT computation
171#ifdef TTK_ENABLE_OPENMP
172#pragma omp atomic capture
173#endif
174 vertToAdd = posSegm[sa]++;
175 mt_data_.segments_[sa][vertToAdd] = vert;
176 }
177
178 } // end is arc
179 } // end for
180 } // end task
181 }
182#ifdef TTK_ENABLE_OPENMP
183#pragma omp taskwait
184#endif
185
186 printTime(segmentsSet, "segmentation set vertices", 4);
187
188 if(mt_data_.trunkSegments.empty()) {
189 // sort arc that have been filled by the trunk
190 // only for MT
191 Timer segmentsSortTime;
192 for(idSuperArc a = 0; a < nbArcs; ++a) {
193 if(posSegm[a]) {
194#ifdef TTK_ENABLE_OPENMP
195#pragma omp task firstprivate(a) OPTIONAL_PRIORITY(isPrior())
196#endif
197 mt_data_.segments_[a].sort(scalars_.get());
198 }
199 }
200#ifdef TTK_ENABLE_OPENMP
201#pragma omp taskwait
202#endif
203 printTime(segmentsSortTime, "segmentation sort vertices", 4);
204 } else {
205 // Contour tree: we create the arc segmentation for arcs in the trunk
206 Timer segmentsArcTime;
207 for(idSuperArc a = 0; a < nbArcs; ++a) {
208 // CT computation, we have already the vert list
209 if(!mt_data_.trunkSegments[a].empty()) {
210#ifdef TTK_ENABLE_OPENMP
211#pragma omp task firstprivate(a) OPTIONAL_PRIORITY(isPrior())
212#endif
213 mt_data_.segments_[a].createFromList(
216 }
217 }
218#ifdef TTK_ENABLE_OPENMP
219#pragma omp taskwait
220#endif
221
222 printTime(segmentsArcTime, "segmentation arcs lists", 4);
223 }
224
225 // Update SuperArc region
226
227 // ST have a segmentation which is in the reverse-order of its build
228 // ST have a segmentation sorted in ascending order as JT
229 for(idSuperArc arcChunkId = 0; arcChunkId < arcChunkNb; ++arcChunkId) {
230#ifdef TTK_ENABLE_OPENMP
231#pragma omp task firstprivate(arcChunkId) OPTIONAL_PRIORITY(isPrior())
232#endif
233 {
234 const idSuperArc lowerBound = arcChunkId * arcChunkSize;
235 const idSuperArc upperBound
236 = min(nbArcs, (arcChunkId + 1) * arcChunkSize);
237 for(idSuperArc a = lowerBound; a < upperBound; ++a) {
238 // avoid empty region
239 if(mt_data_.segments_[a].size()) {
240 (*mt_data_.superArcs)[a].concat(
241 mt_data_.segments_[a].begin(), mt_data_.segments_[a].end());
242 }
243 }
244 }
245 }
246#ifdef TTK_ENABLE_OPENMP
247#pragma omp taskwait
248#endif
249}
250
251std::shared_ptr<FTMTree_MT> FTMTree_MT::clone() const {
252 auto newMT
253 = std::make_shared<FTMTree_MT>(params_, scalars_, mt_data_.treeType);
254
255 newMT->mt_data_.superArcs = mt_data_.superArcs;
256 newMT->mt_data_.nodes = mt_data_.nodes;
257 newMT->mt_data_.leaves = mt_data_.leaves;
258 newMT->mt_data_.roots = mt_data_.roots;
259 newMT->mt_data_.vert2tree = mt_data_.vert2tree;
260
261 return newMT;
262}
263
264void FTMTree_MT::closeArcsUF(idNode closeNode, UF uf) {
265 for(const auto &sa : uf->find()->getOpenedArcs()) {
266 closeSuperArc(sa, closeNode);
267 }
268 uf->find()->clearOpenedArcs();
269}
270
271void FTMTree_MT::closeSuperArc(idSuperArc superArcId, idNode upNodeId) {
272#ifndef TTK_ENABLE_KAMIKAZE
273
274 if(superArcId >= getNumberOfSuperArcs()) {
275 cout << "[Merge Tree] closeSuperArc on a inexisting arc !" << endl;
276 return;
277 }
278
279 if(upNodeId >= getNumberOfNodes()) {
280 cout << "[Merge Tree] closeOpenedArc on a inexisting node !" << endl;
281 return;
282 }
283
284#endif
285 (*mt_data_.superArcs)[superArcId].setUpNodeId(upNodeId);
286 (*mt_data_.nodes)[upNodeId].addDownSuperArcId(superArcId);
287}
288
290 Node *mainNode = getNode(node);
291
292 if(mainNode->getNumberOfUpSuperArcs() == 0) {
293
294 // Root: No Superarc
295#ifndef TTK_ENABLE_KAMIKAZE
296 if(mainNode->getNumberOfDownSuperArcs() != 1) {
297 // Root with several children: impossible /\ .
298 cout << endl << "[FTMTree_MT]:delNode won't delete ";
299 cout << mainNode->getVertexId() << " (root) with ";
300 cout << static_cast<unsigned>(mainNode->getNumberOfDownSuperArcs())
301 << " down ";
302 cout << static_cast<unsigned>(mainNode->getNumberOfUpSuperArcs())
303 << " up ";
304 return;
305 }
306#endif
307
308 idSuperArc const downArc = mainNode->getDownSuperArcId(0);
309 Node *downNode = getNode((*mt_data_.superArcs)[downArc].getDownNodeId());
310
311 downNode->removeUpSuperArc(downArc);
312 mainNode->clearDownSuperArcs();
313
314 } else if(mainNode->getNumberOfDownSuperArcs() < 2) {
315 // Have one up arc
316
317 // We delete the upArc of this node,
318 // if there is a down arc, we reattach it to the upNode
319
320 idSuperArc const upArc = mainNode->getUpSuperArcId(0);
321 idNode const upId = (*mt_data_.superArcs)[upArc].getUpNodeId();
322 Node *upNode = getNode(upId);
323
324 upNode->removeDownSuperArc(upArc);
325 mainNode->clearUpSuperArcs();
326
327 if(mainNode->getNumberOfDownSuperArcs()) {
328 // Have one down arc
329
330 // Reconnect
331 idSuperArc const downArc = mainNode->getDownSuperArcId(0);
332 (*mt_data_.superArcs)[downArc].setUpNodeId(upId);
333 upNode->addDownSuperArcId(downArc);
334 mainNode->clearDownSuperArcs();
335
336 // Segmentation
337 (*mt_data_.superArcs)[downArc].concat((*mt_data_.superArcs)[upArc]);
338 }
339 }
340#ifndef TTK_ENABLE_KAMIKAZE
341 else
342 cerr << "delete node with multiple childrens " << endl;
343#endif
344}
345
347 for(auto &arc : *mt_data_.superArcs) {
348 arc.createSegmentation(scalars_.get());
349 }
350}
351
352tuple<SimplexId, SimplexId>
353 FTMTree_MT::getBoundsFromVerts(const vector<SimplexId> &trunkVerts) const {
354 SimplexId begin, stop;
355
356 if(isST()) {
357 begin = 0;
358 stop = scalars_->offsets[trunkVerts[0]];
359 } else {
360 begin = scalars_->offsets[trunkVerts[0]];
361 stop = scalars_->size;
362 }
363
364 return make_tuple(begin, stop);
365}
366
368 return &((*mt_data_.nodes)[a->getDownNodeId()]);
369}
370
372 return a->getDownNodeId();
373}
374
376 if(isST())
377 return getUpNode(a);
378
379 return getDownNode(a);
380}
381
383 if(isST())
384 return getUpNodeId(a);
385
386 return getDownNodeId(a);
387}
388
390 return &((*mt_data_.nodes)[a->getUpNodeId()]);
391}
392
394 return a->getUpNodeId();
395}
396
398 if(isST())
399 return getDownNode(a);
400
401 return getUpNode(a);
402}
403
405 if(isST())
406 return getDownNodeId(a);
407
408 return getUpNodeId(a);
409}
410
411idNode FTMTree_MT::getVertInRange(const vector<SimplexId> &range,
412 const SimplexId v,
413 const idNode last) const {
414 idNode idRes = last;
415 const idNode rangeSize = range.size();
416 while(idRes + 1 < rangeSize && comp_.vertLower(range[idRes + 1], v)) {
417 ++idRes;
418 }
419 return idRes;
420}
421
422idSuperArc FTMTree_MT::insertNode(Node *node, const bool segm) {
423 // Normal insert : existing arc stay below inserted (JT example)
424 // * - <- upNodeId
425 // | \ | <- newSA
426 // | * <- newNodeId
427 // | | <- currentSA
428 // - - -
429 // already present
430 if(isCorrespondingNode(node->getVertexId())) {
431 Node *myNode = vertex2Node(node->getVertexId());
432 // If it has been hidden / replaced we need to re-make it
433 idSuperArc const correspondingArcId = myNode->getUpSuperArcId(0);
434 updateCorrespondingArc(myNode->getVertexId(), correspondingArcId);
435 }
436
437 idNode upNodeId, newNodeId;
438 idSuperArc currentSA, newSA;
439 SimplexId origin;
440
441 // Create new node
442 currentSA = getCorrespondingSuperArcId(node->getVertexId());
443 upNodeId = (*mt_data_.superArcs)[currentSA].getUpNodeId();
444 origin = (*mt_data_.nodes)[(*mt_data_.superArcs)[currentSA].getDownNodeId()]
445 .getOrigin();
446 newNodeId = makeNode(node, origin);
447
448 // Connectivity
449 // Insert only node inside the partition : created arc don t cross
450 newSA = makeSuperArc(newNodeId, upNodeId);
451
452 (*mt_data_.superArcs)[currentSA].setUpNodeId(newNodeId);
453 (*mt_data_.nodes)[upNodeId].removeDownSuperArc(currentSA);
454 (*mt_data_.nodes)[newNodeId].addDownSuperArcId(currentSA);
455
456 // cut the vertex list at the node position and
457 // give each arc its part.
458 if(segm) {
460 (*mt_data_.superArcs)[newSA].concat(
461 get<1>((*mt_data_.superArcs)[currentSA].splitBack(
462 node->getVertexId(), scalars_.get())));
463 } else {
464 (*mt_data_.superArcs)[newSA].concat(
465 get<1>((*mt_data_.superArcs)[currentSA].splitFront(
466 node->getVertexId(), scalars_.get())));
467 }
468 }
469
470 return newSA;
471}
472
474#ifndef TTK_ENABLE_KAMIKAZE
475 if(vertexId < 0 || vertexId >= scalars_->size) {
476 this->printMsg({{"make node, wrong vertex :", std::to_string(vertexId)},
477 {" on ", std::to_string(scalars_->size)}},
479 return -1;
480 }
481#endif
482
483 if(isCorrespondingNode(vertexId)) {
484 return getCorrespondingNodeId(vertexId);
485 }
486
487 idNode const newNodeId = mt_data_.nodes->getNext();
488 (*mt_data_.nodes)[newNodeId].setVertexId(vertexId);
489 (*mt_data_.nodes)[newNodeId].setTermination(term);
490 updateCorrespondingNode(vertexId, newNodeId);
491
492 return newNodeId;
493}
494
496 return makeNode(n->getVertexId());
497}
498
500
501{
502 idSuperArc const newSuperArcId = mt_data_.superArcs->getNext();
503 (*mt_data_.superArcs)[newSuperArcId].setDownNodeId(downNodeId);
504 (*mt_data_.superArcs)[newSuperArcId].setUpNodeId(upNodeId);
505
506 (*mt_data_.nodes)[downNodeId].addUpSuperArcId(newSuperArcId);
507 (*mt_data_.nodes)[upNodeId].addDownSuperArcId(newSuperArcId);
508
509 return newSuperArcId;
510}
511
513 // we already have common data
515 mt.mt_data_.superArcs.reset();
517 mt.mt_data_.nodes.reset();
518 mt_data_.leaves = std::move(mt.mt_data_.leaves);
520 mt.mt_data_.roots.reset();
521 mt_data_.vert2tree = std::move(mt.mt_data_.vert2tree);
522}
523
525 Timer normTime;
526 sortLeaves(true);
527 if(this->params_->treeType != TreeType::Contour) {
528 sortNodes();
529 sortArcs();
530 }
531
532 auto getNodeParentArcNb
533 = [&](const idNode curNode, const bool goUp) -> idSuperArc {
534 if(goUp) {
535 return getNode(curNode)->getNumberOfUpSuperArcs();
536 }
537
538 return getNode(curNode)->getNumberOfDownSuperArcs();
539 };
540
541 auto getNodeParentArc
542 = [&](const idNode curNode, const bool goUp, idSuperArc i) -> idSuperArc {
543 if(goUp) {
544 return getNode(curNode)->getUpSuperArcId(i);
545 }
546
547 return getNode(curNode)->getDownSuperArcId(i);
548 };
549
550 auto getArcParentNode
551 = [&](const idSuperArc curArc, const bool goUp) -> idNode {
552 if(goUp) {
553 return getSuperArc(curArc)->getUpNodeId();
554 }
555
556 return getSuperArc(curArc)->getDownNodeId();
557 };
558
559 std::queue<tuple<idNode, bool>> q;
560 std::stack<tuple<idNode, bool>> qr;
561 for(const idNode n : mt_data_.leaves) {
562 bool const goUp = isJT() || isST() || getNode(n)->getNumberOfUpSuperArcs();
563 if(goUp)
564 q.emplace(n, goUp);
565 else
566 qr.emplace(n, goUp);
567 }
568
569 while(!qr.empty()) {
570 q.emplace(qr.top());
571 qr.pop();
572 }
573
574 // Normalized id
575 idSuperArc nIdMin = 0;
576 idSuperArc nIdMax = getNumberOfSuperArcs() - 1;
577
578 vector<bool> seenUp(getNumberOfSuperArcs(), false);
579 vector<bool> seenDown(getNumberOfSuperArcs(), false);
580
581 while(!q.empty()) {
582 bool goUp;
583 idNode curNodeId;
584 tie(curNodeId, goUp) = q.front();
585 q.pop();
586
587 if(goUp)
588 sortUpArcs(curNodeId);
589 else
590 sortDownArcs(curNodeId);
591
592 // Assign arc above
593 const idSuperArc nbArcParent = getNodeParentArcNb(curNodeId, goUp);
594 for(idSuperArc pid = 0; pid < nbArcParent; pid++) {
595 const idSuperArc currentArcId = getNodeParentArc(curNodeId, goUp, pid);
596 if(goUp) {
597 if(getSuperArc(currentArcId)->getNormalizedId() == nullSuperArc) {
598 getSuperArc(currentArcId)->setNormalizeIds(nIdMin++);
599 }
600 if(!seenUp[currentArcId]) {
601 q.emplace(getArcParentNode(currentArcId, goUp), goUp);
602 seenUp[currentArcId] = true;
603 }
604 } else {
605 if(getSuperArc(currentArcId)->getNormalizedId() == nullSuperArc) {
606 getSuperArc(currentArcId)->setNormalizeIds(nIdMax--);
607 }
608 if(!seenDown[currentArcId]) {
609 q.emplace(getArcParentNode(currentArcId, goUp), goUp);
610 seenDown[currentArcId] = true;
611 }
612 }
613 }
614 }
615
616#ifndef TTK_ENABLE_KAMIKAZE
617 if(std::abs((long)nIdMax - (long)nIdMin) > 1) {
618 this->printMsg({"error during normalize, tree compromised: ",
619 std::to_string(nIdMin), " ", std::to_string(nIdMax)},
621 }
622#endif
623
624 printTime(normTime, "normalize ids", 4);
625}
626
628#ifndef TTK_ENABLE_KAMIKAZE
629 if(downNodeId >= getNumberOfNodes()) {
630 this->printErr("openSuperArc on a inexisting node !");
631 return -2;
632 }
633#endif
634
635 idSuperArc const newSuperArcId = mt_data_.superArcs->getNext();
636 (*mt_data_.superArcs)[newSuperArcId].setDownNodeId(downNodeId);
637 (*mt_data_.nodes)[downNodeId].addUpSuperArcId(newSuperArcId);
638
639 return newSuperArcId;
640}
641
643 const SuperArc *sa = getSuperArc(a);
644 stringstream res;
645 const SimplexId dv = getNode(sa->getDownNodeId())->getVertexId();
646 const SimplexId uv = getNode(sa->getUpNodeId())->getVertexId();
647 res << a;
648 res << " : ";
649 if(dv != nullVertex) {
650 res << dv << " -- ";
651 } else {
652 res << "XX -- ";
653 }
654 if(uv != nullVertex) {
655 res << uv;
656 } else {
657 res << "XX";
658 }
659
660 res.seekg(0, ios::end);
661 while(res.tellg() < 25) {
662 res << " ";
663 res.seekg(0, ios::end);
664 }
665 res.seekg(0, ios::beg);
666
667 res << "segm #" << sa->regionSize() << " / " << scalars_->size; // << " -> ";
668
669 res.seekg(0, ios::end);
670
671 while(res.tellg() < 45) {
672 res << " ";
673 res.seekg(0, ios::end);
674 }
675 res.seekg(0, ios::beg);
676
677 res << sa->printReg();
678 return res.str();
679}
680
682 const Node *node = getNode(n);
683 stringstream res;
684 res << n;
685 res << " : (";
686 res << node->getVertexId() << ") \\ ";
687
688 for(idSuperArc i = 0; i < node->getNumberOfDownSuperArcs(); ++i) {
689 res << "+";
690 res << node->getDownSuperArcId(i) << " ";
691 }
692
693 res << " / ";
694
695 for(idSuperArc i = 0; i < node->getNumberOfUpSuperArcs(); ++i) {
696 res << "+";
697 res << node->getUpSuperArcId(i) << " ";
698 }
699
700 return res.str();
701}
702
704 if(debugLevel_ > 1) {
705 if(debugLevel_ > 2) {
707 }
708 this->printMsg("number of threads : " + std::to_string(threadNumber_));
709 if(debugLevel_ > 2) {
710 this->printMsg("* debug lvl : " + std::to_string(debugLevel_));
711 string tt;
712 if(params_->treeType == TreeType::Contour) {
713 tt = "Contour";
714 } else if(params_->treeType == TreeType::Join) {
715 tt = "Join";
716 } else if(params_->treeType == TreeType::Split) {
717 tt = "Split";
718 } else if(params_->treeType == TreeType::Join_Split) {
719 tt = "Join + Split";
720 }
721 this->printMsg("* tree type : " + tt);
723 }
724 }
725}
726
728 const string &s,
729 const int debugLevel) const {
730
731 if(this->debugLevel_ >= debugLevel) {
732 stringstream st;
733
734 for(int i = 3; i < debugLevel; i++)
735 st << "-";
736 st << s;
737
738 this->printMsg(st.str(), 1, t.getElapsedTime(), this->threadNumber_);
739 }
740 return 1;
741}
742
744#ifdef TTK_ENABLE_OPENMP
745#pragma omp critical
746#endif
747 {
748 cout << "Nodes----------" << endl;
749 for(idNode nid = 0; nid < getNumberOfNodes(); nid++) {
750 cout << printNode(nid) << endl;
751 }
752
753 cout << "Arcs-----------" << endl;
754 for(idSuperArc said = 0; said < getNumberOfSuperArcs(); ++said) {
755 cout << printArc(said) << endl;
756 }
757
758 cout << "Leaves" << endl;
759 for(const auto &l : mt_data_.leaves)
760 cout << " " << (*mt_data_.nodes)[l].getVertexId();
761 cout << endl;
762
763 cout << "Roots" << endl;
764 for(const auto &r : *mt_data_.roots)
765 cout << " " << (*mt_data_.nodes)[r].getVertexId();
766 cout << endl;
767 }
768}
769
770void FTMTree_MT::sortLeaves(const bool para) {
771 auto indirect_sort = [&](const idNode a, const idNode b) {
772 return comp_.vertLower(
774 };
775
776 if(para) {
777 TTK_PSORT(this->threadNumber_, mt_data_.leaves.begin(),
778 mt_data_.leaves.end(), indirect_sort);
779 } else {
780 std::sort(mt_data_.leaves.begin(), mt_data_.leaves.end(), indirect_sort);
781 }
782}
783
785 std::vector<idNode> sortedNodes(this->mt_data_.nodes->size());
786 std::iota(sortedNodes.begin(), sortedNodes.end(), 0);
787
788 const auto direct_sort = [&](const Node &a, const Node &b) {
789 // sort according to scalar field
790 return this->comp_.vertLower(a.getVertexId(), b.getVertexId());
791 };
792
793 const auto indirect_sort = [&](const idNode a, const idNode b) {
794 return direct_sort(*this->getNode(a), *this->getNode(b));
795 };
796
797 TTK_PSORT(
798 this->threadNumber_, sortedNodes.begin(), sortedNodes.end(), indirect_sort);
799
800 TTK_PSORT(this->threadNumber_, this->mt_data_.nodes->begin(),
801 this->mt_data_.nodes->end(), direct_sort);
802
803 // reverse sortedNodes
804 std::vector<idNode> revSortedNodes(sortedNodes.size());
805 for(size_t i = 0; i < sortedNodes.size(); ++i) {
806 revSortedNodes[sortedNodes[i]] = i;
807 }
808
809 // update leaves
810 for(auto &leaf : this->mt_data_.leaves) {
811 leaf = revSortedNodes[leaf];
812 }
813
814 // update roots
815 for(auto &root : (*this->mt_data_.roots)) {
816 root = revSortedNodes[root];
817 }
818
819 // update arcs
820 for(auto &arc : (*this->mt_data_.superArcs)) {
821 arc.setDownNodeId(revSortedNodes[arc.getDownNodeId()]);
822 arc.setUpNodeId(revSortedNodes[arc.getUpNodeId()]);
823 }
824
825 // update vert2tree
826 for(size_t i = 0; i < sortedNodes.size(); ++i) {
827 const auto &node{(*this->mt_data_.nodes)[i]};
828 if(this->isCorrespondingNode(node.getVertexId())) {
829 this->updateCorrespondingNode(node.getVertexId(), i);
830 }
831 }
832}
833
835 std::vector<idNode> sortedArcs(this->mt_data_.superArcs->size());
836 std::iota(sortedArcs.begin(), sortedArcs.end(), 0);
837
838 const auto direct_sort = [&](const SuperArc &a, const SuperArc &b) {
839 // sort by NodeId (nodes should be already sorted with sortNodes)
840 const auto adn{a.getDownNodeId()};
841 const auto aun{a.getUpNodeId()};
842 const auto bdn{b.getDownNodeId()};
843 const auto bun{b.getUpNodeId()};
844 return std::tie(adn, aun) < std::tie(bdn, bun);
845 };
846
847 auto indirect_sort = [&](const idSuperArc &a, const idSuperArc &b) {
848 const auto aa{this->getSuperArc(a)};
849 const auto bb{this->getSuperArc(b)};
850 return direct_sort(*aa, *bb);
851 };
852
853 TTK_PSORT(
854 this->threadNumber_, sortedArcs.begin(), sortedArcs.end(), indirect_sort);
855
856 TTK_PSORT(this->threadNumber_, this->mt_data_.superArcs->begin(),
857 this->mt_data_.superArcs->end(), direct_sort);
858
859 // reverse sortedArcs
860 std::vector<idSuperArc> revSortedArcs(sortedArcs.size());
861 for(size_t i = 0; i < sortedArcs.size(); ++i) {
862 revSortedArcs[sortedArcs[i]] = i;
863 }
864
865 // update nodes
866 std::vector<idSuperArc> updatedArcs{};
867 for(auto &node : (*this->mt_data_.nodes)) {
868 {
869 const auto da{node.getNumberOfDownSuperArcs()};
870 updatedArcs.clear();
871 updatedArcs.resize(da);
872 for(idSuperArc i = 0; i < da; ++i) {
873 updatedArcs[i] = revSortedArcs[node.getDownSuperArcId(i)];
874 }
875 node.clearDownSuperArcs();
876 for(const auto &arc : updatedArcs) {
877 node.addDownSuperArcId(arc);
878 }
879 }
880 {
881 const auto ua{node.getNumberOfUpSuperArcs()};
882 updatedArcs.clear();
883 updatedArcs.resize(ua);
884 for(idSuperArc i = 0; i < ua; ++i) {
885 updatedArcs[i] = revSortedArcs[node.getUpSuperArcId(i)];
886 }
887 node.clearUpSuperArcs();
888 for(const auto &arc : updatedArcs) {
889 node.addUpSuperArcId(arc);
890 }
891 }
892 }
893
894 // update vert2tree
895 for(size_t i = 0; i < this->mt_data_.vert2tree.size(); ++i) {
896 if(this->isCorrespondingArc(i)) {
897 this->updateCorrespondingArc(
898 i, revSortedArcs[this->getCorrespondingSuperArcId(i)]);
899 }
900 }
901}
902
903vector<idNode> FTMTree_MT::sortedNodes(const bool para) {
904 vector<idNode> sortedNodes(mt_data_.nodes->size());
905 std::iota(sortedNodes.begin(), sortedNodes.end(), 0);
906
907 auto indirect_sort = [&](const idNode a, const idNode b) {
908 return comp_.vertLower(
910 };
911
912 if(para) {
913 TTK_PSORT(this->threadNumber_, sortedNodes.begin(), sortedNodes.end(),
914 indirect_sort);
915 } else {
916#ifdef TTK_ENABLE_OPENMP
917#pragma omp single
918#endif
919 { std::sort(sortedNodes.begin(), sortedNodes.end(), indirect_sort); }
920 }
921
922 return sortedNodes;
923}
924
925SimplexId FTMTree_MT::trunkCTSegmentation(const vector<SimplexId> &trunkVerts,
926 const SimplexId begin,
927 const SimplexId stop) {
928 const int nbTasksThreads = 40;
929 const auto sizeBackBone = abs(stop - begin);
930 const auto chunkSize = getChunkSize(sizeBackBone, nbTasksThreads);
931 const auto chunkNb = getChunkCount(sizeBackBone, nbTasksThreads);
932 // si pas efficace vecteur de la taille de node ici a la place de acc
933 idNode lastVertInRange = 0;
935 for(SimplexId chunkId = 0; chunkId < chunkNb; ++chunkId) {
936#ifdef TTK_ENABLE_OPENMP
937#pragma omp task firstprivate(chunkId, lastVertInRange) shared(trunkVerts) \
938 OPTIONAL_PRIORITY(isPrior())
939#endif
940 {
941 vector<SimplexId> regularList;
942 if(params_->segm) {
943 regularList.reserve(25);
944 }
945 const SimplexId lowerBound = begin + chunkId * chunkSize;
946 const SimplexId upperBound
947 = min(stop, (begin + (chunkId + 1) * chunkSize));
948 if(lowerBound != upperBound) {
949 const SimplexId pos = isST() ? upperBound - 1 : lowerBound;
950 lastVertInRange
951 = getVertInRange(trunkVerts, scalars_->sortedVertices[pos], 0);
952 }
953 for(SimplexId v = lowerBound; v < upperBound; ++v) {
954 const SimplexId s
955 = isST() ? scalars_->sortedVertices[lowerBound + upperBound - 1 - v]
956 : scalars_->sortedVertices[v];
957 if(isCorrespondingNull(s)) {
958 const idNode oldVertInRange = lastVertInRange;
959 lastVertInRange = getVertInRange(trunkVerts, s, lastVertInRange);
960 const idSuperArc thisArc = upArcFromVert(trunkVerts[lastVertInRange]);
961 updateCorrespondingArc(s, thisArc);
962
963 if(params_->segm) {
964 if(oldVertInRange == lastVertInRange) {
965 regularList.emplace_back(s);
966 } else {
967 // accumulated to have only one atomic update when needed
968 const idSuperArc oldArc
969 = upArcFromVert(trunkVerts[oldVertInRange]);
970 if(regularList.size()) {
971#ifdef TTK_ENABLE_OPENMP
972#pragma omp critical
973#endif
974 {
975 mt_data_.trunkSegments[oldArc].emplace_back(regularList);
976 regularList.clear();
977 }
978 }
979 // hand.vtu, sequential: 28554
980 regularList.emplace_back(s);
981 }
982 }
983 }
984 }
985 // force increment last arc
986 const idNode baseNode
987 = getCorrespondingNodeId(trunkVerts[lastVertInRange]);
988 const idSuperArc upArc = getNode(baseNode)->getUpSuperArcId(0);
989 if(regularList.size()) {
990#ifdef TTK_ENABLE_OPENMP
991#pragma omp critical
992#endif
993 {
994 mt_data_.trunkSegments[upArc].emplace_back(regularList);
995 regularList.clear();
996 }
997 }
998 }
999 }
1000#ifdef TTK_ENABLE_OPENMP
1001#pragma omp taskwait
1002#endif
1003 return {};
1004}
1005
1006SimplexId FTMTree_MT::trunkSegmentation(const vector<SimplexId> &trunkVerts,
1007 const SimplexId begin,
1008 const SimplexId stop) {
1009 // Assign missing vert to the good arc
1010 // and also add the corresponding number for
1011 // futur arc reserve
1012 const int nbTasksThreads = 40;
1013 const auto sizeBackBone = abs(stop - begin);
1014 const auto chunkSize = getChunkSize(sizeBackBone, nbTasksThreads);
1015 const auto chunkNb = getChunkCount(sizeBackBone, nbTasksThreads);
1016 // si pas efficace vecteur de la taille de node ici a la place de acc
1017 for(SimplexId chunkId = 0; chunkId < chunkNb; ++chunkId) {
1018#ifdef TTK_ENABLE_OPENMP
1019#pragma omp task firstprivate(chunkId) shared(trunkVerts) \
1020 OPTIONAL_PRIORITY(isPrior())
1021#endif
1022 {
1023 idNode lastVertInRange = 0;
1024 SimplexId acc = 0;
1025
1026 const SimplexId lowerBound = begin + chunkId * chunkSize;
1027 const SimplexId upperBound
1028 = min(stop, (begin + (chunkId + 1) * chunkSize));
1029 for(SimplexId v = lowerBound; v < upperBound; ++v) {
1030 const SimplexId s
1031 = isST() ? scalars_->sortedVertices[lowerBound + upperBound - 1 - v]
1032 : scalars_->sortedVertices[v];
1033 if(isCorrespondingNull(s)) {
1034 const idNode oldVertInRange = lastVertInRange;
1035 lastVertInRange = getVertInRange(trunkVerts, s, lastVertInRange);
1036 const idSuperArc thisArc = upArcFromVert(trunkVerts[lastVertInRange]);
1037 updateCorrespondingArc(s, thisArc);
1038
1039 if(params_->segm) {
1040 if(oldVertInRange == lastVertInRange) {
1041 ++acc;
1042 } else {
1043 // accumulated to have only one atomic update when needed
1044 const idSuperArc oldArc
1045 = upArcFromVert(trunkVerts[oldVertInRange]);
1046 getSuperArc(oldArc)->atomicIncVisited(acc);
1047 acc = 1;
1048 }
1049 }
1050 }
1051 }
1052 // force increment last arc
1053 const idNode baseNode
1054 = getCorrespondingNodeId(trunkVerts[lastVertInRange]);
1055 const idSuperArc upArc = getNode(baseNode)->getUpSuperArcId(0);
1056 getSuperArc(upArc)->atomicIncVisited(acc);
1057 } // end task
1058 }
1059#ifdef TTK_ENABLE_OPENMP
1060#pragma omp taskwait
1061#endif
1062 return {};
1063}
1064
1065std::ostream &ttk::ftm::operator<<(std::ostream &o,
1066 ttk::ftm::SuperArc const &a) {
1067 o << a.getDownNodeId() << " <>> " << a.getUpNodeId();
1068 return o;
1069}
1070
1071std::ostream &ttk::ftm::operator<<(std::ostream &o, ttk::ftm::Node const &n) {
1072 o << n.getNumberOfDownSuperArcs() << " .-. " << n.getNumberOfUpSuperArcs();
1073 return o;
1074}
#define TTK_PSORT(NTHREADS,...)
Parallel sort macro.
Definition OpenMP.h:46
int debugLevel_
Definition Debug.h:379
void setDebugMsgPrefix(const std::string &prefix)
Definition Debug.h:364
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition Debug.h:149
double getElapsedTime()
Definition Timer.h:15
FTMAtomicVector< idSuperArc > & getOpenedArcs()
Definition FTMAtomicUF.h:86
AtomicUF * find()
Definition FTMAtomicUF.h:38
Node * getUpperNode(const SuperArc *a)
bool isCorrespondingNode(const SimplexId val) const
Definition FTMTree_MT.h:449
std::shared_ptr< Params > params_
Definition FTMTree_MT.h:86
std::string printNode(idNode n)
SuperArc * getSuperArc(idSuperArc i)
Definition FTMTree_MT.h:367
Node * getDownNode(const SuperArc *a)
void printParams() const
idNode getNumberOfNodes() const
Definition FTMTree_MT.h:389
idNode getUpNodeId(const SuperArc *a)
SimplexId getChunkSize(const SimplexId nbVerts=-1, const SimplexId nbtasks=100) const
Definition FTMTree_MT.h:821
bool isCorrespondingNull(const SimplexId val) const
Definition FTMTree_MT.h:453
void delNode(idNode node)
std::string printArc(idSuperArc a)
Node * vertex2Node(const SimplexId vert)
Definition FTMTree_MT.h:487
void move(FTMTree_MT &mt)
idSuperArc upArcFromVert(const SimplexId v)
Definition FTMTree_MT.h:817
std::tuple< SimplexId, SimplexId > getBoundsFromVerts(const std::vector< SimplexId > &nodes) const
Node * getLowerNode(const SuperArc *a)
idNode getDownNodeId(const SuperArc *a)
idSuperArc insertNode(Node *node, const bool segm=true)
bool isCorrespondingArc(const SimplexId val) const
Definition FTMTree_MT.h:445
void closeArcsUF(idNode closeNode, UF uf)
idNode getLowerNodeId(const SuperArc *a)
void sortUpArcs(const idNode nid)
Definition FTMTree_MT.h:840
virtual SimplexId trunkSegmentation(const std::vector< SimplexId > &pendingNodesVerts, const SimplexId begin, const SimplexId stop)
void sortNodes()
Sort tree nodes according to vertex order.
idNode getCorrespondingNodeId(const SimplexId val) const
Definition FTMTree_MT.h:459
idSuperArc getNumberOfSuperArcs() const
Definition FTMTree_MT.h:363
idNode getUpperNodeId(const SuperArc *a)
bool isST() const
Definition FTMTree_MT.h:293
void sortLeaves(const bool parallel=false)
idNode getVertInRange(const std::vector< SimplexId > &range, const SimplexId v, const idNode last=0) const
void sortDownArcs(const idNode nid)
Definition FTMTree_MT.h:849
int printTime(Timer &t, const std::string &s, const int debugLevel=2) const
SimplexId getChunkCount(const SimplexId nbVerts=-1, const SimplexId nbTasks=100) const
Definition FTMTree_MT.h:834
SimplexId trunkCTSegmentation(const std::vector< SimplexId > &pendingNodesVerts, const SimplexId begin, const SimplexId stop)
idSuperArc openSuperArc(idNode downNodeId)
std::shared_ptr< FTMTree_MT > clone() const
idNode makeNode(SimplexId vertexId, SimplexId linked=nullVertex)
void buildSegmentation()
use vert2tree to compute the segmentation of the fresh built merge tree.
idSuperArc getCorrespondingSuperArcId(const SimplexId val) const
Definition FTMTree_MT.h:470
idSuperArc makeSuperArc(idNode downNodeId, idNode upNodeId)
void updateCorrespondingNode(const SimplexId vert, const idNode node)
Definition FTMTree_MT.h:498
Node * getUpNode(const SuperArc *a)
bool isJT() const
Definition FTMTree_MT.h:289
Node * getNode(idNode nodeId)
Definition FTMTree_MT.h:393
void closeSuperArc(idSuperArc superArcId, idNode upNodeId)
std::shared_ptr< Scalars > scalars_
Definition FTMTree_MT.h:87
void sortArcs()
Sort tree arcs.
void updateCorrespondingArc(const SimplexId vert, const idSuperArc arc)
Definition FTMTree_MT.h:493
std::vector< idNode > sortedNodes(const bool parallel=false)
idSuperArc getUpSuperArcId(idSuperArc neighborId) const
Definition FTMNode.h:105
idSuperArc clearUpSuperArcs()
Definition FTMNode.h:133
void addDownSuperArcId(idSuperArc downSuperArcId)
Definition FTMNode.h:119
idSuperArc getNumberOfDownSuperArcs() const
Definition FTMNode.h:82
idSuperArc getDownSuperArcId(idSuperArc neighborId) const
Definition FTMNode.h:94
void removeUpSuperArc(idSuperArc idSa)
Definition FTMNode.h:158
idSuperArc clearDownSuperArcs()
Definition FTMNode.h:127
void removeDownSuperArc(idSuperArc idSa)
Definition FTMNode.h:146
SimplexId getVertexId() const
Definition FTMNode.h:54
idSuperArc getNumberOfUpSuperArcs() const
Definition FTMNode.h:86
idSegment size() const
void resize(const std::vector< SimplexId > &sizes)
void setNormalizeIds(const idSuperArc id)
void atomicIncVisited(const SimplexId nb=1)
Definition FTMSuperArc.h:95
idNode getUpNodeId() const
Definition FTMSuperArc.h:68
SimplexId getNbVertSeen() const
std::string printReg() const
size_t regionSize() const
idNode getDownNodeId() const
Definition FTMSuperArc.h:72
long unsigned int idSuperArc
SuperArc index in vect_superArcs_.
std::ostream & operator<<(std::ostream &o, Node const &n)
unsigned int idNode
Node index in vect_nodes_.
The Topology ToolKit.
int SimplexId
Identifier type for simplices of any dimension.
Definition DataTypes.h:22
T begin(std::pair< T, T > &p)
Definition ripserpy.cpp:479
std::shared_ptr< FTMAtomicVector< SuperArc > > superArcs
Definition FTMTree_MT.h:46
std::shared_ptr< FTMAtomicVector< CurrentState > > states
Definition FTMTree_MT.h:60
std::vector< std::list< std::vector< SimplexId > > > trunkSegments
Definition FTMTree_MT.h:54
std::shared_ptr< FTMAtomicVector< Node > > nodes
Definition FTMTree_MT.h:47
std::vector< SimplexId > visitOrder
Definition FTMTree_MT.h:53
std::vector< UF > propagation
Definition FTMTree_MT.h:59
std::vector< idCorresp > vert2tree
Definition FTMTree_MT.h:52
std::vector< UF > ufs
Definition FTMTree_MT.h:58
std::shared_ptr< FTMAtomicVector< idNode > > roots
Definition FTMTree_MT.h:48
std::vector< valence > valences
Definition FTMTree_MT.h:62
std::vector< idNode > leaves
Definition FTMTree_MT.h:49
std::vector< char > openedNodes
Definition FTMTree_MT.h:64
printMsg(debug::output::BOLD+" | | | | | . \\ | | (__| | / __/| |_| / __/|__ _|"+debug::output::ENDCOLOR, debug::Priority::PERFORMANCE, debug::LineMode::NEW, stream)