# -*- coding: utf-8 -*-
"""
AUDIT 3c — close the kernel: is the 6th deformation direction a GAUGE direction?

Audit 3b: dim ker M = 6, but only 5 independent directions from scalings +
the y^2 z-shear. Sol's full graded gauge group (audit 5) contains FOUR more
generators I never linearized: source shears y -> y + q1*x*z,
z -> z + p1*(xy)*y^2, and target shears Y1 -> Y1 + w1*Y2^2,
Y2 -> Y2 + v1*Y1*Y3. Each preserves the weight grading; to FIRST ORDER in
epsilon each should stay inside the equivariant z-affine ansatz class.

If adding their tangent vectors raises the trivial span rank to 6 = dim ker M,
then ker M == Lie algebra of the graded gauge group, and Alpoge's map is
INFINITESIMALLY RIGID MODULO GRADED GAUGE in the full slope-free class.
"""
import sympy as sp
import json, time

t0 = time.time()
x, y, z = sp.symbols('x y z')
a = x*y
ps  = sp.symbols('p0:4');  s1s = sp.symbols('s10:14')
qs  = sp.symbols('q0:3');  s2s = sp.symbols('s20:23')
rs  = sp.symbols('r0:3');  s3s = sp.symbols('s30:32')
c   = sp.Symbol('c')
TH  = list(ps)+list(s1s)+list(qs)+list(s2s)+list(rs)+list(s3s)+[c]

def poly(cs): return sum(co*a**i for i, co in enumerate(cs))
F1s = y**2*poly(ps) + z*poly(s1s)
F2s = y + x*y**2*poly(qs) + x*z*poly(s2s)
F3s = x*poly(rs) + x**3*z*poly(s3s)
D  = sp.expand(sp.Matrix([F1s, F2s, F3s]).jacobian([x, y, z]).det() - c)
K  = sp.Poly(D, x, y, z).coeffs()

th = dict(zip(ps, [4, 7, 3, 0])); th.update(zip(s1s, [1, 3, 3, 1]))
th.update(zip(qs, [12, 9, 0]));   th.update(zip(s2s, [3, 6, 3]))
th.update(zip(rs, [2, -3, 0]));   th.update(zip(s3s, [-1, 0]))
th[c] = -2

M = sp.Matrix([[sp.diff(k, v).subs(th) for v in TH] for k in K])
ns = M.nullspace()
print(f"dim ker M = {len(ns)}", flush=True)

e = sp.Symbol('epsilon')
F0 = [sp.expand(F1s.subs(th)), sp.expand(F2s.subs(th)), sp.expand(F3s.subs(th))]

def lin(expr):
    """first-order truncation in epsilon"""
    return sp.expand(expr.subs(e, 0) + e*sp.diff(expr, e).subs(e, 0))

def extract1(Fs, cfun):
    """extract theta-tangent from a first-order-in-epsilon map; None if it
    leaves the class AT FIRST ORDER."""
    out = {}
    specs = [(Fs[0], [((j, j+2, 0), ps[j]) for j in range(4)] + [((j, j, 1), s1s[j]) for j in range(4)]),
             (sp.expand(Fs[1]-y), [((j+1, j+2, 0), qs[j]) for j in range(3)] + [((j+1, j, 1), s2s[j]) for j in range(3)]),
             (Fs[2], [((j+1, j, 0), rs[j]) for j in range(3)] + [((j+3, j, 1), s3s[j]) for j in range(2)])]
    for comp, slots in specs:
        p = sp.Poly(sp.expand(comp), x, y, z)
        table = {mo: co for mo, co in zip(p.monoms(), p.coeffs())}
        for mo, sym in slots:
            out[sym] = table.pop(mo, 0)
        if any(sp.expand(v) != 0 for v in table.values()):
            return None
    return [sp.diff(out.get(v_, 0), e).subs(e, 0) if v_ != c else sp.diff(cfun, e).subs(e, 0)
            for v_ in TH]

dirs, names = [], []
def add(Fe, cfun, name):
    v = extract1([lin(f) for f in Fe], cfun)
    if v is None:
        print(f"   [{name}] leaves class at first order", flush=True)
    else:
        dirs.append(v); names.append(name)

lam = 1 + e
add([f.subs({x: lam*x}, simultaneous=True) for f in F0], -2*lam, "pre-x scale")
add([f.subs({z: lam*z}, simultaneous=True) for f in F0], -2*lam, "pre-z scale")
add([lam*F0[0], F0[1], F0[2]], -2*lam, "post-F1 scale")
add([F0[0], F0[1], lam*F0[2]], -2*lam, "post-F3 scale")
add([f.subs({z: z + e*y**2}, simultaneous=True) for f in F0], -2 + 0*e, "src z-shear +e*y^2")
Fy = [f.subs({y: lam*y}, simultaneous=True) for f in F0]
add([Fy[0], sp.expand(Fy[1]/lam), Fy[2]], -2 + 0*e, "combined pre-y/post-F2")
# ---- the four gauge shears never linearized before ----
add([f.subs({y: y + e*x*z}, simultaneous=True) for f in F0], -2 + 0*e, "src y-shear +e*xz")
add([f.subs({z: z + e*a*y**2}, simultaneous=True) for f in F0], -2 + 0*e, "src z-shear +e*xy*y^2")
add([sp.expand(F0[0] + e*F0[1]**2), F0[1], F0[2]], -2 + 0*e, "tgt Y1 += e*Y2^2")
add([F0[0], sp.expand(F0[1] + e*F0[0]*F0[2]), F0[2]], -2 + 0*e, "tgt Y2 += e*Y1*Y3")

Tspan = sp.Matrix(dirs).T
rankT = Tspan.rank()
okT = all(all(sp.simplify(val) == 0 for val in M*sp.Matrix(d)) for d in dirs)
print(f"gauge directions in-class: {len(dirs)} ({', '.join(names)})", flush=True)
print(f"rank of gauge span = {rankT};  all lie in ker M: {okT}", flush=True)

Kmat = sp.Matrix([list(n) for n in ns]).T
joint = Kmat.row_join(Tspan)
closed = (joint.rank() == rankT) and (rankT == len(ns))
print(f"dim ker M = {len(ns)}  vs  gauge rank = {rankT}  ; ker == gauge span: {closed}", flush=True)
verdict = "RIGID-MOD-GRADED-GAUGE" if closed else "DEFORMATION-BEYOND-GAUGE"
print(f"\n[{verdict}]", flush=True)
if closed:
    print("THEOREM (Infinitesimal Rigidity, final form): in the FULL torus-equivariant",
          "\nz-affine class (slope free, 19 coefficients + c), the kernel of the Keller",
          "\nlinearization at Alpoge's map equals the tangent space of the graded gauge",
          "\ngroup. Every infinitesimal Keller deformation is a gauge motion: the map is",
          "\nrigid modulo graded equivalence — matching Sol's global finding (audits 4-6)",
          "\nthat all apparent new specimens are gauge orbits of the original.", flush=True)
json.dump({"dim_ker": len(ns), "gauge_rank": rankT, "all_gauge_in_ker": bool(okT),
           "ker_equals_gauge": bool(closed), "directions": names,
           "verdict": verdict, "elapsed": round(time.time()-t0, 2)},
          open('post-jc-program/audit/audit3c.json', 'w'), indent=2)
print(f"artifact: audit3c.json  ({round(time.time()-t0,2)}s)", flush=True)