abacus: bignum: add exponentiation

src/bignum/bignum.cc

@@ -129,6 +129,49 @@ std::pair<digits_type, digits_type> do_div_mod(digits_type const& lhs,
     return std::make_pair(quotient, remainder);
 }
 
+// More optimised than full-on div_mod
+digits_type do_halve(digits_type num) {
+    assert(num.size() != 0);
+
+    int carry = 0;
+    for (auto i = num.rbegin(); i != num.rend(); ++i) {
+        auto const was_odd = (*i % 2) == 1;
+        *i /= 2;
+        *i += carry;
+        if (was_odd) {
+            carry = BASE / 2;
+        }
+    }
+
+    do_trim_leading_zeros(num);
+
+    return num;
+}
+
+bool is_odd(digits_type const& num) {
+    if (num.size() == 0) {
+        return false;
+    }
+
+    return (num.front() % 2) == 1;
+}
+
+digits_type do_pow(digits_type lhs, digits_type rhs) {
+    assert(rhs.size() != 0);
+
+    auto original = lhs;
+
+    while (rhs.size() != 0 && !(rhs.size() == 1 && rhs.front() == 1)) {
+        lhs = do_multiplication(lhs, lhs);
+        if (is_odd(rhs)) {
+            lhs = do_multiplication(lhs, original);
+        }
+        rhs = do_halve(rhs);
+    }
+
+    return lhs;
+}
+
 } // namespace
 
 BigNum::BigNum(std::int64_t number) {
@@ -362,4 +405,23 @@ std::pair<BigNum, BigNum> div_mod(BigNum const& lhs, BigNum const& rhs) {
     return std::make_pair(quotient, remainder);
 }
 
+BigNum pow(BigNum const& lhs, BigNum const& rhs) {
+    assert(lhs.is_canonicalized());
+    assert(rhs.is_canonicalized());
+
+    if (rhs.is_zero()) {
+        return BigNum(1);
+    } else if (rhs.is_negative()) {
+        return BigNum();
+    }
+
+    BigNum res;
+    res.digits_ = do_pow(lhs.digits_, rhs.digits_);
+
+    res.sign_ = is_odd(rhs.digits_) ? lhs.sign_ : 1;
+    res.canonicalize();
+
+    return res;
+}
+
 } // namespace abacus::bignum

src/bignum/bignum.hh

@@ -87,6 +87,8 @@ public:
     friend std::pair<BigNum, BigNum> div_mod(BigNum const& lhs,
                                              BigNum const& rhs);
 
+    friend BigNum pow(BigNum const& lhs, BigNum const& rhs);
+
     friend bool operator==(BigNum const& lhs, BigNum const& rhs) {
         return lhs.equal(rhs);
     }

tests/unit/bignum.cc

@@ -267,3 +267,66 @@ TEST(BigNum, div_mod_identity) {
                   + (minus_five % minus_three),
               minus_five);
 }
+
+TEST(BigNum, pow_negative_exponent) {
+    auto const zero = BigNum(0);
+    auto const minus_one = BigNum(-1);
+    auto const three = BigNum(3);
+
+    EXPECT_EQ(pow(three, minus_one), zero);
+}
+
+TEST(BigNum, pow_zero) {
+    auto const zero = BigNum(0);
+    auto const one = BigNum(1);
+    auto const three = BigNum(3);
+
+    EXPECT_EQ(pow(three, zero), one);
+    EXPECT_EQ(pow(zero, three), zero);
+    EXPECT_EQ(pow(zero, zero), one); // unclear mathematically
+}
+
+TEST(BigNum, pow_one) {
+    auto const one = BigNum(1);
+    auto const three = BigNum(3);
+
+    EXPECT_EQ(pow(one, one), one);
+    EXPECT_EQ(pow(one, three), one);
+    EXPECT_EQ(pow(three, one), three);
+}
+
+TEST(BigNum, pow) {
+    auto const one = BigNum(1);
+    auto const two = BigNum(2);
+    auto const three = BigNum(3);
+    auto const four = BigNum(4);
+    auto const eight = BigNum(8);
+    auto const nine = BigNum(9);
+    auto const twenty_seven = BigNum(27);
+    auto const eighty_one = BigNum(81);
+
+    EXPECT_EQ(pow(two, two), four);
+    EXPECT_EQ(pow(two, three), eight);
+    EXPECT_EQ(pow(three, two), nine);
+    EXPECT_EQ(pow(three, three), twenty_seven);
+    EXPECT_EQ(pow(three, four), eighty_one);
+}
+
+TEST(BigNum, pow_negative) {
+    auto const one = BigNum(1);
+    auto const two = BigNum(2);
+    auto const minus_two = BigNum(-2);
+    auto const three = BigNum(3);
+    auto const minus_three = BigNum(-3);
+    auto const four = BigNum(4);
+    auto const minus_eight = BigNum(-8);
+    auto const nine = BigNum(9);
+    auto const minus_twenty_seven = BigNum(-27);
+    auto const eighty_one = BigNum(81);
+
+    EXPECT_EQ(pow(minus_two, two), four);
+    EXPECT_EQ(pow(minus_two, three), minus_eight);
+    EXPECT_EQ(pow(minus_three, two), nine);
+    EXPECT_EQ(pow(minus_three, three), minus_twenty_seven);
+    EXPECT_EQ(pow(three, four), eighty_one);
+}