C++ Arbitrary Precision Whole Number Class

Question

The class below implements an arbitrarily sized whole number C++ class.

It is wrapped around a bitvector class, which does most the binary heavy lifting.

The basic mathematical operators are implemented using whole word manipulation, inspired from Modern Computer Arithmetic - Paul Zimmermann.

I was unable to understand the pseudo-code and math notation for division from that book. So, I came up with my own algorithm. Feedback on that is very welcomed. I would be surprised if it couldn’t be optimized in some way.

A design change was made from the previous bitvector class, to remove the std::string and std::stringstream aliases as recommended here. After thinking about it if I make the changes now and cascade them up as I revise the rest of the library, it will be more manageable, than later going top to bottom.

This class is designed to abstract the arithmetic for arbitrary precision integrals, which will using this class in their design. Hence the boolean error variable in one of the constructors.

whole_number.h

#include "../support_files/string_support.h"
#include "../BitVector/bitvector.h"

namespace Olly {
    namespace MPA {

        /********************************************************************************************/
        //
        //                              'whole_number' Class Declaration
        //
        /********************************************************************************************/

        class whole_number {

            using Word        = bitvector::word_t;
            using Double_Word = bitvector::double_word_t;

        public:

            whole_number();
            whole_number(Word value);
            whole_number(const std::string& value, Word base = 10);
            whole_number(const std::string& value, Word base, bool& error);
            virtual ~whole_number();

            whole_number(whole_number&& obj)                 = default;
            whole_number(const whole_number& obj)            = default;
            whole_number& operator=(whole_number&& obj)      = default;
            whole_number& operator=(const whole_number& obj) = default;

            friend void swap(whole_number& first, whole_number& second);

            operator bool() const;

            bool is() const;

            bool is_odd()  const;
            bool is_even() const;

            bool operator==(const whole_number& b) const;
            bool operator!=(const whole_number& b) const;

            std::partial_ordering operator<=>(const whole_number& b) const;

            whole_number& operator&=(const whole_number& other);
            whole_number& operator|=(const whole_number& other);
            whole_number& operator^=(const whole_number& other);

            whole_number& operator<<=(std::size_t index);
            whole_number& operator>>=(std::size_t index);

            whole_number bin_comp() const;

            whole_number operator&(const whole_number& b) const;
            whole_number operator|(const whole_number& b) const;
            whole_number operator^(const whole_number& b) const;
            whole_number operator~() const;

            whole_number operator<<(std::size_t index) const;
            whole_number operator>>(std::size_t index) const;

            whole_number& operator+=(const whole_number& other);
            whole_number& operator-=(const whole_number& other);
            whole_number& operator*=(const whole_number& other);
            whole_number& operator/=(const whole_number& other);
            whole_number& operator%=(const whole_number& other);

            whole_number operator+(const whole_number& b) const;
            whole_number operator-(const whole_number& b) const;
            whole_number operator*(const whole_number& b) const;
            whole_number operator/(const whole_number& b) const;
            whole_number operator%(const whole_number& b) const;

            friend whole_number operator&(whole_number&& a, const whole_number& b);
            friend whole_number operator|(whole_number&& a, const whole_number& b);
            friend whole_number operator^(whole_number&& a, const whole_number& b);

            friend whole_number operator<<(whole_number&& a, std::size_t index);
            friend whole_number operator>>(whole_number&& a, std::size_t index);

            friend whole_number operator+(whole_number&& a, const whole_number& b);
            friend whole_number operator-(whole_number&& a, const whole_number& b);
            friend whole_number operator*(whole_number&& a, const whole_number& b);
            friend whole_number operator/(whole_number&& a, const whole_number& b);
            friend whole_number operator%(whole_number&& a, const whole_number& b);

            whole_number& operator++();
            whole_number  operator++(int);

            whole_number& operator--();
            whole_number  operator--(int);

            void div_rem(const whole_number& other, whole_number& qot, whole_number& rem) const;

            whole_number pow(std::size_t b) const;

            whole_number sqrt() const;
            whole_number root(const whole_number& b) const;

            std::string to_string()                 const;
            std::string to_string(std::size_t base) const;

            template<typename N>
            N to_integral() const;

            const bitvector& get_bitvector() const;

        private:

            bitvector _bitvec;

            whole_number(const bitvector& bitvec);

            void trim();

            void divide_remainder(Word d, whole_number& q, whole_number& r, int stop) const;

            bool set_numeric_value(const std::string& text, const Word& base);
        };

        template<typename N>
        inline N whole_number::to_integral() const {
            return _bitvec.to_integral<N>();
        }
    }
}

whole_number.cpp

#include "whole_number.h"

namespace Olly {
    namespace MPA {

        whole_number::whole_number() : _bitvec() {
        }

        whole_number::whole_number(Word value) : _bitvec(1, value) {
        }

        whole_number::whole_number(const std::string& value, Word base) : _bitvec() {
            set_numeric_value(value, base);
        }

        whole_number::whole_number(const std::string& value, Word base, bool& error) : _bitvec() {
            error = set_numeric_value(value, base);
        }

        whole_number::whole_number(const bitvector& bitvec) : _bitvec(bitvec) {
        }

        whole_number::~whole_number() {
        }

        void swap(whole_number& left, whole_number& right) {

            whole_number temp = std::move(left);

            left = std::move(right);
            right = std::move(temp);
        }

        whole_number::operator bool() const {
            return is();
        }

        bool whole_number::is() const {
            return _bitvec.is();
        }

        bool whole_number::is_odd() const {
            return _bitvec.last_word() & 1;
        }

        bool whole_number::is_even() const {
            return !is_odd();
        }

        bool whole_number::operator==(const whole_number& b) const {
            return operator<=>(b) == std::partial_ordering::equivalent;
        }

        bool whole_number::operator!=(const whole_number& b) const {
            return operator<=>(b) != std::partial_ordering::equivalent;
        }

        std::partial_ordering whole_number::operator<=>(const whole_number& b) const {
            return _bitvec <=> b._bitvec;
        }

        whole_number& whole_number::operator&=(const whole_number& other) {
            _bitvec &= other._bitvec;

            trim();

            return *this;
        }

        whole_number& whole_number::operator|=(const whole_number& other) {
            _bitvec |= other._bitvec;

            trim();

            return *this;
        }

        whole_number& whole_number::operator^=(const whole_number& other) {
            _bitvec ^= other._bitvec;

            trim();

            return *this;
        }

        whole_number& whole_number::operator<<=(std::size_t index) {
            _bitvec <<= index;

            trim();

            return *this;
        }

        whole_number& whole_number::operator>>=(std::size_t index) {
            _bitvec >>= index;

            trim();

            return *this;
        }

        whole_number whole_number::bin_comp() const {

            whole_number n;

            n._bitvec = _bitvec.bin_comp();

            return n;
        }

        whole_number whole_number::operator&(const whole_number& b) const {

            whole_number a = *this;

            a &= b;

            return a;
        }

        whole_number whole_number::operator|(const whole_number& b) const {

            whole_number a = *this;

            a |= b;

            return a;
        }

        whole_number whole_number::operator^(const whole_number& b) const {

            whole_number a = *this;

            a ^= b;

            return a;
        }

        whole_number whole_number::operator~() const {

            whole_number a;

            a._bitvec = ~_bitvec;

            return a;
        }

        whole_number whole_number::operator<<(std::size_t index) const {

            whole_number a = *this;

            a <<= index;

            return a;
        }

        whole_number whole_number::operator>>(std::size_t index) const {

            whole_number a = *this;

            a >>= index;

            return a;
        }

        whole_number& whole_number::operator+=(const whole_number& other) {

            std::size_t limit = _bitvec.size_words() > other._bitvec.size_words() ? _bitvec.size_words() : other._bitvec.size_words();

            Double_Word n = 0;

            for (std::size_t i = 0; i < limit; i += 1) {

                n = n + _bitvec.at_word(i) + other._bitvec.at_word(i);

                _bitvec.at_word(i) = static_cast<Word>(n);

                n >>= _bitvec.value_type;
            }

            if (n != 0) {

                _bitvec.at_word(limit) = static_cast<Word>(n);
            }

            trim();

            return *this;
        }

        whole_number& whole_number::operator-=(const whole_number& other) {

            if (other > *this) {
                *this = whole_number();
                return *this;
            }

            std::size_t limit = _bitvec.size_words() > other._bitvec.size_words() ? _bitvec.size_words() : other._bitvec.size_words();

            Double_Word n = 0;

            for (std::size_t i = 0; i < limit; i += 1) {

                n = n + _bitvec.at_word(i) - other._bitvec.at_word(i);

                _bitvec.at_word(i) = static_cast<Word>(n);

                n = ((n >> _bitvec.value_type) ? -1 : 0);
            }

            trim();

            return *this;
        }

        whole_number& whole_number::operator*=(const whole_number& other) {

            *this = *this * other;

            return *this;
        }

        whole_number& whole_number::operator/=(const whole_number& other) {

            *this = *this / other;

            return *this;
        }

        whole_number& whole_number::operator%=(const whole_number& other) {

            *this = *this % other;

            return *this;
        }

        whole_number whole_number::operator+(const whole_number& b) const {

            whole_number a = *this;

            a += b;

            return a;
        }

        whole_number whole_number::operator-(const whole_number& b) const {

            whole_number a = *this;

            a -= b;

            return a;
        }

        whole_number whole_number::operator*(const whole_number& b) const {

            std::size_t size_a = _bitvec.size_words();
            std::size_t size_b = b._bitvec.size_words();

            bitvector r((size_a + size_b + 1), 0);

            for (std::size_t j = 0; j < size_b; j += 1) {

                Double_Word n = 0;

                for (std::size_t i = 0; i < size_a; i += 1) {

                    std::size_t k = i + j;

                    n += static_cast<Double_Word>(_bitvec.at_word(i)) * b._bitvec.at_word(j) + r.at_word(k);

                    r.at_word(k) = static_cast<Word>(n);

                    n >>= _bitvec.value_type;
                }
                r.at_word(j + size_a) = static_cast<Word>(n);
            }

            r.trim();

            return r;
        }

        whole_number whole_number::operator/(const whole_number& b) const {

            whole_number q;
            whole_number r;

            div_rem(b, q, r);

            return q;
        }

        whole_number whole_number::operator%(const whole_number& b) const {

            whole_number q;
            whole_number r;

            div_rem(b, q, r);

            return r;
        }

        whole_number operator&(whole_number&& a, const whole_number& b) {
            return a &= b;
        }

        whole_number operator|(whole_number&& a, const whole_number& b) {
            return a |= b;
        }

        whole_number operator^(whole_number&& a, const whole_number& b) {
            return a ^= b;
        }

        whole_number operator<<(whole_number&& a, std::size_t index) {
            return a <<= index;
        }

        whole_number operator>>(whole_number&& a, std::size_t index) {
            return a >>= index;
        }

        whole_number operator+(whole_number&& a, const whole_number& b) {
            return a += b;
        }

        whole_number operator-(whole_number&& a, const whole_number& b) {
            return a -= b;
        }

        whole_number operator*(whole_number&& a, const whole_number& b) {
            return a *= b;
        }

        whole_number operator/(whole_number&& a, const whole_number& b) {
            return a /= b;
        }

        whole_number operator%(whole_number&& a, const whole_number& b) {
            return a %= b;
        }

        whole_number& whole_number::operator++() {

            ++_bitvec;

            trim();

            return *this;
        }

        whole_number whole_number::operator++(int) {

            whole_number a(*this);

            operator++();

            return a;
        }

        whole_number& whole_number::operator--() {

            --_bitvec;

            trim();

            return *this;
        }

        whole_number whole_number::operator--(int) {

            whole_number a(*this);

            operator--();

            return a;
        }

        void whole_number::div_rem(const whole_number& other, whole_number& qot, whole_number& rem) const {

            // Ensure both the qotient and remander are initalized to zero.
            qot = whole_number();
            rem = whole_number();

            if (!other.is()) {
                // Division by zero.
                return;
            }

            if (other > *this) {
                // Division by a greater value.
                qot = whole_number();
                rem = *this;
                return;
            }

            if (*this > other) {

                if (other._bitvec.size_words() == 1) {

                    divide_remainder(other._bitvec.at_word(0), qot, rem, 0);

                    return;
                }

                int stop = static_cast<int>(other._bitvec.size_words() - 1);

                auto d = other._bitvec.at_word(stop);

                for (int i = stop - 1; i >= 0; i -= 1) {
                    // Add one to 'd' if any other digits are defined.
                    if (other._bitvec.at_word(i) != 0) {
                        d += 1;
                        break;
                    }
                }

                // Perform long division.
                whole_number n = *this;

                whole_number q = whole_number();

                whole_number guard = whole_number();

                while (n >= other && n != guard) {

                    n.divide_remainder(d, q, rem, stop);

                    qot += q;

                    guard = n;

                    n -= (other * q);

                    q = whole_number();
                }

                // Confirm the qoutent is correct.
                q = other * qot;

                while (q < *this) {
                    qot += whole_number(1);
                    q += other;
                }

                while (q > *this) {
                    qot -= whole_number(1);
                    q -= other;
                }

                // Determine the remainder.
                rem = *this - q;
                return;
            }

            // Return  division by two equal values.
            qot = whole_number(1);
            rem = whole_number();
        }

        whole_number whole_number::pow(std::size_t b) const {

            if (b == 2) {
                whole_number a = *this;
                return a * a;
            }

            if (b == 1) {
                return *this;
            }

            if (b == 0) {
                return 1;
            }

            whole_number a = *this;
            whole_number res = 1;

            while (b) {

                if (b & 1) {

                    res *= a;
                }

                b >>= 1;

                if (b) {
                    a *= a;
                }
            }

            return res;
        }

        whole_number whole_number::sqrt() const {
            return static_cast<Word>(_bitvec.lead_bit() - 1);
        }

        whole_number whole_number::root(const whole_number& b) const {

            std::size_t n = b.to_integral<std::size_t>();

            whole_number low = 0;
            whole_number high = 1;

            whole_number ONE = 1;
            whole_number TWO = 2;

            while (high.pow(n) <= *this) {
                low = high;
                high *= TWO;
            }

            while (low != high - ONE) {

                whole_number step = (high - low) / TWO;

                whole_number candidate = low + step;

                auto value = candidate.pow(n);

                if (value == *this) {
                    return candidate;
                }
                if (value < *this) {
                    low = candidate;
                }
                else {
                    high = candidate;
                }
            }

            return low;
        }

        std::string whole_number::to_string() const {
            return to_string(10);
        }

        std::string whole_number::to_string(std::size_t base) const {

            if (base == 10 || base == 0 || base == 2 || base == 8 || base == 16) {

                std::string result = "";

                if (!is()) {
                    return "0";
                }

                whole_number radix = static_cast<Word>(base != 0 ? base : 10);
                whole_number n = *this;

                std::stringstream stream;

                int count = 0;

                while (n.is()) {

                    whole_number q;
                    whole_number r = n;

                    n.div_rem(radix, q, r);

                    n = q;

                    if (base == 16) {
                        stream << std::hex << r._bitvec.at_word(0);
                    }
                    else if (base == 8) {
                        stream << std::oct << r._bitvec.at_word(0);
                    }
                    else {
                        stream << r._bitvec.at_word(0);
                    }

                    if (base == 10) {
                        count += 1;

                        if (count == 3) {
                            stream << ',';
                            count = 0;
                        }
                    }
                }

                std::string res = stream.str();

                if (res.back() == ',') {
                    res.pop_back();
                }

                for (auto i = res.crbegin(); i != res.crend(); ++i) {

                    result += *i;
                }

                return result;
            }

            return "";
        }

        const bitvector& whole_number::get_bitvector() const {
            return _bitvec;
        }

        void whole_number::divide_remainder(Word d, whole_number& q, whole_number& r, int stop) const {

            Double_Word n(0);

            for (int i = static_cast<int>(_bitvec.size_words() - 1); i >= stop; i -= 1) {

                n += _bitvec.at_word(i);

                q._bitvec.at_word(static_cast<std::size_t>(i) - stop) = static_cast<Word>(n / d);

                n %= d;
                n <<= _bitvec.value_type;
            }

            r._bitvec.at_word(0) = n >> _bitvec.value_type;

            return;
        }

        bool whole_number::set_numeric_value(const std::string& text, const Word& base) {

            if (base == 10) {  // Parse a decimal number.

                Word x = 0;

                for (const auto n : text) {

                    if (!std::isspace(n) && n != ',') {

                        x = n - '0';

                        if (x >= 0 && x < base) {

                            operator*=(base);
                            operator+=(x);
                        }
                        else {
                            _bitvec = bitvector();
                            return true;
                        }
                    }
                }
            }

            else if (base == 16) {  // Parse a hexidecimal number.

                Word x;

                for (const auto n : text) {

                    if (!std::isspace(n)) {

                        switch (n) {

                        case('a'):
                        case('A'):
                            x = 10;
                            break;
                        case('b'):
                        case('B'):
                            x = 11;
                            break;
                        case('c'):
                        case('C'):
                            x = 12;
                            break;
                        case('d'):
                        case('D'):
                            x = 13;
                            break;
                        case('e'):
                        case('E'):
                            x = 14;
                            break;
                        case('f'):
                        case('F'):
                            x = 15;
                            break;
                        default:
                            x = (n - '0');
                        }

                        if (x >= 0 && x < base) {

                            operator*=(base);
                            operator+=(x);
                        }
                        else {
                            _bitvec = bitvector();
                            return true;
                        }
                    }
                }
            }

            else if (base == 8) {  // Parse an octal number.

                Word x;

                for (const auto n : text) {

                    if (!std::isspace(n)) {

                        x = n - '0';

                        if (x >= 0 && x < base) {

                            operator*=(base);
                            operator+=(x);
                        }
                        else {
                            _bitvec = bitvector();
                            return true;
                        }
                    }
                }
            }

            else if (base == 2) {  // Parse a binary number.

                for (const auto n : text) {

                    operator<<=(1);

                    if (n == '1') {
                        operator+=(1);
                    }
                    else if (n != '0') {
                        _bitvec = bitvector();
                        return true;
                    }
                }
            }

            return false;
        }

        void whole_number::trim() {
            _bitvec.trim();
        }
    }
}

Answer 1

Interface

I think the constructor that can be invoked with a std::string should be explicit. We should have implicit conversion only from integer types. And consider accepting std::string_view to avoid converting other string-like types to std::string.

The only virtual function is the destructor. That's a likely sign that we're not intending to derive from this class and own it as a pointer-to-base, so we can make it non-virtual and eliminate the need for a vtable.

operator!= can be defaulted.

Function names is() and bin_func() are not self-explanatory - probably need comments.

We don't need the member versions of binary operations, given the friend functions that implement them. Just pass the first argument by value.

Is operator~() really useful with an arbitrary-width type?

Consider naming the square-root function isqrt() - that's more explicit that we find only the integer part, and is a commonly-used name for this.

Consider providing character stream << and >> operators, rather than requiring construction of std::string objects for this.

to_string() could be a single function with default base.

Do we really need get_bitvector? If we eliminate this, we would then be free to change the internal representation.

Consider using the platform's native word size for arithmetic - GCC has __builtin_*_overflow() functions to help with this, as do many other platforms. Fall back to the usual overflow computation for unknown targets.

Implementation

It's hard to review the implementation without the bitvector interface.

There's a lot more blank lines than I would use - it almost looks double-spaced.

The member swap() exactly duplicates std::swap. It would be better implemented as std::swap(left.bitvec, right.bitvec) - or perhaps bitvector provides a better swap() that could be used instead?

Binary operators implemented in terms of assignment operators can be simplified, by passing first argument by value. That's done correctly in the friend functions, but not in the member versions which I've already said should be removed.

It's not clear why +, ++, *, etc. require trim() to be called - isn't that necessary only for functions which might reduce the length?

Division by zero gives results that are indistinguishable from division of zero. We need a better way to signal this error.

sqrt() implementation looks strange - unless bitvector::lead_bit() does something surprising. I'd like to see the tests for this one.

root() should probably just accept a std::size_t, rather than misleading the user. The values ONE and TWO should probably be const (and lower-case, since we normally use all-caps to draw attention to macro expansions).

to_string could do the error return first, and remove an indentation for most of the function. We should be inserting , (or .) only if the locale specifies a thousands separator. What's the return type of bitvector::at_word()? It's not clear that streaming from this will emit a single digit. In particular, how does stream << r._bitvec.at_word(0); work for both binary and decimal? That's another case where seeing the unit tests would be helpful.

div_rem() is expensive for the exact power-of-two bases - but sufficiently general that we can support more than just the specified set. I think we should use >>= for bases 2, 8 and 16.

set_numeric_value() passes plain char to std::isspace() - that's UB for negative values. Octal and hex input could use <<= instead of *= for better performance.

Unit tests

... are completely absent. I can't approve this code without them.

Answer 2

I gave some feedback on an earlier version of this, but there was one thing that I didn’t go into detail about.

The Double_Word Type

Currently, you use this in four places: addition, subtraction, multiplication and division. You’ve moved the definition (which previously had some problems) into a different file that you didn’t show us.

Addition and subtraction don’t technically need this type, although it’s likely faster than the alternatives. If the result of an ddition is larger than either of its addends, it overflowed, and you can carry one bit.

Division doesn’t need it either; you can do unsigned long division using only / and %, although some architectures might have speed-ups.

So that leaves multiplication, where you need some way to get the high word of the product. You have a few imperfect options:

• In practical, real-world use, all compilers provide uint32_t and uint64_t. This is an especially good option if you want a simple learning exercise and will not use this in production. However,
• It would be more efficient to store data in chunks of the native word size of the target machine. In the real world, size_t holds an unsigned machine word and operations on them will be fast. A good way to test how wide this is in the preprocessor is the constant SIZE_MAX in <stdint.h>.
• When size_t is larger than 32 bits (#if SIZE_MAX > 0xFFFFFFFFUL), there isn’t a standard type guaranteed to be large enough to hold the product of two of them. However, most compilers offer extensions that solve this:
  • uintmax_t is theoretically supposed to be the widest supported unsigned integral type (although this is frequently not true), so if any type can hold the square of SIZE_MAX, it might be this. And you can check portably.
  • GCC, Clang and ICX have an unsigned __int128 type.
  • MSVC for 64-bit targets has a __umulh intrinsic function in <intrin.h> that returns the upper word of the product of two unsigned 64-bit numbers, and others for add-with carry. You could use these as non-portable alternatives to a double word.
  • Some compilers have uint128_t as an extension, which might be added to the Standard in the future.