diff --git a/src/Simplex_tree/concept/FiltrationValue.h b/src/Simplex_tree/concept/FiltrationValue.h index 6cf314fa9d..ffd67bdb27 100644 --- a/src/Simplex_tree/concept/FiltrationValue.h +++ b/src/Simplex_tree/concept/FiltrationValue.h @@ -5,19 +5,102 @@ * Copyright (C) 2014 Inria * * Modification(s): + * - 2024/08 Hannah Schreiber: Update of the concept after several additions to the Simplex tree. * - YYYY/MM Author: Description of the modification */ /** \brief Value type for a filtration function on a cell complex. - * - * A filtration of a cell complex (see FilteredComplex) is - * a function \f$f:\mathbf{K} \rightarrow \mathbb{R}\f$ satisfying \f$f(\tau)\leq - * f(\sigma)\f$ whenever \f$\tau \subseteq \sigma\f$. Ordering the simplices - * by increasing filtration values (breaking ties so as a simplex appears after - * its subsimplices of same filtration value) provides an indexing scheme - * (see IndexingTag). - */ - struct FiltrationValue { - /** \brief Operator < is a StrictWeakOrdering. */ - bool operator<(FiltrationValue f1, FiltrationValue f2); - }; + * + * Needs to implement `std::numeric_limits::has_infinity`, + * `std::numeric_limits::infinity()` and `std::numeric_limits::max()`. + * But when `std::numeric_limits::has_infinity` returns `true`, + * `std::numeric_limits::max()` can simply throw when called, as well as, + * `std::numeric_limits::infinity()` if `std::numeric_limits::has_infinity` + * returns `false`. + * + * A filtration of a cell complex (see FilteredComplex) is + * a function \f$f:\mathbf{K} \rightarrow \mathbb{R}\f$ satisfying \f$f(\tau)\leq + * f(\sigma)\f$ whenever \f$\tau \subseteq \sigma\f$. Ordering the simplices + * by increasing filtration values (breaking ties so as a simplex appears after + * its subsimplices of same filtration value) provides an indexing scheme + * (see IndexingTag). + */ +struct FiltrationValue { + /** + * @brief Has to construct the default value of FiltrationValue. + */ + FiltrationValue(); + /** + * @brief Has to be able to take `0` as input at construction and in this case, construct an empty object. + * E.g., 0 for a numerical value or {} for a vector. + */ + FiltrationValue(Any_arithmetic_type v); + + // only for default ordering of filtration_vect_ in initialize_filtration and for prune_above_filtration + /** + * @brief Strictly smaller operator. If the filtration values are totally ordered, should be a StrictWeakOrdering. + */ + friend bool operator<(const FiltrationValue& f1, const FiltrationValue& f2); + /** + * @brief Equality operator + */ + friend bool operator==(const FiltrationValue& f1, const FiltrationValue& f2); + + /** + * @brief Given two filtration values at which a simplex exists, computes the minimal union of births generating + * a lifetime including those two values. The result is stored in the first parameter. + * The overload for arithmetic types like `double` or `int` is already implemented as the minimum of the + * two given values and can also be used for non native arithmetic types like `CGAL::Gmpq` as long as it has an + * `operator<`. The overload is available with @ref Gudhi::unify_lifetimes "". + * + * For a k-critical filtration, FiltrationValue should be able to store an union of values (corresponding to the + * different births of a same simplex) and this method adds the values of @p f2 in @p f1 and removes the values + * from @p f1 which are comparable and greater than other values. + * In the special case of 1-critical filtration, as the union should not contain more than one birth element, + * this method is expected to throw if the two given elements in the filtration values are not comparable. + * If they are comparable, the union is simply the minimum of both. + * + * @return True if and only if the values in @p f1 were actually modified. + */ + friend bool unify_lifetimes(FiltrationValue& f1, const FiltrationValue& f2); + + /** + * @brief Given two filtration values, stores in the first value the lowest common upper bound of the two values. + * The overload for arithmetic types like `double` or `int` is already implemented as the maximum of the two + * given values and can also be used for non native arithmetic types like `CGAL::Gmpq` as long as it has an + * `operator<`. The overload is available with @ref Gudhi::intersect_lifetimes "". + * + * @return True if and only if the values in @p f1 were actually modified. + */ + friend bool intersect_lifetimes(FiltrationValue& f1, const FiltrationValue& f2); + + /** + * @brief Only necessary when serializing the simplex tree. Serialize the given value and insert it at start position. + * Overloads for native arithmetic types or other simple types are already implemented with + * @ref Gudhi::simplex_tree::serialize_trivial "". + * + * @param value The value to serialize. + * @param start Start position where the value is serialized. + * @return The new position in the array of char for the next serialization. + */ + friend char* serialize_trivial(const FiltrationValue& value, char* start); + + /** + * @brief Only necessary when deserializing the simplex tree. Deserialize at the start position in an array of char + * and sets the value with it. + * Overloads for native arithmetic types or other simple types are already implemented with + * @ref Gudhi::simplex_tree::deserialize_trivial "". + * + * @param value The value where to deserialize based on its type. + * @param start Start position where the value is serialized. + * @return The new position in the array of char for the next deserialization. + */ + friend const char* deserialize_trivial(FiltrationValue& value, const char* start); + + /** + * @brief Only necessary when (de)serializing the simplex tree. Returns the serialization size of the given object. + * Overloads for native arithmetic types or other simple types are already implemented with + * @ref Gudhi::simplex_tree::get_serialization_size_of "". + */ + friend std::size_t get_serialization_size_of(const FiltrationValue& value); +}; diff --git a/src/Simplex_tree/concept/SimplexTreeOptions.h b/src/Simplex_tree/concept/SimplexTreeOptions.h index 70f936a2f6..96dfcade73 100644 --- a/src/Simplex_tree/concept/SimplexTreeOptions.h +++ b/src/Simplex_tree/concept/SimplexTreeOptions.h @@ -19,7 +19,7 @@ struct SimplexTreeOptions { typedef IndexingTag Indexing_tag; /** @brief Must be a signed integer type. It admits a total order <. */ typedef VertexHandle Vertex_handle; - /** @brief Must be comparable with operator<. */ + /** @brief Filtration value type which should implement the @ref FiltrationValue concept. */ typedef FiltrationValue Filtration_value; /** @brief Must be an integer type. */ typedef SimplexKey Simplex_key; diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree.h b/src/Simplex_tree/include/gudhi/Simplex_tree.h index 45a955dd87..bc85900fcc 100644 --- a/src/Simplex_tree/include/gudhi/Simplex_tree.h +++ b/src/Simplex_tree/include/gudhi/Simplex_tree.h @@ -10,6 +10,7 @@ * - 2023/05 Clément Maria: Edge insertion method for flag complexes * - 2023/05 Hannah Schreiber: Factorization of expansion methods * - 2023/08 Hannah Schreiber (& Clément Maria): Add possibility of stable simplex handles. + * - 2024/08 Hannah Schreiber: Generalization of the notion of filtration values. * - 2024/08 Hannah Schreiber: Addition of customizable copy constructor. * - YYYY/MM Author: Description of the modification */ @@ -99,7 +100,7 @@ class Simplex_tree { typedef typename Options::Indexing_tag Indexing_tag; /** \brief Type for the value of the filtration function. * - * Must be comparable with <. */ + * Should implement the @ref FiltrationValue concept. */ typedef typename Options::Filtration_value Filtration_value; /** \brief Key associated to each simplex. * @@ -129,8 +130,6 @@ class Simplex_tree { /** \brief Set of nodes sharing a same parent in the simplex tree. */ typedef Simplex_tree_siblings Siblings; - - struct Key_simplex_base_real { Key_simplex_base_real() : key_(-1) {} Key_simplex_base_real(Simplex_key k) : key_(k) {} @@ -149,8 +148,11 @@ class Simplex_tree { struct Extended_filtration_data { Filtration_value minval; Filtration_value maxval; - Extended_filtration_data(){} - Extended_filtration_data(Filtration_value vmin, Filtration_value vmax): minval(vmin), maxval(vmax) {} + Extended_filtration_data() {} + + Extended_filtration_data(const Filtration_value& vmin, const Filtration_value& vmax) + : minval(vmin), maxval(vmax) + {} }; typedef typename std::conditional::type Key_simplex_base; @@ -158,20 +160,33 @@ class Simplex_tree { struct Filtration_simplex_base_real { Filtration_simplex_base_real() : filt_(0) {} Filtration_simplex_base_real(Filtration_value f) : filt_(f) {} - void assign_filtration(Filtration_value f) { filt_ = f; } - Filtration_value filtration() const { return filt_; } + void assign_filtration(const Filtration_value& f) { filt_ = f; } + const Filtration_value& filtration() const { return filt_; } + Filtration_value& filtration() { return filt_; } + + static const Filtration_value& get_infinity() { return inf_; } + private: Filtration_value filt_; + + inline static const Filtration_value inf_ = std::numeric_limits::has_infinity + ? std::numeric_limits::infinity() + : std::numeric_limits::max(); }; struct Filtration_simplex_base_dummy { Filtration_simplex_base_dummy() {} Filtration_simplex_base_dummy(Filtration_value GUDHI_CHECK_code(f)) { - GUDHI_CHECK(f == 0, "filtration value specified for a complex that does not store them"); + GUDHI_CHECK(f == Filtration_value(), + "filtration value specified in the constructor for a complex that does not store them"); } - void assign_filtration(Filtration_value GUDHI_CHECK_code(f)) { - GUDHI_CHECK(f == 0, "filtration value specified for a complex that does not store them"); + + void assign_filtration(const Filtration_value& GUDHI_CHECK_code(f)) { + GUDHI_CHECK(f == Filtration_value(), "filtration value assigned for a complex that does not store them"); } - Filtration_value filtration() const { return 0; } + const Filtration_value& filtration() const { return null_; } + + private: + static constexpr const Filtration_value null_ = Filtration_value(); }; typedef typename std::conditional::type Filtration_simplex_base; @@ -637,7 +652,7 @@ class Simplex_tree { for (auto sh1 = s1->members().begin(); (sh1 != s1->members().end() && sh2 != s2->members().end()); ++sh1, ++sh2) { - if (sh1->first != sh2->first || sh1->second.filtration() != sh2->second.filtration()) + if (sh1->first != sh2->first || !(sh1->second.filtration() == sh2->second.filtration())) return false; if (has_children(sh1) != has_children(sh2)) return false; @@ -653,7 +668,7 @@ class Simplex_tree { * * Same as `filtration()`, but does not handle `null_simplex()`. */ - static Filtration_value filtration_(Simplex_handle sh) { + static const Filtration_value& filtration_(Simplex_handle sh) { GUDHI_CHECK (sh != null_simplex(), "null simplex"); return sh->second.filtration(); } @@ -681,18 +696,18 @@ class Simplex_tree { * Called on the null_simplex, it returns infinity. * If SimplexTreeOptions::store_filtration is false, returns 0. */ - static Filtration_value filtration(Simplex_handle sh) { + static const Filtration_value& filtration(Simplex_handle sh) { if (sh != null_simplex()) { return sh->second.filtration(); } else { - return std::numeric_limits::infinity(); + return Filtration_simplex_base_real::get_infinity(); } } /** \brief Sets the filtration value of a simplex. * \exception std::invalid_argument In debug mode, if sh is a null_simplex. */ - void assign_filtration(Simplex_handle sh, Filtration_value fv) { + void assign_filtration(Simplex_handle sh, const Filtration_value& fv) { GUDHI_CHECK(sh != null_simplex(), std::invalid_argument("Simplex_tree::assign_filtration - cannot assign filtration on null_simplex")); sh->second.assign_filtration(fv); @@ -905,6 +920,49 @@ class Simplex_tree { return true; } + private: + /** + * @brief Inserts a Node in the set of siblings nodes. Calls `update_simplex_tree_after_node_insertion` + * if the insertion succeeded. + * + * @tparam update_fil If true, as well as Options::store_filtration, and the node is already present, assigns + * the "union" of the input filtration_value and the value already present in the node as filtration value. + * @tparam update_children If true and the node has no children, create a child with sib as uncle. + * @tparam set_to_null If true and the node is already present, sets the returned simplex handle to null if the + * filtration value of the simplex was not modified. + * @tparam Filt Filtration value type. + * @param sib Sibling in where the node has to be inserted. + * @param v Label of the node. + * @param filtration_value Filtration value stored in the node. + * @return Pair of the iterator to the new node and a boolean indicating if the insertion succeeded. + */ + template + std::pair insert_node_(Siblings *sib, Vertex_handle v, Filt&& filtration_value) { + std::pair ins = sib->members_.try_emplace(v, sib, std::forward(filtration_value)); + + if constexpr (update_children){ + if (!(has_children(ins.first))) { + ins.first->second.assign_children(new Siblings(sib, v)); + } + } + + if (ins.second){ + // Only required when insertion is successful + update_simplex_tree_after_node_insertion(ins.first); + return ins; + } + + if constexpr (Options::store_filtration && update_fil){ + if (unify_lifetimes(ins.first->second.filtration(), filtration_value)) return ins; + } + + if constexpr (set_to_null){ + ins.first = null_simplex(); + } + + return ins; + } + protected: /** \brief Inserts a simplex represented by a range of vertex. * @param[in] simplex range of Vertex_handles, representing the vertices of the new simplex. The range must be @@ -922,42 +980,25 @@ class Simplex_tree { */ template > std::pair insert_simplex_raw(const RandomVertexHandleRange& simplex, - Filtration_value filtration) { + const Filtration_value& filtration) { Siblings * curr_sib = &root_; std::pair res_insert; auto vi = simplex.begin(); + for (; vi != std::prev(simplex.end()); ++vi) { GUDHI_CHECK(*vi != null_vertex(), "cannot use the dummy null_vertex() as a real vertex"); - res_insert = curr_sib->members_.emplace(*vi, Node(curr_sib, filtration)); - if (res_insert.second) { - // Only required when insertion is successful - update_simplex_tree_after_node_insertion(res_insert.first); - } - if (!(has_children(res_insert.first))) { - res_insert.first->second.assign_children(new Siblings(curr_sib, *vi)); - } + res_insert = insert_node_(curr_sib, *vi, filtration); curr_sib = res_insert.first->second.children(); } + GUDHI_CHECK(*vi != null_vertex(), "cannot use the dummy null_vertex() as a real vertex"); - res_insert = curr_sib->members_.emplace(*vi, Node(curr_sib, filtration)); - if (!res_insert.second) { - // if already in the complex - if (res_insert.first->second.filtration() > filtration) { - // if filtration value modified - res_insert.first->second.assign_filtration(filtration); - return res_insert; + res_insert = insert_node_(curr_sib, *vi, filtration); + if (res_insert.second){ + int dim = static_cast(boost::size(simplex)) - 1; + if (dim > dimension_) { + // Update dimension if needed + dimension_ = dim; } - // if filtration value unchanged - return std::pair(null_simplex(), false); - } else { - // Only required when insertion is successful - update_simplex_tree_after_node_insertion(res_insert.first); - } - // otherwise the insertion has succeeded - size is a size_type - int dim = static_cast(boost::size(simplex)) - 1; - if (dim > dimension_) { - // Update dimension if needed - dimension_ = dim; } return res_insert; } @@ -988,7 +1029,7 @@ class Simplex_tree { * .end() return input iterators, with 'value_type' Vertex_handle. */ template> std::pair insert_simplex(const InputVertexRange & simplex, - Filtration_value filtration = 0) { + const Filtration_value& filtration = Filtration_value()) { auto first = std::begin(simplex); auto last = std::end(simplex); @@ -1015,9 +1056,11 @@ class Simplex_tree { * output pair to the Simplex_handle of the simplex. Otherwise, we set the Simplex_handle part to * null_simplex. */ - template> - std::pair insert_simplex_and_subfaces(const InputVertexRange& Nsimplex, - Filtration_value filtration = 0) { + template > + std::pair insert_simplex_and_subfaces( + const InputVertexRange& Nsimplex, + const Filtration_value& filtration = Filtration_value()) + { auto first = std::begin(Nsimplex); auto last = std::end(Nsimplex); @@ -1042,39 +1085,24 @@ class Simplex_tree { private: // To insert {1,2,3,4}, we insert {2,3,4} twice, once at the root, and once below 1. - template + template std::pair rec_insert_simplex_and_subfaces_sorted(Siblings* sib, - ForwardVertexIterator first, - ForwardVertexIterator last, - Filtration_value filt) { + ForwardVertexIterator first, + ForwardVertexIterator last, + const Filtration_value& filt) + { // An alternative strategy would be: // - try to find the complete simplex, if found (and low filtration) exit // - insert all the vertices at once in sib // - loop over those (new or not) simplices, with a recursive call(++first, last) Vertex_handle vertex_one = *first; - auto&& dict = sib->members(); - auto insertion_result = dict.emplace(vertex_one, Node(sib, filt)); - // update extra data structures in the insertion is successful - if (insertion_result.second) { - // Only required when insertion is successful - update_simplex_tree_after_node_insertion(insertion_result.first); - } - Simplex_handle simplex_one = insertion_result.first; - bool one_is_new = insertion_result.second; - if (!one_is_new) { - if (filtration(simplex_one) > filt) { - assign_filtration(simplex_one, filt); - } else { - // FIXME: this interface makes no sense, and it doesn't seem to be tested. - insertion_result.first = null_simplex(); - } - } - if (++first == last) return insertion_result; - if (!has_children(simplex_one)) - // TODO: have special code here, we know we are building the whole subtree from scratch. - simplex_one->second.assign_children(new Siblings(sib, vertex_one)); - auto res = rec_insert_simplex_and_subfaces_sorted(simplex_one->second.children(), first, last, filt); + if (++first == last) return insert_node_(sib, vertex_one, filt); + + // TODO: have special code here, we know we are building the whole subtree from scratch. + auto insertion_result = insert_node_(sib, vertex_one, filt); + + auto res = rec_insert_simplex_and_subfaces_sorted(insertion_result.first->second.children(), first, last, filt); // No need to continue if the full simplex was already there with a low enough filtration value. if (res.first != null_simplex()) rec_insert_simplex_and_subfaces_sorted(sib, first, last, filt); return res; @@ -1121,39 +1149,123 @@ class Simplex_tree { dimension_ = dimension; } + private: + /** \brief Returns true iff the list of vertices of sh1 + * is smaller than the list of vertices of sh2 w.r.t. + * lexicographic order on the lists read in reverse. + * + * It defines a StrictWeakOrdering on simplices. The Simplex_vertex_iterators + * must traverse the Vertex_handle in decreasing order. Reverse lexicographic order satisfy + * the property that a subsimplex of a simplex is always strictly smaller with this order. */ + bool reverse_lexicographic_order(Simplex_handle sh1, Simplex_handle sh2) { + Simplex_vertex_range rg1 = simplex_vertex_range(sh1); + Simplex_vertex_range rg2 = simplex_vertex_range(sh2); + Simplex_vertex_iterator it1 = rg1.begin(); + Simplex_vertex_iterator it2 = rg2.begin(); + while (it1 != rg1.end() && it2 != rg2.end()) { + if (*it1 == *it2) { + ++it1; + ++it2; + } else { + return *it1 < *it2; + } + } + return ((it1 == rg1.end()) && (it2 != rg2.end())); + } + + /** \brief StrictWeakOrdering, for the simplices, defined by the filtration. + * + * It corresponds to the partial order + * induced by the filtration values, with ties resolved using reverse lexicographic order. + * Reverse lexicographic order has the property to always consider the subsimplex of a simplex + * to be smaller. The filtration function must be monotonic. */ + struct is_before_in_totally_ordered_filtration { + explicit is_before_in_totally_ordered_filtration(Simplex_tree * st) + : st_(st) { } + + bool operator()(const Simplex_handle sh1, const Simplex_handle sh2) const { + // Not using st_->filtration(sh1) because it uselessly tests for null_simplex. + if (!(sh1->second.filtration() == sh2->second.filtration())) { + return sh1->second.filtration() < sh2->second.filtration(); + } + // is sh1 a proper subface of sh2 + return st_->reverse_lexicographic_order(sh1, sh2); + } + + Simplex_tree * st_; + }; + public: /** \brief Initializes the filtration cache, i.e. sorts the * simplices according to their order in the filtration. + * Assumes that the filtration values have a total order. + * If not, please use @ref initialize_filtration(Comparator&&, Ignorer&&) instead. + * Two simplices with same filtration value are ordered by a reverse lexicographic order. * * It always recomputes the cache, even if one already exists. * * Any insertion, deletion or change of filtration value invalidates this cache, - * which can be cleared with clear_filtration(). */ + * which can be cleared with clear_filtration(). + */ void initialize_filtration(bool ignore_infinite_values = false) { + if (ignore_infinite_values){ + initialize_filtration(is_before_in_totally_ordered_filtration(this), [&](Simplex_handle sh) -> bool { + return filtration(sh) == Filtration_simplex_base_real::get_infinity(); + }); + } else { + initialize_filtration(is_before_in_totally_ordered_filtration(this), [](Simplex_handle) -> bool { + return false; + }); + } + } + + /** + * @brief Initializes the filtration cache, i.e. sorts the simplices according to the specified order. + * The given order has to totally order all simplices in the simplex tree such that the resulting order is a valid + * 1-parameter filtration. + * For example, if the stored filtration is a bi-filtration (or 2-parameter filtration), the given method can + * indicate, that the simplices have to be ordered by the first parameter and in case of equality by reverse + * lexicographic order. This would generate a valid 1-parameter filtration which can than be interpreted as a line + * in the 2-parameter module and used with @ref filtration_simplex_range "". Also, the persistence along this line + * can than be computed with @ref Gudhi::persistent_cohomology::Persistent_cohomology "Persistent_cohomology". + * Just note that it is important to call @ref initialize_filtration(Comparator&&, Ignorer&&) explicitly before call + * @ref Gudhi::persistent_cohomology::Persistent_cohomology::compute_persistent_cohomology + * "Persistent_cohomology::compute_persistent_cohomology" in that case, otherwise, the computation of persistence + * will fail (if no call to @ref initialize_filtration() or @ref initialize_filtration(Comparator&&, Ignorer&&) + * was made, it will call it it-self with no arguments and the simplices will be wrongly sorted). + * + * It always recomputes the cache, even if one already exists. + * + * Any insertion, deletion or change of filtration value invalidates this cache, + * which can be cleared with @ref clear_filtration(). + * + * @tparam Comparator Method type taking two Simplex_handle as input and returns a bool. + * @tparam Ignorer Method type taking one Simplex_handle as input and returns a bool. + * @param is_before_in_filtration Method used to compare two simplices with respect to their position in the + * wanted filtration. Takes two Simplex_handle as input and returns true if and only if the first simplex appears + * strictly before the second simplex in the resulting 1-parameter filtration. + * @param ignore_simplex Method taking a simplex handle as input and returns true if and only if the corresponding + * simplex should be ignored in the filtration and therefore not be cached. But note that it is important that + * the resulting filtration is still a valid filtration, that is all proper faces of a not ignored simplex + * should also not be ignored and have to appear before this simplex. There will be no tests to ensure that, so + * it is of the users responsibility. + */ + template + void initialize_filtration(Comparator&& is_before_in_filtration, Ignorer&& ignore_simplex) { filtration_vect_.clear(); filtration_vect_.reserve(num_simplices()); for (Simplex_handle sh : complex_simplex_range()) { - if (ignore_infinite_values && - std::numeric_limits::has_infinity && - filtration(sh) == std::numeric_limits::infinity()) continue; + if (ignore_simplex(sh)) continue; filtration_vect_.push_back(sh); } - /* We use stable_sort here because with libstdc++ it is faster than sort. - * is_before_in_filtration is now a total order, but we used to call - * stable_sort for the following heuristic: - * The use of a depth-first traversal of the simplex tree, provided by - * complex_simplex_range(), combined with a stable sort is meant to - * optimize the order of simplices with same filtration value. The - * heuristic consists in inserting the cofaces of a simplex as soon as - * possible. - */ #ifdef GUDHI_USE_TBB - tbb::parallel_sort(filtration_vect_.begin(), filtration_vect_.end(), is_before_in_filtration(this)); + tbb::parallel_sort(filtration_vect_.begin(), filtration_vect_.end(), is_before_in_filtration); #else - std::stable_sort(filtration_vect_.begin(), filtration_vect_.end(), is_before_in_filtration(this)); + std::stable_sort(filtration_vect_.begin(), filtration_vect_.end(), is_before_in_filtration); #endif } + /** \brief Initializes the filtration cache if it isn't initialized yet. * * Automatically called by filtration_simplex_range(). */ @@ -1276,52 +1388,6 @@ class Simplex_tree { } } - private: - /** \brief Returns true iff the list of vertices of sh1 - * is smaller than the list of vertices of sh2 w.r.t. - * lexicographic order on the lists read in reverse. - * - * It defines a StrictWeakOrdering on simplices. The Simplex_vertex_iterators - * must traverse the Vertex_handle in decreasing order. Reverse lexicographic order satisfy - * the property that a subsimplex of a simplex is always strictly smaller with this order. */ - bool reverse_lexicographic_order(Simplex_handle sh1, Simplex_handle sh2) { - Simplex_vertex_range rg1 = simplex_vertex_range(sh1); - Simplex_vertex_range rg2 = simplex_vertex_range(sh2); - Simplex_vertex_iterator it1 = rg1.begin(); - Simplex_vertex_iterator it2 = rg2.begin(); - while (it1 != rg1.end() && it2 != rg2.end()) { - if (*it1 == *it2) { - ++it1; - ++it2; - } else { - return *it1 < *it2; - } - } - return ((it1 == rg1.end()) && (it2 != rg2.end())); - } - - /** \brief StrictWeakOrdering, for the simplices, defined by the filtration. - * - * It corresponds to the partial order - * induced by the filtration values, with ties resolved using reverse lexicographic order. - * Reverse lexicographic order has the property to always consider the subsimplex of a simplex - * to be smaller. The filtration function must be monotonic. */ - struct is_before_in_filtration { - explicit is_before_in_filtration(Simplex_tree * st) - : st_(st) { } - - bool operator()(const Simplex_handle sh1, const Simplex_handle sh2) const { - // Not using st_->filtration(sh1) because it uselessly tests for null_simplex. - if (sh1->second.filtration() != sh2->second.filtration()) { - return sh1->second.filtration() < sh2->second.filtration(); - } - // is sh1 a proper subface of sh2 - return st_->reverse_lexicographic_order(sh1, sh2); - } - - Simplex_tree * st_; - }; - public: /** \brief Inserts a 1-skeleton in an empty Simplex_tree. * @@ -1396,9 +1462,7 @@ class Simplex_tree { sh->second.assign_children(new Siblings(&root_, sh->first)); } - auto insertion_res = sh->second.children()->members().emplace( - v, Node(sh->second.children(), get(edge_filtration_t(), skel_graph, edge))); - if (insertion_res.second) update_simplex_tree_after_node_insertion(insertion_res.first); + insert_node_(sh->second.children(), v, get(edge_filtration_t(), skel_graph, edge)); } } @@ -1410,7 +1474,7 @@ class Simplex_tree { * The complex does not need to be empty before calling this function. However, if a vertex is * already present, its filtration value is not modified, unlike with other insertion functions. */ template - void insert_batch_vertices(VertexRange const& vertices, Filtration_value filt = 0) { + void insert_batch_vertices(VertexRange const& vertices, const Filtration_value& filt = Filtration_value()) { auto verts = vertices | boost::adaptors::transformed([&](auto v){ return Dit_value_t(v, Node(&root_, filt)); }); root_.members_.insert(boost::begin(verts), boost::end(verts)); @@ -1484,7 +1548,7 @@ class Simplex_tree { */ void insert_edge_as_flag( Vertex_handle u , Vertex_handle v - , Filtration_value fil + , const Filtration_value& fil , int dim_max , std::vector& added_simplices) { @@ -1510,10 +1574,9 @@ class Simplex_tree { static_assert(Options::link_nodes_by_label, "Options::link_nodes_by_label must be true"); if (u == v) { // Are we inserting a vertex? - auto res_ins = root_.members().emplace(u, Node(&root_,fil)); + auto res_ins = insert_node_(&root_, u, fil); if (res_ins.second) { //if the vertex is not in the complex, insert it added_simplices.push_back(res_ins.first); //no more insert in root_.members() - update_simplex_tree_after_node_insertion(res_ins.first); if (dimension_ == -1) dimension_ = 0; } return; //because the vertex is isolated, no more insertions. @@ -1598,7 +1661,7 @@ class Simplex_tree { */ void compute_punctual_expansion( Vertex_handle v , Siblings * sib - , Filtration_value fil + , const Filtration_value& fil , int k //k == dim_max - dimension simplices in sib , std::vector& added_simplices ) { //insertion always succeeds because the edge {u,v} used to not be here. @@ -1645,10 +1708,10 @@ class Simplex_tree { * k must be > 0 */ void create_local_expansion( - Simplex_handle sh_v //Node with label v which has just been inserted - , Siblings * curr_sib //Siblings containing the node sh_v - , Filtration_value fil_uv //Fil value of the edge uv in the zz filtration - , int k //Stopping condition for recursion based on max dim + Simplex_handle sh_v //Node with label v which has just been inserted + , Siblings *curr_sib //Siblings containing the node sh_v + , const Filtration_value& fil_uv //Fil value of the edge uv in the zz filtration + , int k //Stopping condition for recursion based on max dim , std::vector &added_simplices) //range of all new simplices { //pick N^+(v) //intersect N^+(v) with labels y > v in curr_sib @@ -1672,9 +1735,9 @@ class Simplex_tree { * Only called in the case of `void insert_edge_as_flag(...)`. */ void siblings_expansion( - Siblings * siblings // must contain elements - , Filtration_value fil - , int k // == max_dim expansion - dimension curr siblings + Siblings * siblings // must contain elements + , const Filtration_value& fil + , int k // == max_dim expansion - dimension curr siblings , std::vector & added_simplices ) { if (dimension_ > k) { @@ -1717,7 +1780,7 @@ class Simplex_tree { void create_expansion(Siblings * siblings, Dictionary_it& s_h, Dictionary_it& next, - Filtration_value fil, + const Filtration_value& fil, int k, std::vector* added_simplices = nullptr) { @@ -1764,7 +1827,7 @@ class Simplex_tree { static void intersection(std::vector >& intersection, Dictionary_it begin1, Dictionary_it end1, Dictionary_it begin2, Dictionary_it end2, - Filtration_value filtration_) { + const Filtration_value& filtration_) { if (begin1 == end1 || begin2 == end2) return; // ----->> while (true) { @@ -1772,7 +1835,9 @@ class Simplex_tree { if constexpr (force_filtration_value){ intersection.emplace_back(begin1->first, Node(nullptr, filtration_)); } else { - Filtration_value filt = (std::max)({begin1->second.filtration(), begin2->second.filtration(), filtration_}); + Filtration_value filt = begin1->second.filtration(); + intersect_lifetimes(filt, begin2->second.filtration()); + intersect_lifetimes(filt, filtration_); intersection.emplace_back(begin1->first, Node(nullptr, filt)); } if (++begin1 == end1 || ++begin2 == end2) @@ -1841,7 +1906,7 @@ class Simplex_tree { to_be_inserted=false; break; } - filt = (std::max)(filt, filtration(border_child)); + intersect_lifetimes(filt, filtration(border_child)); } if (to_be_inserted) { intersection.emplace_back(next->first, Node(nullptr, filt)); @@ -1974,26 +2039,20 @@ class Simplex_tree { bool modified = false; auto fun = [&modified, this](Simplex_handle sh, int dim) -> void { if (dim == 0) return; - // Find the maximum filtration value in the border - Boundary_simplex_range&& boundary = boundary_simplex_range(sh); - Boundary_simplex_iterator max_border = std::max_element(std::begin(boundary), std::end(boundary), - [](Simplex_handle sh1, Simplex_handle sh2) { - return filtration(sh1) < filtration(sh2); - }); - - Filtration_value max_filt_border_value = filtration(*max_border); - // Replacing if(f=max)) would mean that if f is NaN, we replace it with the max of the children. - // That seems more useful than keeping NaN. - if (!(sh->second.filtration() >= max_filt_border_value)) { - // Store the filtration modification information - modified = true; - sh->second.assign_filtration(max_filt_border_value); + + Filtration_value& current_filt = sh->second.filtration(); + + // Find the maximum filtration value in the border and assigns it if it is greater than the current + for (Simplex_handle b : boundary_simplex_range(sh)) { + // considers NaN as the lowest possible value. + // stores the filtration modification information + modified |= intersect_lifetimes(current_filt, b->second.filtration()); } }; // Loop must be from the end to the beginning, as higher dimension simplex are always on the left part of the tree for_each_simplex(fun); - if(modified) + if (modified) clear_filtration(); // Drop the cache. return modified; } @@ -2009,15 +2068,19 @@ class Simplex_tree { nodes_label_to_list_.clear(); } - /** \brief Prune above filtration value given as parameter. + /** \brief Prune above filtration value given as parameter. That is: if \f$ f \f$ is the given filtration value + * and \f$ f_s \f$ the filtration value associated to the simplex \f$ s \f$ contained in the tree, \f$ s \f$ is + * removed if and only if \f$ f < f_s \f$ or \f$ f_s = NaN \f$, where \f$ < \f$ is the + * @ref FiltrationValue::operator< "operator<" defined for the type + * @ref SimplexTreeOptions::Filtration_value "Simplex_tree::Filtration_value". * @param[in] filtration Maximum threshold value. * @return True if any simplex was removed, false if all simplices already had a value below the threshold. * \post Note that the dimension of the simplicial complex may be lower after calling `prune_above_filtration()` * than it was before. However, `upper_bound_dimension()` will return the old value, which remains a valid upper * bound. If you care, you can call `dimension()` to recompute the exact dimension. */ - bool prune_above_filtration(Filtration_value filtration) { - if (std::numeric_limits::has_infinity && filtration == std::numeric_limits::infinity()) + bool prune_above_filtration(const Filtration_value& filtration) { + if (filtration == Filtration_simplex_base_real::get_infinity()) return false; // ---->> bool modified = rec_prune_above_filtration(root(), filtration); if(modified) @@ -2026,18 +2089,22 @@ class Simplex_tree { } private: - bool rec_prune_above_filtration(Siblings* sib, Filtration_value filt) { + bool rec_prune_above_filtration(Siblings* sib, const Filtration_value& filt) { auto&& list = sib->members(); bool modified = false; bool emptied = false; Simplex_handle last; auto to_remove = [this, filt](Dit_value_t& simplex) { - if (simplex.second.filtration() <= filt) return false; - if (has_children(&simplex)) rec_delete(simplex.second.children()); - // dimension may need to be lowered - dimension_to_be_lowered_ = true; - return true; + //if filt and simplex.second.filtration() are non comparable, should return false. + //if simplex.second.filtration() is NaN, should return true. + if (filt < simplex.second.filtration() || !(simplex.second.filtration() == simplex.second.filtration())) { + if (has_children(&simplex)) rec_delete(simplex.second.children()); + // dimension may need to be lowered + dimension_to_be_lowered_ = true; + return true; + } + return false; }; //TODO: `if constexpr` replaceable by `std::erase_if` in C++20? Has a risk of additional runtime, @@ -2201,10 +2268,12 @@ class Simplex_tree { * @param[in] efd Structure containing the minimum and maximum values of the original filtration. This the output of `extend_filtration()`. * @return A pair containing the original filtration value of the simplex as well as the simplex type. */ - std::pair decode_extended_filtration(Filtration_value f, const Extended_filtration_data& efd){ + std::pair decode_extended_filtration(const Filtration_value& f, + const Extended_filtration_data& efd) + { std::pair p; - Filtration_value minval = efd.minval; - Filtration_value maxval = efd.maxval; + const Filtration_value& minval = efd.minval; + const Filtration_value& maxval = efd.maxval; if (f >= -2 && f <= -1) { p.first = minval + (maxval-minval)*(f + 2); p.second = Extended_simplex_type::UP; } else if (f >= 1 && f <= 2) { @@ -2215,6 +2284,7 @@ class Simplex_tree { return p; }; + //TODO: externalize this method and `decode_extended_filtration` /** \brief Extend filtration for computing extended persistence. * This function only uses the filtration values at the 0-dimensional simplices, * and computes the extended persistence diagram induced by the lower-star filtration @@ -2222,6 +2292,10 @@ class Simplex_tree { * \post Note that after calling this function, the filtration * values are actually modified. The function `decode_extended_filtration()` * retrieves the original values and outputs the extended simplex type. + * + * @warning Currently only works for @ref SimplexTreeOptions::Filtration_value which are + * float types like `float` or `double`. + * * @exception std::invalid_argument In debug mode if the Simplex tree contains a vertex with the largest * Vertex_handle, as this method requires to create an extra vertex internally. * @return A data structure containing the maximum and minimum values of the original filtration. @@ -2233,10 +2307,10 @@ class Simplex_tree { // Compute maximum and minimum of filtration values Vertex_handle maxvert = std::numeric_limits::min(); - Filtration_value minval = std::numeric_limits::infinity(); - Filtration_value maxval = -std::numeric_limits::infinity(); + Filtration_value minval = Filtration_simplex_base_real::get_infinity(); + Filtration_value maxval = -Filtration_simplex_base_real::get_infinity(); for (auto sh = root_.members().begin(); sh != root_.members().end(); ++sh) { - Filtration_value f = this->filtration(sh); + const Filtration_value& f = this->filtration(sh); minval = std::min(minval, f); maxval = std::max(maxval, f); maxvert = std::max(sh->first, maxvert); @@ -2264,7 +2338,7 @@ class Simplex_tree { // Create cone on simplex vr.push_back(maxvert); if (this->dimension(sh) == 0) { - Filtration_value v = this->filtration(sh); + const Filtration_value& v = this->filtration(sh); Filtration_value scaled_v = (v - minval) * scale; // Assign ascending value between -2 and -1 to vertex this->assign_filtration(sh, -2 + scaled_v); @@ -2458,7 +2532,7 @@ class Simplex_tree { * @param[in] filt_value The new filtration value. * @param[in] min_dim The minimal dimension. Default value is 0. */ - void reset_filtration(Filtration_value filt_value, int min_dim = 0) { + void reset_filtration(const Filtration_value& filt_value, int min_dim = 0) { rec_reset_filtration(&root_, filt_value, min_dim); clear_filtration(); // Drop the cache. } @@ -2469,7 +2543,7 @@ class Simplex_tree { * @param[in] filt_value The new filtration value. * @param[in] min_depth The minimal depth. */ - void rec_reset_filtration(Siblings * sib, Filtration_value filt_value, int min_depth) { + void rec_reset_filtration(Siblings * sib, const Filtration_value& filt_value, int min_depth) { for (auto sh = sib->members().begin(); sh != sib->members().end(); ++sh) { if (min_depth <= 0) { sh->second.assign_filtration(filt_value); @@ -2480,6 +2554,22 @@ class Simplex_tree { } } + std::size_t num_simplices_and_filtration_size(Siblings* sib, std::size_t& fv_byte_size) const { + using namespace Gudhi::simplex_tree; + + auto sib_begin = sib->members().begin(); + auto sib_end = sib->members().end(); + size_t simplices_number = sib->members().size(); + for (auto sh = sib_begin; sh != sib_end; ++sh) { + if constexpr (SimplexTreeOptions::store_filtration) + fv_byte_size += get_serialization_size_of(sh->second.filtration()); + if (has_children(sh)) { + simplices_number += num_simplices_and_filtration_size(sh->second.children(), fv_byte_size); + } + } + return simplices_number; + } + public: /** @private @brief Returns the serialization required buffer size. * @@ -2489,9 +2579,12 @@ class Simplex_tree { * architecture. */ std::size_t get_serialization_size() { + using namespace Gudhi::simplex_tree; + const std::size_t vh_byte_size = sizeof(Vertex_handle); - const std::size_t fv_byte_size = SimplexTreeOptions::store_filtration ? sizeof(Filtration_value) : 0; - const std::size_t buffer_byte_size = vh_byte_size + num_simplices() * (fv_byte_size + 2 * vh_byte_size); + std::size_t fv_byte_size = 0; + const std::size_t tree_size = num_simplices_and_filtration_size(&root_, fv_byte_size); + const std::size_t buffer_byte_size = vh_byte_size + fv_byte_size + tree_size * 2 * vh_byte_size; #ifdef DEBUG_TRACES std::clog << "Gudhi::simplex_tree::get_serialization_size - buffer size = " << buffer_byte_size << std::endl; #endif // DEBUG_TRACES @@ -2536,15 +2629,16 @@ class Simplex_tree { private: /** \brief Serialize each element of the sibling and recursively call serialization. */ char* rec_serialize(Siblings *sib, char* buffer) { + using namespace Gudhi::simplex_tree; char* ptr = buffer; - ptr = Gudhi::simplex_tree::serialize_trivial(static_cast(sib->members().size()), ptr); + ptr = serialize_trivial(static_cast(sib->members().size()), ptr); #ifdef DEBUG_TRACES std::clog << "\n" << sib->members().size() << " : "; #endif // DEBUG_TRACES for (auto& map_el : sib->members()) { - ptr = Gudhi::simplex_tree::serialize_trivial(map_el.first, ptr); // Vertex + ptr = serialize_trivial(map_el.first, ptr); // Vertex if (Options::store_filtration) - ptr = Gudhi::simplex_tree::serialize_trivial(map_el.second.filtration(), ptr); // Filtration + ptr = serialize_trivial(map_el.second.filtration(), ptr); // Filtration #ifdef DEBUG_TRACES std::clog << " [ " << map_el.first << " | " << map_el.second.filtration() << " ] "; #endif // DEBUG_TRACES @@ -2553,7 +2647,7 @@ class Simplex_tree { if (has_children(&map_el)) { ptr = rec_serialize(map_el.second.children(), ptr); } else { - ptr = Gudhi::simplex_tree::serialize_trivial(static_cast(0), ptr); + ptr = serialize_trivial(static_cast(0), ptr); #ifdef DEBUG_TRACES std::cout << "\n0 : "; #endif // DEBUG_TRACES @@ -2577,12 +2671,40 @@ class Simplex_tree { * */ void deserialize(const char* buffer, const std::size_t buffer_size) { + deserialize(buffer, buffer_size, [](Filtration_value& filtration, const char* ptr) { + using namespace Gudhi::simplex_tree; + return deserialize_trivial(filtration, ptr); + }); + } + + /** @private @brief Deserialize the array of char (flatten version of the tree) to initialize a Simplex tree. + * It is the user's responsibility to provide an 'empty' Simplex_tree, there is no guarantee otherwise. + * + * @tparam F Method taking an @ref Filtration_value and a `const char*` as input and returning a + * `const char*`. + * @param[in] buffer A pointer on a buffer that contains a serialized Simplex_tree. + * @param[in] buffer_size The size of the buffer. + * @param[in] deserialize_filtration_value To provide if the type of @ref Filtration_value is not trivially + * convertible from the filtration value type of the serialized simplex tree. Takes the filtration value to fill and + * a pointer to the current position in the buffer as arguments and returns the new position of the pointer after + * reading the filtration value and transforming it into a element of the host's filtration value type. + * + * @exception std::invalid_argument In case the deserialization does not finish at the correct buffer_size. + * @exception std::logic_error In debug mode, if the Simplex_tree is not 'empty'. + * + * @warning Serialize/Deserialize is not portable. It is meant to be read in a Simplex_tree with the same + * SimplexTreeOptions (except for the @ref Filtration_value type) and on a computer with the same architecture. + * + */ + template + void deserialize(const char* buffer, const std::size_t buffer_size, F&& deserialize_filtration_value) { + using namespace Gudhi::simplex_tree; GUDHI_CHECK(num_vertices() == 0, std::logic_error("Simplex_tree::deserialize - Simplex_tree must be empty")); const char* ptr = buffer; // Needs to read size before recursivity to manage new siblings for children Vertex_handle members_size; - ptr = Gudhi::simplex_tree::deserialize_trivial(members_size, ptr); - ptr = rec_deserialize(&root_, members_size, ptr, 0); + ptr = deserialize_trivial(members_size, ptr); + ptr = rec_deserialize(&root_, members_size, ptr, 0, deserialize_filtration_value); if (static_cast(ptr - buffer) != buffer_size) { throw std::invalid_argument("Deserialization does not match end of buffer"); } @@ -2590,31 +2712,38 @@ class Simplex_tree { private: /** \brief Serialize each element of the sibling and recursively call serialization. */ - const char* rec_deserialize(Siblings *sib, Vertex_handle members_size, const char* ptr, int dim) { + template + const char* rec_deserialize(Siblings* sib, + Vertex_handle members_size, + const char* ptr, + int dim, + [[maybe_unused]] F&& deserialize_filtration_value) + { + using namespace Gudhi::simplex_tree; // In case buffer is just a 0 char if (members_size > 0) { if constexpr (!Options::stable_simplex_handles) sib->members_.reserve(members_size); Vertex_handle vertex; - Filtration_value filtration; + // Set explicitly to zero to avoid false positive error raising in debug mode when store_filtration == false + // and to force array like Filtration_value value to be empty. + Filtration_value filtration(0); for (Vertex_handle idx = 0; idx < members_size; idx++) { - ptr = Gudhi::simplex_tree::deserialize_trivial(vertex, ptr); - if (Options::store_filtration) { - ptr = Gudhi::simplex_tree::deserialize_trivial(filtration, ptr); - // Default is no children - sib->members_.emplace_hint(sib->members_.end(), vertex, Node(sib, filtration)); - } else { - // Default is no children - sib->members_.emplace_hint(sib->members_.end(), vertex, Node(sib)); + ptr = deserialize_trivial(vertex, ptr); + if constexpr (Options::store_filtration) { + ptr = deserialize_filtration_value(filtration, ptr); } + // Default is no children + // If store_filtration is false, `filtration` is ignored. + sib->members_.emplace_hint(sib->members_.end(), vertex, Node(sib, filtration)); } Vertex_handle child_size; for (auto sh = sib->members().begin(); sh != sib->members().end(); ++sh) { update_simplex_tree_after_node_insertion(sh); - ptr = Gudhi::simplex_tree::deserialize_trivial(child_size, ptr); + ptr = deserialize_trivial(child_size, ptr); if (child_size > 0) { Siblings* child = new Siblings(sib, sh->first); sh->second.assign_children(child); - ptr = rec_deserialize(child, child_size, ptr, dim + 1); + ptr = rec_deserialize(child, child_size, ptr, dim + 1, deserialize_filtration_value); } } if (dim > dimension_) { diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_node_explicit_storage.h b/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_node_explicit_storage.h index 738f13b817..6c5eb0f15e 100644 --- a/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_node_explicit_storage.h +++ b/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_node_explicit_storage.h @@ -43,7 +43,13 @@ struct GUDHI_EMPTY_BASE_CLASS_OPTIMIZATION Simplex_tree_node_explicit_storage typedef typename SimplexTree::Simplex_key Simplex_key; typedef typename SimplexTree::Simplex_data Simplex_data; - Simplex_tree_node_explicit_storage(Siblings* sib = nullptr, Filtration_value filtration = 0) + // Simplex_tree_node_explicit_storage() : children_(nullptr) {} + + // Simplex_tree_node_explicit_storage(Siblings* sib, const Filtration_value& filtration) : children_(sib) + // { + // this->assign_filtration(filtration); + // } + Simplex_tree_node_explicit_storage(Siblings* sib = nullptr, const Filtration_value& filtration = Filtration_value()) : SimplexTree::Filtration_simplex_base(filtration), children_(sib) {} diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_siblings.h b/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_siblings.h index d849eeba07..c0334ff014 100644 --- a/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_siblings.h +++ b/src/Simplex_tree/include/gudhi/Simplex_tree/Simplex_tree_siblings.h @@ -11,13 +11,11 @@ #ifndef SIMPLEX_TREE_SIMPLEX_TREE_SIBLINGS_H_ #define SIMPLEX_TREE_SIMPLEX_TREE_SIBLINGS_H_ +#include #include #include -#include -#include - namespace Gudhi { /** \addtogroup simplex_tree @@ -72,18 +70,6 @@ class Simplex_tree_siblings { } } - /** \brief Inserts a Node in the set of siblings nodes. - * - * If already present, assigns the minimal filtration value - * between input filtration_value and the value already - * present in the node. - */ - void insert(Vertex_handle v, Filtration_value filtration_value) { - auto ins = members_.emplace(v, Node(this, filtration_value)); - if (!ins.second && filtration(ins.first) > filtration_value) - ins.first->second.assign_filtration(filtration_value); - } - Dictionary_it find(Vertex_handle v) { return members_.find(v); } diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree/serialization_utils.h b/src/Simplex_tree/include/gudhi/Simplex_tree/serialization_utils.h index c0bb045c45..1892be9a44 100644 --- a/src/Simplex_tree/include/gudhi/Simplex_tree/serialization_utils.h +++ b/src/Simplex_tree/include/gudhi/Simplex_tree/serialization_utils.h @@ -12,13 +12,14 @@ #define SIMPLEX_TREE_SERIALIZATION_UTILS_H_ #include // for memcpy and std::size_t -#include namespace Gudhi { namespace simplex_tree { -/** \brief Serialize the given value and insert it at start position. +/** + * @ingroup simplex_tree + * @brief Serialize the given value and insert it at start position. * * @param[in] value The value to serialize. * @param[in] start Start position where the value is serialized. @@ -33,7 +34,9 @@ char* serialize_trivial(ArgumentType value, char* start) { return start + arg_size; } -/** \brief Deserialize at the start position in an array of char and sets the value with it. +/** + * @ingroup simplex_tree + * \brief Deserialize at the start position in an array of char and sets the value with it. * * @param[in] value The value where to deserialize based on its type. * @param[in] start Start position where the value is serialized. @@ -48,6 +51,15 @@ const char* deserialize_trivial(ArgumentType& value, const char* start) { return (start + arg_size); } +/** + * @ingroup simplex_tree + * @brief Returns the size of the serialization of the given object. + */ +template +constexpr std::size_t get_serialization_size_of([[maybe_unused]] ArgumentType value) { + return sizeof(ArgumentType); +} + } // namespace simplex_tree } // namespace Gudhi diff --git a/src/Simplex_tree/include/gudhi/Simplex_tree/simplex_tree_options.h b/src/Simplex_tree/include/gudhi/Simplex_tree/simplex_tree_options.h index d8db37d244..ba8ae5eb68 100644 --- a/src/Simplex_tree/include/gudhi/Simplex_tree/simplex_tree_options.h +++ b/src/Simplex_tree/include/gudhi/Simplex_tree/simplex_tree_options.h @@ -14,7 +14,9 @@ #include #include -#include // void_t +#include +#include // void_t +#include // std::isnan namespace Gudhi { @@ -99,6 +101,54 @@ template struct Get_simplex_data_type { typedef No_simplex_d template struct Get_simplex_data_type> { typedef typename O::Simplex_data type; }; +/** + * @brief Given two filtration values at which a simplex exists, stores in the first value the minimal union of births + * generating a lifetime including those two values. + * This is the overload for when @ref FiltrationValue is an arithmetic type, like double, int etc. + * Because the filtration values are totally ordered then, the union is simply the minimum of the two values. + * + * NaN values are not supported. + */ +template +bool unify_lifetimes(Arithmetic_filtration_value& f1, Arithmetic_filtration_value f2) +{ + if (f2 < f1){ + f1 = f2; + return true; + } + return false; +} + +/** + * @brief Given two filtration values, stores in the first value the lowest common upper bound of the two values. + * If a filtration value has value `NaN`, it should be considered as the lowest value possible. + * This is the overload for when @ref FiltrationValue is an arithmetic type, like double, float, int etc. + * Because the filtration values are totally ordered then, the upper bound is always the maximum of the two values. + */ +template +bool intersect_lifetimes(Arithmetic_filtration_value& f1, Arithmetic_filtration_value f2) +{ + if constexpr (std::numeric_limits::has_quiet_NaN) { + if (std::isnan(f1)) { + f1 = f2; + return !std::isnan(f2); + } + + // Computes the max while handling NaN as lowest value. + if (!(f1 < f2)) return false; + + f1 = f2; + return true; + } else { + // NaN not possible. + if (f1 < f2){ + f1 = f2; + return true; + } + return false; + } +} + } // namespace Gudhi #endif // SIMPLEX_TREE_SIMPLEX_TREE_OPTIONS_H_ diff --git a/src/Simplex_tree/test/simplex_tree_ctor_and_move_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_ctor_and_move_unit_test.cpp index 1d370adfef..8f2f979b6f 100644 --- a/src/Simplex_tree/test/simplex_tree_ctor_and_move_unit_test.cpp +++ b/src/Simplex_tree/test/simplex_tree_ctor_and_move_unit_test.cpp @@ -24,6 +24,8 @@ // /!\ Nothing else from Simplex_tree shall be included to test includes are well defined. #include "gudhi/Simplex_tree.h" +#include "test_vector_filtration_simplex_tree.h" + using namespace Gudhi; typedef boost::mpl::list, @@ -177,42 +179,6 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_copy_constructor, Simplex_tree, list_of_te } -struct Simplex_tree_options_custom_fil_values_default { - typedef linear_indexing_tag Indexing_tag; - typedef std::int16_t Vertex_handle; - typedef std::vector Filtration_value; - typedef std::int32_t Simplex_key; - static const bool store_key = false; - static const bool store_filtration = true; - static const bool contiguous_vertices = false; - static const bool link_nodes_by_label = false; - static const bool stable_simplex_handles = false; -}; - -struct Simplex_tree_options_custom_fil_values_fast_persistence { - typedef linear_indexing_tag Indexing_tag; - typedef std::int16_t Vertex_handle; - typedef std::vector Filtration_value; - typedef std::int32_t Simplex_key; - static const bool store_key = true; - static const bool store_filtration = true; - static const bool contiguous_vertices = true; - static const bool link_nodes_by_label = false; - static const bool stable_simplex_handles = false; -}; - -struct Simplex_tree_options_custom_fil_values_full_featured { - typedef linear_indexing_tag Indexing_tag; - typedef std::int16_t Vertex_handle; - typedef std::vector Filtration_value; - typedef std::int32_t Simplex_key; - static const bool store_key = true; - static const bool store_filtration = true; - static const bool contiguous_vertices = false; - static const bool link_nodes_by_label = true; - static const bool stable_simplex_handles = true; -}; - typedef boost::mpl::list, Simplex_tree, Simplex_tree > diff --git a/src/Simplex_tree/test/simplex_tree_extended_filtration_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_extended_filtration_unit_test.cpp index 39ec3e5ee3..6eb8ab5d84 100644 --- a/src/Simplex_tree/test/simplex_tree_extended_filtration_unit_test.cpp +++ b/src/Simplex_tree/test/simplex_tree_extended_filtration_unit_test.cpp @@ -10,6 +10,7 @@ #include #include // for std::uint8_t +#include // std::isnan #define BOOST_TEST_DYN_LINK #define BOOST_TEST_MODULE "simplex_tree_extended_filtration" diff --git a/src/Simplex_tree/test/simplex_tree_serialization_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_serialization_unit_test.cpp index a8450cfa72..c2ea8e24ed 100644 --- a/src/Simplex_tree/test/simplex_tree_serialization_unit_test.cpp +++ b/src/Simplex_tree/test/simplex_tree_serialization_unit_test.cpp @@ -11,8 +11,7 @@ #include #include // for std::size_t and strncmp #include -#include // for std::distance -#include +#include #include // for std::uint8_t #include // for std::setfill, setw #include // for std::hex, uppercase @@ -26,6 +25,8 @@ #include // for de/serialize_trivial #include // for GUDHI_TEST_FLOAT_EQUALITY_CHECK +#include "test_vector_filtration_simplex_tree.h" + using namespace Gudhi; using namespace Gudhi::simplex_tree; @@ -49,17 +50,57 @@ typedef boost::mpl::list, Simplex_tree, Simplex_tree, Simplex_tree, - Simplex_tree > list_of_tested_variants; - -template -Filtration_type random_filtration(Filtration_type lower_bound = 0, Filtration_type upper_bound = 1) { + Simplex_tree, + Simplex_tree, + Simplex_tree, + Simplex_tree > list_of_tested_variants; + +template +Filtration_type random_filtration_ar(Filtration_type lower_bound = 0, + Filtration_type upper_bound = 1) +{ std::uniform_real_distribution unif(lower_bound, upper_bound); std::random_device rand_dev; std::mt19937 rand_engine(rand_dev()); - + return unif(rand_engine); } +template +Filtration_type random_filtration_vec(typename Filtration_type::value_type lower_bound = 0, + typename Filtration_type::value_type upper_bound = 10, + unsigned int number_of_parameters = 2) +{ + std::uniform_int_distribution unif(lower_bound, upper_bound); + std::random_device rand_dev; + std::mt19937 rand_engine(rand_dev()); + + Filtration_type res(number_of_parameters); + for (unsigned int i = 0; i < number_of_parameters; ++i) res[i] = unif(rand_engine); + + return res; +} + +template +Filtration_type random_filtration() +{ + if constexpr (std::is_arithmetic_v) { + return random_filtration_ar(); + } else { + return random_filtration_vec(); + } +} + +template +void test_equality(const Filtration_type& filt1, const Filtration_type& filt2) +{ + if constexpr (std::is_arithmetic_v) { + GUDHI_TEST_FLOAT_EQUALITY_CHECK(filt1, filt2); + } else { + BOOST_CHECK(filt1 == filt2); + } +} + BOOST_AUTO_TEST_CASE_TEMPLATE(basic_simplex_tree_serialization, Stree, list_of_tested_variants) { std::clog << "********************************************************************" << std::endl; std::clog << "BASIC SIMPLEX TREE SERIALIZATION/DESERIALIZATION" << std::endl; @@ -111,12 +152,13 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(basic_simplex_tree_serialization, Stree, list_of_t std::clog << "Serialization size in bytes = " << buffer_size << std::endl; // Sizes are expressed in bytes const std::size_t vertex_size = sizeof(Vertex_type); - const std::size_t filtration_size = Stree::Options::store_filtration ? sizeof(Filtration_type) : 0; + const std::size_t filtration_size = + Stree::Options::store_filtration ? get_serialization_size_of(random_filtration()) : 0; const std::size_t serialization_size = vertex_size + st.num_simplices() * (2 * vertex_size + filtration_size); - BOOST_CHECK(serialization_size == buffer_size); + BOOST_CHECK_EQUAL(serialization_size, buffer_size); Vertex_type vertex = 0; - Filtration_type filtration = 0; + Filtration_type filtration(0); // Reset position pointer at start const char* c_ptr = buffer; // 3 simplices ({0}, {1}, {2}) and its filtration values @@ -126,19 +168,19 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(basic_simplex_tree_serialization, Stree, list_of_t BOOST_CHECK(vertex == 0); if (Stree::Options::store_filtration) { c_ptr = deserialize_trivial(filtration, c_ptr); - GUDHI_TEST_FLOAT_EQUALITY_CHECK(filtration, st.filtration(st.find({0}))); + test_equality(filtration, st.filtration(st.find({0}))); } c_ptr = deserialize_trivial(vertex, c_ptr); BOOST_CHECK(vertex == 1); if (Stree::Options::store_filtration) { c_ptr = deserialize_trivial(filtration, c_ptr); - GUDHI_TEST_FLOAT_EQUALITY_CHECK(filtration, st.filtration(st.find({1}))); + test_equality(filtration, st.filtration(st.find({1}))); } c_ptr = deserialize_trivial(vertex, c_ptr); BOOST_CHECK(vertex == 2); if (Stree::Options::store_filtration) { c_ptr = deserialize_trivial(filtration, c_ptr); - GUDHI_TEST_FLOAT_EQUALITY_CHECK(filtration, st.filtration(st.find({2}))); + test_equality(filtration, st.filtration(st.find({2}))); } // 1 simplex (2) from {0, 2} and its filtration values c_ptr = deserialize_trivial(vertex, c_ptr); @@ -147,7 +189,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(basic_simplex_tree_serialization, Stree, list_of_t BOOST_CHECK(vertex == 2); if (Stree::Options::store_filtration) { c_ptr = deserialize_trivial(filtration, c_ptr); - GUDHI_TEST_FLOAT_EQUALITY_CHECK(filtration, st.filtration(st.find({0, 2}))); + test_equality(filtration, st.filtration(st.find({0, 2}))); } c_ptr = deserialize_trivial(vertex, c_ptr); // (0, 2) end of leaf BOOST_CHECK(vertex == 0); @@ -320,3 +362,84 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(simplex_tree_serialization_and_cofaces, Stree, lis BOOST_CHECK(num_stars == 5); } + +struct Simplex_tree_vec_fil : Simplex_tree_options_default { + typedef Vector_filtration_value Filtration_value; +}; + +BOOST_AUTO_TEST_CASE(simplex_tree_custom_deserialization_vec_to_double) { + using source_st_type = Simplex_tree; + using target_st_type = Simplex_tree; + + source_st_type st; + target_st_type st_copy; + target_st_type st_witness; + + st.insert_simplex_and_subfaces({2, 1, 0}, {3, 1}); + st.insert_simplex_and_subfaces({0, 1, 6, 7}, {4, 1}); + st.insert_simplex_and_subfaces({3, 0}, {2, 1}); + st.insert_simplex_and_subfaces({3, 4, 5}, {3, 1}); + st.insert_simplex_and_subfaces({8}, {1, 1}); + + st_witness.insert_simplex_and_subfaces({2, 1, 0}, 3.0); + st_witness.insert_simplex_and_subfaces({0, 1, 6, 7}, 4.0); + st_witness.insert_simplex_and_subfaces({3, 0}, 2.0); + st_witness.insert_simplex_and_subfaces({3, 4, 5}, 3.0); + st_witness.insert_simplex_and_subfaces({8}, 1.0); + + auto deserialize_filtration_value = [](typename target_st_type::Filtration_value& fil, + const char* start) -> const char* { + typename source_st_type::Filtration_value origin_fil; + const char* ptr = deserialize_trivial(origin_fil, start); //naive way to do it, but fine enough for a test + fil = origin_fil[0]; + return ptr; + }; + + const std::size_t stree_buffer_size = st.get_serialization_size(); + char* stree_buffer = new char[stree_buffer_size]; + st.serialize(stree_buffer, stree_buffer_size); + + st_copy.deserialize(stree_buffer, stree_buffer_size, deserialize_filtration_value); + delete[] stree_buffer; + + BOOST_CHECK(st_witness == st_copy); +} + +BOOST_AUTO_TEST_CASE(simplex_tree_custom_deserialization_double_to_vec) { + using source_st_type = Simplex_tree; + using target_st_type = Simplex_tree; + + source_st_type st; + target_st_type st_copy; + target_st_type st_witness; + + st_witness.insert_simplex_and_subfaces({2, 1, 0}, {3, 1}); + st_witness.insert_simplex_and_subfaces({0, 1, 6, 7}, {4, 1}); + st_witness.insert_simplex_and_subfaces({3, 0}, {2, 1}); + st_witness.insert_simplex_and_subfaces({3, 4, 5}, {3, 1}); + st_witness.insert_simplex_and_subfaces({8}, {1, 1}); + + st.insert_simplex_and_subfaces({2, 1, 0}, 3.0); + st.insert_simplex_and_subfaces({0, 1, 6, 7}, 4.0); + st.insert_simplex_and_subfaces({3, 0}, 2.0); + st.insert_simplex_and_subfaces({3, 4, 5}, 3.0); + st.insert_simplex_and_subfaces({8}, 1.0); + + auto deserialize_filtration_value = [](typename target_st_type::Filtration_value& fil, + const char* start) -> const char* { + typename source_st_type::Filtration_value origin_fil; + const char* ptr = deserialize_trivial(origin_fil, start); + fil.resize(2, 1); + fil[0] = origin_fil; + return ptr; + }; + + const std::size_t stree_buffer_size = st.get_serialization_size(); + char* stree_buffer = new char[stree_buffer_size]; + st.serialize(stree_buffer, stree_buffer_size); + + st_copy.deserialize(stree_buffer, stree_buffer_size, deserialize_filtration_value); + delete[] stree_buffer; + + BOOST_CHECK(st_witness == st_copy); +} diff --git a/src/Simplex_tree/test/simplex_tree_unit_test.cpp b/src/Simplex_tree/test/simplex_tree_unit_test.cpp index d0dc0140ec..e4a168591b 100644 --- a/src/Simplex_tree/test/simplex_tree_unit_test.cpp +++ b/src/Simplex_tree/test/simplex_tree_unit_test.cpp @@ -29,6 +29,8 @@ // /!\ Nothing else from Simplex_tree shall be included to test includes are well defined. #include "gudhi/Simplex_tree.h" +#include "test_vector_filtration_simplex_tree.h" + using namespace Gudhi; typedef boost::mpl::list, @@ -1242,3 +1244,214 @@ BOOST_AUTO_TEST_CASE(simplex_data) { st.simplex_data(st.find({0, 1, 2})) = 4; BOOST_CHECK(st.simplex_data(st.find({0, 1})) == 5); } + +typedef boost::mpl::list, + Simplex_tree, + Simplex_tree > + list_of_custom_fil_variants; + +BOOST_AUTO_TEST_CASE_TEMPLATE(initialize_filtration_std, typeST, list_of_tested_variants) { + typeST st; + using Vertex_handle = typename typeST::Vertex_handle; + using Filtration_value = typename typeST::Filtration_value; + + auto inf = typeST::Filtration_simplex_base::get_infinity(); + + st.insert_simplex_and_subfaces({2, 1, 0}, 3.); + st.insert_simplex_and_subfaces({2, 3, 4}, 2.); + st.insert_simplex_and_subfaces({2, 3, 1}, inf); + + std::vector > simplices{ + {2}, + {3}, + {3, 2}, + {4}, + {4, 2}, + {4, 3}, + {4, 3, 2}, + {0}, + {1}, + {1, 0}, + {2, 0}, + {2, 1}, + {2, 1, 0}, + {3, 1}, + {3, 2, 1} + }; + std::vector fil{ + 2., + 2., + 2., + 2., + 2., + 2., + 2., + 3., + 3., + 3., + 3., + 3., + 3., + inf, + inf + }; + + st.initialize_filtration(); + unsigned int i = 0; + for (auto sh : st.filtration_simplex_range()){ + unsigned int j = 0; + for (auto v : st.simplex_vertex_range(sh)){ + BOOST_CHECK_EQUAL(v, simplices[i][j]); + ++j; + } + BOOST_CHECK_EQUAL(st.filtration(sh), fil[i]); + ++i; + } + BOOST_CHECK_EQUAL(i, 15); + + st.initialize_filtration(true); + i = 0; + for (auto sh : st.filtration_simplex_range()){ + unsigned int j = 0; + for (auto v : st.simplex_vertex_range(sh)){ + BOOST_CHECK_EQUAL(v, simplices[i][j]); + ++j; + } + BOOST_CHECK_EQUAL(st.filtration(sh), fil[i]); + ++i; + } + BOOST_CHECK_EQUAL(i, 13); +} + +BOOST_AUTO_TEST_CASE_TEMPLATE(initialize_filtration_vector, typeST, list_of_custom_fil_variants) { + typeST st; + using Vertex_handle = typename typeST::Vertex_handle; + using Filtration_value = typename typeST::Filtration_value; + using Simplex_handle = typename typeST::Simplex_handle; + + auto inf = typeST::Filtration_simplex_base::get_infinity(); + auto max = std::numeric_limits::max(); + + st.insert_simplex_and_subfaces({2, 1, 0}, {3, 1}); + st.insert_simplex_and_subfaces({2, 3, 4}, {2, max}); + st.insert_simplex_and_subfaces({2, 3, 1}, inf); + + std::vector > simplices1{ + {2}, + {3}, + {3, 2}, + {4}, + {4, 2}, + {4, 3}, + {4, 3, 2}, + {0}, + {1}, + {1, 0}, + {2, 0}, + {2, 1}, + {2, 1, 0}, + {3, 1}, + {3, 2, 1} + }; + std::vector fil1{ + {2, 1}, + {2, max}, + {2, max}, + {2, max}, + {2, max}, + {2, max}, + {2, max}, + {3, 1}, + {3, 1}, + {3, 1}, + {3, 1}, + {3, 1}, + {3, 1}, + inf, + inf + }; + + std::vector > simplices2{ + {0}, + {1}, + {1, 0}, + {2}, + {2, 0}, + {2, 1}, + {2, 1, 0} + }; + std::vector fil2{ + {3, 1}, + {3, 1}, + {3, 1}, + {2, 1}, + {3, 1}, + {3, 1}, + {3, 1} + }; + + auto lex_order = [&st](Simplex_handle sh1, Simplex_handle sh2) -> bool { + auto rg1 = st.simplex_vertex_range(sh1); + auto rg2 = st.simplex_vertex_range(sh2); + auto it1 = rg1.begin(); + auto it2 = rg2.begin(); + while (it1 != rg1.end() && it2 != rg2.end()) { + if (*it1 == *it2) { + ++it1; + ++it2; + } else { + return *it1 < *it2; + } + } + return ((it1 == rg1.end()) && (it2 != rg2.end())); + }; + + st.initialize_filtration( + [&](Simplex_handle sh1, Simplex_handle sh2) { + auto f1 = st.filtration(sh1); + auto f2 = st.filtration(sh2); + if (f1[0] == f2[0]) { + return lex_order(sh1, sh2); + } + return st.filtration(sh1)[0] < st.filtration(sh2)[0]; + }, + [](Simplex_handle sh) { return false; }); + unsigned int i = 0; + for (auto sh : st.filtration_simplex_range()){ + unsigned int j = 0; + for (auto v : st.simplex_vertex_range(sh)){ + BOOST_CHECK_EQUAL(v, simplices1[i][j]); + ++j; + } + BOOST_CHECK(st.filtration(sh) == fil1[i]); + ++i; + } + BOOST_CHECK_EQUAL(i, 15); + + st.initialize_filtration( + [&](Simplex_handle sh1, Simplex_handle sh2) { + auto f1 = st.filtration(sh1); + auto f2 = st.filtration(sh2); + auto v1 = f1.size() == 1 ? f1[0] : f1[1]; + auto v2 = f2.size() == 1 ? f2[0] : f2[1]; + if (v1 == v2) { + return lex_order(sh1, sh2); + } + return v1 < v2; + }, + [&st, max](Simplex_handle sh) { + auto f = st.filtration(sh); + return (f.size() == 1 ? f[0] : f[1]) == max; + }); + i = 0; + for (auto sh : st.filtration_simplex_range()){ + unsigned int j = 0; + for (auto v : st.simplex_vertex_range(sh)){ + BOOST_CHECK_EQUAL(v, simplices2[i][j]); + ++j; + } + BOOST_CHECK(st.filtration(sh) == fil2[i]); + ++i; + } + BOOST_CHECK_EQUAL(i, 7); +} diff --git a/src/Simplex_tree/test/test_vector_filtration_simplex_tree.h b/src/Simplex_tree/test/test_vector_filtration_simplex_tree.h new file mode 100644 index 0000000000..c826cce49f --- /dev/null +++ b/src/Simplex_tree/test/test_vector_filtration_simplex_tree.h @@ -0,0 +1,169 @@ +/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT. + * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details. + * Author(s): Hannah Schreiber + * + * Copyright (C) 2024 Inria + * + * Modification(s): + * - YYYY/MM Author: Description of the modification + */ + +#ifndef SIMPLEX_TREE_TEST_CUSTOM_SIMPLEX_TREE_H_ +#define SIMPLEX_TREE_TEST_CUSTOM_SIMPLEX_TREE_H_ + +#include +#include +#include +#include +#include + +#include + +namespace Gudhi { + +class Vector_filtration_value : public std::vector +{ + using Base = std::vector; + + public: + using const_iterator = Base::const_iterator; + + Vector_filtration_value() : Base() {} + Vector_filtration_value(std::size_t count) : Base(count) {} + Vector_filtration_value(std::initializer_list init) : Base(init) {} + Vector_filtration_value(const_iterator start, const_iterator end) : Base(start, end) {} + + friend bool unify_lifetimes(Vector_filtration_value& f1, const Vector_filtration_value& f2) { + int max = std::numeric_limits::max(); + bool f1_is_inf = f1.size() == 1 && f1[0] == max; + bool f2_is_inf = f2.size() == 1 && f2[0] == max; + if (f1_is_inf){ + if (f2_is_inf){ + return false; + } + f1 = f2; + return true; + } + if (f2_is_inf){ + return false; + } + // same size otherwise + unsigned int i = 0; + bool modified = false; + for (int v : f2){ + if (v < f1[i]){ + f1[i] = v; + modified = true; + } + ++i; + } + + return modified; + } + + friend bool intersect_lifetimes(Vector_filtration_value& f1, const Vector_filtration_value& f2) { + if (f1 < f2) { + f1 = f2; + return true; + } + return false; + } + + friend char* serialize_trivial(const Vector_filtration_value& value, char* start) + { + const auto length = value.size(); + const std::size_t arg_size = sizeof(int) * length; + const std::size_t type_size = sizeof(Vector_filtration_value::size_type); + memcpy(start, &length, type_size); + memcpy(start + type_size, value.data(), arg_size); + return start + arg_size + type_size; + } + + friend const char* deserialize_trivial(Vector_filtration_value& value, const char* start) + { + const std::size_t type_size = sizeof(Vector_filtration_value::size_type); + Vector_filtration_value::size_type length; + memcpy(&length, start, type_size); + std::size_t arg_size = sizeof(int) * length; + value.resize(length); + memcpy(value.data(), start + type_size, arg_size); + return start + arg_size + type_size; + } + + friend std::size_t get_serialization_size_of(const Vector_filtration_value& value) { + return sizeof(Vector_filtration_value::size_type) + sizeof(int) * value.size(); + } + + friend std::ostream &operator<<(std::ostream &stream, const Vector_filtration_value &f) { + if (f.empty()) { + stream << "[]"; + return stream; + } + stream << "["; + for (std::size_t i = 0; i < f.size() - 1; i++) { + stream << f[i] << ", "; + } + stream << f.back(); + stream << "]"; + return stream; + } +}; + +struct Simplex_tree_options_custom_fil_values_default { + typedef linear_indexing_tag Indexing_tag; + typedef std::int16_t Vertex_handle; + typedef Vector_filtration_value Filtration_value; + typedef std::int32_t Simplex_key; + static const bool store_key = false; + static const bool store_filtration = true; + static const bool contiguous_vertices = false; + static const bool link_nodes_by_label = false; + static const bool stable_simplex_handles = false; +}; + +struct Simplex_tree_options_custom_fil_values_fast_persistence { + typedef linear_indexing_tag Indexing_tag; + typedef std::int16_t Vertex_handle; + typedef Vector_filtration_value Filtration_value; + typedef std::int32_t Simplex_key; + static const bool store_key = true; + static const bool store_filtration = true; + static const bool contiguous_vertices = true; + static const bool link_nodes_by_label = false; + static const bool stable_simplex_handles = false; +}; + +struct Simplex_tree_options_custom_fil_values_full_featured { + typedef linear_indexing_tag Indexing_tag; + typedef std::int16_t Vertex_handle; + typedef Vector_filtration_value Filtration_value; + typedef std::int32_t Simplex_key; + static const bool store_key = true; + static const bool store_filtration = true; + static const bool contiguous_vertices = false; + static const bool link_nodes_by_label = true; + static const bool stable_simplex_handles = true; +}; + +} // namespace Gudhi + +namespace std { + +template<> +class numeric_limits { + public: + static constexpr bool has_infinity = true; + + static Gudhi::Vector_filtration_value infinity() noexcept { + return {std::numeric_limits::max()}; + }; + + static Gudhi::Vector_filtration_value max() noexcept(false) { + throw std::logic_error( + "The maximal value cannot be represented with no finite numbers of parameters."); + }; +}; + +} // namespace std + +#endif // SIMPLEX_TREE_TEST_CUSTOM_SIMPLEX_TREE_H_