Return to Answer

Here's the least ugly form of this micro optimisation that I could come up with:

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = exponent & 1 ? x : 1;

   while (exponent >>= 1)
   {
      x *= x;

      if (exponent & 1)
      {
         result *= x;
      }
   }
       
   return result;
}

Here's a simple demo implementation of this algorithm that uses std::vector<> as a faux big integer:

typedef std::vector<uint32_t> fake_bigint;

fake_bigint &nth_power (fake_bigint &result, uint32_t base, unsigned exponent)
{
   fake_bigint current_power(1, base);

   result.resize(1);
   result[0] = exponent & 1 ? base : 1;

   while (exponent >>= 1)
   {
      current_power *= current_power;

      if (exponent & 1)
      {
         result *= current_power;
      }
   } 

   return result;
}

Almost exactly the same as simple version that used big integers... That's why I like C++.

Here's the implementation for operator *=:

void operator *= (fake_bigint &multiplicand, fake_bigint const &multiplier)
{
   fake_bigint result(multiplicand.size() * 2);

   for (unsigned n = unsigned(multiplier.size()), i = 0; i < n; ++i)
   {
      addmul_1(result, i, multiplicand, multiplier[i]);
   }

   while (result.size() > 1 && result.back() == 0)
   {
      result.pop_back();
   }

   multiplicand.swap(result);
}

The real meat is in the addmul_1 helper function that multiplies a big integer by a single word and adds the product to a big integer, at a certain offset from its lower end to reflect the implied power of two for the multiplicand.

void addmul_1 (fake_bigint &result, unsigned offset, fake_bigint const &multiplicand, fake_bigint::value_type multiplier)
{
   result.resize(std::max(result.size(), offset + multiplicand.size()));

   uint64_t accu = 0;

   for (fake_bigint::size_type n = multiplicand.size(), i = 0; i < n; ++i)
   {
      accu += result[offset + i] + uint64_t(multiplicand[i]) * multiplier;
      result[offset + i] = uint32_t(accu);
      accu >>= 32;
   }

   for (fake_bigint::size_type k = offset + multiplicand.size(); accu; ++k)
   {
      if (k < result.size())
         result[k] = uint32_t(accu += result[k]);
      else
         result.push_back(uint32_t(accu));
      accu >>= 32;
   }
}

On my laptop this computes 9999^9999 in 6.4 ms; I put the number in my pastebin (39996 digits).

In case Martin R's excellent solution should not be fast enough: there is a trick that can drastically reduce the number of multiplications needed for raising a number to the nth power. Instead of doing (n - 1) multiplications with x, you do only a small number of repeated squarings and multiplications.

The basic idea is easily demonstrated:

$$x^4 = x * x * x * x = (x^2)^2$$

Or, as a more drastic example:

$$x^{16} = x * x * x * x * x * x * x * x * x * x * x * x * x * x * x * x = (((x^2)^2)^2)^2$$

Exponents that are not a power of two can also be handled efficiently:

$$x^5 = (x^2)^2 * x$$

How to decide which squaring of x to multiply into the result and which not? Easy: if the binary representation of the exponent has bit i set, then the ith squaring of x needs to be multiplied into the result. Although this is not a very mathematical way of putting things, it happens to be very practical.

The principles can be explored by using a builtin type as a stand-in for a big integer:

typedef uint64_t bigint_t;

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = 1;
   
   for ( ; exponent; exponent >>= 1)
   {
      if (exponent & 1)
      {
         result *= x;
      }
      
      x *= x;
   }
   
   return result;
}

This requires k squarings for an exponent where the most significant bit is k, and as many multiplications into the result as there are '1' bits in the exponent:$$x^{1025}$$ requires only ten squarings and two separate multiplications, instead of 1024 for simple repeated multiplication.

A few more multiplications can he shaved off by handling beginning and end more efficiently, but that is an uglifying optimisation that obscures the basic idea. A tiny bit more speed, but the code gets more complicated and less readable. Don't uglify your code unless you have incontrovertible proof that you have to.

Here's the least ugly form of this micro optimisation that I could come up with:

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = exponent & 1 ? x : 1;

   for ( ; exponent >>= 1; )
   {
      x *= x;

      if (exponent & 1)
      {
         result *= x;
      }
   }
       
   return result;
}

In case Martin R's excellent solution should not be fast enough: there is a trick that can drastically reduce the number of multiplications needed for raising a number to the nth power. Instead of doing (n - 1) multiplications with x, you do only a small number of repeated squarings and multiplications.

The basic idea is easily demonstrated:

$$x^4 = x * x * x * x = (x^2)^2$$

Or, as a more drastic example:

$$x^{16} = x * x * x * x * x * x * x * x * x * x * x * x * x * x * x * x = (((x^2)^2)^2)^2$$

Exponents that are not a power of two can also be handled efficiently:

$$x^5 = (x^2)^2 * x$$

How to decide which squaring of x to multiply into the result and which not? Easy: if the binary representation of the exponent has bit i set, then the ith squaring of x needs to be multiplied into the result. Although this is not a very mathematical way of putting things, it happens to be very practical.

The principles can be explored by using a builtin type as a stand-in for a big integer:

typedef uint64_t bigint_t;

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = 1;
   
   for ( ; exponent; exponent >>= 1)
   {
      if (exponent & 1)
      {
         result *= x;
      }
      
      x *= x;
   }
   
   return result;
}

This requires k squarings for an exponent where the most significant bit is k, and as many multiplications into the result as there are '1' bits in the exponent:$$x^{1025}$$ requires only ten squarings and two separate multiplications, instead of 1024 for simple repeated multiplication.

A few more multiplications can he shaved off by handling beginning and end more efficiently, but that is an uglifying optimisation that obscures the basic idea. A tiny bit more speed, but the code gets more complicated and less readable. Don't uglify your code unless you have incontrovertible proof that you have to.

Here's the least ugly form of this micro optimisation that I could come up with:

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = exponent & 1 ? x : 1;

   while (exponent >>= 1)
   {
      x *= x;

      if (exponent & 1)
      {
         result *= x;
      }
   }
       
   return result;
}

In case Martin R's excellent solution should not be fast enough: there is a trick that can drastically reduce the number of multiplications needed for raising a number to the nth power. Instead of doing (n - 1) multiplications with x, you do only a small number of repeated squarings and multiplications.

The basic idea is easily demonstrated:

$$x^4 = x * x * x * x = (x^2)^2$$

Or, as a more drastic example:

$$x^{16} = x * x * x * x * x * x * x * x * x * x * x * x * x * x * x * x = (((x^2)^2)^2)^2$$

Exponents that are not a power of two can also be handled efficiently:

$$x^5 = (x^2)^2 * x$$

How to decide which squaring of x to multiply into the result and which not? Easy: if the binary representation of the exponent has bit i set, then the ith squaring of x needs to be multiplied into the result. Although this is not a very mathematical way of putting things, it happens to be very practical.

The principles can be explored by using a builtin type as a stand-in for a big integer:

typedef uint64_t bigint_t;

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = 1;
   
   for ( ; exponent; exponent >>= 1)
   {
      if (exponent & 1)
      {
         result *= x;
      }
      
      x *= x;
   }
   
   return result;
}

This requires k squarings for an exponent where the most significant bit is k, and as many multiplications into the result as there are '1' bits in the exponent:$$x^{1025}$$ requires only ten squarings and two separate multiplications, instead of 1024 for simple repeated multiplication.

A few more multiplications can he shaved off by handling beginning and end more efficiently, but that is an uglifying optimisation that obscures the basic idea. A tiny bit more speed, but the code gets more complicated and less readable. Don't uglify your code unless you have incontrovertible proof that you have to.

In case Martin R's excellent solution should not be fast enough: there is a trick that can drastically reduce the number of multiplications needed for raising a number to the nth power. Instead of doing (n - 1) multiplications with x, you do only a small number of repeated squarings and multiplications.

The basic idea is easily demonstrated:

$$x^4 = x * x * x * x = (x^2)^2$$

Or, as a more drastic example:

$$x^{16} = x * x * x * x * x * x * x * x * x * x * x * x * x * x * x * x = (((x^2)^2)^2)^2$$

Exponents that are not a power of two can also be handled efficiently:

$$x^5 = (x^2)^2 * x$$

How to decide which squaring of x to multiply into the result and which not? Easy: if the binary representation of the exponent has bit i set, then the ith squaring of x needs to be multiplied into the result. Although this is not a very mathematical way of putting things, it happens to be very practical.

The principles can be explored by using a builtin type as a stand-in for a big integer:

typedef uint64_t bigint_t;

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = 1;
   
   for ( ; exponent; exponent >>= 1)
   {
      if (exponent & 1)
      {
         result *= x;
      }
      
      x *= x;
   }
   
   return result;
}

This requires k squarings for an exponent where the most significant bit is k, and as many multiplications into the result as there are '1' bits in the exponent:$$x^{1025}$$ requires only ten squarings and two separate multiplications, instead of 1024 for simple repeated multiplication.

A few more multiplications can he shaved off by handling beginning and end more efficiently, but that is an uglifying optimisation that obscures the basic idea. A tiny bit more speed, but the code gets more complicated and less readable. Don't uglify your code unless you have incontrovertible proof that you have to.

Here's the least ugly form of this micro optimisation that I could come up with:

bigint_t nth_power (bigint_t x, unsigned exponent)
{
   bigint_t result = exponent & 1 ? x : 1;

   for ( ; exponent >>= 1; )
   {
      x *= x;

      if (exponent & 1)
      {
         result *= x;
      }
   }
       
   return result;
}