Detailed mathematical approach

Problem Summary

Define

$P(x)=1+x+x^2+x^3+x^4.$

For each nonnegative integer \(k\), let \(Q(k)\) be the number of odd coefficients in the expansion of \(P(x)^k\). The problem asks for

$\sum_{t=1}^{18} Q(10^t).$

A direct expansion is hopeless when \(k\) is as large as \(10^{18}\). The implementations therefore never build the polynomial itself. They only track coefficient parity, because the question is about which coefficients are odd.

Mathematical Approach

Write the binary expansion of \(k\) as

$k=\sum_{b\ge 0} k_b 2^b,\qquad k_b\in\{0,1\}.$

The key idea is that, modulo \(2\), the problem becomes a finite carry automaton with only five possible carry values and therefore only \(2^5=32\) parity states.

Step 1: Reduce the problem modulo \(2\)

If \(a_m(k)\) denotes the coefficient of \(x^m\) in \(P(x)^k\), then

$Q(k)=\#\{m\ge 0: a_m(k)\equiv 1 \pmod{2}\}.$

So the problem is not asking for exact coefficients. It only asks which coefficients survive after reduction in \(\mathbb{F}_2\).

Over \(\mathbb{F}_2\), squaring obeys

$f(x)^2=f(x^2).$

As a consequence, \(Q(2n)=Q(n)\). This observation is useful conceptually, even though the implementations evaluate the requested powers directly.

Step 2: Split the exponent into binary layers

Using the binary digits of \(k\), we can rewrite

$P(x)^k=\prod_{b\ge 0} P(x)^{k_b 2^b}\equiv \prod_{b:k_b=1} P(x^{2^b}) \pmod{2}.$

For one active bit \(b\), the factor is

$P(x^{2^b})=1+x^{2^b}+x^{2\cdot 2^b}+x^{3\cdot 2^b}+x^{4\cdot 2^b}.$

Choosing one term from each active factor is equivalent to choosing a digit \(d_b\). If \(k_b=0\), then necessarily \(d_b=0\). If \(k_b=1\), then

$d_b\in\{0,1,2,3,4\}.$

Every such choice contributes an exponent

$m=\sum_{b\ge 0} d_b 2^b.$

Hence the coefficient of \(x^m\) modulo \(2\) is the parity of the number of digit assignments that produce that same integer \(m\).

Step 3: Interpret the digit sum through carries

The digits \(d_b\) are not binary digits, so the representation above is not unique. The right model is binary addition with carries. Process bits from least significant to most significant. If the carry entering bit \(b\) is \(c\), then the local sum is

$s_b=c+d_b.$

The output bit of \(m\) at this position is

$m_b\equiv s_b \pmod{2},$

and the carry sent to the next position is

$c'=\left\lfloor\frac{s_b}{2}\right\rfloor.$

Because both \(c\) and \(d_b\) lie between \(0\) and \(4\), the next carry also lies in \(\{0,1,2,3,4\}\). So only five carry values ever occur.

Step 4: Store only parity information

For a fixed prefix of the output bits of \(m\), we do not need the exact number of partial digit assignments leading to each carry. We only need whether that number is odd or even.

So the state is a subset of \(\{0,1,2,3,4\}\): carry \(c\) is present exactly when the number of partial assignments reaching \(c\) is odd. Equivalently, the state is a 5-bit mask, so there are only \(32\) possible states.

Now fix \(t\in\{0,1\}\), where \(t=k_b\), and fix the output bit \(r\in\{0,1\}\). The allowed digits are

$D_t=\begin{cases} \{0\}, & t=0,\\ \{0,1,2,3,4\}, & t=1. \end{cases}$

From each active carry \(c\), inspect every \(d\in D_t\) satisfying \(c+d\equiv r\pmod{2}\). Each such pair produces the next carry \(c'=\lfloor(c+d)/2\rfloor\). If the same next carry appears twice, the two contributions cancel modulo \(2\), so the transition is naturally expressed by toggling bits in the next mask.

Step 5: Count all exponents by dynamic programming on masks

After the transition rule is known, we count output bit strings rather than coefficient assignments. Initially, before any bits are processed, the only active carry is \(0\). At each bit position there are two choices for the next output bit, namely \(0\) or \(1\), and each choice sends the current mask to one deterministic next mask.

If the highest set bit of \(k\) is \(B\), then the degree of \(P(x)^k\) is \(4k\), so only \(O(\log k)\) bit positions matter. After the significant bits of \(k\) are consumed, three extra zero-digit positions are enough to flush any remaining carry, because

$4\mapsto 2\mapsto 1\mapsto 0.$

Once all relevant bits are processed, an exponent contributes to \(Q(k)\) exactly when the final parity attached to carry \(0\) is odd. Therefore \(Q(k)\) is the total number of final masks in which carry \(0\) is active, counted with multiplicity from the outer dynamic program.

Step 6: Worked Example with \(k=3\)

Because \(3=11_2\), only bits \(0\) and \(1\) are active, so

$P(x)^3\equiv P(x)\,P(x^2)\pmod{2}.$

Here we choose digits \(d_0,d_1\in\{0,1,2,3,4\}\), and the exponent is

$m=d_0+2d_1.$

Some exponents appear an odd number of times and survive. Others appear an even number of times and disappear. For example, \(m=4\) occurs three times, namely from \((4,0)\), \((2,1)\), and \((0,2)\), so its coefficient is odd. But \(m=5\) occurs twice, from \((3,1)\) and \((1,2)\), so it cancels modulo \(2\).

Carrying this through for all pairs gives

$P(x)^3\equiv 1+x+x^4+x^6+x^8+x^{11}+x^{12}\pmod{2},$

hence

$Q(3)=7.$

This matches the small reference value used by the implementations. Two further checkpoints are \(Q(10)=17\) and \(Q(100)=35\).

How the Code Works

The C++, Python, and Java implementations all use the same structure. First they precompute the transition for every triple consisting of the current bit of \(k\), the requested output bit, and the current 5-bit carry-parity mask. Since there are only \(2\times 2\times 32\) such triples, this preprocessing is tiny.

For one concrete input \(k\), the implementation keeps an array of length \(32\). Entry \(M\) stores how many output prefixes currently lead to mask \(M\). The initial array has value \(1\) only at the mask representing carry \(0\) and value \(0\) everywhere else.

The program then scans bit positions from low to high. At each position it branches on the next output bit being \(0\) or \(1\), uses the precomputed transition to obtain the next mask, and adds the current count to the corresponding destination state. In other words, there is an outer count over output prefixes and an inner parity transition over carries.

The number of processed positions is the bit length of \(k\) plus three extra steps, which is sufficient to remove all remaining carry. After that, the implementation sums all state counts whose mask has the carry-\(0\) bit active. That sum is exactly \(Q(k)\). Repeating the same procedure for \(10,10^2,\dots,10^{18}\) and adding the eighteen values yields the required answer.

Complexity Analysis

The transition table has constant size. For one value of \(k\), if

$L=\lfloor\log_2 k\rfloor+3,$

then the dynamic program processes \(32\) masks at each of the \(L\) positions. Therefore the running time is \(O(32L)=O(\log k)\), and the working memory is \(O(32)=O(1)\).

The full Project Euler task requires only eighteen inputs, so the total cost stays extremely small. The important point is that the method depends on the bit length of \(k\), not on the huge degree \(4k\) of the expanded polynomial.

Footnotes and References

1. Problem page: https://projecteuler.net/problem=588
2. Finite field: Wikipedia - Finite field
3. Freshman's dream: Wikipedia - Freshman's dream
4. Lucas's theorem: Wikipedia - Lucas's theorem
5. Generating function: Wikipedia - Generating function