Problem 458: Permutations of Project
View on Project EulerProject Euler Problem 458 Solution
EulerSolve provides an optimized solution for Project Euler Problem 458, Permutations of Project, with C++, Python, Java, and a step-by-step mathematical explanation.
Problem Summary Let \(T(n)\) be the number of length-\(n\) words over the seven distinct letters of “project” such that no contiguous block of length \(7\) is a permutation of those seven letters. We must compute \(T(10^{12}) \bmod 10^9\). Because the alphabet has size \(7\), a forbidden block of length \(7\) is exactly a block whose seven letters are all different. So the problem is to count words in which every length-\(7\) window contains at least one repeated letter. Mathematical Approach Step 1: Only the newest 7-window can become invalid Suppose we already have a valid word and append one more letter. Every old window of length \(7\) stays unchanged, so the only new constraint comes from the final window that ends at the appended position. Therefore, to know which letters may be appended, we only need information about the previous \(6\) characters. This gives the problem a finite-state structure: the future depends only on the suffix of length \(6\), not on the entire earlier prefix. Step 2: Replace concrete letters by a canonical equality pattern The actual names of the letters do not matter; only the equality pattern inside the last six positions matters. Two suffixes are equivalent if equal positions match equal positions after consistently renaming letters....
Detailed mathematical approach
Problem Summary
Let \(T(n)\) be the number of length-\(n\) words over the seven distinct letters of “project” such that no contiguous block of length \(7\) is a permutation of those seven letters. We must compute \(T(10^{12}) \bmod 10^9\).
Because the alphabet has size \(7\), a forbidden block of length \(7\) is exactly a block whose seven letters are all different. So the problem is to count words in which every length-\(7\) window contains at least one repeated letter.
Mathematical Approach
Step 1: Only the newest 7-window can become invalid
Suppose we already have a valid word and append one more letter. Every old window of length \(7\) stays unchanged, so the only new constraint comes from the final window that ends at the appended position. Therefore, to know which letters may be appended, we only need information about the previous \(6\) characters.
This gives the problem a finite-state structure: the future depends only on the suffix of length \(6\), not on the entire earlier prefix.
Step 2: Replace concrete letters by a canonical equality pattern
The actual names of the letters do not matter; only the equality pattern inside the last six positions matters. Two suffixes are equivalent if equal positions match equal positions after consistently renaming letters. For example, the suffix \(a\,b\,a\,c\,b\,d\) becomes the canonical pattern
$$\left(0,1,0,2,1,3\right).$$
Reading left to right, the first new letter receives label \(0\), the next unseen letter receives label \(1\), and so on. These canonical patterns are exactly the set partitions of a 6-element set, so the number of states is the Bell number
$$B_6 = 203.$$
That is why the implementations work with a \(203\)-state transfer matrix instead of all \(7^6\) concrete suffixes.
Step 3: Weighted transitions between states
Take a canonical suffix pattern using \(k\) distinct labels. This means the current suffix contains exactly \(k\) distinct letters.
If we append a letter already present in the suffix, the new 7-window still contains at most \(k \le 6\) distinct letters, so it is always valid. There are exactly \(k\) such concrete choices. After dropping the oldest character and canonicalizing the new suffix again, each of those choices contributes weight \(1\) to some next state.
If we append a letter not present in the suffix, then the new 7-window contains \(k+1\) distinct letters. This is allowed only when \(k \le 5\). There are \(7-k\) such letters, and after canonical renaming they all lead to the same next pattern, so that transition receives weight \(7-k\).
When \(k=6\), no new letter may be appended, because that would create a window with all seven letters distinct, which is exactly the forbidden case.
Therefore, if \(v_m\) is the column vector whose entries count valid words of length \(m\) ending in each canonical state, then
$$v_{m+1}=A\,v_m \pmod{10^9},$$
where \(A\) is the \(203\times203\) weighted transition matrix.
Step 4: Initial vector at length 6
For length \(6\), every word is valid because no length-\(7\) window exists yet. If a canonical pattern uses \(k\) labels, then we must assign \(k\) distinct actual letters to those labels. The number of such assignments is the falling factorial
$$ (7)_k = 7\cdot 6\cdot 5\cdots(7-k+1)=\frac{7!}{(7-k)!}. $$
So the entry of \(v_6\) corresponding to a state with \(k\) labels is exactly \((7)_k\).
Step 5: Fast exponentiation in the word length
For \(n \le 6\), we simply have
$$T(n)=7^n.$$
For \(n \ge 6\), repeated application of the transition matrix gives
$$v_n=A^{n-6}v_6 \pmod{10^9},\qquad T(n)=\sum_{s}(v_n)_s \pmod{10^9}.$$
The huge exponent \(n-6\) is handled by exponentiation by squaring, so the dependence on \(n\) is only logarithmic.
Checkpoint Example
For \(n=7\), the only invalid words are the permutations of the seven distinct letters, one for each ordering of the alphabet. Hence
$$T(7)=7^7-7!=823543-5040=818503.$$
Applying one more transfer step yields
$$T(8)=5699281,$$
which matches the checkpoints verified by the implementations.
How the Code Works
The C++, Python, and Java implementations generate all 203 canonical suffix patterns, index them, and build the weighted transition matrix described above. They then create the length-\(6\) start vector using the falling-factorial count \((7)_k\), raise the matrix to the power \(n-6\) by repeated squaring modulo \(10^9\), and finally sum the entries of the resulting state vector.
Complexity Analysis
Let \(S=203\) be the number of states. State generation and transition construction are fixed-size preprocessing for this problem. The dominant cost is matrix exponentiation, which requires \(O(S^3\log n)\) time and \(O(S^2)\) memory. Since \(S\) is small and fixed, this easily handles \(n=10^{12}\).
Footnotes and References
- Problem page: https://projecteuler.net/problem=458
- Bell numbers and set partitions: Wikipedia — Bell number
- Falling factorials: Wikipedia — Falling and rising factorials
- Transfer-matrix method: Wikipedia — Transfer-matrix method
- Exponentiation by squaring: Wikipedia — Exponentiation by squaring
Problem 458 source code
C++
#include <array>
#include <cstdint>
#include <iostream>
#include <string>
#include <unordered_map>
#include <vector>
namespace {
using u32 = std::uint32_t;
using u64 = std::uint64_t;
using u128 = __uint128_t;
constexpr u32 kMod = 1'000'000'000U;
constexpr int kLen = 6;
constexpr int kAlphabet = 7;
struct Options {
u64 n = 1'000'000'000'000ULL;
bool run_checkpoints = true;
};
bool parse_u64_after_prefix(const std::string& arg, const std::string& prefix, u64& out) {
if (arg.rfind(prefix, 0U) != 0U) {
return false;
}
const std::string tail = arg.substr(prefix.size());
if (tail.empty()) {
return false;
}
try {
out = static_cast<u64>(std::stoull(tail));
} catch (...) {
return false;
}
return true;
}
bool parse_arguments(int argc, char** argv, Options& options) {
for (int i = 1; i < argc; ++i) {
const std::string arg(argv[i]);
if (arg == "--skip-checkpoints") {
options.run_checkpoints = false;
continue;
}
if (parse_u64_after_prefix(arg, "--n=", options.n)) {
continue;
}
std::cerr << "Unknown argument: " << arg << '\n';
return false;
}
return true;
}
u32 mod_pow(u64 base, u64 exp) {
u64 result = 1ULL;
base %= kMod;
while (exp > 0ULL) {
if ((exp & 1ULL) != 0ULL) {
result = (result * base) % kMod;
}
base = (base * base) % kMod;
exp >>= 1ULL;
}
return static_cast<u32>(result);
}
using Pattern = std::array<std::uint8_t, kLen>;
u32 encode_pattern(const Pattern& p) {
u32 key = 0U;
for (int i = 0; i < kLen; ++i) {
key |= static_cast<u32>(p[static_cast<std::size_t>(i)]) << (3 * i);
}
return key;
}
Pattern canonicalize(const Pattern& p) {
std::array<int, kAlphabet + 1> remap{};
remap.fill(-1);
int next = 0;
Pattern out{};
for (int i = 0; i < kLen; ++i) {
const int v = p[static_cast<std::size_t>(i)];
if (remap[static_cast<std::size_t>(v)] == -1) {
remap[static_cast<std::size_t>(v)] = next++;
}
out[static_cast<std::size_t>(i)] =
static_cast<std::uint8_t>(remap[static_cast<std::size_t>(v)]);
}
return out;
}
void generate_patterns_recursive(int pos,
int current_max,
Pattern& current,
std::vector<Pattern>& out) {
if (pos == kLen) {
out.push_back(current);
return;
}
for (int v = 0; v <= current_max + 1; ++v) {
current[static_cast<std::size_t>(pos)] = static_cast<std::uint8_t>(v);
generate_patterns_recursive(pos + 1, (v > current_max ? v : current_max), current, out);
}
}
using Matrix = std::vector<std::vector<u32>>;
Matrix multiply_matrix(const Matrix& a, const Matrix& b) {
const std::size_t n = a.size();
Matrix c(n, std::vector<u32>(n, 0U));
for (std::size_t i = 0; i < n; ++i) {
for (std::size_t j = 0; j < n; ++j) {
u128 sum = 0;
for (std::size_t k = 0; k < n; ++k) {
sum += static_cast<u128>(a[i][k]) * b[k][j];
}
c[i][j] = static_cast<u32>(sum % kMod);
}
}
return c;
}
std::vector<u32> multiply_mat_vec(const Matrix& a, const std::vector<u32>& v) {
const std::size_t n = a.size();
std::vector<u32> out(n, 0U);
for (std::size_t i = 0; i < n; ++i) {
u128 sum = 0;
for (std::size_t j = 0; j < n; ++j) {
sum += static_cast<u128>(a[i][j]) * v[j];
}
out[i] = static_cast<u32>(sum % kMod);
}
return out;
}
u32 solve(const u64 n) {
if (n == 0ULL) {
return 1U;
}
if (n <= 6ULL) {
return mod_pow(7U, n);
}
std::vector<Pattern> patterns;
patterns.reserve(203);
Pattern seed{};
seed.fill(0U);
seed[0] = 0U;
generate_patterns_recursive(1, 0, seed, patterns);
std::unordered_map<u32, u32> id_by_key;
id_by_key.reserve(patterns.size() * 2U);
for (u32 id = 0U; id < static_cast<u32>(patterns.size()); ++id) {
id_by_key.emplace(encode_pattern(patterns[static_cast<std::size_t>(id)]), id);
}
const std::size_t state_count = patterns.size();
Matrix trans(state_count, std::vector<u32>(state_count, 0U));
for (u32 id = 0U; id < static_cast<u32>(state_count); ++id) {
const Pattern& p = patterns[static_cast<std::size_t>(id)];
int k = 0;
for (int i = 0; i < kLen; ++i) {
const int label = p[static_cast<std::size_t>(i)];
if (label + 1 > k) {
k = label + 1;
}
}
for (int label = 0; label < k; ++label) {
Pattern next{};
for (int t = 0; t < kLen - 1; ++t) {
next[static_cast<std::size_t>(t)] = p[static_cast<std::size_t>(t + 1)];
}
next[static_cast<std::size_t>(kLen - 1)] = static_cast<std::uint8_t>(label);
next = canonicalize(next);
const u32 nid = id_by_key[encode_pattern(next)];
trans[static_cast<std::size_t>(nid)][static_cast<std::size_t>(id)] =
(trans[static_cast<std::size_t>(nid)][static_cast<std::size_t>(id)] + 1U) % kMod;
}
const int new_count = kAlphabet - k;
if (new_count > 0 && k < 6) {
Pattern next{};
for (int t = 0; t < kLen - 1; ++t) {
next[static_cast<std::size_t>(t)] = p[static_cast<std::size_t>(t + 1)];
}
next[static_cast<std::size_t>(kLen - 1)] = static_cast<std::uint8_t>(k);
next = canonicalize(next);
const u32 nid = id_by_key[encode_pattern(next)];
trans[static_cast<std::size_t>(nid)][static_cast<std::size_t>(id)] =
(trans[static_cast<std::size_t>(nid)][static_cast<std::size_t>(id)] +
static_cast<u32>(new_count)) %
kMod;
}
}
std::vector<u32> vec(state_count, 0U);
for (u32 id = 0U; id < static_cast<u32>(state_count); ++id) {
const Pattern& p = patterns[static_cast<std::size_t>(id)];
int k = 0;
for (int i = 0; i < kLen; ++i) {
const int label = p[static_cast<std::size_t>(i)];
if (label + 1 > k) {
k = label + 1;
}
}
u32 ways = 1U;
for (int t = 0; t < k; ++t) {
ways = static_cast<u32>((static_cast<u64>(ways) * (kAlphabet - t)) % kMod);
}
vec[static_cast<std::size_t>(id)] = ways;
}
u64 steps = n - 6ULL;
Matrix power = trans;
while (steps > 0ULL) {
if ((steps & 1ULL) != 0ULL) {
vec = multiply_mat_vec(power, vec);
}
steps >>= 1ULL;
if (steps > 0ULL) {
power = multiply_matrix(power, power);
}
}
u64 answer = 0ULL;
for (const u32 v : vec) {
answer += v;
answer %= kMod;
}
return static_cast<u32>(answer);
}
bool run_checkpoints() {
if (solve(7ULL) != 818'503U) {
std::cerr << "Checkpoint failed: T(7)\n";
return false;
}
if (solve(8ULL) != 5'699'281U) {
std::cerr << "Checkpoint failed: T(8)\n";
return false;
}
return true;
}
} // namespace
int main(int argc, char** argv) {
Options options;
if (!parse_arguments(argc, argv, options)) {
return 1;
}
if (options.run_checkpoints && !run_checkpoints()) {
return 2;
}
std::cout << solve(options.n) << '\n';
return 0;
}
Python
def solve():
MOD = 1000000000
kLen = 6
kAlphabet = 7
n = 1000000000000
def mod_pow(base, exp):
r = 1; base %= MOD
while exp > 0:
if exp & 1: r = r * base % MOD
base = base * base % MOD
exp >>= 1
return r
if n <= 6: return str(mod_pow(7, n))
def canonicalize(p):
remap = {}; nxt = 0; out = []
for v in p:
if v not in remap: remap[v] = nxt; nxt += 1
out.append(remap[v])
return tuple(out)
def gen_patterns(pos, mx, cur, out):
if pos == kLen: out.append(tuple(cur)); return
for v in range(mx + 2):
cur.append(v)
gen_patterns(pos+1, max(mx,v), cur, out)
cur.pop()
patterns = []
gen_patterns(0, -1, [], patterns)
id_map = {p: i for i, p in enumerate(patterns)}
sc = len(patterns)
# Build transition matrix
trans = [[0]*sc for _ in range(sc)]
for idx, p in enumerate(patterns):
k = max(p) + 1
for label in range(k):
nxt = canonicalize(p[1:] + (label,))
nid = id_map[nxt]
trans[nid][idx] = (trans[nid][idx] + 1) % MOD
new_count = kAlphabet - k
if new_count > 0 and k < 6:
nxt = canonicalize(p[1:] + (k,))
nid = id_map[nxt]
trans[nid][idx] = (trans[nid][idx] + new_count) % MOD
# Initial vector
vec = [0] * sc
for idx, p in enumerate(patterns):
k = max(p) + 1; w = 1
for t in range(k): w = w * (kAlphabet - t) % MOD
vec[idx] = w
# Matrix-vector power
def mat_vec_mul(M, v):
r = [0]*len(v)
for i in range(len(v)):
s = 0
for j in range(len(v)):
s += M[i][j] * v[j]
r[i] = s % MOD
return r
def mat_mul(A, B):
sz = len(A)
C = [[0]*sz for _ in range(sz)]
for i in range(sz):
for k in range(sz):
if A[i][k] == 0: continue
for j in range(sz):
C[i][j] += A[i][k] * B[k][j]
for j in range(sz):
C[i][j] %= MOD
return C
steps = n - 6
power = trans
while steps > 0:
if steps & 1: vec = mat_vec_mul(power, vec)
steps >>= 1
if steps > 0: power = mat_mul(power, power)
return str(sum(vec) % MOD)
if __name__ == '__main__':
print(solve())
Java
import java.util.ArrayList;
import java.util.Arrays;
import java.util.HashMap;
import java.util.List;
import java.util.Map;
public class Euler458 {
static final long MOD = 1000000000L;
static final int kLen = 6;
static final int kAlphabet = 7;
static long modPow(long base, long exp) {
long result = 1;
base %= MOD;
while (exp > 0) {
if ((exp & 1) != 0) {
result = (result * base) % MOD;
}
base = (base * base) % MOD;
exp >>= 1;
}
return result;
}
static class Pattern {
byte[] p = new byte[kLen];
Pattern() {
}
Pattern(byte[] p) {
System.arraycopy(p, 0, this.p, 0, kLen);
}
@Override
public int hashCode() {
return Arrays.hashCode(p);
}
@Override
public boolean equals(Object obj) {
return Arrays.equals(p, ((Pattern) obj).p);
}
}
static Pattern canonicalize(Pattern p) {
int[] remap = new int[kAlphabet + 1];
Arrays.fill(remap, -1);
int next = 0;
Pattern out = new Pattern();
for (int i = 0; i < kLen; i++) {
int v = p.p[i];
if (remap[v] == -1) {
remap[v] = next++;
}
out.p[i] = (byte) remap[v];
}
return out;
}
static void generatePatterns(int pos, int currentMax, Pattern current, List<Pattern> out) {
if (pos == kLen) {
out.add(new Pattern(current.p));
return;
}
for (int v = 0; v <= currentMax + 1; v++) {
current.p[pos] = (byte) v;
generatePatterns(pos + 1, Math.max(v, currentMax), current, out);
}
}
static long[][] multiplyMatrix(long[][] a, long[][] b) {
int n = a.length;
long[][] c = new long[n][n];
for (int i = 0; i < n; i++) {
for (int j = 0; j < n; j++) {
long sum = 0;
for (int k = 0; k < n; k++) {
sum = (sum + a[i][k] * b[k][j]) % MOD;
}
c[i][j] = sum;
}
}
return c;
}
static long[] multiplyMatVec(long[][] a, long[] v) {
int n = a.length;
long[] out = new long[n];
for (int i = 0; i < n; i++) {
long sum = 0;
for (int j = 0; j < n; j++) {
sum = (sum + a[i][j] * v[j]) % MOD;
}
out[i] = sum;
}
return out;
}
public static String solve() {
long n = 1000000000000L;
if (n == 0)
return "1";
if (n <= 6)
return Long.toString(modPow(7, n));
List<Pattern> patterns = new ArrayList<>();
Pattern seed = new Pattern();
generatePatterns(1, 0, seed, patterns);
Map<Pattern, Integer> idByKey = new HashMap<>();
for (int i = 0; i < patterns.size(); i++) {
idByKey.put(patterns.get(i), i);
}
int stateCount = patterns.size();
long[][] trans = new long[stateCount][stateCount];
for (int i = 0; i < stateCount; i++) {
Pattern p = patterns.get(i);
int k = 0;
for (byte b : p.p) {
k = Math.max(k, b + 1);
}
for (int label = 0; label < k; label++) {
Pattern next = new Pattern();
for (int t = 0; t < kLen - 1; t++) {
next.p[t] = p.p[t + 1];
}
next.p[kLen - 1] = (byte) label;
next = canonicalize(next);
int nid = idByKey.get(next);
trans[nid][i] = (trans[nid][i] + 1) % MOD;
}
int newCount = kAlphabet - k;
if (newCount > 0 && k < 6) {
Pattern next = new Pattern();
for (int t = 0; t < kLen - 1; t++) {
next.p[t] = p.p[t + 1];
}
next.p[kLen - 1] = (byte) k;
next = canonicalize(next);
int nid = idByKey.get(next);
trans[nid][i] = (trans[nid][i] + newCount) % MOD;
}
}
long[] vec = new long[stateCount];
for (int i = 0; i < stateCount; i++) {
Pattern p = patterns.get(i);
int k = 0;
for (byte b : p.p) {
k = Math.max(k, b + 1);
}
long ways = 1;
for (int t = 0; t < k; t++) {
ways = (ways * (kAlphabet - t)) % MOD;
}
vec[i] = ways;
}
long steps = n - 6;
long[][] power = trans;
while (steps > 0) {
if ((steps & 1) != 0) {
vec = multiplyMatVec(power, vec);
}
steps >>= 1;
if (steps > 0) {
power = multiplyMatrix(power, power);
}
}
long answer = 0;
for (long v : vec) {
answer = (answer + v) % MOD;
}
return Long.toString(answer);
}
public static void main(String[] args) {
System.out.println(solve());
}
}