X-Git-Url: https://www.ginac.de/ginac.git//ginac.git?p=ginac.git;a=blobdiff_plain;f=ginac%2Futils.h;h=7c67085d2d2aab537fdf2d1b1783bf311b207ebf;hp=ae17f53ae34df211c07fd1a0f02c001a219b8e1e;hb=dcb735feae4c3a6962a7824e03b0ce4496e3ffe7;hpb=22abfbe8c78e339188096a5bf749a7c2d4f0a368 diff --git a/ginac/utils.h b/ginac/utils.h index ae17f53a..7c67085d 100644 --- a/ginac/utils.h +++ b/ginac/utils.h @@ -4,7 +4,7 @@ * of any interest to the user of the library. */ /* - * GiNaC Copyright (C) 1999-2005 Johannes Gutenberg University Mainz, Germany + * GiNaC Copyright (C) 1999-2021 Johannes Gutenberg University Mainz, Germany * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by @@ -21,15 +21,14 @@ * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA */ -#ifndef __GINAC_UTILS_H__ -#define __GINAC_UTILS_H__ +#ifndef GINAC_UTILS_H +#define GINAC_UTILS_H -#include "config.h" +#include "assertion.h" -#include #include - -#include "assertion.h" +#include // for uintptr_t +#include namespace GiNaC { @@ -44,6 +43,13 @@ class dunno {}; unsigned log2(unsigned n); +/** Rotate bits of unsigned value by one bit to the left. + * This can be necessary if the user wants to define its own hashes. */ +inline unsigned rotate_left(unsigned n) +{ + return (n & 0x80000000U) ? (n << 1 | 0x00000001U) : (n << 1); +} + /** Compare two pointers (just to establish some sort of canonical order). * @return -1, 0, or 1 */ template @@ -58,35 +64,10 @@ inline int compare_pointers(const T * a, const T * b) return 0; } -/** Rotate bits of unsigned value by one bit to the left. */ -inline unsigned rotate_left(unsigned n) -{ - return (n & 0x80000000U) ? (n << 1 | 0x00000001U) : (n << 1); -} - /** Truncated multiplication with golden ratio, for computing hash values. */ -inline unsigned golden_ratio_hash(unsigned n) +inline unsigned golden_ratio_hash(uintptr_t n) { - // This function works much better when fast arithmetic with at - // least 64 significant bits is available. -#if SIZEOF_LONG >= 8 - // So 'long' has 64 bits. Excellent! We prefer it because it might be - // more efficient than 'long long'. - unsigned long l = n * 0x4f1bbcddUL; - return (unsigned)l; -#elif SIZEOF_LONG_LONG >= 8 - // This requires 'long long' (or an equivalent 64 bit type)---which is, - // unfortunately, not ANSI-C++-compliant. - // (Yet C99 demands it, which is reason for hope.) - unsigned long long l = n * 0x4f1bbcddULL; - return (unsigned)l; -#else - // Without a type with 64 significant bits do the multiplication manually - // by splitting n up into the lower and upper two bytes. - const unsigned n0 = (n & 0x0000ffffU); - const unsigned n1 = (n & 0xffff0000U) >> 16; - return (n0 * 0x0000bcddU) + ((n1 * 0x0000bcddU + n0 * 0x00004f1bU) << 16); -#endif + return n * UINT64_C(0x4f1bbcdd); } /* Compute the sign of a permutation of a container, with and without an @@ -95,6 +76,7 @@ inline unsigned golden_ratio_hash(unsigned n) template int permutation_sign(It first, It last) { + using std::swap; if (first == last) return 0; --last; @@ -109,13 +91,15 @@ int permutation_sign(It first, It last) bool swapped = false; while (i != first) { if (*i < *other) { - std::iter_swap(other, i); + swap(*other, *i); flag = other; swapped = true; sign = -sign; } else if (!(*other < *i)) return 0; - --i; --other; + --i; + if (i != first) + --other; } if (!swapped) return sign; @@ -128,13 +112,15 @@ int permutation_sign(It first, It last) swapped = false; while (i != last) { if (*other < *i) { - std::iter_swap(i, other); + swap(*i, *other); flag = other; swapped = true; sign = -sign; } else if (!(*i < *other)) return 0; - ++i; ++other; + ++i; + if (i != last) + ++other; } if (!swapped) return sign; @@ -168,7 +154,9 @@ int permutation_sign(It first, It last, Cmp comp, Swap swapit) sign = -sign; } else if (!comp(*other, *i)) return 0; - --i; --other; + --i; + if (i != first) + --other; } if (!swapped) return sign; @@ -187,7 +175,9 @@ int permutation_sign(It first, It last, Cmp comp, Swap swapit) sign = -sign; } else if (!comp(*i, *other)) return 0; - ++i; ++other; + ++i; + if (i != last) + ++other; } if (!swapped) return sign; @@ -219,7 +209,9 @@ void shaker_sort(It first, It last, Cmp comp, Swap swapit) flag = other; swapped = true; } - --i; --other; + --i; + if (i != first) + --other; } if (!swapped) return; @@ -236,7 +228,9 @@ void shaker_sort(It first, It last, Cmp comp, Swap swapit) flag = other; swapped = true; } - ++i; ++other; + ++i; + if (i != last) + ++other; } if (!swapped) return; @@ -278,6 +272,225 @@ again: } } +/** Base class for generating all bounded combinatorial partitions of an integer + * n with exactly m parts in non-decreasing order. + */ +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; + unsigned n; // n>0 + unsigned m; // 0 partition; // current partition + mutable bool current_updated; // whether partition vector has been updated +public: + 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& get() const + { + 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]; + + current_updated = true; + } + return partition; + } + bool next() + { + current_updated = false; + if (!mpgen.next_partition()) { + if (mpgen.m == m || mpgen.m == mpgen.n) + return false; // current is last + // increment number of parts + mpgen = mpartition2(mpgen.n, mpgen.m + 1); + } + return true; + } +}; + +/** 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]; + + current_updated = true; + } + 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. + */ +class composition_generator { +private: + // Generates all distinct permutations of a multiset. + // (Based on Aaron Williams' algorithm 1 from "Loopless Generation of + // Multiset Permutations using a Constant Number of Variables by Prefix + // Shifts." ) + struct coolmulti { + // element of singly linked list + struct element { + unsigned value; + element* next; + element(unsigned val, element* n) + : value(val), next(n) {} + ~element() + { // recurses down to the end of the singly linked list + delete next; + } + }; + element *head, *i, *after_i; + // NB: Partition must be sorted in non-decreasing order. + explicit coolmulti(const std::vector& partition) + : head(nullptr), i(nullptr), after_i(nullptr) + { + for (unsigned n = 0; n < partition.size(); ++n) { + head = new element(partition[n], head); + if (n <= 1) + i = head; + } + after_i = i->next; + } + ~coolmulti() + { // deletes singly linked list + delete head; + } + void next_permutation() + { + element *before_k; + if (after_i->next != nullptr && i->value >= after_i->next->value) + before_k = after_i; + else + before_k = i; + element *k = before_k->next; + before_k->next = k->next; + k->next = head; + if (k->value < head->value) + i = k; + after_i = i->next; + head = k; + } + bool finished() const + { + return after_i->next == nullptr && after_i->value >= head->value; + } + } 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 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()), current_updated(false) + { + for (unsigned i=1; i& get() const + { + 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; + } + bool next() + { + // This ugly contortion is needed because the original coolmulti + // algorithm requires code duplication of the payload procedure, + // one before the loop and one inside it. + if (trivial || atend) + return false; + cmgen.next_permutation(); + current_updated = false; + atend = cmgen.finished(); + return true; + } +}; + +/** Compute the multinomial coefficient n!/(p1!*p2!*...*pk!) where + * n = p1+p2+...+pk, i.e. p is a partition of n. + */ +const numeric +multinomial_coefficient(const std::vector & p); + // Collection of `construct on first use' wrappers for safely avoiding // internal object replication without running into the `static @@ -391,21 +604,7 @@ extern const ex _ex120; // Helper macros for class implementations (mostly useful for trivial classes) #define DEFAULT_CTOR(classname) \ -classname::classname() : inherited(TINFO_##classname) { setflag(status_flags::evaluated | status_flags::expanded); } - -#define DEFAULT_UNARCHIVE(classname) \ -ex classname::unarchive(const archive_node &n, lst &sym_lst) \ -{ \ - return (new classname(n, sym_lst))->setflag(status_flags::dynallocated); \ -} - -#define DEFAULT_ARCHIVING(classname) \ -classname::classname(const archive_node &n, lst &sym_lst) : inherited(n, sym_lst) {} \ -DEFAULT_UNARCHIVE(classname) \ -void classname::archive(archive_node &n) const \ -{ \ - inherited::archive(n); \ -} +classname::classname() { setflag(status_flags::evaluated | status_flags::expanded); } #define DEFAULT_COMPARE(classname) \ int classname::compare_same_type(const basic & other) const \ @@ -429,5 +628,4 @@ void classname::do_print_latex(const print_latex & c, unsigned level) const \ } // namespace GiNaC - -#endif // ndef __GINAC_UTILS_H__ +#endif // ndef GINAC_UTILS_H