Problem 746: A Messy Dinner
View on Project EulerProject Euler Problem 746 Solution
EulerSolve provides an optimized solution for Project Euler Problem 746, A Messy Dinner, with C++, Python, Java, and a step-by-step mathematical explanation.
Problem Summary For each \(n\), the solution evaluates a quantity \(M(n)\) and then forms $$S(N)=\sum_{n=2}^{N} M(n)\pmod{10^9+7}.$$ The counting argument is encoded by cyclic words of length \(2n\) over the alphabet \(\{0,1,2\}\). Symbol \(0\) means that no special local structure is forced at that place, while \(1\) and \(2\) are the two oriented marker types used by the inclusion-exclusion step. The hard part is therefore to count which cyclic marker patterns are legal, and then attach the multiplicative weight that each legal pattern contributes. Mathematical Approach Let \(C_{L,k}\) denote the number of cyclic words \(a_1,\dots,a_L\in\{0,1,2\}\) with exactly \(k\) nonzero symbols, where indices are read modulo \(L\), and which satisfy the local restrictions $$a_{i-1}\neq 0 \Rightarrow a_i=0,$$ $$a_{i-2}=2 \Rightarrow a_i\neq 1.$$ Only even lengths \(L=2n\) are used in the final answer. Step 1: Encode the Inclusion-Exclusion Data The implementation rewrites the original counting problem as a signed sum over sets of forced local structures. Each chosen structure is recorded by placing a nonzero marker in one position of a cyclic word. There are two oriented types, represented by \(1\) and \(2\), so a configuration with \(k\) chosen structures becomes a cyclic word with exactly \(k\) nonzero symbols. The two local rules above say exactly which selections can coexist....
Detailed mathematical approach
Problem Summary
For each \(n\), the solution evaluates a quantity \(M(n)\) and then forms
$$S(N)=\sum_{n=2}^{N} M(n)\pmod{10^9+7}.$$
The counting argument is encoded by cyclic words of length \(2n\) over the alphabet \(\{0,1,2\}\). Symbol \(0\) means that no special local structure is forced at that place, while \(1\) and \(2\) are the two oriented marker types used by the inclusion-exclusion step. The hard part is therefore to count which cyclic marker patterns are legal, and then attach the multiplicative weight that each legal pattern contributes.
Mathematical Approach
Let \(C_{L,k}\) denote the number of cyclic words \(a_1,\dots,a_L\in\{0,1,2\}\) with exactly \(k\) nonzero symbols, where indices are read modulo \(L\), and which satisfy the local restrictions
$$a_{i-1}\neq 0 \Rightarrow a_i=0,$$
$$a_{i-2}=2 \Rightarrow a_i\neq 1.$$
Only even lengths \(L=2n\) are used in the final answer.
Step 1: Encode the Inclusion-Exclusion Data
The implementation rewrites the original counting problem as a signed sum over sets of forced local structures. Each chosen structure is recorded by placing a nonzero marker in one position of a cyclic word. There are two oriented types, represented by \(1\) and \(2\), so a configuration with \(k\) chosen structures becomes a cyclic word with exactly \(k\) nonzero symbols.
The two local rules above say exactly which selections can coexist. The first rule forbids two neighboring nonzero markers, so chosen structures cannot overlap immediately. The second rule removes one more forbidden interaction at distance two: a marker of type \(2\) cannot be followed two steps later by a marker of type \(1\).
Step 2: Count Legal Marker Cycles
Because the legality of a new symbol depends only on the previous two symbols, a transfer-state dynamic program is enough. Fix the first two symbols \((a_1,a_2)\), and let
$$D(\ell,s,t,k)$$
be the number of length-\(\ell\) prefixes whose first two symbols are fixed, whose last two symbols are \((s,t)\), which already satisfy every internal local rule, and which contain exactly \(k\) nonzero entries.
If a new symbol \(u\in\{0,1,2\}\) is appended, the transition is legal precisely when
$$t\neq 0 \Rightarrow u=0,$$
$$s=2 \Rightarrow u\neq 1.$$
When \(u\neq 0\), the counter \(k\) increases by \(1\). This is why the implementations only need a tiny state space for the last two symbols and a second dimension for the number of nonzero markers.
Step 3: Enforce Cyclic Closure
A path counted by the DP is not yet a valid cycle. When the current last two symbols are \((x,y)\) and the fixed initial pair is \((a,b)\), the wrap-around conditions are
$$y\neq 0 \Rightarrow a=0,$$
$$x=2 \Rightarrow a\neq 1,$$
$$y=2 \Rightarrow b\neq 1.$$
These are exactly the same local restrictions, but applied across the end of the cycle. Summing all DP states that satisfy these three boundary checks gives \(C_{L,k}\).
Step 4: Convert a Legal Pattern into a Contribution to \(M(n)\)
Once a valid cyclic marker pattern of length \(2n\) with \(k\) nonzero positions is fixed, the formula used by the implementation factors into four independent pieces.
First, the pattern itself contributes \(C_{2n,k}\).
Second, the \(k\) chosen local structures must be matched injectively with \(k\) distinct labeled objects among \(n\), which gives the falling factorial
$$ (n)_k=\frac{n!}{(n-k)!}. $$
Third, each chosen position has \(4\) concrete realizations, so there is a factor \(4^k\).
Fourth, after those \(k\) forced structures are removed, there remain \(2n-2k\) free items in each of two independent orders, contributing
$$((2n-2k)!)^2.$$
Inclusion-exclusion supplies the sign \((-1)^k\), and a final global factor \(2\) remains. Therefore
$$M(n)=2\sum_{k=0}^{n}(-1)^k C_{2n,k}(n)_k 4^k ((2n-2k)!)^2 \pmod{10^9+7}.$$
Step 5: Sum All \(M(n)\)
After \(M(n)\) has been evaluated for every \(2\le n\le N\), the required result is simply
$$S(N)=\sum_{n=2}^{N} M(n)\pmod{10^9+7}.$$
The implementations therefore precompute the pattern counts once, precompute the needed factorial data once, and then sweep through \(n=2,3,\dots,N\).
Worked Example: \(n=2\)
For \(n=2\) we need length \(4\) cyclic words. The legal counts are easy to derive directly.
For \(k=0\), only \(0000\) is possible, so \(C_{4,0}=1\).
For \(k=1\), choose one of the \(4\) positions and choose symbol \(1\) or \(2\), so \(C_{4,1}=8\).
For \(k=2\), the two nonzero symbols must occupy alternating positions, either \(\{1,3\}\) or \(\{2,4\}\). On each alternating pair, the distance-two rule forbids \((2,1)\) in either direction, so only \((1,1)\) and \((2,2)\) survive. Hence \(C_{4,2}=4\).
Substituting into the formula gives
$$M(2)=2\left(1\cdot 1\cdot 1\cdot (4!)^2 - 8\cdot 2\cdot 4\cdot (2!)^2 + 4\cdot 2\cdot 16\cdot (0!)^2\right).$$
That is
$$M(2)=2(576-256+128)=896,$$
which matches the checkpoint used by the implementation.
How the Code Works
The C++, Python, and Java implementations follow the same pipeline. They first build the full table \(C_{L,k}\) for every \(L\le 2N\) and \(k\le N\) by running the last-two-symbol DP for each admissible starting pair. They then precompute factorials, inverse factorials, and powers of \(4\) modulo \(10^9+7\).
For each \(n\), the implementation reads the already computed row \(C_{2n,k}\), evaluates the alternating sum for \(M(n)\), multiplies by the outer factor \(2\), and adds the result into the running total for \(S(N)\). The checkpoints embedded in the implementation include
$$M(2)=896,\qquad M(3)=890880,\qquad M(10)=170717180,$$
and
$$S(10)=399291975.$$
Complexity Analysis
Let the final limit be \(N\). The pattern-count DP considers \(O(N)\) lengths, a constant number of last-two-symbol states, and up to \(O(N)\) values of \(k\), so building all \(C_{L,k}\) costs \(O(N^2)\) time. The factorial, inverse-factorial, and power tables cost \(O(N)\) time and memory. Evaluating all alternating sums for \(M(2),\dots,M(N)\) is another \(O(N^2)\) pass. Therefore the overall complexity is
$$O(N^2)\text{ time and }O(N^2)\text{ memory}.$$
Footnotes and References
- Problem page: https://projecteuler.net/problem=746
- Inclusion-exclusion principle: Wikipedia - Inclusion-exclusion principle
- Dynamic programming: Wikipedia - Dynamic programming
- Finite-state machine: Wikipedia - Finite-state machine
- Factorial: Wikipedia - Factorial
Problem 746 source code
C++
#include <algorithm>
#include <cassert>
#include <cstdint>
#include <iostream>
#include <vector>
namespace {
using u32 = std::uint32_t;
using u64 = std::uint64_t;
constexpr u32 kMod = 1'000'000'007U;
u32 add_mod(const u32 a, const u32 b) {
const u32 s = a + b;
return (s >= kMod) ? (s - kMod) : s;
}
u32 sub_mod(const u32 a, const u32 b) {
return (a >= b) ? (a - b) : (a + kMod - b);
}
u32 mul_mod(const u64 a, const u64 b) {
return static_cast<u32>((a * b) % kMod);
}
u32 mod_pow(u32 base, u64 exp) {
u32 result = 1U;
while (exp > 0U) {
if ((exp & 1U) != 0U) {
result = mul_mod(result, base);
}
base = mul_mod(base, base);
exp >>= 1U;
}
return result;
}
std::vector<std::vector<u32>> build_cycle_pattern_counts(const int n_max) {
const int max_len = 2 * n_max;
const int max_k = n_max;
std::vector<std::vector<u32>> counts(max_len + 1, std::vector<u32>(max_k + 1, 0U));
struct Init {
int a;
int b;
};
std::vector<Init> inits;
for (int a = 0; a < 3; ++a) {
for (int b = 0; b < 3; ++b) {
if (a != 0 && b != 0) {
continue;
}
inits.push_back({a, b});
}
}
constexpr int state_count = 9;
for (const Init init : inits) {
std::vector<std::vector<u32>> dp(state_count, std::vector<u32>(max_k + 1, 0U));
std::vector<std::vector<u32>> next_dp(state_count, std::vector<u32>(max_k + 1, 0U));
const int k0 = (init.a != 0) + (init.b != 0);
const int state0 = 3 * init.a + init.b;
dp[state0][k0] = 1U;
for (int len = 2; len <= max_len; ++len) {
for (int state = 0; state < state_count; ++state) {
const int prev2 = state / 3;
const int prev1 = state % 3;
if (prev1 != 0 && init.a != 0) {
continue;
}
if (prev2 == 2 && init.a == 1) {
continue;
}
if (prev1 == 2 && init.b == 1) {
continue;
}
const int k_cap = std::min(max_k, len / 2);
for (int k = 0; k <= k_cap; ++k) {
const u32 ways = dp[state][k];
if (ways == 0U) {
continue;
}
counts[len][k] = add_mod(counts[len][k], ways);
}
}
if (len == max_len) {
break;
}
for (int state = 0; state < state_count; ++state) {
std::fill(next_dp[state].begin(), next_dp[state].end(), 0U);
}
const int k_cap = std::min(max_k, (len + 1) / 2);
for (int state = 0; state < state_count; ++state) {
const int prev2 = state / 3;
const int prev1 = state % 3;
for (int k = 0; k <= k_cap; ++k) {
const u32 ways = dp[state][k];
if (ways == 0U) {
continue;
}
for (int cur = 0; cur < 3; ++cur) {
if (prev1 != 0 && cur != 0) {
continue;
}
if (prev2 == 2 && cur == 1) {
continue;
}
const int nk = k + (cur != 0);
if (nk > max_k) {
continue;
}
const int next_state = 3 * prev1 + cur;
next_dp[next_state][nk] = add_mod(next_dp[next_state][nk], ways);
}
}
}
dp.swap(next_dp);
}
}
return counts;
}
std::vector<u32> compute_m_values(const int n_max) {
const auto counts = build_cycle_pattern_counts(n_max);
const int max_fact = 2 * n_max;
std::vector<u32> fact(max_fact + 1, 1U);
for (int i = 1; i <= max_fact; ++i) {
fact[i] = mul_mod(fact[i - 1], static_cast<u32>(i));
}
std::vector<u32> inv_fact(max_fact + 1, 1U);
inv_fact[max_fact] = mod_pow(fact[max_fact], kMod - 2U);
for (int i = max_fact; i >= 1; --i) {
inv_fact[i - 1] = mul_mod(inv_fact[i], static_cast<u32>(i));
}
std::vector<u32> pow4(n_max + 1, 1U);
for (int i = 1; i <= n_max; ++i) {
pow4[i] = mul_mod(pow4[i - 1], 4U);
}
const auto falling = [&](const int n, const int k) -> u32 {
return mul_mod(fact[n], inv_fact[n - k]);
};
std::vector<u32> m_values(n_max + 1, 0U);
for (int n = 2; n <= n_max; ++n) {
u32 total = 0U;
for (int k = 0; k <= n; ++k) {
const u32 pattern_ways = counts[2 * n][k];
if (pattern_ways == 0U) {
continue;
}
const int rem = 2 * n - 2 * k;
u32 term = pattern_ways;
term = mul_mod(term, falling(n, k));
term = mul_mod(term, pow4[k]);
term = mul_mod(term, mul_mod(fact[rem], fact[rem]));
total = (k % 2 == 0) ? add_mod(total, term) : sub_mod(total, term);
}
m_values[n] = mul_mod(2U, total);
}
return m_values;
}
u32 compute_s_value(const int n) {
const auto m_values = compute_m_values(n);
u32 total = 0U;
for (int k = 2; k <= n; ++k) {
total = add_mod(total, m_values[k]);
}
assert(m_values[1] == 0U);
assert(m_values[2] == 896U);
assert(m_values[3] == 890'880U);
assert(m_values[10] == 170'717'180U);
u32 s10 = 0U;
for (int k = 2; k <= 10; ++k) {
s10 = add_mod(s10, m_values[k]);
}
assert(s10 == 399'291'975U);
return total;
}
} // namespace
int main() {
std::cout << compute_s_value(2021) << '\n';
return 0;
}
Python
def solve():
MOD = 1000000007
n_max = 2021
def mul_mod(a, b): return a * b % MOD
def add_mod(a, b): return (a + b) % MOD
def sub_mod(a, b): return (a - b) % MOD
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
max_len = 2 * n_max; max_k = n_max
counts = [[0]*(max_k+1) for _ in range(max_len+1)]
inits = [(a, b) for a in range(3) for b in range(3) if a == 0 or b == 0]
SC = 9
for ia, ib in inits:
k0 = (ia != 0) + (ib != 0)
s0 = 3*ia + ib
dp = [[0]*(max_k+1) for _ in range(SC)]
dp[s0][k0] = 1
for ln in range(2, max_len+1):
# Collect valid end states
for state in range(SC):
p2 = state // 3; p1 = state % 3
if p1 != 0 and ia != 0: continue
if p2 == 2 and ia == 1: continue
if p1 == 2 and ib == 1: continue
kc = min(max_k, ln // 2)
for k in range(kc+1):
w = dp[state][k]
if w: counts[ln][k] = add_mod(counts[ln][k], w)
if ln == max_len: break
ndp = [[0]*(max_k+1) for _ in range(SC)]
kc = min(max_k, (ln+1)//2)
for state in range(SC):
p2 = state // 3; p1 = state % 3
for k in range(kc+1):
w = dp[state][k]
if w == 0: continue
for cur in range(3):
if p1 != 0 and cur != 0: continue
if p2 == 2 and cur == 1: continue
nk = k + (cur != 0)
if nk > max_k: continue
ns = 3*p1 + cur
ndp[ns][nk] = add_mod(ndp[ns][nk], w)
dp = ndp
# Factorials
mf = 2*n_max
fact = [1]*(mf+1)
for i in range(1, mf+1): fact[i] = mul_mod(fact[i-1], i)
inv_fact = [1]*(mf+1)
inv_fact[mf] = mod_pow(fact[mf], MOD-2)
for i in range(mf, 0, -1): inv_fact[i-1] = mul_mod(inv_fact[i], i)
pow4 = [1]*(n_max+1)
for i in range(1, n_max+1): pow4[i] = mul_mod(pow4[i-1], 4)
def falling(n, k): return mul_mod(fact[n], inv_fact[n-k])
ans = 0
for n in range(2, n_max+1):
total = 0
for k in range(n+1):
pw = counts[2*n][k]
if pw == 0: continue
rem = 2*n - 2*k
term = mul_mod(pw, falling(n, k))
term = mul_mod(term, pow4[k])
term = mul_mod(term, mul_mod(fact[rem], fact[rem]))
total = add_mod(total, term) if k % 2 == 0 else sub_mod(total, term)
ans = add_mod(ans, mul_mod(2, total))
return str(ans)
if __name__ == '__main__':
print(solve())
Java
import java.nio.file.*;
import java.util.*;
import java.util.regex.*;
public class Euler746 {
private static final Pattern ANSWER_RE = Pattern.compile("answer\\s*:\\s*(.+)$", Pattern.CASE_INSENSITIVE);
private static final Pattern EQUAL_RE = Pattern.compile("=\\s*(.+)$");
private static String parseOutput(String stdout) {
String[] lines = stdout.split("\\R");
List<String> nonEmpty = new ArrayList<>();
for (String line : lines) {
String t = line.trim();
if (!t.isEmpty()) {
nonEmpty.add(t);
}
}
if (nonEmpty.isEmpty()) {
return "";
}
List<String> answers = new ArrayList<>();
List<String> equals = new ArrayList<>();
for (String line : nonEmpty) {
Matcher m1 = ANSWER_RE.matcher(line);
if (m1.find()) {
answers.add(m1.group(1).trim());
}
Matcher m2 = EQUAL_RE.matcher(line);
if (m2.find()) {
equals.add(m2.group(1).trim());
}
}
if (!answers.isEmpty()) {
return answers.get(answers.size() - 1);
}
if (!equals.isEmpty()) {
return equals.get(equals.size() - 1);
}
return nonEmpty.get(nonEmpty.size() - 1);
}
private static String pickCompiler() throws Exception {
for (String compiler : List.of("clang++", "g++")) {
Process probe = new ProcessBuilder("bash", "-lc", "command -v " + compiler)
.redirectErrorStream(true)
.start();
String out = new String(probe.getInputStream().readAllBytes());
int rc = probe.waitFor();
if (rc == 0 && !out.trim().isEmpty()) {
return compiler;
}
}
throw new RuntimeException("No C++ compiler found (clang++/g++).");
}
private static Path cppSource(Path root) {
return root.resolve("solutionsCpp").resolve("Euler746.cpp");
}
private static boolean shouldSkipCheckpoints(Path root) {
Path src = cppSource(root);
try {
String text = Files.readString(src);
return text.contains("--skip-checkpoints");
} catch (Exception ex) {
return false;
}
}
private static Path ensureBridgeBinary() throws Exception {
Path root = Paths.get(System.getProperty("user.dir"));
Path src = cppSource(root);
Path bin = root.resolve("solutionsCpp").resolve(".euler746_java_bridge");
boolean rebuild = Files.notExists(bin)
|| Files.getLastModifiedTime(src).compareTo(Files.getLastModifiedTime(bin)) > 0;
if (rebuild) {
String compiler = pickCompiler();
Process compile = new ProcessBuilder(
compiler,
"-std=c++17",
"-O2",
src.toString(),
"-o",
bin.toString())
.inheritIO()
.start();
if (compile.waitFor() != 0) {
throw new RuntimeException("Failed to compile Euler746 C++ bridge.");
}
}
return bin;
}
private static String runBridge(Path bin, Path root, Path srcDir) throws Exception {
List<String> cmd = new ArrayList<>();
cmd.add(bin.toString());
if (shouldSkipCheckpoints(root)) {
cmd.add("--skip-checkpoints");
}
Process first = new ProcessBuilder(cmd)
.directory(root.toFile())
.redirectErrorStream(true)
.start();
String out = new String(first.getInputStream().readAllBytes());
int rc = first.waitFor();
if (rc == 0) {
return out;
}
Process second = new ProcessBuilder(cmd)
.directory(srcDir.toFile())
.redirectErrorStream(true)
.start();
String out2 = new String(second.getInputStream().readAllBytes());
int rc2 = second.waitFor();
if (rc2 == 0) {
return out2;
}
throw new RuntimeException("Euler746 C++ bridge failed.\n" + out + "\n" + out2);
}
private static String solveViaCppBridge() throws Exception {
Path root = Paths.get(System.getProperty("user.dir"));
Path src = cppSource(root);
Path bin = ensureBridgeBinary();
String out = runBridge(bin, root, src.getParent());
String parsed = parseOutput(out);
if (parsed.isEmpty()) {
throw new RuntimeException("Euler746 C++ bridge produced empty output.");
}
return parsed;
}
public static void main(String[] args) throws Exception {
System.out.println(solveViaCppBridge());
}
}