Rhapsody 🏳️‍⚧️ Θ∆

Exponentiation by squaring is a useful trick that comes up pretty often if you need to multiply something by itself a large number of times. I mentioned this trick in passing in my blog post on Blazing fast Fibonacci numbers using monoids but I thought I'd go into more depth about how it actually works.

Let's imagine that we want to compute n⁶⁴ in as few multiplications as possible. The inefficient way would be to do something like this:

n × n × n × … × n
─────────────────
    64 times

… because we have to use the multiplication operation 63 times. We can substantially reduce the number of multiplications by instead doing something like this:

n²  = n   × n
n⁴  = n²  × n²
n⁸  = n⁴  × n⁴
n¹⁶ = n⁸  × n⁸
n³² = n¹⁶ × n¹⁶
n⁶⁴ = n³² × n³²

result = n⁶⁴

That only takes 6 multiplications, which is a lot fewer than 63.

However, this is the best case scenario where the exponent is a power of two. What would we do if we wanted to compute n⁶⁵? Well, we only need one additional multiplication at the very end:

n²  = n   × n
n⁴  = n²  × n²
n⁸  = n⁴  × n⁴
n¹⁶ = n⁸  × n⁸
n³² = n¹⁶ × n¹⁶
n⁶⁴ = n³² × n³²

result = n⁶⁴ × n

Still only 7 multiplications, which is way better than 64 multiplications.

Okay, maybe that was still a bit too easy. What if we need to compute n⁹⁹? Then we do this:

n²  = n   × n
n⁴  = n²  × n²
n⁸  = n⁴  × n⁴
n¹⁶ = n⁸  × n⁸
n³² = n¹⁶ × n¹⁶
n⁶⁴ = n³² × n³²

result = n⁶⁴ × n³² × n² × n

Only 9 multiplications, which is still way better than 98 multiplications.

In fact, for any whole exponent up to n¹²⁷ you can compute it as the product of some subset of n, n², n⁴, n⁸, n¹⁶, n³², n⁶⁴ (which I'll refer to as the "squares of n"). Or in other words, you can compute the result as a product at most 7 of those squares of n, using each one no more than once. The reason this works is because every natural number from 0 to 127 can be represented as a binary number with up to 7 digits (bits), which we'll denote using 0bNNNNNNN:

  0: 0b0000000
  1: 0b0000001
  2: 0b0000010
  3: 0b0000011
  4: 0b0000100
  …
 64: 0b1000000
 65: 0b1000001
  …
 99: 0b1100011
  …
127: 0b1111111

… and we can use the binary representation of each number to figure out which of n, n², n⁴, n⁸, n¹⁶, n³², n⁶⁴ to multiply into our final result.

To illustrate how, I'll start with the concrete example of 99, which has a binary representation of:

0b1100011

… and we can pick out which squares of n to combine to get n⁹⁹ by doing this:

0b1100011
  ↓↓↓↓↓↓↓
  1        --       include n⁶⁴ (because this digit is 1)
   1       --       include n³² (because this digit is 1)
    0      -- don't include n¹⁶ (because this digit is 0)
     0     -- don't include n⁸  (because this digit is 0)
      0    -- don't include n⁴  (because this digit is 0)
       1   --       include n²  (because this digit is 1)
        1  --       include n   (because this digit is 1)

result = n⁶⁴ × n³² × n² × n

… and that works because n⁹⁹ = n^0b1100011 = n^{(64 + 32 + 2 + 1)} = n⁶⁴ × n³² × n² × n¹.

This means that we only need up to 12 multiplications to compute all whole powers of n up to n¹²⁷:

• a most 6 multiplications to compute the squares of n we need (i.e. n², n⁴, n⁸, n¹⁶, n³², n⁶⁴)
• at most 6 multiplications to combine those squares of n into the final result

… with the worst case being n¹²⁷ (which requires all 12 multiplications).

In fact, we can generalize this result to handle arbitrary whole powers of n. To compute n^m we need:

• at most ⌊log₂(m)⌋ multiplications to compute the squares of n we need
• at most ⌊log₂(m)⌋ multiplications to combine those squares of n into the final result

… and that's because the binary representation of the exponent (m) requires (⌊log₂(m)⌋ + 1) bits.

Computer scientists usually don't need the time complexity to be that precise, so we can simplify that by just saying that we need at most O(log(n)) multiplications, which is a really good time complexity. For example, this trick lets me compute something like n^1000000000 in only 41 multiplications.

That's essentially the entirety of how exponentiation by squaring works. To compute n^m you:

• compute as many squares of n as you need
• use the binary representation of m to select which squares of n to combine into the final result (based on which digits in the binary representation are 1)

If you want to see a particularly neat application of this trick, my blog post on fibonacci numbers that I mentioned above has a worked example that builds on this trick.