#!/usr/bin/env python3
"""M10: Tseitin-Formeln (XOR ohne Redundanz) vs. 3-SAT-Kerne gleicher Variablenzahl.

Tseitin auf einem zusammenhaengenden 3-regulaeren Graphen mit V Knoten: Variable = Kante (n = 3V/2),
je Knoten eine Gleichung  XOR(3 Kanten) = Ladung,  Summe der Ladungen ungerade -> unloesbar.
Auf zusammenhaengenden Graphen ist die GANZE Formel minimal unloesbar (jede Klausel noetig): kein kleiner
Widerspruch, keine Redundanz -> direkt vergleichbar mit 3-SAT-Kernen (MUS) aus M9.
Graphen:
  expander  zufaellig 3-regulaer (Konfigurationsmodell, zusammenhaengend) -> fuer Resolution beweisbar exponentiell
            (Urquhart 1987; Ben-Sasson/Wigderson 2001)
  prisma    Prismengraph C_k x K2 (schmale Baumweite) -> Positivkontrolle: Resolution/DP sollte hier flach bleiben
Je Formel:  l2cnf, weite (DP minres wie ebene1_dp), aufbl = l2cnf - log2(Klauselzahl), ns (NS-Grad GF(2), erwartet 3)
Vergleich:  3-SAT-Kerne aus m9u.csv (oder --m9), nach Kern-Variablenzahl gebinnt, gleiche Aufblaehungs-Definition.
Nutzung:  python3 m10_tseitin.py --vs 8,10,12,14,16 --inst 60 --tmax 60 --jobs 15 --m9 m9u.csv
"""
import argparse, csv, math, sys, time
import numpy as np
from multiprocessing import Pool

def graph_expander(V, rng):
    while True:
        pts = np.repeat(np.arange(V), 3); rng.shuffle(pts)
        E = [tuple(sorted((int(pts[2 * j]), int(pts[2 * j + 1])))) for j in range(len(pts) // 2)]
        if any(a == b for a, b in E) or len(set(E)) != len(E): continue
        adj = {v: [] for v in range(V)}
        for a, b in E: adj[a].append(b); adj[b].append(a)
        seen, st = {0}, [0]
        while st:
            for w in adj[st.pop()]:
                if w not in seen: seen.add(w); st.append(w)
        if len(seen) == V: return E

def graph_prisma(V, rng):
    k = V // 2; E = []
    for i in range(k):
        E += [(i, (i + 1) % k), (k + i, k + (i + 1) % k), (i, k + i)]
    return [tuple(sorted(e)) for e in E]

def tseitin(V, E, rng):
    perm = rng.permutation(len(E))                      # zufaellige Variablennummerierung
    var = {E[p]: int(j) + 1 for j, p in enumerate(perm)}
    inc = {v: [] for v in range(V)}
    for e in E:
        inc[e[0]].append(var[e]); inc[e[1]].append(var[e])
    ch = [int(x) for x in rng.integers(2, size=V)]
    if sum(ch) % 2 == 0: ch[-1] ^= 1                   # ungerade Gesamtladung -> UNSAT
    cls = []
    for v in range(V):
        xs3 = inc[v]
        for bits in range(8):
            xs = [(bits >> j) & 1 for j in range(3)]
            if (xs[0] ^ xs[1] ^ xs[2]) != ch[v]:
                cls.append(tuple((-x if b else x) for x, b in zip(xs3, xs)))
    return len(E), cls

G = {}
def _init(a): G.update(a)

def job(key):
    typ, V, i, seed = key
    import ebene1_dp, m8_nsgrad as N
    rng = np.random.default_rng(seed)
    E = graph_expander(V, rng) if typ == "expander" else graph_prisma(V, rng)
    n, cls = tseitin(V, E, rng)
    minimal = float("nan")
    if G["pruefen"] and n <= 21:
        import m7_bewegung as M, m9_kern as K
        verl = M.verletzt_sat3(M.bits(n), cls); k = K.mus(n, verl, range(len(cls)))
        minimal = float(k is not None and len(k) == len(cls))
    mx, capped, _, w = ebene1_dp.dp(n, [frozenset(c) for c in cls], "minres", G["cap"], G["tmax"])
    l2 = math.log2(max(mx, 1))
    ns, nsg, _ = N.ns_grad(n, cls, G["dmax"], G["ns_tmax"])
    return dict(typ=typ, V=V, i=i, n=n, m=len(cls), l2cnf=l2, aufbl=l2 - math.log2(len(cls)), weite=w,
                gek=int(capped), ns=ns, nsgek=nsg, minimal=minimal)

def med(v): return float(np.median(v)) if v else float("nan")
def q90(v): return float(np.quantile(v, 0.9)) if v else float("nan")

def main():
    ap = argparse.ArgumentParser()
    ap.add_argument("--vs", default="8,10,12,14,16", help="Knotenzahlen (gerade); Variablen = 1.5*V")
    ap.add_argument("--inst", type=int, default=60)
    ap.add_argument("--inst-prisma", type=int, default=10)
    ap.add_argument("--jobs", type=int, default=15)
    ap.add_argument("--seed", type=int, default=11)
    ap.add_argument("--cap", type=int, default=400000)
    ap.add_argument("--tmax", type=float, default=60.0)
    ap.add_argument("--dmax", type=int, default=4)
    ap.add_argument("--ns-tmax", type=float, default=120.0)
    ap.add_argument("--pruefen", action="store_true", help="Minimalitaet per Brute Force pruefen (n <= 21)")
    ap.add_argument("--m9", default="m9u.csv", help="M9-CSV mit 3-SAT-Kernen zum Vergleich ('' = ohne)")
    ap.add_argument("--csv", default="m10_tseitin.csv")
    ap.add_argument("--report", default="report_m10.txt")
    a = ap.parse_args()
    vs = [int(x) for x in a.vs.split(",")]
    keys = [(t, V, i, a.seed * 1000003 + V * 7919 + i * 104729 + (0 if t == "expander" else 77777))
            for V in vs for t, k in (("expander", a.inst), ("prisma", a.inst_prisma)) for i in range(k)]
    keys.sort(key=lambda k: -k[1])
    R, t0 = [], time.time()
    with Pool(a.jobs, initializer=_init, initargs=(dict(cap=a.cap, tmax=a.tmax, dmax=a.dmax, ns_tmax=a.ns_tmax, pruefen=a.pruefen),)) as p:
        for k, r in enumerate(p.imap_unordered(job, keys), 1):
            R.append(r)
            if k % max(1, len(keys) // 40) == 0 or k == len(keys):
                el = time.time() - t0
                print(f"[{k}/{len(keys)}] {el/60:6.1f} min, Rest ~{el/k*(len(keys)-k)/60:6.1f} min", file=sys.stderr, flush=True)
    with open(a.csv, "w", newline="") as f:
        w = csv.DictWriter(f, fieldnames=list(R[0].keys())); w.writeheader(); w.writerows(R)

    L = [f"M10  seed={a.seed}  inst={a.inst}/{a.inst_prisma}  tmax={a.tmax}s", "",
         "Tseitin (ganze Formel = Kern):",
         "typ        V   n   m | aufbl med  q90 | weite med | ns (Verteilung) | gek | minimal"]
    for t in ("expander", "prisma"):
        for V in vs:
            rr = [r for r in R if r["typ"] == t and r["V"] == V]
            if not rr: continue
            ok = [r for r in rr if not r["gek"]]
            nsv = {}
            for r in rr: nsv[r["ns"]] = nsv.get(r["ns"], 0) + 1
            mn = [r["minimal"] for r in rr if not math.isnan(r["minimal"])]
            L.append(f"{t:9s} {V:3d} {rr[0]['n']:3d} {rr[0]['m']:3d} | {med([r['aufbl'] for r in ok]):5.2f} {q90([r['aufbl'] for r in ok]):5.2f} | "
                     f"{med([r['weite'] for r in ok]):5.1f}     | {' '.join(f'{k}:{v}' for k, v in sorted(nsv.items())):15s} | "
                     f"{sum(r['gek'] for r in rr):3d} | {(f'{sum(mn)/len(mn):.2f}' if mn else '-')}")
        L.append("")
    if a.m9:
        try:
            S = [r for r in csv.DictReader(open(a.m9)) if r["typ"] == "sat3" and r["gek_k"] == "0"]
        except FileNotFoundError:
            S = []; L.append(f"({a.m9} nicht gefunden - kein 3-SAT-Vergleich)")
        if S:
            L.append(f"Vergleich bei gleicher Variablenzahl: 3-SAT-Kerne aus {a.m9} (alle alpha) vs Tseitin-Expander")
            L.append("Variablen | 3-SAT-Kern: N  aufbl med  q90  weite | Tseitin-Exp.: N  aufbl med  q90  weite")
            for V in vs:
                nv = 3 * V // 2
                ss = [r for r in S if abs(float(r["kvar"]) - nv) <= 1]
                a3 = [float(r["l2cnf_k"]) - math.log2(float(r["mus"])) for r in ss]
                tt = [r for r in R if r["typ"] == "expander" and r["V"] == V and not r["gek"]]
                L.append(f"  {nv:3d}+-1  | {len(ss):5d}  {med(a3):5.2f} {q90(a3):5.2f}  {med([float(r['weite_k']) for r in ss]):4.1f} | "
                         f"{len(tt):4d}  {med([r['aufbl'] for r in tt]):5.2f} {q90([r['aufbl'] for r in tt]):5.2f}  {med([r['weite'] for r in tt]):4.1f}")
            L.append("")
    L += ["Lesart:",
          " - minimal = Anteil Formeln, deren Klausel-MUS die ganze Formel ist (Brute Force, --pruefen). Soll 1.00 sein.",
          " - ns = 3 ueberall erwartet: Gauss in GF(2) braucht Grad 3 -> Ebene 2 vorhanden.",
          " - Expander-Tseitin blaeht bei gleicher Variablenzahl mindestens so stark auf wie 3-SAT-Kerne -> XOR ist auch",
          "   redundanzfrei fuer Resolution nicht leichter; Leichtigkeit liegt ganz auf Ebene 2.",
          " - Prisma bleibt flach, Expander nicht -> die Aufblaehung misst die Graphstruktur (Expansion), nicht die Relation.",
          " - gek > 0 -> aufbl nur ueber nicht gekappte Laeufe; viele Kappungen = Aufblaehung jenseits des Deckels (selbst ein Befund).",
          " - 3-SAT-Kerne mit > 20 Variablen gibt es in M9 nicht (n <= 20); groesste V daher ohne Vergleichszeile."]
    txt = "\n".join(L); print(txt); open(a.report, "w").write(txt + "\n")

if __name__ == "__main__":
    main()