Project Euler Problem 29

Consider all integer combinations of $a^b$ for $2 \le a \le 5$ and $2 \le b \le 5$:

$$\begin{array}{rrrr} 2^2=4, &2^3=8, &2^4=16, &2^5=32\\ 3^2=9, &3^3=27, &3^4=81, &3^5=243\\ 4^2=16, &4^3=64, &4^4=256, &4^5=1024\\ 5^2=25, &5^3=125, &5^4=625, &5^5=3125 \end{array}$$

If they are then placed in numerical order, with any repeats removed, we get the following sequence of $15$ distinct terms:

$$4, 8, 9, 16, 25, 27, 32, 64, 81, 125, 243, 256, 625, 1024, 3125.$$

How many distinct terms are in the sequence generated by $a^b$ for $2 \le a \le 100$ and $2 \le b \le 100$?

We'll solve the more general problem of counting distinct values of $a^b$ for $2 \le a, b \le N$.

Computing $a^b$ directly requires arbitrary precision integers, and this is definitely fine for the original problem. But going to numbers like $1000^{1000}$ and larger we'll want something more efficient.

I tried using logarithms (since $a^b = c^d \Leftrightarrow b \log a = d \log c$) but floating-point precision is an issue and I could not find a good similarity threshold to classify two floats as corresponding to the same value of $a^b$, although a threshold of $10^{-10}$ worked for the original problem.

So instead of computing powers, let's track exponents and reduce this down to a counting problem. The idea is that $4^3 = 2^6$ and $8^2 = 2^6$ are duplicates because they're both $2^6$. More generally, if $a = r^k$ for some primitive root $r$ (a number that isn't itself a perfect power), then $a^b = r^{kb}$.

Looking at powers of $2$ to start and assuming $N = 100$ the bases ${2, 4, 8, 16, 32, 64}$ all reduce to powers of 2. So the number of unique powers with root 2 is the cardinality of the set

where $m$ is the number of powers to include. $m = 6$ for a base root of 2 since we have 6 powers of 2 below $N = 100$.

We don't have to look at the base anymore since all the numbers have a base power of $2$, so we just look at the values of the exponents.

Writing it a bit more generally, we can say

So for $N = 100$ since we have 6 numbers that all share the primitive root 2, we have $S_6$ unique powers of 2. Similarly, we will have $S_4$ unique powers of 3 (from bases 3, 9, 27, and 81), $S_2$ unique powers of 5 (from bases 5 and 25), $S_2$ unique powers of 6, and so on. Powers of 4 are counted in powers of 2, so are powers of 8, powers of 9 are counted in powers of 3, and so on. The maximum value of $m$ we need to go up to is for the powers of 2, for which we only need to go up to $\lfloor \log_2 N \rfloor$.

function count_distinct_powers(N)
    max_log = floor(Int, log2(N))

    # Precompute S_m for m = 1 to log₂(N)
    unique_exponent_counts = Vector{Int}(undef, max_log)
    for m in 1:max_log
        seen = Set{Int}()
        for j in 1:m
            for b in 2:N
                push!(seen, j * b)
            end
        end
        unique_exponent_counts[m] = length(seen)
    end

    # Find primitive roots and sum their contributions
    is_perfect_power = falses(N)
    result = 0

    for base in 2:N
        if is_perfect_power[base]
            continue
        end

        power_count = 1
        val = base
        while true
            next_val = val * base
            if next_val > N
                break
            end
            val = next_val
            power_count += 1
            is_perfect_power[val] = true
        end

        result += unique_exponent_counts[power_count]
    end

    return result
end

This runs pretty fast even up to $N = 10^5$ as tabulated below.

The precomputation loops over $m$ from 1 to $\log_2 N$, and for each $m$ iterates over $j \in [1, m]$ and $b \in [2, N]$, giving $N (1 + 2 + \cdots + \log_2 N) = \mathcal{O}(N \log^2 N)$ set insertions. Since the main loop is $\mathcal{O}(N)$, the whole thing runs in $\mathcal{O}(N \log^2 N)$.