From 490fed367bb97a96b90caa6ef04265c063d91df1 Mon Sep 17 00:00:00 2001 From: Arnur Nigmetov Date: Sat, 7 Mar 2020 08:44:20 +0100 Subject: Fix multiple bugs in matching distance for modules, add example. --- README.txt | 15 +++++-- bottleneck/README | 2 +- matching/CMakeLists.txt | 4 ++ matching/README.md | 72 +++++++++++++++++++++++++++++++-- matching/example/matching_dist.cpp | 3 -- matching/include/matching_distance.h | 2 +- matching/include/matching_distance.hpp | 1 + matching/include/persistence_module.h | 2 +- matching/include/persistence_module.hpp | 19 ++++++--- 9 files changed, 101 insertions(+), 19 deletions(-) diff --git a/README.txt b/README.txt index ae3f6a0..48a77c0 100644 --- a/README.txt +++ b/README.txt @@ -1,5 +1,6 @@ This repository contains software to compute bottleneck and Wasserstein -distances between persistence diagrams. +distances between persistence diagrams, and matching distance +between 2-parameter persistence modules and (1-critical) bi-filtrations. The software is licensed under BSD license, see license.txt file. @@ -9,13 +10,14 @@ you probably do not need to worry about that. See README files in subdirectories for usage and building. If you use Hera in your project, we would appreciate if you -cite the corresponding paper: +cite the corresponding paper. + +Bottleneck or Wasserstein distance: + Michael Kerber, Dmitriy Morozov, and Arnur Nigmetov, "Geometry Helps to Compare Persistence Diagrams.", Journal of Experimental Algorithmics, vol. 22, 2017, pp. 1--20. (conference version: ALENEX 2016). -The BibTeX is below: - @article{jea_hera, title={Geometry helps to compare persistence diagrams}, @@ -26,3 +28,8 @@ The BibTeX is below: year={2017}, publisher={ACM New York, NY, USA} } + +Matching distance: + +Michael Kerber, Arnur Nigmetov, "Efficient Approximation of the Matching +Distance for 2-parameter persistence.", SoCG 2020 diff --git a/bottleneck/README b/bottleneck/README index 8b368af..c04390b 100644 --- a/bottleneck/README +++ b/bottleneck/README @@ -1,6 +1,6 @@ Accompanying paper: M. Kerber, D. Morozov, A. Nigmetov. Geometry Helps To Compare Persistence Diagrams (ALENEX 2016, http://www.geometrie.tugraz.at/nigmetov/geom_dist.pdf) -Bug reports can be sent to "nigmetov EMAIL SIGN tugraz DOT at". +Bug reports can be sent to "anigmetov EMAIL SIGN lbl DOT gov". # Dependencies diff --git a/matching/CMakeLists.txt b/matching/CMakeLists.txt index 3ee0f6b..a391d84 100644 --- a/matching/CMakeLists.txt +++ b/matching/CMakeLists.txt @@ -48,5 +48,9 @@ endif() add_executable(matching_dist "example/matching_dist.cpp" ${MD_HEADERS} ${BT_HEADERS} ) target_link_libraries(matching_dist PUBLIC ${libraries}) +add_executable(module_example "example/module_example.cpp" ${MD_HEADERS} ${BT_HEADERS} ) +target_link_libraries(module_example PUBLIC ${libraries}) + + add_executable(matching_distance_test ${SRC_TEST_FILES} ${BT_HEADERS} ${MD_HEADERS}) target_link_libraries(matching_distance_test PUBLIC ${libraries}) diff --git a/matching/README.md b/matching/README.md index 7dc1874..fc97441 100644 --- a/matching/README.md +++ b/matching/README.md @@ -1,5 +1,69 @@ -Matching distance between bifiltrations. +# Matching distance between bifiltrations and 2-persistence modules. -Currently supports only 1-critical bi-filtrations -in PHAT-like format: boundary matrix + critical values -of a simplex in each row +## Accompanying paper +M. Kerber, A. Nigmetov, +Efficient Approximation of the Matching Distance for 2-parameter persistence. +SoCG 2020. + +Bug reports can be sent to "anigmetov EMAIL SIGN lbl DOT gov". + +## Dependencies + +* Your compiler must support C++11. +* Boost. + +## Usage: + +1. To use a standalone command-line utility matching_dist: + +`./matching_dist -d dimension -e relative_error file1 file2` + +If `relative_error` is not specified, the default 0.1 is used. +If `dimension` is not specified, the default value is 0. +Run `./matching_dist` without parameters to see other options. + +The output is an approximation of the exact distance (which is assumed to be non-zero). +Precisely: if *d_exact* is the true distance and *d_approx* is the output, then + +> | d_exact - d_approx | / d_exact < relative_error. + +Files file1 and file2 must contain 1-critical bi-filtrations in a plain text format which is similar to PHAT. The first line of a file must say *bifiltration_phat_like*. The second line contains the total number of simplices *N*. The next *N* lines contain simplices in the format *dim x y boundary*. +* *dim*: the dimension of the simplex +* *x, y*: coordinates of the critical value +* *boundary*: indices of simplices forming the boundary of the current simplex. Indices are separated by space. +* Simplices are indexed starting from 0. + +For example, the bi-filtration of a segment with vertices appearing at (0,0) and the 1-segment appearing at (3,4) shall be written as: + +> bifiltration_phat_like +> 3 +> \# lines starting with \# are ignored +> \# vertex A has dimension 0, hence no boundary, its index is 0 +> 0 0 0 +> \# vertex B has index 1 +> 0 0 0 +> \# 1-dimensional simplex {A, B} +> 1 3 4 0 1 + +2. To use from your code. + +Here you can compute the matching distance either between bi-filtrations or between persistence modules. +First, you need to include `#include "matching_distance.h"` Practically every class you need is parameterized by Real type, which should be either float or double. The header provides two functions called `matching_distance.` +See `example/module_example.cpp` for additional details. + +## License + +See `licence.txt` in the repository root folder. + +## Building + +CMakeLists.txt can be used to build the command-line utility in the standard +way. On Linux/Mac/Windows with Cygwin: +> `mkdir build` +> `cd build` +> `cmake ..` +> `make` + +On Windows with Visual Studio: use `cmake-gui` to create the solution in build directory and build it with VS. + +The library itself is header-only and does not require separate compilation. diff --git a/matching/example/matching_dist.cpp b/matching/example/matching_dist.cpp index 13e5c6d..734f6cc 100644 --- a/matching/example/matching_dist.cpp +++ b/matching/example/matching_dist.cpp @@ -13,9 +13,6 @@ #include "spdlog/spdlog.h" #include "spdlog/fmt/ostr.h" -//#include "persistence_module.h" -#include "bifiltration.h" -#include "box.h" #include "matching_distance.h" using Real = double; diff --git a/matching/include/matching_distance.h b/matching/include/matching_distance.h index 276cc3a..e82a97c 100644 --- a/matching/include/matching_distance.h +++ b/matching/include/matching_distance.h @@ -152,7 +152,7 @@ namespace md { int max_depth {6}; // maximal number of refinenemnts int initialization_depth {3}; int dim {0}; // in which dim to calculate the distance; use ALL_DIMENSIONS to get max over all dims - BoundStrategy bound_strategy {BoundStrategy::bruteforce}; + BoundStrategy bound_strategy {BoundStrategy::local_combined}; TraverseStrategy traverse_strategy {TraverseStrategy::breadth_first}; bool tolerate_max_iter_exceeded {true}; Real actual_error {std::numeric_limits::max()}; diff --git a/matching/include/matching_distance.hpp b/matching/include/matching_distance.hpp index 48c8464..7cff073 100644 --- a/matching/include/matching_distance.hpp +++ b/matching/include/matching_distance.hpp @@ -193,6 +193,7 @@ namespace md { // make all coordinates non-negative auto min_coord = std::min(module_a_.minimal_coordinate(), module_b_.minimal_coordinate()); + spd::debug("in DistanceCalculator ctor, min_coord = {}", min_coord); if (min_coord < 0) { module_a_.translate(-min_coord); module_b_.translate(-min_coord); diff --git a/matching/include/persistence_module.h b/matching/include/persistence_module.h index e99771f..4a261bb 100644 --- a/matching/include/persistence_module.h +++ b/matching/include/persistence_module.h @@ -47,7 +47,7 @@ namespace md { IndexVec components_; Relation() {} - Relation(const Point& _pos, const IndexVec& _components); + Relation(const Point& _pos, const IndexVec& _components) : position_(_pos), components_(_components) {} Real get_x() const { return position_.x; } Real get_y() const { return position_.y; } diff --git a/matching/include/persistence_module.hpp b/matching/include/persistence_module.hpp index 6e49b2e..128fed9 100644 --- a/matching/include/persistence_module.hpp +++ b/matching/include/persistence_module.hpp @@ -34,10 +34,10 @@ namespace md { template void ModulePresentation::init_boundaries() { - max_x_ = std::numeric_limits::max(); - max_y_ = std::numeric_limits::max(); - min_x_ = -std::numeric_limits::max(); - min_y_ = -std::numeric_limits::max(); + max_x_ = -std::numeric_limits::max(); + max_y_ = -std::numeric_limits::max(); + min_x_ = std::numeric_limits::max(); + min_y_ = std::numeric_limits::max(); for(const auto& gen : positions_) { min_x_ = std::min(gen.x, min_x_); @@ -55,6 +55,7 @@ namespace md { generators_(_generators), relations_(_relations) { + positions_ = concat_gen_and_rel_positions(generators_, relations_); init_boundaries(); } @@ -85,6 +86,7 @@ namespace md { void ModulePresentation::project_generators(const DualPoint& slice, IndexVec& sorted_indices, RealVec& projections) const { + spd::debug("Enter project_generators, slice = {}", slice); size_t num_gens = generators_.size(); RealVec gen_values; @@ -97,6 +99,7 @@ namespace md { projections.reserve(num_gens); for(auto i : sorted_indices) { projections.push_back(gen_values[i]); + spd::debug("added push = {}", gen_values[i]); } } @@ -104,6 +107,8 @@ namespace md { void ModulePresentation::project_relations(const DualPoint& slice, IndexVec& sorted_rel_indices, RealVec& projections) const { + + spd::debug("Enter project_relations, slice = {}", slice); size_t num_rels = relations_.size(); RealVec rel_values; @@ -116,12 +121,14 @@ namespace md { projections.reserve(num_rels); for(auto i : sorted_rel_indices) { projections.push_back(rel_values[i]); + spd::debug("added push = {}", rel_values[i]); } } template Diagram ModulePresentation::weighted_slice_diagram(const DualPoint& slice) const { + spd::debug("Enter weighted_slice_diagram, slice = {}", slice); IndexVec sorted_gen_indices, sorted_rel_indices; RealVec gen_projections, rel_projections; @@ -138,7 +145,8 @@ namespace md { j = sorted_gen_indices[j]; } std::sort(current_relation.begin(), current_relation.end()); - phat_matrix.set_dim(i, current_relation.size()); + // modules do not have dimension, set all to 0 + phat_matrix.set_dim(i, 0); phat_matrix.set_col(i, current_relation); } @@ -154,6 +162,7 @@ namespace md { bool is_finite_pair = new_pair.second != phat::k_infinity_index; Real birth = gen_projections.at(new_pair.first); Real death = is_finite_pair ? rel_projections.at(new_pair.second) : real_inf; + spd::debug("i = {}, birth = {}, death = {}", i, new_pair.first, new_pair.second); if (birth != death) { dgm.emplace_back(birth, death); } -- cgit v1.2.3