"""Independent exact re-verification of an antinorm certificate (does not use the stored
convex combinations).  The antinorm is ell(P) = min { n . Tc(P) : n >= 0, n . Tc(z) >= 1 for z in Z },
Tc = cumulative coordinates.  For Tc(P) >= 0 the minimum is attained at a vertex of the dual
polyhedron; we enumerate ALL vertices exactly (intersections of 3 of the |Z|+3 planes,
rational Gaussian elimination, feasibility check) and evaluate.
Checks: ell(z_i o z_j) >= 1 for all i, j; ell(Init(m)) >= e1 e2^m for m <= m0;
ell((y/2,0,0)) >= e1 e2^(m0+1); and reports min_i fin(z_i)."""
import itertools, json, sys
from fractions import Fraction as Fr

def mrg(P, Q):
    # parent P = (in X, out dominated-by-child, out undominated), child Q
    A, B, C = P; a, b, c = Q
    return (A * (2 * a + b + c), B * (a + b + c / 2) + C * a, C * (b + c / 2))

def cum(P):
    return (P[0], P[0] + P[1], P[0] + P[1] + P[2])

def solve3(M, r):
    M = [list(row) + [rv] for row, rv in zip(M, r)]
    for col in range(3):
        piv = next((i for i in range(col, 3) if M[i][col] != 0), None)
        if piv is None:
            return None
        M[col], M[piv] = M[piv], M[col]
        for i in range(3):
            if i != col and M[i][col] != 0:
                f = M[i][col] / M[col][col]
                M[i] = [x - f * yv for x, yv in zip(M[i], M[col])]
    return [M[i][3] / M[i][i] for i in range(3)]

def dual_vertices(TZ):
    planes = [(tuple(tz), Fr(1)) for tz in TZ] + [((Fr(1), Fr(0), Fr(0)), Fr(0)), ((Fr(0), Fr(1), Fr(0)), Fr(0)), ((Fr(0), Fr(0), Fr(1)), Fr(0))]
    verts = set()
    for (p1, p2, p3) in itertools.combinations(planes, 3):
        sol = solve3([p1[0], p2[0], p3[0]], [p1[1], p2[1], p3[1]])
        if sol is None:
            continue
        if any(x < 0 for x in sol):
            continue
        if all(sum(a * b for a, b in zip(tz, sol)) >= 1 for tz in TZ):
            verts.add(tuple(sol))
    return list(verts)

def ell(P, V):
    TP = cum(P)
    assert all(x >= 0 for x in TP)
    return min(sum(a * b for a, b in zip(v, TP)) for v in V)

if __name__ == '__main__':
    for f in sys.argv[1:]:
        cert = json.load(open(f))
        y = Fr(*cert['y']); e1 = Fr(*cert['e1']); e2 = Fr(*cert['e2']); m0 = cert['m0']
        Z = [tuple(Fr(*v) for v in z) for z in cert['Z']]
        assert all(min(z) >= 0 for z in Z) and 0 < e2 < 1 and e1 > 0
        V = dual_vertices([cum(z) for z in Z])
        assert all(abs(ell(z, V)) >= 1 for z in Z)
        worst = min(ell(mrg(zi, zj), V) for zi in Z for zj in Z)
        assert worst >= 1, worst
        for m in range(m0 + 1):
            assert ell((y / 2, Fr(0), Fr(1, 2 ** m)), V) >= e1 * e2 ** m, m
        assert ell((y / 2, Fr(0), Fr(0)), V) >= e1 * e2 ** (m0 + 1)
        c = min(z[0] + z[1] + z[2] / 2 for z in Z)
        print(f"{f}: INDEPENDENTLY VERIFIED ({len(Z)} points, {len(V)} dual vertices); min ell(z_i o z_j) = {float(worst):.8f}; c = {float(c):.6f}")
