verify-phase-free-nahm.py

#!/usr/bin/env python3
"""Symbolic degree identities and a certified odd-numerator control.

Standard library only. No slope census, floating-point arithmetic, or
saturation of the native edge/filling image is used.
"""
from fractions import Fraction as F
from functools import lru_cache
from math import isqrt
from pathlib import Path
import json

ROOT = Path(__file__).resolve().parent
checks, tails = [], []
def check(name, value):
    assert value, name
    checks.append(name)

# Exact Laurent polynomials in p,q,t,a1,a2. The Bezout relation is imposed
# by r=(1+p*t)/q; q is nonzero throughout the theorem.
NV = 5
def clean(p): return {m:c for m,c in p.items() if c}
def pc(c): return clean({(0,)*NV:F(c)})
def var(j, power=1):
    v=[0]*NV; v[j]=power
    return {tuple(v):F(1)}
def pa(*ps):
    out={}
    for p in ps:
        for m,c in p.items():out[m]=out.get(m,0)+c
    return clean(out)
def pm(p,q):
    out={}
    for m,c in p.items():
        for n,d in q.items():
            k=tuple(a+b for a,b in zip(m,n))
            out[k]=out.get(k,0)+c*d
    return clean(out)
def scale(c,p):return pm(pc(c),p)
def neg(p):return scale(-1,p)
def sub(p,q):return pa(p,neg(q))
def pdot(x,y):return pa(*(pm(a,b) for a,b in zip(x,y)))
p,q,t,a1,a2=[var(j) for j in range(NV)]
iq=var(1,-1)
r=pm(pa(pc(1),pm(p,t)),iq)
check('symbolic Bezout relation qr-pt=1',sub(sub(pm(q,r),pm(p,t)),pc(1))=={})
diag=sub(pc(F(1,2)),scale(F(1,8),pm(p,iq)))
off=scale(F(1,8),pm(p,iq))
def T(x):return (pa(pm(diag,x[0]),pm(off,x[1])),pa(pm(off,x[0]),pm(diag,x[1])))
def Q(x):return sub(pdot(x,x),pdot(x,T(x)))
v=(sub(pc(1),iq),pa(pc(1),iq))
lin0=(pa(pc(1),iq),sub(pc(1),iq))
alpha=(sub(r,scale(4,t)),r)
beta=(sub(r,scale(2,t)),pc(0))
rho=pa(neg(r),scale(2,t))
delta=tuple(sub(b,c) for b,c in zip(beta,T(alpha)))
check('longitude fractional defect first component',delta[0]==scale(F(1,2),iq))
check('longitude fractional defect second component',delta[1]==scale(F(-1,2),iq))
kappa=pa(scale(2,pdot(alpha,delta)),pdot(v,alpha),scale(2,rho))
check('native scalar cancels the complete square constant',kappa=={})
check('native scalar cannot be omitted',pa(scale(2,pdot(alpha,delta)),pdot(v,alpha))==scale(-2,rho))
avec=(a1,a2)
x=tuple(sub(a,b) for a,b in zip(avec,alpha))
raw=pa(scale(2,Q(x)),scale(2,pa(*x)),neg(pdot(v,x)),
       scale(4,pdot(tuple(sub(a,b) for a,b in zip(alpha,beta)),x)),
       scale(2,pa(pdot(alpha,alpha),pa(*alpha),neg(pdot(alpha,beta)))),scale(2,rho))
expected=pa(scale(2,Q(avec)),pm(sub(pc(1),iq),a1),pm(pa(pc(1),iq),a2))
check('all-parameter native longitude-square degree identity',raw==expected)
unit=pa(scale(2,Q(avec)),scale(2,pa(*avec)),neg(pdot(v,avec)))
check('all-parameter unit degree identity',unit==pa(scale(2,Q(avec)),pdot(lin0,avec)))
check('copositivity decomposition at the boundary',
      Q(avec)==pa(pm(a1,a2),scale(F(1,8),pm(pa(pc(4),pm(p,iq)),pm(sub(a1,a2),sub(a1,a2))))))

def dot(x,y):return sum(a*b for a,b in zip(x,y))
def addv(x,y):return tuple(a+b for a,b in zip(x,y))
def sv(c,x):return tuple(c*a for a in x)
def omega(x,y):return dot(x[2:],y[:2])-dot(x[:2],y[2:])
def collapse(R):return tuple(r-s for r,s,t in R)+tuple(t-s for r,s,t in R)
E1=((2,0,1),(2,0,1)); E2=((0,2,1),(0,2,1))
MU=((0,-1,0),(1,0,0)); LONG=((-2,2,0),(0,0,0))
h=collapse(E1); mu=collapse(MU); longitude=collapse(LONG)
check('retained actual raw rows',h==(2,2,1,1) and collapse(E2)==sv(-1,h)
      and mu==(1,1,1,0) and longitude==(-4,0,-2,0))
check('retained actual scalar exponents',[sum(x[1] for x in R) for R in (E1,E2,MU,LONG)]==[0,4,-1,2])

# A single deliberately excluded-by-the-old-phase-rule control: p=q=1,
# r=1,t=0, h_s=(-3,1;-1,0), c_s=iA, l=-iA^-1 W(mu).
# Gaussian integer coefficients are pairs of Python integers, never floats.
def ga(a,b):return (a[0]+b[0],a[1]+b[1])
def gm(a,b):return (a[0]*b[0]-a[1]*b[1],a[0]*b[1]+a[1]*b[0])
def ipow(k):return ((1,0),(0,1),(-1,0),(0,-1))[k%4]
def sclean(s):return {e:c for e,c in s.items() if c!=(0,0)}
def sa(*ss):
    out={}
    for s in ss:
        for e,c in s.items():out[e]=ga(out.get(e,(0,0)),c)
    return sclean(out)
def shift(s,p,phase=(1,0),cap=None):
    return sclean({e+p:gm(phase,c) for e,c in s.items() if cap is None or e+p<=cap})
@lru_cache(None)
def invpoch(indices,cap):
    if cap<0:return {}
    v=[0]*(cap+1);v[0]=1
    for n in indices:
        for j in range(1,n+1):
            for e in range(4*j,cap+1):v[e]+=v[e-4*j]
    return {e:(c,0) for e,c in enumerate(v) if c}
def T1(x):return (F(3,8)*x[0]+F(1,8)*x[1],F(1,8)*x[0]+F(3,8)*x[1])
def Q1(x):return dot(x,x)-dot(x,T1(x))
def H(u,cap):
    alpha=tuple(F(i,2) for i in u[:2]); beta=tuple(F(i,2) for i in u[2:])
    delta=tuple(a-b for a,b in zip(beta,T1(alpha)))
    L=(2-4*delta[0],-4*delta[1])
    const=2*dot(alpha,delta)+2*alpha[1]
    # For p=q=1, 2Q(a)>=||a||^2 on all of R^2. Hence
    # D >= ||a||^2/2 - ||L||^2/2 + const.
    bound=max(F(0),2*(cap-const)+dot(L,L))
    R=isqrt(bound.numerator//bound.denominator)+1
    tails.append({'input':u,'cap_A':cap,'radius_squared_bound':str(bound),'square_cutoff':R})
    out={}
    for a1 in range(R+1):
        for a2 in range(R+1):
            a=(a1,a2); x=tuple(a-b for a,b in zip(a,alpha))
            n=(x[1]-x[0])/2; e=x[1]-n/2
            if n.denominator!=1 or e.denominator!=1:continue
            n,e=int(n),int(e)
            lam=addv(sv(e,h),sv(n,(-3,1,-1,0)))
            check_graph=tuple(F(i,2) for i in lam[2:])==T1(x)
            assert check_graph and tuple(F(i,2) for i in lam[:2])==x
            ds=tuple(a+b for a,b in zip(T1(x),beta))
            if any(v.denominator!=1 for v in ds):continue
            # chi(eh+nh_s)=(-1)^e i^n A^(2e+n).
            deg=2*(dot(a,a)+sum(a)-dot(a,ds))+F(omega(lam,u),2)-(2*e+n)
            expanded=2*Q1(a)+dot(L,a)+const
            assert deg==expanded and deg.denominator==1
            deg=int(deg)
            if deg>cap:continue
            phase=ipow(2*sum(a)-2*e-n)
            out=sa(out,shift(invpoch(a,cap-deg),deg,phase,cap))
    return out

CAP=48
hs=(-3,1,-1,0); b=mu
unit_series=H((0,0,0,0),CAP)
long_series=shift(H(sv(2,b),CAP+2),-2,(-1,0),CAP)
vacuum=sa(unit_series,shift(long_series,0,(-1,0)))
check('odd-numerator unit has exactly one constant',unit_series.get(0)==(1,0) and min(unit_series)==0)
check('odd-numerator longitude-square has strictly positive order',long_series and min(long_series)>0)
check('odd-numerator filled vacuum constant one',vacuum.get(0)==(1,0) and min(vacuum)==0)
check('odd-numerator unit has a negative native summand',
      addv(sv(4,h),sv(2,hs))==(2,10,2,4)
      and ipow(2*(1+5)-2*4-2)==(-1,0))
check('odd-numerator commutator retains orientation',omega(hs,b)==2)

for u in [(0,0,0,0),hs,(0,0,1,1),(2,0,-2,2)]:
    for j in range(2):
        uz=list(u);uz[j]-=2
        uw=list(u);uw[j+2]+=2
        first=shift(H(tuple(uz),CAP+u[j+2]),-u[j+2],cap=CAP)
        second=shift(H(tuple(uw),CAP+u[j]),-u[j],cap=CAP)
        check(f'left tetrahedron multiplier, j={j}, input={u}',
              sa(first,second,shift(H(u,CAP),0,(-1,0)))=={})
    for row,power0,phase,name in [(h,0,(1,0),'edge1'),(sv(-1,h),4,(1,0),'edge2'),(hs,1,(0,1),'filling')]:
        power=F(omega(row,u),2)+power0
        assert power.denominator==1
        power=int(power)
        check(f'right multiplier, {name}, input={u}',
              sa(shift(H(addv(row,u),CAP-power),power,phase,CAP),shift(H(u,CAP-2),2,cap=CAP))=={})
check('odd-root supported input has the native imaginary phase',
      H(hs,CAP)==shift(H((0,0,0,0),CAP-1),1,(0,1),CAP))
check('odd-root supported input is nonzero',H(hs,CAP)!={})
check('fractional endpoint defect is not integer at q=1',F(1,2).denominator!=1)

result={'status':'passed','checks':checks,'check_count':len(checks),
        'theorem_family':'gcd(p,q)=1, q!=0, p/q>-4; no numerator congruence restriction',
        'symbolic_scope':'Laurent polynomial identities for arbitrary p,q,t with r=(1+pt)/q',
        'control_slope':[1,1],
        'unit_coefficients_A':{str(e):list(c) for e,c in sorted(unit_series.items())},
        'longitude_square_coefficients_A':{str(e):list(c) for e,c in sorted(long_series.items())},
        'vacuum_coefficients_A':{str(e):list(c) for e,c in sorted(vacuum.items())},
        'coefficient_encoding':'[real,imaginary] Gaussian integers',
        'quadratic_tail_certificates':tails,
        'limits':'One fixed actual triangulation; dependencies and proof in note; not an unrestricted conjecture resolution.'}
(ROOT/'exact-results.json').write_text(json.dumps(result,indent=2)+'\n')
print(f'PASS: {len(checks)} exact checks, with certified infinite-series tails.')
print('Unit:',sorted(unit_series.items()))
print('Longitude square:',sorted(long_series.items()))
print('Filled vacuum:',sorted(vacuum.items()))