From e05b10616f73a1ee2adcf249d4804c3931657855 Mon Sep 17 00:00:00 2001 From: hschreiber Date: Wed, 17 Apr 2024 15:58:58 +0200 Subject: [PATCH] doc + revisting fields --- .../include/gudhi/Fields/Multi_field.h | 139 ++-- .../gudhi/Fields/Multi_field_operators.h | 2 +- .../include/gudhi/Fields/Multi_field_shared.h | 205 ++--- .../include/gudhi/Fields/Multi_field_small.h | 560 ++++++------- .../Fields/Multi_field_small_operators.h | 2 +- .../gudhi/Fields/Multi_field_small_shared.h | 590 ++++++-------- .../include/gudhi/Fields/Z2_field.h | 82 +- .../include/gudhi/Fields/Zp_field.h | 748 +++++++++--------- .../include/gudhi/Fields/Zp_field_operators.h | 432 ++++++---- .../include/gudhi/Fields/Zp_field_shared.h | 699 ++++++++-------- .../test/Persistence_matrix_field_tests.cpp | 66 +- src/common/doc/main_page.md | 152 ++-- 12 files changed, 1791 insertions(+), 1886 deletions(-) diff --git a/src/Persistence_matrix/include/gudhi/Fields/Multi_field.h b/src/Persistence_matrix/include/gudhi/Fields/Multi_field.h index 6d0e3fd4b3..c1762d1989 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Multi_field.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Multi_field.h @@ -10,7 +10,7 @@ /** * @file Multi_field.h - * @author Hannah Schreiber + * @author Hannah Schreiber, Clément Maria * @brief Contains the @ref Multi_field_element class. */ @@ -36,10 +36,10 @@ namespace persistence_fields { * @tparam maximum Interval closed upper bound. */ template -class Multi_field_element -{ +class Multi_field_element { public: - using element_type = mpz_class; /**< Type for the elements in the field. */ + using element_type = mpz_class; /**< Type for the elements in the field. */ + using characteristic_type = element_type; /**< Type for the field characteristic. */ /** * @brief Default constructor. Sets the element to 0. @@ -47,19 +47,19 @@ class Multi_field_element Multi_field_element(); /** * @brief Constructor setting the element to the given value. - * + * * @param element Value of the element. */ - Multi_field_element(mpz_class element); + Multi_field_element(element_type element); /** * @brief Copy constructor. - * + * * @param toCopy Element to copy. */ Multi_field_element(const Multi_field_element& toCopy); /** * @brief Move constructor. - * + * * @param toMove Element to move. */ Multi_field_element(Multi_field_element&& toMove) noexcept; @@ -81,22 +81,22 @@ class Multi_field_element /** * @brief operator+= */ - friend void operator+=(Multi_field_element& f, mpz_class const v) { + friend void operator+=(Multi_field_element& f, element_type const v) { mpz_add(f.element_.get_mpz_t(), f.element_.get_mpz_t(), v.get_mpz_t()); mpz_mod(f.element_.get_mpz_t(), f.element_.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); } /** * @brief operator+ */ - friend Multi_field_element operator+(Multi_field_element f, mpz_class const v) { + friend Multi_field_element operator+(Multi_field_element f, element_type const v) { f += v; return f; } /** * @brief operator+ */ - friend mpz_class operator+(mpz_class v, Multi_field_element const& f) { - mpz_class e(v); + friend element_type operator+(element_type v, Multi_field_element const& f) { + element_type e(v); mpz_add(e.get_mpz_t(), e.get_mpz_t(), f.element_.get_mpz_t()); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e; @@ -119,22 +119,22 @@ class Multi_field_element /** * @brief operator-= */ - friend void operator-=(Multi_field_element& f, mpz_class const v) { + friend void operator-=(Multi_field_element& f, element_type const v) { mpz_sub(f.element_.get_mpz_t(), f.element_.get_mpz_t(), v.get_mpz_t()); mpz_mod(f.element_.get_mpz_t(), f.element_.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); } /** * @brief operator- */ - friend Multi_field_element operator-(Multi_field_element f, mpz_class const v) { + friend Multi_field_element operator-(Multi_field_element f, element_type const v) { f -= v; return f; } /** * @brief operator- */ - friend mpz_class operator-(mpz_class v, Multi_field_element const& f) { - mpz_class e(v); + friend element_type operator-(element_type v, Multi_field_element const& f) { + element_type e(v); if (e >= productOfAllCharacteristics_) mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); if (f.element_ > e) mpz_add(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); @@ -159,22 +159,22 @@ class Multi_field_element /** * @brief operator*= */ - friend void operator*=(Multi_field_element& f, mpz_class const v) { + friend void operator*=(Multi_field_element& f, element_type const v) { mpz_mul(f.element_.get_mpz_t(), f.element_.get_mpz_t(), v.get_mpz_t()); mpz_mod(f.element_.get_mpz_t(), f.element_.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); } /** * @brief operator* */ - friend Multi_field_element operator*(Multi_field_element f, mpz_class const v) { + friend Multi_field_element operator*(Multi_field_element f, element_type const v) { f *= v; return f; } /** * @brief operator* */ - friend mpz_class operator*(mpz_class v, Multi_field_element const& f) { - mpz_class e(v); + friend element_type operator*(element_type v, Multi_field_element const& f) { + element_type e(v); mpz_mul(e.get_mpz_t(), e.get_mpz_t(), f.element_.get_mpz_t()); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e; @@ -189,36 +189,18 @@ class Multi_field_element /** * @brief operator== */ - friend bool operator==(const mpz_class v, const Multi_field_element& f) { + friend bool operator==(const element_type v, const Multi_field_element& f) { if (v < productOfAllCharacteristics_) return v == f.element_; - mpz_class e(v); + element_type e(v); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e == f.element_; } /** * @brief operator== */ - friend bool operator==(const Multi_field_element& f, const mpz_class v) { + friend bool operator==(const Multi_field_element& f, const element_type v) { if (v < productOfAllCharacteristics_) return v == f.element_; - mpz_class e(v); - mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); - return e == f.element_; - } - /** - * @brief operator== - */ - friend bool operator==(const unsigned int v, const Multi_field_element& f) { - mpz_class e(v); - if (e < productOfAllCharacteristics_) return e == f.element_; - mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); - return e == f.element_; - } - /** - * @brief operator== - */ - friend bool operator==(const Multi_field_element& f, const unsigned int v) { - mpz_class e(v); - if (e < productOfAllCharacteristics_) return e == f.element_; + element_type e(v); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e == f.element_; } @@ -229,19 +211,11 @@ class Multi_field_element /** * @brief operator!= */ - friend bool operator!=(const mpz_class v, const Multi_field_element& f) { return !(v == f); } - /** - * @brief operator!= - */ - friend bool operator!=(const Multi_field_element& f, const mpz_class v) { return !(v == f); } + friend bool operator!=(const element_type v, const Multi_field_element& f) { return !(v == f); } /** * @brief operator!= */ - friend bool operator!=(const unsigned int v, const Multi_field_element& f) { return !(v == f); } - /** - * @brief operator!= - */ - friend bool operator!=(const Multi_field_element& f, const unsigned int v) { return !(v == f); } + friend bool operator!=(const Multi_field_element& f, const element_type v) { return !(v == f); } /** * @brief Assign operator. @@ -250,7 +224,7 @@ class Multi_field_element /** * @brief Assign operator. */ - Multi_field_element& operator=(const mpz_class value); + Multi_field_element& operator=(const element_type value); /** * @brief Swap operator. */ @@ -267,57 +241,58 @@ class Multi_field_element /** * @brief Returns the inverse of the element in the multi-field, see @cite boissonnat:hal-00922572. - * + * * @return The inverse. */ Multi_field_element get_inverse() const; /** * @brief Returns the inverse of the element with respect to a sub-product of the characteristics in the multi-field, * see @cite boissonnat:hal-00922572. - * + * * @param productOfCharacteristics Sub-product of the characteristics. * @return Pair of the inverse and the characteristic the inverse corresponds to. */ - std::pair get_partial_inverse(const mpz_class& productOfCharacteristics) const; + std::pair get_partial_inverse( + const characteristic_type& productOfCharacteristics) const; /** * @brief Returns the additive identity of a field. - * + * * @return The additive identity of a field. */ static Multi_field_element get_additive_identity(); /** * @brief Returns the multiplicative identity of a field. - * + * * @return The multiplicative identity of a field. */ static Multi_field_element get_multiplicative_identity(); /** * @brief Returns the partial multiplicative identity of the multi-field from the given product. * See @cite boissonnat:hal-00922572 for more details. - * + * * @param productOfCharacteristics Product of the different characteristics to take into account in the multi-field. * @return The partial multiplicative identity of the multi-field. */ - static Multi_field_element get_partial_multiplicative_identity(const mpz_class& productOfCharacteristics); + static Multi_field_element get_partial_multiplicative_identity(const characteristic_type& productOfCharacteristics); /** * @brief Returns the product of all characteristics. - * + * * @return The product of all characteristics. */ - static mpz_class get_characteristic(); + static characteristic_type get_characteristic(); /** * @brief Returns the value of the element. - * + * * @return Value of the element. */ - mpz_class get_value() const; + element_type get_value() const; // static constexpr bool handles_only_z2() { return false; } private: - mpz_class element_; + element_type element_; static inline const std::vector primes_ = []() { std::vector res; @@ -341,16 +316,16 @@ class Multi_field_element return res; }(); - static inline const mpz_class productOfAllCharacteristics_ = []() { - mpz_class res = 1; + static inline const characteristic_type productOfAllCharacteristics_ = []() { + characteristic_type res = 1; for (const auto p : primes_) { res *= p; } return res; }(); - static inline const std::vector partials_ = []() { - std::vector res; + static inline const std::vector partials_ = []() { + std::vector res; if (productOfAllCharacteristics_ == 1) return res; @@ -364,7 +339,7 @@ class Multi_field_element }(); // If I understood the paper well, multiplicativeID_ always equals to 1. But in Clement's code, // multiplicativeID_ is computed (see commented lambda function below). TODO: verify with Clement. - static inline const mpz_class multiplicativeID_ = 1; /*[](){ + static inline const element_type multiplicativeID_ = 1; /*[](){ mpz_class res = 0; for (unsigned int i = 0; i < partials_.size(); ++i){ res = (res + partials_[i]) % productOfAllCharacteristics_; @@ -387,7 +362,7 @@ inline Multi_field_element::Multi_field_element() : element_(0 } template -inline Multi_field_element::Multi_field_element(mpz_class element) : element_(element) { +inline Multi_field_element::Multi_field_element(element_type element) : element_(element) { static_assert(maximum >= 2, "Characteristics has to be positive."); static_assert(minimum <= maximum, "The given interval is not valid."); static_assert(minimum != maximum || _is_prime(minimum), "The given interval does not contain a prime number."); @@ -415,7 +390,8 @@ inline Multi_field_element& Multi_field_element -inline Multi_field_element& Multi_field_element::operator=(mpz_class const value) { +inline Multi_field_element& Multi_field_element::operator=( + element_type const value) { mpz_mod(element_.get_mpz_t(), value.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return *this; } @@ -436,16 +412,17 @@ inline Multi_field_element Multi_field_element -inline std::pair, mpz_class> -Multi_field_element::get_partial_inverse(const mpz_class& productOfCharacteristics) const { - mpz_class QR; +inline std::pair, + typename Multi_field_element::characteristic_type> +Multi_field_element::get_partial_inverse(const characteristic_type& productOfCharacteristics) const { + characteristic_type QR; mpz_gcd(QR.get_mpz_t(), element_.get_mpz_t(), productOfCharacteristics.get_mpz_t()); // QR <- gcd(x,QS) if (QR == productOfCharacteristics) return {Multi_field_element(), multiplicativeID_}; // partial inverse is 0 - mpz_class QT = productOfCharacteristics / QR; + characteristic_type QT = productOfCharacteristics / QR; - mpz_class inv_qt; + characteristic_type inv_qt; mpz_invert(inv_qt.get_mpz_t(), element_.get_mpz_t(), QT.get_mpz_t()); auto res = get_partial_multiplicative_identity(QT); @@ -466,7 +443,7 @@ inline Multi_field_element Multi_field_element inline Multi_field_element Multi_field_element::get_partial_multiplicative_identity( - const mpz_class& productOfCharacteristics) { + const characteristic_type& productOfCharacteristics) { if (productOfCharacteristics == 0) { return Multi_field_element(multiplicativeID_); } @@ -480,12 +457,14 @@ inline Multi_field_element Multi_field_element -inline mpz_class Multi_field_element::get_characteristic() { +inline typename Multi_field_element::characteristic_type +Multi_field_element::get_characteristic() { return productOfAllCharacteristics_; } template -inline mpz_class Multi_field_element::get_value() const { +inline typename Multi_field_element::element_type Multi_field_element::get_value() + const { return element_; } diff --git a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_operators.h b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_operators.h index 58a6e9473a..3053f5c388 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_operators.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_operators.h @@ -10,7 +10,7 @@ /** * @file Multi_field_operators.h - * @author Hannah Schreiber + * @author Hannah Schreiber, Clément Maria * @brief Contains the @ref Multi_field_operators class. */ diff --git a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_shared.h b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_shared.h index 8a6c772cde..ff352f7ed3 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_shared.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_shared.h @@ -10,7 +10,7 @@ /** * @file Multi_field_shared.h - * @author Hannah Schreiber + * @author Hannah Schreiber, Clément Maria * @brief Contains the @ref Shared_multi_field_element class. */ @@ -37,7 +37,8 @@ namespace persistence_fields { class Shared_multi_field_element { public: - using element_type = mpz_class; /**< Type for the elements in the field. */ + using element_type = mpz_class; /**< Type for the elements in the field. */ + using characteristic_type = element_type; /**< Type for the field characteristic. */ /** * @brief Default constructor. Sets the element to 0. @@ -48,7 +49,7 @@ class Shared_multi_field_element * * @param element Value of the element. */ - Shared_multi_field_element(mpz_class element); + Shared_multi_field_element(element_type element); /** * @brief Copy constructor. * @@ -73,9 +74,6 @@ class Shared_multi_field_element /** * @brief operator+= - * - * @param f1 Shared_multi_field_element - * @param f2 Shared_multi_field_element */ friend void operator+=(Shared_multi_field_element& f1, Shared_multi_field_element const& f2) { mpz_add(f1.element_.get_mpz_t(), f1.element_.get_mpz_t(), f2.element_.get_mpz_t()); @@ -83,10 +81,6 @@ class Shared_multi_field_element } /** * @brief operator+ - * - * @param f1 Shared_multi_field_element - * @param f2 const Shared_multi_field_element & - * @return Shared_multi_field_element */ friend Shared_multi_field_element operator+(Shared_multi_field_element f1, Shared_multi_field_element const& f2) { f1 += f2; @@ -94,34 +88,23 @@ class Shared_multi_field_element } /** * @brief operator+= - * - * @param f Shared_multi_field_element& - * @param v const mpz_class */ - friend void operator+=(Shared_multi_field_element& f, mpz_class const v) { + friend void operator+=(Shared_multi_field_element& f, element_type const v) { mpz_add(f.element_.get_mpz_t(), f.element_.get_mpz_t(), v.get_mpz_t()); mpz_mod(f.element_.get_mpz_t(), f.element_.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); } /** * @brief operator+ - * - * @param f Shared_multi_field_element - * @param v const mpz_class - * @return Shared_multi_field_element */ - friend Shared_multi_field_element operator+(Shared_multi_field_element f, mpz_class const v) { + friend Shared_multi_field_element operator+(Shared_multi_field_element f, element_type const v) { f += v; return f; } /** * @brief operator+ - * - * @param v mpz_class - * @param f const Shared_multi_field_element & - * @return mpz_class */ - friend mpz_class operator+(mpz_class v, Shared_multi_field_element const& f) { - mpz_class e(v); + friend element_type operator+(element_type v, Shared_multi_field_element const& f) { + element_type e(v); mpz_add(e.get_mpz_t(), e.get_mpz_t(), f.element_.get_mpz_t()); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e; @@ -129,9 +112,6 @@ class Shared_multi_field_element /** * @brief operator-= - * - * @param f1 Shared_multi_field_element& - * @param f2 const Shared_multi_field_element & */ friend void operator-=(Shared_multi_field_element& f1, Shared_multi_field_element const& f2) { mpz_sub(f1.element_.get_mpz_t(), f1.element_.get_mpz_t(), f2.element_.get_mpz_t()); @@ -139,10 +119,6 @@ class Shared_multi_field_element } /** * @brief operator- - * - * @param f1 Shared_multi_field_element - * @param f2 const Shared_multi_field_element & - * @return Shared_multi_field_element */ friend Shared_multi_field_element operator-(Shared_multi_field_element f1, Shared_multi_field_element const& f2) { f1 -= f2; @@ -150,34 +126,23 @@ class Shared_multi_field_element } /** * @brief operator-= - * - * @param f Shared_multi_field_element& - * @param v const mpz_class */ - friend void operator-=(Shared_multi_field_element& f, mpz_class const v) { + friend void operator-=(Shared_multi_field_element& f, element_type const v) { mpz_sub(f.element_.get_mpz_t(), f.element_.get_mpz_t(), v.get_mpz_t()); mpz_mod(f.element_.get_mpz_t(), f.element_.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); } /** * @brief operator- - * - * @param f Shared_multi_field_element - * @param v const mpz_class - * @return Shared_multi_field_element */ - friend Shared_multi_field_element operator-(Shared_multi_field_element f, mpz_class const v) { + friend Shared_multi_field_element operator-(Shared_multi_field_element f, element_type const v) { f -= v; return f; } /** * @brief operator- - * - * @param v mpz_class - * @param f const Shared_multi_field_element & - * @return mpz_class */ - friend mpz_class operator-(mpz_class v, Shared_multi_field_element const& f) { - mpz_class e(v); + friend element_type operator-(element_type v, Shared_multi_field_element const& f) { + element_type e(v); if (e >= productOfAllCharacteristics_) mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); if (f.element_ > e) mpz_add(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); @@ -187,9 +152,6 @@ class Shared_multi_field_element /** * @brief operator*= - * - * @param f1 Shared_multi_field_element& - * @param f2 const Shared_multi_field_element & */ friend void operator*=(Shared_multi_field_element& f1, Shared_multi_field_element const& f2) { mpz_mul(f1.element_.get_mpz_t(), f1.element_.get_mpz_t(), f2.element_.get_mpz_t()); @@ -197,10 +159,6 @@ class Shared_multi_field_element } /** * @brief operator* - * - * @param f1 Shared_multi_field_element - * @param f2 const Shared_multi_field_element & - * @return Shared_multi_field_element */ friend Shared_multi_field_element operator*(Shared_multi_field_element f1, Shared_multi_field_element const& f2) { f1 *= f2; @@ -208,34 +166,23 @@ class Shared_multi_field_element } /** * @brief operator*= - * - * @param f Shared_multi_field_element& - * @param v const mpz_class */ - friend void operator*=(Shared_multi_field_element& f, mpz_class const v) { + friend void operator*=(Shared_multi_field_element& f, element_type const v) { mpz_mul(f.element_.get_mpz_t(), f.element_.get_mpz_t(), v.get_mpz_t()); mpz_mod(f.element_.get_mpz_t(), f.element_.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); } /** * @brief operator* - * - * @param f Shared_multi_field_element - * @param v const mpz_class - * @return Shared_multi_field_element */ - friend Shared_multi_field_element operator*(Shared_multi_field_element f, mpz_class const v) { + friend Shared_multi_field_element operator*(Shared_multi_field_element f, element_type const v) { f *= v; return f; } /** * @brief operator* - * - * @param v mpz_class - * @param f const Shared_multi_field_element & - * @return mpz_class */ - friend mpz_class operator*(mpz_class v, Shared_multi_field_element const& f) { - mpz_class e(v); + friend element_type operator*(element_type v, Shared_multi_field_element const& f) { + element_type e(v); mpz_mul(e.get_mpz_t(), e.get_mpz_t(), f.element_.get_mpz_t()); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e; @@ -243,108 +190,42 @@ class Shared_multi_field_element /** * @brief operator== - * - * @param f1 const Shared_multi_field_element& - * @param f2 const Shared_multi_field_element& - * @return bool */ friend bool operator==(const Shared_multi_field_element& f1, const Shared_multi_field_element& f2) { return f1.element_ == f2.element_; } /** * @brief operator== - * - * @param f1 const Shared_multi_field_element& - * @param f2 const Shared_multi_field_element& - * @return bool */ - friend bool operator==(const mpz_class v, const Shared_multi_field_element& f) { + friend bool operator==(const element_type v, const Shared_multi_field_element& f) { if (v < productOfAllCharacteristics_) return v == f.element_; - mpz_class e(v); + element_type e(v); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e == f.element_; } /** * @brief operator== - * - * @param f const Shared_multi_field_element& - * @param v const mpz_class - * @return bool */ - friend bool operator==(const Shared_multi_field_element& f, const mpz_class v) { + friend bool operator==(const Shared_multi_field_element& f, const element_type v) { if (v < productOfAllCharacteristics_) return v == f.element_; - mpz_class e(v); - mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); - return e == f.element_; - } - /** - * @brief operator== - * - * @param v const unsigned int - * @param f const Shared_multi_field_element& - * @return bool - */ - friend bool operator==(const unsigned int v, const Shared_multi_field_element& f) { - mpz_class e(v); - if (e < productOfAllCharacteristics_) return e == f.element_; - mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); - return e == f.element_; - } - /** - * @brief operator== - * - * @param f const Shared_multi_field_element& - * @param v const unsigned int - * @return bool - */ - friend bool operator==(const Shared_multi_field_element& f, const unsigned int v) { - mpz_class e(v); - if (e < productOfAllCharacteristics_) return e == f.element_; + element_type e(v); mpz_mod(e.get_mpz_t(), e.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return e == f.element_; } /** * @brief operator!= - * - * @param f1 const Shared_multi_field_element& - * @param f2 const Shared_multi_field_element& - * @return bool */ friend bool operator!=(const Shared_multi_field_element& f1, const Shared_multi_field_element& f2) { return !(f1 == f2); } /** * @brief operator!= - * - * @param f1 const Shared_multi_field_element& - * @param f2 const Shared_multi_field_element& - * @return bool - */ - friend bool operator!=(const mpz_class v, const Shared_multi_field_element& f) { return !(v == f); } - /** - * @brief operator!= - * - * @param f const Shared_multi_field_element& - * @param v const mpz_class - * @return bool */ - friend bool operator!=(const Shared_multi_field_element& f, const mpz_class v) { return !(v == f); } + friend bool operator!=(const element_type v, const Shared_multi_field_element& f) { return !(v == f); } /** * @brief operator!= - * - * @param v const unsigned int - * @param f const Shared_multi_field_element& - * @return bool */ - friend bool operator!=(const unsigned int v, const Shared_multi_field_element& f) { return !(v == f); } - /** - * @brief operator!= - * - * @param f const Shared_multi_field_element& - * @param v const unsigned int - * @return bool - */ - friend bool operator!=(const Shared_multi_field_element& f, const unsigned int v) { return !(v == f); } + friend bool operator!=(const Shared_multi_field_element& f, const element_type v) { return !(v == f); } /** * @brief Assign operator. @@ -353,7 +234,7 @@ class Shared_multi_field_element /** * @brief Assign operator. */ - Shared_multi_field_element& operator=(const mpz_class value); + Shared_multi_field_element& operator=(const element_type value); /** * @brief Swap operator. */ @@ -383,7 +264,8 @@ class Shared_multi_field_element * @param productOfCharacteristics Sub-product of the characteristics. * @return Pair of the inverse and the characteristic the inverse corresponds to. */ - std::pair get_partial_inverse(const mpz_class& productOfCharacteristics) const; + std::pair get_partial_inverse( + const characteristic_type& productOfCharacteristics) const; /** * @brief Returns the additive identity of a field. @@ -404,36 +286,37 @@ class Shared_multi_field_element * @param productOfCharacteristics Product of the different characteristics to take into account in the multi-field. * @return The partial multiplicative identity of the multi-field. */ - static Shared_multi_field_element get_partial_multiplicative_identity(const mpz_class& productOfCharacteristics); + static Shared_multi_field_element get_partial_multiplicative_identity( + const characteristic_type& productOfCharacteristics); /** * @brief Returns the product of all characteristics. * * @return The product of all characteristics. */ - static mpz_class get_characteristic(); + static characteristic_type get_characteristic(); /** * @brief Returns the value of the element. * * @return Value of the element. */ - mpz_class get_value() const; + element_type get_value() const; // static constexpr bool handles_only_z2() { return false; } private: - mpz_class element_; /**< Element. */ + element_type element_; /**< Element. */ static inline std::vector primes_; /**< All characteristics. */ - static inline mpz_class productOfAllCharacteristics_ = 0; /**< Product of all characteristics. */ - static inline std::vector partials_; /**< Partial products of the characteristics. */ - static inline const mpz_class multiplicativeID_ = 1; /**< Multiplicative identity. */ + static inline characteristic_type productOfAllCharacteristics_ = 0; /**< Product of all characteristics. */ + static inline std::vector partials_; /**< Partial products of the characteristics. */ + static inline const element_type multiplicativeID_ = 1; /**< Multiplicative identity. */ static constexpr bool _is_prime(const int p); }; inline Shared_multi_field_element::Shared_multi_field_element() : element_(0) {} -inline Shared_multi_field_element::Shared_multi_field_element(mpz_class element) : element_(element) { +inline Shared_multi_field_element::Shared_multi_field_element(element_type element) : element_(element) { mpz_mod(element_.get_mpz_t(), element_.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); } @@ -494,7 +377,7 @@ inline Shared_multi_field_element& Shared_multi_field_element::operator=(Shared_ return *this; } -inline Shared_multi_field_element& Shared_multi_field_element::operator=(mpz_class const value) { +inline Shared_multi_field_element& Shared_multi_field_element::operator=(element_type const value) { mpz_mod(element_.get_mpz_t(), value.get_mpz_t(), productOfAllCharacteristics_.get_mpz_t()); return *this; } @@ -507,16 +390,16 @@ inline Shared_multi_field_element Shared_multi_field_element::get_inverse() cons return get_partial_inverse(productOfAllCharacteristics_).first; } -inline std::pair Shared_multi_field_element::get_partial_inverse( - const mpz_class& productOfCharacteristics) const { - mpz_class QR; +inline std::pair +Shared_multi_field_element::get_partial_inverse(const characteristic_type& productOfCharacteristics) const { + element_type QR; mpz_gcd(QR.get_mpz_t(), element_.get_mpz_t(), productOfCharacteristics.get_mpz_t()); // QR <- gcd(x,QS) if (QR == productOfCharacteristics) return {Shared_multi_field_element(), multiplicativeID_}; // partial inverse is 0 - mpz_class QT = productOfCharacteristics / QR; + element_type QT = productOfCharacteristics / QR; - mpz_class inv_qt; + element_type inv_qt; mpz_invert(inv_qt.get_mpz_t(), element_.get_mpz_t(), QT.get_mpz_t()); auto res = get_partial_multiplicative_identity(QT); @@ -534,7 +417,7 @@ inline Shared_multi_field_element Shared_multi_field_element::get_multiplicative } inline Shared_multi_field_element Shared_multi_field_element::get_partial_multiplicative_identity( - const mpz_class& productOfCharacteristics) { + const characteristic_type& productOfCharacteristics) { if (productOfCharacteristics == 0) { return Shared_multi_field_element(multiplicativeID_); } @@ -547,9 +430,11 @@ inline Shared_multi_field_element Shared_multi_field_element::get_partial_multip return mult; } -inline mpz_class Shared_multi_field_element::get_characteristic() { return productOfAllCharacteristics_; } +inline Shared_multi_field_element::characteristic_type Shared_multi_field_element::get_characteristic() { + return productOfAllCharacteristics_; +} -inline mpz_class Shared_multi_field_element::get_value() const { return element_; } +inline Shared_multi_field_element::element_type Shared_multi_field_element::get_value() const { return element_; } inline constexpr bool Shared_multi_field_element::_is_prime(const int p) { if (p <= 1) return false; diff --git a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small.h b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small.h index 8d15a2cc93..db8594b318 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small.h @@ -10,7 +10,7 @@ /** * @file Multi_field_small.h - * @author Hannah Schreiber + * @author Hannah Schreiber, Clément Maria * @brief Contains the @ref Multi_field_element_with_small_characteristics class. */ @@ -31,46 +31,59 @@ namespace persistence_fields { * @ingroup persistence_fields * * @brief Class representing an element of a multi-field, such that the product of all characteristics fits into - * an unsigned int. The characteristics will corresponds to all prime numbers in the interval given as template. - * + * the given @p Unsigned_integer_type template argument. The characteristics will corresponds to all prime numbers + * in the interval given as other template arguments. + * * @tparam minimum Interval closed lower bound. * @tparam maximum Interval closed upper bound. + * @tparam Unsigned_integer_type A native unsigned integer type: unsigned int, long unsigned int, etc. + * Will be used as the field element type. */ -template +template > > class Multi_field_element_with_small_characteristics { public: - using element_type = unsigned int; /**< Type for the elements in the field. */ + using element_type = Unsigned_integer_type; /**< Type for the elements in the field. */ + using characteristic_type = element_type; /**< Type for the field characteristic. */ template using isInteger = std::enable_if_t >; /** * @brief Default constructor. Sets the element to 0. */ - Multi_field_element_with_small_characteristics(); - /** - * @brief Constructor setting the element to the given value. - * - * @param element Value of the element. - */ - Multi_field_element_with_small_characteristics(unsigned int element); + Multi_field_element_with_small_characteristics() : element_(0) { + static_assert(maximum >= 2, "Characteristics have to be positive."); + static_assert(minimum <= maximum, "The given interval is not valid."); + static_assert(minimum != maximum || _is_prime(minimum), "The given interval does not contain a prime number."); + static_assert(productOfAllCharacteristics_ != 1, "The given interval does not contain a prime number."); + } /** * @brief Constructor setting the element to the given value. - * + * * @param element Value of the element. */ - Multi_field_element_with_small_characteristics(int element); + template > + Multi_field_element_with_small_characteristics(Integer_type element) + : element_(_get_value(element)) { + static_assert(maximum >= 2, "Characteristics has to be positive."); + static_assert(minimum <= maximum, "The given interval is not valid."); + static_assert(minimum != maximum || _is_prime(minimum), "The given interval does not contain a prime number."); + static_assert(productOfAllCharacteristics_ != 1, "The given interval does not contain a prime number."); + } /** * @brief Copy constructor. - * + * * @param toCopy Element to copy. */ - Multi_field_element_with_small_characteristics(const Multi_field_element_with_small_characteristics& toCopy); + Multi_field_element_with_small_characteristics(const Multi_field_element_with_small_characteristics& toCopy) + : element_(toCopy.element_) {} /** * @brief Move constructor. - * + * * @param toMove Element to move. */ - Multi_field_element_with_small_characteristics(Multi_field_element_with_small_characteristics&& toMove) noexcept; + Multi_field_element_with_small_characteristics(Multi_field_element_with_small_characteristics&& toMove) noexcept + : element_(std::exchange(toMove.element_, 0)) {} /** * @brief operator+= @@ -89,29 +102,33 @@ class Multi_field_element_with_small_characteristics { } /** * @brief operator+= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - friend void operator+=(Multi_field_element_with_small_characteristics& f, unsigned int const v) { - f.element_ = _add(f.element_, v < productOfAllCharacteristics_ ? v : v % productOfAllCharacteristics_); + template > + friend void operator+=(Multi_field_element_with_small_characteristics& f, const Integer_type& v) { + f.element_ = _add(f.element_, _get_value(v)); } /** * @brief operator+ - * - * @warning @p v is assumed to be positive and will be casted into an unsigned int. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend Multi_field_element_with_small_characteristics operator+(Multi_field_element_with_small_characteristics f, - const Integer_type v) { + const Integer_type& v) { f += v; return f; } /** * @brief operator+ + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend Integer_type operator+(Integer_type v, Multi_field_element_with_small_characteristics const& f) { - v += f.element_; - v %= productOfAllCharacteristics_; - return v; + friend Integer_type operator+(const Integer_type& v, Multi_field_element_with_small_characteristics f) { + f += v; + return f.element_; } /** @@ -131,32 +148,32 @@ class Multi_field_element_with_small_characteristics { } /** * @brief operator-= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - friend void operator-=(Multi_field_element_with_small_characteristics& f, unsigned int const v) { - f.element_ = _substract(f.element_, v < productOfAllCharacteristics_ ? v : v % productOfAllCharacteristics_); + template > + friend void operator-=(Multi_field_element_with_small_characteristics& f, const Integer_type& v) { + f.element_ = _substract(f.element_, _get_value(v)); } /** * @brief operator- - * - * @warning @p v is assumed to be positive and will be casted into an unsigned int. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend Multi_field_element_with_small_characteristics operator-(Multi_field_element_with_small_characteristics f, - const Integer_type v) { + const Integer_type& v) { f -= v; return f; } /** * @brief operator- - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend Integer_type operator-(Integer_type v, Multi_field_element_with_small_characteristics const& f) { - if (v >= productOfAllCharacteristics_) v %= productOfAllCharacteristics_; - if (f.element_ > v) v += productOfAllCharacteristics_; - v -= f.element_; - return v; + friend Integer_type operator-(const Integer_type& v, const Multi_field_element_with_small_characteristics& f) { + return _substract(_get_value(v), f.element_); } /** @@ -176,45 +193,33 @@ class Multi_field_element_with_small_characteristics { } /** * @brief operator*= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - friend void operator*=(Multi_field_element_with_small_characteristics& f, unsigned int const v) { - f.element_ = _multiply(f.element_, v < productOfAllCharacteristics_ ? v : v % productOfAllCharacteristics_); + template > + friend void operator*=(Multi_field_element_with_small_characteristics& f, const Integer_type& v) { + f.element_ = _multiply(f.element_, _get_value(v)); } /** * @brief operator* - * - * @warning @p v is assumed to be positive and will be casted into an unsigned int. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend Multi_field_element_with_small_characteristics operator*(Multi_field_element_with_small_characteristics f, - const Integer_type v) { + const Integer_type& v) { f *= v; return f; } /** * @brief operator* - * - * @warning Uses bitwise operations on @p v, so be carefull with signed integers. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend Integer_type operator*(Integer_type v, Multi_field_element_with_small_characteristics const& f) { - unsigned int b = f.element_; - unsigned int res = 0; - unsigned int temp_b; - - while (v != 0) { - if (v & 1) { - if (b >= productOfAllCharacteristics_ - res) res -= productOfAllCharacteristics_; - res += b; - } - v >>= 1; - - temp_b = b; - if (b >= productOfAllCharacteristics_ - b) temp_b -= productOfAllCharacteristics_; - b += temp_b; - } - - return res; + friend Integer_type operator*(const Integer_type& v, Multi_field_element_with_small_characteristics f) { + f *= v; + return f.element_; } /** @@ -226,23 +231,21 @@ class Multi_field_element_with_small_characteristics { } /** * @brief operator== - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend bool operator==(const Integer_type v, const Multi_field_element_with_small_characteristics& f) { - if (v < productOfAllCharacteristics_) return v == f.element_; - return (v % productOfAllCharacteristics_) == f.element_; + return _get_value(v) == f.element_; } /** * @brief operator== - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend bool operator==(const Multi_field_element_with_small_characteristics& f, const Integer_type v) { - if (v < productOfAllCharacteristics_) return v == f.element_; - return (v % productOfAllCharacteristics_) == f.element_; + return _get_value(v) == f.element_; } /** * @brief operator!= @@ -253,8 +256,8 @@ class Multi_field_element_with_small_characteristics { } /** * @brief operator!= - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend bool operator!=(const Integer_type v, const Multi_field_element_with_small_characteristics& f) { @@ -262,8 +265,8 @@ class Multi_field_element_with_small_characteristics { } /** * @brief operator!= - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend bool operator!=(const Multi_field_element_with_small_characteristics& f, const Integer_type v) { @@ -273,11 +276,20 @@ class Multi_field_element_with_small_characteristics { /** * @brief Assign operator. */ - Multi_field_element_with_small_characteristics& operator=(Multi_field_element_with_small_characteristics other); + Multi_field_element_with_small_characteristics& operator=(Multi_field_element_with_small_characteristics other) { + std::swap(element_, other.element_); + return *this; + } /** * @brief Assign operator. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - Multi_field_element_with_small_characteristics& operator=(const unsigned int value); + template > + Multi_field_element_with_small_characteristics& operator=(const Integer_type& value) { + element_ = _get_value(value); + return *this; + } /** * @brief Swap operator. */ @@ -289,96 +301,206 @@ class Multi_field_element_with_small_characteristics { /** * @brief Casts the element into an unsigned int. */ - operator unsigned int() const; + operator unsigned int() const { return element_; } /** * @brief Returns the inverse of the element in the multi-field, see @cite boissonnat:hal-00922572. - * + * * @return The inverse. */ - Multi_field_element_with_small_characteristics get_inverse() const; + Multi_field_element_with_small_characteristics get_inverse() const { + return get_partial_inverse(productOfAllCharacteristics_).first; + } /** * @brief Returns the inverse of the element with respect to a sub-product of the characteristics in the multi-field, * see @cite boissonnat:hal-00922572. - * + * * @param productOfCharacteristics Sub-product of the characteristics. * @return Pair of the inverse and the characteristic the inverse corresponds to. */ - std::pair get_partial_inverse( - unsigned int productOfCharacteristics) const; + std::pair get_partial_inverse( + characteristic_type productOfCharacteristics) const { + characteristic_type gcd = std::gcd(element_, productOfAllCharacteristics_); + + if (gcd == productOfCharacteristics) + return {Multi_field_element_with_small_characteristics(), multiplicativeID_}; // partial inverse is 0 + + characteristic_type QT = productOfCharacteristics / gcd; + + const element_type inv_qt = _get_inverse(element_, QT); + + auto res = get_partial_multiplicative_identity(QT); + res *= inv_qt; + + return {res, QT}; + } /** * @brief Returns the additive identity of a field. - * + * * @return The additive identity of a field. */ - static Multi_field_element_with_small_characteristics get_additive_identity(); + static Multi_field_element_with_small_characteristics get_additive_identity() { + return Multi_field_element_with_small_characteristics(); + } /** * @brief Returns the multiplicative identity of a field. - * + * * @return The multiplicative identity of a field. */ - static Multi_field_element_with_small_characteristics get_multiplicative_identity(); + static Multi_field_element_with_small_characteristics get_multiplicative_identity() { + return Multi_field_element_with_small_characteristics(multiplicativeID_); + } /** * @brief Returns the partial multiplicative identity of the multi-field from the given product. * See @cite boissonnat:hal-00922572 for more details. - * + * * @param productOfCharacteristics Product of the different characteristics to take into account in the multi-field. * @return The partial multiplicative identity of the multi-field. */ static Multi_field_element_with_small_characteristics get_partial_multiplicative_identity( - const mpz_class& productOfCharacteristics); + const mpz_class& productOfCharacteristics) { + if (productOfCharacteristics == 0) { + return Multi_field_element_with_small_characteristics(multiplicativeID_); + } + Multi_field_element_with_small_characteristics mult; + for (characteristic_type idx = 0; idx < primes_.size(); ++idx) { + if ((productOfCharacteristics % primes_[idx]) == 0) { + mult += partials_[idx]; + } + } + return mult; + } /** * @brief Returns the product of all characteristics. - * + * * @return The product of all characteristics. */ - static constexpr unsigned int get_characteristic(); + static constexpr characteristic_type get_characteristic() { return productOfAllCharacteristics_; } /** * @brief Returns the value of the element. - * + * * @return Value of the element. */ - unsigned int get_value() const; + element_type get_value() const { return element_; } // static constexpr bool handles_only_z2() { return false; } private: - static constexpr bool _is_prime(const int p); - static constexpr unsigned int _multiply(unsigned int a, unsigned int b); - static constexpr unsigned int _add(unsigned int element, unsigned int v); - static constexpr unsigned int _substract(unsigned int element, unsigned int v); - static constexpr int _get_inverse(unsigned int element, const unsigned int mod); - - unsigned int element_; - static inline const std::vector primes_ = []() { - std::vector res; - for (unsigned int i = minimum; i <= maximum; ++i) { + static constexpr bool _is_prime(const unsigned int p) { + if (p <= 1) return false; + if (p <= 3) return true; + if (p % 2 == 0 || p % 3 == 0) return false; + + for (long i = 5; i * i <= p; i = i + 6) + if (p % i == 0 || p % (i + 2) == 0) return false; + + return true; + } + static constexpr element_type _multiply(element_type a, element_type b) { + element_type res = 0; + element_type temp_b = 0; + + if (b < a) std::swap(a, b); + + while (a != 0) { + if (a & 1) { + /* Add b to res, modulo m, without overflow */ + if (b >= productOfAllCharacteristics_ - res) res -= productOfAllCharacteristics_; + res += b; + } + a >>= 1; + + /* Double b, modulo m */ + temp_b = b; + if (b >= productOfAllCharacteristics_ - b) temp_b -= productOfAllCharacteristics_; + b += temp_b; + } + return res; + } + static constexpr element_type _add(element_type element, element_type v) { + if (UINT_MAX - element < v) { + // automatic unsigned integer overflow behaviour will make it work + element += v; + element -= productOfAllCharacteristics_; + return element; + } + + element += v; + if (element >= productOfAllCharacteristics_) element -= productOfAllCharacteristics_; + + return element; + } + static constexpr element_type _substract(element_type element, element_type v) { + if (element < v) { + element += productOfAllCharacteristics_; + } + element -= v; + + return element; + } + static constexpr int _get_inverse(element_type element, const element_type mod) { + // to solve: Ax + My = 1 + int M = mod; + int A = element; + int y = 0, x = 1; + // extended euclidien division + while (A > 1) { + int quotient = A / M; + int temp = M; + + M = A % M, A = temp; + temp = y; + + y = x - quotient * y; + x = temp; + } + + if (x < 0) x += mod; + + return x; + } + + template > + static constexpr element_type _get_value(Integer_type e) { + if constexpr (std::is_signed_v){ + if (e < -static_cast(productOfAllCharacteristics_)) e = e % productOfAllCharacteristics_; + if (e < 0) return e += productOfAllCharacteristics_; + return e < static_cast(productOfAllCharacteristics_) ? e : e % productOfAllCharacteristics_; + } else { + return e < productOfAllCharacteristics_ ? e : e % productOfAllCharacteristics_; + } + } + + element_type element_; + static inline const std::vector primes_ = []() { + std::vector res; + for (characteristic_type i = minimum; i <= maximum; ++i) { if (_is_prime(i)) { res.push_back(i); } } return res; }(); - static inline constexpr unsigned int productOfAllCharacteristics_ = []() { - unsigned int res = 1; - for (unsigned int i = minimum; i <= maximum; ++i) { + static inline constexpr characteristic_type productOfAllCharacteristics_ = []() { + characteristic_type res = 1; + for (characteristic_type i = minimum; i <= maximum; ++i) { if (_is_prime(i)) { res *= i; } } return res; }(); - static inline const std::vector partials_ = []() { - std::vector res; + static inline const std::vector partials_ = []() { + std::vector res; if (productOfAllCharacteristics_ == 1) return res; - for (unsigned int i = 0; i < primes_.size(); ++i) { - unsigned int p = primes_[i]; - unsigned int base = productOfAllCharacteristics_ / p; - unsigned int exp = p - 1; + for (characteristic_type i = 0; i < primes_.size(); ++i) { + characteristic_type p = primes_[i]; + characteristic_type base = productOfAllCharacteristics_ / p; + characteristic_type exp = p - 1; res.push_back(1); while (exp > 0) { @@ -394,7 +516,7 @@ class Multi_field_element_with_small_characteristics { }(); // If I understood the paper well, multiplicativeID_ always equals to 1. But in Clement's code, // multiplicativeID_ is computed (see commented lambda function below). TODO: verify with Clement. - static inline constexpr unsigned int multiplicativeID_ = 1; /*= [](){ + static inline constexpr element_type multiplicativeID_ = 1; /*= [](){ unsigned int res = 0; for (unsigned int i = 0; i < partials_.size(); ++i){ res = (res + partials_[i]) % productOfAllCharacteristics_; @@ -404,216 +526,6 @@ class Multi_field_element_with_small_characteristics { }();*/ }; -template -inline Multi_field_element_with_small_characteristics::Multi_field_element_with_small_characteristics() - : element_(0) { - static_assert(maximum >= 2, "Characteristics have to be positive."); - static_assert(minimum <= maximum, "The given interval is not valid."); - static_assert(minimum != maximum || _is_prime(minimum), "The given interval does not contain a prime number."); - static_assert(productOfAllCharacteristics_ != 1, "The given interval does not contain a prime number."); -} - -template -inline Multi_field_element_with_small_characteristics::Multi_field_element_with_small_characteristics( - unsigned int element) - : element_(element % productOfAllCharacteristics_) { - static_assert(maximum >= 2, "Characteristics has to be positive."); - static_assert(minimum <= maximum, "The given interval is not valid."); - static_assert(minimum != maximum || _is_prime(minimum), "The given interval does not contain a prime number."); - static_assert(productOfAllCharacteristics_ != 1, "The given interval does not contain a prime number."); -} - -template -inline Multi_field_element_with_small_characteristics::Multi_field_element_with_small_characteristics( - int element) - : element_(element % productOfAllCharacteristics_) { - static_assert(maximum >= 2, "Characteristics has to be positive."); - static_assert(minimum <= maximum, "The given interval is not valid."); - static_assert(minimum != maximum || _is_prime(minimum), "The given interval does not contain a prime number."); - static_assert(productOfAllCharacteristics_ != 1, "The given interval does not contain a prime number."); -} - -template -inline Multi_field_element_with_small_characteristics::Multi_field_element_with_small_characteristics( - const Multi_field_element_with_small_characteristics& toCopy) - : element_(toCopy.element_) {} - -template -inline Multi_field_element_with_small_characteristics::Multi_field_element_with_small_characteristics( - Multi_field_element_with_small_characteristics&& toMove) noexcept - : element_(std::exchange(toMove.element_, 0)) {} - -template -inline Multi_field_element_with_small_characteristics& -Multi_field_element_with_small_characteristics::operator=( - Multi_field_element_with_small_characteristics other) { - std::swap(element_, other.element_); - return *this; -} - -template -inline Multi_field_element_with_small_characteristics& -Multi_field_element_with_small_characteristics::operator=(unsigned int const value) { - element_ = value % productOfAllCharacteristics_; - return *this; -} - -template -inline Multi_field_element_with_small_characteristics::operator unsigned int() const { - return element_; -} - -template -inline Multi_field_element_with_small_characteristics -Multi_field_element_with_small_characteristics::get_inverse() const { - return get_partial_inverse(productOfAllCharacteristics_).first; -} - -template -inline std::pair, unsigned int> -Multi_field_element_with_small_characteristics::get_partial_inverse( - unsigned int productOfCharacteristics) const { - unsigned int gcd = std::gcd(element_, productOfAllCharacteristics_); - - if (gcd == productOfCharacteristics) - return {Multi_field_element_with_small_characteristics(), multiplicativeID_}; // partial inverse is 0 - - unsigned int QT = productOfCharacteristics / gcd; - - const unsigned int inv_qt = _get_inverse(element_, QT); - - auto res = get_partial_multiplicative_identity(QT); - res *= inv_qt; - - return {res, QT}; -} - -template -inline Multi_field_element_with_small_characteristics -Multi_field_element_with_small_characteristics::get_additive_identity() { - return Multi_field_element_with_small_characteristics(); -} - -template -inline Multi_field_element_with_small_characteristics -Multi_field_element_with_small_characteristics::get_multiplicative_identity() { - return Multi_field_element_with_small_characteristics(multiplicativeID_); -} - -template -inline Multi_field_element_with_small_characteristics -Multi_field_element_with_small_characteristics::get_partial_multiplicative_identity( - const mpz_class& productOfCharacteristics) { - if (productOfCharacteristics == 0) { - return Multi_field_element_with_small_characteristics(multiplicativeID_); - } - Multi_field_element_with_small_characteristics mult; - for (unsigned int idx = 0; idx < primes_.size(); ++idx) { - if ((productOfCharacteristics % primes_[idx]) == 0) { - mult += partials_[idx]; - } - } - return mult; -} - -template -inline constexpr unsigned int Multi_field_element_with_small_characteristics::get_characteristic() { - return productOfAllCharacteristics_; -} - -template -inline unsigned int Multi_field_element_with_small_characteristics::get_value() const { - return element_; -} - -template -inline constexpr unsigned int Multi_field_element_with_small_characteristics::_add( - unsigned int element, unsigned int v) { - if (UINT_MAX - element < v) { - // automatic unsigned integer overflow behaviour will make it work - element += v; - element -= productOfAllCharacteristics_; - return element; - } - - element += v; - if (element >= productOfAllCharacteristics_) element -= productOfAllCharacteristics_; - - return element; -} - -template -inline constexpr unsigned int Multi_field_element_with_small_characteristics::_substract( - unsigned int element, unsigned int v) { - if (element < v) { - element += productOfAllCharacteristics_; - } - element -= v; - - return element; -} - -template -inline constexpr unsigned int Multi_field_element_with_small_characteristics::_multiply( - unsigned int a, unsigned int b) { - unsigned int res = 0; - unsigned int temp_b = 0; - - if (b < a) std::swap(a, b); - - while (a != 0) { - if (a & 1) { - /* Add b to res, modulo m, without overflow */ - if (b >= productOfAllCharacteristics_ - res) res -= productOfAllCharacteristics_; - res += b; - } - a >>= 1; - - /* Double b, modulo m */ - temp_b = b; - if (b >= productOfAllCharacteristics_ - b) temp_b -= productOfAllCharacteristics_; - b += temp_b; - } - return res; -} - -template -inline constexpr int Multi_field_element_with_small_characteristics::_get_inverse( - unsigned int element, const unsigned int mod) { - // to solve: Ax + My = 1 - int M = mod; - int A = element; - int y = 0, x = 1; - // extended euclidien division - while (A > 1) { - int quotient = A / M; - int temp = M; - - M = A % M, A = temp; - temp = y; - - y = x - quotient * y; - x = temp; - } - - if (x < 0) x += mod; - - return x; -} - -template -inline constexpr bool Multi_field_element_with_small_characteristics::_is_prime(const int p) { - if (p <= 1) return false; - if (p <= 3) return true; - if (p % 2 == 0 || p % 3 == 0) return false; - - for (long i = 5; i * i <= p; i = i + 6) - if (p % i == 0 || p % (i + 2) == 0) return false; - - return true; -} - } // namespace persistence_fields } // namespace Gudhi diff --git a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_operators.h b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_operators.h index d3b9024505..2608870f2f 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_operators.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_operators.h @@ -10,7 +10,7 @@ /** * @file Multi_field_small_operators.h - * @author Hannah Schreiber + * @author Hannah Schreiber, Clément Maria * @brief Contains the @ref Multi_field_operators_with_small_characteristics class. */ diff --git a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_shared.h b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_shared.h index d769ef4061..f3e1748fae 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_shared.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Multi_field_small_shared.h @@ -10,7 +10,7 @@ /** * @file Multi_field_small_shared.h - * @author Hannah Schreiber + * @author Hannah Schreiber, Clément Maria * @brief Contains the @ref Shared_multi_field_element_with_small_characteristics class. */ @@ -32,58 +32,99 @@ namespace persistence_fields { * @ingroup persistence_fields * * @brief Class representing an element of a multi-field, such that `productOfAllCharacteristics ^ 2` fits into - * an unsigned int. If each instanciation of the class can represent another element, they all share the same - * characteritics. That is if the characteristics are set for one, they will be set for all the others. - * The characteristics can be set before instianciating the elements with the static - * @ref Shared_multi_field_element::initialize method. + * the given @p Unsigned_integer_type template argument. If each instanciation of the class can represent another + * element, they all share the same characteritics. That is if the characteristics are set for one, they will be + * set for all the others. The characteristics can be set before instanciating the elements with the static + * @ref Shared_multi_field_element_with_small_characteristics::initialize method. + * + * @tparam Unsigned_integer_type A native unsigned integer type: unsigned int, long unsigned int, etc. + * Will be used as the field element type. */ -class Shared_multi_field_element_with_small_characteristics -{ +template > > +class Shared_multi_field_element_with_small_characteristics { public: - using element_type = unsigned int; /**< Type for the elements in the field. */ + using element_type = Unsigned_integer_type; /**< Type for the elements in the field. */ + using characteristic_type = element_type; /**< Type for the field characteristic. */ template using isInteger = std::enable_if_t >; /** * @brief Default constructor. Sets the element to 0. */ - Shared_multi_field_element_with_small_characteristics(); - /** - * @brief Constructor setting the element to the given value. - * - * @param element Value of the element. - */ - Shared_multi_field_element_with_small_characteristics(unsigned int element); + Shared_multi_field_element_with_small_characteristics() : element_(0) {} /** * @brief Constructor setting the element to the given value. - * + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if negative. * @param element Value of the element. */ - Shared_multi_field_element_with_small_characteristics(int element); + template > + Shared_multi_field_element_with_small_characteristics(Integer_type element) : element_(_get_value(element)) {} /** * @brief Copy constructor. - * + * * @param toCopy Element to copy. */ Shared_multi_field_element_with_small_characteristics( - const Shared_multi_field_element_with_small_characteristics& toCopy); + const Shared_multi_field_element_with_small_characteristics& toCopy) + : element_(toCopy.element_) {} /** * @brief Move constructor. - * + * * @param toMove Element to move. */ Shared_multi_field_element_with_small_characteristics( - Shared_multi_field_element_with_small_characteristics&& toMove) noexcept; + Shared_multi_field_element_with_small_characteristics&& toMove) noexcept + : element_(std::exchange(toMove.element_, 0)) {} /** * @brief Initialize the multi-field to the characteristics (primes) contained in the given interval. * Should be called first before constructing the field elements. * The characteristics must be small enough such that `productOfAllCharacteristics ^ 2` fits into an unsigned int. - * + * * @param minimum Lowest value in the interval. * @param maximum Highest value in the interval. */ - static void initialize(unsigned int minimum, unsigned int maximum); + static void initialize(unsigned int minimum, unsigned int maximum) { + if (maximum < 2) throw std::invalid_argument("Characteristic must be strictly positive"); + if (minimum > maximum) throw std::invalid_argument("The given interval is not valid."); + if (minimum == maximum && !_is_prime(minimum)) + throw std::invalid_argument("The given interval does not contain a prime number."); + + productOfAllCharacteristics_ = 1; + primes_.clear(); + for (unsigned int i = minimum; i <= maximum; ++i) { + if (_is_prime(i)) { + primes_.push_back(i); + productOfAllCharacteristics_ *= i; + } + } + + if (primes_.empty()) throw std::invalid_argument("The given interval does not contain a prime number."); + + partials_.resize(primes_.size()); + for (characteristic_type i = 0; i < primes_.size(); ++i) { + characteristic_type p = primes_[i]; + characteristic_type base = productOfAllCharacteristics_ / p; + characteristic_type exp = p - 1; + partials_[i] = 1; + + while (exp > 0) { + // If exp is odd, multiply with result + if (exp & 1) partials_[i] = _multiply(partials_[i], base); + // y must be even now + exp = exp >> 1; // y = y/2 + base = _multiply(base, base); + } + } + + // If I understood the paper well, multiplicativeID_ always equals to 1. But in Clement's code, + // multiplicativeID_ is computed (see commented loop below). TODO: verify with Clement. + // for (unsigned int i = 0; i < partials_.size(); ++i){ + // multiplicativeID_ = (multiplicativeID_ + partials_[i]) % productOfAllCharacteristics_; + // } + } /** * @brief operator+= @@ -103,29 +144,33 @@ class Shared_multi_field_element_with_small_characteristics } /** * @brief operator+= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - friend void operator+=(Shared_multi_field_element_with_small_characteristics& f, unsigned int const v) { - f.element_ = _add(f.element_, v < productOfAllCharacteristics_ ? v : v % productOfAllCharacteristics_); + template > + friend void operator+=(Shared_multi_field_element_with_small_characteristics& f, const Integer_type& v) { + f.element_ = _add(f.element_, _get_value(v)); } /** * @brief operator+ - * - * @warning @p v is assumed to be positive and will be casted into an unsigned int. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend Shared_multi_field_element_with_small_characteristics operator+( - Shared_multi_field_element_with_small_characteristics f, const Integer_type v) { + Shared_multi_field_element_with_small_characteristics f, const Integer_type& v) { f += v; return f; } /** * @brief operator+ + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend Integer_type operator+(Integer_type v, Shared_multi_field_element_with_small_characteristics const& f) { - v += f.element_; - v %= productOfAllCharacteristics_; - return v; + friend Integer_type operator+(const Integer_type& v, Shared_multi_field_element_with_small_characteristics f) { + f += v; + return f.element_; } /** @@ -146,32 +191,32 @@ class Shared_multi_field_element_with_small_characteristics } /** * @brief operator-= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - friend void operator-=(Shared_multi_field_element_with_small_characteristics& f, unsigned int const v) { - f.element_ = _substract(f.element_, v < productOfAllCharacteristics_ ? v : v % productOfAllCharacteristics_); + template > + friend void operator-=(Shared_multi_field_element_with_small_characteristics& f, const Integer_type& v) { + f.element_ = _substract(f.element_, _get_value(v)); } /** * @brief operator- - * - * @warning @p v is assumed to be positive and will be casted into an unsigned int. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend Shared_multi_field_element_with_small_characteristics operator-( - Shared_multi_field_element_with_small_characteristics f, const Integer_type v) { + Shared_multi_field_element_with_small_characteristics f, const Integer_type& v) { f -= v; return f; } /** * @brief operator- - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend Integer_type operator-(Integer_type v, Shared_multi_field_element_with_small_characteristics const& f) { - if (v >= productOfAllCharacteristics_) v %= productOfAllCharacteristics_; - if (f.element_ > v) v += productOfAllCharacteristics_; - v -= f.element_; - return v; + friend Integer_type operator-(const Integer_type& v, const Shared_multi_field_element_with_small_characteristics& f) { + return _substract(_get_value(v), f.element_); } /** @@ -192,45 +237,33 @@ class Shared_multi_field_element_with_small_characteristics } /** * @brief operator*= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - friend void operator*=(Shared_multi_field_element_with_small_characteristics& f, unsigned int const v) { - f.element_ = _multiply(f.element_, v < productOfAllCharacteristics_ ? v : v % productOfAllCharacteristics_); + template > + friend void operator*=(Shared_multi_field_element_with_small_characteristics& f, const Integer_type& v) { + f.element_ = _multiply(f.element_, _get_value(v)); } /** * @brief operator* - * - * @warning @p v is assumed to be positive and will be casted into an unsigned int. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend Shared_multi_field_element_with_small_characteristics operator*( - Shared_multi_field_element_with_small_characteristics f, const Integer_type v) { + Shared_multi_field_element_with_small_characteristics f, const Integer_type& v) { f *= v; return f; } /** * @brief operator* - * - * @warning Uses bitwise operations on @p v, so be carefull with signed integers. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend Integer_type operator*(Integer_type v, Shared_multi_field_element_with_small_characteristics const& f) { - unsigned int b = f.element_; - unsigned int res = 0; - unsigned int temp_b; - - while (v != 0) { - if (v & 1) { - if (b >= productOfAllCharacteristics_ - res) res -= productOfAllCharacteristics_; - res += b; - } - v >>= 1; - - temp_b = b; - if (b >= productOfAllCharacteristics_ - b) temp_b -= productOfAllCharacteristics_; - b += temp_b; - } - - return res; + friend Integer_type operator*(const Integer_type& v, Shared_multi_field_element_with_small_characteristics f) { + f *= v; + return f.element_; } /** @@ -242,23 +275,21 @@ class Shared_multi_field_element_with_small_characteristics } /** * @brief operator== - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend bool operator==(const Integer_type v, const Shared_multi_field_element_with_small_characteristics& f) { - if (v < productOfAllCharacteristics_) return v == f.element_; - return (v % productOfAllCharacteristics_) == f.element_; + friend bool operator==(const Integer_type& v, const Shared_multi_field_element_with_small_characteristics& f) { + return _get_value(v) == f.element_; } /** * @brief operator== - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > - friend bool operator==(const Shared_multi_field_element_with_small_characteristics& f, const Integer_type v) { - if (v < productOfAllCharacteristics_) return v == f.element_; - return (v % productOfAllCharacteristics_) == f.element_; + friend bool operator==(const Shared_multi_field_element_with_small_characteristics& f, const Integer_type& v) { + return _get_value(v) == f.element_; } /** * @brief operator!= @@ -269,8 +300,8 @@ class Shared_multi_field_element_with_small_characteristics } /** * @brief operator!= - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend bool operator!=(const Integer_type v, const Shared_multi_field_element_with_small_characteristics& f) { @@ -278,8 +309,8 @@ class Shared_multi_field_element_with_small_characteristics } /** * @brief operator!= - * - * @warning @p v is assumed to be positive. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ template > friend bool operator!=(const Shared_multi_field_element_with_small_characteristics& f, const Integer_type v) { @@ -290,11 +321,20 @@ class Shared_multi_field_element_with_small_characteristics * @brief Assign operator. */ Shared_multi_field_element_with_small_characteristics& operator=( - Shared_multi_field_element_with_small_characteristics other); + Shared_multi_field_element_with_small_characteristics other) { + std::swap(element_, other.element_); + return *this; + } /** * @brief Assign operator. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. */ - Shared_multi_field_element_with_small_characteristics& operator=(const unsigned int value); + template > + Shared_multi_field_element_with_small_characteristics& operator=(const Integer_type& value) { + element_ = _get_value(value); + return *this; + } /** * @brief Swap operator. */ @@ -306,328 +346,184 @@ class Shared_multi_field_element_with_small_characteristics /** * @brief Casts the element into an unsigned int. */ - operator unsigned int() const; + operator unsigned int() const { return element_; } /** * @brief Returns the inverse of the element in the multi-field, see @cite boissonnat:hal-00922572. - * + * * @return The inverse. */ - Shared_multi_field_element_with_small_characteristics get_inverse() const; + Shared_multi_field_element_with_small_characteristics get_inverse() const { + return get_partial_inverse(productOfAllCharacteristics_).first; + } /** * @brief Returns the inverse of the element with respect to a sub-product of the characteristics in the multi-field, * see @cite boissonnat:hal-00922572. - * + * * @param productOfCharacteristics Sub-product of the characteristics. * @return Pair of the inverse and the characteristic the inverse corresponds to. */ - std::pair get_partial_inverse( - unsigned int productOfCharacteristics) const; + std::pair get_partial_inverse( + characteristic_type productOfCharacteristics) const { + characteristic_type gcd = std::gcd(element_, productOfAllCharacteristics_); + + if (gcd == productOfCharacteristics) + return {Shared_multi_field_element_with_small_characteristics(), multiplicativeID_}; // partial inverse is 0 + + characteristic_type QT = productOfCharacteristics / gcd; + + const element_type inv_qt = _get_inverse(element_, QT); + + auto res = get_partial_multiplicative_identity(QT); + res *= inv_qt; + + return {res, QT}; + } /** * @brief Returns the additive identity of a field. - * + * * @return The additive identity of a field. */ - static Shared_multi_field_element_with_small_characteristics get_additive_identity(); + static Shared_multi_field_element_with_small_characteristics get_additive_identity() { + return Shared_multi_field_element_with_small_characteristics(); + } /** * @brief Returns the multiplicative identity of a field. - * + * * @return The multiplicative identity of a field. */ - static Shared_multi_field_element_with_small_characteristics get_multiplicative_identity(); + static Shared_multi_field_element_with_small_characteristics get_multiplicative_identity() { + return Shared_multi_field_element_with_small_characteristics(multiplicativeID_); + } /** * @brief Returns the partial multiplicative identity of the multi-field from the given product. * See @cite boissonnat:hal-00922572 for more details. - * + * * @param productOfCharacteristics Product of the different characteristics to take into account in the multi-field. * @return The partial multiplicative identity of the multi-field. */ static Shared_multi_field_element_with_small_characteristics get_partial_multiplicative_identity( - const mpz_class& productOfCharacteristics); + const mpz_class& productOfCharacteristics) { + if (productOfCharacteristics == 0) { + return Shared_multi_field_element_with_small_characteristics(multiplicativeID_); + } + Shared_multi_field_element_with_small_characteristics mult; + for (characteristic_type idx = 0; idx < primes_.size(); ++idx) { + if ((productOfCharacteristics % primes_[idx]) == 0) { + mult += partials_[idx]; + } + } + return mult; + } /** * @brief Returns the product of all characteristics. - * + * * @return The product of all characteristics. */ - static unsigned int get_characteristic(); + static characteristic_type get_characteristic() { return productOfAllCharacteristics_; } /** * @brief Returns the value of the element. - * + * * @return Value of the element. */ - unsigned int get_value() const; + element_type get_value() const { return element_; } // static constexpr bool handles_only_z2() { return false; } private: - static constexpr bool _is_prime(const int p); - static unsigned int _multiply(unsigned int a, unsigned int b); - static unsigned int _add(unsigned int element, unsigned int v); - static unsigned int _substract(unsigned int element, unsigned int v); - static constexpr int _get_inverse(unsigned int element, const unsigned int mod); - - unsigned int element_; /**< Element. */ - static inline std::vector primes_; /**< All characteristics. */ - static inline unsigned int productOfAllCharacteristics_; /**< Product of all characteristics. */ - static inline std::vector partials_; /**< Partial products of the characteristics. */ - static inline constexpr unsigned int multiplicativeID_ = 1; /**< Multiplicative identity. */ -}; + static constexpr bool _is_prime(const unsigned int p) { + if (p <= 1) return false; + if (p <= 3) return true; + if (p % 2 == 0 || p % 3 == 0) return false; -inline Shared_multi_field_element_with_small_characteristics::Shared_multi_field_element_with_small_characteristics() - : element_(0) {} - -inline Shared_multi_field_element_with_small_characteristics::Shared_multi_field_element_with_small_characteristics( - unsigned int element) - : element_(element % productOfAllCharacteristics_) {} - -inline Shared_multi_field_element_with_small_characteristics::Shared_multi_field_element_with_small_characteristics( - int element) - : element_(element % productOfAllCharacteristics_) {} - -inline Shared_multi_field_element_with_small_characteristics::Shared_multi_field_element_with_small_characteristics( - const Shared_multi_field_element_with_small_characteristics& toCopy) - : element_(toCopy.element_) {} - -inline Shared_multi_field_element_with_small_characteristics::Shared_multi_field_element_with_small_characteristics( - Shared_multi_field_element_with_small_characteristics&& toMove) noexcept - : element_(std::exchange(toMove.element_, 0)) {} - -inline void Shared_multi_field_element_with_small_characteristics::initialize(unsigned int minimum, - unsigned int maximum) { - if (maximum < 2) throw std::invalid_argument("Characteristic must be strictly positive"); - if (minimum > maximum) throw std::invalid_argument("The given interval is not valid."); - if (minimum == maximum && !_is_prime(minimum)) - throw std::invalid_argument("The given interval does not contain a prime number."); - - productOfAllCharacteristics_ = 1; - primes_.clear(); - for (unsigned int i = minimum; i <= maximum; ++i) { - if (_is_prime(i)) { - primes_.push_back(i); - productOfAllCharacteristics_ *= i; - } + for (long i = 5; i * i <= p; i = i + 6) + if (p % i == 0 || p % (i + 2) == 0) return false; + + return true; } + static element_type _multiply(element_type a, element_type b) { + element_type res = 0; + element_type temp_b = 0; - if (primes_.empty()) throw std::invalid_argument("The given interval does not contain a prime number."); + if (b < a) std::swap(a, b); - partials_.resize(primes_.size()); - for (unsigned int i = 0; i < primes_.size(); ++i) { - unsigned int p = primes_[i]; - unsigned int base = productOfAllCharacteristics_ / p; - unsigned int exp = p - 1; - partials_[i] = 1; + while (a != 0) { + if (a & 1) { + /* Add b to res, modulo m, without overflow */ + if (b >= productOfAllCharacteristics_ - res) res -= productOfAllCharacteristics_; + res += b; + } + a >>= 1; - while (exp > 0) { - // If exp is odd, multiply with result - if (exp & 1) partials_[i] = _multiply(partials_[i], base); - // y must be even now - exp = exp >> 1; // y = y/2 - base = _multiply(base, base); + /* Double b, modulo m */ + temp_b = b; + if (b >= productOfAllCharacteristics_ - b) temp_b -= productOfAllCharacteristics_; + b += temp_b; } + return res; } - - // If I understood the paper well, multiplicativeID_ always equals to 1. But in Clement's code, - // multiplicativeID_ is computed (see commented loop below). TODO: verify with Clement. - // for (unsigned int i = 0; i < partials_.size(); ++i){ - // multiplicativeID_ = (multiplicativeID_ + partials_[i]) % productOfAllCharacteristics_; - // } -} - -// inline Shared_multi_field_element_with_small_characteristics -// &Shared_multi_field_element_with_small_characteristics::operator+=(Shared_multi_field_element_with_small_characteristics -// const &f) -//{ -// _add(f.element_); -// return *this; -// } - -// inline Shared_multi_field_element_with_small_characteristics -// &Shared_multi_field_element_with_small_characteristics::operator+=(unsigned int const v) -//{ -// _add(v % productOfAllCharacteristics_); -// return *this; -// } - -// inline Shared_multi_field_element_with_small_characteristics -// &Shared_multi_field_element_with_small_characteristics::operator-=(Shared_multi_field_element_with_small_characteristics -// const &f) -//{ -// _substract(f.element_); -// return *this; -// } - -// inline Shared_multi_field_element_with_small_characteristics -// &Shared_multi_field_element_with_small_characteristics::operator-=(unsigned int const v) -//{ -// _substract(v % productOfAllCharacteristics_); -// return *this; -// } - -// inline Shared_multi_field_element_with_small_characteristics -// &Shared_multi_field_element_with_small_characteristics::operator*=(Shared_multi_field_element_with_small_characteristics -// const &f) -//{ -// element_ = _multiply(element_, f.element_); -// return *this; -// } - -// inline Shared_multi_field_element_with_small_characteristics -// &Shared_multi_field_element_with_small_characteristics::operator*=(unsigned int const v) -//{ -// element_ = _multiply(element_, v % productOfAllCharacteristics_); -// return *this; -// } - -inline Shared_multi_field_element_with_small_characteristics& -Shared_multi_field_element_with_small_characteristics::operator=( - Shared_multi_field_element_with_small_characteristics other) { - std::swap(element_, other.element_); - return *this; -} - -inline Shared_multi_field_element_with_small_characteristics& -Shared_multi_field_element_with_small_characteristics::operator=(unsigned int const value) { - element_ = value % productOfAllCharacteristics_; - return *this; -} - -inline Shared_multi_field_element_with_small_characteristics::operator unsigned int() const { return element_; } - -inline Shared_multi_field_element_with_small_characteristics -Shared_multi_field_element_with_small_characteristics::get_inverse() const { - return get_partial_inverse(productOfAllCharacteristics_).first; -} - -inline std::pair -Shared_multi_field_element_with_small_characteristics::get_partial_inverse( - unsigned int productOfCharacteristics) const { - unsigned int gcd = std::gcd(element_, productOfAllCharacteristics_); - - if (gcd == productOfCharacteristics) - return {Shared_multi_field_element_with_small_characteristics(), multiplicativeID_}; // partial inverse is 0 - - unsigned int QT = productOfCharacteristics / gcd; - - const unsigned int inv_qt = _get_inverse(element_, QT); - - auto res = get_partial_multiplicative_identity(QT); - res *= inv_qt; - - return {res, QT}; -} - -inline Shared_multi_field_element_with_small_characteristics -Shared_multi_field_element_with_small_characteristics::get_additive_identity() { - return Shared_multi_field_element_with_small_characteristics(); -} - -inline Shared_multi_field_element_with_small_characteristics -Shared_multi_field_element_with_small_characteristics::get_multiplicative_identity() { - return Shared_multi_field_element_with_small_characteristics(multiplicativeID_); -} - -inline Shared_multi_field_element_with_small_characteristics -Shared_multi_field_element_with_small_characteristics::get_partial_multiplicative_identity( - const mpz_class& productOfCharacteristics) { - if (productOfCharacteristics == 0) { - return Shared_multi_field_element_with_small_characteristics(multiplicativeID_); - } - Shared_multi_field_element_with_small_characteristics mult; - for (unsigned int idx = 0; idx < primes_.size(); ++idx) { - if ((productOfCharacteristics % primes_[idx]) == 0) { - mult += partials_[idx]; + static element_type _add(element_type element, element_type v) { + if (UINT_MAX - element < v) { + // automatic unsigned integer overflow behaviour will make it work + element += v; + element -= productOfAllCharacteristics_; + return element; } - } - return mult; -} - -inline unsigned int Shared_multi_field_element_with_small_characteristics::get_characteristic() { - return productOfAllCharacteristics_; -} - -inline unsigned int Shared_multi_field_element_with_small_characteristics::get_value() const { return element_; } -inline unsigned int Shared_multi_field_element_with_small_characteristics::_add(unsigned int element, unsigned int v) { - if (UINT_MAX - element < v) { - // automatic unsigned integer overflow behaviour will make it work element += v; - element -= productOfAllCharacteristics_; - return element; - } - - element += v; - if (element >= productOfAllCharacteristics_) element -= productOfAllCharacteristics_; - - return element; -} + if (element >= productOfAllCharacteristics_) element -= productOfAllCharacteristics_; -inline unsigned int Shared_multi_field_element_with_small_characteristics::_substract(unsigned int element, - unsigned int v) { - if (element < v) { - element += productOfAllCharacteristics_; + return element; } - element -= v; - - return element; -} - -inline unsigned int Shared_multi_field_element_with_small_characteristics::_multiply(unsigned int a, unsigned int b) { - unsigned int res = 0; - unsigned int temp_b = 0; - - if (b < a) std::swap(a, b); - - while (a != 0) { - if (a & 1) { - /* Add b to res, modulo m, without overflow */ - if (b >= productOfAllCharacteristics_ - res) res -= productOfAllCharacteristics_; - res += b; + static element_type _substract(element_type element, element_type v) { + if (element < v) { + element += productOfAllCharacteristics_; } - a >>= 1; + element -= v; - /* Double b, modulo m */ - temp_b = b; - if (b >= productOfAllCharacteristics_ - b) temp_b -= productOfAllCharacteristics_; - b += temp_b; + return element; } - return res; -} - -inline constexpr int Shared_multi_field_element_with_small_characteristics::_get_inverse(unsigned int element, - const unsigned int mod) { - // to solve: Ax + My = 1 - int M = mod; - int A = element; - int y = 0, x = 1; - // extended euclidien division - while (A > 1) { - int quotient = A / M; - int temp = M; + static constexpr int _get_inverse(element_type element, const characteristic_type mod) { + // to solve: Ax + My = 1 + int M = mod; + int A = element; + int y = 0, x = 1; + // extended euclidien division + while (A > 1) { + int quotient = A / M; + int temp = M; + + M = A % M, A = temp; + temp = y; + + y = x - quotient * y; + x = temp; + } - M = A % M, A = temp; - temp = y; + if (x < 0) x += mod; - y = x - quotient * y; - x = temp; + return x; } - if (x < 0) x += mod; - - return x; -} - -inline constexpr bool Shared_multi_field_element_with_small_characteristics::_is_prime(const int p) { - if (p <= 1) return false; - if (p <= 3) return true; - if (p % 2 == 0 || p % 3 == 0) return false; - - for (long i = 5; i * i <= p; i = i + 6) - if (p % i == 0 || p % (i + 2) == 0) return false; + template > + static constexpr element_type _get_value(Integer_type e) { + if constexpr (std::is_signed_v){ + if (e < -static_cast(productOfAllCharacteristics_)) e = e % productOfAllCharacteristics_; + if (e < 0) return e += productOfAllCharacteristics_; + return e < static_cast(productOfAllCharacteristics_) ? e : e % productOfAllCharacteristics_; + } else { + return e < productOfAllCharacteristics_ ? e : e % productOfAllCharacteristics_; + } + } - return true; -} + element_type element_; /**< Element. */ + static inline std::vector primes_; /**< All characteristics. */ + static inline characteristic_type productOfAllCharacteristics_; /**< Product of all characteristics. */ + static inline std::vector partials_; /**< Partial products of the characteristics. */ + static inline constexpr element_type multiplicativeID_ = 1; /**< Multiplicative identity. */ +}; } // namespace persistence_fields } // namespace Gudhi diff --git a/src/Persistence_matrix/include/gudhi/Fields/Z2_field.h b/src/Persistence_matrix/include/gudhi/Fields/Z2_field.h index 40865cf36f..4877fba31e 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Z2_field.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Z2_field.h @@ -2,7 +2,7 @@ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details. * Author(s): Hannah Schreiber * - * Copyright (C) 2022 Inria + * Copyright (C) 2022-24 Inria * * Modification(s): * - YYYY/MM Author: Description of the modification @@ -31,7 +31,7 @@ namespace persistence_fields { class Z2_field_element { public: - using element_type = unsigned int; /**< Type for the elements in the field. */ + using element_type = bool; /**< Type for the elements in the field. */ template using isInteger = std::enable_if_t >; @@ -42,15 +42,11 @@ class Z2_field_element /** * @brief Constructor setting the element to the given value. * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic. * @param element Value of the element. */ - Z2_field_element(unsigned int element); - /** - * @brief Constructor setting the element to the given value. - * - * @param element Value of the element. - */ - Z2_field_element(int element); + template > + Z2_field_element(Integer_type element); /** * @brief Copy constructor. * @@ -79,17 +75,18 @@ class Z2_field_element } /** * @brief operator+= + * + * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ - friend void operator+=(Z2_field_element& f, unsigned int const v) { f.element_ = (f.element_ != (v % 2)); } + template > + friend void operator+=(Z2_field_element& f, const Integer_type& v) { f.element_ = (f.element_ != _get_value(v)); } /** * @brief operator+ * - * @warning @p v will be casted into an unsigned int. - * * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ template > - friend Z2_field_element operator+(Z2_field_element f, const Integer_type v) { + friend Z2_field_element operator+(Z2_field_element f, const Integer_type& v) { f += v; return f; } @@ -99,8 +96,8 @@ class Z2_field_element * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ template > - friend Integer_type operator+(const Integer_type v, Z2_field_element const& f) { - return f.element_ != (v % 2); + friend Integer_type operator+(const Integer_type& v, const Z2_field_element& f) { + return f.element_ != _get_value(v); } /** @@ -118,17 +115,18 @@ class Z2_field_element } /** * @brief operator-= + * + * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ - friend void operator-=(Z2_field_element& f, unsigned int const v) { f.element_ = (f.element_ != (v % 2)); } + template > + friend void operator-=(Z2_field_element& f, const Integer_type& v) { f.element_ = (f.element_ != _get_value(v)); } /** * @brief operator- * - * @warning @p v will be casted into an unsigned int. - * * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ template > - friend Z2_field_element operator-(Z2_field_element f, const Integer_type v) { + friend Z2_field_element operator-(Z2_field_element f, const Integer_type& v) { f -= v; return f; } @@ -139,7 +137,7 @@ class Z2_field_element */ template > friend Integer_type operator-(const Integer_type v, Z2_field_element const& f) { - return f.element_ != (v % 2); + return f.element_ != _get_value(v); } /** @@ -157,18 +155,18 @@ class Z2_field_element } /** * @brief operator*= + * + * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ - friend void operator*=(Z2_field_element& f, unsigned int const v) { f.element_ = (f.element_ && (v % 2)); } - // v will be casted into an unsigned int + template > + friend void operator*=(Z2_field_element& f, const Integer_type& v) { f.element_ = (f.element_ && _get_value(v)); } /** * @brief operator* * - * @warning @p v will be casted into an unsigned int. - * * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ template > - friend Z2_field_element operator*(Z2_field_element f, const Integer_type v) { + friend Z2_field_element operator*(Z2_field_element f, const Integer_type& v) { f *= v; return f; } @@ -178,8 +176,8 @@ class Z2_field_element * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ template > - friend Integer_type operator*(const Integer_type v, Z2_field_element const& f) { - return f.element_ && (v % 2); + friend Integer_type operator*(const Integer_type& v, Z2_field_element const& f) { + return f.element_ && _get_value(v); } /** @@ -192,8 +190,8 @@ class Z2_field_element * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ template > - friend bool operator==(const Integer_type v, const Z2_field_element& f) { - return (v % 2) == f.element_; + friend bool operator==(const Integer_type& v, const Z2_field_element& f) { + return _get_value(v) == f.element_; } /** * @brief operator== @@ -201,8 +199,8 @@ class Z2_field_element * @tparam Integer_type A native integer type: int, unsigned int, long int, bool, etc. */ template > - friend bool operator==(const Z2_field_element& f, const Integer_type v) { - return (v % 2) == f.element_; + friend bool operator==(const Z2_field_element& f, const Integer_type& v) { + return _get_value(v) == f.element_; } /** * @brief operator!= @@ -290,19 +288,27 @@ class Z2_field_element * * @return Value of the element. */ - unsigned int get_value() const; + element_type get_value() const; // static constexpr bool handles_only_z2() { return true; } private: - bool element_; + element_type element_; + + template > + static constexpr element_type _get_value(Integer_type e) { + if constexpr (std::is_same_v) { + return e; + } else { + return e < 2 && e >= 0 ? e : e % 2; // returns bool, so %2 won't be negative and is optimized to & + } + } }; inline Z2_field_element::Z2_field_element() : element_(false) {} -inline Z2_field_element::Z2_field_element(unsigned int element) : element_(element & 1U) {} - -inline Z2_field_element::Z2_field_element(int element) : element_(element & 1U) {} +template +inline Z2_field_element::Z2_field_element(Integer_type element) : element_(_get_value(element)) {} inline Z2_field_element::Z2_field_element(const Z2_field_element& toCopy) : element_(toCopy.element_) {} @@ -315,7 +321,7 @@ inline Z2_field_element& Z2_field_element::operator=(Z2_field_element other) { } inline Z2_field_element& Z2_field_element::operator=(unsigned int const value) { - element_ = value & 1U; + element_ = _get_value(value); return *this; } @@ -341,7 +347,7 @@ inline Z2_field_element Z2_field_element::get_partial_multiplicative_identity( inline constexpr unsigned int Z2_field_element::get_characteristic() { return 2; } -inline unsigned int Z2_field_element::get_value() const { return element_; } +inline Z2_field_element::element_type Z2_field_element::get_value() const { return element_; } } // namespace persistence_fields } // namespace Gudhi diff --git a/src/Persistence_matrix/include/gudhi/Fields/Zp_field.h b/src/Persistence_matrix/include/gudhi/Fields/Zp_field.h index 1fa6385b01..56fb4cf4f3 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Zp_field.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Zp_field.h @@ -2,12 +2,18 @@ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details. * Author(s): Hannah Schreiber * - * Copyright (C) 2022 Inria + * Copyright (C) 2022-24 Inria * * Modification(s): * - YYYY/MM Author: Description of the modification */ +/** + * @file Zp_field.h + * @author Hannah Schreiber + * @brief Contains the @ref Zp_field_element class. + */ + #ifndef MATRIX_FIELD_ZP_H_ #define MATRIX_FIELD_ZP_H_ @@ -19,366 +25,396 @@ namespace Gudhi { namespace persistence_fields { /** + * @class Zp_field_element Zp_field.h gudhi/Fields/Zp_field.h * @ingroup persistence_fields * - * @brief TODO: - * + * @brief Class representing an element of the \f$ \mathbb{F}_p \f$ field for any prime number \f$ p \f$. + * + * @tparam characteristic Value of the characteristic of the field. Has to be a positive prime number. + * @tparam Unsigned_integer_type A native unsigned integer type: unsigned int, long unsigned int, etc. + * Will be used as the field element type. */ -template +template > > class Zp_field_element { -public: - using element_type = unsigned int; - template - using isInteger = std::enable_if_t >; - - Zp_field_element(); - Zp_field_element(unsigned int element); - Zp_field_element(int element); - Zp_field_element(const Zp_field_element& toCopy); - Zp_field_element(Zp_field_element&& toMove) noexcept; - - friend void operator+=(Zp_field_element& f1, Zp_field_element const& f2){ - f1.element_ = Zp_field_element::_add(f1.element_, f2.element_); - } - friend Zp_field_element operator+(Zp_field_element f1, Zp_field_element const& f2){ - f1 += f2; - return f1; - } - friend void operator+=(Zp_field_element& f, const unsigned int v){ - f.element_ = Zp_field_element::_add(f.element_, v < characteristic ? v : v % characteristic); - } - //v is assumed to be positive and will be casted into an unsigned int - template > - friend Zp_field_element operator+(Zp_field_element f, const Integer_type v){ - f += v; - return f; - } - template > - friend Integer_type operator+(Integer_type v, Zp_field_element const& f){ - v += f.element_; - if (v >= characteristic) v %= characteristic; - return v; - } - - friend void operator-=(Zp_field_element& f1, Zp_field_element const& f2){ - f1.element_ = Zp_field_element::_substract(f1.element_, f2.element_); - } - friend Zp_field_element operator-(Zp_field_element f1, Zp_field_element const& f2){ - f1 -= f2; - return f1; - } - friend void operator-=(Zp_field_element& f, unsigned int const v){ - f.element_ = Zp_field_element::_substract(f.element_, v < characteristic ? v : v % characteristic); - } - //v is assumed to be positive and will be casted into an unsigned int - template > - friend Zp_field_element operator-(Zp_field_element f, const Integer_type v){ - f -= v; - return f; - } - //v is assumed to be positive - template > - friend unsigned int operator-(Integer_type v, Zp_field_element const& f){ - if (v >= characteristic) v %= characteristic; - if (f.element_ > v) v += characteristic; - v -= f.element_; - return v; - } - - friend void operator*=(Zp_field_element& f1, Zp_field_element const& f2){ - f1.element_ = Zp_field_element::_multiply(f1.element_, f2.element_); - } - friend Zp_field_element operator*(Zp_field_element f1, Zp_field_element const& f2){ - f1 *= f2; - return f1; - } - friend void operator*=(Zp_field_element& f, unsigned int const v){ - f.element_ = Zp_field_element::_multiply(f.element_, v < characteristic ? v : v % characteristic); - } - //v is assumed to be positive and will be casted into an unsigned int - template > - friend Zp_field_element operator*(Zp_field_element f, const Integer_type v){ - f *= v; - return f; - } - //uses bitwise operations on v, so be carefull with signed integers - template > - friend unsigned int operator*(Integer_type v, Zp_field_element const& f){ - unsigned int b = f.element_; - unsigned int res = 0; - unsigned int temp_b; - - while (v != 0) { - if (v & 1) { - if (b >= characteristic - res) - res -= characteristic; - res += b; - } - v >>= 1; - - temp_b = b; - if (b >= characteristic - b) - temp_b -= characteristic; - b += temp_b; - } - - return res; - } - - friend bool operator==(const Zp_field_element& f1, const Zp_field_element& f2){ - return f1.element_ == f2.element_; - } - //v is assumed to be positive - template > - friend bool operator==(const Integer_type v, const Zp_field_element& f){ - if (v < characteristic) return v == f.element_; - return (v % characteristic) == f.element_; - } - //v is assumed to be positive - template > - friend bool operator==(const Zp_field_element& f, const Integer_type v){ - if (v < characteristic) return v == f.element_; - return (v % characteristic) == f.element_; - } - friend bool operator!=(const Zp_field_element& f1, const Zp_field_element& f2){ - return !(f1 == f2); - } - //v is assumed to be positive - template > - friend bool operator!=(const Integer_type v, const Zp_field_element& f){ - return !(v == f); - } - //v is assumed to be positive - template > - friend bool operator!=(const Zp_field_element& f, const Integer_type v){ - return !(v == f); - } - - Zp_field_element& operator=(Zp_field_element other); - Zp_field_element& operator=(const unsigned int value); - operator unsigned int() const; - friend void swap(Zp_field_element& f1, Zp_field_element& f2){ - std::swap(f1.element_, f2.element_); - } - - Zp_field_element get_inverse() const; - std::pair get_partial_inverse(unsigned int productOfCharacteristics) const; - - static Zp_field_element get_additive_identity(); - static Zp_field_element get_multiplicative_identity(); - static Zp_field_element get_partial_multiplicative_identity([[maybe_unused]] unsigned int productOfCharacteristics); - static constexpr unsigned int get_characteristic(); - - unsigned int get_value() const; - - static constexpr bool handles_only_z2(){ - return false; - } - -private: - unsigned int element_; - static inline std::array inverse_; - - static unsigned int _add(unsigned int element, unsigned int v); - static unsigned int _substract(unsigned int element, unsigned int v); - static unsigned int _multiply(unsigned int element, unsigned int v); - static int _get_inverse(unsigned int element); - - static constexpr bool _is_prime(); + public: + using element_type = Unsigned_integer_type; /**< Type for the elements in the field. */ + using characteristic_type = element_type; /**< Type for the field characteristic. */ + template + using isInteger = std::enable_if_t >; + + /** + * @brief Default constructor. Sets the element to 0. + */ + Zp_field_element() : element_(0) { static_assert(_is_prime(), "Characteristic has to be a prime number."); } + /** + * @brief Constructor setting the element to the given value. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + * @param element Value of the element. + */ + template > + Zp_field_element(Integer_type element) : element_(_get_value(element)) { + static_assert(_is_prime(), "Characteristic has to be a prime number."); + } + /** + * @brief Copy constructor. + * + * @param toCopy Element to copy. + */ + Zp_field_element(const Zp_field_element& toCopy) : element_(toCopy.element_) {} + /** + * @brief Move constructor. + * + * @param toMove Element to move. + */ + Zp_field_element(Zp_field_element&& toMove) noexcept : element_(std::exchange(toMove.element_, 0)) {} + + /** + * @brief operator+= + */ + friend void operator+=(Zp_field_element& f1, const Zp_field_element& f2) { + f1.element_ = Zp_field_element::_add(f1.element_, f2.element_); + } + /** + * @brief operator+ + */ + friend Zp_field_element operator+(Zp_field_element f1, const Zp_field_element& f2) { + f1 += f2; + return f1; + } + /** + * @brief operator+= + */ + template > + friend void operator+=(Zp_field_element& f, const Integer_type& v) { + f.element_ = Zp_field_element::_add(f.element_, _get_value(v)); + } + /** + * @brief operator+ + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Zp_field_element operator+(Zp_field_element f, const Integer_type& v) { + f += v; + return f; + } + /** + * @brief operator+ + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Integer_type operator+(const Integer_type& v, Zp_field_element f) { + f += v; + return f.element_; + } + + /** + * @brief operator-= + */ + friend void operator-=(Zp_field_element& f1, const Zp_field_element& f2) { + f1.element_ = Zp_field_element::_substract(f1.element_, f2.element_); + } + /** + * @brief operator- + */ + friend Zp_field_element operator-(Zp_field_element f1, const Zp_field_element& f2) { + f1 -= f2; + return f1; + } + /** + * @brief operator-= + */ + template > + friend void operator-=(Zp_field_element& f, const Integer_type& v) { + f.element_ = Zp_field_element::_substract(f.element_, _get_value(v)); + } + /** + * @brief operator- + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Zp_field_element operator-(Zp_field_element f, const Integer_type& v) { + f -= v; + return f; + } + /** + * @brief operator- + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Integer_type operator-(const Integer_type& v, const Zp_field_element& f) { + return Zp_field_element::_substract(_get_value(v), f.element_); + } + + /** + * @brief operator*= + */ + friend void operator*=(Zp_field_element& f1, const Zp_field_element& f2) { + f1.element_ = Zp_field_element::_multiply(f1.element_, f2.element_); + } + /** + * @brief operator* + */ + friend Zp_field_element operator*(Zp_field_element f1, const Zp_field_element& f2) { + f1 *= f2; + return f1; + } + /** + * @brief operator*= + */ + template > + friend void operator*=(Zp_field_element& f, const Integer_type& v) { + f.element_ = Zp_field_element::_multiply(f.element_, _get_value(v)); + } + /** + * @brief operator* + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Zp_field_element operator*(Zp_field_element f, const Integer_type& v) { + f *= v; + return f; + } + /** + * @brief operator* + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Integer_type operator*(const Integer_type& v, Zp_field_element f) { + f *= v; + return f.element_; + } + + /** + * @brief operator== + */ + friend bool operator==(const Zp_field_element& f1, const Zp_field_element& f2) { + return f1.element_ == f2.element_; + } + /** + * @brief operator== + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator==(const Integer_type& v, const Zp_field_element& f) { + return Zp_field_element::_get_value(v) == f.element_; + } + /** + * @brief operator== + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator==(const Zp_field_element& f, const Integer_type& v) { + return Zp_field_element::_get_value(v) == f.element_; + } + /** + * @brief operator!= + */ + friend bool operator!=(const Zp_field_element& f1, const Zp_field_element& f2) { return !(f1 == f2); } + /** + * @brief operator!= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator!=(const Integer_type& v, const Zp_field_element& f) { + return !(v == f); + } + /** + * @brief operator!= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator!=(const Zp_field_element& f, const Integer_type& v) { + return !(v == f); + } + + /** + * @brief Assign operator. + */ + Zp_field_element& operator=(Zp_field_element other) { + std::swap(element_, other.element_); + return *this; + } + /** + * @brief Assign operator. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + Zp_field_element& operator=(const Integer_type& value) { + element_ = _get_value(value); + return *this; + } + /** + * @brief Swap operator. + */ + friend void swap(Zp_field_element& f1, Zp_field_element& f2) { std::swap(f1.element_, f2.element_); } + + /** + * @brief Casts the element into an unsigned int. + */ + operator unsigned int() const { return element_; } + + /** + * @brief Returns the inverse of the element in the field. + * + * @return The inverse. + */ + Zp_field_element get_inverse() const { + if (element_ != 0 && inverse_[element_] == 0) { // initialize everything at instanciation instead? + inverse_[element_] = _get_inverse(element_); + } + + return Zp_field_element(inverse_[element_]); + } + /** + * @brief For interface purposes with multi-fields. Returns the inverse together with the argument. + * + * @param productOfCharacteristics Some value. + * @return Pair whose first element is the inverse and the second element is @p productOfCharacteristics. + */ + std::pair get_partial_inverse(unsigned int productOfCharacteristics) const { + return {get_inverse(), productOfCharacteristics}; + } + + /** + * @brief Returns the additive identity of the field. + * + * @return 0. + */ + static Zp_field_element get_additive_identity() { return Zp_field_element(); } + /** + * @brief Returns the multiplicative identity of the field. + * + * @return 1. + */ + static Zp_field_element get_multiplicative_identity() { return Zp_field_element(1); } + /** + * @brief For interface purposes with multi-fields. Returns the multiplicative identity of the field. + * + * @param productOfCharacteristics Some value. + * @return 1. + */ + static Zp_field_element get_partial_multiplicative_identity([[maybe_unused]] unsigned int productOfCharacteristics) { + return Zp_field_element(1); + } + /** + * @brief Returns the current characteristic. + * + * @return The value of the current characteristic. + */ + static constexpr unsigned int get_characteristic() { return characteristic; } + + /** + * @brief Returns the value of the element. + * + * @return Value of the element. + */ + element_type get_value() const { return element_; } + + // static constexpr bool handles_only_z2() { return false; } + + private: + element_type element_; /**< Field element. */ + static inline std::array inverse_; /**< All inverse elements. */ + + static element_type _add(element_type element, element_type v) { + if (UINT_MAX - element < v) { + // automatic unsigned integer overflow behaviour will make it work + element += v; + element -= characteristic; + return element; + } + + element += v; + if (element >= characteristic) element -= characteristic; + + return element; + } + static element_type _substract(element_type element, element_type v) { + if (element < v) { + element += characteristic; + } + element -= v; + + return element; + } + static element_type _multiply(element_type element, element_type v) { + element_type a = element; + element = 0; + element_type temp_b; + + while (a != 0) { + if (a & 1) { + if (v >= characteristic - element) element -= characteristic; + element += v; + } + a >>= 1; + + temp_b = v; + if (v >= characteristic - v) temp_b -= characteristic; + v += temp_b; + } + + return element; + } + static int _get_inverse(element_type element) { + // to solve: Ax + My = 1 + int M = characteristic; + int A = element; + int y = 0, x = 1; + // extended euclidien division + while (A > 1) { + int quotient = A / M; + int temp = M; + + M = A % M, A = temp; + temp = y; + + y = x - quotient * y; + x = temp; + } + + if (x < 0) x += characteristic; + + return x; + } + + template > + static constexpr element_type _get_value(Integer_type e) { + if constexpr (std::is_signed_v) { + if (e < -static_cast(characteristic)) e = e % characteristic; + if (e < 0) return e += characteristic; + return e < static_cast(characteristic) ? e : e % characteristic; + } else { + return e < characteristic ? e : e % characteristic; + } + } + + static constexpr bool _is_prime() { + if (characteristic <= 1) return false; + if (characteristic <= 3) return true; + if (characteristic % 2 == 0 || characteristic % 3 == 0) return false; + + for (long i = 5; i * i <= characteristic; i = i + 6) + if (characteristic % i == 0 || characteristic % (i + 2) == 0) return false; + + return true; + } }; -template -inline Zp_field_element::Zp_field_element() - : element_(0) -{ - static_assert(_is_prime(), "Characteristic has to be a prime number."); -} - -template -inline Zp_field_element::Zp_field_element(unsigned int element) - : element_(element < characteristic ? element : element % characteristic) -{ - static_assert(_is_prime(), "Characteristic has to be a prime number."); -} - -template -inline Zp_field_element::Zp_field_element(int element) -{ - int res = element < static_cast(characteristic) ? element : element % static_cast(characteristic); - if (res < 0) res += characteristic; - element_ = res; - - static_assert(_is_prime(), "Characteristic has to be a prime number."); -} - -template -inline Zp_field_element::Zp_field_element(const Zp_field_element &toCopy) - : element_(toCopy.element_) -{} - -template -inline Zp_field_element::Zp_field_element(Zp_field_element &&toMove) noexcept - : element_(std::exchange(toMove.element_, 0)) -{} - -template -inline Zp_field_element &Zp_field_element::operator=(Zp_field_element other) -{ - std::swap(element_, other.element_); - return *this; -} - -template -inline Zp_field_element &Zp_field_element::operator=(unsigned int const value) -{ - element_ = value < characteristic ? value : value % characteristic; - return *this; -} - -template -inline Zp_field_element::operator unsigned int() const -{ - return element_; -} - -template -inline Zp_field_element Zp_field_element::get_inverse() const -{ - if (element_ != 0 && inverse_[element_] == 0) { //initialize everything at instanciation instead? - inverse_[element_] = _get_inverse(element_); - } - - return Zp_field_element(inverse_[element_]); -} - -template -inline std::pair, unsigned int> -Zp_field_element::get_partial_inverse(unsigned int productOfCharacteristics) const -{ - return {get_inverse(), productOfCharacteristics}; -} - -template -inline Zp_field_element Zp_field_element::get_additive_identity() -{ - return Zp_field_element(); -} - -template -inline Zp_field_element Zp_field_element::get_multiplicative_identity() -{ - return Zp_field_element(1); -} - -template -inline Zp_field_element Zp_field_element::get_partial_multiplicative_identity([[maybe_unused]] unsigned int productOfCharacteristics) -{ - return Zp_field_element(1); -} - -template -inline constexpr unsigned int Zp_field_element::get_characteristic() -{ - return characteristic; -} - -template -inline unsigned int Zp_field_element::get_value() const -{ - return element_; -} - -template -inline unsigned int Zp_field_element::_add(unsigned int element, unsigned int v) -{ - if (UINT_MAX - element < v) { - //automatic unsigned integer overflow behaviour will make it work - element += v; - element -= characteristic; - return element; - } - - element += v; - if (element >= characteristic) element -= characteristic; - - return element; -} - -template -inline unsigned int Zp_field_element::_substract(unsigned int element, unsigned int v) -{ - if (element < v){ - element += characteristic; - } - element -= v; - - return element; -} - -template -inline unsigned int Zp_field_element::_multiply(unsigned int element, unsigned int v) -{ - unsigned int a = element; - element = 0; - unsigned int temp_b; - - while (a != 0) { - if (a & 1) { - if (v >= characteristic - element) - element -= characteristic; - element += v; - } - a >>= 1; - - temp_b = v; - if (v >= characteristic - v) - temp_b -= characteristic; - v += temp_b; - } - - return element; -} - -template -inline int Zp_field_element::_get_inverse(unsigned int element) -{ - //to solve: Ax + My = 1 - int M = characteristic; - int A = element; - int y = 0, x = 1; - //extended euclidien division - while (A > 1) { - int quotient = A / M; - int temp = M; - - M = A % M, A = temp; - temp = y; - - y = x - quotient * y; - x = temp; - } - - if (x < 0) - x += characteristic; - - return x; -} - -template -inline constexpr bool Zp_field_element::_is_prime() -{ - if (characteristic <= 1) return false; - if (characteristic <= 3) return true; - if (characteristic % 2 == 0 || characteristic % 3 == 0) return false; - - for (long i = 5; i * i <= characteristic; i = i + 6) - if (characteristic % i == 0 || characteristic % (i + 2) == 0) - return false; - - return true; -} - -} //namespace persistence_fields -} //namespace Gudhi +} // namespace persistence_fields +} // namespace Gudhi #endif // MATRIX_FIELD_ZP_H_ diff --git a/src/Persistence_matrix/include/gudhi/Fields/Zp_field_operators.h b/src/Persistence_matrix/include/gudhi/Fields/Zp_field_operators.h index dc7a0deb0e..5fa1d7c19f 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Zp_field_operators.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Zp_field_operators.h @@ -8,6 +8,12 @@ * - YYYY/MM Author: Description of the modification */ +/** + * @file Zp_field_operators.h + * @author Hannah Schreiber + * @brief Contains the @ref Zp_field_operators class. + */ + #ifndef MATRIX_FIELD_ZP_OPERATOR_H_ #define MATRIX_FIELD_ZP_OPERATOR_H_ @@ -21,169 +27,277 @@ namespace Gudhi { namespace persistence_fields { /** + * @class Zp_field_operators Zp_field_operators.h gudhi/Fields/Zp_field_operators.h * @ingroup persistence_fields * - * @brief TODO: - * + * @brief Class defining operators for the \f$ \mathbb{F}_p \f$ field for any prime number \f$ p \f$. + * + * @tparam Unsigned_integer_type A native unsigned integer type: unsigned int, long unsigned int, etc. + * Will be used as the field element type. */ -template > > -class Zp_field_operators { -public: - using element_type = Unsigned_integer_type; - using characteristic_type = element_type; - template - using isSignedInteger = std::enable_if_t >; - - Zp_field_operators(characteristic_type characteristic = 0): characteristic_(0) { - if (characteristic != 0) set_characteristic(characteristic); - } - Zp_field_operators(const Zp_field_operators& toCopy) - : characteristic_(toCopy.characteristic_), inverse_(toCopy.inverse_) - {} - Zp_field_operators(Zp_field_operators&& toMove) noexcept - : characteristic_(std::exchange(toMove.characteristic_, 0)), inverse_(std::move(toMove.inverse_)) - {} - - void set_characteristic(characteristic_type characteristic){ - if (characteristic <= 1) - throw std::invalid_argument("Characteristic must be strictly positive and a prime number."); - - inverse_.resize(characteristic); - inverse_[0] = 0; - for (unsigned int i = 1; i < characteristic; ++i) { - unsigned int inv = 1; - unsigned int mult = inv * i; - while ((mult % characteristic) != 1) { - ++inv; - if (mult == characteristic) - throw std::invalid_argument("Characteristic must be a prime number."); - mult = inv * i; - } - inverse_[i] = inv; - } - - characteristic_ = characteristic; - } - characteristic_type get_characteristic() const{ - return characteristic_; - } - - element_type get_value(element_type e) const{ - return e < characteristic_ ? e : e % characteristic_; - } - - template > - element_type get_value(Signed_integer_type e) const{ - if (e < -static_cast(characteristic_)) e = e % characteristic_; - if (e < 0) return e += characteristic_; - return e < static_cast(characteristic_) ? e : e % characteristic_; - } - - //r = e1 + e2 - element_type add(element_type e1, element_type e2) const{ - return _add(get_value(e1), get_value(e2), characteristic_); - } - - //r = e1 - e2 - element_type substract(element_type e1, element_type e2) const{ - return _substract(get_value(e1), get_value(e2), characteristic_); - } - - //r = e1 * e2 - element_type multiply(element_type e1, element_type e2) const{ - return _multiply(get_value(e1), get_value(e2), characteristic_); - } - - //r = e * m + a - //WARNING: not overflow safe. - element_type multiply_and_add(element_type e, element_type m, element_type a) const{ - return get_value(e * m + a); - } - - //r = (e + a) * m - //WARNING: not overflow safe. - element_type add_and_multiply(element_type e, element_type a, element_type m) const{ - return get_value((e + a) * m); - } - - bool are_equal(element_type e1, element_type e2) const{ - return get_value(e1) == get_value(e2); - } - - element_type get_inverse(element_type e) const{ - return inverse_[get_value(e)]; - } - std::pair get_partial_inverse(element_type e, characteristic_type productOfCharacteristics) const{ - return {get_inverse(e), productOfCharacteristics}; - } - - static constexpr element_type get_additive_identity(){ return 0; } - static constexpr element_type get_multiplicative_identity(){ return 1; } - static constexpr element_type get_partial_multiplicative_identity([[maybe_unused]] characteristic_type productOfCharacteristics){ return 1; } - - static constexpr bool handles_only_z2(){ - return false; - } - - Zp_field_operators& operator=(Zp_field_operators other){ - std::swap(characteristic_, other.characteristic_); - inverse_.swap(other.inverse_); - return *this; - } - friend void swap(Zp_field_operators& f1, Zp_field_operators& f2){ - std::swap(f1.characteristic_, f2.characteristic_); - f1.inverse_.swap(f2.inverse_); - } - -private: - characteristic_type characteristic_; - std::vector inverse_; - - static element_type _add(element_type e1, element_type e2, characteristic_type characteristic){ - if (UINT_MAX - e1 < e2) { - //automatic unsigned integer overflow behaviour will make it work - e1 += e2; - e1 -= characteristic; - return e1; - } - - e1 += e2; - if (e1 >= characteristic) e1 -= characteristic; - - return e1; - } - static element_type _substract(element_type e1, element_type e2, characteristic_type characteristic){ - if (e1 < e2){ - e1 += characteristic; - } - e1 -= e2; - - return e1; - } - static element_type _multiply(element_type e1, element_type e2, characteristic_type characteristic){ - unsigned int a = e1; - e1 = 0; - unsigned int temp_b; - - while (a != 0) { - if (a & 1) { - if (e2 >= characteristic - e1) - e1 -= characteristic; - e1 += e2; - } - a >>= 1; - - temp_b = e2; - if (e2 >= characteristic - e2) - temp_b -= characteristic; - e2 += temp_b; - } - - return e1; - } +template > > +class Zp_field_operators +{ + public: + using element_type = Unsigned_integer_type; /**< Type for the elements in the field. */ + using characteristic_type = element_type; /**< Type for the field characteristic. */ + template + using isSignedInteger = std::enable_if_t >; + + /** + * @brief Default constructor. If a non-zero characteristic is given, initializes the field with it. + * The characteristic can later be changed again or initialized with @ref set_characteristic. + * + * @param characteristic Prime number corresponding to the desired characteristic of the field. + */ + Zp_field_operators(characteristic_type characteristic = 0) : characteristic_(0) { + if (characteristic != 0) set_characteristic(characteristic); + } + /** + * @brief Copy constructor. + * + * @param toCopy Operators to copy. + */ + Zp_field_operators(const Zp_field_operators& toCopy) + : characteristic_(toCopy.characteristic_), inverse_(toCopy.inverse_) {} + /** + * @brief Move constructor. + * + * @param toMove Operators to move. + */ + Zp_field_operators(Zp_field_operators&& toMove) noexcept + : characteristic_(std::exchange(toMove.characteristic_, 0)), inverse_(std::move(toMove.inverse_)) {} + + /** + * @brief Sets the characteristic of the field. + * + * @param characteristic Prime number corresponding to the desired characteristic of the field. + */ + void set_characteristic(characteristic_type characteristic) { + if (characteristic <= 1) + throw std::invalid_argument("Characteristic must be strictly positive and a prime number."); + + inverse_.resize(characteristic); + inverse_[0] = 0; + for (unsigned int i = 1; i < characteristic; ++i) { + unsigned int inv = 1; + unsigned int mult = inv * i; + while ((mult % characteristic) != 1) { + ++inv; + if (mult == characteristic) throw std::invalid_argument("Characteristic must be a prime number."); + mult = inv * i; + } + inverse_[i] = inv; + } + + characteristic_ = characteristic; + } + /** + * @brief Returns the current characteristic. + * + * @return The value of the current characteristic. + */ + characteristic_type get_characteristic() const { return characteristic_; } + + /** + * @brief Returns the value of an integer in the field. + * That is the positive value of the integer modulo the current characteristic. + * + * @param e Unsigned integer to return the value from. + * @return @p e modulo the current characteristic, such that the result is positive. + */ + element_type get_value(element_type e) const { return e < characteristic_ ? e : e % characteristic_; } + /** + * @brief Returns the value of an integer in the field. + * That is the positive value of the integer modulo the current characteristic. + * + * @tparam Signed_integer_type A native signed integer type: int, long int, etc. + * @param e Integer to return the value from. + * @return @p e modulo the current characteristic, such that the result is positive. + */ + template > + element_type get_value(Signed_integer_type e) const { + if (e < -static_cast(characteristic_)) e = e % characteristic_; + if (e < 0) return e += characteristic_; + return e < static_cast(characteristic_) ? e : e % characteristic_; + } + + /** + * @brief Returns the sum of two elements in the field. + * + * @param e1 First element. + * @param e2 Second element. + * @return `(e1 + e2) % characteristic`, such that the result is positive. + */ + element_type add(element_type e1, element_type e2) const { + return _add(get_value(e1), get_value(e2), characteristic_); + } + + /** + * @brief Returns the substraction in the field of the first element by the second element. + * + * @param e1 First element. + * @param e2 Second element. + * @return `(e1 - e2) % characteristic`, such that the result is positive. + */ + element_type substract(element_type e1, element_type e2) const { + return _substract(get_value(e1), get_value(e2), characteristic_); + } + + /** + * @brief Returns the multiplication of two elements in the field. + * + * @param e1 First element. + * @param e2 Second element. + * @return `(e1 * e2) % characteristic`, such that the result is positive. + */ + element_type multiply(element_type e1, element_type e2) const { + return _multiply(get_value(e1), get_value(e2), characteristic_); + } + + /** + * @brief Multiplies the first element with the second one and adds the third one. Returns the result in the field. + * + * @warning Not overflow safe. + * + * @param e First element. + * @param m Second element. + * @param a Third element. + * @return `(e * m + a) % characteristic`, such that the result is positive. + */ + element_type multiply_and_add(element_type e, element_type m, element_type a) const { return get_value(e * m + a); } + + /** + * @brief Adds the first element to the second one and multiplies the third one with it. + * Returns the result in the field. + * + * @warning Not overflow safe. + * + * @param e First element. + * @param a Second element. + * @param m Third element. + * @return `((e + a) * m) % characteristic`, such that the result is positive. + */ + element_type add_and_multiply(element_type e, element_type a, element_type m) const { return get_value((e + a) * m); } + + /** + * @brief Returns true if the two given elements are equal in the field, false otherwise. + * + * @param e1 First element to compare. + * @param e2 Second element to compare. + * @return true If `e1 % characteristic == e2 % characteristic`. + * @return false Otherwise. + */ + bool are_equal(element_type e1, element_type e2) const { return get_value(e1) == get_value(e2); } + + /** + * @brief Returns the inverse of the given element in the field. + * + * @param e Element to get the inverse from. + * @return Inverse in the current field of `e % characteristic`. + */ + element_type get_inverse(element_type e) const { return inverse_[get_value(e)]; } + /** + * @brief For interface purposes with multi-fields. Returns the inverse together with the second argument. + * + * @param e Element to get the inverse from. + * @param productOfCharacteristics Some value. + * @return Pair whose first element is the inverse of @p e and the second element is @p productOfCharacteristics. + */ + std::pair get_partial_inverse(element_type e, + characteristic_type productOfCharacteristics) const { + return {get_inverse(e), productOfCharacteristics}; + } + + /** + * @brief Returns the additive identity of the field. + * + * @return 0. + */ + static constexpr element_type get_additive_identity() { return 0; } + /** + * @brief Returns the multiplicative identity of the field. + * + * @return 1. + */ + static constexpr element_type get_multiplicative_identity() { return 1; } + /** + * @brief For interface purposes with multi-fields. Returns the multiplicative identity of the field. + * + * @param productOfCharacteristics Some value. + * @return 1. + */ + static constexpr element_type get_partial_multiplicative_identity( + [[maybe_unused]] characteristic_type productOfCharacteristics) { + return 1; + } + + // static constexpr bool handles_only_z2() { return false; } + + /** + * @brief Assign operator. + */ + Zp_field_operators& operator=(Zp_field_operators other) { + std::swap(characteristic_, other.characteristic_); + inverse_.swap(other.inverse_); + return *this; + } + /** + * @brief Swap operator. + */ + friend void swap(Zp_field_operators& f1, Zp_field_operators& f2) { + std::swap(f1.characteristic_, f2.characteristic_); + f1.inverse_.swap(f2.inverse_); + } + + private: + characteristic_type characteristic_; /**< Current characteristic of the field. */ + std::vector inverse_; /**< All inverse elements. */ + + static element_type _add(element_type e1, element_type e2, characteristic_type characteristic) { + if (UINT_MAX - e1 < e2) { + // automatic unsigned integer overflow behaviour will make it work + e1 += e2; + e1 -= characteristic; + return e1; + } + + e1 += e2; + if (e1 >= characteristic) e1 -= characteristic; + + return e1; + } + static element_type _substract(element_type e1, element_type e2, characteristic_type characteristic) { + if (e1 < e2) { + e1 += characteristic; + } + e1 -= e2; + + return e1; + } + static element_type _multiply(element_type e1, element_type e2, characteristic_type characteristic) { + unsigned int a = e1; + e1 = 0; + unsigned int temp_b; + + while (a != 0) { + if (a & 1) { + if (e2 >= characteristic - e1) e1 -= characteristic; + e1 += e2; + } + a >>= 1; + + temp_b = e2; + if (e2 >= characteristic - e2) temp_b -= characteristic; + e2 += temp_b; + } + + return e1; + } }; -} //namespace persistence_fields -} //namespace Gudhi +} // namespace persistence_fields +} // namespace Gudhi #endif // MATRIX_FIELD_ZP_OPERATOR_H_ diff --git a/src/Persistence_matrix/include/gudhi/Fields/Zp_field_shared.h b/src/Persistence_matrix/include/gudhi/Fields/Zp_field_shared.h index 63d8ad9d9e..d553559405 100644 --- a/src/Persistence_matrix/include/gudhi/Fields/Zp_field_shared.h +++ b/src/Persistence_matrix/include/gudhi/Fields/Zp_field_shared.h @@ -8,6 +8,12 @@ * - YYYY/MM Author: Description of the modification */ +/** + * @file Zp_field_shared.h + * @author Hannah Schreiber + * @brief Contains the @ref Shared_Zp_field_element class. + */ + #ifndef MATRIX_FIELD_ZP_VAR_H_ #define MATRIX_FIELD_ZP_VAR_H_ @@ -20,322 +26,393 @@ namespace Gudhi { namespace persistence_fields { /** + * @class Shared_Zp_field_element Zp_field_shared.h gudhi/Fields/Zp_field_shared.h * @ingroup persistence_fields * - * @brief TODO: - * + * @brief Class representing an element of the \f$ \mathbb{F}_p \f$ field for any prime number \f$ p \f$. + * If each instanciation of the class can represent another element, they all share the same characteritics. + * That is if the characteristics are set for one, they will be set for all the others. The characteristics can + * be set before instianciating the elements with the static @ref Shared_Zp_field_element::initialize method. + * + * @tparam Unsigned_integer_type A native unsigned integer type: unsigned int, long unsigned int, etc. + * Will be used as the field element type. */ +template > > class Shared_Zp_field_element { -public: - using element_type = unsigned int; - template - using isInteger = std::enable_if_t >; - - Shared_Zp_field_element(); - Shared_Zp_field_element(unsigned int element); - Shared_Zp_field_element(int element); //only works if characteristic can be contained in an int - Shared_Zp_field_element(const Shared_Zp_field_element& toCopy); - Shared_Zp_field_element(Shared_Zp_field_element&& toMove) noexcept; - - static void initialize(unsigned int characteristic); - - friend void operator+=(Shared_Zp_field_element& f1, Shared_Zp_field_element const &f2){ - f1.element_ = Shared_Zp_field_element::_add(f1.element_, f2.element_); - } - friend Shared_Zp_field_element operator+(Shared_Zp_field_element f1, Shared_Zp_field_element const& f2){ - f1 += f2; - return f1; - } - friend void operator+=(Shared_Zp_field_element& f, unsigned int const v){ - f.element_ = Shared_Zp_field_element::_add(f.element_, v < characteristic_ ? v : v % characteristic_); - } - //v is assumed to be positive and will be casted into an unsigned int - template > - friend Shared_Zp_field_element operator+(Shared_Zp_field_element f, const Integer_type v){ - f += v; - return f; - } - //v is assumed to be positive - template > - friend Integer_type operator+(Integer_type v, Shared_Zp_field_element const& f){ - v += f.element_; - if (v >= characteristic_) v %= characteristic_; - return v; - } - - friend void operator-=(Shared_Zp_field_element& f1, Shared_Zp_field_element const &f2){ - f1.element_ = Shared_Zp_field_element::_substract(f1.element_, f2.element_); - } - friend Shared_Zp_field_element operator-(Shared_Zp_field_element f1, Shared_Zp_field_element const& f2){ - f1 -= f2; - return f1; - } - friend void operator-=(Shared_Zp_field_element& f, unsigned int const v){ - f.element_ = Shared_Zp_field_element::_substract(f.element_, v < characteristic_ ? v : v % characteristic_); - } - //v is assumed to be positive and will be casted into an unsigned int - template > - friend Shared_Zp_field_element operator-(Shared_Zp_field_element f, const Integer_type v){ - f -= v; - return f; - } - //v is assumed to be positive - template > - friend Integer_type operator-(Integer_type v, Shared_Zp_field_element const& f){ - if (v >= characteristic_) v %= characteristic_; - if (f.element_ > v) v += characteristic_; - v -= f.element_; - return v; - } - - friend void operator*=(Shared_Zp_field_element& f1, Shared_Zp_field_element const &f2){ - f1.element_ = Shared_Zp_field_element::_multiply(f1.element_, f2.element_); - } - friend Shared_Zp_field_element operator*(Shared_Zp_field_element f1, Shared_Zp_field_element const& f2){ - f1 *= f2; - return f1; - } - friend void operator*=(Shared_Zp_field_element& f, unsigned int const v){ - f.element_ = Shared_Zp_field_element::_multiply(f.element_, v < characteristic_ ? v : v % characteristic_); - } - //v is assumed to be positive and will be casted into an unsigned int - template > - friend Shared_Zp_field_element operator*(Shared_Zp_field_element f, const Integer_type v){ - f *= v; - return f; - } - //uses bitwise operations on v, so be carefull with signed integers - template > - friend Integer_type operator*(Integer_type v, Shared_Zp_field_element const& f){ - unsigned int b = f.element_; - unsigned int res = 0; - unsigned int temp_b; - - while (v != 0) { - if (v & 1) { - if (b >= characteristic_ - res) - res -= characteristic_; - res += b; - } - v >>= 1; - - temp_b = b; - if (b >= characteristic_ - b) - temp_b -= characteristic_; - b += temp_b; - } - - return res; - } - - friend bool operator==(const Shared_Zp_field_element& f1, const Shared_Zp_field_element& f2){ - return f1.element_ == f2.element_; - } - //v is assumed to be positive - template > - friend bool operator==(const Integer_type v, const Shared_Zp_field_element& f){ - if (v < characteristic_) return v == f.element_; - return (v % characteristic_) == f.element_; - } - //v is assumed to be positive - template > - friend bool operator==(const Shared_Zp_field_element& f, const Integer_type v){ - if (v < characteristic_) return v == f.element_; - return (v % characteristic_) == f.element_; - } - friend bool operator!=(const Shared_Zp_field_element& f1, const Shared_Zp_field_element& f2){ - return !(f1 == f2); - } - //v is assumed to be positive - template > - friend bool operator!=(const Integer_type v, const Shared_Zp_field_element& f){ - return !(v == f); - } - //v is assumed to be positive - template > - friend bool operator!=(const Shared_Zp_field_element& f, const Integer_type v){ - return !(v == f); - } - - Shared_Zp_field_element& operator=(Shared_Zp_field_element other); - Shared_Zp_field_element& operator=(const unsigned int value); - operator unsigned int() const; - friend void swap(Shared_Zp_field_element& f1, Shared_Zp_field_element& f2){ - std::swap(f1.element_, f2.element_); - } - - Shared_Zp_field_element get_inverse() const; - std::pair get_partial_inverse(unsigned int productOfCharacteristics) const; - - static Shared_Zp_field_element get_additive_identity(); - static Shared_Zp_field_element get_multiplicative_identity(); - static Shared_Zp_field_element get_partial_multiplicative_identity([[maybe_unused]] unsigned int productOfCharacteristics); - static unsigned int get_characteristic(); - - unsigned int get_value() const; - - static constexpr bool handles_only_z2(){ - return false; - } - -private: - unsigned int element_; - static inline unsigned int characteristic_; - static inline std::vector inverse_; - - static unsigned int _add(unsigned int element, unsigned int v); - static unsigned int _substract(unsigned int element, unsigned int v); - static unsigned int _multiply(unsigned int element, unsigned int v); + public: + using element_type = Unsigned_integer_type; /**< Type for the elements in the field. */ + using characteristic_type = element_type; /**< Type for the field characteristic. */ + template + using isInteger = std::enable_if_t >; + + /** + * @brief Default constructor. Sets the element to 0. + */ + Shared_Zp_field_element() : element_(0) {} + /** + * @brief Constructor setting the element to the given value. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + * @param element Value of the element. + */ + template > + Shared_Zp_field_element(Integer_type element) : element_(_get_value(element)) {} + /** + * @brief Copy constructor. + * + * @param toCopy Element to copy. + */ + Shared_Zp_field_element(const Shared_Zp_field_element& toCopy) : element_(toCopy.element_) {} + /** + * @brief Move constructor. + * + * @param toMove Element to move. + */ + Shared_Zp_field_element(Shared_Zp_field_element&& toMove) noexcept : element_(std::exchange(toMove.element_, 0)) {} + + /** + * @brief Initialize the field to the given characteristic. + * Should be called first before constructing the field elements. + * + * @param characteristic Characteristic of the field. A positive prime number. + */ + static void initialize(characteristic_type characteristic) { + if (characteristic <= 1) + throw std::invalid_argument("Characteristic must be strictly positive and a prime number."); + + inverse_.resize(characteristic); + inverse_[0] = 0; + for (element_type i = 1; i < characteristic; ++i) { + element_type inv = 1; + element_type mult = inv * i; + while ((mult % characteristic) != 1) { + ++inv; + if (mult == characteristic) throw std::invalid_argument("Characteristic must be a prime number."); + mult = inv * i; + } + inverse_[i] = inv; + } + + characteristic_ = characteristic; + } + + /** + * @brief Returns the value of the element. + * + * @return Value of the element. + */ + element_type get_value() const { return element_; } + + /** + * @brief operator+= + */ + friend void operator+=(Shared_Zp_field_element& f1, const Shared_Zp_field_element& f2) { + f1.element_ = Shared_Zp_field_element::_add(f1.element_, f2.element_); + } + /** + * @brief operator+ + */ + friend Shared_Zp_field_element operator+(Shared_Zp_field_element f1, const Shared_Zp_field_element& f2) { + f1 += f2; + return f1; + } + /** + * @brief operator+= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend void operator+=(Shared_Zp_field_element& f, const Integer_type& v) { + f.element_ = Shared_Zp_field_element::_add(f.element_, _get_value(v)); + } + /** + * @brief operator+ + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Shared_Zp_field_element operator+(Shared_Zp_field_element f, const Integer_type& v) { + f += v; + return f; + } + /** + * @brief operator+ + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Integer_type operator+(const Integer_type& v, Shared_Zp_field_element f) { + f += v; + return f.element_; + } + + /** + * @brief operator-= + */ + friend void operator-=(Shared_Zp_field_element& f1, const Shared_Zp_field_element& f2) { + f1.element_ = Shared_Zp_field_element::_substract(f1.element_, f2.element_); + } + /** + * @brief operator- + */ + friend Shared_Zp_field_element operator-(Shared_Zp_field_element f1, const Shared_Zp_field_element& f2) { + f1 -= f2; + return f1; + } + /** + * @brief operator-= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend void operator-=(Shared_Zp_field_element& f, const Integer_type& v) { + f.element_ = Shared_Zp_field_element::_substract(f.element_, _get_value(v)); + } + /** + * @brief operator- + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Shared_Zp_field_element operator-(Shared_Zp_field_element f, const Integer_type& v) { + f -= v; + return f; + } + /** + * @brief operator- + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Integer_type operator-(const Integer_type& v, const Shared_Zp_field_element& f) { + return Shared_Zp_field_element::_substract(_get_value(v), f.element_); + } + + /** + * @brief operator*= + */ + friend void operator*=(Shared_Zp_field_element& f1, const Shared_Zp_field_element& f2) { + f1.element_ = Shared_Zp_field_element::_multiply(f1.element_, f2.element_); + } + /** + * @brief operator* + */ + friend Shared_Zp_field_element operator*(Shared_Zp_field_element f1, const Shared_Zp_field_element& f2) { + f1 *= f2; + return f1; + } + /** + * @brief operator*= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend void operator*=(Shared_Zp_field_element& f, const Integer_type& v) { + f.element_ = Shared_Zp_field_element::_multiply(f.element_, _get_value(v)); + } + /** + * @brief operator* + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend Shared_Zp_field_element operator*(Shared_Zp_field_element f, const Integer_type& v) { + f *= v; + return f; + } + /** + * @brief operator* + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic. if signed. + */ + template > + friend Integer_type operator*(const Integer_type& v, Shared_Zp_field_element f) { + f *= v; + return f.element_; + } + + /** + * @brief operator== + */ + friend bool operator==(const Shared_Zp_field_element& f1, const Shared_Zp_field_element& f2) { + return f1.element_ == f2.element_; + } + /** + * @brief operator== + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator==(const Integer_type& v, const Shared_Zp_field_element& f) { + return Shared_Zp_field_element::_get_value(v) == f.element_; + } + /** + * @brief operator== + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator==(const Shared_Zp_field_element& f, const Integer_type& v) { + return Shared_Zp_field_element::_get_value(v) == f.element_; + } + /** + * @brief operator!= + */ + friend bool operator!=(const Shared_Zp_field_element& f1, const Shared_Zp_field_element& f2) { return !(f1 == f2); } + /** + * @brief operator!= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator!=(const Integer_type& v, const Shared_Zp_field_element& f) { + return !(v == f); + } + /** + * @brief operator!= + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + friend bool operator!=(const Shared_Zp_field_element& f, const Integer_type& v) { + return !(v == f); + } + + /** + * @brief Assign operator. + */ + Shared_Zp_field_element& operator=(Shared_Zp_field_element other) { + std::swap(element_, other.element_); + return *this; + } + /** + * @brief Assign operator. + * + * @tparam Integer_type A native integer type. Should be able to contain the characteristic if signed. + */ + template > + Shared_Zp_field_element& operator=(const Integer_type& value) { + element_ = Shared_Zp_field_element::_get_value(value); + return *this; + } + /** + * @brief Swap operator. + */ + friend void swap(Shared_Zp_field_element& f1, Shared_Zp_field_element& f2) { std::swap(f1.element_, f2.element_); } + + /** + * @brief Casts the element into an unsigned int. + */ + operator unsigned int() const { return element_; } + + /** + * @brief Returns the inverse of the element in the field. + * + * @return The inverse. + */ + Shared_Zp_field_element get_inverse() const { return Shared_Zp_field_element(inverse_[element_]); } + /** + * @brief For interface purposes with multi-fields. Returns the inverse together with the argument. + * + * @param productOfCharacteristics Some value. + * @return Pair whose first element is the inverse and the second element is @p productOfCharacteristics. + */ + std::pair get_partial_inverse( + characteristic_type productOfCharacteristics) const { + return {get_inverse(), productOfCharacteristics}; + } + + /** + * @brief Returns the additive identity of the field. + * + * @return 0. + */ + static Shared_Zp_field_element get_additive_identity() { return Shared_Zp_field_element(); } + /** + * @brief Returns the multiplicative identity of the field. + * + * @return 1. + */ + static Shared_Zp_field_element get_multiplicative_identity() { return Shared_Zp_field_element(1); } + /** + * @brief For interface purposes with multi-fields. Returns the multiplicative identity of the field. + * + * @param productOfCharacteristics Some value. + * @return 1. + */ + static Shared_Zp_field_element get_partial_multiplicative_identity( + [[maybe_unused]] characteristic_type productOfCharacteristics) { + return Shared_Zp_field_element(1); + } + /** + * @brief Returns the current characteristic. + * + * @return The value of the current characteristic. + */ + static characteristic_type get_characteristic() { return characteristic_; } + + // static constexpr bool handles_only_z2() { return false; } + + private: + element_type element_; /**< Field element. */ + static inline characteristic_type characteristic_; /**< Current characteristic of the field. */ + static inline std::vector inverse_; /**< All inverse elements. */ + + static element_type _add(element_type element, element_type v) { + if (UINT_MAX - element < v) { + // automatic unsigned integer overflow behaviour will make it work + element += v; + element -= characteristic_; + return element; + } + + element += v; + if (element >= characteristic_) element -= characteristic_; + + return element; + } + static element_type _substract(element_type element, element_type v) { + if (element < v) { + element += characteristic_; + } + element -= v; + + return element; + } + static element_type _multiply(element_type element, element_type v) { + element_type a = element; + element = 0; + element_type temp_b; + + while (a != 0) { + if (a & 1) { + if (v >= characteristic_ - element) element -= characteristic_; + element += v; + } + a >>= 1; + + temp_b = v; + if (v >= characteristic_ - v) temp_b -= characteristic_; + v += temp_b; + } + + return element; + } + + template > + static constexpr element_type _get_value(Integer_type e) { + if constexpr (std::is_signed_v){ + if (e < -static_cast(characteristic_)) e = e % characteristic_; + if (e < 0) return e += characteristic_; + return e < static_cast(characteristic_) ? e : e % characteristic_; + } else { + return e < characteristic_ ? e : e % characteristic_; + } + } }; -inline Shared_Zp_field_element::Shared_Zp_field_element() - : element_(0) -{} - -inline Shared_Zp_field_element::Shared_Zp_field_element(unsigned int element) - : element_(element < characteristic_ ? element : element % characteristic_) -{} - -inline Shared_Zp_field_element::Shared_Zp_field_element(int element) - : element_() -{ - int res = element < static_cast(characteristic_) ? element : element % static_cast(characteristic_); - if (res < 0) res += characteristic_; - element_ = res; -} - -inline Shared_Zp_field_element::Shared_Zp_field_element(const Shared_Zp_field_element &toCopy) - : element_(toCopy.element_) -{} - -inline Shared_Zp_field_element::Shared_Zp_field_element(Shared_Zp_field_element &&toMove) noexcept - : element_(std::exchange(toMove.element_, 0)) -{} - -inline void Shared_Zp_field_element::initialize(unsigned int characteristic) -{ - if (characteristic <= 1) - throw std::invalid_argument("Characteristic must be strictly positive and a prime number."); - - inverse_.resize(characteristic); - inverse_[0] = 0; - for (unsigned int i = 1; i < characteristic; ++i) { - unsigned int inv = 1; - unsigned int mult = inv * i; - while ((mult % characteristic) != 1) { - ++inv; - if (mult == characteristic) - throw std::invalid_argument("Characteristic must be a prime number."); - mult = inv * i; - } - inverse_[i] = inv; - } - - characteristic_ = characteristic; -} - -inline Shared_Zp_field_element &Shared_Zp_field_element::operator=(Shared_Zp_field_element other) -{ - std::swap(element_, other.element_); - return *this; -} - -inline Shared_Zp_field_element &Shared_Zp_field_element::operator=(unsigned int const value) -{ - element_ = value < characteristic_ ? value : value % characteristic_; - return *this; -} - -inline Shared_Zp_field_element::operator unsigned int() const -{ - return element_; -} - -inline Shared_Zp_field_element Shared_Zp_field_element::get_inverse() const -{ - return Shared_Zp_field_element(inverse_[element_]); -} - -inline std::pair -Shared_Zp_field_element::get_partial_inverse(unsigned int productOfCharacteristics) const -{ - return {get_inverse(), productOfCharacteristics}; -} - -inline Shared_Zp_field_element Shared_Zp_field_element::get_additive_identity() -{ - return Shared_Zp_field_element(); -} - -inline Shared_Zp_field_element Shared_Zp_field_element::get_multiplicative_identity() -{ - return Shared_Zp_field_element(1); -} - -inline Shared_Zp_field_element Shared_Zp_field_element::get_partial_multiplicative_identity([[maybe_unused]] unsigned int productOfCharacteristics) -{ - return Shared_Zp_field_element(1); -} - -inline unsigned int Shared_Zp_field_element::get_characteristic() -{ - return characteristic_; -} - -inline unsigned int Shared_Zp_field_element::get_value() const -{ - return element_; -} - -inline unsigned int Shared_Zp_field_element::_add(unsigned int element, unsigned int v) -{ - if (UINT_MAX - element < v) { - //automatic unsigned integer overflow behaviour will make it work - element += v; - element -= characteristic_; - return element; - } - - element += v; - if (element >= characteristic_) element -= characteristic_; - - return element; -} - -inline unsigned int Shared_Zp_field_element::_substract(unsigned int element, unsigned int v) -{ - if (element < v){ - element += characteristic_; - } - element -= v; - - return element; -} - -inline unsigned int Shared_Zp_field_element::_multiply(unsigned int element, unsigned int v) -{ - unsigned int a = element; - element = 0; - unsigned int temp_b; - - while (a != 0) { - if (a & 1) { - if (v >= characteristic_ - element) - element -= characteristic_; - element += v; - } - a >>= 1; - - temp_b = v; - if (v >= characteristic_ - v) - temp_b -= characteristic_; - v += temp_b; - } - - return element; -} - -} //namespace persistence_fields -} //namespace Gudhi +} // namespace persistence_fields +} // namespace Gudhi #endif // MATRIX_FIELD_ZP_VAR_H_ diff --git a/src/Persistence_matrix/test/Persistence_matrix_field_tests.cpp b/src/Persistence_matrix/test/Persistence_matrix_field_tests.cpp index 1c7131c7f2..3874ed70d7 100644 --- a/src/Persistence_matrix/test/Persistence_matrix_field_tests.cpp +++ b/src/Persistence_matrix/test/Persistence_matrix_field_tests.cpp @@ -326,35 +326,35 @@ BOOST_AUTO_TEST_CASE(Field_properties) BOOST_AUTO_TEST_CASE(Shared_Field_constructors) { - Shared_Zp_field_element::initialize(2); - test_z2_standart_field_constructors(); + Shared_Zp_field_element<>::initialize(2); + test_z2_standart_field_constructors >(); - Shared_Zp_field_element::initialize(5); - test_z5_standart_field_constructors(); + Shared_Zp_field_element<>::initialize(5); + test_z5_standart_field_constructors >(); - Shared_Zp_field_element::initialize(13); - test_z13_standart_field_constructors(); + Shared_Zp_field_element<>::initialize(13); + test_z13_standart_field_constructors >(); } BOOST_AUTO_TEST_CASE(Shared_Field_operators) { - Shared_Zp_field_element::initialize(2); - test_z2_standart_field_operators(); + Shared_Zp_field_element<>::initialize(2); + test_z2_standart_field_operators >(); - Shared_Zp_field_element::initialize(5); - test_z5_standart_field_operators(); + Shared_Zp_field_element<>::initialize(5); + test_z5_standart_field_operators >(); } BOOST_AUTO_TEST_CASE(Shared_Field_properties) { - Shared_Zp_field_element::initialize(2); - test_z2_standart_field_properties(); + Shared_Zp_field_element<>::initialize(2); + test_z2_standart_field_properties >(); - Shared_Zp_field_element::initialize(5); - test_z5_standart_field_properties(); + Shared_Zp_field_element<>::initialize(5); + test_z5_standart_field_properties >(); - Shared_Zp_field_element::initialize(7); - test_z7_standart_field_properties(); + Shared_Zp_field_element<>::initialize(7); + test_z7_standart_field_properties >(); } template @@ -362,29 +362,29 @@ void test_multi_field_constructors(){ using T = typename MF::element_type; //default constructor - Multi_field_element<5,13> m_d; + MF m_d; BOOST_CHECK_EQUAL(m_d, T(0)); //value constructor - Multi_field_element<5,13> m_v(5006); + MF m_v(5006); BOOST_CHECK_EQUAL(m_v, T(1)); //copy constructor - Multi_field_element<5,13> m_c1(5006); - Multi_field_element<5,13> m_c2 = m_c1; + MF m_c1(5006); + MF m_c2 = m_c1; BOOST_CHECK_EQUAL(m_c2, T(1)); - Multi_field_element<5,13> m_c3(m_c2); + MF m_c3(m_c2); BOOST_CHECK_EQUAL(m_c3, T(1)); //move constructor - Multi_field_element<5,13> m_m1(5006); - Multi_field_element<5,13> m_m2(std::move(m_m1)); + MF m_m1(5006); + MF m_m2(std::move(m_m1)); BOOST_CHECK_EQUAL(m_m2, T(1)); BOOST_CHECK_EQUAL(m_m1, T(0)); //swap - Multi_field_element<5,13> m_s1(5006); - Multi_field_element<5,13> m_s2(5005); + MF m_s1(5006); + MF m_s2(5005); swap(m_s1, m_s2); BOOST_CHECK_EQUAL(m_s2, T(1)); BOOST_CHECK_EQUAL(m_s1, T(0)); @@ -499,8 +499,8 @@ BOOST_AUTO_TEST_CASE(Shared_Multi_Field_constructors) Shared_multi_field_element::initialize(5, 13); test_multi_field_constructors(); - Shared_multi_field_element_with_small_characteristics::initialize(5, 13); - test_multi_field_constructors(); + Shared_multi_field_element_with_small_characteristics<>::initialize(5, 13); + test_multi_field_constructors >(); } BOOST_AUTO_TEST_CASE(Shared_Multi_Field_operators) @@ -508,8 +508,8 @@ BOOST_AUTO_TEST_CASE(Shared_Multi_Field_operators) Shared_multi_field_element::initialize(5, 13); test_multi_field_operators(); - Shared_multi_field_element_with_small_characteristics::initialize(5, 13); - test_multi_field_operators(); + Shared_multi_field_element_with_small_characteristics<>::initialize(5, 13); + test_multi_field_operators >(); } BOOST_AUTO_TEST_CASE(Shared_Multi_Field_properties) @@ -517,13 +517,13 @@ BOOST_AUTO_TEST_CASE(Shared_Multi_Field_properties) Shared_multi_field_element::initialize(5, 13); test_multi_field_properties(); - Shared_multi_field_element_with_small_characteristics::initialize(5, 13); - test_multi_field_properties(); + Shared_multi_field_element_with_small_characteristics<>::initialize(5, 13); + test_multi_field_properties >(); Shared_multi_field_element::initialize(3, 30); - Shared_multi_field_element_with_small_characteristics::initialize(3, 30); + Shared_multi_field_element_with_small_characteristics<>::initialize(3, 30); Shared_multi_field_element mb1(2); - Shared_multi_field_element_with_small_characteristics mb2(2); + Shared_multi_field_element_with_small_characteristics<> mb2(2); BOOST_CHECK_EQUAL(mb1.get_characteristic(), mb2.get_characteristic()); // == 3234846615 BOOST_CHECK_EQUAL(mb1.get_partial_inverse(35).first.get_value(), mb2.get_partial_inverse(35).first.get_value()); // == 2033332158 diff --git a/src/common/doc/main_page.md b/src/common/doc/main_page.md index d32990d8fe..943edbbf44 100644 --- a/src/common/doc/main_page.md +++ b/src/common/doc/main_page.md @@ -6,7 +6,7 @@ ## Data structures for cell complexes {#Complexes} ### Cubical complexes - + - - + + - +
@@ -21,12 +21,12 @@ Introduced in: GUDHI 1.3.0
Copyright: MIT
- User manual: \ref cubical_complex + User manual: \ref cubical_complex
### Simplicial complexes @@ -39,20 +39,20 @@ The simplex tree is an efficient and flexible - data structure for representing general (filtered) simplicial complexes. The data structure - is described in \cite boissonnatmariasimplextreealgorithmica . + data structure for representing general (filtered) simplicial complexes. The data structure + is described in \cite boissonnatmariasimplextreealgorithmica . Author: Clément Maria
Introduced in: GUDHI 1.0.0
Copyright: MIT
- - + + - User manual: \ref simplex_tree + User manual: \ref simplex_tree - + #### Toplex Map @@ -72,12 +72,12 @@ Introduced in: GUDHI 2.1.0
Copyright: MIT
- - + + - User manual: \ref toplex_map + User manual: \ref toplex_map - + #### Skeleton blocker @@ -100,12 +100,12 @@ Introduced in: GUDHI 1.1.0
Copyright: MIT
- - + + - User manual: \ref skbl + User manual: \ref skbl - + #### Basic operation: contraction @@ -127,12 +127,12 @@ Copyright: MIT [(LGPL v3)](../../licensing/)
Requires: \ref cgal ≥ 4.11.0 - - + + - User manual: \ref contr + User manual: \ref contr - + ## Filtrations @@ -158,18 +158,18 @@ Copyright: MIT [(GPL v3)](../../licensing/)
Requires: \ref eigen ≥ 3.1.0 and \ref cgal ≥ 4.11.0 - - + + - User manual: \ref alpha_complex + User manual: \ref alpha_complex - + ### Čech complex - + @@ -183,12 +183,12 @@ Copyright: MIT [(LGPL v3)](../../licensing/)
Requires: \ref cgal - - + + - +
\image html "cech_complex_representation.png"
- User manual: \ref cech_complex + User manual: \ref cech_complex
### Rips complex @@ -209,12 +209,12 @@ Introduced in: GUDHI 2.0.0
Copyright: MIT
- - + + - User manual: \ref rips_complex + User manual: \ref rips_complex - + ### Edge collapse @@ -238,12 +238,12 @@ Introduced in: GUDHI 3.3.0
Copyright: MIT - - + + - User manual: \ref edge_collapse + User manual: \ref edge_collapse - + ### Witness complex @@ -263,12 +263,12 @@ Copyright: MIT ([GPL v3](../../licensing/) for Euclidean version)
Euclidean version requires: \ref eigen ≥ 3.1.0 and \ref cgal ≥ 4.11.0 - - + + - User manual: \ref witness_complex + User manual: \ref witness_complex - + ### Cover Complexes @@ -289,12 +289,12 @@ Copyright: MIT [(GPL v3)](../../licensing/)
Requires: \ref cgal ≥ 4.11.0 - - + + - User manual: \ref cover_complex + User manual: \ref cover_complex - + ## Manifold reconstructions @@ -315,12 +315,12 @@ Copyright: MIT [(LGPL v3)](../../licensing/)
Requires: \ref eigen ≥ 3.1.0 - - + + User manual: \ref coxeter_triangulation - + ### Tangential complex @@ -343,12 +343,12 @@ Copyright: MIT [(GPL v3)](../../licensing/)
Requires: \ref eigen ≥ 3.1.0 and \ref cgal ≥ 4.11.0 - - + + - User manual: \ref tangential_complex + User manual: \ref tangential_complex - + ## Topological descriptors computation {#TopologicalDescriptorsComputation} @@ -375,16 +375,17 @@ Introduced in: GUDHI 1.0.0
Copyright: MIT
- - + + - User manual: \ref persistent_cohomology + User manual: \ref persistent_cohomology - + ### Persistence Matrix +[//]: # (TODO: picture) - - + + - +
@@ -393,19 +394,18 @@ Matrix structure for filtered complexes with multiple functionnalities related to persistence homology, such as representative cycles computation or vineyards. - TODO: picture Author: Hannah Schreiber
Introduced in: GUDHI 1.0.0
Copyright: MIT
- User manual: \ref persistence_matrix and \ref persistence_fields + User manual: \ref persistence_matrix and \ref persistence_fields
## Topological descriptors tools {#TopologicalDescriptorsTools} @@ -431,12 +431,12 @@ Copyright: MIT [(GPL v3)](../../licensing/)
Requires: \ref cgal ≥ 4.11.0 - - + + - User manual: \ref bottleneck_distance + User manual: \ref bottleneck_distance - + ### Persistence representations @@ -456,12 +456,12 @@ Introduced in: GUDHI 2.1.0
Copyright: MIT
- - + + - User manual: \ref Persistence_representations + User manual: \ref Persistence_representations - + ## Point cloud utilities {#PointCloudUtils} @@ -479,10 +479,10 @@ Introduced in: GUDHI 1.3.0
Copyright: MIT [(GPL v3)](../../licensing/)
- - + + - Manuals: \ref spatial_searching, \ref subsampling + Manuals: \ref spatial_searching, \ref subsampling - +