"""Independent exact verification of the Weyl and Poisson lifts.

Uses hand-rolled rational multivariate-polynomial arithmetic and shares no
code with a computer algebra system.
"""

# Independent cross-check of all identities using hand-rolled exact polynomial
# arithmetic over Q (dict-based, fractions.Fraction). Shares no code with SymPy.
from fractions import Fraction as Fr

# polynomial = dict {(ex,ey,ez): Fraction coeff}
def pclean(p): return {m: c for m, c in p.items() if c != 0}
def padd(p, q):
    r = dict(p)
    for m, c in q.items(): r[m] = r.get(m, Fr(0)) + c
    return pclean(r)
def pscale(p, s): return pclean({m: c*Fr(s) for m, c in p.items()})
def pmul(p, q):
    r = {}
    for m1, c1 in p.items():
        for m2, c2 in q.items():
            m = (m1[0]+m2[0], m1[1]+m2[1], m1[2]+m2[2])
            r[m] = r.get(m, Fr(0)) + c1*c2
    return pclean(r)
def pdiff(p, i):
    r = {}
    for m, c in p.items():
        if m[i] > 0:
            mm = list(m); e = mm[i]; mm[i] -= 1
            r[tuple(mm)] = r.get(tuple(mm), Fr(0)) + c*e
    return pclean(r)
def peval(p, pt):
    return sum(c * pt[0]**m[0] * pt[1]**m[1] * pt[2]**m[2] for m, c in p.items())

ONE = {(0,0,0): Fr(1)}; X = {(1,0,0): Fr(1)}; Y = {(0,1,0): Fr(1)}; Z = {(0,0,1): Fr(1)}
def lin(*terms):  # sum of (coeff, poly)
    r = {}
    for c, p in terms: r = padd(r, pscale(p, c))
    return r

U  = padd(ONE, pmul(X, Y))                       # 1+xy
U2 = pmul(U, U); U3 = pmul(U2, U)
W  = padd(pscale(ONE, 4), pscale(pmul(X, Y), 3)) # 4+3xy
F1 = padd(pmul(U3, Z), pmul(pmul(Y, Y), pmul(U, W)))
F2 = padd(Y, padd(pscale(pmul(X, pmul(U2, Z)), 3), pscale(pmul(pmul(X, pmul(Y, Y)), W), 3)))
F3 = lin((2, X), (-3, pmul(pmul(X, X), Y)), (-1, pmul(pmul(pmul(X, X), X), Z)))
F = [F1, F2, F3]

J = [[pdiff(F[i], k) for k in range(3)] for i in range(3)]

def det3(M):
    d = {}
    for (a,b,c,s) in [(0,1,2,1),(1,2,0,1),(2,0,1,1),(2,1,0,-1),(0,2,1,-1),(1,0,2,-1)]:
        d = padd(d, pscale(pmul(M[0][a], pmul(M[1][b], M[2][c])), s))
    return d

D = det3(J)
assert D == {(0,0,0): Fr(-2)}, D
print("check 1: det JF == -2 identically  OK")

pts = [(Fr(0),Fr(0),Fr(-1,4)), (Fr(1),Fr(-3,2),Fr(13,2)), (Fr(-1),Fr(3,2),Fr(13,2))]
for pt in pts:
    v = tuple(peval(Fi, pt) for Fi in F)
    assert v == (Fr(-1,4), Fr(0), Fr(0)), (pt, v)
print("check 2: three-point collision at (-1/4,0,0)  OK")

# adjugate via cofactors, A = adj(J) * (-1/2)
def minor(M, i, k):
    rows = [r for r in range(3) if r != i]; cols = [c for c in range(3) if c != k]
    a, b = rows; c, d = cols
    return padd(pmul(M[a][c], M[b][d]), pscale(pmul(M[a][d], M[b][c]), -1))
ADJ = [[pscale(minor(J, k, i), (-1)**(i+k)) for k in range(3)] for i in range(3)]  # adj = C^T
A = [[pscale(ADJ[i][k], Fr(-1,2)) for k in range(3)] for i in range(3)]

for i in range(3):
    for j in range(3):
        s = {}
        for k in range(3): s = padd(s, pmul(J[i][k], A[k][j]))
        assert s == ({(0,0,0): Fr(1)} if i == j else {}), (i, j)
        s2 = {}
        for k in range(3): s2 = padd(s2, pmul(A[i][k], J[k][j]))
        assert s2 == ({(0,0,0): Fr(1)} if i == j else {}), (i, j)
print("check 3: J A == A J == I with polynomial A  OK")

# commutator coefficients b_l^{(ij)} = sum_k ( A[k][i] d_k A[l][j] - A[k][j] d_k A[l][i] )
for i in range(3):
    for j in range(i+1, 3):
        for l in range(3):
            b = {}
            for k in range(3):
                b = padd(b, pmul(A[k][i], pdiff(A[l][j], k)))
                b = padd(b, pscale(pmul(A[k][j], pdiff(A[l][i], k)), -1))
            assert b == {}, (i, j, l, b)
print("check 4: all nine commuting-vector-field coefficients vanish  OK")
print("check 5: the same identities verify the Weyl and Poisson relations  OK")

print("\nAll Weyl and Poisson lift identities independently verified with exact arithmetic (no CAS).")