"""Thresholds for the y-branching of the lazy 2-separator rule (Appendix C.9 of the paper).

Lemma 13 bounds the decrease in the branch G - y from below by
    L = w_{d(y)} + sum_{c in C+} w_{d(c)} + zeta + o_y dw_Delta,
with zeta = w_{d(z)} - w_{o_z} if o_z >= 2 and zeta = w_{d(z)} + dw_Delta if o_z = 1, and the branching is valid when
L - [section 1] psi >= 1.  This script minimizes L exactly (fractions) over the parameters of a lazy configuration:
the admissible type of x, whether y and z are neighbors of x, d(y), d(z) in [3, Delta], o_y, o_z >= 1 with
o_y <= d(y) - 1 and o_z <= d(z) - 1 (y and z have a neighbor in C), and e further vertices of C+ outside N[x] of
degree >= 2, under the hypotheses of Lemma 15 (e >= e0, or, in section 1, a vertex of degree >= 3 in C+ \\ {x}).
Part C evaluates e and L for the 47 lazy pairs of G(80) (imported from g2_enumeration.py) as test cases.
Standard library only.  Run from the repository root:
    python3 -I experiments/algorithmic-consequences/g2_thresholds.py
"""
from fractions import Fraction as Fr
from itertools import combinations_with_replacement as cwr
import sys

SECTIONS = {
    1: ((Fr("0.0052182"), Fr("0.176755"), Fr("0.3360004")), Fr("0.0061458")),
    2: ((Fr("0.0536473"), Fr("0.176755"), Fr("0.2875713")), Fr(0)),
    3: ((Fr("0.0708019"), Fr("0.1796141"), Fr("0.2761349")), Fr(0)),
    4: ((Fr("0.089188"), Fr("0.1857428"), Fr("0.2700062")), Fr(0)),
    5: ((Fr("0.10766765"), Fr("0.19366265"), Fr("0.26736625")), Fr(0)),
}
GROUP = {(2, 2, 2, 2): 1, (2, 2, 2, 3): 1, (2, 2, 2, 4): 1, (2, 2, 3, 3): 1, (2, 2, 3, 4): 1, (2, 2, 4, 4): 1,
         (2, 3, 3, 3): 2, (2, 3, 3, 4): 3, (2, 3, 4, 4): 3, (2, 4, 4, 4): 3, (3, 3, 3, 3): 3, (3, 3, 3, 4): 3,
         (3, 3, 4, 4): 4, (3, 4, 4, 4): 4, (4, 4, 4, 4): 5}
HIGH = {
    "poly (eq:is-general-weights)": [Fr(0), Fr(0), Fr("0.1136655"), Fr("0.1986263"), Fr("0.2673669"), Fr("0.295498"), Fr("0.3034362")],
    "exp (eq:is-exp-weights)": [Fr(0), Fr(0), Fr("0.1103289"), Fr("0.1932144"), Fr("0.2500011"), Fr("0.2867833"), Fr("0.296828"), Fr("0.2970601")],
}
fails = []


def min_L(w, Delta, admissible, extra_degrees, need_deg3=False):
    """Minimum of L (and a minimizer) over the parameters; extra_degrees = degrees of the further vertices."""
    dw = w[Delta] - w[Delta - 1]
    best = None
    for T in cwr(range(2, Delta + 1), Delta):
        if not admissible(T):
            continue
        # positions of y, z in N(x): None (not adjacent) or an index with T[i] >= 3; dedup by the degree value
        choices = [None] + sorted({T[i] for i in range(Delta) if T[i] >= 3})
        for dy_in in choices:
            for dz_in in choices:
                core = list(T)
                if dy_in is not None:
                    core.remove(dy_in)
                if dz_in is not None:
                    if dz_in not in core:
                        continue
                    core.remove(dz_in)
                if need_deg3 and max(core + list(extra_degrees), default=2) < 3:
                    continue
                for dy in ([dy_in] if dy_in is not None else range(3, Delta + 1)):
                    for dz in ([dz_in] if dz_in is not None else range(3, Delta + 1)):
                        for oy in range(1, dy):
                            for oz in range(1, dz):
                                zeta = w[dz] - w[oz] if oz >= 2 else w[dz] + dw
                                L = w[dy] + w[Delta] + sum(w[d] for d in core) + sum(w[d] for d in extra_degrees) + zeta + oy * dw
                                desc = (T, "y in N(x)" if dy_in else "y not in N(x)", "z in N(x)" if dz_in else "z not in N(x)", dy, dz, oy, oz)
                                if best is None or L < best[0]:
                                    best = (L, desc)
    return best


def report(label, w, Delta, admissible, psi, hyp, extra_degrees, need_deg3=False, expect=True):
    L, desc = min_L(w, Delta, admissible, extra_degrees, need_deg3)
    val = L - psi
    ok = val >= 1
    print("   %s, %s: min L%s = %.7f (%s) at %s" % (label, hyp, " - psi" if psi else "", float(val), "valid" if ok else "NOT valid", desc))
    if ok != expect:
        fails.append((label, hyp))


print("Lemma 15: minima of the closed-form bound L of Lemma 13 (exact rationals)")
print("A. Maximum degree 4, sections 1 to 5 (D = 4)")
for k, ((w2, w3, w4), psi) in SECTIONS.items():
    w = [Fr(0), Fr(0), w2, w3, w4]
    adm = lambda t, k=k: GROUP[t] >= k
    lab = "section %d" % k
    if k == 1:
        report(lab, w, 4, adm, psi, "a vertex of degree >= 3 in C+ \\ {x} (core)", [], need_deg3=True)
        report(lab, w, 4, adm, psi, "a vertex of degree >= 3 in C+ \\ {x} (one further vertex of degree 3)", [3], need_deg3=True)
        for e in range(0, 4):
            report(lab, w, 4, adm, psi, "e = %d further vertices of degree 2, no degree-3 vertex" % e, [2] * e, expect=False)
    elif k in (2, 3):
        report(lab, w, 4, adm, psi, "e = 0", [], expect=False)
        report(lab, w, 4, adm, psi, "e = 1 (degree 2)", [2], expect=False)
        report(lab, w, 4, adm, psi, "e = 2 (degree 2)", [2, 2])
        report(lab, w, 4, adm, psi, "e = 1 (degree 3)", [3])
    else:
        report(lab, w, 4, adm, psi, "e = 0", [])
print("B. Maximum degree 5 and 6 (D = Delta), non-regular types")
for name, w in HIGH.items():
    for Delta in (5, 6):
        adm = lambda t, Delta=Delta: any(d < Delta for d in t)
        lab = "%s, Delta = %d" % (name, Delta)
        report(lab, w, Delta, adm, Fr(0), "e = 0", [], expect=False)
        e1_ok = not (name.startswith("exp") and Delta == 5)
        report(lab, w, Delta, adm, Fr(0), "e = 1 (degree 2)", [2], expect=e1_ok)
        if not e1_ok:
            report(lab, w, Delta, adm, Fr(0), "e = 2 (degree 2)", [2, 2])
print("C. Test cases: the lazy configurations of G(80) (section 2), e = |C+ \\ N[x]| and the closed-form bound L")
import os
sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)))
from g2_enumeration import counterexample, comps, mu_sections                      # noqa: E402
from itertools import combinations                                                   # noqa: E402
adj = counterexample(80, 0)
V = set(adj)
K = 65
w = [Fr(0), Fr(0)] + list(SECTIONS[2][0])
dw4 = w[4] - w[3]
nx = adj["x"] | {"x"}
by_e, minL = {}, {}
for p, q in combinations(sorted((u for u in V if len(adj[u]) >= 3), key=str), 2):
    cc = comps(adj, V - {p, q})
    small = [c for c in cc if len(c) <= K]
    if len(cc) < 2 or not any((nx - {p, q}) <= c for c in small):
        continue
    Cp = set().union(*small)
    e = len(Cp - nx)
    for a, b in ((p, q), (q, p)):                # branch on a, fold at b
        oa = len(adj[a] - Cp - {b})
        ob = len(adj[b] - Cp - {a})
        da, db = len(adj[a]), len(adj[b])
        zeta = w[db] - w[ob] if ob >= 2 else w[db] + dw4
        L = w[da] + sum(w[len(adj[c])] for c in Cp) + zeta + oa * dw4
        by_e[e] = by_e.get(e, 0) + 1
        if e not in minL or L < minL[e][0]:
            minL[e] = (L, a, b)
for e in sorted(by_e):
    print("   e = %2d: %2d oriented branchings, smallest L = %.7f (branch on %s, fold at %s)%s" % (e, by_e[e], float(minL[e][0]), minL[e][1], minL[e][2], "" if (e >= 2) == (minL[e][0] >= 1) or e < 2 else "  UNEXPECTED"))
    if e >= 2 and minL[e][0] < 1:
        fails.append(("G(80)", e))
print("ALL THRESHOLDS AS STATED" if not fails else "UNEXPECTED: %s" % fails)
sys.exit(1 if fails else 0)