#48 - Self-powers

The series, $1^1+2^2+3^3+\dots+10^{10}=10405071317$.

Find the last ten digits of the series, $1^1+2^2+3^3+\dots+1000^{1000}$.


All we need is a brute force loop. We then take modulus $10^{10}$ to get the last ten digits.

# file: "problem048.py"
nums = list(range(1, 1001))
print(sum(nums ** nums) % 10 ** 10)

Running this short program results in,

9110846700
0.010996122083266172 seconds.

Therefore, the last ten digits is 9110846700.

Bonus

If you have Python, then solving this problem is trivial. However, most other languages have a cap on the size of integers that can be dealt with. Is there a way to quickly calculate this sum without needing to store a large number such as $1000^{1000}$?

First, notice that we took a modulus to get the last ten digits. Modulus is distributive among addition, so that means

\((a + b)\mod c\equiv (a\mod c + b\mod c)\mod c\) and we can reduce the size of the numbers we are adding. The modulus also allows us to quickly calculate the large raw values needed for the sum.

Just like with addition, modulus is also distributive with multiplication:

\((a\times b)\mod c \equiv ((a\mod c)\times(b\mod c))\mod c\) Performing this repeatedly allows us to calculate $a^k\mod c$ where $k$ is a power of two. For example, $3^{16}\mod 15$:

\(\begin{aligned} 3^2\equiv 9\mod 15 \\ 3^4 = (3^2)^2\equiv ((3^2\mod 15)\times(3^2\mod 15))\mod 15\equiv(9\times9)\mod 15\equiv5\mod15 \\ 3^8=(3^4)^2\equiv((3^4\mod 15)\times(3^4\mod 15))\mod 15\equiv (5\times5)\mod 15\equiv10\mod15 \\ 3^{16}\equiv((3^8\mod15)\times(3^8\mod15))\mod15\equiv(10\times10)\mod 15\equiv10\mod 15 \end{aligned}\) How would we do this for all $k$ though? We can use the binary expansion, as that tells which powers of 2 are used to make the number. For example, let’s do $4^{49}\mod 11$. First convert $49$ to binary. This is $110001_2=2^5+2^4+2^0$. So this means we are trying to find $4^{2^5+2^4+2^0}\mod 11 = \left(4^{2^5}\times4^{2^4}\times4^{2^0}\right)\mod 11$. In that case, we can use the same method and find the appropriate powers of two:

\(\begin{aligned} 4\equiv 4\mod 11 \\ 4^2\equiv (4\times4)\mod 11 \equiv16\mod11\equiv 5\mod 11 \\ 4^4\equiv (4^2\times4^2)\mod 11\equiv 25\mod 11\equiv 3\mod 11 \\ 4^8\equiv (4^4\times4^4)\mod 11\equiv9\mod 11\equiv 9\mod 11 \\ 4^{16}\equiv (4^8\times 4^8)\mod 11\equiv 81\mod 11\equiv 4\mod 11 \\ 4^{32}\equiv (4^{16}\times 4^{16})\mod 11\equiv 16\mod 11\equiv 5\mod 11 \end{aligned}\) Now we pick the powers of 2 we need, plug them in, and calculate:

\(\left(4^{2^5}\times4^{2^4}\times4^{2^0}\right)\mod 11\equiv(5\times4\times4)\mod11\equiv80\mod 11\equiv\boxed{3\mod 11}\) And that’s it! Thun number of times we multiply within the parethesis is proportional to the log base 2 of the power, since we are writing in binary.

In terms of code, we first generate all mods of powers of 2, and then we multiply the appropriate results using the binary representation.

# file: "problem048.py"
def largeMod(a, b, c):
    # Calculate mod c for all power of 2 <= b.
    # Index is the power. so the array looks
    # like [a^(2^0) mod c, a^(2^1) mod c, ...]
    largestPower = int(math.log2(b))
    # First element is a^(2^0) = a^1 = a
    modsOfTwo = [a % c]
    for _ in range(largestPower):
        modsOfTwo.append(modsOfTwo[-1] ** 2 % c)
    # Convert the exponent to binary
    # to see which mods to use...
    bbinary = bin(b)[2:][::-1]
    prod = 1
    for i, bit in enumerate(bbinary, 0):
        if bit == '1':
            prod *= modsOfTwo[i]
    return prod % c

Now solving the problem is just a simple loop.

# file: "problem048.py"
limit = 1000
s = 0
for i in range(1, limit + 1):
    s += largeMod(i, i, 10 ** 10)

print(s % 10 ** 10)

Running the code results in,

9110846700
0.008301412421343397 seconds.

We get the same answer, and slightly faster at that. The function to calculate a large power modulus will be extremely useful, so let’s bookmark this for the future.