Project Euler Problem 651
An infinitely long cylinder has its curved surface fully covered with different coloured but otherwise identical rectang
Solution
Answer: 448233151
Let the stickers form an $a\times b$ rectangular grid wrapped onto a cylinder:
- $a$ = period along the axis,
- $b$ = number of stickers around the circumference.
Because the pattern is periodic in both directions, the surface is effectively an $a\times b$ torus.
We must count colourings using exactly $m$ colours, modulo all geometric symmetries generated by
- translations along the axis,
- rotations around the cylinder,
- reflections.
The symmetry group is therefore
$$G \cong D_a \times D_b,$$
with size
$$|G| = 4ab.$$
The problem is a direct application of Burnside’s Lemma.
1. Burnside’s lemma
For a group element $g$, let $c(g)$ be the number of cycles of the induced permutation on the $ab$ sticker positions.
A colouring fixed by $g$ must be constant on each cycle, so with $k$ available colours there are
$$k^{c(g)}$$
fixed colourings.
To force the use of exactly $m$ colours, use inclusion–exclusion:
$$E_m(c) = \sum_{j=0}^{m} (-1)^{m-j} \binom{m}{j} j^c.$$
This equals $m!,S(c,m)$, where $S(c,m)$ is a Stirling number of the second kind.
Hence
$$f(m,a,b) = \frac1{4ab} \sum_{g\in G} E_m(c(g)).$$
So the whole task reduces to computing the cycle counts $c(g)$.
2. Cycle structures
A dihedral element is either:
- a rotation,
- or a reflection.
For a rotation by $r$ on a cycle of length $n$:
$$\gcd(n,r)$$
cycles occur, each of length $n/\gcd(n,r)$.
For reflections:
-
if $n$ is odd:
-
one fixed point,
-
$(n-1)/2$ transpositions;
-
if $n$ is even:
-
half the reflections have two fixed points,
-
half have none.
For product actions, if one factor has cycles of lengths $p_i$ and the other has lengths $q_j$, then the product permutation contributes
$$\sum_{i,j}\gcd(p_i,q_j)$$
cycles.
This yields efficient divisor-sum formulas and makes the computation practical even for
$$(F_{39},F_{40})=(63245986,102334155).$$
3. Python implementation
from math import gcd, isqrt, comb
MOD = 1_000_000_007
# ------------------------------------------------------------
# Euler phi
# ------------------------------------------------------------
phi_cache = {1: 1}
def phi(n):
if n in phi_cache:
return phi_cache[n]
x = n
result = n
p = 2
while p * p <= x:
if x % p == 0:
while x % p == 0:
x //= p
result -= result // p
p += 1
if x > 1:
result -= result // x
phi_cache[n] = result
return result
# ------------------------------------------------------------
# Divisors
# ------------------------------------------------------------
def divisors(n):
ds = []
for d in range(1, isqrt(n) + 1):
if n % d == 0:
ds.append(d)
if d * d != n:
ds.append(n // d)
return ds
# ------------------------------------------------------------
# Number of surjections onto exactly m colours
# ------------------------------------------------------------
def exact_colourings(cycles, m):
total = 0
for j in range(m + 1):
total += (
(-1) ** (m - j)
* comb(m, j)
* pow(j, cycles, MOD)
)
return total % MOD
# ------------------------------------------------------------
# Burnside count
# ------------------------------------------------------------
def f(m, a, b):
total = 0
# --------------------------------------------------------
# rotation / rotation
# --------------------------------------------------------
for d1 in divisors(a):
for d2 in divisors(b):
count = phi(a // d1) * phi(b // d2)
cycles = d1 * d2
total += count * exact_colourings(cycles, m)
# --------------------------------------------------------
# rotation / reflection
# --------------------------------------------------------
for d1 in divisors(a):
A = a // d1
count = phi(A)
if b % 2 == 1:
cycles = d1 * (
1 + ((b - 1) // 2) * gcd(A, 2)
)
total += b * count * exact_colourings(cycles, m)
else:
cycles_even = d1 * (b // 2 + 1)
cycles_odd = d1 * (b // 2)
total += (
(b // 2)
* count
* (
exact_colourings(cycles_even, m)
+ exact_colourings(cycles_odd, m)
)
)
# --------------------------------------------------------
# reflection / rotation
# --------------------------------------------------------
for d2 in divisors(b):
B = b // d2
count = phi(B)
if a % 2 == 1:
cycles = d2 * (
1 + ((a - 1) // 2) * gcd(B, 2)
)
total += a * count * exact_colourings(cycles, m)
else:
cycles_even = d2 * (a // 2 + 1)
cycles_odd = d2 * (a // 2)
total += (
(a // 2)
* count
* (
exact_colourings(cycles_even, m)
+ exact_colourings(cycles_odd, m)
)
)
# --------------------------------------------------------
# reflection / reflection
# --------------------------------------------------------
def reflection_types(n):
if n % 2 == 1:
# (fixed points, transpositions), multiplicity
return [((1, (n - 1) // 2), n)]
return [
((2, (n - 2) // 2), n // 2),
((0, n // 2), n // 2),
]
for (f1, t1), mult1 in reflection_types(a):
for (f2, t2), mult2 in reflection_types(b):
cycles = (
f1 * f2
+ f1 * t2
+ f2 * t1
+ 2 * t1 * t2
)
total += (
mult1
* mult2
* exact_colourings(cycles, m)
)
total %= MOD
return total * pow(4 * a * b, MOD - 2, MOD) % MOD
# ------------------------------------------------------------
# Fibonacci numbers
# ------------------------------------------------------------
F = [0, 1]
for _ in range(2, 41):
F.append(F[-1] + F[-2])
# ------------------------------------------------------------
# Required sum
# ------------------------------------------------------------
ans = 0
for i in range(4, 41):
ans += f(i, F[i - 1], F[i])
ans %= MOD
print(ans)
4. Verification on the given examples
The program reproduces all supplied checks:
$$f(2,2,3)=11$$
$$f(3,2,3)=56$$
$$f(2,3,4)=156$$
$$f(8,13,21)\equiv 49718354 \pmod{10^9+7}$$
$$f(13,144,233)\equiv 907081451 \pmod{10^9+7}$$
all exactly as stated.
5. Final computation
Evaluating
$$\sum_{i=4}^{40} f(i,F_{i-1},F_i) \pmod{1,000,000,007}$$
gives:
$$448233151$$
Answer: 448233151