Detailed mathematical approach

Problem Summary

The task is to count rotationally distinct valid colourings of a two-row loop of length \(n\) using \(k\) colours, modulo \(M=1000004321\). A direct search over all full loops is hopeless for the target input, so the implementation scans the loop one column at a time and remembers only the local boundary information needed to continue a legal colouring.

The admissible local patterns are exactly those encoded by the transition system: a colour region may continue horizontally on the top row or the bottom row for at most three columns, or a column may be a single vertical monochromatic domino. Ordinary columns have different top and bottom colours, while a vertical domino has equal top and bottom colours. A simultaneous restart of both rows is allowed only immediately after a vertical domino, which is how the automaton excludes forbidden four-corner junctions.

Mathematical Approach

Let \(F_k(n)\) denote the number of rotational classes of valid colourings of length \(n\). The solution has three layers: encode local legality with a finite automaton, turn rooted colourings into matrix traces, and then quotient by rotation with Burnside's lemma.

Step 1: Encode the Boundary State

Walk around the loop column by column and record a state

$\bigl(r_t,r_b,p,c_t,c_b\bigr).$

Here \(r_t,r_b\in\{0,1,2\}\) are the numbers of future columns still occupied by the current top and bottom horizontal runs after the present column, \(c_t,c_b\) are the current colours, and \(p\in\{0,1\}\) marks whether the current column is a vertical monochromatic domino.

The validity conditions are

$p=0 \implies c_t\ne c_b,$

$p=1 \implies r_t=r_b=0 \land c_t=c_b.$

So the state space is finite, with

$S=9k(k-1)+k=k(9k-8)$

valid states.

Step 2: Build the Transition Matrix

From each state, the next column is determined by simple local continuation rules.

If \(r_t>0\) and \(r_b>0\), both horizontal runs continue, so the next state is obtained by decrementing both remainders. If exactly one row has ended, the unfinished row keeps its colour while the finished row starts a new run of length \(1\), \(2\), or \(3\), with a colour different from both its left neighbour and its vertical neighbour in the new column.

If both rows have ended, a fresh colour may create a vertical domino. Two independent new horizontal runs may also start, but only when the previous column was vertical; this is the local rule that removes four-corner meetings. These legal moves form a directed graph. Let \(A\) be its adjacency matrix.

Step 3: Count Rooted Colourings by Matrix Traces

A rooted colouring of length \(m\) that closes consistently after exactly \(m\) columns corresponds to a closed walk of length \(m\) in the state graph. Therefore the rooted count is

$T_m=\operatorname{tr}(A^m).$

This trace counts closed walks because the diagonal entries of \(A^m\) count returns to the same state after \(m\) steps. Since \(A\) is a finite matrix, the sequence \((T_m)\) satisfies a linear recurrence modulo \(M\).

Step 4: Recover the Linear Recurrence

The implementation first computes a prefix

$T_0,T_1,T_2,\dots$

by repeatedly multiplying the current matrix power by the sparse transition graph and taking the trace after each step. Once enough terms are available, Berlekamp-Massey reconstructs the shortest recurrence

$T_{m+L}=c_1T_{m+L-1}+c_2T_{m+L-2}+\cdots+c_LT_m \pmod M.$

After that, very large indices can be evaluated quickly using binary exponentiation in the recurrence algebra, so the huge target value of \(n\) is reduced to a logarithmic-time query in \(n\).

Step 5: Quotient by Rotation with Burnside's Lemma

The loop is cyclic, so rooted colourings that differ only by a rotation represent the same final object. Burnside's lemma gives

$F_k(n)=\frac{1}{n}\sum_{s=0}^{n-1}\operatorname{Fix}(s),$

where \(\operatorname{Fix}(s)\) is the number of colourings fixed by rotation through \(s\) columns. A colouring fixed by that rotation has period \(\gcd(n,s)\), hence

$\operatorname{Fix}(s)=T_{\gcd(n,s)}.$

Grouping rotations by their order yields the divisor form used in the program:

$F_k(n)=\frac{1}{n}\sum_{d\mid n}\varphi(d)\,T_{n/d}\pmod M.$

Worked Example: The Case \(n=3\)

For a loop of length \(3\), there are exactly three rotations. The identity fixes all rooted colourings of period \(3\), contributing \(T_3\). The other two rotations both force period \(1\), so each contributes \(T_1\). Therefore

$F_k(3)=\frac{T_3+2T_1}{3}.$

This is a clean small example of the divisor formula, because \(3\) is prime. For \(k=4\), the implementation checks that

$F_4(3)=104,$

which is exactly the expected Burnside average for that case.

How the Code Works

The C++, Python, and Java implementations all follow the same pipeline. First they enumerate every valid boundary state and build the directed adjacency list of legal next states. This adjacency list is the sparse representation of the transfer matrix.

Next they start from the identity matrix, repeatedly apply one more transition step, and record the traces \(T_0,T_1,\dots\). Once the trace sequence is long enough, they recover a minimal linear recurrence modulo \(M\) and keep the initial values needed to seed that recurrence.

Finally they factor \(n\), enumerate all divisors together with their Euler totient weights, evaluate \(T_{n/d}\) from the recurrence for each divisor \(d\), sum

$\varphi(d)\,T_{n/d},$

and multiply by the modular inverse of \(n\). The published answer is that Burnside average modulo \(M\).

Complexity Analysis

Let \(S\) be the number of states, \(E\) the number of directed transitions, and \(L\) the recovered recurrence order. Building the automaton takes \(O(S+E)\) time and memory.

If a prefix of \(m\) trace values is generated, the repeated dense-by-sparse matrix updates cost roughly \(O(mSE)\) time and \(O(S^2+E)\) memory, because the current matrix power is stored explicitly while the transition graph is sparse.

Recovering a recurrence from a prefix of length \(O(L)\) is quadratic in \(L\) in the standard Berlekamp-Massey formulation. After that, each large-index query \(T_x\) costs \(O(L^2\log x)\), and the final Burnside sum needs one such query for each divisor of \(n\). Thus the cyclic count after recurrence discovery is

$O\bigl(\tau(n)L^2\log n\bigr),$

plus the comparatively small cost of factoring \(n\) and enumerating its divisors.

Footnotes and References

1. Problem page: Project Euler 671
2. Transfer-matrix method: Wikipedia - Transfer-matrix method
3. Burnside's lemma: Wikipedia - Burnside's lemma
4. Euler's totient function: Wikipedia - Euler's totient function
5. Berlekamp-Massey algorithm: Wikipedia - Berlekamp-Massey algorithm
6. Cayley-Hamilton theorem: Wikipedia - Cayley-Hamilton theorem