From: Richard Kreckel Date: Wed, 15 Feb 2017 11:58:48 +0000 (+0100) Subject: Clean up combinatorial helpers. X-Git-Tag: release_1-7-3~22 X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=commitdiff_plain;h=a528cac19001a5beba9d54832b341ed84dc5bd53;ds=sidebyside Clean up combinatorial helpers. * Convert all signatures from int to unsiged. * Renamed partition_generator to partition_with_zero_parts_generator... * ...and add class partition_generator which does not include zero parts. * Rename ::current() to ::get() and make them compute vectors on first use. --- diff --git a/ginac/power.cpp b/ginac/power.cpp index 3c718027..26194928 100644 --- a/ginac/power.cpp +++ b/ginac/power.cpp @@ -979,9 +979,9 @@ ex power::expand_add(const add & a, long n, unsigned options) // Multinomial expansion of power(+(x,...,z;0),k)*c^(n-k): // Iterate over all partitions of k with exactly as many parts as // there are symbolic terms in the basis (including zero parts). - partition_generator partitions(k, a.seq.size()); + partition_with_zero_parts_generator partitions(k, a.seq.size()); do { - const std::vector& partition = partitions.current(); + const std::vector& partition = partitions.get(); // All monomials of this partition have the same number of terms and the same coefficient. const unsigned msize = std::count_if(partition.begin(), partition.end(), [](int i) { return i > 0; }); const numeric coeff = multinomial_coefficient(partition) * binomial_coefficient; @@ -989,7 +989,7 @@ ex power::expand_add(const add & a, long n, unsigned options) // Iterate over all compositions of the current partition. composition_generator compositions(partition); do { - const std::vector& exponent = compositions.current(); + const std::vector& exponent = compositions.get(); epvector monomial; monomial.reserve(msize); numeric factor = coeff; diff --git a/ginac/utils.cpp b/ginac/utils.cpp index 89217a93..8085ba3b 100644 --- a/ginac/utils.cpp +++ b/ginac/utils.cpp @@ -23,7 +23,6 @@ #include "ex.h" #include "numeric.h" -#include "operators.h" #include "utils.h" #include "version.h" @@ -58,14 +57,14 @@ unsigned log2(unsigned n) * n = p1+p2+...+pk, i.e. p is a partition of n. */ const numeric -multinomial_coefficient(const std::vector & p) +multinomial_coefficient(const std::vector & p) { numeric n = 0, d = 1; for (auto & it : p) { - n += numeric(it); - d *= factorial(numeric(it)); + n = n.add(numeric(it)); + d = d.mul(factorial(numeric(it))); } - return factorial(n) / d; + return factorial(n).div(d); } ////////// diff --git a/ginac/utils.h b/ginac/utils.h index 816679d7..8d3f2633 100644 --- a/ginac/utils.h +++ b/ginac/utils.h @@ -272,32 +272,32 @@ again: } } -/** Generate all bounded combinatorial partitions of an integer n with exactly - * m parts (including zero parts) in non-decreasing order. +/** Base class for generating all bounded combinatorial partitions of an integer + * n with exactly m parts in non-decreasing order. */ -class partition_generator { -private: +class basic_partition_generator { +protected: // Partitions n into m parts, not including zero parts. // (Cf. OEIS sequence A008284; implementation adapted from Jörg Arndt's // FXT library) struct mpartition2 { // partition: x[1] + x[2] + ... + x[m] = n and sentinel x[0] == 0 - std::vector x; - int n; // n>0 - int m; // 0 x; + unsigned n; // n>0 + unsigned m; // 0 partition; // current partition + }; + mpartition2 mpgen; + basic_partition_generator(unsigned n_, unsigned m_) + : mpgen(n_, m_) + { } +}; + +/** Generate all bounded combinatorial partitions of an integer n with exactly + * m parts (including zero parts) in non-decreasing order. + */ +class partition_with_zero_parts_generator : public basic_partition_generator { +private: + unsigned m; // number of parts 0 partition; // current partition + mutable bool current_updated; // whether partition vector has been updated public: - partition_generator(unsigned n_, unsigned m_) - : mpgen(n_, 1), m(m_), partition(m_) + partition_with_zero_parts_generator(unsigned n_, unsigned m_) + : basic_partition_generator(n_, 1), m(m_), partition(m_), current_updated(false) { } // returns current partition in non-decreasing order, padded with zeros - const std::vector& current() const + const std::vector& get() const { - for (int i = 0; i < m - mpgen.m; ++i) - partition[i] = 0; // pad with zeros - - for (int i = m - mpgen.m; i < m; ++i) - partition[i] = mpgen.x[i - m + mpgen.m + 1]; + if (!current_updated) { + for (unsigned i = 0; i < m - mpgen.m; ++i) + partition[i] = 0; // pad with zeros + for (unsigned i = m - mpgen.m; i < m; ++i) + partition[i] = mpgen.x[i - m + mpgen.m + 1]; + } return partition; } bool next() { + current_updated = false; if (!mpgen.next_partition()) { if (mpgen.m == m || mpgen.m == mpgen.n) return false; // current is last @@ -344,6 +358,33 @@ public: } }; +/** Generate all bounded combinatorial partitions of an integer n with exactly + * m parts (not including zero parts) in non-decreasing order. + */ +class partition_generator : public basic_partition_generator { +private: + mutable std::vector partition; // current partition + mutable bool current_updated; // whether partition vector has been updated +public: + partition_generator(unsigned n_, unsigned m_) + : basic_partition_generator(n_, m_), partition(m_), current_updated(false) + { } + // returns current partition in non-decreasing order, padded with zeros + const std::vector& get() const + { + if (!current_updated) { + for (unsigned i = 0; i < mpgen.m; ++i) + partition[i] = mpgen.x[i + 1]; + } + return partition; + } + bool next() + { + current_updated = false; + return mpgen.next_partition(); + } +}; + /** Generate all compositions of a partition of an integer n, starting with the * compositions which has non-decreasing order. */ @@ -356,9 +397,9 @@ private: struct coolmulti { // element of singly linked list struct element { - int value; + unsigned value; element* next; - element(int val, element* n) + element(unsigned val, element* n) : value(val), next(n) {} ~element() { // recurses down to the end of the singly linked list @@ -367,7 +408,7 @@ private: }; element *head, *i, *after_i; // NB: Partition must be sorted in non-decreasing order. - explicit coolmulti(const std::vector& partition) + explicit coolmulti(const std::vector& partition) : head(nullptr), i(nullptr), after_i(nullptr) { for (unsigned n = 0; n < partition.size(); ++n) { @@ -403,22 +444,26 @@ private: } cmgen; bool atend; // needed for simplifying iteration over permutations bool trivial; // likewise, true if all elements are equal - mutable std::vector composition; // current compositions + mutable std::vector composition; // current compositions + mutable bool current_updated; // whether composition vector has been updated public: - explicit composition_generator(const std::vector& partition) - : cmgen(partition), atend(false), trivial(true), composition(partition.size()) + explicit composition_generator(const std::vector& partition) + : cmgen(partition), atend(false), trivial(true), composition(partition.size()), current_updated(false) { for (unsigned i=1; i& current() const + const std::vector& get() const { - coolmulti::element* it = cmgen.head; - size_t i = 0; - while (it != nullptr) { - composition[i] = it->value; - it = it->next; - ++i; + if (!current_updated) { + coolmulti::element* it = cmgen.head; + size_t i = 0; + while (it != nullptr) { + composition[i] = it->value; + it = it->next; + ++i; + } + current_updated = true; } return composition; } @@ -430,6 +475,7 @@ public: if (trivial || atend) return false; cmgen.next_permutation(); + current_updated = false; atend = cmgen.finished(); return true; } @@ -439,7 +485,7 @@ public: * n = p1+p2+...+pk, i.e. p is a partition of n. */ const numeric -multinomial_coefficient(const std::vector & p); +multinomial_coefficient(const std::vector & p); // Collection of `construct on first use' wrappers for safely avoiding