使用基2^32 ( up 256)跟踪的巨型整数类

Question

再次感谢大家对大型整数类我几天前在这里贴的的非常有用的建议。我已经做了许多建议的修改，也许是执行中最大的改动。特别是，我现在使用一个uint32_t类型的数组来存储数字，每个成员代表一个基本的2^{32}数字。这取代了以前版本中使用的uint8_t。基数补码现在是基- 2^{32}补码.这种类型的改变大大提高了速度性能。例如，在以前的版本中，计算1000!在我的笔记本上花费了大约43秒。现在还不到1秒就可以完成了！我没有像建议的那样让类成为模板，但这很容易做到。我没有让成千上万的分离器坚持当地的做法。我承认这一点的重要性，应该这样做，但现在，我只需要一些对数字人来说有意义的东西。

其他变动包括：

1. 该类位于自己的命名空间中。
2. 非类内联实用程序函数现在位于匿名名称空间中，以限制对翻译单元(文件)的访问。
3. 转换构造函数现在从long long int创建/转换，而不是从long int创建/转换，以使代码对其他操作系统更加友好。
4. 从C字符串中清除构造函数的代码，后者现在验证所提供的字符串。(我没有将类型更改为std::string，因为提供的字符串有99%是文字的。)
5. 修正了isZero()和toRawString()中未定义的(好的，简单的)行为。感谢那些指出这一点的人。
6. 添加静态成员函数getMinimum()和getMaximum()，分别返回基2^{32}方案中可表示的整数的最小值和最大值。
7. 添加了一个成员函数numDecimalDigits()，它返回*this的十进制表示形式中的小数位数。
8. 在算术运算符+、-等的内循环中增加显式转换，以确保运算数在执行加法、乘法等操作数之前被提升为uint64_t，从而正确传播。
9. 更新了驱动程序，它现在运行得更快，计算更多的数字。
10. 其他轻微的代码美学变化。

默认实现将numDigits硬编码为300，这将产生2890小数位数。如果你的硬件能处理的话，你可以比这个高得多。但是，到long double的转换将悄然中断(如果您确实需要的话)输出NaNs，除非您采取步骤将CPU置于对这些错误有错误的状态。

再次感谢所有对这些变化作出贡献的人。这个项目产生了一个我认为非常有用的工人阶级，我对此感到非常鼓舞。

新代码如下。我很感谢你对改进的意见和建议。谢谢您抽时间见我。

HugeInt.h

/*
 * HugeInt.h
 * 
 * Definition of the huge integer class
 * February, 2020
 * 
 * RADIX 2^32 VERSION
 *
 * Huge integers are represented as N-digit arrays of uint32_t types, where
 * each uint32_t value represents a base-2^32 digit. By default N = 300, which 
 * corresponds to a maximum of 2890 decimal digits. Each uint32_t contains 
 * a single base-2^32 digit in the range 0 <= digit <= 2^32 - 1. If `index' 
 * represents the index of the array of uint32_t digits[N], 
 * i.e., 0 <= index <= N - 1, and 'value' represents the power of 2^32 
 * corresponding to the radix 2^32 digit at 'index', then we have the following 
 * correspondence:
 *
 * index  |...... |     4    |     3    |     2    |     1    |     0    |
 * -----------------------------------------------------------------------
 * value  |...... | (2^32)^4 | (2^32)^3 | (2^32)^2 | (2^32)^1 | (2^32)^0 |
 *
 * The physical layout of the uint32_t array in memory is:
 *
 * uint32_t digits[N] = {digits[0], digits[1], digits[2], digits[3], ... }
 *
 * which means that the units (2^32)^0 appear first in memory, while the power 
 * (2^32)^(N-1) appears last. This LITTLE ENDIAN storage represents the 
 * number in memory in the REVERSE order of the way we write decimal numbers, 
 * but is convenient.
 *
 * Negative integers are represented by their radix complement. With the 
 * base 2^32 implementation here, we represent negative integers by their base 
 * 2^32 complement. With this convention the range of 
 * non-negative integers is:
 *                      0 <= x <= (2^32)^N/2 - 1
 * The range of base 2^32 integers CORRESPONDING to negative values in the
 * base 2^32 complement scheme is:
 *                      (2^32)^N/2 <= x <= (2^32)^N - 1 
 * So -1 corresponds to (2^32)^N - 1, -2 corresponds to (2^32)^N - 2, and so on.
 * 
 * The complete range of integers represented by a HugeInt using radix 
 * complement is:
 * 
 *                     -(2^32)^N/2 <= x <= (2^32)^N/2 - 1
 */


#ifndef HUGEINT_H
#define HUGEINT_H

#include <string>
#include <iostream>

namespace iota {

class HugeInt {
public:
    HugeInt() = default;
    HugeInt(const long long int);  // conversion constructor from long long int
    explicit HugeInt(const char *const); // conversion constructor from C string
    HugeInt(const HugeInt&);    // copy/conversion constructor

    // assignment operator
    const HugeInt& operator=(const HugeInt&);

    // unary minus operator
    HugeInt operator-() const;

    // conversion to long double
    explicit operator long double() const;

    // basic arithmetic
    friend HugeInt operator+(const HugeInt&, const HugeInt&);
    friend HugeInt operator-(const HugeInt&, const HugeInt&);
    friend HugeInt operator*(const HugeInt&, const HugeInt&);
  //  friend HugeInt operator/(const HugeInt&, const HugeInt&); // TODO:

    // increment and decrement operators
    HugeInt& operator+=(const HugeInt&);
    HugeInt& operator-=(const HugeInt&);
    HugeInt& operator*=(const HugeInt&);
   // HugeInt& operator/=(const HugeInt&); TODO:
    HugeInt& operator++();     // prefix
    HugeInt  operator++(int);  // postfix
    HugeInt& operator--();     // prefix
    HugeInt  operator--(int);  // postfix

    // relational operators
    friend bool operator==(const HugeInt&, const HugeInt&);
    friend bool operator!=(const HugeInt&, const HugeInt&);
    friend bool operator<(const HugeInt&, const HugeInt&);
    friend bool operator>(const HugeInt&, const HugeInt&);
    friend bool operator<=(const HugeInt&, const HugeInt&);
    friend bool operator>=(const HugeInt&, const HugeInt&);

    // input/output 
    std::string toRawString() const;
    std::string toDecimalString() const;
    friend std::ostream& operator<<(std::ostream&, const HugeInt&);
    friend std::istream& operator>>(std::istream&, HugeInt&);

    // informational
    int numDecimalDigits() const;
    static HugeInt getMinimum();
    static HugeInt getMaximum();

private:
    static const int      numDigits{300};       // max. no. radix 2^32 digits
    static const uint64_t twoPow32{4294967296}; // 2^32, for convenience
    uint32_t              digits[numDigits]{0}; // radix 2^32 digits; default 0

    // private utility functions
    bool     isZero() const;
    bool     isNegative() const;
    HugeInt& radixComplementSelf();  
    HugeInt  radixComplement() const;
    HugeInt  shortDivide(uint32_t) const;
    uint32_t shortModulo(uint32_t) const;
    HugeInt  shortMultiply(uint32_t) const;
    HugeInt& shiftLeftDigits(int);
};

} /* namespace iota */

#endif /* HUGEINT_H */

HugeInt.cpp

/*
 * HugeInt.cpp
 *
 * Implementation of the HugeInt class. See comments in HugeInt.h for
 * details of representation, etc.
 *
 * February, 2020
 *
 * RADIX 2^32 VERSION
 * 
 */

#include <cstdlib>   // for abs(), labs(), etc.
#include <iostream>
#include <iomanip>
#include <sstream>
#include <cstring>
#include <stdexcept>
#include <cmath>
#include "HugeInt.h"

/*
 * Non-member utility functions (in anonymous namespace -- file scope only).
 * 
 */

namespace { /* anonymous namespace */

/*
 * Simple function to check for non-digit characters in a C string.
 * 
 * Returns true if string contains all digit characters; otherwise
 * false.
 * 
 */

inline bool validate_digits(const char *const str) {
    bool retval = true;

    for (size_t i = 0; i < std::strlen(str); ++i) {
        if (std::isdigit(static_cast<unsigned char>(str[i])) == 0) {
            retval = false;
            break;
        } 
    }

    return retval;
}

/**
 * get_carry32
 *
 * Return the high 32 bits of a 64-bit uint64_t.
 * Return this 32-bit value as a uint32_t.
 * 
 * @param value
 * @return 
 */

inline uint32_t get_carry32(uint64_t value) {
    return static_cast<uint32_t>(value >> 32 & 0xffffffff);
}

/**
 * get_digit32
 * 
 * Return the low 32 bits of the two-byte word stored as an int.
 * Return this 32-bit value as a uint32_t.
 * 
 * @param value
 * @return 
 */

inline uint32_t get_digit32(uint64_t value) {
    return static_cast<uint32_t>(value & 0xffffffff);
}

} /* anonymous namespace */



/*
 * Member functions in namespace iota
 *
 */

namespace iota {

/**
 * Constructor (conversion constructor)
 *
 * Construct a HugeInt from a long long int.
 *
 */ 

HugeInt::HugeInt(const long long int x) {
    if (x == 0LL) {
        return;
    }

    long long int xp{std::llabs(x)};
    int i{0};

    // Successively determine units, 2^32's, (2^32)^2's, (2^32)^3's, etc.
    // storing them in digits[0], digits[1], digits[2], ...,
    // respectively. That is units = digits[0], 2^32's = digits[1], etc.
    while (xp > 0LL) {
        digits[i++] = xp % twoPow32;
        xp /= twoPow32;
    }

    if (x < 0LL) {
        radixComplementSelf();
    }
}

/**
 * Constructor (conversion constructor)
 *
 * Construct a HugeInt from a null-terminated C string representing the
 * base 10 representation of the number. The string is assumed to have 
 * the form "[+/-]31415926", including an optional '+' or '-' sign. 
 *
 * WARNING: No spaces are allowed in the decimal string containing numerals.
 * 
 * 
 * @param str
 */

HugeInt::HugeInt(const char *const str) {
    const int len{static_cast<int>(std::strlen(str))};

    if (len == 0) {
        throw std::invalid_argument{"empty decimal string in constructor."};
    }

    // Check for explicit positive and negative signs and adjust accordingly.
    // If negative, we flag the case and perform a radix complement at the end.
    bool flagNegative{false};
    int  numDecimalDigits{len};
    int  offset{0};

    if (str[0] == '+') {
        --numDecimalDigits;
        ++offset;
    } 

    if (str[0] == '-') {
        flagNegative = true;
        --numDecimalDigits;
        ++offset;
    }

    // validate the string of numerals
    if (!validate_digits(str + offset)) {
        throw std::invalid_argument{"string contains non-digit in constructor."};
    }

    // Loop (backwards) through each decimal digit, digit[i], in the string, 
    // adding its numerical contribution, digit[i]*10^i, to theNumber. Here i 
    // runs upwards from zero, starting at the right-most digit of the string 
    // of decimal digits.
    uint32_t digitValue{0};
    HugeInt  theNumber{0LL};
    HugeInt  powerOfTen{1LL}; // initially 10^0 = 1

    for (int i = 0; i < numDecimalDigits; ++i) {
        digitValue = static_cast<uint32_t>(str[len - 1 - i]) - '0';
        theNumber += powerOfTen.shortMultiply(digitValue);
        powerOfTen = powerOfTen.shortMultiply(10);
    }

    if (flagNegative) {
        theNumber.radixComplementSelf();
    }

    for (int i = 0; i < numDigits; ++i) {
        digits[i] = theNumber.digits[i];
    }
}

/**
 * Copy constructor (could be defaulted)
 * 
 * @param rhs
 */

HugeInt::HugeInt(const HugeInt& rhs) {
    // TODO: perhaps call copy assignment?
    for (int i = 0; i < numDigits; ++i)
        digits[i] = rhs.digits[i];
}

/**
 * Assignment operator
 * 
 * @param rhs
 * @return 
 */

const HugeInt& HugeInt::operator=(const HugeInt& rhs) {
    if (&rhs != this) {
        for (int i = 0; i < numDigits; ++i) {
            digits[i] = rhs.digits[i]; 
        }
    }

    return *this;
}

/**
 * Unary minus operator
 * 
 * @return 
 */

HugeInt HugeInt::operator-() const {
    return radixComplement();
}

/**
 * operator long double() 
 *
 * Use with static_cast<long double>(hugeint) to convert hugeint to its
 * approximate (long double) floating point value.
 * 
 */
HugeInt::operator long double() const {
    long double sign{1.0L};
    HugeInt     copy{*this};

    if (copy.isNegative()) {
        copy.radixComplementSelf();
        sign = -1.0L;
    }

    long double retval{0.0L};
    long double pwrOfBase{1.0L}; // Base = 2^32; (2^32)^0 initially

    for (int i = 0; i < numDigits; ++i) {
        retval += copy.digits[i] * pwrOfBase;
        pwrOfBase *= twoPow32;
    }

    return retval*sign;
}

/**
 * Operator +=
 *
 * NOTE: With the conversion constructors provided, also
 *       provides operator+=(long int) and
 *                operator+=(const char *const)
 * 
 * @param increment
 * @return 
 */

HugeInt& HugeInt::operator+=(const HugeInt& increment) {
    *this = *this + increment;
    return *this;
}

/**
 * Operator -=
 * 
 * NOTE: With the conversion constructors provided, also
 *       provides operator-=(long int) and
 *                operator-=(const char *const)
 * 
 * 
 * @param decrement
 * @return 
 */

HugeInt& HugeInt::operator-=(const HugeInt& decrement) {
    *this = *this - decrement;
    return *this;
}

/**
 * Operator *=
 * 
 * NOTE: With the conversion constructors provided, also
 *       provides operator*=(long int) and
 *                operator*=(const char *const)
 * 
 * @param multiplier
 * @return 
 */

HugeInt& HugeInt::operator*=(const HugeInt& multiplier) {
    *this = *this * multiplier;
    return *this;
}

/**
 * Operator ++ (prefix)
 * 
 * @return 
 */

HugeInt& HugeInt::operator++() {
    *this = *this + 1LL;
    return *this;
}

/**
 * Operator ++ (postfix)
 * 
 * @param 
 * @return 
 */

HugeInt HugeInt::operator++(int) {
   HugeInt retval{*this};
   ++(*this);

   return retval;
}

/**
 * Operator -- (prefix)
 * 
 * @return 
 */

HugeInt& HugeInt::operator--() {
   *this = *this - 1LL;
   return *this;
}

/**
 * Operator -- (postfix)
 * 
 * @param 
 * @return 
 */

HugeInt HugeInt::operator--(int) {
   HugeInt retval{*this};
   --(*this);

   return retval;
}


////////////////////////////////////////////////////////////////////////////
// Input/Output                                                           //
////////////////////////////////////////////////////////////////////////////

/**
 * toRawString()
 * 
 * Format a HugeInt as string in raw internal format, i.e., as a sequence 
 * of base-2^32 digits (each in decimal form, 0 <= digit <= 2^32 - 1).
 *  
 * @return 
 */

std::string HugeInt::toRawString() const {
    int istart{numDigits - 1};

    for ( ; istart >= 0; --istart) {
        if (digits[istart] != 0) {
            break;
        }
    }

    std::ostringstream oss;

    if (istart == -1) // the number is zero
    {
        oss << digits[0];
    } else {
        for (int i = istart; i >= 0; --i) {
            oss << std::setw(10) << std::setfill('0') << digits[i] << " ";
        }
    }

    return oss.str();
}

/**
 * toDecimalString()
 * 
 * Format HugeInt as a string of decimal digits. The length of the decimal 
 * string is estimated (roughly) by solving for x:
 *
 *     (2^32)^N = 10^x    ==>    x = N log_10(2^32) = N * 9.63296 (approx)
 *
 * where N is the number of base 2^32 digits. A safety margin of 5 is added
 * for good measure.
 * 
 * @return 
 */

std::string HugeInt::toDecimalString() const {
    std::ostringstream oss;

    // Special case HugeInt == 0 is easy
    if (isZero()) {
        oss << "0";
        return oss.str();
    }

    // set copy to the absolute value of *this
    // for use in shortDivide and shortModulo
    HugeInt tmp;

    if (isNegative()) {
        oss << "-";
        tmp = this->radixComplement();
    } else {
        tmp = *this;
    }

    // determine the decimal digits of the absolute value 
    int       i{0};
    const int numDecimal{static_cast<int>(numDigits * 9.63296) + 5};
    uint32_t  decimalDigits[numDecimal]{0};

    while (!tmp.isZero()) {
        decimalDigits[i++] = tmp.shortModulo(10);
        tmp = tmp.shortDivide(10);
    }

    // output the decimal digits
    for (int j = i - 1; j >= 0; --j) {
        if (j < i - 1) {
            if ((j + 1) % 3 == 0) {
                oss << ','; // thousands separator
            }
        }

        oss << decimalDigits[j];
    }

    return oss.str();
}

/////////////////////////////////////////////////////////////////////////////
// Useful informational member functions                                   //
/////////////////////////////////////////////////////////////////////////////

/**
 * getMinimum()
 * 
 * Return the minimum representable value for a HugeInt. Static member
 * function.
 * 
 * @return 
 */

HugeInt HugeInt::getMinimum() {
    HugeInt retval;

    retval.digits[numDigits - 1] = 2147483648;

    return retval;
}

/**
 * getMaximum()
 * 
 * Return the maximum representable value for a HugeInt. Static member 
 * function.
 * 
 * @return 
 */

HugeInt HugeInt::getMaximum() {
    HugeInt retval;

    retval.digits[numDigits - 1] = 2147483648;

    --retval;

    return retval;
}

/**
 * numDecimalDigits()
 * 
 * Return the number of decimal digits this HugeInt has. 
 * 
 * We use a simple algorithm using base-10 logarithms. Consider, e.g., 457, 
 * which we can write as 4.57 * 10^2. Taking base-10 logs: 
 *         
 *          log10(4.57 * 10^2) = log10(4.57) + 2.
 * 
 * Since 0 < log10(4.57) < log10(10) = 1, we need to round up (always) to get
 * the extra digit, corresponding to the fractional part in the eq. above. 
 * Hence the use of ceil below. Values of x in the range -10 < x < 10 are dealt
 * with as a special case.
 * 
 * @return 
 */

int HugeInt::numDecimalDigits() const {

    if (-10 < *this && *this < 10) {
        return 1;
    }
    else {
        long double approx = static_cast<long double>(*this);
        return static_cast<int>(std::ceil(std::log10(std::fabs(approx))));
    }
}

////////////////////////////////////////////////////////////////////////////
// friend functions                                                       //
////////////////////////////////////////////////////////////////////////////

/**
 * friend binary operator +
 *
 * Add two HugeInts a and b and return c = a + b.
 *
 * Note: since we provide a conversion constructor for long long int's, this 
 *       function, in effect, also provides the following functionality by 
 *       implicit conversion of long long int's to HugeInt
 *
 *       c = a + <some long long int>    e.g.  c = a + 2412356LL
 *       c = <some long long int> + a    e.g.  c = 2412356LL + a
 * 
 * @param a
 * @param b
 * @return 
 */

HugeInt operator+(const HugeInt& a, const HugeInt& b) {
    uint32_t carry{0};
    uint64_t partial{0};
    HugeInt sum;

    for (int i = 0; i < HugeInt::numDigits; ++i) {
        partial =   static_cast<uint64_t>(a.digits[i]) 
                  + static_cast<uint64_t>(b.digits[i]) 
                  + static_cast<uint64_t>(carry);
        carry = get_carry32(partial);
        sum.digits[i] = get_digit32(partial);
    }

    return sum;
}

/**
 * friend binary operator-
 *
 * Subtract HugeInt a from HugeInt a and return the value c = a - b.
 *
 * Note: since we provide a conversion constructor for long long int's, this 
 *       function, in effect, also provides the following functionality by 
 *       implicit conversion of long long int's to HugeInt:
 *
 *       c = a - <some long long int>    e.g.  c = a - 2412356LL
 *       c = <some long long int> - a    e.g.  c = 2412356LL - a
 * 
 * @param a
 * @param b
 * @return 
 */

HugeInt operator-(const HugeInt& a, const HugeInt& b) {
    return a + (-b);
}

/**
 * friend binary operator *
 *
 * Multiply two HugeInt numbers. Uses standard long multipication algorithm
 * adapted to base 2^32. See comments on implicit conversion before 
 * HugeInt operator+(const HugeInt&, const HugeInt& ) above.
 * 
 * @param a
 * @param b
 * @return 
 */

HugeInt operator*(const HugeInt& a, const HugeInt& b) {
    HugeInt product{0LL};
    HugeInt partial;

    for (int i = 0; i < HugeInt::numDigits; ++i) {
        partial = a.shortMultiply(b.digits[i]);
        product += partial.shiftLeftDigits(i);
    }

    return product;
}

////////////////////////////////////////////////////////////////////////////
// Relational operators (friends)                                         //
////////////////////////////////////////////////////////////////////////////

/**
 * Operator ==
 * 
 * @param lhs
 * @param rhs
 * @return 
 */

bool operator==(const HugeInt& lhs, const HugeInt& rhs) {
   HugeInt diff{rhs - lhs};

   return diff.isZero();
}

/**
 * Operator !=
 * 
 * @param lhs
 * @param rhs
 * @return 
 */

bool operator!=(const HugeInt& lhs, const HugeInt& rhs) {
   return !(rhs == lhs);
}

/**
 * Operator <
 * 
 * @param lhs
 * @param rhs
 * @return 
 */

bool operator<(const HugeInt& lhs, const HugeInt& rhs) {
   HugeInt diff{lhs - rhs};

   return diff.isNegative();
}

/**
 * Operator >
 * 
 * @param lhs
 * @param rhs
 * @return 
 */

bool operator>(const HugeInt& lhs, const HugeInt& rhs) {
    return rhs < lhs;
}

/**
 * Operator <=
 * 
 * @param lhs
 * @param rhs
 * @return 
 */

bool operator<=(const HugeInt& lhs, const HugeInt& rhs) {
    return !(lhs > rhs);
}

/**
 * Operator >=
 * 
 * @param lhs
 * @param rhs
 * @return 
 */

bool operator>=(const HugeInt& lhs, const HugeInt& rhs) {
    return !(lhs < rhs);
}


////////////////////////////////////////////////////////////////////////////
// Private utility functions                                              //
////////////////////////////////////////////////////////////////////////////

/**
 * isZero()
 * 
 * Return true if the HugeInt is zero, otherwise false.
 * 
 * @return 
 */

bool HugeInt::isZero() const {
    int i{numDigits - 1};

    for ( ; i >= 0; --i) {
        if (digits[i] != 0) {
            break;
        }
    }

    return i == -1;
}

/**
 * isNegative()
 * 
 * Return true if a number x is negative (x < 0). If x >=0, then
 * return false.
 * 
 * NOTE: In the radix-2^32 complement convention, negative numbers, x, are 
 *       represented by the range of values: (2^32)^N/2 <= x <=(2^32)^N - 1.
 *       Since (2^32)^N/2 = (2^32/2)*(2^32)^(N-1) = 2147483648*(2^32)^(N-1), 
 *       we need only check whether the (N - 1)'th base 2^32 digit is at 
 *       least 2147483648. 
 * 
 * @return 
 */

bool HugeInt::isNegative() const {
    return digits[numDigits - 1] >= 2147483648;
}

/**
 * shortDivide:
 * 
 * Return the result of a base 2^32 short division by divisor, where 
 * 0 < divisor <= 2^32 - 1, using the usual primary school algorithm 
 * adapted to radix 2^32.
 *
 * WARNING: assumes both HugeInt and the divisor are POSITIVE.
 * 
 * @param divisor
 * @return 
 */

HugeInt HugeInt::shortDivide(uint32_t divisor) const {
    uint64_t j;
    uint64_t remainder{0};
    HugeInt quotient;

    for (int i = numDigits - 1; i >= 0; --i) {
        j = twoPow32 * remainder + static_cast<uint64_t>(digits[i]);
        quotient.digits[i] = static_cast<uint32_t>(j / divisor);
        remainder = j % divisor;
    }

    return quotient;
}

/**
 * shortModulo
 * 
 * Return the remainder of a base 2^32 short division by divisor, where 
 * 0 < divisor <= 2^32 - 1.
 *
 * WARNING: assumes both HugeInt and the divisor are POSITIVE.
 * 
 * @param divisor
 * @return 
 */

uint32_t HugeInt::shortModulo(uint32_t divisor) const {
    uint64_t j;
    uint64_t remainder{0};

    for (int i = numDigits - 1; i >= 0; --i) {
        j = twoPow32 * remainder + static_cast<uint64_t>(digits[i]);
        remainder = j % divisor;
    }

    return static_cast<uint32_t>(remainder);
}

/**
 * shortMultiply
 * 
 * Return the result of a base 2^32 short multiplication by multiplier, where
 * 0 <= multiplier <= 2^32 - 1.
 *
 * WARNING: assumes both HugeInt and multiplier are POSITIVE.
 * 
 * @param multiplier
 * @return 
 */

HugeInt HugeInt::shortMultiply(uint32_t multiplier) const {
    uint32_t carry{0};
    uint64_t tmp;
    HugeInt product;

    for (int i = 0; i < numDigits; ++i) {
        tmp = static_cast<uint64_t>(digits[i]) * multiplier + carry;
        carry = get_carry32(tmp);
        product.digits[i] = get_digit32(tmp);
    }

    return product;
}

/**
 * shiftLeftDigits
 *
 * Shift this HugeInt's radix-2^32 digits left by num places, filling
 * with zeroes from the right.
 * 
 * @param num
 * @return 
 */

HugeInt& HugeInt::shiftLeftDigits(int num) {
    if (num == 0) {
        return *this;
    }

    for (int i = numDigits - num - 1; i >= 0; --i) {
        digits[i + num] = digits[i];
    }

    for (int i = 0; i < num; ++i) {
        digits[i] = 0;
    }

    return *this;
}

/**
 * radixComplementSelf()
 *
 * Perform a radix complement on the object in place (changes object).
 * 
 * @return 
 */

HugeInt& HugeInt::radixComplementSelf() {
    if (!isZero()) {
        uint64_t sum{0};
        uint32_t carry{1};

        for (int i = 0; i < numDigits; ++i) {
            sum =   static_cast<uint64_t>(twoPow32 - 1) 
                  - static_cast<uint64_t>(digits[i]) 
                  + static_cast<uint64_t>(carry);
            carry = get_carry32(sum);
            digits[i] = get_digit32(sum);
        }
    }

    return *this;
}

/**
 * radixComplement()
 * 
 * Return the radix-2^32 (base-2^32) complement of HugeInt.
 * 
 * @return 
 */

HugeInt HugeInt::radixComplement() const {
    HugeInt result{*this};

    return result.radixComplementSelf();
}

/**
 * operator<<
 * 
 * Overloaded stream insertion for HugeInt.
 * 
 * @param output
 * @param x
 * @return 
 */

std::ostream& operator<<(std::ostream& output, const HugeInt& x) {
    output << x.toDecimalString();

    return output;
}

/**
 * operator >>
 * 
 * Overloaded stream extraction for HugeInt.
 * 
 * @param input
 * @param x
 * @return 
 */

std::istream& operator>>(std::istream& input, HugeInt& x) {
    std::string str;

    input >> str;
    x = HugeInt(str.c_str());

    return input;
}

} /* namespace iota */

示例驱动程序代码：

/*
 * Simple driver to test a few features of the HugeInt class.
 * 
 * Improved version.
 * 
 */

#include <iostream>
#include <iomanip>
#include <limits>
#include "HugeInt.h"

iota::HugeInt read_bounded_hugeint(const iota::HugeInt&, const iota::HugeInt&);
iota::HugeInt factorial_recursive(const iota::HugeInt&);
iota::HugeInt factorial_iterative(const iota::HugeInt&);
iota::HugeInt fibonacci_recursive(const iota::HugeInt&);
iota::HugeInt fibonacci_iterative(const iota::HugeInt&);
void preamble();

// limits to avoid overflow
const iota::HugeInt FACTORIAL_LIMIT{1100LL};
const iota::HugeInt FIBONACCI_LIMIT{13000LL};

int main() {

    preamble(); // blah

    iota::HugeInt nfac = read_bounded_hugeint(0LL, FACTORIAL_LIMIT);

    iota::HugeInt factorial = factorial_iterative(nfac);
    long double factorial_dec = static_cast<long double>(factorial);

    std::cout << "\nThe value of " << nfac << "! is:\n";
    std::cout << factorial << '\n';
    std::cout << "\nThis value has " << factorial.numDecimalDigits() 
              << " decimal digits.\n";
    std::cout.precision(std::numeric_limits<long double>::digits10);
    std::cout << "\nIts decimal approximation is: " << factorial_dec << "\n\n";


    iota::HugeInt nfib = read_bounded_hugeint(0LL, FIBONACCI_LIMIT);

    iota::HugeInt fibonacci = fibonacci_iterative(nfib);
    long double fibonacci_dec = static_cast<long double>(fibonacci);

    std::cout << "\nThe " << nfib << "th Fibonacci number is:\n";
    std::cout << fibonacci << '\n';
    std::cout << "\nThis value has " << fibonacci.numDecimalDigits() 
              << " decimal digits.\n";
    std::cout << "\nIts decimal approximation is: " << fibonacci_dec << '\n';

    std::cout << "\nComparing these two values we observe that ";
    if (factorial == fibonacci) {
        std::cout << nfac << "! == Fibonacci_{" << nfib << "}\n";
    }

    if (factorial < fibonacci) {
        std::cout << nfac << "! < Fibonacci_{" << nfib << "}\n";
    }

    if (factorial > fibonacci) {
        std::cout << nfac << "! > Fibonacci_{" << nfib << "}\n";
    }

    iota::HugeInt sum = factorial + fibonacci;
    iota::HugeInt diff = factorial - fibonacci;

    std::cout << "\nTheir sum (factorial + fibonacci) is:\n";
    std::cout << sum << '\n';
    std::cout << "\n\twhich is approximately " << static_cast<long double>(sum);
    std::cout << '\n';

    std::cout << "\nTheir difference (factorial - fibonacci) is:\n";
    std::cout << diff << '\n';
    std::cout << "\n\twhich is approximately " << static_cast<long double>(diff);
    std::cout << '\n';

    iota::HugeInt x{"-80538738812075974"};
    iota::HugeInt y{"80435758145817515"};
    iota::HugeInt z{"12602123297335631"};

    iota::HugeInt k = x*x*x + y*y*y + z*z*z;

    std::cout << "\nDid you know that, with:\n";
    std::cout << "\tx = " << x << '\n';
    std::cout << "\ty = " << y << '\n';
    std::cout << "\tz = " << z << '\n';
    std::cout << "\nx^3 + y^3 + z^3 = " << k << '\n';

    return 0;
}


/**
 * read_bounded_hugeint
 * 
 * Read and return a value in the range min <= value <= max.
 * Dies after 5 retries.
 * 
 * @param min
 * @param max
 * @return 
 */
iota::HugeInt read_bounded_hugeint(const iota::HugeInt& min, 
                                   const iota::HugeInt& max) {
    iota::HugeInt value;
    bool fail;
    int retries = 0;

    do {
        try {
            std::cout << "Enter an integer (" << min << " - " << max << "): "; 
            std::cin >> value;

            if (value < min || value > max) {
                fail = true;
                ++retries;
            }
            else {
                fail = false;
            }
        }
        catch (std::invalid_argument& error) {
            std::cout << "You entered an invalid HugeInt value.";
            std::cout << " Please use, e.g., [+/-]1234567876376763.\n";
            //std::cout << "Exception: " << error.what() << '\n';
            fail = true;
            ++retries;
        }  
    } while (fail && retries < 5);

    if (retries == 5) {
        std::cerr << "Giving up...\n";
        exit(EXIT_FAILURE);
    }

    return value;
}

/**
 * factorial_recursive:
 * 
 * Recursive factorial function using HugeInt. Not too slow.
 * 
 * @param n
 * @return 
 */

iota::HugeInt factorial_recursive(const iota::HugeInt& n) {
    const iota::HugeInt one{1LL};

    if (n <= one) {
        return one;
    } else {
        return n * factorial_recursive(n - one);
    }
}

iota::HugeInt factorial_iterative(const iota::HugeInt& n) {
    iota::HugeInt result{1LL};

    if (n == 0LL) {
        return result;
    }

    for (iota::HugeInt i = n; i >= 1; --i) {
        result *= i;
    }

    return result;
}

/**
 * fibonacci_recursive:
 * 
 * Recursively calculate the n'th Fibonacci number, where n>=0.
 * 
 * WARNING: S l o w . . .
 * 
 * @param n
 * @return 
 */
iota::HugeInt fibonacci_recursive(const iota::HugeInt& n) {
    const iota::HugeInt zero;
    const iota::HugeInt one{1LL};

    if ((n == zero) || (n == one)) {
        return n;
    } 
    else {
        return fibonacci_recursive(n - 1LL) + fibonacci_recursive(n - 2LL);
    }  
}

iota::HugeInt fibonacci_iterative(const iota::HugeInt& n) {
    const iota::HugeInt zero;
    const iota::HugeInt one{1LL};

    if ((n == zero) || (n == one)) {
        return n;
    }

    iota::HugeInt retval;
    iota::HugeInt fib_nm1 = one;
    iota::HugeInt fib_nm2 = zero;

    for (iota::HugeInt i = 2; i <= n; ++i) {
        retval = fib_nm1 + fib_nm2;
        fib_nm2 = fib_nm1;
        fib_nm1 = retval;
    }

    return retval;
}

void preamble() {
    long double min = static_cast<long double>(iota::HugeInt::getMinimum());
    long double max = static_cast<long double>(iota::HugeInt::getMaximum());

    std::cout.precision(std::numeric_limits<long double>::digits10);
    std::cout <<"**************************************************************"
              <<"*************\n\n";
    std::cout << "The range of integers, x, that can be represented in " 
              << "the default HugeInt\nconfiguration is, approximately\n"
              << "      " << min << " <= x <= " << max << '\n';

    std::cout << "\nThe precise values of the upper and lower limits can be "
              << "found using\nHugeInt::getMinimum()/HugeInt::getMaximum().\n";

    std::cout << "\nThe maximum number of decimal digits of an integer "
              << "representable with\na HugeInt is: " 
              << iota::HugeInt::getMaximum().numDecimalDigits()
              << "\n\n";   
    std::cout <<"**************************************************************"
              <<"*************\n\n";
}

Answer 1

我建议首先包括内部标头，然后再包含任何标准标头。这有助于暴露任何意外的依赖关系，使您很难在其他程序中使用您的类型。

所以在HugeInt.cpp里，我会写

#include "HugeInt.h"

#include <cstdlib>   // for abs(), labs(), etc.
#include <iostream>
#include <iomanip>
#include <sstream>
#include <cstring>
#include <stdexcept>
#include <cmath>

不要在标题中包括<iostream>；这太过分了。相反，我们可以包含<iosfwd>，这要轻得多，只需给出我们需要在签名中指定std::istream&和std::ostream&的前向声明，而不需要引入<iostream>中的完整模板定义。这意味着不执行I/O操作的翻译单元不承担开销。

拼写(通篇)：uint32_t和uint64_t位于std命名空间中。您的编译器也可以在全局命名空间中声明它们，但不需要这样做，因此您有一个可移植性错误。

也许我们应该提供一个divmod()函数，在一个运算中给出商和余数，然后用它来实现除法和模？这个“短”版本当然是有用的，并且在执行toDecimalString()时可以节省重复。

我们使用字符串流来实现toDecimalString()，然后使用字符串流实现operator<<()。我认为我们应该采取相反的方法:使用operator<<()实现toDecimalString()。然后流输出就不需要创建临时字符串了。如果您希望保留我的救生器的操作程序状态，请考虑使用它。

在toDecimalString()中，我们有一个循环，每次生成一个数字：

while (!tmp.isZero()) {
    decimalDigits[i++] = tmp.shortModulo(10);
    tmp = tmp.shortDivide(10);
}

我们可以减少九倍的除数，因为shortModulo()接受std::uint32_t。如果我们仔细考虑如何打印第一个“块”，我们可以除以10亿，而不是10：

// determine the decimal digits of the absolute value
constexpr int numDecimal{static_cast<int>(numDigits * 1.07032) + 1};
std::array<std::uint32_t, numDecimal> decimal;

int i{0};
while (!tmp.isZero()) {
    decimal[i++] = tmp.shortModulo(1'000'000'000);
    tmp = tmp.shortDivide(1'000'000'000);
}

// output the decimal digits, in threes
oss << std::setw(0) << std::setfill('0');
--i;
// first digits
auto const d0 = decimal[i];
if (d0 > 1'000'000) {
    oss << (d0 / 1'000'000) << ',' << std::setw(3);
}
if (d0 > 1'000) {
    oss << (d0 / 1'000 % 1'000) << ',' << std::setw(3);
}
oss << (d0 % 1'000);
// subsequent digits
while (i--) {
    auto const d = decimal[i];
    oss << ',' << std::setw(3) << (d / 1'000'000)
        << ',' << std::setw(3) << (d / 1'000 % 1'000)
        << ',' << std::setw(3) << (d % 1'000);
}

Answer 2

即使没有不可移植的优化，乘法方法也可以得到改进(例如，使用_mulx_u64)。遗憾的是，它会使乘法的代码变得不那么漂亮，现在看起来很好，按照下面的方法，它看起来就不那么漂亮了。

这里发生的一件不太好的事情是，计算的顺序迫使创建一个又一个的大的部分乘积，每个乘积都有一个明确的移位和一个完整的BigInt加法。对于阶乘基准，这不是一个问题，因为一个操作数总是很小。例如，对于强调“平衡”乘法的基准，可以将Fibonacci测试扩展到对2x2矩阵进行幂运算。

另一种安排是按权重的顺序生成小的部分乘积(乘数的一个分支乘以被乘数的一个分支)，这样一组权重相等的部分乘积(和一个进位)可以立即相加，而不需要大量临时存储，然后得到最终结果的一个分支。搬运行李有点棘手。

为了清晰起见，这里有一个排序图：

(来源:加密硬件和嵌入式系统)