From 1912b703496607999c3ce09ef24e680304a6a093 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Thu, 29 Sep 2016 15:32:19 +0200 Subject: suppress output of 0-persistence pairs in dim 0 --- ripser.cpp | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 3d4eb29..41600b8 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -926,7 +926,8 @@ int main(int argc, char** argv) { if (u != v) { #ifdef PRINT_PERSISTENCE_PAIRS - std::cout << " [0," << get_diameter(e) << ")" << std::endl; + if (get_diameter(e) > 0) + std::cout << " [0," << get_diameter(e) << ")" << std::endl; #endif dset.link(u, v); } else @@ -954,4 +955,4 @@ int main(int argc, char** argv) { assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, comp, dim, n, threshold, binomial_coeff); } } -} \ No newline at end of file +} -- cgit v1.2.3 From 57079b5a02472f4e9e81481d87925b08094a8b53 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Fri, 30 Sep 2016 18:38:40 +0200 Subject: code restructuring --- ripser.cpp | 187 ++++++++++++++++++++------------------- ripser.xcodeproj/project.pbxproj | 17 ++-- 2 files changed, 105 insertions(+), 99 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 41600b8..39a5594 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -98,32 +98,34 @@ std::vector multiplicative_inverse_vector(const coefficient_t m) return inverse; } +index_t get_next_vertex(index_t& v, const index_t idx, const index_t k, const binomial_coeff_table& binomial_coeff) { + if (binomial_coeff(v, k) > idx) { + index_t count = v; + while (count > 0) { + index_t i = v; + index_t step = count >> 1; + i -= step; + if (binomial_coeff(i, k) > idx) { + v = --i; + count -= step + 1; + } else + count = step; + } + } + assert(binomial_coeff(v, k) <= idx); + assert(binomial_coeff(v + 1, k) > idx); + return v; +} + template -OutputIterator get_simplex_vertices(index_t idx, const index_t dim, index_t n, +OutputIterator get_simplex_vertices(index_t idx, const index_t dim, index_t v, const binomial_coeff_table& binomial_coeff, OutputIterator out) { - --n; - + --v; for (index_t k = dim + 1; k > 0; --k) { - if (binomial_coeff(n, k) > idx) { - index_t count = n; - while (count > 0) { - index_t i = n; - index_t step = count >> 1; - i -= step; - if (binomial_coeff(i, k) > idx) { - n = --i; - count -= step + 1; - } else - count = step; - } - } - assert(binomial_coeff(n, k) <= idx); - assert(binomial_coeff(n + 1, k) > idx); - - *out++ = n; - idx -= binomial_coeff(n, k); + get_next_vertex(v, idx, k, binomial_coeff); + *out++ = v; + idx -= binomial_coeff(v, k); } - return out; } @@ -165,7 +167,7 @@ typedef index_t entry_t; const index_t get_index(entry_t i) { return i; } index_t get_coefficient(entry_t i) { return 1; } entry_t make_entry(index_t _index, coefficient_t _value) { return entry_t(_index); } -void set_coefficient(index_t& e, const coefficient_t c) { e = c; } +void set_coefficient(index_t& e, const coefficient_t c) { } #endif @@ -175,7 +177,11 @@ template struct smaller_index { bool operator()(const Entry& a, const Entry& b) { return get_index(a) < get_index(b); } }; -typedef std::pair diameter_index_t; +class diameter_index_t: public std::pair { +public: + diameter_index_t() : std::pair() {} + diameter_index_t(std::pair p) : std::pair(p) {} +}; value_t get_diameter(diameter_index_t i) { return i.first; } index_t get_index(diameter_index_t i) { return i.second; } @@ -184,6 +190,12 @@ public: diameter_entry_t(std::pair p) : std::pair(p) {} diameter_entry_t(entry_t e) : std::pair(0, e) {} diameter_entry_t() : diameter_entry_t(0) {} + diameter_entry_t(value_t _diameter, index_t _index, coefficient_t _coefficient) + : std::pair(_diameter, make_entry(_index, _coefficient)) {} + diameter_entry_t(diameter_index_t _diameter_index, coefficient_t _coefficient) + : std::pair(get_diameter(_diameter_index), + make_entry(get_index(_diameter_index), _coefficient)) {} + diameter_entry_t(diameter_index_t _diameter_index) : diameter_entry_t(_diameter_index, 1) {} }; const entry_t& get_entry(const diameter_entry_t& p) { return p.second; } @@ -192,12 +204,6 @@ const index_t get_index(const diameter_entry_t& p) { return get_index(get_entry( const coefficient_t get_coefficient(const diameter_entry_t& p) { return get_coefficient(get_entry(p)); } const value_t& get_diameter(const diameter_entry_t& p) { return p.first; } void set_coefficient(diameter_entry_t& p, const coefficient_t c) { set_coefficient(get_entry(p), c); } -diameter_entry_t make_diameter_entry(value_t _diameter, index_t _index, coefficient_t _coefficient) { - return std::make_pair(_diameter, make_entry(_index, _coefficient)); -} -diameter_entry_t make_diameter_entry(diameter_index_t _diameter_index, coefficient_t _coefficient) { - return std::make_pair(get_diameter(_diameter_index), make_entry(get_index(_diameter_index), _coefficient)); -} template struct greater_diameter_or_smaller_index { bool operator()(const Entry& a, const Entry& b) { @@ -242,14 +248,23 @@ public: } }; -class simplex_coboundary_enumerator { +template class simplex_coboundary_enumerator { private: index_t idx_below, idx_above, v, k; + std::vector vertices; + const diameter_entry_t simplex; + const coefficient_t modulus; + const DistanceMatrix& dist; const binomial_coeff_table& binomial_coeff; public: - simplex_coboundary_enumerator(index_t _idx, index_t _dim, index_t _n, const binomial_coeff_table& _binomial_coeff) - : idx_below(_idx), idx_above(0), v(_n - 1), k(_dim + 1), binomial_coeff(_binomial_coeff) {} + simplex_coboundary_enumerator(const diameter_entry_t _simplex, index_t _dim, index_t _n, + const coefficient_t _modulus, const DistanceMatrix& _dist, + const binomial_coeff_table& _binomial_coeff) + : simplex(_simplex), idx_below(get_index(_simplex)), idx_above(0), v(_n - 1), k(_dim + 1), modulus(_modulus), + binomial_coeff(_binomial_coeff), dist(_dist), vertices(_dim + 1) { + get_simplex_vertices(get_index(_simplex), _dim, _n, binomial_coeff, vertices.begin()); + } bool has_next() { while ((v != -1) && (binomial_coeff(v, k) <= idx_below)) { @@ -263,10 +278,13 @@ public: return v != -1; } - std::pair next() { - auto result = std::make_pair(make_entry(idx_above + binomial_coeff(v, k + 1) + idx_below, k & 1 ? -1 : 1), v); - --v; - return result; + index_t next_index() { return idx_above + binomial_coeff(v--, k + 1) + idx_below; } + + diameter_entry_t next() { + value_t coface_diameter = get_diameter(simplex); + for (index_t w : vertices) coface_diameter = std::max(coface_diameter, dist(v, w)); + coefficient_t coface_coefficient = (k & 1 ? -1 + modulus : 1) * get_coefficient(simplex) % modulus; + return diameter_entry_t(coface_diameter, idx_above + binomial_coeff(v--, k + 1) + idx_below, coface_coefficient); } }; @@ -484,6 +502,11 @@ void assemble_columns_to_reduce(std::vector& columns_to_reduce if (pivot_column_index.find(index) == pivot_column_index.end()) { value_t diameter = comp.diameter(index); if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); +#ifdef INDICATE_PROGRESS + if ((index + 1) % 1000 == 0) + std::cout << "\033[K" + << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) << "/" << num_simplices << " columns" << std::flush << "\r"; +#endif } } @@ -501,9 +524,10 @@ void assemble_columns_to_reduce(std::vector& columns_to_reduce template void compute_pairs(std::vector& columns_to_reduce, hash_map& pivot_column_index, - const DistanceMatrix& dist, const ComparatorCofaces& comp, const Comparator& comp_prev, index_t dim, - index_t n, value_t threshold, coefficient_t modulus, + index_t dim, index_t n, value_t threshold, coefficient_t modulus, const std::vector& multiplicative_inverse, + const DistanceMatrix& dist, + const ComparatorCofaces& comp, const Comparator& comp_prev, const binomial_coeff_table& binomial_coeff) { #ifdef PRINT_PERSISTENCE_PAIRS @@ -511,7 +535,7 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map reduction_matrix; + compressed_sparse_matrix reduction_coefficients; #else #ifdef USE_COEFFICIENTS std::vector reduction_coefficients; @@ -519,7 +543,6 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map coface_entries; - std::vector vertices; for (index_t i = 0; i < columns_to_reduce.size(); ++i) { auto column_to_reduce = columns_to_reduce[i]; @@ -546,15 +569,15 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map& columns_to_reduce, hash_map cofaces(simplex, dim, n, modulus, dist, binomial_coeff); while (cofaces.has_next()) { - auto coface_descriptor = cofaces.next(); - entry_t coface = coface_descriptor.first; - index_t covertex = coface_descriptor.second; - index_t coface_index = get_index(coface); - value_t coface_diameter = get_diameter(simplex); - for (index_t v : vertices) { coface_diameter = std::max(coface_diameter, dist(v, covertex)); } - assert(comp.diameter(coface_index) == coface_diameter); - - if (coface_diameter <= threshold) { - coefficient_t coface_coefficient = - (get_coefficient(coface) + modulus) * simplex_coefficient % modulus; - assert(coface_coefficient >= 0); - - diameter_entry_t coface_entry = - make_diameter_entry(coface_diameter, coface_index, coface_coefficient); - coface_entries.push_back(coface_entry); - - if (might_be_apparent_pair && (get_diameter(simplex) == coface_diameter)) { - if (pivot_column_index.find(coface_index) == pivot_column_index.end()) { - pivot = coface_entry; + diameter_entry_t coface = cofaces.next(); + if (get_diameter(coface) <= threshold) { + coface_entries.push_back(coface); + if (might_be_apparent_pair && (get_diameter(simplex) == get_diameter(coface))) { + if (pivot_column_index.find(get_index(coface)) == pivot_column_index.end()) { + pivot = coface; goto found_persistence_pair; } might_be_apparent_pair = false; @@ -656,25 +660,22 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map 0); -#else - const coefficient_t coefficient = 1; + set_coefficient(e, inverse * get_coefficient(e) % modulus); + assert(get_coefficient(e) > 0); #endif - reduction_matrix.push_back(make_diameter_entry(get_diameter(e), index, coefficient)); + reduction_coefficients.push_back(e); } #else #ifdef USE_COEFFICIENTS reduction_coefficients.pop_back(); - reduction_coefficients.push_back(make_diameter_entry(column_to_reduce, inverse)); + reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, inverse)); #endif #endif break; @@ -910,7 +911,7 @@ int main(int argc, char** argv) { rips_filtration_comparator comp(dist, 1, binomial_coeff); for (index_t index = binomial_coeff(n, 2); index-- > 0;) { value_t diameter = comp.diameter(index); - if (diameter <= threshold) edges.push_back(diameter_index_t(diameter, index)); + if (diameter <= threshold) edges.push_back(std::make_pair(diameter, index)); } std::sort(edges.rbegin(), edges.rend(), greater_diameter_or_smaller_index()); @@ -948,8 +949,8 @@ int main(int argc, char** argv) { hash_map pivot_column_index; pivot_column_index.reserve(columns_to_reduce.size()); - compute_pairs(columns_to_reduce, pivot_column_index, dist, comp, comp_prev, dim, n, threshold, modulus, - multiplicative_inverse, binomial_coeff); + compute_pairs(columns_to_reduce, pivot_column_index, dim, n, threshold, modulus, multiplicative_inverse, + dist, comp, comp_prev, binomial_coeff); if (dim < dim_max) { assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, comp, dim, n, threshold, binomial_coeff); diff --git a/ripser.xcodeproj/project.pbxproj b/ripser.xcodeproj/project.pbxproj index 04bc3a6..24c7d8b 100644 --- a/ripser.xcodeproj/project.pbxproj +++ b/ripser.xcodeproj/project.pbxproj @@ -88,7 +88,7 @@ 551018401BD63CB300990BFF /* Project object */ = { isa = PBXProject; attributes = { - LastUpgradeCheck = 0700; + LastUpgradeCheck = 0800; ORGANIZATIONNAME = "ulrich-bauer.org"; TargetAttributes = { 551018471BD63CB300990BFF = { @@ -129,7 +129,7 @@ isa = XCBuildConfiguration; buildSettings = { ALWAYS_SEARCH_USER_PATHS = NO; - CLANG_CXX_LANGUAGE_STANDARD = "gnu++0x"; + CLANG_CXX_LANGUAGE_STANDARD = "c++0x"; CLANG_CXX_LIBRARY = "libc++"; CLANG_ENABLE_MODULES = YES; CLANG_ENABLE_OBJC_ARC = YES; @@ -138,15 +138,17 @@ CLANG_WARN_DIRECT_OBJC_ISA_USAGE = YES_ERROR; CLANG_WARN_EMPTY_BODY = YES; CLANG_WARN_ENUM_CONVERSION = YES; + CLANG_WARN_INFINITE_RECURSION = YES; CLANG_WARN_INT_CONVERSION = YES; CLANG_WARN_OBJC_ROOT_CLASS = YES_ERROR; + CLANG_WARN_SUSPICIOUS_MOVE = YES; CLANG_WARN_UNREACHABLE_CODE = YES; CLANG_WARN__DUPLICATE_METHOD_MATCH = YES; COPY_PHASE_STRIP = NO; DEBUG_INFORMATION_FORMAT = dwarf; ENABLE_STRICT_OBJC_MSGSEND = YES; ENABLE_TESTABILITY = YES; - GCC_C_LANGUAGE_STANDARD = gnu99; + GCC_C_LANGUAGE_STANDARD = c11; GCC_DYNAMIC_NO_PIC = NO; GCC_NO_COMMON_BLOCKS = YES; GCC_OPTIMIZATION_LEVEL = 0; @@ -172,7 +174,7 @@ isa = XCBuildConfiguration; buildSettings = { ALWAYS_SEARCH_USER_PATHS = NO; - CLANG_CXX_LANGUAGE_STANDARD = "gnu++0x"; + CLANG_CXX_LANGUAGE_STANDARD = "c++0x"; CLANG_CXX_LIBRARY = "libc++"; CLANG_ENABLE_MODULES = YES; CLANG_ENABLE_OBJC_ARC = YES; @@ -181,15 +183,17 @@ CLANG_WARN_DIRECT_OBJC_ISA_USAGE = YES_ERROR; CLANG_WARN_EMPTY_BODY = YES; CLANG_WARN_ENUM_CONVERSION = YES; + CLANG_WARN_INFINITE_RECURSION = YES; CLANG_WARN_INT_CONVERSION = YES; CLANG_WARN_OBJC_ROOT_CLASS = YES_ERROR; + CLANG_WARN_SUSPICIOUS_MOVE = YES; CLANG_WARN_UNREACHABLE_CODE = YES; CLANG_WARN__DUPLICATE_METHOD_MATCH = YES; COPY_PHASE_STRIP = NO; DEBUG_INFORMATION_FORMAT = "dwarf-with-dsym"; ENABLE_NS_ASSERTIONS = NO; ENABLE_STRICT_OBJC_MSGSEND = YES; - GCC_C_LANGUAGE_STANDARD = gnu99; + GCC_C_LANGUAGE_STANDARD = c11; GCC_NO_COMMON_BLOCKS = YES; GCC_WARN_64_TO_32_BIT_CONVERSION = YES; GCC_WARN_ABOUT_RETURN_TYPE = YES_ERROR; @@ -212,7 +216,8 @@ _FILE_FORMAT_DISTANCE_MATRIX, _FILE_FORMAT_LOWER_DISTANCE_MATRIX, _FILE_FORMAT_POINT_CLOUD, - USE_COEFFICIENTS, + _USE_COEFFICIENTS, + INDICATE_PROGRESS, ); PRODUCT_NAME = "$(TARGET_NAME)"; }; -- cgit v1.2.3 From dea2c3d290ffbab6441f0e03864db34fcd3ef9d5 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Fri, 30 Sep 2016 18:43:25 +0200 Subject: pretty-print --- ripser.cpp | 28 +++++++++++++--------------- 1 file changed, 13 insertions(+), 15 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 39a5594..9d57fa9 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -18,7 +18,6 @@ with this program. If not, see . */ - //#define ASSEMBLE_REDUCTION_MATRIX //#define USE_COEFFICIENTS @@ -167,7 +166,7 @@ typedef index_t entry_t; const index_t get_index(entry_t i) { return i; } index_t get_coefficient(entry_t i) { return 1; } entry_t make_entry(index_t _index, coefficient_t _value) { return entry_t(_index); } -void set_coefficient(index_t& e, const coefficient_t c) { } +void set_coefficient(index_t& e, const coefficient_t c) {} #endif @@ -177,7 +176,7 @@ template struct smaller_index { bool operator()(const Entry& a, const Entry& b) { return get_index(a) < get_index(b); } }; -class diameter_index_t: public std::pair { +class diameter_index_t : public std::pair { public: diameter_index_t() : std::pair() {} diameter_index_t(std::pair p) : std::pair(p) {} @@ -284,7 +283,8 @@ public: value_t coface_diameter = get_diameter(simplex); for (index_t w : vertices) coface_diameter = std::max(coface_diameter, dist(v, w)); coefficient_t coface_coefficient = (k & 1 ? -1 + modulus : 1) * get_coefficient(simplex) % modulus; - return diameter_entry_t(coface_diameter, idx_above + binomial_coeff(v--, k + 1) + idx_below, coface_coefficient); + return diameter_entry_t(coface_diameter, idx_above + binomial_coeff(v--, k + 1) + idx_below, + coface_coefficient); } }; @@ -503,9 +503,10 @@ void assemble_columns_to_reduce(std::vector& columns_to_reduce value_t diameter = comp.diameter(index); if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); #ifdef INDICATE_PROGRESS - if ((index + 1) % 1000 == 0) - std::cout << "\033[K" - << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) << "/" << num_simplices << " columns" << std::flush << "\r"; + if ((index + 1) % 1000 == 0) + std::cout << "\033[K" + << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) << "/" + << num_simplices << " columns" << std::flush << "\r"; #endif } } @@ -525,8 +526,7 @@ void assemble_columns_to_reduce(std::vector& columns_to_reduce template void compute_pairs(std::vector& columns_to_reduce, hash_map& pivot_column_index, index_t dim, index_t n, value_t threshold, coefficient_t modulus, - const std::vector& multiplicative_inverse, - const DistanceMatrix& dist, + const std::vector& multiplicative_inverse, const DistanceMatrix& dist, const ComparatorCofaces& comp, const Comparator& comp_prev, const binomial_coeff_table& binomial_coeff) { @@ -596,8 +596,7 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map 0) - std::cout << " [0," << get_diameter(e) << ")" << std::endl; + if (get_diameter(e) > 0) std::cout << " [0," << get_diameter(e) << ")" << std::endl; #endif dset.link(u, v); } else @@ -949,8 +947,8 @@ int main(int argc, char** argv) { hash_map pivot_column_index; pivot_column_index.reserve(columns_to_reduce.size()); - compute_pairs(columns_to_reduce, pivot_column_index, dim, n, threshold, modulus, multiplicative_inverse, - dist, comp, comp_prev, binomial_coeff); + compute_pairs(columns_to_reduce, pivot_column_index, dim, n, threshold, modulus, multiplicative_inverse, dist, + comp, comp_prev, binomial_coeff); if (dim < dim_max) { assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, comp, dim, n, threshold, binomial_coeff); -- cgit v1.2.3 From f46a9148452480916fc1fda4b501e2af7e5bd570 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Mon, 24 Oct 2016 00:23:07 +0200 Subject: updated readme with reference to live.ripser.org --- README.md | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/README.md b/README.md index 0a3f304..0c5a2e6 100644 --- a/README.md +++ b/README.md @@ -7,6 +7,8 @@ Copyright © 2015–2016 [Ulrich Bauer]. Ripser is a lean C++ code for the computation of Vietoris–Rips persistence barcodes. It can do just this one thing, but does it extremely well. +To see a live demo of Ripser's capabilities, go to [live.ripser.org]. The computation happens inside the browser (using [PNaCl] on Chrome and JavaScript via [Emscripten] on other browsers). + The main features of Ripser: - time- and memory-efficient @@ -14,7 +16,7 @@ The main features of Ripser: - support for coefficients in prime finite fields - no external dependencies (optional support for Google's [sparsehash]) -Currently, Ripser outperforms other codes ([Dionysus], [DIPHA], [GUDHI], [Perseus], [PHAT]) by a factor of more than 40 in computation time and a factor of more than 15 in memory efficiency. (Note that [PHAT] does not contain code for generating Vietoris–Rips filtrations). +Currently, Ripser outperforms other codes ([Dionysus], [DIPHA], [GUDHI], [Perseus], [PHAT]) by a factor of more than 40 in computation time and a factor of more than 15 in memory efficiency (for the example linked at [live.ripser.org]). (Note that [PHAT] does not contain code for generating Vietoris–Rips filtrations). Input formats currently supported by Ripser: @@ -57,7 +59,7 @@ Ripser supports several compile-time options. They are switched on by defining t - `USE_COEFFICIENTS`: enable support for coefficients in a prime field - `INDICATE_PROGRESS`: indicate the current progress in the console - `PRINT_PERSISTENCE_PAIRS`: output the computed persistence pairs (enabled by default in the code; comment out to disable) - - `USE_GOOGLE_HASHMAP`: enable support for Google's [sparsehash] data structure; may further reducue memory footprint + - `USE_GOOGLE_HASHMAP`: enable support for Google's [sparsehash] data structure; may further reduce memory footprint For example, to build Ripser with support for coefficients: @@ -89,13 +91,17 @@ The following features are currently planned for future versions: - computation of representative cycles for persistent homology (currenly only *co*cycles are computed) - support for sparse distance matrices +Prototype implementations are already avaliable; please contact the author if one of these features might be relevant for your research. + ### License Ripser is licensed under the [LGPL] 3.0. Please contact the author if you want to use Ripser in your software under a different license. - [Ulrich Bauer]: +[live.ripser.org]: +[PNaCl]: +[Emscripten]: [latest-release]: [Dionysus]: [DIPHA]: -- cgit v1.2.3 From 370d31554bf52353cfc2af542242faf621066588 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Mon, 7 Nov 2016 10:55:29 +0100 Subject: javaplex benchmarks --- benchmarks/benchmarks.txt | 1 + benchmarks/javaplex.txt | 19 +++++++++++++++++++ 2 files changed, 20 insertions(+) create mode 100644 benchmarks/javaplex.txt diff --git a/benchmarks/benchmarks.txt b/benchmarks/benchmarks.txt index c5fc91f..bb34d08 100644 --- a/benchmarks/benchmarks.txt +++ b/benchmarks/benchmarks.txt @@ -8,3 +8,4 @@ c++ -std=c++11 ripser.cpp -o ripser -Ofast -D NDEBUG -D FILE_FORMAT_DIPHA -D PRINT_PERSISTENCE_PAIRS && /usr/bin/time -l ./ripser --dim 2 ~/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex +/usr/bin/time -l java -Xmx16G -cp Bitbucket/phat-paper/benchmark:Source/javaplex/dist/javaplex-4.2.1.jar RipsFiltration ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2 3 \ No newline at end of file diff --git a/benchmarks/javaplex.txt b/benchmarks/javaplex.txt new file mode 100644 index 0000000..a39e804 --- /dev/null +++ b/benchmarks/javaplex.txt @@ -0,0 +1,19 @@ +geometry74:~ uli$ /usr/bin/time -l java -Xmx16G -cp Bitbucket/phat-paper/benchmark:Source/javaplex/dist/javaplex-4.2.1.jar RipsFiltration ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2 3 +Constructing Rips filtration in 425.57800000000003 s +Computing persistence (default) in 2848.477 s + 3275.48 real 3900.79 user 36.35 sys +12636958720 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 4791446 page reclaims + 2 page faults + 0 swaps + 5 block input operations + 5 block output operations + 0 messages sent + 0 messages received + 5 signals received + 9 voluntary context switches + 9898205 involuntary context switches + -- cgit v1.2.3 From 5c55dd581f101b68941150ab803cb79eef6a79e0 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Thu, 10 Nov 2016 08:43:53 +0100 Subject: renamed variable --- ripser.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ripser.cpp b/ripser.cpp index 9d57fa9..8bd5a83 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -619,7 +619,7 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map Date: Thu, 10 Nov 2016 18:21:18 +0100 Subject: Ripser binary file format --- ripser.cpp | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/ripser.cpp b/ripser.cpp index 8bd5a83..3e4cf1f 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -686,7 +686,7 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map T read(std::istream& s) { T result; @@ -787,6 +787,12 @@ compressed_lower_distance_matrix read_dipha(std::istream& input_stream) { return compressed_lower_distance_matrix(std::move(distances)); } +compressed_lower_distance_matrix read_ripser(std::istream& input_stream) { + std::vector distances; + while (!input_stream.eof()) distances.push_back(read(input_stream)); + return compressed_lower_distance_matrix(std::move(distances)); +} + compressed_lower_distance_matrix read_file(std::istream& input_stream, file_format format) { switch (format) { case LOWER_DISTANCE_MATRIX: @@ -799,6 +805,8 @@ compressed_lower_distance_matrix read_file(std::istream& input_stream, file_form return read_point_cloud(input_stream); case DIPHA: return read_dipha(input_stream); + case RIPSER: + return read_ripser(input_stream); } } @@ -816,6 +824,7 @@ void print_usage_and_exit(int exit_code) { << " distance (full distance matrix)" << std::endl << " point-cloud (point cloud in Euclidean space)" << std::endl << " dipha (distance matrix in DIPHA file format)" << std::endl + << " ripser (distance matrix in Ripser binary file format)" << std::endl << " --dim compute persistent homology up to dimension " << std::endl << " --threshold compute Rips complexes up to diameter " << std::endl #ifdef USE_COEFFICIENTS @@ -867,6 +876,8 @@ int main(int argc, char** argv) { format = POINT_CLOUD; else if (parameter == "dipha") format = DIPHA; + else if (parameter == "ripser") + format = RIPSER; else print_usage_and_exit(-1); #ifdef USE_COEFFICIENTS -- cgit v1.2.3 From 5cfc80f27c19e2db52373a155b8dbf1ff627ac26 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Tue, 22 Nov 2016 17:18:54 -0500 Subject: coeff reduction in makefile --- Makefile | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/Makefile b/Makefile index 73e09fb..893d776 100644 --- a/Makefile +++ b/Makefile @@ -1,7 +1,7 @@ build: ripser -all: ripser ripser-coeff ripser-reduction ripser-debug +all: ripser ripser-coeff ripser-reduction ripser-coeff-reduction ripser-debug ripser: ripser.cpp @@ -13,9 +13,12 @@ ripser-coeff: ripser.cpp ripser-reduction: ripser.cpp c++ -std=c++11 ripser.cpp -o ripser-reduction -Ofast -D NDEBUG -D ASSEMBLE_REDUCTION_MATRIX +ripser-coeff-reduction: ripser.cpp + c++ -std=c++11 ripser.cpp -o ripser-coeff-reduction -Ofast -D NDEBUG -D USE_COEFFICIENTS -D ASSEMBLE_REDUCTION_MATRIX + ripser-debug: ripser.cpp c++ -std=c++11 ripser.cpp -o ripser-debug -g clean: - rm -f ripser ripser-coeff ripser-reduction ripser-debug + rm -f ripser ripser-coeff ripser-reduction ripser-coeff-reduction ripser-debug -- cgit v1.2.3 From c0574ec97f6d413b597447636a39b550cb55ee5d Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Tue, 22 Nov 2016 17:19:27 -0500 Subject: don’t store diagonal in reduction matrix when not using coefficients MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- ripser.cpp | 20 +++++++++++++++----- 1 file changed, 15 insertions(+), 5 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 3e4cf1f..022b5fe 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -548,7 +548,7 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map, smaller_index> + std::priority_queue, greater_diameter_or_smaller_index> reduction_column; #endif @@ -574,20 +574,25 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map coeffs(0); + coeffs.push_back(columns_to_reduce[j]); + for (auto it = reduction_coefficients.cbegin(j); it != reduction_coefficients.cend(j); ++it) coeffs.push_back(*it); + auto coeffs_begin = coeffs.begin(), coeffs_end = coeffs.end(); +#endif #else #ifdef USE_COEFFICIENTS auto coeffs_begin = &reduction_coefficients[j], coeffs_end = &reduction_coefficients[j] + 1; @@ -661,7 +666,12 @@ void compute_pairs(std::vector& columns_to_reduce, hash_map Date: Tue, 22 Nov 2016 18:26:01 -0500 Subject: moved attribute packed --- ripser.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 022b5fe..8f6ec89 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -136,13 +136,13 @@ std::vector vertices_of_simplex(const index_t simplex_index, const inde } #ifdef USE_COEFFICIENTS -struct entry_t { +struct __attribute__((packed)) entry_t { index_t index : 8 * (sizeof(index_t) - sizeof(coefficient_t)); coefficient_t coefficient; entry_t(index_t _index, coefficient_t _coefficient) : index(_index), coefficient(_coefficient) {} entry_t(index_t _index) : index(_index), coefficient(1) {} entry_t() : index(0), coefficient(1) {} -} __attribute__((packed)); +}; static_assert(sizeof(entry_t) == sizeof(index_t), "size of entry_t is not the same as index_t"); -- cgit v1.2.3 From 6f36a1288f8de8c9cddc8b90ca8256935291ce28 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Sat, 26 Nov 2016 17:20:41 -0500 Subject: code reorganization --- ripser.cpp | 433 ++++++++++++++++++++++++++++--------------------------------- 1 file changed, 198 insertions(+), 235 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 8f6ec89..dc6b1bd 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -47,8 +47,6 @@ template class hash_map : public std::unordered_map #endif typedef float value_t; -// typedef uint16_t value_t; - typedef int64_t index_t; typedef int16_t coefficient_t; @@ -74,8 +72,7 @@ public: } index_t operator()(index_t n, index_t k) const { - assert(n <= n_max); - assert(k <= k_max); + assert(n <= n_max && k <= k_max); return B[n][k]; } }; @@ -111,8 +108,7 @@ index_t get_next_vertex(index_t& v, const index_t idx, const index_t k, const bi count = step; } } - assert(binomial_coeff(v, k) <= idx); - assert(binomial_coeff(v + 1, k) > idx); + assert(binomial_coeff(v, k) <= idx && binomial_coeff(v + 1, k) > idx); return v; } @@ -128,13 +124,6 @@ OutputIterator get_simplex_vertices(index_t idx, const index_t dim, index_t v, return out; } -std::vector vertices_of_simplex(const index_t simplex_index, const index_t dim, const index_t n, - const binomial_coeff_table& binomial_coeff) { - std::vector vertices; - get_simplex_vertices(simplex_index, dim, n, binomial_coeff, std::back_inserter(vertices)); - return vertices; -} - #ifdef USE_COEFFICIENTS struct __attribute__((packed)) entry_t { index_t index : 8 * (sizeof(index_t) - sizeof(coefficient_t)); @@ -172,10 +161,6 @@ void set_coefficient(index_t& e, const coefficient_t c) {} const entry_t& get_entry(const entry_t& e) { return e; } -template struct smaller_index { - bool operator()(const Entry& a, const Entry& b) { return get_index(a) < get_index(b); } -}; - class diameter_index_t : public std::pair { public: diameter_index_t() : std::pair() {} @@ -211,42 +196,6 @@ template struct greater_diameter_or_smaller_index { } }; -template class rips_filtration_comparator { -public: - const DistanceMatrix& dist; - const index_t dim; - -private: - mutable std::vector vertices; - const binomial_coeff_table& binomial_coeff; - -public: - rips_filtration_comparator(const DistanceMatrix& _dist, const index_t _dim, - const binomial_coeff_table& _binomial_coeff) - : dist(_dist), dim(_dim), vertices(_dim + 1), binomial_coeff(_binomial_coeff){}; - - value_t diameter(const index_t index) const { - value_t diam = 0; - get_simplex_vertices(index, dim, dist.size(), binomial_coeff, vertices.begin()); - - for (index_t i = 0; i <= dim; ++i) - for (index_t j = 0; j < i; ++j) { diam = std::max(diam, dist(vertices[i], vertices[j])); } - return diam; - } - - bool operator()(const index_t a, const index_t b) const { - assert(a < binomial_coeff(dist.size(), dim + 1)); - assert(b < binomial_coeff(dist.size(), dim + 1)); - - return greater_diameter_or_smaller_index()(diameter_index_t(diameter(a), a), - diameter_index_t(diameter(b), b)); - } - - template bool operator()(const Entry& a, const Entry& b) const { - return operator()(get_index(a), get_index(b)); - } -}; - template class simplex_coboundary_enumerator { private: index_t idx_below, idx_above, v, k; @@ -269,7 +218,6 @@ public: while ((v != -1) && (binomial_coeff(v, k) <= idx_below)) { idx_below -= binomial_coeff(v, k); idx_above += binomial_coeff(v, k + 1); - --v; --k; assert(k != -1); @@ -485,216 +433,289 @@ template void push_entry(Heap& column, index_t i, coefficient_t column.push(std::make_pair(diameter, e)); } -template -void assemble_columns_to_reduce(std::vector& columns_to_reduce, - hash_map& pivot_column_index, const Comparator& comp, index_t dim, - index_t n, value_t threshold, const binomial_coeff_table& binomial_coeff) { - index_t num_simplices = binomial_coeff(n, dim + 2); +class ripser { + index_t dim_max, n; + value_t threshold; + const binomial_coeff_table binomial_coeff; + std::vector multiplicative_inverse; + coefficient_t modulus; + compressed_lower_distance_matrix dist; + mutable std::vector vertices; + +public: + ripser(compressed_lower_distance_matrix&& _dist, index_t _dim_max, value_t _threshold, coefficient_t _modulus) + : dist(_dist), n(_dist.size()), dim_max(std::min(_dim_max, index_t(_dist.size() - 2))), threshold(_threshold), + modulus(_modulus), binomial_coeff(n, dim_max + 2), + multiplicative_inverse(multiplicative_inverse_vector(_modulus)) {} + + value_t compute_diameter(const index_t index, index_t dim) const { + value_t diam = 0; + + vertices.clear(); + get_simplex_vertices(index, dim, dist.size(), binomial_coeff, std::back_inserter(vertices)); - columns_to_reduce.clear(); + for (index_t i = 0; i <= dim; ++i) + for (index_t j = 0; j < i; ++j) { diam = std::max(diam, dist(vertices[i], vertices[j])); } + return diam; + } + + void assemble_columns_to_reduce(std::vector& columns_to_reduce, + hash_map& pivot_column_index, index_t dim) { + index_t num_simplices = binomial_coeff(n, dim + 1); + + columns_to_reduce.clear(); #ifdef INDICATE_PROGRESS - std::cout << "\033[K" - << "assembling " << num_simplices << " columns" << std::flush << "\r"; + std::cout << "\033[K" + << "assembling " << num_simplices << " columns" << std::flush << "\r"; #endif - for (index_t index = 0; index < num_simplices; ++index) { - if (pivot_column_index.find(index) == pivot_column_index.end()) { - value_t diameter = comp.diameter(index); - if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); + for (index_t index = 0; index < num_simplices; ++index) { + if (pivot_column_index.find(index) == pivot_column_index.end()) { + value_t diameter = compute_diameter(index, dim); + if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); #ifdef INDICATE_PROGRESS - if ((index + 1) % 1000 == 0) - std::cout << "\033[K" - << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) << "/" - << num_simplices << " columns" << std::flush << "\r"; + if ((index + 1) % 1000 == 0) + std::cout << "\033[K" + << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) << "/" + << num_simplices << " columns" << std::flush << "\r"; #endif + } } - } #ifdef INDICATE_PROGRESS - std::cout << "\033[K" - << "sorting " << num_simplices << " columns" << std::flush << "\r"; + std::cout << "\033[K" + << "sorting " << num_simplices << " columns" << std::flush << "\r"; #endif - std::sort(columns_to_reduce.begin(), columns_to_reduce.end(), - greater_diameter_or_smaller_index()); + std::sort(columns_to_reduce.begin(), columns_to_reduce.end(), + greater_diameter_or_smaller_index()); #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif -} + } -template -void compute_pairs(std::vector& columns_to_reduce, hash_map& pivot_column_index, - index_t dim, index_t n, value_t threshold, coefficient_t modulus, - const std::vector& multiplicative_inverse, const DistanceMatrix& dist, - const ComparatorCofaces& comp, const Comparator& comp_prev, - const binomial_coeff_table& binomial_coeff) { + void compute_pairs(std::vector& columns_to_reduce, hash_map& pivot_column_index, + index_t dim) { #ifdef PRINT_PERSISTENCE_PAIRS - std::cout << "persistence intervals in dim " << dim << ":" << std::endl; + std::cout << "persistence intervals in dim " << dim << ":" << std::endl; #endif #ifdef ASSEMBLE_REDUCTION_MATRIX - compressed_sparse_matrix reduction_coefficients; + compressed_sparse_matrix reduction_coefficients; #else #ifdef USE_COEFFICIENTS - std::vector reduction_coefficients; + std::vector reduction_coefficients; #endif #endif - std::vector coface_entries; + std::vector coface_entries; - for (index_t i = 0; i < columns_to_reduce.size(); ++i) { - auto column_to_reduce = columns_to_reduce[i]; + for (index_t i = 0; i < columns_to_reduce.size(); ++i) { + auto column_to_reduce = columns_to_reduce[i]; #ifdef ASSEMBLE_REDUCTION_MATRIX - std::priority_queue, greater_diameter_or_smaller_index> - reduction_column; + std::priority_queue, + greater_diameter_or_smaller_index> + reduction_column; #endif - std::priority_queue, - greater_diameter_or_smaller_index> - working_coboundary; + std::priority_queue, + greater_diameter_or_smaller_index> + working_coboundary; - value_t diameter = get_diameter(column_to_reduce); + value_t diameter = get_diameter(column_to_reduce); #ifdef INDICATE_PROGRESS - if ((i + 1) % 1000 == 0) - std::cout << "\033[K" - << "reducing column " << i + 1 << "/" << columns_to_reduce.size() << " (diameter " << diameter - << ")" << std::flush << "\r"; + if ((i + 1) % 1000 == 0) + std::cout << "\033[K" + << "reducing column " << i + 1 << "/" << columns_to_reduce.size() << " (diameter " << diameter + << ")" << std::flush << "\r"; #endif - index_t j = i; + index_t j = i; - // start with a dummy pivot entry with coefficient -1 in order to initialize - // working_coboundary with the coboundary of the simplex with index column_to_reduce - diameter_entry_t pivot(0, -1, -1 + modulus); + // start with a dummy pivot entry with coefficient -1 in order to initialize + // working_coboundary with the coboundary of the simplex with index column_to_reduce + diameter_entry_t pivot(0, -1, -1 + modulus); #ifdef ASSEMBLE_REDUCTION_MATRIX - // initialize reduction_coefficients as identity matrix - reduction_coefficients.append_column(); + // initialize reduction_coefficients as identity matrix + reduction_coefficients.append_column(); #endif #ifdef USE_COEFFICIENTS - reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, 1)); + reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, 1)); #endif - - bool might_be_apparent_pair = (i == j); - do { - const coefficient_t factor = modulus - get_coefficient(pivot); + bool might_be_apparent_pair = (i == j); + + do { + const coefficient_t factor = modulus - get_coefficient(pivot); #ifdef ASSEMBLE_REDUCTION_MATRIX #ifdef USE_COEFFICIENTS - auto coeffs_begin = reduction_coefficients.cbegin(j), coeffs_end = reduction_coefficients.cend(j); + auto coeffs_begin = reduction_coefficients.cbegin(j), coeffs_end = reduction_coefficients.cend(j); #else - std::vector coeffs(0); - coeffs.push_back(columns_to_reduce[j]); - for (auto it = reduction_coefficients.cbegin(j); it != reduction_coefficients.cend(j); ++it) coeffs.push_back(*it); - auto coeffs_begin = coeffs.begin(), coeffs_end = coeffs.end(); + std::vector coeffs(0); + coeffs.push_back(columns_to_reduce[j]); + for (auto it = reduction_coefficients.cbegin(j); it != reduction_coefficients.cend(j); ++it) + coeffs.push_back(*it); + auto coeffs_begin = coeffs.begin(), coeffs_end = coeffs.end(); #endif #else #ifdef USE_COEFFICIENTS - auto coeffs_begin = &reduction_coefficients[j], coeffs_end = &reduction_coefficients[j] + 1; + auto coeffs_begin = &reduction_coefficients[j], coeffs_end = &reduction_coefficients[j] + 1; #else - auto coeffs_begin = &columns_to_reduce[j], coeffs_end = &columns_to_reduce[j] + 1; + auto coeffs_begin = &columns_to_reduce[j], coeffs_end = &columns_to_reduce[j] + 1; #endif #endif - for (auto it = coeffs_begin; it != coeffs_end; ++it) { - diameter_entry_t simplex = *it; - set_coefficient(simplex, get_coefficient(simplex) * factor % modulus); + for (auto it = coeffs_begin; it != coeffs_end; ++it) { + diameter_entry_t simplex = *it; + set_coefficient(simplex, get_coefficient(simplex) * factor % modulus); #ifdef ASSEMBLE_REDUCTION_MATRIX - reduction_column.push(simplex); -#endif - - coface_entries.clear(); - simplex_coboundary_enumerator cofaces(simplex, dim, n, modulus, dist, binomial_coeff); - while (cofaces.has_next()) { - diameter_entry_t coface = cofaces.next(); - if (get_diameter(coface) <= threshold) { - coface_entries.push_back(coface); - if (might_be_apparent_pair && (get_diameter(simplex) == get_diameter(coface))) { - if (pivot_column_index.find(get_index(coface)) == pivot_column_index.end()) { - pivot = coface; - goto found_persistence_pair; + reduction_column.push(simplex); +#endif + + coface_entries.clear(); + simplex_coboundary_enumerator cofaces(simplex, dim, n, modulus, dist, + binomial_coeff); + while (cofaces.has_next()) { + diameter_entry_t coface = cofaces.next(); + if (get_diameter(coface) <= threshold) { + coface_entries.push_back(coface); + if (might_be_apparent_pair && (get_diameter(simplex) == get_diameter(coface))) { + if (pivot_column_index.find(get_index(coface)) == pivot_column_index.end()) { + pivot = coface; + goto found_persistence_pair; + } + might_be_apparent_pair = false; } - might_be_apparent_pair = false; } } + for (auto coface : coface_entries) working_coboundary.push(coface); } - for (auto coface : coface_entries) working_coboundary.push(coface); - } - pivot = get_pivot(working_coboundary, modulus); + pivot = get_pivot(working_coboundary, modulus); - if (get_index(pivot) != -1) { - auto pair = pivot_column_index.find(get_index(pivot)); + if (get_index(pivot) != -1) { + auto pair = pivot_column_index.find(get_index(pivot)); - if (pair != pivot_column_index.end()) { - j = pair->second; - continue; - } - } else { + if (pair != pivot_column_index.end()) { + j = pair->second; + continue; + } + } else { #ifdef PRINT_PERSISTENCE_PAIRS #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif - std::cout << " [" << diameter << ", )" << std::endl << std::flush; + std::cout << " [" << diameter << ", )" << std::endl << std::flush; #endif - break; - } + break; + } - found_persistence_pair: + found_persistence_pair: #ifdef PRINT_PERSISTENCE_PAIRS - value_t death = get_diameter(pivot); - if (diameter != death) { + value_t death = get_diameter(pivot); + if (diameter != death) { #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif - std::cout << " [" << diameter << "," << death << ")" << std::endl << std::flush; - } + std::cout << " [" << diameter << "," << death << ")" << std::endl << std::flush; + } #endif - pivot_column_index.insert(std::make_pair(get_index(pivot), i)); + pivot_column_index.insert(std::make_pair(get_index(pivot), i)); #ifdef USE_COEFFICIENTS - const coefficient_t inverse = multiplicative_inverse[get_coefficient(pivot)]; + const coefficient_t inverse = multiplicative_inverse[get_coefficient(pivot)]; #endif #ifdef ASSEMBLE_REDUCTION_MATRIX - // replace current column of reduction_coefficients (with a single diagonal 1 entry) - // by reduction_column (possibly with a different entry on the diagonal) +// replace current column of reduction_coefficients (with a single diagonal 1 entry) +// by reduction_column (possibly with a different entry on the diagonal) #ifdef USE_COEFFICIENTS - reduction_coefficients.pop_back(); + reduction_coefficients.pop_back(); #else - pop_pivot(reduction_column, modulus); + pop_pivot(reduction_column, modulus); #endif - - while (true) { - diameter_entry_t e = pop_pivot(reduction_column, modulus); - if (get_index(e) == -1) break; + + while (true) { + diameter_entry_t e = pop_pivot(reduction_column, modulus); + if (get_index(e) == -1) break; #ifdef USE_COEFFICIENTS - set_coefficient(e, inverse * get_coefficient(e) % modulus); - assert(get_coefficient(e) > 0); + set_coefficient(e, inverse * get_coefficient(e) % modulus); + assert(get_coefficient(e) > 0); #endif - reduction_coefficients.push_back(e); - } + reduction_coefficients.push_back(e); + } #else #ifdef USE_COEFFICIENTS - reduction_coefficients.pop_back(); - reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, inverse)); + reduction_coefficients.pop_back(); + reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, inverse)); #endif #endif - break; - } while (true); - } + break; + } while (true); + } #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif -} + } + + void compute_barcodes() { + + std::vector columns_to_reduce; + + { + union_find dset(n); + std::vector edges; + for (index_t index = binomial_coeff(n, 2); index-- > 0;) { + value_t diameter = compute_diameter(index, 1); + if (diameter <= threshold) edges.push_back(std::make_pair(diameter, index)); + } + std::sort(edges.rbegin(), edges.rend(), greater_diameter_or_smaller_index()); + +#ifdef PRINT_PERSISTENCE_PAIRS + std::cout << "persistence intervals in dim 0:" << std::endl; +#endif + + std::vector vertices_of_edge(2); + for (auto e : edges) { + vertices_of_edge.clear(); + get_simplex_vertices(get_index(e), 1, n, binomial_coeff, std::back_inserter(vertices_of_edge)); + index_t u = dset.find(vertices_of_edge[0]), v = dset.find(vertices_of_edge[1]); + + if (u != v) { +#ifdef PRINT_PERSISTENCE_PAIRS + if (get_diameter(e) > 0) std::cout << " [0," << get_diameter(e) << ")" << std::endl; +#endif + dset.link(u, v); + } else + columns_to_reduce.push_back(e); + } + std::reverse(columns_to_reduce.begin(), columns_to_reduce.end()); + +#ifdef PRINT_PERSISTENCE_PAIRS + for (index_t i = 0; i < n; ++i) + if (dset.find(i) == i) std::cout << " [0, )" << std::endl << std::flush; +#endif + } + + for (index_t dim = 1; dim <= dim_max; ++dim) { + hash_map pivot_column_index; + pivot_column_index.reserve(columns_to_reduce.size()); + + compute_pairs(columns_to_reduce, pivot_column_index, dim); + + if (dim < dim_max) { assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, dim + 1); } + } + } +}; enum file_format { LOWER_DISTANCE_MATRIX, UPPER_DISTANCE_MATRIX, DISTANCE_MATRIX, POINT_CLOUD, DIPHA, RIPSER }; @@ -911,68 +932,10 @@ int main(int argc, char** argv) { compressed_lower_distance_matrix dist = read_file(filename ? file_stream : std::cin, format); - index_t n = dist.size(); - - std::cout << "distance matrix with " << n << " points" << std::endl; + std::cout << "distance matrix with " << dist.size() << " points" << std::endl; auto value_range = std::minmax_element(dist.distances.begin(), dist.distances.end()); std::cout << "value range: [" << *value_range.first << "," << *value_range.second << "]" << std::endl; - dim_max = std::min(dim_max, n - 2); - - binomial_coeff_table binomial_coeff(n, dim_max + 2); - std::vector multiplicative_inverse(multiplicative_inverse_vector(modulus)); - - std::vector columns_to_reduce; - - { - union_find dset(n); - std::vector edges; - rips_filtration_comparator comp(dist, 1, binomial_coeff); - for (index_t index = binomial_coeff(n, 2); index-- > 0;) { - value_t diameter = comp.diameter(index); - if (diameter <= threshold) edges.push_back(std::make_pair(diameter, index)); - } - std::sort(edges.rbegin(), edges.rend(), greater_diameter_or_smaller_index()); - -#ifdef PRINT_PERSISTENCE_PAIRS - std::cout << "persistence intervals in dim 0:" << std::endl; -#endif - - std::vector vertices_of_edge(2); - for (auto e : edges) { - vertices_of_edge.clear(); - get_simplex_vertices(get_index(e), 1, n, binomial_coeff, std::back_inserter(vertices_of_edge)); - index_t u = dset.find(vertices_of_edge[0]), v = dset.find(vertices_of_edge[1]); - - if (u != v) { -#ifdef PRINT_PERSISTENCE_PAIRS - if (get_diameter(e) > 0) std::cout << " [0," << get_diameter(e) << ")" << std::endl; -#endif - dset.link(u, v); - } else - columns_to_reduce.push_back(e); - } - std::reverse(columns_to_reduce.begin(), columns_to_reduce.end()); - -#ifdef PRINT_PERSISTENCE_PAIRS - for (index_t i = 0; i < n; ++i) - if (dset.find(i) == i) std::cout << " [0, )" << std::endl << std::flush; -#endif - } - - for (index_t dim = 1; dim <= dim_max; ++dim) { - rips_filtration_comparator comp(dist, dim + 1, binomial_coeff); - rips_filtration_comparator comp_prev(dist, dim, binomial_coeff); - - hash_map pivot_column_index; - pivot_column_index.reserve(columns_to_reduce.size()); - - compute_pairs(columns_to_reduce, pivot_column_index, dim, n, threshold, modulus, multiplicative_inverse, dist, - comp, comp_prev, binomial_coeff); - - if (dim < dim_max) { - assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, comp, dim, n, threshold, binomial_coeff); - } - } + ripser(std::move(dist), dim_max, threshold, modulus).compute_barcodes(); } -- cgit v1.2.3 From f177efa1b4ceac3cc72a3702ebc1dd508559dcd5 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Sat, 26 Nov 2016 17:20:59 -0500 Subject: more code reorganization --- ripser.cpp | 145 ++++++++++++++++++++++++++++++------------------------------- 1 file changed, 71 insertions(+), 74 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index dc6b1bd..6209730 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -94,35 +94,6 @@ std::vector multiplicative_inverse_vector(const coefficient_t m) return inverse; } -index_t get_next_vertex(index_t& v, const index_t idx, const index_t k, const binomial_coeff_table& binomial_coeff) { - if (binomial_coeff(v, k) > idx) { - index_t count = v; - while (count > 0) { - index_t i = v; - index_t step = count >> 1; - i -= step; - if (binomial_coeff(i, k) > idx) { - v = --i; - count -= step + 1; - } else - count = step; - } - } - assert(binomial_coeff(v, k) <= idx && binomial_coeff(v + 1, k) > idx); - return v; -} - -template -OutputIterator get_simplex_vertices(index_t idx, const index_t dim, index_t v, - const binomial_coeff_table& binomial_coeff, OutputIterator out) { - --v; - for (index_t k = dim + 1; k > 0; --k) { - get_next_vertex(v, idx, k, binomial_coeff); - *out++ = v; - idx -= binomial_coeff(v, k); - } - return out; -} #ifdef USE_COEFFICIENTS struct __attribute__((packed)) entry_t { @@ -196,46 +167,6 @@ template struct greater_diameter_or_smaller_index { } }; -template class simplex_coboundary_enumerator { -private: - index_t idx_below, idx_above, v, k; - std::vector vertices; - const diameter_entry_t simplex; - const coefficient_t modulus; - const DistanceMatrix& dist; - const binomial_coeff_table& binomial_coeff; - -public: - simplex_coboundary_enumerator(const diameter_entry_t _simplex, index_t _dim, index_t _n, - const coefficient_t _modulus, const DistanceMatrix& _dist, - const binomial_coeff_table& _binomial_coeff) - : simplex(_simplex), idx_below(get_index(_simplex)), idx_above(0), v(_n - 1), k(_dim + 1), modulus(_modulus), - binomial_coeff(_binomial_coeff), dist(_dist), vertices(_dim + 1) { - get_simplex_vertices(get_index(_simplex), _dim, _n, binomial_coeff, vertices.begin()); - } - - bool has_next() { - while ((v != -1) && (binomial_coeff(v, k) <= idx_below)) { - idx_below -= binomial_coeff(v, k); - idx_above += binomial_coeff(v, k + 1); - --v; - --k; - assert(k != -1); - } - return v != -1; - } - - index_t next_index() { return idx_above + binomial_coeff(v--, k + 1) + idx_below; } - - diameter_entry_t next() { - value_t coface_diameter = get_diameter(simplex); - for (index_t w : vertices) coface_diameter = std::max(coface_diameter, dist(v, w)); - coefficient_t coface_coefficient = (k & 1 ? -1 + modulus : 1) * get_coefficient(simplex) % modulus; - return diameter_entry_t(coface_diameter, idx_above + binomial_coeff(v--, k + 1) + idx_below, - coface_coefficient); - } -}; - enum compressed_matrix_layout { LOWER_TRIANGULAR, UPPER_TRIANGULAR }; template class compressed_distance_matrix { @@ -447,17 +378,84 @@ public: : dist(_dist), n(_dist.size()), dim_max(std::min(_dim_max, index_t(_dist.size() - 2))), threshold(_threshold), modulus(_modulus), binomial_coeff(n, dim_max + 2), multiplicative_inverse(multiplicative_inverse_vector(_modulus)) {} - + + index_t get_next_vertex(index_t& v, const index_t idx, const index_t k) const { + if (binomial_coeff(v, k) > idx) { + index_t count = v; + while (count > 0) { + index_t i = v; + index_t step = count >> 1; + i -= step; + if (binomial_coeff(i, k) > idx) { + v = --i; + count -= step + 1; + } else + count = step; + } + } + assert(binomial_coeff(v, k) <= idx && binomial_coeff(v + 1, k) > idx); + return v; + } + + template + OutputIterator get_simplex_vertices(index_t idx, const index_t dim, index_t v, + OutputIterator out) const { + --v; + for (index_t k = dim + 1; k > 0; --k) { + get_next_vertex(v, idx, k); + *out++ = v; + idx -= binomial_coeff(v, k); + } + return out; + } + value_t compute_diameter(const index_t index, index_t dim) const { value_t diam = 0; vertices.clear(); - get_simplex_vertices(index, dim, dist.size(), binomial_coeff, std::back_inserter(vertices)); + get_simplex_vertices(index, dim, dist.size(), std::back_inserter(vertices)); for (index_t i = 0; i <= dim; ++i) for (index_t j = 0; j < i; ++j) { diam = std::max(diam, dist(vertices[i], vertices[j])); } return diam; } + + class simplex_coboundary_enumerator { + private: + index_t idx_below, idx_above, v, k; + std::vector vertices; + const diameter_entry_t simplex; + const coefficient_t modulus; + const compressed_lower_distance_matrix& dist; + const binomial_coeff_table& binomial_coeff; + + public: + simplex_coboundary_enumerator(const diameter_entry_t _simplex, index_t _dim, const ripser& parent) + : simplex(_simplex), idx_below(get_index(_simplex)), idx_above(0), v(parent.n - 1), k(_dim + 1), modulus(parent.modulus), + binomial_coeff(parent.binomial_coeff), dist(parent.dist), vertices(_dim + 1) { + parent.get_simplex_vertices(get_index(_simplex), _dim, parent.n, vertices.begin()); + } + + bool has_next() { + while ((v != -1) && (binomial_coeff(v, k) <= idx_below)) { + idx_below -= binomial_coeff(v, k); + idx_above += binomial_coeff(v, k + 1); + --v; + --k; + assert(k != -1); + } + return v != -1; + } + + diameter_entry_t next() { + value_t coface_diameter = get_diameter(simplex); + for (index_t w : vertices) coface_diameter = std::max(coface_diameter, dist(v, w)); + index_t coface_index = idx_above + binomial_coeff(v--, k + 1) + idx_below; + coefficient_t coface_coefficient = (k & 1 ? -1 + modulus : 1) * get_coefficient(simplex) % modulus; + return diameter_entry_t(coface_diameter, coface_index, coface_coefficient); + } + }; + void assemble_columns_to_reduce(std::vector& columns_to_reduce, hash_map& pivot_column_index, index_t dim) { @@ -580,8 +578,7 @@ public: #endif coface_entries.clear(); - simplex_coboundary_enumerator cofaces(simplex, dim, n, modulus, dist, - binomial_coeff); + simplex_coboundary_enumerator cofaces(simplex, dim, *this); while (cofaces.has_next()) { diameter_entry_t coface = cofaces.next(); if (get_diameter(coface) <= threshold) { @@ -687,7 +684,7 @@ public: std::vector vertices_of_edge(2); for (auto e : edges) { vertices_of_edge.clear(); - get_simplex_vertices(get_index(e), 1, n, binomial_coeff, std::back_inserter(vertices_of_edge)); + get_simplex_vertices(get_index(e), 1, n, std::back_inserter(vertices_of_edge)); index_t u = dset.find(vertices_of_edge[0]), v = dset.find(vertices_of_edge[1]); if (u != v) { -- cgit v1.2.3 From eb8c26dc3bad09678290c58d3637cb8336a37704 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Sat, 26 Nov 2016 18:16:27 -0500 Subject: cleanup binomial coefficients --- ripser.cpp | 13 +++---------- 1 file changed, 3 insertions(+), 10 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 6209730..72e784c 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -52,27 +52,20 @@ typedef int16_t coefficient_t; class binomial_coeff_table { std::vector> B; - index_t n_max, k_max; - public: - binomial_coeff_table(index_t n, index_t k) { - n_max = n; - k_max = k; - - B.resize(n + 1); + binomial_coeff_table(index_t n, index_t k) : B(n + 1) { for (index_t i = 0; i <= n; i++) { B[i].resize(k + 1); - for (index_t j = 0; j <= std::min(i, k); j++) { + for (index_t j = 0; j <= std::min(i, k); j++) if (j == 0 || j == i) B[i][j] = 1; else B[i][j] = B[i - 1][j - 1] + B[i - 1][j]; - } } } index_t operator()(index_t n, index_t k) const { - assert(n <= n_max && k <= k_max); + assert(n < B.size() && k < B[n].size()); return B[n][k]; } }; -- cgit v1.2.3 From 8521d24c4584a6eb1149a604309ffa987862414f Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Mon, 28 Nov 2016 14:13:46 -0500 Subject: allow for negative edge filtration values --- ripser.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 72e784c..efa672e 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -403,7 +403,7 @@ public: } value_t compute_diameter(const index_t index, index_t dim) const { - value_t diam = 0; + value_t diam = -std::numeric_limits::infinity(); vertices.clear(); get_simplex_vertices(index, dim, dist.size(), std::back_inserter(vertices)); @@ -682,7 +682,7 @@ public: if (u != v) { #ifdef PRINT_PERSISTENCE_PAIRS - if (get_diameter(e) > 0) std::cout << " [0," << get_diameter(e) << ")" << std::endl; + if (get_diameter(e) != 0) std::cout << " [0," << get_diameter(e) << ")" << std::endl; #endif dset.link(u, v); } else -- cgit v1.2.3 From c9596839245779792ddf7c8235ab96feecc8aeb8 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Wed, 30 Nov 2016 14:42:47 -0500 Subject: fix for missing constructor of diameter_entry_t --- ripser.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ripser.cpp b/ripser.cpp index efa672e..403516b 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -548,7 +548,7 @@ public: #ifdef USE_COEFFICIENTS auto coeffs_begin = reduction_coefficients.cbegin(j), coeffs_end = reduction_coefficients.cend(j); #else - std::vector coeffs(0); + std::vector coeffs; coeffs.push_back(columns_to_reduce[j]); for (auto it = reduction_coefficients.cbegin(j); it != reduction_coefficients.cend(j); ++it) coeffs.push_back(*it); -- cgit v1.2.3 From d32402d37998dd0d2c24be8e5a139ebc525a0a56 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Mon, 5 Dec 2016 10:51:16 +0100 Subject: fixed move constructor --- ripser.cpp | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 403516b..aad196d 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -170,7 +170,7 @@ public: void init_rows(); compressed_distance_matrix(std::vector&& _distances) - : distances(_distances), rows((1 + std::sqrt(1 + 8 * distances.size())) / 2) { + : distances(std::move(_distances)), rows((1 + std::sqrt(1 + 8 * distances.size())) / 2) { assert(distances.size() == size() * (size() - 1) / 2); init_rows(); } @@ -222,7 +222,7 @@ class euclidean_distance_matrix { public: std::vector> points; - euclidean_distance_matrix(std::vector>&& _points) : points(_points) {} + euclidean_distance_matrix(std::vector>&& _points) : points(std::move(_points)) {} value_t operator()(const index_t i, const index_t j) const { return std::sqrt(std::inner_product(points[i].begin(), points[i].end(), points[j].begin(), value_t(), @@ -358,17 +358,17 @@ template void push_entry(Heap& column, index_t i, coefficient_t } class ripser { + compressed_lower_distance_matrix dist; index_t dim_max, n; value_t threshold; + coefficient_t modulus; const binomial_coeff_table binomial_coeff; std::vector multiplicative_inverse; - coefficient_t modulus; - compressed_lower_distance_matrix dist; mutable std::vector vertices; public: ripser(compressed_lower_distance_matrix&& _dist, index_t _dim_max, value_t _threshold, coefficient_t _modulus) - : dist(_dist), n(_dist.size()), dim_max(std::min(_dim_max, index_t(_dist.size() - 2))), threshold(_threshold), + : dist(std::move(_dist)), dim_max(std::min(_dim_max, index_t(dist.size() - 2))), n(dist.size()), threshold(_threshold), modulus(_modulus), binomial_coeff(n, dim_max + 2), multiplicative_inverse(multiplicative_inverse_vector(_modulus)) {} -- cgit v1.2.3 From 8b42fb2a7885b0e371ff6b0e1bbe31332ae63adb Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Sat, 17 Dec 2016 14:05:20 +0100 Subject: cleanup compressed distance matrix access --- ripser.cpp | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index aad196d..0f9d02f 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -206,13 +206,11 @@ template <> void compressed_distance_matrix::init_rows() { } template <> value_t compressed_distance_matrix::operator()(index_t i, index_t j) const { - if (i > j) std::swap(i, j); - return i == j ? 0 : rows[i][j]; + return i == j ? 0 : i > j ? rows[j][i] : rows[i][j]; } template <> value_t compressed_distance_matrix::operator()(index_t i, index_t j) const { - if (i > j) std::swap(i, j); - return i == j ? 0 : rows[j][i]; + return i == j ? 0 : i < j ? rows[j][i] : rows[i][j]; } typedef compressed_distance_matrix compressed_lower_distance_matrix; -- cgit v1.2.3 From 2158b4a884d1239ede305b7ce91c9f9beeb0668b Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Thu, 29 Dec 2016 16:23:53 +0100 Subject: reordered initializers --- ripser.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 0f9d02f..052c07d 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -422,8 +422,8 @@ public: public: simplex_coboundary_enumerator(const diameter_entry_t _simplex, index_t _dim, const ripser& parent) - : simplex(_simplex), idx_below(get_index(_simplex)), idx_above(0), v(parent.n - 1), k(_dim + 1), modulus(parent.modulus), - binomial_coeff(parent.binomial_coeff), dist(parent.dist), vertices(_dim + 1) { + : idx_below(get_index(_simplex)), idx_above(0), v(parent.n - 1), k(_dim + 1), vertices(_dim + 1), simplex(_simplex), modulus(parent.modulus), dist(parent.dist), + binomial_coeff(parent.binomial_coeff) { parent.get_simplex_vertices(get_index(_simplex), _dim, parent.n, vertices.begin()); } -- cgit v1.2.3 From 11dfa016218793929b6a8b7f8c93a3669f3229d5 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Thu, 29 Dec 2016 17:52:49 +0100 Subject: code reformatting --- .clang-format | 2 +- ripser.cpp | 469 ++++++++++++++++++++++++++++++++-------------------------- 2 files changed, 259 insertions(+), 212 deletions(-) diff --git a/.clang-format b/.clang-format index 598d7c2..1975b19 100644 --- a/.clang-format +++ b/.clang-format @@ -1,7 +1,7 @@ BasedOnStyle: LLVM IndentWidth: 4 TabWidth: 4 -ColumnLimit: 120 +ColumnLimit: 100 AccessModifierOffset: -4 AllowShortIfStatementsOnASingleLine: true AllowShortLoopsOnASingleLine: true diff --git a/ripser.cpp b/ripser.cpp index 052c07d..02bfac9 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -52,6 +52,7 @@ typedef int16_t coefficient_t; class binomial_coeff_table { std::vector> B; + public: binomial_coeff_table(index_t n, index_t k) : B(n + 1) { for (index_t i = 0; i <= n; i++) { @@ -87,19 +88,21 @@ std::vector multiplicative_inverse_vector(const coefficient_t m) return inverse; } - #ifdef USE_COEFFICIENTS struct __attribute__((packed)) entry_t { index_t index : 8 * (sizeof(index_t) - sizeof(coefficient_t)); coefficient_t coefficient; - entry_t(index_t _index, coefficient_t _coefficient) : index(_index), coefficient(_coefficient) {} + entry_t(index_t _index, coefficient_t _coefficient) + : index(_index), coefficient(_coefficient) {} entry_t(index_t _index) : index(_index), coefficient(1) {} entry_t() : index(0), coefficient(1) {} }; static_assert(sizeof(entry_t) == sizeof(index_t), "size of entry_t is not the same as index_t"); -entry_t make_entry(index_t _index, coefficient_t _coefficient) { return entry_t(_index, _coefficient); } +entry_t make_entry(index_t _index, coefficient_t _coefficient) { + return entry_t(_index, _coefficient); +} index_t get_index(entry_t e) { return e.index; } index_t get_coefficient(entry_t e) { return e.coefficient; } void set_coefficient(entry_t& e, const coefficient_t c) { e.coefficient = c; } @@ -149,9 +152,13 @@ public: const entry_t& get_entry(const diameter_entry_t& p) { return p.second; } entry_t& get_entry(diameter_entry_t& p) { return p.second; } const index_t get_index(const diameter_entry_t& p) { return get_index(get_entry(p)); } -const coefficient_t get_coefficient(const diameter_entry_t& p) { return get_coefficient(get_entry(p)); } +const coefficient_t get_coefficient(const diameter_entry_t& p) { + return get_coefficient(get_entry(p)); +} const value_t& get_diameter(const diameter_entry_t& p) { return p.first; } -void set_coefficient(diameter_entry_t& p, const coefficient_t c) { set_coefficient(get_entry(p), c); } +void set_coefficient(diameter_entry_t& p, const coefficient_t c) { + set_coefficient(get_entry(p), c); +} template struct greater_diameter_or_smaller_index { bool operator()(const Entry& a, const Entry& b) { @@ -205,11 +212,13 @@ template <> void compressed_distance_matrix::init_rows() { } } -template <> value_t compressed_distance_matrix::operator()(index_t i, index_t j) const { +template <> +value_t compressed_distance_matrix::operator()(index_t i, index_t j) const { return i == j ? 0 : i > j ? rows[j][i] : rows[i][j]; } -template <> value_t compressed_distance_matrix::operator()(index_t i, index_t j) const { +template <> +value_t compressed_distance_matrix::operator()(index_t i, index_t j) const { return i == j ? 0 : i < j ? rows[j][i] : rows[i][j]; } @@ -220,12 +229,13 @@ class euclidean_distance_matrix { public: std::vector> points; - euclidean_distance_matrix(std::vector>&& _points) : points(std::move(_points)) {} + euclidean_distance_matrix(std::vector>&& _points) + : points(std::move(_points)) {} value_t operator()(const index_t i, const index_t j) const { - return std::sqrt(std::inner_product(points[i].begin(), points[i].end(), points[j].begin(), value_t(), - std::plus(), - [](value_t u, value_t v) { return (u - v) * (u - v); })); + return std::sqrt(std::inner_product( + points[i].begin(), points[i].end(), points[j].begin(), value_t(), std::plus(), + [](value_t u, value_t v) { return (u - v) * (u - v); })); } size_t size() const { return points.size(); } @@ -350,7 +360,8 @@ public: } }; -template void push_entry(Heap& column, index_t i, coefficient_t c, value_t diameter) { +template +void push_entry(Heap& column, index_t i, coefficient_t c, value_t diameter) { entry_t e = make_entry(i, c); column.push(std::make_pair(diameter, e)); } @@ -365,11 +376,12 @@ class ripser { mutable std::vector vertices; public: - ripser(compressed_lower_distance_matrix&& _dist, index_t _dim_max, value_t _threshold, coefficient_t _modulus) - : dist(std::move(_dist)), dim_max(std::min(_dim_max, index_t(dist.size() - 2))), n(dist.size()), threshold(_threshold), - modulus(_modulus), binomial_coeff(n, dim_max + 2), + ripser(compressed_lower_distance_matrix&& _dist, index_t _dim_max, value_t _threshold, + coefficient_t _modulus) + : dist(std::move(_dist)), dim_max(std::min(_dim_max, index_t(dist.size() - 2))), + n(dist.size()), threshold(_threshold), modulus(_modulus), binomial_coeff(n, dim_max + 2), multiplicative_inverse(multiplicative_inverse_vector(_modulus)) {} - + index_t get_next_vertex(index_t& v, const index_t idx, const index_t k) const { if (binomial_coeff(v, k) > idx) { index_t count = v; @@ -387,10 +399,10 @@ public: assert(binomial_coeff(v, k) <= idx && binomial_coeff(v + 1, k) > idx); return v; } - + template OutputIterator get_simplex_vertices(index_t idx, const index_t dim, index_t v, - OutputIterator out) const { + OutputIterator out) const { --v; for (index_t k = dim + 1; k > 0; --k) { get_next_vertex(v, idx, k); @@ -399,7 +411,7 @@ public: } return out; } - + value_t compute_diameter(const index_t index, index_t dim) const { value_t diam = -std::numeric_limits::infinity(); @@ -407,10 +419,12 @@ public: get_simplex_vertices(index, dim, dist.size(), std::back_inserter(vertices)); for (index_t i = 0; i <= dim; ++i) - for (index_t j = 0; j < i; ++j) { diam = std::max(diam, dist(vertices[i], vertices[j])); } + for (index_t j = 0; j < i; ++j) { + diam = std::max(diam, dist(vertices[i], vertices[j])); + } return diam; } - + class simplex_coboundary_enumerator { private: index_t idx_below, idx_above, v, k; @@ -419,14 +433,16 @@ public: const coefficient_t modulus; const compressed_lower_distance_matrix& dist; const binomial_coeff_table& binomial_coeff; - + public: - simplex_coboundary_enumerator(const diameter_entry_t _simplex, index_t _dim, const ripser& parent) - : idx_below(get_index(_simplex)), idx_above(0), v(parent.n - 1), k(_dim + 1), vertices(_dim + 1), simplex(_simplex), modulus(parent.modulus), dist(parent.dist), - binomial_coeff(parent.binomial_coeff) { + simplex_coboundary_enumerator(const diameter_entry_t _simplex, index_t _dim, + const ripser& parent) + : idx_below(get_index(_simplex)), idx_above(0), v(parent.n - 1), k(_dim + 1), + vertices(_dim + 1), simplex(_simplex), modulus(parent.modulus), dist(parent.dist), + binomial_coeff(parent.binomial_coeff) { parent.get_simplex_vertices(get_index(_simplex), _dim, parent.n, vertices.begin()); } - + bool has_next() { while ((v != -1) && (binomial_coeff(v, k) <= idx_below)) { idx_below -= binomial_coeff(v, k); @@ -437,276 +453,248 @@ public: } return v != -1; } - + diameter_entry_t next() { value_t coface_diameter = get_diameter(simplex); for (index_t w : vertices) coface_diameter = std::max(coface_diameter, dist(v, w)); index_t coface_index = idx_above + binomial_coeff(v--, k + 1) + idx_below; - coefficient_t coface_coefficient = (k & 1 ? -1 + modulus : 1) * get_coefficient(simplex) % modulus; + coefficient_t coface_coefficient = + (k & 1 ? -1 + modulus : 1) * get_coefficient(simplex) % modulus; return diameter_entry_t(coface_diameter, coface_index, coface_coefficient); } }; - void assemble_columns_to_reduce(std::vector& columns_to_reduce, - hash_map& pivot_column_index, index_t dim) { - index_t num_simplices = binomial_coeff(n, dim + 1); + hash_map& pivot_column_index, index_t dim); + + void compute_pairs(std::vector& columns_to_reduce, + hash_map& pivot_column_index, index_t dim); + + void compute_barcodes(); +}; - columns_to_reduce.clear(); +void ripser::assemble_columns_to_reduce(std::vector& columns_to_reduce, + hash_map& pivot_column_index, + index_t dim) { + index_t num_simplices = binomial_coeff(n, dim + 1); + + columns_to_reduce.clear(); #ifdef INDICATE_PROGRESS - std::cout << "\033[K" - << "assembling " << num_simplices << " columns" << std::flush << "\r"; + std::cout << "\033[K" + << "assembling " << num_simplices << " columns" << std::flush << "\r"; #endif - for (index_t index = 0; index < num_simplices; ++index) { - if (pivot_column_index.find(index) == pivot_column_index.end()) { - value_t diameter = compute_diameter(index, dim); - if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); + for (index_t index = 0; index < num_simplices; ++index) { + if (pivot_column_index.find(index) == pivot_column_index.end()) { + value_t diameter = compute_diameter(index, dim); + if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); #ifdef INDICATE_PROGRESS - if ((index + 1) % 1000 == 0) - std::cout << "\033[K" - << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) << "/" - << num_simplices << " columns" << std::flush << "\r"; + if ((index + 1) % 1000 == 0) + std::cout << "\033[K" + << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) + << "/" << num_simplices << " columns" << std::flush << "\r"; #endif - } } + } #ifdef INDICATE_PROGRESS - std::cout << "\033[K" - << "sorting " << num_simplices << " columns" << std::flush << "\r"; + std::cout << "\033[K" + << "sorting " << num_simplices << " columns" << std::flush << "\r"; #endif - std::sort(columns_to_reduce.begin(), columns_to_reduce.end(), - greater_diameter_or_smaller_index()); + std::sort(columns_to_reduce.begin(), columns_to_reduce.end(), + greater_diameter_or_smaller_index()); #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif - } +} - void compute_pairs(std::vector& columns_to_reduce, hash_map& pivot_column_index, - index_t dim) { +void ripser::compute_pairs(std::vector& columns_to_reduce, + hash_map& pivot_column_index, index_t dim) { #ifdef PRINT_PERSISTENCE_PAIRS - std::cout << "persistence intervals in dim " << dim << ":" << std::endl; + std::cout << "persistence intervals in dim " << dim << ":" << std::endl; #endif #ifdef ASSEMBLE_REDUCTION_MATRIX - compressed_sparse_matrix reduction_coefficients; + compressed_sparse_matrix reduction_coefficients; #else #ifdef USE_COEFFICIENTS - std::vector reduction_coefficients; + std::vector reduction_coefficients; #endif #endif - std::vector coface_entries; + std::vector coface_entries; - for (index_t i = 0; i < columns_to_reduce.size(); ++i) { - auto column_to_reduce = columns_to_reduce[i]; + for (index_t i = 0; i < columns_to_reduce.size(); ++i) { + auto column_to_reduce = columns_to_reduce[i]; #ifdef ASSEMBLE_REDUCTION_MATRIX - std::priority_queue, - greater_diameter_or_smaller_index> - reduction_column; + std::priority_queue, + greater_diameter_or_smaller_index> + reduction_column; #endif - std::priority_queue, - greater_diameter_or_smaller_index> - working_coboundary; + std::priority_queue, + greater_diameter_or_smaller_index> + working_coboundary; - value_t diameter = get_diameter(column_to_reduce); + value_t diameter = get_diameter(column_to_reduce); #ifdef INDICATE_PROGRESS - if ((i + 1) % 1000 == 0) - std::cout << "\033[K" - << "reducing column " << i + 1 << "/" << columns_to_reduce.size() << " (diameter " << diameter - << ")" << std::flush << "\r"; + if ((i + 1) % 1000 == 0) + std::cout << "\033[K" + << "reducing column " << i + 1 << "/" << columns_to_reduce.size() + << " (diameter " << diameter << ")" << std::flush << "\r"; #endif - index_t j = i; + index_t j = i; - // start with a dummy pivot entry with coefficient -1 in order to initialize - // working_coboundary with the coboundary of the simplex with index column_to_reduce - diameter_entry_t pivot(0, -1, -1 + modulus); + // start with a dummy pivot entry with coefficient -1 in order to initialize + // working_coboundary with the coboundary of the simplex with index column_to_reduce + diameter_entry_t pivot(0, -1, -1 + modulus); #ifdef ASSEMBLE_REDUCTION_MATRIX - // initialize reduction_coefficients as identity matrix - reduction_coefficients.append_column(); + // initialize reduction_coefficients as identity matrix + reduction_coefficients.append_column(); #endif #ifdef USE_COEFFICIENTS - reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, 1)); + reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, 1)); #endif - bool might_be_apparent_pair = (i == j); + bool might_be_apparent_pair = (i == j); - do { - const coefficient_t factor = modulus - get_coefficient(pivot); + do { + const coefficient_t factor = modulus - get_coefficient(pivot); #ifdef ASSEMBLE_REDUCTION_MATRIX #ifdef USE_COEFFICIENTS - auto coeffs_begin = reduction_coefficients.cbegin(j), coeffs_end = reduction_coefficients.cend(j); + auto coeffs_begin = reduction_coefficients.cbegin(j), + coeffs_end = reduction_coefficients.cend(j); #else - std::vector coeffs; - coeffs.push_back(columns_to_reduce[j]); - for (auto it = reduction_coefficients.cbegin(j); it != reduction_coefficients.cend(j); ++it) - coeffs.push_back(*it); - auto coeffs_begin = coeffs.begin(), coeffs_end = coeffs.end(); + std::vector coeffs; + coeffs.push_back(columns_to_reduce[j]); + for (auto it = reduction_coefficients.cbegin(j); it != reduction_coefficients.cend(j); + ++it) + coeffs.push_back(*it); + auto coeffs_begin = coeffs.begin(), coeffs_end = coeffs.end(); #endif #else #ifdef USE_COEFFICIENTS - auto coeffs_begin = &reduction_coefficients[j], coeffs_end = &reduction_coefficients[j] + 1; + auto coeffs_begin = &reduction_coefficients[j], + coeffs_end = &reduction_coefficients[j] + 1; #else - auto coeffs_begin = &columns_to_reduce[j], coeffs_end = &columns_to_reduce[j] + 1; + auto coeffs_begin = &columns_to_reduce[j], coeffs_end = &columns_to_reduce[j] + 1; #endif #endif - for (auto it = coeffs_begin; it != coeffs_end; ++it) { - diameter_entry_t simplex = *it; - set_coefficient(simplex, get_coefficient(simplex) * factor % modulus); + for (auto it = coeffs_begin; it != coeffs_end; ++it) { + diameter_entry_t simplex = *it; + set_coefficient(simplex, get_coefficient(simplex) * factor % modulus); #ifdef ASSEMBLE_REDUCTION_MATRIX - reduction_column.push(simplex); -#endif - - coface_entries.clear(); - simplex_coboundary_enumerator cofaces(simplex, dim, *this); - while (cofaces.has_next()) { - diameter_entry_t coface = cofaces.next(); - if (get_diameter(coface) <= threshold) { - coface_entries.push_back(coface); - if (might_be_apparent_pair && (get_diameter(simplex) == get_diameter(coface))) { - if (pivot_column_index.find(get_index(coface)) == pivot_column_index.end()) { - pivot = coface; - goto found_persistence_pair; - } - might_be_apparent_pair = false; + reduction_column.push(simplex); +#endif + + coface_entries.clear(); + simplex_coboundary_enumerator cofaces(simplex, dim, *this); + while (cofaces.has_next()) { + diameter_entry_t coface = cofaces.next(); + if (get_diameter(coface) <= threshold) { + coface_entries.push_back(coface); + if (might_be_apparent_pair && + (get_diameter(simplex) == get_diameter(coface))) { + if (pivot_column_index.find(get_index(coface)) == + pivot_column_index.end()) { + pivot = coface; + goto found_persistence_pair; } + might_be_apparent_pair = false; } } - for (auto coface : coface_entries) working_coboundary.push(coface); } + for (auto coface : coface_entries) working_coboundary.push(coface); + } - pivot = get_pivot(working_coboundary, modulus); + pivot = get_pivot(working_coboundary, modulus); - if (get_index(pivot) != -1) { - auto pair = pivot_column_index.find(get_index(pivot)); + if (get_index(pivot) != -1) { + auto pair = pivot_column_index.find(get_index(pivot)); - if (pair != pivot_column_index.end()) { - j = pair->second; - continue; - } - } else { + if (pair != pivot_column_index.end()) { + j = pair->second; + continue; + } + } else { #ifdef PRINT_PERSISTENCE_PAIRS #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif - std::cout << " [" << diameter << ", )" << std::endl << std::flush; + std::cout << " [" << diameter << ", )" << std::endl << std::flush; #endif - break; - } + break; + } - found_persistence_pair: + found_persistence_pair: #ifdef PRINT_PERSISTENCE_PAIRS - value_t death = get_diameter(pivot); - if (diameter != death) { + value_t death = get_diameter(pivot); + if (diameter != death) { #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif - std::cout << " [" << diameter << "," << death << ")" << std::endl << std::flush; - } + std::cout << " [" << diameter << "," << death << ")" << std::endl << std::flush; + } #endif - pivot_column_index.insert(std::make_pair(get_index(pivot), i)); + pivot_column_index.insert(std::make_pair(get_index(pivot), i)); #ifdef USE_COEFFICIENTS - const coefficient_t inverse = multiplicative_inverse[get_coefficient(pivot)]; + const coefficient_t inverse = multiplicative_inverse[get_coefficient(pivot)]; #endif #ifdef ASSEMBLE_REDUCTION_MATRIX // replace current column of reduction_coefficients (with a single diagonal 1 entry) // by reduction_column (possibly with a different entry on the diagonal) #ifdef USE_COEFFICIENTS - reduction_coefficients.pop_back(); + reduction_coefficients.pop_back(); #else - pop_pivot(reduction_column, modulus); + pop_pivot(reduction_column, modulus); #endif - while (true) { - diameter_entry_t e = pop_pivot(reduction_column, modulus); - if (get_index(e) == -1) break; + while (true) { + diameter_entry_t e = pop_pivot(reduction_column, modulus); + if (get_index(e) == -1) break; #ifdef USE_COEFFICIENTS - set_coefficient(e, inverse * get_coefficient(e) % modulus); - assert(get_coefficient(e) > 0); + set_coefficient(e, inverse * get_coefficient(e) % modulus); + assert(get_coefficient(e) > 0); #endif - reduction_coefficients.push_back(e); - } + reduction_coefficients.push_back(e); + } #else #ifdef USE_COEFFICIENTS - reduction_coefficients.pop_back(); - reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, inverse)); + reduction_coefficients.pop_back(); + reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, inverse)); #endif #endif - break; - } while (true); - } - -#ifdef INDICATE_PROGRESS - std::cout << "\033[K"; -#endif + break; + } while (true); } - void compute_barcodes() { - - std::vector columns_to_reduce; - - { - union_find dset(n); - std::vector edges; - for (index_t index = binomial_coeff(n, 2); index-- > 0;) { - value_t diameter = compute_diameter(index, 1); - if (diameter <= threshold) edges.push_back(std::make_pair(diameter, index)); - } - std::sort(edges.rbegin(), edges.rend(), greater_diameter_or_smaller_index()); - -#ifdef PRINT_PERSISTENCE_PAIRS - std::cout << "persistence intervals in dim 0:" << std::endl; -#endif - - std::vector vertices_of_edge(2); - for (auto e : edges) { - vertices_of_edge.clear(); - get_simplex_vertices(get_index(e), 1, n, std::back_inserter(vertices_of_edge)); - index_t u = dset.find(vertices_of_edge[0]), v = dset.find(vertices_of_edge[1]); - - if (u != v) { -#ifdef PRINT_PERSISTENCE_PAIRS - if (get_diameter(e) != 0) std::cout << " [0," << get_diameter(e) << ")" << std::endl; -#endif - dset.link(u, v); - } else - columns_to_reduce.push_back(e); - } - std::reverse(columns_to_reduce.begin(), columns_to_reduce.end()); - -#ifdef PRINT_PERSISTENCE_PAIRS - for (index_t i = 0; i < n; ++i) - if (dset.find(i) == i) std::cout << " [0, )" << std::endl << std::flush; +#ifdef INDICATE_PROGRESS + std::cout << "\033[K"; #endif - } - - for (index_t dim = 1; dim <= dim_max; ++dim) { - hash_map pivot_column_index; - pivot_column_index.reserve(columns_to_reduce.size()); - - compute_pairs(columns_to_reduce, pivot_column_index, dim); +} - if (dim < dim_max) { assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, dim + 1); } - } - } +enum file_format { + LOWER_DISTANCE_MATRIX, + UPPER_DISTANCE_MATRIX, + DISTANCE_MATRIX, + POINT_CLOUD, + DIPHA, + RIPSER }; -enum file_format { LOWER_DISTANCE_MATRIX, UPPER_DISTANCE_MATRIX, DISTANCE_MATRIX, POINT_CLOUD, DIPHA, RIPSER }; - template T read(std::istream& s) { T result; s.read(reinterpret_cast(&result), sizeof(T)); @@ -733,7 +721,8 @@ compressed_lower_distance_matrix read_point_cloud(std::istream& input_stream) { index_t n = eucl_dist.size(); - std::cout << "point cloud with " << n << " points in dimension " << eucl_dist.points.front().size() << std::endl; + std::cout << "point cloud with " << n << " points in dimension " + << eucl_dist.points.front().size() << std::endl; std::vector distances; @@ -830,26 +819,30 @@ compressed_lower_distance_matrix read_file(std::istream& input_stream, file_form } void print_usage_and_exit(int exit_code) { - std::cerr << "Usage: " - << "ripser " - << "[options] [filename]" << std::endl - << std::endl - << "Options:" << std::endl - << std::endl - << " --help print this screen" << std::endl - << " --format use the specified file format for the input. Options are:" << std::endl - << " lower-distance (lower triangular distance matrix; default)" << std::endl - << " upper-distance (upper triangular distance matrix)" << std::endl - << " distance (full distance matrix)" << std::endl - << " point-cloud (point cloud in Euclidean space)" << std::endl - << " dipha (distance matrix in DIPHA file format)" << std::endl - << " ripser (distance matrix in Ripser binary file format)" << std::endl - << " --dim compute persistent homology up to dimension " << std::endl - << " --threshold compute Rips complexes up to diameter " << std::endl + std::cerr + << "Usage: " + << "ripser " + << "[options] [filename]" << std::endl + << std::endl + << "Options:" << std::endl + << std::endl + << " --help print this screen" << std::endl + << " --format use the specified file format for the input. Options are:" + << std::endl + << " lower-distance (lower triangular distance matrix; default)" + << std::endl + << " upper-distance (upper triangular distance matrix)" << std::endl + << " distance (full distance matrix)" << std::endl + << " point-cloud (point cloud in Euclidean space)" << std::endl + << " dipha (distance matrix in DIPHA file format)" << std::endl + << " ripser (distance matrix in Ripser binary file format)" + << std::endl + << " --dim compute persistent homology up to dimension " << std::endl + << " --threshold compute Rips complexes up to diameter " << std::endl #ifdef USE_COEFFICIENTS - << " --modulus

compute homology with coefficients in the prime field Z/

Z" + << " --modulus

compute homology with coefficients in the prime field Z/

Z" #endif - << std::endl; + << std::endl; exit(exit_code); } @@ -923,7 +916,61 @@ int main(int argc, char** argv) { std::cout << "distance matrix with " << dist.size() << " points" << std::endl; auto value_range = std::minmax_element(dist.distances.begin(), dist.distances.end()); - std::cout << "value range: [" << *value_range.first << "," << *value_range.second << "]" << std::endl; + std::cout << "value range: [" << *value_range.first << "," << *value_range.second << "]" + << std::endl; ripser(std::move(dist), dim_max, threshold, modulus).compute_barcodes(); } + +void ripser::compute_barcodes() { + + std::vector columns_to_reduce; + + { + union_find dset(n); + std::vector edges; + for (index_t index = binomial_coeff(n, 2); index-- > 0;) { + value_t diameter = compute_diameter(index, 1); + if (diameter <= threshold) edges.push_back(std::make_pair(diameter, index)); + } + std::sort(edges.rbegin(), edges.rend(), + greater_diameter_or_smaller_index()); + +#ifdef PRINT_PERSISTENCE_PAIRS + std::cout << "persistence intervals in dim 0:" << std::endl; +#endif + + std::vector vertices_of_edge(2); + for (auto e : edges) { + vertices_of_edge.clear(); + get_simplex_vertices(get_index(e), 1, n, std::back_inserter(vertices_of_edge)); + index_t u = dset.find(vertices_of_edge[0]), v = dset.find(vertices_of_edge[1]); + + if (u != v) { +#ifdef PRINT_PERSISTENCE_PAIRS + if (get_diameter(e) != 0) + std::cout << " [0," << get_diameter(e) << ")" << std::endl; +#endif + dset.link(u, v); + } else + columns_to_reduce.push_back(e); + } + std::reverse(columns_to_reduce.begin(), columns_to_reduce.end()); + +#ifdef PRINT_PERSISTENCE_PAIRS + for (index_t i = 0; i < n; ++i) + if (dset.find(i) == i) std::cout << " [0, )" << std::endl << std::flush; +#endif + } + + for (index_t dim = 1; dim <= dim_max; ++dim) { + hash_map pivot_column_index; + pivot_column_index.reserve(columns_to_reduce.size()); + + compute_pairs(columns_to_reduce, pivot_column_index, dim); + + if (dim < dim_max) { + assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, dim + 1); + } + } +} -- cgit v1.2.3 From 2bed1bcf7f3d654b0f8332e003bb3ec602b84681 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Thu, 29 Dec 2016 17:59:03 +0100 Subject: code reformatting --- ripser.cpp | 247 ++++++++++++++++++++++++++++++------------------------------- 1 file changed, 121 insertions(+), 126 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 02bfac9..e17bc6a 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -464,227 +464,222 @@ public: } }; - void assemble_columns_to_reduce(std::vector& columns_to_reduce, - hash_map& pivot_column_index, index_t dim); - - void compute_pairs(std::vector& columns_to_reduce, - hash_map& pivot_column_index, index_t dim); - void compute_barcodes(); -}; -void ripser::assemble_columns_to_reduce(std::vector& columns_to_reduce, - hash_map& pivot_column_index, - index_t dim) { - index_t num_simplices = binomial_coeff(n, dim + 1); + void assemble_columns_to_reduce(std::vector& columns_to_reduce, + hash_map& pivot_column_index, index_t dim) { + index_t num_simplices = binomial_coeff(n, dim + 1); - columns_to_reduce.clear(); + columns_to_reduce.clear(); #ifdef INDICATE_PROGRESS - std::cout << "\033[K" - << "assembling " << num_simplices << " columns" << std::flush << "\r"; + std::cout << "\033[K" + << "assembling " << num_simplices << " columns" << std::flush << "\r"; #endif - for (index_t index = 0; index < num_simplices; ++index) { - if (pivot_column_index.find(index) == pivot_column_index.end()) { - value_t diameter = compute_diameter(index, dim); - if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); + for (index_t index = 0; index < num_simplices; ++index) { + if (pivot_column_index.find(index) == pivot_column_index.end()) { + value_t diameter = compute_diameter(index, dim); + if (diameter <= threshold) + columns_to_reduce.push_back(std::make_pair(diameter, index)); #ifdef INDICATE_PROGRESS - if ((index + 1) % 1000 == 0) - std::cout << "\033[K" - << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) - << "/" << num_simplices << " columns" << std::flush << "\r"; + if ((index + 1) % 1000 == 0) + std::cout << "\033[K" + << "assembled " << columns_to_reduce.size() << " out of " + << (index + 1) << "/" << num_simplices << " columns" << std::flush + << "\r"; #endif + } } - } #ifdef INDICATE_PROGRESS - std::cout << "\033[K" - << "sorting " << num_simplices << " columns" << std::flush << "\r"; + std::cout << "\033[K" + << "sorting " << num_simplices << " columns" << std::flush << "\r"; #endif - std::sort(columns_to_reduce.begin(), columns_to_reduce.end(), - greater_diameter_or_smaller_index()); + std::sort(columns_to_reduce.begin(), columns_to_reduce.end(), + greater_diameter_or_smaller_index()); #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif -} + } -void ripser::compute_pairs(std::vector& columns_to_reduce, - hash_map& pivot_column_index, index_t dim) { + void compute_pairs(std::vector& columns_to_reduce, + hash_map& pivot_column_index, index_t dim) { #ifdef PRINT_PERSISTENCE_PAIRS - std::cout << "persistence intervals in dim " << dim << ":" << std::endl; + std::cout << "persistence intervals in dim " << dim << ":" << std::endl; #endif #ifdef ASSEMBLE_REDUCTION_MATRIX - compressed_sparse_matrix reduction_coefficients; + compressed_sparse_matrix reduction_coefficients; #else #ifdef USE_COEFFICIENTS - std::vector reduction_coefficients; + std::vector reduction_coefficients; #endif #endif - std::vector coface_entries; + std::vector coface_entries; - for (index_t i = 0; i < columns_to_reduce.size(); ++i) { - auto column_to_reduce = columns_to_reduce[i]; + for (index_t i = 0; i < columns_to_reduce.size(); ++i) { + auto column_to_reduce = columns_to_reduce[i]; #ifdef ASSEMBLE_REDUCTION_MATRIX - std::priority_queue, - greater_diameter_or_smaller_index> - reduction_column; + std::priority_queue, + greater_diameter_or_smaller_index> + reduction_column; #endif - std::priority_queue, - greater_diameter_or_smaller_index> - working_coboundary; + std::priority_queue, + greater_diameter_or_smaller_index> + working_coboundary; - value_t diameter = get_diameter(column_to_reduce); + value_t diameter = get_diameter(column_to_reduce); #ifdef INDICATE_PROGRESS - if ((i + 1) % 1000 == 0) - std::cout << "\033[K" - << "reducing column " << i + 1 << "/" << columns_to_reduce.size() - << " (diameter " << diameter << ")" << std::flush << "\r"; + if ((i + 1) % 1000 == 0) + std::cout << "\033[K" + << "reducing column " << i + 1 << "/" << columns_to_reduce.size() + << " (diameter " << diameter << ")" << std::flush << "\r"; #endif - index_t j = i; + index_t j = i; - // start with a dummy pivot entry with coefficient -1 in order to initialize - // working_coboundary with the coboundary of the simplex with index column_to_reduce - diameter_entry_t pivot(0, -1, -1 + modulus); + // start with a dummy pivot entry with coefficient -1 in order to initialize + // working_coboundary with the coboundary of the simplex with index column_to_reduce + diameter_entry_t pivot(0, -1, -1 + modulus); #ifdef ASSEMBLE_REDUCTION_MATRIX - // initialize reduction_coefficients as identity matrix - reduction_coefficients.append_column(); + // initialize reduction_coefficients as identity matrix + reduction_coefficients.append_column(); #endif #ifdef USE_COEFFICIENTS - reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, 1)); + reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, 1)); #endif - bool might_be_apparent_pair = (i == j); + bool might_be_apparent_pair = (i == j); - do { - const coefficient_t factor = modulus - get_coefficient(pivot); + do { + const coefficient_t factor = modulus - get_coefficient(pivot); #ifdef ASSEMBLE_REDUCTION_MATRIX #ifdef USE_COEFFICIENTS - auto coeffs_begin = reduction_coefficients.cbegin(j), - coeffs_end = reduction_coefficients.cend(j); + auto coeffs_begin = reduction_coefficients.cbegin(j), + coeffs_end = reduction_coefficients.cend(j); #else - std::vector coeffs; - coeffs.push_back(columns_to_reduce[j]); - for (auto it = reduction_coefficients.cbegin(j); it != reduction_coefficients.cend(j); - ++it) - coeffs.push_back(*it); - auto coeffs_begin = coeffs.begin(), coeffs_end = coeffs.end(); + std::vector coeffs; + coeffs.push_back(columns_to_reduce[j]); + for (auto it = reduction_coefficients.cbegin(j); + it != reduction_coefficients.cend(j); ++it) + coeffs.push_back(*it); + auto coeffs_begin = coeffs.begin(), coeffs_end = coeffs.end(); #endif #else #ifdef USE_COEFFICIENTS - auto coeffs_begin = &reduction_coefficients[j], - coeffs_end = &reduction_coefficients[j] + 1; + auto coeffs_begin = &reduction_coefficients[j], + coeffs_end = &reduction_coefficients[j] + 1; #else - auto coeffs_begin = &columns_to_reduce[j], coeffs_end = &columns_to_reduce[j] + 1; + auto coeffs_begin = &columns_to_reduce[j], coeffs_end = &columns_to_reduce[j] + 1; #endif #endif - for (auto it = coeffs_begin; it != coeffs_end; ++it) { - diameter_entry_t simplex = *it; - set_coefficient(simplex, get_coefficient(simplex) * factor % modulus); + for (auto it = coeffs_begin; it != coeffs_end; ++it) { + diameter_entry_t simplex = *it; + set_coefficient(simplex, get_coefficient(simplex) * factor % modulus); #ifdef ASSEMBLE_REDUCTION_MATRIX - reduction_column.push(simplex); -#endif - - coface_entries.clear(); - simplex_coboundary_enumerator cofaces(simplex, dim, *this); - while (cofaces.has_next()) { - diameter_entry_t coface = cofaces.next(); - if (get_diameter(coface) <= threshold) { - coface_entries.push_back(coface); - if (might_be_apparent_pair && - (get_diameter(simplex) == get_diameter(coface))) { - if (pivot_column_index.find(get_index(coface)) == - pivot_column_index.end()) { - pivot = coface; - goto found_persistence_pair; + reduction_column.push(simplex); +#endif + + coface_entries.clear(); + simplex_coboundary_enumerator cofaces(simplex, dim, *this); + while (cofaces.has_next()) { + diameter_entry_t coface = cofaces.next(); + if (get_diameter(coface) <= threshold) { + coface_entries.push_back(coface); + if (might_be_apparent_pair && + (get_diameter(simplex) == get_diameter(coface))) { + if (pivot_column_index.find(get_index(coface)) == + pivot_column_index.end()) { + pivot = coface; + goto found_persistence_pair; + } + might_be_apparent_pair = false; } - might_be_apparent_pair = false; } } + for (auto coface : coface_entries) working_coboundary.push(coface); } - for (auto coface : coface_entries) working_coboundary.push(coface); - } - pivot = get_pivot(working_coboundary, modulus); + pivot = get_pivot(working_coboundary, modulus); - if (get_index(pivot) != -1) { - auto pair = pivot_column_index.find(get_index(pivot)); + if (get_index(pivot) != -1) { + auto pair = pivot_column_index.find(get_index(pivot)); - if (pair != pivot_column_index.end()) { - j = pair->second; - continue; - } - } else { + if (pair != pivot_column_index.end()) { + j = pair->second; + continue; + } + } else { #ifdef PRINT_PERSISTENCE_PAIRS #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif - std::cout << " [" << diameter << ", )" << std::endl << std::flush; + std::cout << " [" << diameter << ", )" << std::endl << std::flush; #endif - break; - } + break; + } - found_persistence_pair: + found_persistence_pair: #ifdef PRINT_PERSISTENCE_PAIRS - value_t death = get_diameter(pivot); - if (diameter != death) { + value_t death = get_diameter(pivot); + if (diameter != death) { #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif - std::cout << " [" << diameter << "," << death << ")" << std::endl << std::flush; - } + std::cout << " [" << diameter << "," << death << ")" << std::endl << std::flush; + } #endif - pivot_column_index.insert(std::make_pair(get_index(pivot), i)); + pivot_column_index.insert(std::make_pair(get_index(pivot), i)); #ifdef USE_COEFFICIENTS - const coefficient_t inverse = multiplicative_inverse[get_coefficient(pivot)]; + const coefficient_t inverse = multiplicative_inverse[get_coefficient(pivot)]; #endif #ifdef ASSEMBLE_REDUCTION_MATRIX // replace current column of reduction_coefficients (with a single diagonal 1 entry) // by reduction_column (possibly with a different entry on the diagonal) #ifdef USE_COEFFICIENTS - reduction_coefficients.pop_back(); + reduction_coefficients.pop_back(); #else - pop_pivot(reduction_column, modulus); + pop_pivot(reduction_column, modulus); #endif - while (true) { - diameter_entry_t e = pop_pivot(reduction_column, modulus); - if (get_index(e) == -1) break; + while (true) { + diameter_entry_t e = pop_pivot(reduction_column, modulus); + if (get_index(e) == -1) break; #ifdef USE_COEFFICIENTS - set_coefficient(e, inverse * get_coefficient(e) % modulus); - assert(get_coefficient(e) > 0); + set_coefficient(e, inverse * get_coefficient(e) % modulus); + assert(get_coefficient(e) > 0); #endif - reduction_coefficients.push_back(e); - } + reduction_coefficients.push_back(e); + } #else #ifdef USE_COEFFICIENTS - reduction_coefficients.pop_back(); - reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, inverse)); + reduction_coefficients.pop_back(); + reduction_coefficients.push_back(diameter_entry_t(column_to_reduce, inverse)); #endif #endif - break; - } while (true); - } + break; + } while (true); + } #ifdef INDICATE_PROGRESS - std::cout << "\033[K"; + std::cout << "\033[K"; #endif -} + } +}; enum file_format { LOWER_DISTANCE_MATRIX, -- cgit v1.2.3 From 82b1e5ec6c9ee4aedc7c3e6701c0dec9f66bd56f Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Sat, 7 Jan 2017 00:24:40 +0100 Subject: constructors fixed --- ripser.cpp | 34 ++++++++++++++++++++-------------- 1 file changed, 20 insertions(+), 14 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index e17bc6a..53fa4fa 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -103,8 +103,8 @@ static_assert(sizeof(entry_t) == sizeof(index_t), "size of entry_t is not the sa entry_t make_entry(index_t _index, coefficient_t _coefficient) { return entry_t(_index, _coefficient); } -index_t get_index(entry_t e) { return e.index; } -index_t get_coefficient(entry_t e) { return e.coefficient; } +index_t get_index(const entry_t& e) { return e.index; } +index_t get_coefficient(const entry_t& e) { return e.coefficient; } void set_coefficient(entry_t& e, const coefficient_t c) { e.coefficient = c; } bool operator==(const entry_t& e1, const entry_t& e2) { @@ -119,10 +119,10 @@ std::ostream& operator<<(std::ostream& stream, const entry_t& e) { #else typedef index_t entry_t; -const index_t get_index(entry_t i) { return i; } -index_t get_coefficient(entry_t i) { return 1; } +const index_t get_index(const entry_t& i) { return i; } +index_t get_coefficient(const entry_t& i) { return 1; } entry_t make_entry(index_t _index, coefficient_t _value) { return entry_t(_index); } -void set_coefficient(index_t& e, const coefficient_t c) {} +void set_coefficient(entry_t& e, const coefficient_t c) {} #endif @@ -131,22 +131,22 @@ const entry_t& get_entry(const entry_t& e) { return e; } class diameter_index_t : public std::pair { public: diameter_index_t() : std::pair() {} - diameter_index_t(std::pair p) : std::pair(p) {} + diameter_index_t(std::pair&& p) : std::pair(std::move(p)) {} }; -value_t get_diameter(diameter_index_t i) { return i.first; } -index_t get_index(diameter_index_t i) { return i.second; } +value_t get_diameter(const diameter_index_t& i) { return i.first; } +index_t get_index(const diameter_index_t& i) { return i.second; } class diameter_entry_t : public std::pair { public: diameter_entry_t(std::pair p) : std::pair(p) {} - diameter_entry_t(entry_t e) : std::pair(0, e) {} - diameter_entry_t() : diameter_entry_t(0) {} + diameter_entry_t(entry_t&& e) : std::pair(0, std::move(e)) {} + diameter_entry_t() : diameter_entry_t(entry_t()) {} diameter_entry_t(value_t _diameter, index_t _index, coefficient_t _coefficient) : std::pair(_diameter, make_entry(_index, _coefficient)) {} - diameter_entry_t(diameter_index_t _diameter_index, coefficient_t _coefficient) + diameter_entry_t(const diameter_index_t& _diameter_index, coefficient_t _coefficient) : std::pair(get_diameter(_diameter_index), make_entry(get_index(_diameter_index), _coefficient)) {} - diameter_entry_t(diameter_index_t _diameter_index) : diameter_entry_t(_diameter_index, 1) {} + diameter_entry_t(const diameter_index_t& _diameter_index) : diameter_entry_t(_diameter_index, 1) {} }; const entry_t& get_entry(const diameter_entry_t& p) { return p.second; } @@ -230,9 +230,15 @@ public: std::vector> points; euclidean_distance_matrix(std::vector>&& _points) - : points(std::move(_points)) {} + : points(std::move(_points)) { + for (auto p: points) { + assert(p.size() == points.front().size()); + } + } value_t operator()(const index_t i, const index_t j) const { + assert(i < points.size()); + assert(j < points.size()); return std::sqrt(std::inner_product( points[i].begin(), points[i].end(), points[j].begin(), value_t(), std::plus(), [](value_t u, value_t v) { return (u - v) * (u - v); })); @@ -537,7 +543,7 @@ public: value_t diameter = get_diameter(column_to_reduce); #ifdef INDICATE_PROGRESS - if ((i + 1) % 1000 == 0) + //if ((i + 1) % 1000 == 0) std::cout << "\033[K" << "reducing column " << i + 1 << "/" << columns_to_reduce.size() << " (diameter " << diameter << ")" << std::flush << "\r"; -- cgit v1.2.3 From 0dc8dc9c9dd444951928dcf6c59967f9f617f488 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Mon, 9 Jan 2017 11:56:45 +0100 Subject: benchmark H^1 --- benchmarks/sphere_3_192/run_benchmarks.sh | 4 ++++ benchmarks/sphere_3_192/run_dionysus.sh | 1 + benchmarks/sphere_3_192/run_dipha.sh | 1 + benchmarks/sphere_3_192/run_gudhi.sh | 1 + benchmarks/sphere_3_192/run_ripser.sh | 1 + 5 files changed, 8 insertions(+) create mode 100644 benchmarks/sphere_3_192/run_benchmarks.sh create mode 100644 benchmarks/sphere_3_192/run_dionysus.sh create mode 100644 benchmarks/sphere_3_192/run_dipha.sh create mode 100644 benchmarks/sphere_3_192/run_gudhi.sh create mode 100644 benchmarks/sphere_3_192/run_ripser.sh diff --git a/benchmarks/sphere_3_192/run_benchmarks.sh b/benchmarks/sphere_3_192/run_benchmarks.sh new file mode 100644 index 0000000..14cac6c --- /dev/null +++ b/benchmarks/sphere_3_192/run_benchmarks.sh @@ -0,0 +1,4 @@ +. run_dionysus.sh +. run_dipha.sh +. run_gudhi.sh +. run_ripser.sh diff --git a/benchmarks/sphere_3_192/run_dionysus.sh b/benchmarks/sphere_3_192/run_dionysus.sh new file mode 100644 index 0000000..1e5c7d3 --- /dev/null +++ b/benchmarks/sphere_3_192/run_dionysus.sh @@ -0,0 +1 @@ +/usr/bin/time -l ~/Source/Dionysus/examples/cohomology/rips-pairwise-cohomology -s2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee dionysus.out.txt diff --git a/benchmarks/sphere_3_192/run_dipha.sh b/benchmarks/sphere_3_192/run_dipha.sh new file mode 100644 index 0000000..312bcd7 --- /dev/null +++ b/benchmarks/sphere_3_192/run_dipha.sh @@ -0,0 +1 @@ +/usr/bin/time -l ~/Bitbucket/dipha/dipha --benchmark --upper_dim 2 --dual ~/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex /dev/null 2>&1 | tee dipha.out.txt diff --git a/benchmarks/sphere_3_192/run_gudhi.sh b/benchmarks/sphere_3_192/run_gudhi.sh new file mode 100644 index 0000000..a4c6c51 --- /dev/null +++ b/benchmarks/sphere_3_192/run_gudhi.sh @@ -0,0 +1 @@ +/usr/bin/time -l ~/Source/Gudhi_library_1.3.1/example/Persistent_cohomology/rips_persistence -d2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat -o/dev/null 2>&1 | tee gudhi.out.txt diff --git a/benchmarks/sphere_3_192/run_ripser.sh b/benchmarks/sphere_3_192/run_ripser.sh new file mode 100644 index 0000000..209c410 --- /dev/null +++ b/benchmarks/sphere_3_192/run_ripser.sh @@ -0,0 +1 @@ +/usr/bin/time -l ~/Bitbucket/ripser/ripser ~/Bitbucket/ripser/examples/sphere_3_192.lower_distance_matrix 2>&1 | tee ripser.out.txt -- cgit v1.2.3 From 3c3f3e74e26a1534202cf0235956c8a7924085ed Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Wed, 11 Jan 2017 13:16:58 +0100 Subject: gudhi benchmark output --- benchmarks/sphere_3_192/run_gudhi.sh | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/sphere_3_192/run_gudhi.sh b/benchmarks/sphere_3_192/run_gudhi.sh index a4c6c51..80ae403 100644 --- a/benchmarks/sphere_3_192/run_gudhi.sh +++ b/benchmarks/sphere_3_192/run_gudhi.sh @@ -1 +1 @@ -/usr/bin/time -l ~/Source/Gudhi_library_1.3.1/example/Persistent_cohomology/rips_persistence -d2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat -o/dev/null 2>&1 | tee gudhi.out.txt +/usr/bin/time -l ~/Source/Gudhi_library_1.3.1/example/Persistent_cohomology/rips_persistence -d2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee gudhi.out.txt -- cgit v1.2.3 From f84042ed1b0b1de7f66e57d9a516cf60736af91a Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Wed, 11 Jan 2017 13:22:21 +0100 Subject: benchmark results --- benchmarks/sphere_3_192/dionysus.out.txt | 27 ++++ benchmarks/sphere_3_192/dipha.out.txt | 55 +++++++ benchmarks/sphere_3_192/gudhi.out.txt | 262 ++++++++++++++++++++++++++++++ benchmarks/sphere_3_192/ripser.out.txt | 264 +++++++++++++++++++++++++++++++ 4 files changed, 608 insertions(+) create mode 100644 benchmarks/sphere_3_192/dionysus.out.txt create mode 100644 benchmarks/sphere_3_192/dipha.out.txt create mode 100644 benchmarks/sphere_3_192/gudhi.out.txt create mode 100644 benchmarks/sphere_3_192/ripser.out.txt diff --git a/benchmarks/sphere_3_192/dionysus.out.txt b/benchmarks/sphere_3_192/dionysus.out.txt new file mode 100644 index 0000000..3191f98 --- /dev/null +++ b/benchmarks/sphere_3_192/dionysus.out.txt @@ -0,0 +1,27 @@ +Boundary matrix: +Cocycles: *.ccl +Vertices: +Diagram: +Simplex vector generated, size: 1179808 + +0% 10 20 30 40 50 60 70 80 90 100% +|----|----|----|----|----|----|----|----|----|----| +*************************************************** +Rips timer : Elapsed time [4.39] seconds +Persistence timer : Elapsed time [2.55] seconds +Total timer : Elapsed time [7.88] seconds + 7.99 real 7.88 user 0.08 sys + 141664256 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 34606 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 0 voluntary context switches + 2629 involuntary context switches diff --git a/benchmarks/sphere_3_192/dipha.out.txt b/benchmarks/sphere_3_192/dipha.out.txt new file mode 100644 index 0000000..bae5fc7 --- /dev/null +++ b/benchmarks/sphere_3_192/dipha.out.txt @@ -0,0 +1,55 @@ + +Input filename: +/Users/uli/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex + +upper_dim: 2 + +Number of processes used: +1 + +Detailed information for rank 0: + time prior mem peak mem bytes recv + 0.0s 3 MB 4 MB 0 MB complex.load_binary( input_filename, upper_dim ); + +Number of cells in input: +1179808 + 0.2s 4 MB 40 MB 0 MB get_filtration_to_cell_map( complex, dualize, filtration_to_cell_map ); + 0.1s 40 MB 127 MB 18 MB get_cell_to_filtration_map( complex.get_num_cells(), filtration_to_cell_map, cell_to_filtration_map ); + 0.0s 154 MB 156 MB 0 MB generate_unreduced_columns( complex, filtration_to_cell_map, cell_to_filtration_map, cur_dim, dualize, unreduced_columns ); + 0.0s 153 MB 156 MB 0 MB reduction_kernel( complex.get_num_cells(), unreduced_columns, reduced_columns ); + 0.3s 151 MB 165 MB 53 MB generate_unreduced_columns( complex, filtration_to_cell_map, cell_to_filtration_map, cur_dim, dualize, unreduced_columns ); + 0.1s 136 MB 169 MB 0 MB reduction_kernel( complex.get_num_cells(), unreduced_columns, reduced_columns ); + 0.1s 167 MB 169 MB 1 MB dipha::outputs::save_persistence_diagram( output_filename, complex, filtration_to_cell_map, reduced_columns, dualize, upper_dim ); + +Overall running time in seconds: +0.9 + +Reduction kernel running time in seconds: +0.1 + +Overall peak mem in GB of all ranks: +0.2 + +Individual peak mem in GB of per rank: +0.2 + +Maximal communication traffic (without sorting) in GB between any pair of nodes: +0.1 + +Total communication traffic (without sorting) in GB between all pairs of nodes: +0.1 + 0.94 real 0.75 user 0.16 sys + 177598464 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 78159 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 5 voluntary context switches + 333 involuntary context switches diff --git a/benchmarks/sphere_3_192/gudhi.out.txt b/benchmarks/sphere_3_192/gudhi.out.txt new file mode 100644 index 0000000..1a0cddf --- /dev/null +++ b/benchmarks/sphere_3_192/gudhi.out.txt @@ -0,0 +1,262 @@ +The complex contains 1179808 simplices + and has dimension 2 +2 0 0 inf +2 1 0.316534 0.682559 +2 1 0.304761 0.647803 +2 0 0 0.332695 +2 0 0 0.323964 +2 0 0 0.319493 +2 0 0 0.316023 +2 1 0.34298 0.65758 +2 0 0 0.304062 +2 0 0 0.297933 +2 0 0 0.297053 +2 0 0 0.293466 +2 1 0.347806 0.638039 +2 0 0 0.289189 +2 0 0 0.288591 +2 0 0 0.287316 +2 0 0 0.285416 +2 0 0 0.283881 +2 0 0 0.283562 +2 0 0 0.282169 +2 0 0 0.280629 +2 0 0 0.277582 +2 0 0 0.277191 +2 0 0 0.276431 +2 0 0 0.274611 +2 0 0 0.270996 +2 0 0 0.268806 +2 0 0 0.266635 +2 0 0 0.266592 +2 0 0 0.264483 +2 0 0 0.264334 +2 0 0 0.257997 +2 0 0 0.25725 +2 0 0 0.253227 +2 0 0 0.253163 +2 0 0 0.249345 +2 0 0 0.249339 +2 0 0 0.248879 +2 0 0 0.247042 +2 0 0 0.246491 +2 0 0 0.24509 +2 0 0 0.24348 +2 0 0 0.240632 +2 0 0 0.23879 +2 0 0 0.237571 +2 0 0 0.234969 +2 0 0 0.234112 +2 0 0 0.231675 +2 0 0 0.231179 +2 1 0.35058 0.580726 +2 0 0 0.228917 +2 0 0 0.2286 +2 0 0 0.227933 +2 0 0 0.227524 +2 0 0 0.223232 +2 0 0 0.218561 +2 0 0 0.218526 +2 0 0 0.215759 +2 0 0 0.214018 +2 1 0.347388 0.559978 +2 0 0 0.210372 +2 0 0 0.208851 +2 0 0 0.208063 +2 0 0 0.207714 +2 0 0 0.207261 +2 0 0 0.207016 +2 0 0 0.205389 +2 0 0 0.204544 +2 0 0 0.204208 +2 0 0 0.203812 +2 0 0 0.202564 +2 0 0 0.201961 +2 0 0 0.199673 +2 0 0 0.197377 +2 0 0 0.195722 +2 0 0 0.193804 +2 0 0 0.193791 +2 0 0 0.192733 +2 0 0 0.192406 +2 0 0 0.192356 +2 0 0 0.191661 +2 0 0 0.190642 +2 1 0.352356 0.541947 +2 0 0 0.189106 +2 0 0 0.184757 +2 0 0 0.18428 +2 0 0 0.179609 +2 1 0.301986 0.478334 +2 0 0 0.176306 +2 0 0 0.175016 +2 0 0 0.17493 +2 0 0 0.174172 +2 0 0 0.173684 +2 0 0 0.173586 +2 0 0 0.173412 +2 0 0 0.168587 +2 0 0 0.167455 +2 0 0 0.165907 +2 0 0 0.165592 +2 0 0 0.165494 +2 0 0 0.165314 +2 0 0 0.165202 +2 0 0 0.164815 +2 0 0 0.164809 +2 0 0 0.162464 +2 1 0.322683 0.482958 +2 0 0 0.159282 +2 0 0 0.158413 +2 0 0 0.157373 +2 0 0 0.157293 +2 0 0 0.157251 +2 0 0 0.1569 +2 0 0 0.15587 +2 0 0 0.154095 +2 0 0 0.15244 +2 0 0 0.151643 +2 0 0 0.150747 +2 0 0 0.149226 +2 0 0 0.149017 +2 0 0 0.148506 +2 0 0 0.147427 +2 0 0 0.147289 +2 0 0 0.147203 +2 0 0 0.145725 +2 0 0 0.145018 +2 0 0 0.143151 +2 0 0 0.143085 +2 0 0 0.142834 +2 0 0 0.142537 +2 0 0 0.142152 +2 0 0 0.141823 +2 0 0 0.138674 +2 1 0.377729 0.514046 +2 0 0 0.135402 +2 0 0 0.134886 +2 0 0 0.133722 +2 0 0 0.133525 +2 0 0 0.133033 +2 0 0 0.131638 +2 0 0 0.131438 +2 0 0 0.129741 +2 0 0 0.127864 +2 0 0 0.127383 +2 0 0 0.126595 +2 0 0 0.124774 +2 1 0.339394 0.463666 +2 0 0 0.123527 +2 0 0 0.122089 +2 0 0 0.121955 +2 0 0 0.119187 +2 1 0.2576 0.373749 +2 0 0 0.115939 +2 1 0.370051 0.483155 +2 0 0 0.111798 +2 0 0 0.111131 +2 0 0 0.109526 +2 0 0 0.108768 +2 1 0.386278 0.4922 +2 1 0.443911 0.54761 +2 1 0.374531 0.477153 +2 1 0.354747 0.454956 +2 0 0 0.0994522 +2 1 0.318032 0.417462 +2 0 0 0.0988752 +2 0 0 0.0970981 +2 0 0 0.0959033 +2 0 0 0.0884133 +2 0 0 0.0870489 +2 0 0 0.0857401 +2 0 0 0.084341 +2 0 0 0.0831521 +2 0 0 0.0831358 +2 0 0 0.0814685 +2 0 0 0.0802451 +2 0 0 0.0800348 +2 0 0 0.0792212 +2 0 0 0.0782455 +2 0 0 0.076909 +2 0 0 0.0753955 +2 1 0.324159 0.399094 +2 0 0 0.0746054 +2 0 0 0.0743236 +2 0 0 0.0741239 +2 1 0.413789 0.487379 +2 0 0 0.0728647 +2 0 0 0.0714516 +2 0 0 0.0702075 +2 0 0 0.0692368 +2 1 0.409968 0.478461 +2 0 0 0.0672609 +2 0 0 0.0651112 +2 0 0 0.0637053 +2 0 0 0.0633425 +2 1 0.298341 0.360737 +2 0 0 0.0623736 +2 1 0.309419 0.37151 +2 1 0.411549 0.471715 +2 0 0 0.0593071 +2 0 0 0.057488 +2 1 0.377019 0.434385 +2 1 0.329066 0.386148 +2 0 0 0.0560443 +2 0 0 0.0559886 +2 0 0 0.0522132 +2 0 0 0.0515157 +2 1 0.412572 0.46308 +2 0 0 0.0496187 +2 0 0 0.048826 +2 0 0 0.0487061 +2 1 0.240869 0.288747 +2 0 0 0.0471186 +2 1 0.531636 0.578093 +2 1 0.530723 0.576869 +2 1 0.33364 0.3797 +2 1 0.431628 0.477277 +2 0 0 0.0433426 +2 1 0.344005 0.387193 +2 1 0.463389 0.505345 +2 1 0.361715 0.403181 +2 1 0.381084 0.421374 +2 1 0.253337 0.292286 +2 0 0 0.0387451 +2 0 0 0.037954 +2 1 0.377147 0.414788 +2 1 0.308779 0.346297 +2 1 0.317521 0.349622 +2 0 0 0.0311057 +2 1 0.331259 0.356381 +2 1 0.279016 0.301331 +2 0 0 0.0217315 +2 1 0.350913 0.369543 +2 1 0.324571 0.343083 +2 0 0 0.0173134 +2 1 0.33392 0.35111 +2 1 0.542696 0.558863 +2 1 0.33836 0.350411 +2 1 0.29997 0.311993 +2 0 0 0.0120098 +2 1 0.189652 0.197515 +2 0 0 0.00391946 +2 0 0 0.0036753 +2 1 0.445398 0.448892 +2 1 0.24546 0.248903 +2 1 0.301045 0.303943 +2 1 0.309336 0.310216 + 0.59 real 0.56 user 0.02 sys + 69210112 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 16916 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 10 voluntary context switches + 125 involuntary context switches diff --git a/benchmarks/sphere_3_192/ripser.out.txt b/benchmarks/sphere_3_192/ripser.out.txt new file mode 100644 index 0000000..bc741ae --- /dev/null +++ b/benchmarks/sphere_3_192/ripser.out.txt @@ -0,0 +1,264 @@ +distance matrix with 192 points +value range: [0.00367531,2] +persistence intervals in dim 0: + [0,0.00367531) + [0,0.00391946) + [0,0.0120098) + [0,0.0173134) + [0,0.0217315) + [0,0.0311057) + [0,0.037954) + [0,0.0387451) + [0,0.0433426) + [0,0.0471186) + [0,0.0487061) + [0,0.048826) + [0,0.0496187) + [0,0.0515157) + [0,0.0522132) + [0,0.0559886) + [0,0.0560443) + [0,0.057488) + [0,0.0593071) + [0,0.0623736) + [0,0.0633425) + [0,0.0637053) + [0,0.0651112) + [0,0.0672609) + [0,0.0692368) + [0,0.0702075) + [0,0.0714516) + [0,0.0728647) + [0,0.0741239) + [0,0.0743236) + [0,0.0746054) + [0,0.0753955) + [0,0.076909) + [0,0.0782455) + [0,0.0792212) + [0,0.0800348) + [0,0.0802451) + [0,0.0814685) + [0,0.0831358) + [0,0.0831521) + [0,0.084341) + [0,0.0857401) + [0,0.0870489) + [0,0.0884133) + [0,0.0959033) + [0,0.0970981) + [0,0.0988752) + [0,0.0994522) + [0,0.108768) + [0,0.109526) + [0,0.111131) + [0,0.111798) + [0,0.115939) + [0,0.119187) + [0,0.121955) + [0,0.122089) + [0,0.123527) + [0,0.124774) + [0,0.126595) + [0,0.127383) + [0,0.127864) + [0,0.129741) + [0,0.131438) + [0,0.131638) + [0,0.133033) + [0,0.133525) + [0,0.133722) + [0,0.134886) + [0,0.135402) + [0,0.138674) + [0,0.141823) + [0,0.142152) + [0,0.142537) + [0,0.142834) + [0,0.143085) + [0,0.143151) + [0,0.145018) + [0,0.145725) + [0,0.147203) + [0,0.147289) + [0,0.147427) + [0,0.148506) + [0,0.149017) + [0,0.149226) + [0,0.150747) + [0,0.151643) + [0,0.15244) + [0,0.154095) + [0,0.15587) + [0,0.1569) + [0,0.157251) + [0,0.157293) + [0,0.157373) + [0,0.158413) + [0,0.159282) + [0,0.162464) + [0,0.164809) + [0,0.164815) + [0,0.165202) + [0,0.165314) + [0,0.165494) + [0,0.165592) + [0,0.165907) + [0,0.167455) + [0,0.168587) + [0,0.173412) + [0,0.173586) + [0,0.173684) + [0,0.174172) + [0,0.17493) + [0,0.175016) + [0,0.176306) + [0,0.179609) + [0,0.18428) + [0,0.184757) + [0,0.189106) + [0,0.190642) + [0,0.191661) + [0,0.192356) + [0,0.192406) + [0,0.192733) + [0,0.193791) + [0,0.193804) + [0,0.195722) + [0,0.197377) + [0,0.199673) + [0,0.201961) + [0,0.202564) + [0,0.203812) + [0,0.204208) + [0,0.204544) + [0,0.205389) + [0,0.207016) + [0,0.207261) + [0,0.207714) + [0,0.208063) + [0,0.208851) + [0,0.210372) + [0,0.214018) + [0,0.215759) + [0,0.218526) + [0,0.218561) + [0,0.223232) + [0,0.227524) + [0,0.227933) + [0,0.2286) + [0,0.228917) + [0,0.231179) + [0,0.231675) + [0,0.234112) + [0,0.234969) + [0,0.237571) + [0,0.23879) + [0,0.240632) + [0,0.24348) + [0,0.24509) + [0,0.246491) + [0,0.247042) + [0,0.248879) + [0,0.249339) + [0,0.249345) + [0,0.253163) + [0,0.253227) + [0,0.25725) + [0,0.257997) + [0,0.264334) + [0,0.264483) + [0,0.266592) + [0,0.266635) + [0,0.268806) + [0,0.270996) + [0,0.274611) + [0,0.276431) + [0,0.277191) + [0,0.277582) + [0,0.280629) + [0,0.282169) + [0,0.283562) + [0,0.283881) + [0,0.285416) + [0,0.287316) + [0,0.288591) + [0,0.289189) + [0,0.293466) + [0,0.297053) + [0,0.297933) + [0,0.304062) + [0,0.316023) + [0,0.319493) + [0,0.323964) + [0,0.332695) + [0, ) +persistence intervals in dim 1: + [0.542696,0.558863) + [0.531636,0.578093) + [0.530723,0.576869) + [0.463389,0.505345) + [0.445398,0.448892) + [0.443911,0.54761) + [0.431628,0.477277) + [0.413789,0.487379) + [0.412572,0.46308) + [0.411549,0.471715) + [0.409968,0.478461) + [0.386278,0.4922) + [0.381084,0.421374) + [0.377729,0.514046) + [0.377147,0.414788) + [0.377019,0.434385) + [0.374531,0.477153) + [0.370051,0.483155) + [0.361715,0.403181) + [0.354747,0.454956) + [0.352356,0.541947) + [0.350913,0.369543) + [0.35058,0.580726) + [0.347806,0.638039) + [0.347388,0.559978) + [0.344005,0.387193) + [0.34298,0.65758) + [0.339394,0.463666) + [0.33836,0.350411) + [0.33392,0.35111) + [0.33364,0.3797) + [0.331259,0.356381) + [0.329066,0.386148) + [0.324571,0.343083) + [0.324159,0.399094) + [0.322683,0.482958) + [0.318032,0.417462) + [0.317521,0.349622) + [0.316534,0.682559) + [0.309419,0.37151) + [0.309336,0.310216) + [0.308779,0.346297) + [0.304761,0.647803) + [0.301986,0.478334) + [0.301045,0.303943) + [0.29997,0.311993) + [0.298341,0.360737) + [0.279016,0.301331) + [0.2576,0.373749) + [0.253337,0.292286) + [0.24546,0.248903) + [0.240869,0.288747) + [0.189652,0.197515) + 0.03 real 0.03 user 0.00 sys + 3276800 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 819 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 7 voluntary context switches + 119 involuntary context switches -- cgit v1.2.3 From eb8a0a7d92c1b3260ff33e8f213b7369dd038d0b Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Wed, 22 Feb 2017 12:55:07 +0100 Subject: typo in readme --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index 0c5a2e6..4aae03d 100644 --- a/README.md +++ b/README.md @@ -74,7 +74,7 @@ The input is given either in a file whose name is passed as an argument, or thro - `--format`: use the specified file format for the input. The following formats are supported: - `lower-distance` (default if no format is specified): lower triangular distance matrix; a comma (or whitespace, or other non-numerical character) separated list of the distance matrix entries below the diagonal, sorted lexicographically by row index, then column index - `upper-distance`: upper triangular distance matrix; similar to the previous, but for the entries above the diagonal; suitable for output from the MATLAB functions `pdist` or `seqpdist`, exported to a CSV file - - `distances`: full distance matrix; similar to the above, but for all entries of the distance matrix + - `distance`: full distance matrix; similar to the above, but for all entries of the distance matrix - `dipha`: DIPHA distance matrix as described on the [DIPHA] website - `point-cloud`: point cloud; a comma (or whitespace, or other non-numerical character) separated list of coordinates of the points in some Euclidean space, one point per line - `--dim k`: compute persistent homology up to dimension *k* -- cgit v1.2.3 From a5ada276181dadded36a8cb19d83466d1f098278 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Mon, 20 Mar 2017 15:45:21 +0100 Subject: updated benchmark --- benchmarks/sphere_3_192/dionysus.dim_1.out.txt | 27 +++ benchmarks/sphere_3_192/dionysus.dim_2.out.txt | 27 +++ benchmarks/sphere_3_192/dipha.dim_1.out.txt | 55 +++++ benchmarks/sphere_3_192/dipha.dim_2.out.txt | 57 ++++++ benchmarks/sphere_3_192/gudhi.dim_1.out.txt | 262 ++++++++++++++++++++++++ benchmarks/sphere_3_192/gudhi.dim_2.out.txt | 263 ++++++++++++++++++++++++ benchmarks/sphere_3_192/ripser.dim_1.out.txt | 264 ++++++++++++++++++++++++ benchmarks/sphere_3_192/ripser.dim_2.out.txt | 266 +++++++++++++++++++++++++ benchmarks/sphere_3_192/run_dionysus.sh | 4 +- benchmarks/sphere_3_192/run_dipha.sh | 4 +- benchmarks/sphere_3_192/run_gudhi.sh | 4 +- benchmarks/sphere_3_192/run_ripser.sh | 4 +- 12 files changed, 1233 insertions(+), 4 deletions(-) create mode 100644 benchmarks/sphere_3_192/dionysus.dim_1.out.txt create mode 100644 benchmarks/sphere_3_192/dionysus.dim_2.out.txt create mode 100644 benchmarks/sphere_3_192/dipha.dim_1.out.txt create mode 100644 benchmarks/sphere_3_192/dipha.dim_2.out.txt create mode 100644 benchmarks/sphere_3_192/gudhi.dim_1.out.txt create mode 100644 benchmarks/sphere_3_192/gudhi.dim_2.out.txt create mode 100644 benchmarks/sphere_3_192/ripser.dim_1.out.txt create mode 100644 benchmarks/sphere_3_192/ripser.dim_2.out.txt diff --git a/benchmarks/sphere_3_192/dionysus.dim_1.out.txt b/benchmarks/sphere_3_192/dionysus.dim_1.out.txt new file mode 100644 index 0000000..0e8cd4d --- /dev/null +++ b/benchmarks/sphere_3_192/dionysus.dim_1.out.txt @@ -0,0 +1,27 @@ +Boundary matrix: +Cocycles: *.ccl +Vertices: +Diagram: +Simplex vector generated, size: 1179808 + +0% 10 20 30 40 50 60 70 80 90 100% +|----|----|----|----|----|----|----|----|----|----| +*************************************************** +Rips timer : Elapsed time [4.19] seconds +Persistence timer : Elapsed time [2.38] seconds +Total timer : Elapsed time [7.45] seconds + 7.54 real 7.46 user 0.06 sys + 142344192 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 34700 page reclaims + 72 page faults + 0 swaps + 5 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 7 voluntary context switches + 212 involuntary context switches diff --git a/benchmarks/sphere_3_192/dionysus.dim_2.out.txt b/benchmarks/sphere_3_192/dionysus.dim_2.out.txt new file mode 100644 index 0000000..3495e23 --- /dev/null +++ b/benchmarks/sphere_3_192/dionysus.dim_2.out.txt @@ -0,0 +1,27 @@ +Boundary matrix: +Cocycles: *.ccl +Vertices: +Diagram: +Simplex vector generated, size: 56050288 + +0% 10 20 30 40 50 60 70 80 90 100% +|----|----|----|----|----|----|----|----|----|----| +*************************************************** +Rips timer : Elapsed time [256.35] seconds +Persistence timer : Elapsed time [205.61] seconds +Total timer : Elapsed time [528.91] seconds + 532.99 real 527.96 user 4.95 sys +3376783360 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 1430837 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 0 voluntary context switches + 17326 involuntary context switches diff --git a/benchmarks/sphere_3_192/dipha.dim_1.out.txt b/benchmarks/sphere_3_192/dipha.dim_1.out.txt new file mode 100644 index 0000000..133977f --- /dev/null +++ b/benchmarks/sphere_3_192/dipha.dim_1.out.txt @@ -0,0 +1,55 @@ + +Input filename: +/Users/uli/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex + +upper_dim: 2 + +Number of processes used: +1 + +Detailed information for rank 0: + time prior mem peak mem bytes recv + 0.0s 3 MB 4 MB 0 MB complex.load_binary( input_filename, upper_dim ); + +Number of cells in input: +1179808 + 0.2s 4 MB 40 MB 0 MB get_filtration_to_cell_map( complex, dualize, filtration_to_cell_map ); + 0.1s 40 MB 127 MB 18 MB get_cell_to_filtration_map( complex.get_num_cells(), filtration_to_cell_map, cell_to_filtration_map ); + 0.0s 154 MB 156 MB 0 MB generate_unreduced_columns( complex, filtration_to_cell_map, cell_to_filtration_map, cur_dim, dualize, unreduced_columns ); + 0.0s 153 MB 156 MB 0 MB reduction_kernel( complex.get_num_cells(), unreduced_columns, reduced_columns ); + 0.3s 150 MB 165 MB 53 MB generate_unreduced_columns( complex, filtration_to_cell_map, cell_to_filtration_map, cur_dim, dualize, unreduced_columns ); + 0.1s 136 MB 169 MB 0 MB reduction_kernel( complex.get_num_cells(), unreduced_columns, reduced_columns ); + 0.1s 167 MB 169 MB 1 MB dipha::outputs::save_persistence_diagram( output_filename, complex, filtration_to_cell_map, reduced_columns, dualize, upper_dim ); + +Overall running time in seconds: +0.8 + +Reduction kernel running time in seconds: +0.1 + +Overall peak mem in GB of all ranks: +0.2 + +Individual peak mem in GB of per rank: +0.2 + +Maximal communication traffic (without sorting) in GB between any pair of nodes: +0.1 + +Total communication traffic (without sorting) in GB between all pairs of nodes: +0.1 + 0.90 real 0.71 user 0.15 sys + 177610752 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 78019 page reclaims + 164 page faults + 0 swaps + 10 block input operations + 8 block output operations + 0 messages sent + 0 messages received + 0 signals received + 48 voluntary context switches + 367 involuntary context switches diff --git a/benchmarks/sphere_3_192/dipha.dim_2.out.txt b/benchmarks/sphere_3_192/dipha.dim_2.out.txt new file mode 100644 index 0000000..1a7b30b --- /dev/null +++ b/benchmarks/sphere_3_192/dipha.dim_2.out.txt @@ -0,0 +1,57 @@ + +Input filename: +/Users/uli/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex + +upper_dim: 3 + +Number of processes used: +1 + +Detailed information for rank 0: + time prior mem peak mem bytes recv + 0.0s 3 MB 4 MB 0 MB complex.load_binary( input_filename, upper_dim ); + +Number of cells in input: +56050288 + 12.3s 4 MB 1714 MB 0 MB get_filtration_to_cell_map( complex, dualize, filtration_to_cell_map ); + 5.5s 431 MB 4416 MB 855 MB get_cell_to_filtration_map( complex.get_num_cells(), filtration_to_cell_map, cell_to_filtration_map ); + 0.3s 2278 MB 4416 MB 0 MB generate_unreduced_columns( complex, filtration_to_cell_map, cell_to_filtration_map, cur_dim, dualize, unreduced_columns ); + 0.0s 2279 MB 4416 MB 0 MB reduction_kernel( complex.get_num_cells(), unreduced_columns, reduced_columns ); + 0.6s 2279 MB 4416 MB 53 MB generate_unreduced_columns( complex, filtration_to_cell_map, cell_to_filtration_map, cur_dim, dualize, unreduced_columns ); + 0.0s 2240 MB 4416 MB 0 MB reduction_kernel( complex.get_num_cells(), unreduced_columns, reduced_columns ); + 18.7s 2252 MB 5285 MB 3349 MB generate_unreduced_columns( complex, filtration_to_cell_map, cell_to_filtration_map, cur_dim, dualize, unreduced_columns ); + 5.9s 3930 MB 5715 MB 0 MB reduction_kernel( complex.get_num_cells(), unreduced_columns, reduced_columns ); + 6.8s 3716 MB 5715 MB 106 MB dipha::outputs::save_persistence_diagram( output_filename, complex, filtration_to_cell_map, reduced_columns, dualize, upper_dim ); + +Overall running time in seconds: +50.9 + +Reduction kernel running time in seconds: +5.9 + +Overall peak mem in GB of all ranks: +5.6 + +Individual peak mem in GB of per rank: +5.6 + +Maximal communication traffic (without sorting) in GB between any pair of nodes: +4.3 + +Total communication traffic (without sorting) in GB between all pairs of nodes: +4.3 + 50.97 real 42.99 user 7.93 sys +5993607168 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 5448458 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 6 block output operations + 0 messages sent + 0 messages received + 0 signals received + 6 voluntary context switches + 6625 involuntary context switches diff --git a/benchmarks/sphere_3_192/gudhi.dim_1.out.txt b/benchmarks/sphere_3_192/gudhi.dim_1.out.txt new file mode 100644 index 0000000..94e7f6e --- /dev/null +++ b/benchmarks/sphere_3_192/gudhi.dim_1.out.txt @@ -0,0 +1,262 @@ +The complex contains 1179808 simplices + and has dimension 2 +2 0 0 inf +2 1 0.316534 0.682559 +2 1 0.304761 0.647803 +2 0 0 0.332695 +2 0 0 0.323964 +2 0 0 0.319493 +2 0 0 0.316023 +2 1 0.34298 0.65758 +2 0 0 0.304062 +2 0 0 0.297933 +2 0 0 0.297053 +2 0 0 0.293466 +2 1 0.347806 0.638039 +2 0 0 0.289189 +2 0 0 0.288591 +2 0 0 0.287316 +2 0 0 0.285416 +2 0 0 0.283881 +2 0 0 0.283562 +2 0 0 0.282169 +2 0 0 0.280629 +2 0 0 0.277582 +2 0 0 0.277191 +2 0 0 0.276431 +2 0 0 0.274611 +2 0 0 0.270996 +2 0 0 0.268806 +2 0 0 0.266635 +2 0 0 0.266592 +2 0 0 0.264483 +2 0 0 0.264334 +2 0 0 0.257997 +2 0 0 0.25725 +2 0 0 0.253227 +2 0 0 0.253163 +2 0 0 0.249345 +2 0 0 0.249339 +2 0 0 0.248879 +2 0 0 0.247042 +2 0 0 0.246491 +2 0 0 0.24509 +2 0 0 0.24348 +2 0 0 0.240632 +2 0 0 0.23879 +2 0 0 0.237571 +2 0 0 0.234969 +2 0 0 0.234112 +2 0 0 0.231675 +2 0 0 0.231179 +2 1 0.35058 0.580726 +2 0 0 0.228917 +2 0 0 0.2286 +2 0 0 0.227933 +2 0 0 0.227524 +2 0 0 0.223232 +2 0 0 0.218561 +2 0 0 0.218526 +2 0 0 0.215759 +2 0 0 0.214018 +2 1 0.347388 0.559978 +2 0 0 0.210372 +2 0 0 0.208851 +2 0 0 0.208063 +2 0 0 0.207714 +2 0 0 0.207261 +2 0 0 0.207016 +2 0 0 0.205389 +2 0 0 0.204544 +2 0 0 0.204208 +2 0 0 0.203812 +2 0 0 0.202564 +2 0 0 0.201961 +2 0 0 0.199673 +2 0 0 0.197377 +2 0 0 0.195722 +2 0 0 0.193804 +2 0 0 0.193791 +2 0 0 0.192733 +2 0 0 0.192406 +2 0 0 0.192356 +2 0 0 0.191661 +2 0 0 0.190642 +2 1 0.352356 0.541947 +2 0 0 0.189106 +2 0 0 0.184757 +2 0 0 0.18428 +2 0 0 0.179609 +2 1 0.301986 0.478334 +2 0 0 0.176306 +2 0 0 0.175016 +2 0 0 0.17493 +2 0 0 0.174172 +2 0 0 0.173684 +2 0 0 0.173586 +2 0 0 0.173412 +2 0 0 0.168587 +2 0 0 0.167455 +2 0 0 0.165907 +2 0 0 0.165592 +2 0 0 0.165494 +2 0 0 0.165314 +2 0 0 0.165202 +2 0 0 0.164815 +2 0 0 0.164809 +2 0 0 0.162464 +2 1 0.322683 0.482958 +2 0 0 0.159282 +2 0 0 0.158413 +2 0 0 0.157373 +2 0 0 0.157293 +2 0 0 0.157251 +2 0 0 0.1569 +2 0 0 0.15587 +2 0 0 0.154095 +2 0 0 0.15244 +2 0 0 0.151643 +2 0 0 0.150747 +2 0 0 0.149226 +2 0 0 0.149017 +2 0 0 0.148506 +2 0 0 0.147427 +2 0 0 0.147289 +2 0 0 0.147203 +2 0 0 0.145725 +2 0 0 0.145018 +2 0 0 0.143151 +2 0 0 0.143085 +2 0 0 0.142834 +2 0 0 0.142537 +2 0 0 0.142152 +2 0 0 0.141823 +2 0 0 0.138674 +2 1 0.377729 0.514046 +2 0 0 0.135402 +2 0 0 0.134886 +2 0 0 0.133722 +2 0 0 0.133525 +2 0 0 0.133033 +2 0 0 0.131638 +2 0 0 0.131438 +2 0 0 0.129741 +2 0 0 0.127864 +2 0 0 0.127383 +2 0 0 0.126595 +2 0 0 0.124774 +2 1 0.339394 0.463666 +2 0 0 0.123527 +2 0 0 0.122089 +2 0 0 0.121955 +2 0 0 0.119187 +2 1 0.2576 0.373749 +2 0 0 0.115939 +2 1 0.370051 0.483155 +2 0 0 0.111798 +2 0 0 0.111131 +2 0 0 0.109526 +2 0 0 0.108768 +2 1 0.386278 0.4922 +2 1 0.443911 0.54761 +2 1 0.374531 0.477153 +2 1 0.354747 0.454956 +2 0 0 0.0994522 +2 1 0.318032 0.417462 +2 0 0 0.0988752 +2 0 0 0.0970981 +2 0 0 0.0959033 +2 0 0 0.0884133 +2 0 0 0.0870489 +2 0 0 0.0857401 +2 0 0 0.084341 +2 0 0 0.0831521 +2 0 0 0.0831358 +2 0 0 0.0814685 +2 0 0 0.0802451 +2 0 0 0.0800348 +2 0 0 0.0792212 +2 0 0 0.0782455 +2 0 0 0.076909 +2 0 0 0.0753955 +2 1 0.324159 0.399094 +2 0 0 0.0746054 +2 0 0 0.0743236 +2 0 0 0.0741239 +2 1 0.413789 0.487379 +2 0 0 0.0728647 +2 0 0 0.0714516 +2 0 0 0.0702075 +2 0 0 0.0692368 +2 1 0.409968 0.478461 +2 0 0 0.0672609 +2 0 0 0.0651112 +2 0 0 0.0637053 +2 0 0 0.0633425 +2 1 0.298341 0.360737 +2 0 0 0.0623736 +2 1 0.309419 0.37151 +2 1 0.411549 0.471715 +2 0 0 0.0593071 +2 0 0 0.057488 +2 1 0.377019 0.434385 +2 1 0.329066 0.386148 +2 0 0 0.0560443 +2 0 0 0.0559886 +2 0 0 0.0522132 +2 0 0 0.0515157 +2 1 0.412572 0.46308 +2 0 0 0.0496187 +2 0 0 0.048826 +2 0 0 0.0487061 +2 1 0.240869 0.288747 +2 0 0 0.0471186 +2 1 0.531636 0.578093 +2 1 0.530723 0.576869 +2 1 0.33364 0.3797 +2 1 0.431628 0.477277 +2 0 0 0.0433426 +2 1 0.344005 0.387193 +2 1 0.463389 0.505345 +2 1 0.361715 0.403181 +2 1 0.381084 0.421374 +2 1 0.253337 0.292286 +2 0 0 0.0387451 +2 0 0 0.037954 +2 1 0.377147 0.414788 +2 1 0.308779 0.346297 +2 1 0.317521 0.349622 +2 0 0 0.0311057 +2 1 0.331259 0.356381 +2 1 0.279016 0.301331 +2 0 0 0.0217315 +2 1 0.350913 0.369543 +2 1 0.324571 0.343083 +2 0 0 0.0173134 +2 1 0.33392 0.35111 +2 1 0.542696 0.558863 +2 1 0.33836 0.350411 +2 1 0.29997 0.311993 +2 0 0 0.0120098 +2 1 0.189652 0.197515 +2 0 0 0.00391946 +2 0 0 0.0036753 +2 1 0.445398 0.448892 +2 1 0.24546 0.248903 +2 1 0.301045 0.303943 +2 1 0.309336 0.310216 + 0.56 real 0.52 user 0.02 sys + 69439488 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 16943 page reclaims + 29 page faults + 0 swaps + 2 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 10 voluntary context switches + 73 involuntary context switches diff --git a/benchmarks/sphere_3_192/gudhi.dim_2.out.txt b/benchmarks/sphere_3_192/gudhi.dim_2.out.txt new file mode 100644 index 0000000..36740c8 --- /dev/null +++ b/benchmarks/sphere_3_192/gudhi.dim_2.out.txt @@ -0,0 +1,263 @@ +The complex contains 56050288 simplices + and has dimension 3 +2 0 0 inf +2 2 0.720484 1.65562 +2 1 0.316534 0.682559 +2 1 0.304761 0.647803 +2 0 0 0.332695 +2 0 0 0.323964 +2 0 0 0.319493 +2 0 0 0.316023 +2 1 0.34298 0.65758 +2 0 0 0.304062 +2 0 0 0.297933 +2 0 0 0.297053 +2 0 0 0.293466 +2 1 0.347806 0.638039 +2 0 0 0.289189 +2 0 0 0.288591 +2 0 0 0.287316 +2 0 0 0.285416 +2 0 0 0.283881 +2 0 0 0.283562 +2 0 0 0.282169 +2 0 0 0.280629 +2 0 0 0.277582 +2 0 0 0.277191 +2 0 0 0.276431 +2 0 0 0.274611 +2 0 0 0.270996 +2 0 0 0.268806 +2 0 0 0.266635 +2 0 0 0.266592 +2 0 0 0.264483 +2 0 0 0.264334 +2 0 0 0.257997 +2 0 0 0.25725 +2 0 0 0.253227 +2 0 0 0.253163 +2 0 0 0.249345 +2 0 0 0.249339 +2 0 0 0.248879 +2 0 0 0.247042 +2 0 0 0.246491 +2 0 0 0.24509 +2 0 0 0.24348 +2 0 0 0.240632 +2 0 0 0.23879 +2 0 0 0.237571 +2 0 0 0.234969 +2 0 0 0.234112 +2 0 0 0.231675 +2 0 0 0.231179 +2 1 0.35058 0.580726 +2 0 0 0.228917 +2 0 0 0.2286 +2 0 0 0.227933 +2 0 0 0.227524 +2 0 0 0.223232 +2 0 0 0.218561 +2 0 0 0.218526 +2 0 0 0.215759 +2 0 0 0.214018 +2 1 0.347388 0.559978 +2 0 0 0.210372 +2 0 0 0.208851 +2 0 0 0.208063 +2 0 0 0.207714 +2 0 0 0.207261 +2 0 0 0.207016 +2 0 0 0.205389 +2 0 0 0.204544 +2 0 0 0.204208 +2 0 0 0.203812 +2 0 0 0.202564 +2 0 0 0.201961 +2 0 0 0.199673 +2 0 0 0.197377 +2 0 0 0.195722 +2 0 0 0.193804 +2 0 0 0.193791 +2 0 0 0.192733 +2 0 0 0.192406 +2 0 0 0.192356 +2 0 0 0.191661 +2 0 0 0.190642 +2 1 0.352356 0.541947 +2 0 0 0.189106 +2 0 0 0.184757 +2 0 0 0.18428 +2 0 0 0.179609 +2 1 0.301986 0.478334 +2 0 0 0.176306 +2 0 0 0.175016 +2 0 0 0.17493 +2 0 0 0.174172 +2 0 0 0.173684 +2 0 0 0.173586 +2 0 0 0.173412 +2 0 0 0.168587 +2 0 0 0.167455 +2 0 0 0.165907 +2 0 0 0.165592 +2 0 0 0.165494 +2 0 0 0.165314 +2 0 0 0.165202 +2 0 0 0.164815 +2 0 0 0.164809 +2 0 0 0.162464 +2 1 0.322683 0.482958 +2 0 0 0.159282 +2 0 0 0.158413 +2 0 0 0.157373 +2 0 0 0.157293 +2 0 0 0.157251 +2 0 0 0.1569 +2 0 0 0.15587 +2 0 0 0.154095 +2 0 0 0.15244 +2 0 0 0.151643 +2 0 0 0.150747 +2 0 0 0.149226 +2 0 0 0.149017 +2 0 0 0.148506 +2 0 0 0.147427 +2 0 0 0.147289 +2 0 0 0.147203 +2 0 0 0.145725 +2 0 0 0.145018 +2 0 0 0.143151 +2 0 0 0.143085 +2 0 0 0.142834 +2 0 0 0.142537 +2 0 0 0.142152 +2 0 0 0.141823 +2 0 0 0.138674 +2 1 0.377729 0.514046 +2 0 0 0.135402 +2 0 0 0.134886 +2 0 0 0.133722 +2 0 0 0.133525 +2 0 0 0.133033 +2 0 0 0.131638 +2 0 0 0.131438 +2 0 0 0.129741 +2 0 0 0.127864 +2 0 0 0.127383 +2 0 0 0.126595 +2 0 0 0.124774 +2 1 0.339394 0.463666 +2 0 0 0.123527 +2 0 0 0.122089 +2 0 0 0.121955 +2 0 0 0.119187 +2 1 0.2576 0.373749 +2 0 0 0.115939 +2 1 0.370051 0.483155 +2 0 0 0.111798 +2 0 0 0.111131 +2 0 0 0.109526 +2 0 0 0.108768 +2 1 0.386278 0.4922 +2 1 0.443911 0.54761 +2 1 0.374531 0.477153 +2 1 0.354747 0.454956 +2 0 0 0.0994522 +2 1 0.318032 0.417462 +2 0 0 0.0988752 +2 0 0 0.0970981 +2 0 0 0.0959033 +2 0 0 0.0884133 +2 0 0 0.0870489 +2 0 0 0.0857401 +2 0 0 0.084341 +2 0 0 0.0831521 +2 0 0 0.0831358 +2 0 0 0.0814685 +2 0 0 0.0802451 +2 0 0 0.0800348 +2 0 0 0.0792212 +2 0 0 0.0782455 +2 0 0 0.076909 +2 0 0 0.0753955 +2 1 0.324159 0.399094 +2 0 0 0.0746054 +2 0 0 0.0743236 +2 0 0 0.0741239 +2 1 0.413789 0.487379 +2 0 0 0.0728647 +2 0 0 0.0714516 +2 0 0 0.0702075 +2 0 0 0.0692368 +2 1 0.409968 0.478461 +2 0 0 0.0672609 +2 0 0 0.0651112 +2 0 0 0.0637053 +2 0 0 0.0633425 +2 1 0.298341 0.360737 +2 0 0 0.0623736 +2 1 0.309419 0.37151 +2 1 0.411549 0.471715 +2 0 0 0.0593071 +2 0 0 0.057488 +2 1 0.377019 0.434385 +2 1 0.329066 0.386148 +2 0 0 0.0560443 +2 0 0 0.0559886 +2 0 0 0.0522132 +2 0 0 0.0515157 +2 1 0.412572 0.46308 +2 0 0 0.0496187 +2 0 0 0.048826 +2 0 0 0.0487061 +2 1 0.240869 0.288747 +2 0 0 0.0471186 +2 1 0.531636 0.578093 +2 1 0.530723 0.576869 +2 1 0.33364 0.3797 +2 1 0.431628 0.477277 +2 0 0 0.0433426 +2 1 0.344005 0.387193 +2 1 0.463389 0.505345 +2 1 0.361715 0.403181 +2 1 0.381084 0.421374 +2 1 0.253337 0.292286 +2 0 0 0.0387451 +2 0 0 0.037954 +2 1 0.377147 0.414788 +2 1 0.308779 0.346297 +2 1 0.317521 0.349622 +2 0 0 0.0311057 +2 1 0.331259 0.356381 +2 1 0.279016 0.301331 +2 0 0 0.0217315 +2 1 0.350913 0.369543 +2 1 0.324571 0.343083 +2 0 0 0.0173134 +2 1 0.33392 0.35111 +2 1 0.542696 0.558863 +2 1 0.33836 0.350411 +2 1 0.29997 0.311993 +2 0 0 0.0120098 +2 1 0.189652 0.197515 +2 0 0 0.00391946 +2 0 0 0.0036753 +2 1 0.445398 0.448892 +2 1 0.24546 0.248903 +2 1 0.301045 0.303943 +2 1 0.309336 0.310216 + 75.58 real 74.18 user 1.35 sys +2911698944 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 820365 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 11 voluntary context switches + 9107 involuntary context switches diff --git a/benchmarks/sphere_3_192/ripser.dim_1.out.txt b/benchmarks/sphere_3_192/ripser.dim_1.out.txt new file mode 100644 index 0000000..f7f631f --- /dev/null +++ b/benchmarks/sphere_3_192/ripser.dim_1.out.txt @@ -0,0 +1,264 @@ +distance matrix with 192 points +value range: [0.00367531,2] +persistence intervals in dim 0: + [0,0.00367531) + [0,0.00391946) + [0,0.0120098) + [0,0.0173134) + [0,0.0217315) + [0,0.0311057) + [0,0.037954) + [0,0.0387451) + [0,0.0433426) + [0,0.0471186) + [0,0.0487061) + [0,0.048826) + [0,0.0496187) + [0,0.0515157) + [0,0.0522132) + [0,0.0559886) + [0,0.0560443) + [0,0.057488) + [0,0.0593071) + [0,0.0623736) + [0,0.0633425) + [0,0.0637053) + [0,0.0651112) + [0,0.0672609) + [0,0.0692368) + [0,0.0702075) + [0,0.0714516) + [0,0.0728647) + [0,0.0741239) + [0,0.0743236) + [0,0.0746054) + [0,0.0753955) + [0,0.076909) + [0,0.0782455) + [0,0.0792212) + [0,0.0800348) + [0,0.0802451) + [0,0.0814685) + [0,0.0831358) + [0,0.0831521) + [0,0.084341) + [0,0.0857401) + [0,0.0870489) + [0,0.0884133) + [0,0.0959033) + [0,0.0970981) + [0,0.0988752) + [0,0.0994522) + [0,0.108768) + [0,0.109526) + [0,0.111131) + [0,0.111798) + [0,0.115939) + [0,0.119187) + [0,0.121955) + [0,0.122089) + [0,0.123527) + [0,0.124774) + [0,0.126595) + [0,0.127383) + [0,0.127864) + [0,0.129741) + [0,0.131438) + [0,0.131638) + [0,0.133033) + [0,0.133525) + [0,0.133722) + [0,0.134886) + [0,0.135402) + [0,0.138674) + [0,0.141823) + [0,0.142152) + [0,0.142537) + [0,0.142834) + [0,0.143085) + [0,0.143151) + [0,0.145018) + [0,0.145725) + [0,0.147203) + [0,0.147289) + [0,0.147427) + [0,0.148506) + [0,0.149017) + [0,0.149226) + [0,0.150747) + [0,0.151643) + [0,0.15244) + [0,0.154095) + [0,0.15587) + [0,0.1569) + [0,0.157251) + [0,0.157293) + [0,0.157373) + [0,0.158413) + [0,0.159282) + [0,0.162464) + [0,0.164809) + [0,0.164815) + [0,0.165202) + [0,0.165314) + [0,0.165494) + [0,0.165592) + [0,0.165907) + [0,0.167455) + [0,0.168587) + [0,0.173412) + [0,0.173586) + [0,0.173684) + [0,0.174172) + [0,0.17493) + [0,0.175016) + [0,0.176306) + [0,0.179609) + [0,0.18428) + [0,0.184757) + [0,0.189106) + [0,0.190642) + [0,0.191661) + [0,0.192356) + [0,0.192406) + [0,0.192733) + [0,0.193791) + [0,0.193804) + [0,0.195722) + [0,0.197377) + [0,0.199673) + [0,0.201961) + [0,0.202564) + [0,0.203812) + [0,0.204208) + [0,0.204544) + [0,0.205389) + [0,0.207016) + [0,0.207261) + [0,0.207714) + [0,0.208063) + [0,0.208851) + [0,0.210372) + [0,0.214018) + [0,0.215759) + [0,0.218526) + [0,0.218561) + [0,0.223232) + [0,0.227524) + [0,0.227933) + [0,0.2286) + [0,0.228917) + [0,0.231179) + [0,0.231675) + [0,0.234112) + [0,0.234969) + [0,0.237571) + [0,0.23879) + [0,0.240632) + [0,0.24348) + [0,0.24509) + [0,0.246491) + [0,0.247042) + [0,0.248879) + [0,0.249339) + [0,0.249345) + [0,0.253163) + [0,0.253227) + [0,0.25725) + [0,0.257997) + [0,0.264334) + [0,0.264483) + [0,0.266592) + [0,0.266635) + [0,0.268806) + [0,0.270996) + [0,0.274611) + [0,0.276431) + [0,0.277191) + [0,0.277582) + [0,0.280629) + [0,0.282169) + [0,0.283562) + [0,0.283881) + [0,0.285416) + [0,0.287316) + [0,0.288591) + [0,0.289189) + [0,0.293466) + [0,0.297053) + [0,0.297933) + [0,0.304062) + [0,0.316023) + [0,0.319493) + [0,0.323964) + [0,0.332695) + [0, ) +persistence intervals in dim 1: + [0.542696,0.558863) + [0.531636,0.578093) + [0.530723,0.576869) + [0.463389,0.505345) + [0.445398,0.448892) + [0.443911,0.54761) + [0.431628,0.477277) + [0.413789,0.487379) + [0.412572,0.46308) + [0.411549,0.471715) + [0.409968,0.478461) + [0.386278,0.4922) + [0.381084,0.421374) + [0.377729,0.514046) + [0.377147,0.414788) + [0.377019,0.434385) + [0.374531,0.477153) + [0.370051,0.483155) + [0.361715,0.403181) + [0.354747,0.454956) + [0.352356,0.541947) + [0.350913,0.369543) + [0.35058,0.580726) + [0.347806,0.638039) + [0.347388,0.559978) + [0.344005,0.387193) + [0.34298,0.65758) + [0.339394,0.463666) + [0.33836,0.350411) + [0.33392,0.35111) + [0.33364,0.3797) + [0.331259,0.356381) + [0.329066,0.386148) + [0.324571,0.343083) + [0.324159,0.399094) + [0.322683,0.482958) + [0.318032,0.417462) + [0.317521,0.349622) + [0.316534,0.682559) + [0.309419,0.37151) + [0.309336,0.310216) + [0.308779,0.346297) + [0.304761,0.647803) + [0.301986,0.478334) + [0.301045,0.303943) + [0.29997,0.311993) + [0.298341,0.360737) + [0.279016,0.301331) + [0.2576,0.373749) + [0.253337,0.292286) + [0.24546,0.248903) + [0.240869,0.288747) + [0.189652,0.197515) + 0.03 real 0.02 user 0.00 sys + 3379200 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 844 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 4 voluntary context switches + 21 involuntary context switches diff --git a/benchmarks/sphere_3_192/ripser.dim_2.out.txt b/benchmarks/sphere_3_192/ripser.dim_2.out.txt new file mode 100644 index 0000000..b3283e7 --- /dev/null +++ b/benchmarks/sphere_3_192/ripser.dim_2.out.txt @@ -0,0 +1,266 @@ +distance matrix with 192 points +value range: [0.00367531,2] +persistence intervals in dim 0: + [0,0.00367531) + [0,0.00391946) + [0,0.0120098) + [0,0.0173134) + [0,0.0217315) + [0,0.0311057) + [0,0.037954) + [0,0.0387451) + [0,0.0433426) + [0,0.0471186) + [0,0.0487061) + [0,0.048826) + [0,0.0496187) + [0,0.0515157) + [0,0.0522132) + [0,0.0559886) + [0,0.0560443) + [0,0.057488) + [0,0.0593071) + [0,0.0623736) + [0,0.0633425) + [0,0.0637053) + [0,0.0651112) + [0,0.0672609) + [0,0.0692368) + [0,0.0702075) + [0,0.0714516) + [0,0.0728647) + [0,0.0741239) + [0,0.0743236) + [0,0.0746054) + [0,0.0753955) + [0,0.076909) + [0,0.0782455) + [0,0.0792212) + [0,0.0800348) + [0,0.0802451) + [0,0.0814685) + [0,0.0831358) + [0,0.0831521) + [0,0.084341) + [0,0.0857401) + [0,0.0870489) + [0,0.0884133) + [0,0.0959033) + [0,0.0970981) + [0,0.0988752) + [0,0.0994522) + [0,0.108768) + [0,0.109526) + [0,0.111131) + [0,0.111798) + [0,0.115939) + [0,0.119187) + [0,0.121955) + [0,0.122089) + [0,0.123527) + [0,0.124774) + [0,0.126595) + [0,0.127383) + [0,0.127864) + [0,0.129741) + [0,0.131438) + [0,0.131638) + [0,0.133033) + [0,0.133525) + [0,0.133722) + [0,0.134886) + [0,0.135402) + [0,0.138674) + [0,0.141823) + [0,0.142152) + [0,0.142537) + [0,0.142834) + [0,0.143085) + [0,0.143151) + [0,0.145018) + [0,0.145725) + [0,0.147203) + [0,0.147289) + [0,0.147427) + [0,0.148506) + [0,0.149017) + [0,0.149226) + [0,0.150747) + [0,0.151643) + [0,0.15244) + [0,0.154095) + [0,0.15587) + [0,0.1569) + [0,0.157251) + [0,0.157293) + [0,0.157373) + [0,0.158413) + [0,0.159282) + [0,0.162464) + [0,0.164809) + [0,0.164815) + [0,0.165202) + [0,0.165314) + [0,0.165494) + [0,0.165592) + [0,0.165907) + [0,0.167455) + [0,0.168587) + [0,0.173412) + [0,0.173586) + [0,0.173684) + [0,0.174172) + [0,0.17493) + [0,0.175016) + [0,0.176306) + [0,0.179609) + [0,0.18428) + [0,0.184757) + [0,0.189106) + [0,0.190642) + [0,0.191661) + [0,0.192356) + [0,0.192406) + [0,0.192733) + [0,0.193791) + [0,0.193804) + [0,0.195722) + [0,0.197377) + [0,0.199673) + [0,0.201961) + [0,0.202564) + [0,0.203812) + [0,0.204208) + [0,0.204544) + [0,0.205389) + [0,0.207016) + [0,0.207261) + [0,0.207714) + [0,0.208063) + [0,0.208851) + [0,0.210372) + [0,0.214018) + [0,0.215759) + [0,0.218526) + [0,0.218561) + [0,0.223232) + [0,0.227524) + [0,0.227933) + [0,0.2286) + [0,0.228917) + [0,0.231179) + [0,0.231675) + [0,0.234112) + [0,0.234969) + [0,0.237571) + [0,0.23879) + [0,0.240632) + [0,0.24348) + [0,0.24509) + [0,0.246491) + [0,0.247042) + [0,0.248879) + [0,0.249339) + [0,0.249345) + [0,0.253163) + [0,0.253227) + [0,0.25725) + [0,0.257997) + [0,0.264334) + [0,0.264483) + [0,0.266592) + [0,0.266635) + [0,0.268806) + [0,0.270996) + [0,0.274611) + [0,0.276431) + [0,0.277191) + [0,0.277582) + [0,0.280629) + [0,0.282169) + [0,0.283562) + [0,0.283881) + [0,0.285416) + [0,0.287316) + [0,0.288591) + [0,0.289189) + [0,0.293466) + [0,0.297053) + [0,0.297933) + [0,0.304062) + [0,0.316023) + [0,0.319493) + [0,0.323964) + [0,0.332695) + [0, ) +persistence intervals in dim 1: + [0.542696,0.558863) + [0.531636,0.578093) + [0.530723,0.576869) + [0.463389,0.505345) + [0.445398,0.448892) + [0.443911,0.54761) + [0.431628,0.477277) + [0.413789,0.487379) + [0.412572,0.46308) + [0.411549,0.471715) + [0.409968,0.478461) + [0.386278,0.4922) + [0.381084,0.421374) + [0.377729,0.514046) + [0.377147,0.414788) + [0.377019,0.434385) + [0.374531,0.477153) + [0.370051,0.483155) + [0.361715,0.403181) + [0.354747,0.454956) + [0.352356,0.541947) + [0.350913,0.369543) + [0.35058,0.580726) + [0.347806,0.638039) + [0.347388,0.559978) + [0.344005,0.387193) + [0.34298,0.65758) + [0.339394,0.463666) + [0.33836,0.350411) + [0.33392,0.35111) + [0.33364,0.3797) + [0.331259,0.356381) + [0.329066,0.386148) + [0.324571,0.343083) + [0.324159,0.399094) + [0.322683,0.482958) + [0.318032,0.417462) + [0.317521,0.349622) + [0.316534,0.682559) + [0.309419,0.37151) + [0.309336,0.310216) + [0.308779,0.346297) + [0.304761,0.647803) + [0.301986,0.478334) + [0.301045,0.303943) + [0.29997,0.311993) + [0.298341,0.360737) + [0.279016,0.301331) + [0.2576,0.373749) + [0.253337,0.292286) + [0.24546,0.248903) + [0.240869,0.288747) + [0.189652,0.197515) +persistence intervals in dim 2: + [0.720484,1.65562) + 1.17 real 1.08 user 0.08 sys + 182210560 maximum resident set size + 0 average shared memory size + 0 average unshared data size + 0 average unshared stack size + 44505 page reclaims + 0 page faults + 0 swaps + 0 block input operations + 0 block output operations + 0 messages sent + 0 messages received + 0 signals received + 2 voluntary context switches + 69 involuntary context switches diff --git a/benchmarks/sphere_3_192/run_dionysus.sh b/benchmarks/sphere_3_192/run_dionysus.sh index 1e5c7d3..7af804f 100644 --- a/benchmarks/sphere_3_192/run_dionysus.sh +++ b/benchmarks/sphere_3_192/run_dionysus.sh @@ -1 +1,3 @@ -/usr/bin/time -l ~/Source/Dionysus/examples/cohomology/rips-pairwise-cohomology -s2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee dionysus.out.txt +/usr/bin/time -l ~/Source/Dionysus/examples/cohomology/rips-pairwise-cohomology -s2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee dionysus.dim_1.out.txt +/usr/bin/time -l ~/Source/Dionysus/examples/cohomology/rips-pairwise-cohomology -s3 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee dionysus.dim_2.out.txt + diff --git a/benchmarks/sphere_3_192/run_dipha.sh b/benchmarks/sphere_3_192/run_dipha.sh index 312bcd7..b700a1a 100644 --- a/benchmarks/sphere_3_192/run_dipha.sh +++ b/benchmarks/sphere_3_192/run_dipha.sh @@ -1 +1,3 @@ -/usr/bin/time -l ~/Bitbucket/dipha/dipha --benchmark --upper_dim 2 --dual ~/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex /dev/null 2>&1 | tee dipha.out.txt +/usr/bin/time -l ~/Bitbucket/dipha/dipha --benchmark --upper_dim 2 --dual ~/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex /dev/null 2>&1 | tee dipha.dim_1.out.txt +/usr/bin/time -l ~/Bitbucket/dipha/dipha --benchmark --upper_dim 3 --dual ~/Bitbucket/phat-paper/benchmark/dipha/sphere_3_192.complex /dev/null 2>&1 | tee dipha.dim_2.out.txt + diff --git a/benchmarks/sphere_3_192/run_gudhi.sh b/benchmarks/sphere_3_192/run_gudhi.sh index 80ae403..7e39da9 100644 --- a/benchmarks/sphere_3_192/run_gudhi.sh +++ b/benchmarks/sphere_3_192/run_gudhi.sh @@ -1 +1,3 @@ -/usr/bin/time -l ~/Source/Gudhi_library_1.3.1/example/Persistent_cohomology/rips_persistence -d2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee gudhi.out.txt +/usr/bin/time -l ~/Source/Gudhi_library_1.3.1/example/Persistent_cohomology/rips_persistence -d2 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee gudhi.dim_1.out.txt +/usr/bin/time -l ~/Source/Gudhi_library_1.3.1/example/Persistent_cohomology/rips_persistence -d3 -p2 ~/Bitbucket/phat-paper/benchmark/point\ cloud/sphere_3_192_points.dat 2>&1 | tee gudhi.dim_2.out.txt + diff --git a/benchmarks/sphere_3_192/run_ripser.sh b/benchmarks/sphere_3_192/run_ripser.sh index 209c410..4cc1c3d 100644 --- a/benchmarks/sphere_3_192/run_ripser.sh +++ b/benchmarks/sphere_3_192/run_ripser.sh @@ -1 +1,3 @@ -/usr/bin/time -l ~/Bitbucket/ripser/ripser ~/Bitbucket/ripser/examples/sphere_3_192.lower_distance_matrix 2>&1 | tee ripser.out.txt +/usr/bin/time -l ~/Bitbucket/ripser/ripser ~/Bitbucket/ripser/examples/sphere_3_192.lower_distance_matrix 2>&1 | tee ripser.dim_1.out.txt +/usr/bin/time -l ~/Bitbucket/ripser/ripser ~/Bitbucket/ripser/examples/sphere_3_192.lower_distance_matrix --dim 2 2>&1 | tee ripser.dim_2.out.txt + -- cgit v1.2.3 From 92536895274a4632eaff452c56227caa2e3e0e66 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Mon, 20 Mar 2017 16:07:33 +0100 Subject: benchmark using sparsehash --- Makefile | 10 +++++----- benchmarks/sphere_3_192/ripser.dim_1.out.txt | 10 +++++----- benchmarks/sphere_3_192/ripser.dim_2.out.txt | 8 ++++---- 3 files changed, 14 insertions(+), 14 deletions(-) diff --git a/Makefile b/Makefile index 893d776..90543db 100644 --- a/Makefile +++ b/Makefile @@ -5,19 +5,19 @@ all: ripser ripser-coeff ripser-reduction ripser-coeff-reduction ripser-debug ripser: ripser.cpp - c++ -std=c++11 ripser.cpp -o ripser -Ofast -D NDEBUG + c++ -std=c++11 ripser.cpp -o ripser -Ofast -D NDEBUG -D USE_GOOGLE_HASHMAP ripser-coeff: ripser.cpp - c++ -std=c++11 ripser.cpp -o ripser-coeff -Ofast -D NDEBUG -D USE_COEFFICIENTS + c++ -std=c++11 ripser.cpp -o ripser-coeff -Ofast -D NDEBUG -D USE_GOOGLE_HASHMAP -D USE_COEFFICIENTS ripser-reduction: ripser.cpp - c++ -std=c++11 ripser.cpp -o ripser-reduction -Ofast -D NDEBUG -D ASSEMBLE_REDUCTION_MATRIX + c++ -std=c++11 ripser.cpp -o ripser-reduction -Ofast -D NDEBUG -D USE_GOOGLE_HASHMAP -D ASSEMBLE_REDUCTION_MATRIX ripser-coeff-reduction: ripser.cpp - c++ -std=c++11 ripser.cpp -o ripser-coeff-reduction -Ofast -D NDEBUG -D USE_COEFFICIENTS -D ASSEMBLE_REDUCTION_MATRIX + c++ -std=c++11 ripser.cpp -o ripser-coeff-reduction -Ofast -D NDEBUG -D USE_GOOGLE_HASHMAP -D USE_COEFFICIENTS -D ASSEMBLE_REDUCTION_MATRIX ripser-debug: ripser.cpp - c++ -std=c++11 ripser.cpp -o ripser-debug -g + c++ -std=c++11 ripser.cpp -o ripser-debug -g -D USE_GOOGLE_HASHMAP clean: diff --git a/benchmarks/sphere_3_192/ripser.dim_1.out.txt b/benchmarks/sphere_3_192/ripser.dim_1.out.txt index f7f631f..4ef9cb2 100644 --- a/benchmarks/sphere_3_192/ripser.dim_1.out.txt +++ b/benchmarks/sphere_3_192/ripser.dim_1.out.txt @@ -247,12 +247,12 @@ persistence intervals in dim 1: [0.24546,0.248903) [0.240869,0.288747) [0.189652,0.197515) - 0.03 real 0.02 user 0.00 sys - 3379200 maximum resident set size + 0.02 real 0.02 user 0.00 sys + 2682880 maximum resident set size 0 average shared memory size 0 average unshared data size 0 average unshared stack size - 844 page reclaims + 674 page reclaims 0 page faults 0 swaps 0 block input operations @@ -260,5 +260,5 @@ persistence intervals in dim 1: 0 messages sent 0 messages received 0 signals received - 4 voluntary context switches - 21 involuntary context switches + 5 voluntary context switches + 4 involuntary context switches diff --git a/benchmarks/sphere_3_192/ripser.dim_2.out.txt b/benchmarks/sphere_3_192/ripser.dim_2.out.txt index b3283e7..7ba2de6 100644 --- a/benchmarks/sphere_3_192/ripser.dim_2.out.txt +++ b/benchmarks/sphere_3_192/ripser.dim_2.out.txt @@ -249,12 +249,12 @@ persistence intervals in dim 1: [0.189652,0.197515) persistence intervals in dim 2: [0.720484,1.65562) - 1.17 real 1.08 user 0.08 sys - 182210560 maximum resident set size + 1.17 real 1.10 user 0.06 sys + 151293952 maximum resident set size 0 average shared memory size 0 average unshared data size 0 average unshared stack size - 44505 page reclaims + 36957 page reclaims 0 page faults 0 swaps 0 block input operations @@ -263,4 +263,4 @@ persistence intervals in dim 2: 0 messages received 0 signals received 2 voluntary context switches - 69 involuntary context switches + 49 involuntary context switches -- cgit v1.2.3 From d898ef27044bc7f10da478d16581915e95edc2c9 Mon Sep 17 00:00:00 2001 From: Ulrich Bauer Date: Thu, 30 Mar 2017 08:36:40 +0200 Subject: progress indicator --- ripser.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/ripser.cpp b/ripser.cpp index 53fa4fa..d3aa615 100644 --- a/ripser.cpp +++ b/ripser.cpp @@ -489,7 +489,7 @@ public: if (diameter <= threshold) columns_to_reduce.push_back(std::make_pair(diameter, index)); #ifdef INDICATE_PROGRESS - if ((index + 1) % 1000 == 0) + if ((index + 1) % 1000000 == 0) std::cout << "\033[K" << "assembled " << columns_to_reduce.size() << " out of " << (index + 1) << "/" << num_simplices << " columns" << std::flush @@ -543,7 +543,7 @@ public: value_t diameter = get_diameter(column_to_reduce); #ifdef INDICATE_PROGRESS - //if ((i + 1) % 1000 == 0) + if ((i + 1) % 1000000 == 0) std::cout << "\033[K" << "reducing column " << i + 1 << "/" << columns_to_reduce.size() << " (diameter " << diameter << ")" << std::flush << "\r"; -- cgit v1.2.3