|
|
@@ -12,9 +12,9 @@
|
|
|
using namespace math::detail;
|
|
|
|
|
|
namespace {
|
|
|
- constexpr int32_t const MAX_SEG{999999999};
|
|
|
- constexpr int32_t const OVER_SEG{1000000000};
|
|
|
- constexpr int32_t const SEG_DIGITS{9};
|
|
|
+constexpr int32_t const MAX_SEG{999999999};
|
|
|
+constexpr int32_t const OVER_SEG{1000000000};
|
|
|
+constexpr int32_t const SEG_DIGITS{9};
|
|
|
}
|
|
|
|
|
|
namespace math { namespace detail {
|
|
|
@@ -23,130 +23,130 @@ namespace math { namespace detail {
|
|
|
return -1; \
|
|
|
else if (lhs lexpr < rhs rexpr) \
|
|
|
return 1
|
|
|
- int compare(data_type const & rhs, data_type const & lhs, size_t offset) {
|
|
|
- IMPL_COMPARE(.size(), .size() + offset);
|
|
|
- for (size_t i = lhs.size(); i > 0; --i) {
|
|
|
- IMPL_COMPARE([i + offset - 1], [i - 1]);
|
|
|
- }
|
|
|
- return 0;
|
|
|
+int compare(data_type const & rhs, data_type const & lhs, size_t offset) {
|
|
|
+ IMPL_COMPARE(.size(), .size() + offset);
|
|
|
+ for (size_t i = lhs.size(); i > 0; --i) {
|
|
|
+ IMPL_COMPARE([i + offset - 1], [i - 1]);
|
|
|
}
|
|
|
+ return 0;
|
|
|
+}
|
|
|
#undef IMPL_COMPARE
|
|
|
|
|
|
- void carry(data_type & rhs, size_t start) {
|
|
|
- auto const ubnd = rhs.size() - 1;
|
|
|
- for (size_t i = start; i < ubnd; ++i) {
|
|
|
- if (rhs[i] > MAX_SEG) {
|
|
|
- rhs[i] -= OVER_SEG;
|
|
|
- rhs[i + 1] += 1;
|
|
|
- }
|
|
|
+void carry(data_type & rhs, size_t start) {
|
|
|
+ auto const ubnd = rhs.size() - 1;
|
|
|
+ for (size_t i = start; i < ubnd; ++i) {
|
|
|
+ if (rhs[i] > MAX_SEG) {
|
|
|
+ rhs[i] -= OVER_SEG;
|
|
|
+ rhs[i + 1] += 1;
|
|
|
}
|
|
|
- if (rhs[ubnd] == 0) { rhs.pop_back(); }
|
|
|
}
|
|
|
+ if (rhs[ubnd] == 0) { rhs.pop_back(); }
|
|
|
+}
|
|
|
|
|
|
- void add(data_type & rhs, data_type const & lhs, size_t offset) {
|
|
|
- rhs.resize(std::max(rhs.size(), lhs.size() + offset) + 1);
|
|
|
- auto const lbnd = lhs.size();
|
|
|
- for (size_t i = 0; i < lbnd; ++i) {
|
|
|
- rhs[i + offset] += lhs[i];
|
|
|
- }
|
|
|
- carry(rhs, offset);
|
|
|
+void add(data_type & rhs, data_type const & lhs, size_t offset) {
|
|
|
+ rhs.resize(std::max(rhs.size(), lhs.size() + offset) + 1);
|
|
|
+ auto const lbnd = lhs.size();
|
|
|
+ for (size_t i = 0; i < lbnd; ++i) {
|
|
|
+ rhs[i + offset] += lhs[i];
|
|
|
}
|
|
|
+ carry(rhs, offset);
|
|
|
+}
|
|
|
|
|
|
- void add(data_type & rhs, int32_t lhs, size_t offset) {
|
|
|
- rhs.resize(std::max(rhs.size(), offset) + 1);
|
|
|
- rhs[offset] += lhs;
|
|
|
- carry(rhs, offset);
|
|
|
- }
|
|
|
+void add(data_type & rhs, int32_t lhs, size_t offset) {
|
|
|
+ rhs.resize(std::max(rhs.size(), offset) + 1);
|
|
|
+ rhs[offset] += lhs;
|
|
|
+ carry(rhs, offset);
|
|
|
+}
|
|
|
|
|
|
- void subtract_nounderflow(data_type & rhs, data_type const & lhs,
|
|
|
- size_t offset) {
|
|
|
- size_t const rbnd = rhs.size(), lbnd = lhs.size();
|
|
|
- // Subtract
|
|
|
- for (size_t i = 0; i < lbnd; ++i) {
|
|
|
- rhs[i + offset] -= lhs[i];
|
|
|
- }
|
|
|
- // Borrow
|
|
|
- for (size_t i = 0; i < rbnd - 1; ++i) {
|
|
|
- if (rhs[i] < 0) {
|
|
|
- rhs[i] += OVER_SEG;
|
|
|
- rhs[i + 1] -= 1;
|
|
|
- }
|
|
|
- }
|
|
|
- if (rhs[rbnd - 1] == 0 && rbnd > 1) { rhs.pop_back(); }
|
|
|
+void subtract_nounderflow(data_type & rhs, data_type const & lhs,
|
|
|
+ size_t offset) {
|
|
|
+ size_t const rbnd = rhs.size(), lbnd = lhs.size();
|
|
|
+ // Subtract
|
|
|
+ for (size_t i = 0; i < lbnd; ++i) {
|
|
|
+ rhs[i + offset] -= lhs[i];
|
|
|
}
|
|
|
-
|
|
|
- void multiply(data_type & rhs, data_type const & lhs) {
|
|
|
- size_t const rbnd = rhs.size(), lbnd = lhs.size();
|
|
|
- size_t const ubnd = rbnd + lbnd;
|
|
|
- rhs.resize(ubnd);
|
|
|
- // Multiply
|
|
|
- for (size_t i = rbnd; i > 0; --i) {
|
|
|
- int32_t const value = rhs[i - 1];
|
|
|
- for (size_t j = lbnd; j > 0; --j) {
|
|
|
- // Max input 999,999,999
|
|
|
- // Max output 999,999,998,000,000,001
|
|
|
- int64_t product = int64_t(value) * lhs[j - 1];
|
|
|
- int64_t overflow = product / OVER_SEG;
|
|
|
- size_t const to = i + j - 2;
|
|
|
- rhs[to] += int32_t(product - (overflow * OVER_SEG));
|
|
|
- rhs[to + 1] += int32_t(overflow);
|
|
|
- }
|
|
|
- rhs[i - 1] -= value;
|
|
|
- }
|
|
|
- // Carry
|
|
|
- for (size_t i = 0; i < ubnd - 1; ++i) {
|
|
|
- if (rhs[i] > MAX_SEG) {
|
|
|
- int32_t overflow = rhs[i] / OVER_SEG;
|
|
|
- rhs[i] -= (overflow * OVER_SEG);
|
|
|
- rhs[i + 1] += overflow;
|
|
|
- }
|
|
|
- }
|
|
|
- while (rhs.back() == 0) {
|
|
|
- rhs.pop_back();
|
|
|
+ // Borrow
|
|
|
+ for (size_t i = 0; i < rbnd - 1; ++i) {
|
|
|
+ if (rhs[i] < 0) {
|
|
|
+ rhs[i] += OVER_SEG;
|
|
|
+ rhs[i + 1] -= 1;
|
|
|
}
|
|
|
}
|
|
|
+ if (rhs[rbnd - 1] == 0 && rbnd > 1) { rhs.pop_back(); }
|
|
|
+}
|
|
|
|
|
|
- data_type shift10(data_type const & data, int32_t places) {
|
|
|
- int32_t shift = places / SEG_DIGITS;
|
|
|
- int64_t const pow = powers[places - (shift * SEG_DIGITS)];
|
|
|
- size_t const bnd = data.size();
|
|
|
- data_type rval(size_t((int32_t)bnd + shift) + 1);
|
|
|
- for (size_t i = 0, o = size_t(std::max(0, shift)); i < bnd; ++i, ++o) {
|
|
|
- int64_t product = int64_t(data[i]) * pow;
|
|
|
+void multiply(data_type & rhs, data_type const & lhs) {
|
|
|
+ size_t const rbnd = rhs.size(), lbnd = lhs.size();
|
|
|
+ size_t const ubnd = rbnd + lbnd;
|
|
|
+ rhs.resize(ubnd);
|
|
|
+ // Multiply
|
|
|
+ for (size_t i = rbnd; i > 0; --i) {
|
|
|
+ int32_t const value = rhs[i - 1];
|
|
|
+ for (size_t j = lbnd; j > 0; --j) {
|
|
|
+ // Max input 999,999,999
|
|
|
+ // Max output 999,999,998,000,000,001
|
|
|
+ int64_t product = int64_t(value) * lhs[j - 1];
|
|
|
int64_t overflow = product / OVER_SEG;
|
|
|
- rval[o] += int32_t(product - (overflow * OVER_SEG));
|
|
|
- rval[o + 1] += int32_t(overflow);
|
|
|
+ size_t const to = i + j - 2;
|
|
|
+ rhs[to] += int32_t(product - (overflow * OVER_SEG));
|
|
|
+ rhs[to + 1] += int32_t(overflow);
|
|
|
+ }
|
|
|
+ rhs[i - 1] -= value;
|
|
|
+ }
|
|
|
+ // Carry
|
|
|
+ for (size_t i = 0; i < ubnd - 1; ++i) {
|
|
|
+ if (rhs[i] > MAX_SEG) {
|
|
|
+ int32_t overflow = rhs[i] / OVER_SEG;
|
|
|
+ rhs[i] -= (overflow * OVER_SEG);
|
|
|
+ rhs[i + 1] += overflow;
|
|
|
}
|
|
|
- if (rval.back() == 0 && rval.size() > 1) { rval.pop_back(); }
|
|
|
- return rval;
|
|
|
}
|
|
|
+ while (rhs.back() == 0) {
|
|
|
+ rhs.pop_back();
|
|
|
+ }
|
|
|
+}
|
|
|
|
|
|
- size_t digits(int32_t val) {
|
|
|
- return val == 0 ? 1 : size_t(floor(log10(val))) + 1;
|
|
|
+data_type shift10(data_type const & data, int32_t places) {
|
|
|
+ int32_t shift = places / SEG_DIGITS;
|
|
|
+ int64_t const pow = powers[places - (shift * SEG_DIGITS)];
|
|
|
+ size_t const bnd = data.size();
|
|
|
+ data_type rval(size_t((int32_t)bnd + shift) + 1);
|
|
|
+ for (size_t i = 0, o = size_t(std::max(0, shift)); i < bnd; ++i, ++o) {
|
|
|
+ int64_t product = int64_t(data[i]) * pow;
|
|
|
+ int64_t overflow = product / OVER_SEG;
|
|
|
+ rval[o] += int32_t(product - (overflow * OVER_SEG));
|
|
|
+ rval[o + 1] += int32_t(overflow);
|
|
|
}
|
|
|
+ if (rval.back() == 0 && rval.size() > 1) { rval.pop_back(); }
|
|
|
+ return rval;
|
|
|
+}
|
|
|
|
|
|
- size_t digits(data_type const & data) {
|
|
|
- size_t segDig = SEG_DIGITS * (data.size() - 1);
|
|
|
- if (data.back() == 0) {
|
|
|
- return segDig - SEG_DIGITS + digits(data[data.size() - 2]);
|
|
|
- } else {
|
|
|
- return segDig + digits(data.back());
|
|
|
- }
|
|
|
+size_t digits(int32_t val) {
|
|
|
+ return val == 0 ? 1 : size_t(floor(log10(val))) + 1;
|
|
|
+}
|
|
|
+
|
|
|
+size_t digits(data_type const & data) {
|
|
|
+ size_t segDig = SEG_DIGITS * (data.size() - 1);
|
|
|
+ if (data.back() == 0) {
|
|
|
+ return segDig - SEG_DIGITS + digits(data[data.size() - 2]);
|
|
|
+ } else {
|
|
|
+ return segDig + digits(data.back());
|
|
|
}
|
|
|
+}
|
|
|
|
|
|
- data_type divide(data_type & remainder, data_type const & divisor) {
|
|
|
- data_type accum{0};
|
|
|
- auto const dig = digits(divisor);
|
|
|
+data_type divide(data_type & remainder, data_type const & divisor) {
|
|
|
+ data_type accum{0};
|
|
|
+ auto const dig = digits(divisor);
|
|
|
+ do {
|
|
|
+ auto const diff = std::max(1UL, digits(remainder) - dig) - 1;
|
|
|
+ auto const shift = diff / SEG_DIGITS;
|
|
|
+ auto const ipow = diff - (shift * SEG_DIGITS);
|
|
|
+ data_type value{shift10(divisor, int32_t(ipow))};
|
|
|
do {
|
|
|
- auto const diff = std::max(1UL, digits(remainder) - dig) - 1;
|
|
|
- auto const shift = diff / SEG_DIGITS;
|
|
|
- auto const ipow = diff - (shift * SEG_DIGITS);
|
|
|
- data_type value{shift10(divisor, int32_t(ipow))};
|
|
|
- do {
|
|
|
- subtract_nounderflow(remainder, value, shift);
|
|
|
- add(accum, powers[ipow], shift);
|
|
|
- } while (compare(remainder, value, shift) >= 0); // Up to 10 times
|
|
|
- } while (compare(remainder, divisor) >= 0);
|
|
|
- return accum;
|
|
|
- }
|
|
|
+ subtract_nounderflow(remainder, value, shift);
|
|
|
+ add(accum, powers[ipow], shift);
|
|
|
+ } while (compare(remainder, value, shift) >= 0); // Up to 10 times
|
|
|
+ } while (compare(remainder, divisor) >= 0);
|
|
|
+ return accum;
|
|
|
+}
|
|
|
}}
|