"""Stage 5: how tight is the boundary lower bound?  Computes the exact
beta_18(n) for 26 <= n <= 40 by solving F(n,18,b) upward from b = 12 with
the project encoder (CaDiCaL); the certificate margin is beta_18(n) - 12.
Models are re-checked by direct counting.  Output: results/robustness/margin.json"""
import json, os, sys, time
ROOT = os.path.dirname(os.path.dirname(os.path.dirname(os.path.abspath(__file__)))); os.chdir(ROOT)
sys.path.insert(0, os.path.join(ROOT, ".pydeps")); sys.path.insert(0, os.path.join(ROOT, "code"))
from pysat.solvers import Solver
from sat_or_ilp_solver import boundary_clauses
def nbrs(n, k=5):
    A = {}
    for i in range(n):
        A[i] = {(i + 1) % n, (i - 1) % n, n + i}; A[n + i] = {n + (i + k) % n, n + (i - k) % n, i}
    return A
out = []
for n in range(26, 41):
    A = nbrs(n); t0 = time.time(); rec = {"n": n}
    for b in range(12, 20):
        cl, x, z, nv = boundary_clauses(n, 5, 18, b)
        with Solver(name="cadical195", bootstrap_with=cl) as sv:
            sat = sv.solve(); model = sv.get_model() if sat else None
        if sat:
            X = [w for w in range(2 * n) if model[x[w] - 1] > 0]; Xs = set(X)
            bd = len({w for q in X for w in A[q] if w not in Xs})
            rec.update({"beta_18": b, "model_boundary_recount": bd, "model_size": len(X), "margin": b - 12,
                        "X": ["%s%d" % ("u" if w < n else "v", w % n) for w in X]})
            break
    rec["seconds"] = round(time.time() - t0, 1); out.append(rec); print(rec, flush=True)
json.dump(out, open("results/robustness/margin.json", "w"), indent=1)