"""Stage 5: exhaustiveness audit of the C oracle (code/brute_force.c) from
its recorded counts alone.

For every row (n, s) of results/tables/brute_force.csv:
  * 'subsets' must equal C(2n, s): the Gosper enumeration (and the m0
    splitting, where used) covered every s-subset;
  * 'orbit_representatives' must equal the number of D_n-orbits on
    s-subsets of V(P(n,5)), computed here by Burnside's lemma from the cycle
    structure of the 2n dihedral permutations (rotations: 2*gcd(n,a) cycles
    of length n/gcd(n,a); reflections i -> a-i: fixed vertices where 2i = a
    mod n, the rest in 2-cycles).  Equality means the canonical-form filter
    tested exactly one subset per orbit, so by P3 nothing was skipped.
Also records, per n, that the exhausted sizes form the full range 1..Z-1
(size monotonicity: a forcing set of size < Z-1 would extend to one of size
Z-1, so exhausting Z-1 alone already suffices; both are recorded).
Output: results/robustness/enumeration_counts_check.json
"""
import csv
import json
import math
import os

ROOT = os.path.dirname(os.path.dirname(os.path.dirname(os.path.abspath(__file__))))
os.chdir(ROOT)


def orbit_count(n, s):
    N = 2 * n
    tot = 0
    for a in range(n):  # rotations
        if a == 0:
            tot += math.comb(N, s)
            continue
        L = n // math.gcd(n, a)
        c = N // L
        tot += math.comb(c, s // L) if s % L == 0 else 0
    for a in range(n):  # reflections i -> a - i on both layers
        fp = 2 * sum(1 for i in range(n) if (2 * i - a) % n == 0)
        two = (N - fp) // 2
        tot += sum(math.comb(fp, j) * math.comb(two, (s - j) // 2) for j in range(0, min(fp, s) + 1) if (s - j) % 2 == 0)
    assert tot % (2 * n) == 0
    return tot // (2 * n)


def main():
    rows = list(csv.DictReader(open("results/tables/brute_force.csv")))
    Z = {int(r["n"]): int(r["Z"]) for r in csv.DictReader(open("results/exact_values.csv"))}
    out = {"rows": [], "mismatches": 0}
    per_n = {}
    for r in rows:
        n, s = int(r["n"]), int(r["s"])
        sub, rep, forcing = int(r["subsets"]), int(r["orbit_representatives"]), int(r["forcing_representatives"])
        exp_sub, exp_rep = math.comb(2 * n, s), orbit_count(n, s)
        ok = sub == exp_sub and rep == exp_rep
        out["rows"].append({"n": n, "s": s, "subsets": sub, "binomial": exp_sub, "representatives": rep, "burnside_orbits": exp_rep, "forcing": forcing, "ok": ok})
        out["mismatches"] += (not ok)
        per_n.setdefault(n, {})[s] = forcing
    out["per_n"] = []
    for n in sorted(per_n):
        sizes = per_n[n]
        z = Z.get(n)
        exhausted = sorted(s for s, f in sizes.items() if f == 0)
        found = sorted(s for s, f in sizes.items() if f > 0)
        out["per_n"].append({"n": n, "Z_table": z, "sizes_with_no_forcing_set": exhausted, "sizes_with_forcing_set": found,
                             "size_Z_minus_1_exhausted": (z - 1) in exhausted if z else None,
                             "all_sizes_1_to_Z_minus_1_exhausted": all(s in exhausted for s in range(1, z)) if z else None,
                             "consistent_with_Z": (z is not None and (z - 1) in exhausted and all(s >= z for s in found))})
    out["all_rows_ok"] = out["mismatches"] == 0
    out["all_n_consistent"] = all(p["consistent_with_Z"] for p in out["per_n"])
    json.dump(out, open("results/robustness/enumeration_counts_check.json", "w"), indent=1)
    print("rows", len(rows), "mismatches", out["mismatches"], "all n consistent with Z:", out["all_n_consistent"])
    for p in out["per_n"]:
        print(p["n"], p["Z_table"], "exhausted", p["sizes_with_no_forcing_set"], "found", p["sizes_with_forcing_set"], "all 1..Z-1:", p["all_sizes_1_to_Z_minus_1_exhausted"])


if __name__ == "__main__":
    main()
