Bignum multiplication of numbers with arbitrary bases, using Karatsuba's algorithm.

bn_karatsuba.cpp

/**
 *  Bignum multiplication of numbers with arbitrary bases, using Karatsuba's algorithm.
 *
 *  g++ -std=c++11 -Wall -Wextra -Werror -pedantic-errors bn_karatsuba.cpp -o bn_karatsuba
 *  
 *  usage: bn_karatsuba.exe 123456789 987654321
 */

#include <algorithm>
#include <cassert>
#include <cstring>
#include <iomanip>
#include <iostream>
#include <iterator>
#include <string>
#include <vector>

#define NUM_BASE_16

#define PRINT_DEBUG std::cout << "DEBUG: " << __FUNCTION__ << " - " << __LINE__ << std::endl

#ifdef NUM_BASE_16
    enum {
        NUM_MIN = 0,
        NUM_MAX = 15,
        NUM_BASE = 16
    };
#elif defined(NUM_BASE_8)
    enum {
        NUM_MIN = 0,
        NUM_MAX = 7,
        NUM_BASE = 8
    };
#elif defined(NUM_BASE_2)
    enum {
        NUM_MIN = 0,
        NUM_MAX = 1,
        NUM_BASE = 2
    };
#else
    enum {
        NUM_MIN = 0,
        NUM_MAX = 9,
        NUM_BASE = 10
    };
#endif

typedef unsigned int bn_t;
typedef unsigned long bn_double_t;
typedef std::vector<bn_t> bignum;

/**
 *  Simple bignum printing operation
 */
std::ostream& operator << (std::ostream& ostr, const bignum& bn) {
    for (unsigned i = bn.size(); i > 0; --i) {
        const unsigned index = i-1;
        bn_t digit = bn[index];

        #ifdef NUM_BASE_16
            ostr << std::hex << digit;
        #else
            ostr << digit;
        #endif
    }

    return ostr;
}

/*
 * 
 */
bool bignum_from_str(bignum& b, const char* const number) {
    int i = strlen(number);
    while (i --> 0) {
        bn_t digit = number[i];

        if (digit >= '0' && digit <= '9') {
            digit -= '0';
        }
#ifdef NUM_BASE_16
        else if (digit >= 'a' && digit <= 'f') {
            digit = 10 + (digit-'a');
        }
        else if (digit >= 'A' && digit <= 'F') {
            digit =  10 + (digit-'A');
        }
#endif
        else {
            std::cerr
                << "Invalid digit \'" << digit
                << "\' in parameter \'" << number << "\'."
                << std::endl;
                b.clear();
                return false;
        }

        b.push_back(digit);
    }

    return true;
}

/*
 * 
 */
bignum bignum_from_num(bn_double_t n) {
    bn_double_t count = n;
    bn_double_t numDigits = 0;
    
    while (count) {
        ++numDigits;
        count /= NUM_BASE;
    }
    
    bignum ret;
    ret.reserve(numDigits);
    
    count = n;
    
    while (count) {
        ret.push_back((bn_t)count % NUM_BASE);
        count /= NUM_BASE;
    }
    
    return ret;
}

/**
 * Convert the first two parameters into numbers
 */
bool tokenize(int argc, char** argv, bignum& out1, bignum& out2){
    bignum first = {};
    bignum second = {};

    if (argc != 3) {
        std::cerr << "Invalid number of arguments." << std::endl;
        return false;
    };

    if (!bignum_from_str(first, argv[1])
    || !bignum_from_str(second, argv[2])
    ) {
        return false;
    }

    out1 = std::move(first);
    out2 = std::move(second);

    return true;
}

/**
 *  Determine what number is larger than the other
 */
bool operator < (const bignum& a, const bignum& b) {
    if (a.size() < b.size()) {
        return true;
    }

    if (a.size() > b.size()) {
        return false;
    }

    // check each digit for its value
    for (unsigned i = 0; i < a.size(); ++i) {
        if (a[i] < b[i]) {
            return true;
        }
    }

    return false;
}

/**
 *  BigNum Addition
 */
bignum operator +(const bignum& x, const bignum& y) {
    // numerical iterators
    typename bignum::size_type outIter = 0;
    typename bignum::size_type inIter = 0;
    const bignum& inNum = y;
    bignum outNum = x;

    // arithmetic remainder
    bn_t remainder = 0;

    while (true) {
        // add a number to this object's container if necessary
        if (outIter == outNum.size()) {
            outNum.push_back(0);
        }

        // get the single digits that are to be added
        const bn_double_t inDigit = (inIter < inNum.size()) ? inNum[inIter] : bn_double_t{0};
        const bn_double_t outDigit = outNum[outIter] + inDigit + remainder;

        // allow integer overflow to handle the resulting value.
        if (NUM_MAX < outDigit) {
            outNum[outIter] = outDigit-bn_double_t{NUM_MAX}-bn_double_t{1};
            remainder = 1;
        }
        else {
            outNum[outIter] = outDigit;
            remainder = 0;
        }

        // onto the next number
        ++outIter;
        ++inIter;

        // determine when to stop counting
        if (outIter == outNum.size() && inIter >= inNum.size() && 0 == remainder) {
            break;
        }
    }

    return outNum;
}

/**
 *  BigNum Subtraction
 *  subtract y from x
 */
bignum operator -(const bignum& x, const bignum& y) {
    const unsigned maxSize = std::max(x.size(), y.size());
    bn_t remainder = 0;
    bn_t nextX = 0;
    bignum sum;

    sum.reserve(maxSize);

    for (unsigned i = 0; i < maxSize; ++i) {

        bn_t a = (i < x.size()) ? x[i]-nextX : 0; // x s ith digit
        bn_t b = (b < y.size()) ? y[i] : 0; // y's ith digit

        if (b > a) {
            a += NUM_MAX;

            if (i < x.size()) {
                nextX = 1;
            }
            else {
                nextX = 0;
            }
        }

        remainder = a-b; // actual sum, use overflow to get the ith result
        sum.push_back(remainder);
        remainder = (remainder > 0) ? NUM_MAX : 0; // remainder, if one exists
    }

    if (remainder != 0) {
        sum.push_back(remainder);
    }

    return sum;
}

/**
 *  Multiply bignums
 */
bignum basic_mul(const bignum& a, const bignum& b) {
    bignum ret{};
    
    const bignum::size_type totalLen = a.size() + b.size();
    
    if (!totalLen) {
        return ret;
    }
    
    ret.resize(totalLen);
    
    for (unsigned i = 0; i < b.size(); ++i) {
        unsigned remainder = 0;
        
        for (unsigned j = 0; j < a.size(); ++j) {
            const bignum::size_type index = i + j;
            ret[index] += remainder + a[j] * b[i];
            remainder = ret[index] / NUM_BASE;
            ret[index] = ret[index] % NUM_BASE;
        }
        
        ret[i+a.size()] += remainder;
    }
    
    return ret;
}

/*
 * 
 */
constexpr bn_double_t const_pow(const bn_double_t& a, const bn_double_t p = 2) {
    return p ? a * const_pow(a, p-1) : 1;
}

/*
 * 
 */
bn_double_t bloat_max(bignum& a, bignum& b) {
    const bn_double_t aLen = a.size();
    const bn_double_t bLen = b.size();

    if (aLen > bLen) {
        b.resize(aLen-bLen-1, 0);
        return aLen;
    }
    else if (bLen > aLen) {
        a.resize(bLen-aLen-1, 0);
        return bLen;
    }
    else {
        return aLen;
    }

    return aLen >= bLen ? aLen : bLen;
}

/*
 * 
 */
bignum split_upper(const bignum& n, bn_double_t splitpoint) {
    return bignum{n.cbegin(), n.cbegin()+splitpoint};
}

/*
 * 
 */
bignum split_lower(const bignum& n, bn_double_t splitpoint) {
    return bignum{n.cend()-splitpoint, n.cend()};
}

/*
 * 
 */
bignum operator * (const bignum& a, const bignum& b) {
    const bignum::size_type aLen = a.size();
    const bignum::size_type bLen = b.size();
    
    const bignum::size_type half = ((aLen >= bLen ? aLen : bLen)+1) / 2;
    
    // seems like a good stopping point.
    if (half <= NUM_BASE*NUM_BASE) {
        return basic_mul(a, b);
    }
    
    const bignum&& aHi = split_upper(a, half);
    const bignum&& aLo = split_lower(a, half);
    const bignum&& bHi = split_upper(b, half);
    const bignum&& bLo = split_lower(b, half);
    
    const bignum&& z0 = aLo * bLo;
    const bignum&& z1 = (aLo+aHi) * (bLo+bHi);
    const bignum&& z2 = aHi * bHi;
    const bignum&& z3 = z1-z2-z0;
    
    const bignum&& p2 = bignum_from_num(NUM_BASE*const_pow(half, 2));
    const bignum&& pm = bignum_from_num(const_pow(NUM_BASE, half));

    return (z2 * p2) + (z3 * pm) + z0;
}

/**
 *  Main()
 */
int main(int, char**) {
    #ifdef NUM_BASE_16
        bignum a = {1,2,3,4,5,6,7,8,9,10,11,12,13,14,15};
        bignum b = {15,14,13,12,11,10,9,8,7,6,5,4,3,2,1};
    #elif defined(NUM_BASE_8)
        bignum a = {0,1,2,3,4,5,6,7,0,1,2,3,4,5,6,7};
        bignum b = {0,7,6,5,4,3,2,1,0,7,6,5,4,3,2,1};
    #elif defined(NUM_BASE_2)
        bignum a = {0,1,1,1,0,1,1,1,0,1,1,1,0,1,1,1,0,1,1,1};
        bignum b = {0,0,0,1,0,0,0,1,0,0,0,1,0,0,0,1,0,0,0,1};
    #else
        bignum a = {1,2,3,4,5,6,7,8,9,0,1,2,3,4,5,6,7,8,9};
        bignum b = {9,8,7,6,5,4,3,2,1,0,9,8,7,6,5,4,3,2,1};
    #endif

    /*
    if (!tokenize(argc, argv, a, b)) {
        return -1;
    }
    */

    std::cout << "Testing bignum multiplication." << std::endl;
    std::cout << a << " * " << b << " = " << a*b << std::endl;

    return 0;
}