1 Preliminaries
This section describes our formal setting and presents an overview of our approach. It contains definitions that we adapted from Tierny et al. [84]. An introduction to topology can be found in [22].
1.1 Background
Input data: Without loss of generality, we assume that the input data is a piecewise linear (PL) scalar field defined on a PL manifold with . Values are given at the vertices of and linearly interpolated on higher dimension simplices.
Critical points: Given an isovalue , the sublevel set of , noted , is the preimage of the open interval onto through : . The surlevel set is symmetrically given by . These two objects serve as segmentation tools in many analysis tasks [8].
The points where the topology of differs from that of are the critical points of and their values are called critical values
. Critical points can be classified with their
index , which is 0 for minima, 1 for 1saddles, for saddles and for maxima, with the dimension of .Persistence diagrams: The distribution of critical points of can visually be represented by a topological abstraction called the persistence diagram [24, 15] (Fig. 1). By applying the Elder Rule [22], critical points can be arranged in a set of pairs, such that each critical point appears in only one pair with and . More precisely, the Elder Rule states that as the value increases, if two topological features of , for instance two connected components, meet at a given saddle of , the youngest of the two (the one with the highest minimal value, ) dies at the advantage of the oldest (the one with the lowest minimal value). Critical points and then form a persistence pair.
The persistence diagram embeds each pair in the plane such that its horizontal coordinate equals , and the vertical coordinate of both and is and , corresponding respectively to the birth and death of the pair. The height of the pair is called the persistence and denotes the lifespan of the topological feature created at and destroyed at . In three dimensions, the persistence of the pairs linking critical points of index , and denotes the lifespan of connected components, voids and noncollapsible cycles of .
In practice, persistence diagrams serve as an important visual representation of the distribution of critical points in a scalar dataset. Small oscillations in the data due to noise are typically represented by critical point pairs with low persistence, in the vicinity of the diagonal. In contrast, topological features that are the most prominent in the data are associated with large vertical bars (Fig. 1). In many applications, persistence diagrams help users as a visual guide to interactively tune simplification thresholds in multiscale data segmentation tasks based on the Reeb graph [58, 12, 53, 85, 30, 83] or the MorseSmale complex [32, 33, 65].
Wasserstein distance: Metrics have been defined to evaluate the distance between two scalar fields . The norm , is a classical example. In the context of topological data analysis, multiple metrics [15, 14] have been introduced in order to compare persistence diagrams. In our context, such metrics are key to identifying zones in the data which are similar to one another. Persistence diagrams can be associated with a pointwise distance, noted inspired by the norm. Given persistence pairs and , is given by Eq. 1:
(1) 
The Wasserstein distance [48, 40], sometimes called the Earth Mover’s Distance [46], noted , between the persistence diagrams and is defined in Eq. 2:
(2) 
where is the set of all possible bijections mapping each critical point of to a critical point of the same index in or to the diagonal, noted – which corresponds to the removal of the corresponding feature from the assignment, with a cost .
1.2 Assignment problem
The assignment problem is the problem of choosing an optimal assignment of workers to jobs , assuming numerical ratings are given for each worker’s performance on each job.
Given ratings are summed up in a cost matrix , finding an optimal assignment means choosing a set of independent entries of the matrix so that the sum of these elements is optimal. Independent means than no two such elements should belong to the same row or column (i.e. no two workers should be assigned to the same job and no worker should be given more than one job). In other words, one must find a map of workers and jobs for which the sum is optimal. There are possible assignments, of which several may be optimal, so that an exhaustive search is impractical as gets large.
Similarly, the unbalanced assignment problem is the problem of finding an optimal assignment of workers to jobs, where some jobs or workers might be left unassigned. This is the case of assignments between sets of persistence pairs; where costs are defined for leaving specific pairs unassigned.
The Hungarian algorithm [43, 51] is the first polynomial algorithm proposed by Kuhn to solve the assignment problem. It is an iterative algorithm based on the following two properties:
Theorem 1
If a number is added or subtracted from all the entries of any one row or column of a cost matrix, then an optimal assignment for the resulting cost matrix is also an optimal assignment for the original cost matrix.
This means that the cost matrix can be replaced with where (resp. ) is an arbitrary number which is fixed for the row (resp. the column).
Theorem 2
If is a matrix and is the maximum number of independent zeros of (i.e. number of entries valued at 0), then there are lines (row or columns) which contain all the zeros of .
This allows to determine whether an optimal assignment has been found and thus constitutes the stop criterion.
The algorithm iteratively performs additions and subtractions on lines and columns of the cost matrix, in a way that globally decreases the matrix cost, until the optimal assignment has been found, that is, until the matrix contains a set of independent zeros.
In the remainder we consider the unbalanced KuhnMunkres algorithm [51, 7], an improvement over Kuhn’s original version which follows the same principles, with an enhanced strategy for finding independent elements. The goal is always to reduce the cost matrix and find a maximal set of independent zeros. These independent zeros are marked with a star: they are candidates for the optimal assignment. Zeros which are candidates for being swapped with a starred zero are marked with a prime. Throughout the algorithm, rows and columns of the matrix are marked as covered to restrict the search for candidate zeros.
1.3 Persistence assignment problem
The assignment problem for persistence pairs is very similar to the standard unbalanced assignment problem, except additional costs are defined for not assigning elements (i.e. matching persistence pairs with the diagonal). The assignment between diagrams and then involves the numerical rating associated with assigning with , along with (resp. ) the numerical rating associated with matching (resp. ) with the diagonal.
If and are sets of persistence pairs such that and , then it is possible to solve the persistence assignment problem using the standard KuhnMunkres algorithm with the cost matrix described by Eq. 3, as proposed in [50]:
(3) 
The first line corresponds to matching pairs from to pairs from ; the second one corresponds to the possibility of matching pairs from to the diagonal; the third one is for matching pairs of to the diagonal and the last one completes the cost matrix. The drawback of this approach is that it requires to solve the assignment problem on a cost matrix (that potentially contains two nonsparse blocks where persistence elements are located, see Fig. 4), though the number of distinct elements is at most . As seen in Sec. 2, our algorithm addresses this issue.
1.4 Overview
This section presents a quick overview of our tracking method, which is illustrated in Fig. 9. The input data is a timevarying PL scalar field defined on a PL manifold with .

First, we compute the persistence diagram of the scalar field for every available timestep.

Next, for each pair of two consecutive timesteps and , we consider the two corresponding persistence diagrams and . For each couple of persistence pairs , we define a distance metric corresponding to the similarity of these pairs: (see Sec. 3).

For each pair of consecutive timesteps, we compute a matching function . Every persistence pair of is associated to to , which is either a persistence pair in or so as to minimize the total distance . Finding the optimal involves solving a variant of the classical Assignment Problem, as presented in Sec. 1.3. Only persistence pairs involving critical points of the same index are taken into account.

We compute tracking trajectories starting from the first timestep. If at timestep the matching associates with , then a segment is traced between and . If , the current trajectory ends. Trajectories are grown following this principle throughout all timesteps. Properties are associated to trajectories (time span, critical index), and to trajectory segments (matching cost, scalar value).

Finally, trajectories are postprocessed to detect feature merging or splitting events with a userdefined geometric threshold.
2 Optimized persistence matching
This section presents our novel extension of the KuhnMunkres algorithm, which has been specifically designed to address the computation time bottleneck described in Sec. 1.2.
2.1 Reduced cost matrix
The classical persistence assignment algorithm based on KuhnMunkres considers , a cost matrix. We propose to work instead with , a reduced matrix defined in Eq. 4, where every zero appearing in the last row is considered independent. This amounts to considering that persistence pairs corresponding to rows are not assigned by default. Fig. 4 summarizes the matrices considered by each assignment method.
(4) 
This last row, emulating the diagonal blocks of Fig. 4b requires a specific handling in the optimization procedure. In particular, it requires the first step of the algorithm to subtract minimum elements from columns (and not rows) so as not to have negative elements in the matrix.
As a reminder, the original algorithm proceeds iteratively in two alternating phases: matrix reduction that makes new zeros appear, and augmenting path that finds a maximal set of independent zeros. At the iteration, the current maximal set of independent zeros is made of starred zeros. After a matrix reduction, new zeros appeared that can potentially belong to the new maximal set of independent zeros. Such candidates are primed. A single augmenting path (as in Fig. 3) replaces a set of starred zeros with primed zeros, forming a new set of independent zeros with one more element. Rows and columns of the matrix are marked as covered to restrict the search for candidates in the augmenting path phase. Blue blocks of Algorithm 1 indicate our extension of the KuhnMunkres algorithm.
In this novel extension, an augmenting path constructed in the corresponding phase can start from a starred zero in the last row (and then potentially find a primed zero in its column), but such a path can never access a starred zero in the last row at another step, for the corresponding column would have been covered prior to this (and thus cannot contain a primed zero, see Algorithm 1). A starred zero in the last row can then never be unstarred.
The KuhnMunkres approach has the property to only increase row values (and only decrease column values). When our algorithm working on the reduced matrix ends, it is therefore not possible that the elements on the topright corner of the corresponding full matrix be negative. Furthermore, given Theorem 1, the resulting matrix corresponds to the same assignment problem.
2.2 Optimality
Working with the reduced matrix , however, does not necessarily yield an optimal assignment. When assignments are found in the bottom row, if there has been additions to the matrix rows, then the corresponding matrix would contain a topright block that is not zero, and a topleft block that is not zero either. Thus, the stop criterion stated by Theorem 2 may not be respected when lines are covered (as the real number of independent zeros in is ). Moreover, in our setup, a starred zero in the last column can never be unstarred; this is allowed in the approach on , owing to the bottomright block, initially filled with zeros.
We therefore use the criterion stated in Eq. 5 to ensure that if, at any given iteration of the algorithm, a zero is starred in the last row of column , the cost of assigning the corresponding persistence pair to any other pair is higher than the cost of leaving both unassigned (0 for the column pair and the residual value for row pairs – see Algorithm 1). This specificity is illustrated in Fig. 5.
(5) 
The Eq. 5 criterion is checked whenever a zero appears on the last row after a subtraction is performed on a column by the algorithm. If it is observed, the corresponding column is removed from the problem and the persistence pair is set unassigned.
If the criterion is not respected, we have to report back the reduced problem onto the full matrix (missing banned columns and with reported found residuals ). For this, we need to keep track of residuals, that is, values that have been added or subtracted from each row and column throughout the course of the algorithm. Once these residuals have been reported onto the full matrix, there can be no negative element, and all of the optimization work has been reported (so that we do not start all over again from the beginning, but we start from the optimized output of the first phase).
In practice for persistence diagrams, we always observed that the first phase is sufficient to find an optimal assignment. Using this approach prevents from working with two potentially large blocks of persistence elements, typically occurring with the complete matrix for and . This property is further motivated by the use of geometrical lifting (Sec. 3). The approach is detailed in Algorithm 1.
2.3 Sparse assignment
In practice, it is often observed that some assignments are not possible, and that reordering columns in the associated cost matrix would enable faster lookups and modifications [18], using sparse matrices. With persistence diagrams, the following simple criterion (Eq. 6) can be used to discard lookups for potential matchings.
(6) 
Working with our version of the KuhnMunkres algorithm then becomes interesting for many assignments verify Eq. 6 (Fig. 6), hence greatly reducing the lookup time for zeros, minimal elements, and the access time for operations performed on rows or columns.
On the contrary, the original fullmatrix version of KuhnMunkres deals with nonsparse blocks which have to be accessed and modified constantly throughout the course of the algorithm.
3 Lifted persistence Wasserstein metric
This section highlights the limitations of the natural Wasserstein metric applied to timevarying persistence diagrams and presents an extension that enhances its geometrical stability. Geometrical considerations are motivated, in terms of accuracy and performance.
Persistence diagrams can be embedded into the geometrical domain (Fig. 1). Doing so, one easily sees how different embeddings can correspond to similar persistence diagrams in the birthdeath space. Working in this 2D space does yield irrelevant matchings: as can be seen in Fig. 7, when only the birthdeath coordinates of persistence pairs are considered, a matching can be optimal even if it happens between geometrically distant zones. As a consequence, the distance metric between persistence pairs must be augmented with geometrical considerations.
(7) 
where , and correspond to geometric distances between the extrema involved in the persistence pairs on a given axis. We process diagonal projections as follows (Eq. 8):
(8) 
where the terms , and correspond to the geometric distance between the critical points of a given pair . Intuitively, it accounts for the distance between the critical points to cancel.
A lifted distance is considered by augmenting the geometric distance with coefficients . This aims at giving more or less importance to the birth, death or some of the coordinates during the matching, depending on applicative contexts. For instance, in practice it is desirable to give less importance to the birth coordinate when dealing with  persistence pairs (in other words, for tracking local maxima, see Fig. 8). For the remainder of the paper and the experiments, we used for maxima and for minima, for normalized geometrical extent and scalar values. We observed that using a lifted metric further increases the cost matrix sparsity, resulting in extra speedups.
4 Feature tracking
This section describes the four main stages of our tracking framework, relying on the discussed theoretical setup. Without loss of generality, we assume that the input data is a timevarying 2D or 3D scalar field defined on a PLmanifold. Topological features are extracted for all timesteps (Fig. 9, ab), then matched (Fig. 9, c); trajectories are built from the successive matchings (Fig. 9, d) and postprocessed to detect merging and splitting events.
4.1 Feature detection
First, we compute persistence diagrams for each timestep. We propose using the algorithm by Gueunet et al. [31], in which only  and  persistence pairs are considered.
When the data is noisy, it is possible to discard pairs of low persistence (typically induced by noise) by applying a simple threshold. In practice, this amounts to only considering the most prominent features. Using such a threshold accelerates the matching process, where for approaches based on overlaps, removing topological noise would require a topological simplification of the domain (for example using the approach in [86]), which is computationally expensive.
4.2 Feature matching
4.3 Trajectory extraction
Trajectories are constructed by simply attaching successively matched segments. For all timesteps , if the feature matching associates with , then a segment is traced between and , and is potentially connected back to the previous segment of ’s trajectory. If , the current trajectory ends. Properties are associated to trajectories (time span, critical point index) and to trajectory segments (matching cost, scalar value, persistence value, embedded volume).
4.4 Handling merging and splitting events
Given a userdefined geometrical threshold , we propose to detect events of merging or splitting along trajectories in the following manner. If are two trajectories spanning throughout and respectively, and if for some , , where is a lifted distance, then an event of merging (or splitting) is detected. We consider that a merging event occurs between and at time , when neither nor start at . We then consider that the oldest trajectory takes over the youngest. For example, and meet (according to the criterion) at the last timestep of , and started before , then we disconnect the remainder of from the trajectory before and we connect it so as to continue until ’s original end. Similarly, a splitting event occurs between and at time , when neither nor end at . The process is illustrated in Fig. 10. It is done separately for distinct critical point types: minima, maxima and saddles are not mixed.
5 Results
This section presents experimental results obtained on a desktop computer with two Xeon CPUs (3.0 GHz, 4 cores each), with 64 GB of RAM. We report experiments on 2D and 3D timevarying datasets, that were either simulated with Gerris [55] (von Kármán Vortex street, Boussinesq flow, starting vortex), or acquired (Sea surface height, Isabel hurricane). Persistence diagrams are computed with the implementation of [31] available in the Topology ToolKit [84]; the tracking is restricted to  and  pairs. We implemented our matching (Sec. 2) and tracking approaches (Sec. 4) in C++ as a Topology ToolKit module.
5.1 Application to simulated and acquired datasets
We applied our tracking framework to both simulated and acquired timevarying datasets to outline specific phenomena.
In Fig. 11, we present the results of the tracking framework applied to an oceanographic dataset. The scalar field (sea surface height) is defined on 365 timesteps on a triangular mesh. We can see interesting trajectories corresponding to wellknown oceanic currents. Drifting (a, b) and turbulent behaviors (c) of local extrema are highlighted. In Fig. 12, tracking is performed on the vorticity of highly unstable Boussinesq flow. Thanks to our analysis, trajectories can be filtered according to their temporal lifespan, revealing clearly different trajectory patterns among the turbulent features. This kind of analysis may be easily performed based on other trajectory attributes, depending on applicative contexts. In Fig. 13, we show our approach on a 3D hurricane dataset whose temporal resolution is such that a method based on overlaps of splittree leaves (see Sec. 5.3) could not extract trajectories. In Fig. 14, our tracking framework correctly follows local extrema of the vorticity field in a simulated vortex street.
5.2 Tracking robustness
In the following two sections, we demonstrate the robustness and performance of our tracking framework.
We compare to the greedy approach based on the overlap of volumes [9, 8, 75, 66] of splittree leaves, which amounts to tracking local maxima. In this approach, for every pair of consecutive timesteps , splittree segmentations and are computed (these are a set of connected regions). Overlap scores are then computed for every pair of regions , as the number of common vertices between and . Scores are sorted and is considered matched to the highest positive scoring such that has not been matched before. Trajectories are extracted by repeating the process for all timesteps.
The robustness of our tracking framework is first assessed on a synthetic dataset consisting of whirling gaussians, on which we applied noise (Fig. 15). Identified trajectories are sensibly the same with a perturbation of 10% of the scalar range. The 75% most important features are still correctly tracked after a 25% random perturbation has been applied to the data.
In Fig. 14, our method is compared to the greedy approach, based on overlaps, while decreasing the temporal resolution. The overlap approach yields trajectories corresponding to noise (Fig. 14e), which can be filtered by applying topological simplification beforehand (this would have a significant computational cost as it requires to modify the original function), or by associating the scalar value of the function to every point in the trajectory and then filtering the trajectory in a postprocess step. In our setup, is is much simpler to discard this noise, by using a threshold for discarding small persistence pairs before the matching (implying a faster matching computation). When downsampling the temporal resolution to only 20% of the timesteps, our approach still gives the correct results (Fig. 14, d vs. e). With 15% of the timesteps, our approach (Fig. 14f) still agrees with the tracking on the fullresolution data (Fig. 14b), until preceding features begin to catch up, resulting in a zigzag pattern. By comparison, the overlap method fails to correctly track meaningful regions from the beginning of the simulation to its end; it is indeed dependent on the geometry of overlaps, which is unstable. It can be argued that the locality captured by overlaps is emulated in our framework by embedding and lifting the Wasserstein metric, when the overlap method does not take persistence into account when matching regions. Also note that if the saddle component of persistence pairs associated to maxima is ignored (i.e. if in Eq. 7 and Eq. 8) during the matching, then the geometrical distance can be insufficient for correctly tracking these persistence pairs (c). Therefore, the problem of matching persistence pairs for tracking topological features cannot be reduced to the (unbalanced) problem of assigning critical points in 4 dimensions (3 for the geometrical extent, one for the scalar value).
Fig. 13 further illustrates the robustness of our approach when downsampling the data temporal resolution. In hurricane datasets, local maxima can be displaced to geometrically distant zones between timesteps if those are taken at multipleday intervals. This unstable behavior and the absence of obvious overlaps makes it particularly difficult to track extrema; nevertheless, our framework managed to track them at a very low timeresolution.
5.3 Tracking performance
We then compare our framework with our implementation of the approach based on overlaps [8] on the ground of performance. Figures are given in Table 1. Note that our approach has the advantage of taking persistence diagrams as inputs, so it can be applied to unstructured or timevarying meshes, for which overlap computations are not trivial. Our approach is also relatively dimensionindependent: though in 3D, computing overlaps is very timeconsuming (Fig. 1Isabel), the complexity of the Wasserstein matcher, which only takes embedded persistence diagrams as inputs, for a given number of persistence pairs is sensibly equivalent. For both Isabel and Sea surface height datasets, we applied a 4% persistence filtering on input persistence diagrams. As the experiments show, our approach is faster in practice than the overlap method with bestmatch search.
Dataset  Preproc (s)  Matching (s)  Postproc  

FTM  [8]  ours  (s)  
Boussinesq  116  75  18  4.7 
Vortex street  45  23  18  2.8 
Isabel (3D)  863  3k  17  162 
Sea height  568  N.A.  277  113 
5.4 Matching performance
Next, we compare the performance of the matching method we introduced in Sec. 2 to two other stateofthe art algorithms.
Dataset  Sizes of diagrams  Time (s)  

[90]  ours  
Starting vortex  68.6  1.26  
Isabel  72.2  3.58  
Boussinesq  11.1k  102  
Sea height  26.5k  155 
We compare it to the reference approach for the exact assignment problem [90] based on the KuhnMunkres algorithm, and to our implementation of the approximate approach based on the auction algorithm [4, 41], on the ground of performance.
Table 2 shows that our new assignment algorithm is up to two orders of magnitude faster than the classical exact approach [90]. In particular, the best speedups occur for the larger datasets which indicates that our approach also benefits from an improved scaling.
It is often useful in practice to discard lowpersistence pairs prior to any topological data analysis as these correspond to noise. Fig. 16 compares the running times of our approach, [90] and [4] as more and more lowpersistence pairs are taken into account. When removing pairs whose persistence is below 5% of the scalar range, which is commonly accepted as a conservative threshold, our approach is faster than all competing alternatives. When considering more lowpersistence features, below 4%, our approach is competitive with the approximated auction approach with 1% error. Below 2%, only noise is typically added in the process. The performance of our algorithm becomes comparable to that of the highprecision auction approximation although our approach guarantees exact results.
5.5 Limitations
As we described, our framework enables the tracking of 01 and  persistence pairs. It would be interesting to extend it to support the tracking of saddlesaddle pairs (in 3D) and see its application to meaningful use cases.
Besides, the lifting coefficients proposed in our metric (Eq. 7) might be seen as supplementary parameters that have to be tuned according to the dataset and applicative domain. Nonetheless, we observed in our experiments that these parameters do not require finetuning to produce meaningful tracking trajectories. The extent to which these can be enhanced by finetuning is left to future work.
The lifted distance can be generalized to take other parameters, such as the geometrical volume, mass, feature speed, into account, and be finetuned to answer the specificity of various scientific domains. Merging and splitting might also be enhanced, or given more flexibility, for instance with additional criteria. We also believe that the performance of the postprocessing phase can be improved.
Additionally, we believe that the approximate auction algorithm can also take the lifted persistence metric into account by performing Wasserstein matchings between persistence pairs in 5 dimensions, and possibly benefit from geometrybased lookup accelerations, as suggested in lower dimension in [41]. It remains to be clarified how the quality of the matchings is affected in practice by using an approximate matching method, and how sparsity can enhance the research phase for the auction algorithm.
We note that the theoretical complexity of our matching method is, as the Munkres method, cubic; however, the two orders of magnitude speedups demonstrated in our experiments allow to study more challenging datasets. For very large case studies, the use of persistence thresholds could prove quite helpful for controlling the computing time of matchings. Among other nontrivial tracking methods, some graph matching methods are based on graphedit distances [26, 1]. Their adaptation to the case of persistence diagrams or other topological structures (such as contour trees and Reeb graphs) may enable an additional structural regularization, this ought to be investigated in future work.
6 Conclusion
In this paper, we presented an original framework for tracking topological features in a robust and efficient way. It is the first approach combining topological data analysis and transport for feature tracking. As the kernel of our approach, we proposed a sparsecompliant extension of the seminal assignment algorithm for the exact matching of persistence diagrams, leveraging in practice important speedups. We introduced a new metric for persistence diagrams that enhances geometrical stability and further improves computation time. Overall, in comparison with overlapbased techniques, our approach displays improved performance and robustness to temporal downsampling, as experiments have shown.
We plan to release the implementation of our tracking framework opensource as a part of TTK [84] in the near future; we hope that it will be useful to the community with an interest for efficient tracking methods. We look forward to adapting it to tracking phenomena in insitu contexts, where the largescale timevarying data is accessed in a streaming fashion. As we are also interested in larger datasets, we are currently carrying out scaling tests on complex physical case studies available at Total S.A., for which one needs specifically adapted rendering techniques [47] to apprehend the resulting graphical complexity of the topology evolution.
We also believe that the application potential of our matching framework can be studied for tasks other than timetracking, for instance, selfpattern matching and symmetry detection
[82], or feature comparisons in ensemble data.Acknowledgements.
This work is partially supported by the Bpifrance grant “AVIDO” (Programme d’Investissements d’Avenir, reference P1120172661376/DOS0021427) and by the French National Association for Research and Technology (ANRT), in the framework of the LIP6  Total SA CIFRE partnership reference 2016/0010. The authors would like to thank the anonymous reviewers for their thoughtful remarks and suggestions.References
 [1] K. Beketayev, D. Yeliussizov, D. Morozov, G. H. Weber, and B. Hamann. Measuring the distance between merge trees. In Topological Methods in Data Analysis and Visualization III, pp. 151–165. Springer, 2014.
 [2] J. C. Bennett, H. Abbasi, P. T. Bremer, R. Grout, A. Gyulassy, T. Jin, S. Klasky, H. Kolla, M. Parashar, V. Pascucci, P. Pebay, D. Thompson, H. Yu, F. Zhang, and J. Chen. Combining insitu and intransit processing to enable extremescale scientific analysis. In High Performance Computing, Networking, Storage and Analysis (SC), 2012 International Conference for, pp. 1–9, Nov 2012.
 [3] D. P. Bertsekas. Network optimization: Continuous and discrete models, 1998.
 [4] D. P. Bertsekas and D. A. Castanon. The auction algorithm for the transportation problem. Ann. Oper. Res., 20(14):67–96, 1989.
 [5] S. Biasotti, D. Giorgio, M. Spagnuolo, and B. Falcidieno. Reeb graphs for shape analysis and applications. TCS, 2008.
 [6] A. Bock, H. Doraiswamy, A. Summers, and C. Silva. Topoangler: Interactive topologybased extraction of fishes. IEEE Transactions on Visualization and Computer Graphics, 24(1):812–821, Jan 2018.
 [7] F. Bourgeois and J.C. Lassalle. An extension of the munkres algorithm for the assignment problem to rectangular matrices. Commun. ACM, 14(12):802–804, Dec. 1971.
 [8] P. Bremer, G. Weber, J. Tierny, V. Pascucci, M. Day, and J. Bell. Interactive exploration and analysis of large scale simulations using topologybased data segmentation. IEEE TVCG, 17(9):1307–1324, 2011.
 [9] P. T. Bremer, G. Weber, V. Pascucci, M. Day, and J. Bell. Analyzing and tracking burning structures in lean premixed hydrogen flames. IEEE Transactions on Visualization and Computer Graphics, 16(2):248–260, March 2010.
 [10] R. Burkard, M. Dell’Amico, and S. Martello. Assignment Problems. Society for Industrial and Applied Mathematics, Philadelphia, PA, USA, 2009.
 [11] H. Carr, J. Snoeyink, and U. Axen. Computing contour trees in all dimensions. Computational Geometry, 24(2):75–94, 2003.
 [12] H. Carr, J. Snoeyink, and M. van de Panne. Simplifying flexible isosurfaces using local geometric measures. In IEEE VIS, pp. 497–504, 2004.
 [13] M. Carrière, M. Cuturi, and S. Oudot. Sliced wasserstein kernel for persistence diagrams. In ICML, 2017.
 [14] F. Chazal, D. CohenSteiner, M. Glisse, L. J. Guibas, and S. Oudot. Proximity of Persistence Modules and their Diagrams. In SoCG, 2009.
 [15] D. CohenSteiner, H. Edelsbrunner, and J. Harer. Stability of persistence diagrams. In Symp. on Comp. Geom., pp. 263–271, 2005.
 [16] D. CohenSteiner, H. Edelsbrunner, and D. Morozov. Vines and vineyards by updating persistence in linear time. In Proceedings of the Twentysecond Annual Symposium on Computational Geometry, SCG ’06, pp. 119–126. ACM, New York, NY, USA, 2006.
 [17] D. A. S. C. A. Committee. Synergistic challenges in dataintensive science and exascale computing. Technical report, DoE Advanced Scientific Computing Advisory Committee, Data Subcommittee, 2013.
 [18] H. Cui, J. Zhang, C. Cui, and Q. Chen. Solving largescale assignment problems by kuhnmunkres algorithm. In International Conference on Advances in Mechanical Engineering and Industrial Informatics (AMEII), 01 2016.
 [19] M. Cuturi. Sinkhorn distances: Lightspeed computation of optimal transport. In C. J. C. Burges, L. Bottou, M. Welling, Z. Ghahramani, and K. Q. Weinberger, eds., Advances in Neural Information Processing Systems 26, pp. 2292–2300. Curran Associates, Inc., 2013.
 [20] L. De Floriani, U. Fugacci, F. Iuricich, and P. Magillo. Morse complexes for shape segmentation and homological analysis: discrete models and algorithms. Comp. Grap. For., 2015.
 [21] H. Edelsbrunner and J. Harer. Jacobi sets of multiple morse functions. In Foundations of Computational Mathematics, 2004.
 [22] H. Edelsbrunner and J. Harer. Computational Topology: An Introduction. American Mathematical Society, 2009.
 [23] H. Edelsbrunner, J. Harer, A. Mascarenhas, and V. Pascucci. Timevarying reeb graphs for continuous spacetime data. In Proceedings of the Twentieth Annual Symposium on Computational Geometry, SCG ’04, pp. 366–372. ACM, 2004.
 [24] H. Edelsbrunner, D. Letscher, and A. Zomorodian. Topological persistence and simplification. Disc. Compu. Geom., 28(4):511–533, 2002.
 [25] G. Favelier, C. Gueunet, and J. Tierny. Visualizing ensembles of viscous fingers. In IEEE SciVis Contest, 2016.
 [26] X. Gao, B. Xiao, D. Tao, and X. Li. A survey of graph edit distance. Pattern Analysis and applications, 13(1):113–129, 2010.
 [27] S. Grottel, G. Reina, J. Vrabec, and T. Ertl. Visual verification and analysis of cluster detection for molecular dynamics. IEEE Transactions on Visualization and Computer Graphics, 13(6):1624–1631, Nov 2007.
 [28] Y. Gu and C. Wang. Transgraph: Hierarchical exploration of transition relationships in timevarying volumetric data. IEEE Transactions on Visualization and Computer Graphics, 17(12):2015–2024, Dec 2011.
 [29] D. Guenther, R. AlvarezBoto, J. ContrerasGarcia, J.P. Piquemal, and J. Tierny. Characterizing molecular interactions in chemical systems. IEEE TVCG, 20(12):2476–2485, 2014.
 [30] C. Gueunet, P. Fortin, J. Jomier, and J. Tierny. Contour forests: Fast multithreaded augmented contour trees. In IEEE LDAV, pp. 85–92, 2016.
 [31] C. Gueunet, P. Fortin, J. Jomier, and J. Tierny. Taskbased augmented merge trees with Fibonacci heaps. In IEEE LDAV, 2017.
 [32] A. Gyulassy, P. T. Bremer, B. Hamann, and V. Pascucci. A practical approach to morsesmale complex computation: Scalability and generality. IEEE TVCG, 14(6):1619–1626, 2008.
 [33] A. Gyulassy, D. Guenther, J. A. Levine, J. Tierny, and V. Pascucci. Conforming morsesmale complexes. IEEE TVCG, 20(12):2595–2603, 2014.
 [34] A. Gyulassy, A. Knoll, K. Lau, B. Wang, P. Bremer, M. Papka, L. A. Curtiss, and V. Pascucci. Interstitial and interlayer ion diffusion geometry extraction in graphitic nanosphere battery materials. IEEE TVCG, 2015.
 [35] S. M. Hannon and J. A. Thomson. Aircraft wake vortex detection and measurement with pulsed solidstate coherent laser radar. Journal of Modern Optics, 41(11):2175–2196, 1994.
 [36] C. Heine, H. Leitte, M. Hlawitschka, F. Iuricich, L. De Floriani, G. Scheuermann, H. Hagen, and C. Garth. A survey of topologybased methods in visualization. Comp. Grap. For., 35(3):643–667, 2016.
 [37] G. Ji and H.W. Shen. Efficient Isosurface Tracking Using Precomputed Correspondence Table. In O. Deussen, C. Hansen, D. Keim, and D. Saupe, eds., Eurographics / IEEE VGTC Symposium on Visualization. The Eurographics Association, 2004.
 [38] G. Ji and H.W. Shen. Feature tracking using earth mover’s distance and global optimization. 2006.
 [39] G. Ji, H.W. Shen, and R. Wenger. Volume tracking using higher dimensional isosurfacing. In Proceedings of the 14th IEEE Visualization 2003 (VIS’03), VIS ’03, pp. 28–, 2003.
 [40] L. Kantorovich. On the translocation of masses. AS USSR, 1942.
 [41] M. Kerber, D. Morozov, and A. Nigmetov. Geometry helps to compare persistence diagrams. J. Exp. Algorithmics, 22:1.4:1–1.4:20, Sept. 2017.
 [42] T. Klein and T. Ertl. Scalespace tracking of critical points in 3d vector fields. In H. Hauser, H. Hagen, and H. Theisel, eds., Topologybased Methods in Visualization, pp. 35–49. Springer Berlin Heidelberg, Berlin, Heidelberg, 2007.
 [43] H. W. Kuhn and B. Yaw. The hungarian method for the assignment problem. Naval Res. Logist. Quart, pp. 83–97, 1955.
 [44] A. Landge, V. Pascucci, A. Gyulassy, J. Bennett, H. Kolla, J. Chen, and T. Bremer. Insitu feature extraction of large scale combustion simulations using segmented merge trees. In SuperComputing, 2014.
 [45] D. E. Laney, P. Bremer, A. Mascarenhas, P. Miller, and V. Pascucci. Understanding the structure of the turbulent mixing layer in hydrodynamic instabilities. IEEE TVCG, 2006.
 [46] E. Levina and P. Bickel. The earthmover’s distance is the mallows distance: some insights from statistics. In IEEE ICCV, vol. 2, pp. 251–256, 2001.
 [47] J. Lukasczyk, G. Weber, R. Maciejewski, C. Garth, and H. Leitte. Nested tracking graphs. Comp. Graph. For., 36(3):12–22, 2017.
 [48] G. Monge. Mémoire sur la théorie des déblais et des remblais. Académie Royale des Sciences de Paris, 1781.
 [49] K. Moreland, R. Oldfield, P. Marion, S. Jourdain, N. Podhorszki, V. Vishwanath, N. Fabian, C. Docan, M. Parashar, M. Hereld, M. E. Papka, and S. Klasky. Examples of in transit visualization. In Proceedings of the 2Nd International Workshop on Petascal Data Analytics: Challenges and Opportunities, PDAC ’11, pp. 1–6. ACM, New York, NY, USA, 2011.
 [50] D. Morozov. Dionysus. http://www.mrzv.org/software/dionysus/, 2010.
 [51] J. Munkres. Algorithms for the assignment and transportation problems, 1957.
 [52] S. Ozer, D. Silver, K. Bemis, and P. Martin. Activity detection in scientific visualization. IEEE Transactions on Visualization and Computer Graphics, 20(3):377–390, March 2014.
 [53] V. Pascucci, G. Scorzelli, P. T. Bremer, and A. Mascarenhas. Robust online computation of Reeb graphs: simplicity and speed. ToG, 26(3):58, 2007.
 [54] V. Pascucci, X. Tricoche, H. Hagen, and J. Tierny. Topological Data Analysis and Visualization: Theory, Algorithms and Applications. Springer, 2010.
 [55] S. Popinet. Gerris: A treebased adaptive solver for the incompressible euler equations in complex geometries. J. Comput. Phys., 190(2):572–600, 2003.
 [56] F. H. Post, B. Vrolijk, H. Hauser, R. S. Laramee, and H. Doleisch. The state of the art in flow visualisation: Feature extraction and tracking. Computer Graphics Forum, 22(4):775–792, 2003.
 [57] M. Rasquin, P. Marion, V. Vishwanath, B. Matthews, M. Hereld, K. Jansen, R. Loy, A. Bauer, M. Zhou, O. Sahni, J. Fu, N. Liu, C. Carothers, M. Shephard, M. Papka, K. Kumaran, and B. Geveci. Electronic poster: Covisualization of full data and in situ data extracts from unstructured grid cfd at 160k cores. In Proceedings of the 2011 Companion on High Performance Computing Networking, Storage and Analysis Companion, SC ’11 Companion, pp. 103–104. ACM, New York, NY, USA, 2011.
 [58] G. Reeb. Sur les points singuliers d’une forme de Pfaff complètement intégrable ou d’une fonction numérique. Acad. des Sci., 1946.
 [59] F. Reinders, F. H. Post, and H. J. Spoelder. Visualization of timedependent data with feature tracking and event detection. The Visual Computer, 17(1):55–71, Feb 2001.
 [60] J. Reininghaus, S. Huber, U. Bauer, and R. Kwitt. A stable multiscale kernel for topological machine learning. In IEEE CVPR, 2015.
 [61] J. Reininghaus, J. Kasten, T. Weinkauf, and I. Hotz. Efficient computation of combinatorial feature flow fields. IEEE Transactions on Visualization and Computer Graphics, 18(9):1563–1573, 2012.
 [62] T. Ringler, M. Petersen, R. L. Higdon, D. Jacobsen, P. W. Jones, and M. Maltrud. A multiresolution approach to global ocean modeling. Ocean Modelling, 69:211 – 232, 2013.
 [63] M. Rivi, L. Calori, G. Muscianisi, and V. Slavnic. Insitu visualization: Stateoftheart and some use cases. 2011.
 [64] K. Robbins, C. Jeffrey, and S. Robbins. Visualization of splitting and merging processes. Journal of Visual Languages and Computing, 11(6):593 – 614, 2000.
 [65] V. Robins, P. Wood, and A. Sheppard. Theory and algorithms for constructing discrete morse complexes from grayscale digital images. IEEE Trans. on Pat. Ana. and Mach. Int., 2011.

[66]
H. Saikia and T. Weinkauf.
Global feature tracking and similarity estimation in timedependent scalar fields.
Comput. Graph. Forum, 36(3):1–11, June 2017.  [67] R. Samtaney, D. Silver, N. Zabusky, and J. Cao. Visualizing features and tracking their evolution. Computer, 27(7):20–27, 1994.
 [68] K. Shi, H. Theisel, H. Hauser, T. Weinkauf, K. Matkovic, H.C. Hege, and H.P. Seidel. Path Line Attributes  an Information Visualization Approach to Analyzing the Dynamic Behavior of 3D TimeDependent Flow Fields, pp. 75–88. Springer, 2009.
 [69] K. Shi, H. Theisel, T. Weinkauf, H. Hauser, H.C. Hege, and H.P. Seidel. Path line oriented topology for periodic 2D timedependent vector fields. In Proc. Eurographics / IEEE VGTC Symposium on Visualization (EuroVis ’06), pp. 139–146. Lisbon, Portugal, May 2006.
 [70] D. Silver. Objectoriented visualization. IEEE Comput. Graph. Appl., 15(3):54–62, May 1995.
 [71] D. Silver and X. Wang. Volume tracking. In Visualization ’96. Proceedings., pp. 157–164, Oct 1996.
 [72] D. Silver and X. Wang. Tracking and visualizing turbulent 3d features. IEEE Transactions on Visualization and Computer Graphics, 3(2):129–141, Apr 1997.
 [73] D. Silver and X. Wang. Tracking scalar features in unstructured data sets. In Visualization ’98. Proceedings, pp. 79–86, Oct 1998.
 [74] D. Silver and X. Wang. Visualizing evolving scalar phenomena. Future Generation Computer Systems, 15(1):99 – 108, 1999.
 [75] B. S. Sohn and C. Bajaj. Timevarying contour topology. IEEE TVCG, 12(1):14–25, 2006.
 [76] J. Solomon, F. de Goes, G. Peyré, M. Cuturi, A. Butscher, A. Nguyen, T. Du, and L. J. Guibas. Convolutional wasserstein distances: efficient optimal transportation on geometric domains. ACM Trans. Graph., 34(4):66:1–66:11, 2015.
 [77] J. Solomon, G. Peyré, V. G. Kim, and S. Sra. Entropic metric alignment for correspondence problems. ACM Trans. Graph., 35(4):72:1–72:13, July 2016.
 [78] T. Sousbie. The persistent cosmic web and its filamentary structure: Theory and implementations. Royal Astronomical Society, 414(1):350–383, 2011.
 [79] H. Theisel and H.P. Seidel. Feature flow fields. In Procs. of the Symp. on Data Visualisation 2003, VISSYM ’03, pp. 141–148. Eurographics Association, AirelaVille, Switzerland, Switzerland, 2003.
 [80] H. Theisel, T. Weinkauf, H. C. Hege, and H. P. Seidel. Stream line and path line oriented topology for 2d timedependent vector fields. In IEEE Visualization 2004, pp. 321–328, Oct 2004.
 [81] H. Theisel, T. Weinkauf, H. C. Hege, and H. P. Seidel. Topological methods for 2d timedependent vector fields based on stream lines and path lines. IEEE Transactions on Visualization and Computer Graphics, 11(4):383–394, July 2005.
 [82] D. M. Thomas and V. Natarajan. Multiscale symmetry detection in scalar fields by clustering contours. IEEE TVCG, 2014.
 [83] J. Tierny and H. Carr. Jacobi fiber surfaces for bivariate Reeb space computation. IEEE TVCG, 23(1):960–969, 2016.
 [84] J. Tierny, G. Favelier, J. A. Levine, C. Gueunet, and M. Michaux. The Topology ToolKit. IEEE TVCG, 24(1):832–842, 2017.
 [85] J. Tierny, A. Gyulassy, E. Simon, and V. Pascucci. Loop surgery for volumetric meshes: Reeb graphs reduced to contour trees. IEEE TVCG, 15(6):1177–1184, 2009.
 [86] J. Tierny and V. Pascucci. Generalized topological simplification of scalar fields on surfaces. IEEE TVCG, 18(12):2005–2013, 2012.
 [87] X. Tricoche, T. Wischgoll, R. G. Scheuermann, and H. Hagen. Topology tracking for the visualization of timedependent twodimensional flows. Computers and Graphics, 26(2):249 – 257, 2002.
 [88] C. Villani. Optimal Transport: Old and New. Springer, 2009 ed., Sept. 2008.
 [89] C. Wang and J. Tao. Graphs in scientific visualization: A survey. Computer Graphics Forum, 36(1):263–287, 2017.
 [90] J. Weaver. An implementation of the kuhn–munkres algorithm. https://github.com/saebyn/munkrescpp, 2013.
 [91] G. Weber, P.T. Bremer, M. Day, J. Bell, and V. Pascucci. Feature Tracking Using Reeb Graphs, pp. 241–253. Springer Berlin Heidelberg, Berlin, Heidelberg, 2011.
 [92] W. Widanagamaachchi, J. Chen, P. Klacansky, V. Pascucci, H. Kolla, A. Bhagatwala, and P. T. Bremer. Tracking features in embedded surfaces: Understanding extinction in turbulent combustion. In 2015 IEEE 5th Symposium on Large Data Analysis and Visualization (LDAV), pp. 9–16, Oct 2015.
 [93] H. Yu, C. Wang, R. W. Grout, J. H. Chen, and K. L. Ma. In situ visualization for largescale combustion simulations. IEEE Computer Graphics and Applications, 30(3):45–57, May 2010.
 [94] F. Zhang, E. Fiorelli, and N. E. Leonard. Exploring scalar fields using multiple sensor platforms: Tracking level curves. In 2007 46th IEEE Conference on Decision and Control, pp. 3579–3584, Dec 2007.
 [95] F. Zhang, S. Lasluisa, T. Jin, I. Rodero, H. Bui, and M. Parashar. Insitu featurebased objects tracking for largescale scientific simulations. In 2012 SC Companion: High Performance Computing, Networking Storage and Analysis, pp. 736–740, Nov 2012.
Comments
There are no comments yet.