verify-normal-sector.py

"""Exact native-index controls. Python standard library only.

Run: python3 verify-normal-sector.py
Series use x=A^(1/2), q=x^8, Gaussian integer coefficients (pairs).
The finite tail cut is certified from the exact piecewise-quadratic degree.
This is not a formalization of the general normal-surface theorem.
"""
from functools import lru_cache
import json

def gadd(a,b): return (a[0]+b[0],a[1]+b[1])
def gmul(a,b): return (a[0]*b[0]-a[1]*b[1],a[0]*b[1]+a[1]*b[0])
def ipow(n): return ((1,0),(0,1),(-1,0),(0,-1))[n%4]
def add(a,b):
    c=dict(a)
    for k,v in b.items(): c[k]=gadd(c.get(k,(0,0)),v)
    return {k:v for k,v in c.items() if v!=(0,0)}
def scale(a,n=0,g=(1,0),cut=64):
    return {k+n:gmul(v,g) for k,v in a.items() if k+n<=cut}
def mul(a,b,cut):
    c={}
    for k,v in a.items():
        for j,w in b.items():
            if k+j<=cut:c[k+j]=gadd(c.get(k+j,(0,0)),gmul(v,w))
    return {k:v for k,v in c.items() if v!=(0,0)}
@lru_cache(None)
def invpoch(n,cut):
    if cut<0:return {}
    c=[1]+[0]*cut
    for j in range(1,n+1):
        for k in range(8*j,cut+1): c[k]+=c[k-8*j]
    return {k:(v,0) for k,v in enumerate(c) if v}

def seed_degree(m,e):
    y=(e,0,-m); low=min(y); z=[a-low for a in y]
    return 4*(z[0]*z[1]+z[1]*z[2]+z[2]*z[0]-low)

@lru_cache(None)
def seed(m,e,cut):
    # Reduce J(e,0,-m), permute to J(p,r,0)=I(-p,r).
    y=(e,0,-m); low=min(y); z=sorted(a-low for a in y)
    p,r=z[1:]; base=4*(p*r-low)
    if base>cut:return {}
    out={}; j=0
    while base+4*j*(j+1)+8*p*j<=cut:
        power=base+4*j*(j+1)+8*p*j
        coeff=mul(invpoch(j,cut-power),invpoch(j+r,cut-power),cut-power)
        out=add(out,scale(coeff,power,ipow(2*(j-low)),cut))
        j+=1
    return out

def omega(a,b): return sum(a[j+1]*b[j]-a[j]*b[j+1] for j in range(0,len(a),2))

DATA={
 'solid': {'rows':[(2,1,2,2,2,2),(0,1,0,0,0,0)],
           'h':(-1,-1,0,0),'chi_phase':1,'chi_A':1,
           'u':(0,2,0,0),'blocks':2},
 'product':{'rows':[(1,2,2,1,1,1,1,1,1),(1,0,0,0,0,0,0,0,0),(0,0,0,1,1,1,1,1,1)],
           'h':(1,0,0,0,0,0),'chi_phase':2,'chi_A':2,
           'u':(0,-2,0,0,0,0),'blocks':3}}

def P(row):return tuple(v for j in range(0,len(row),3) for v in (row[j]-row[j+1],row[j+2]-row[j+1]))

def degree(data,u,n):
    h=data['h']; v=[u[j]+n*h[j] for j in range(len(h))]
    if any(a%2 for a in v):return None
    return (-2*n*data['chi_A']+n*omega(h,u)
            +sum(seed_degree(v[j+1]//2,-v[j]//2) for j in range(0,len(v),2)))

def average(data,u,cut=64):
    h=data['h']; bound=2*(cut+sum(abs(v) for v in u)+12)
    # Outside this bound all local minimum breakpoints have passed.
    # On each supported parity tail the degree is a polynomial of degree <=2.
    for parity in (0,1):
        for sign in (-1,1):
            n=sign*(bound+parity)
            ds=[degree(data,u,n+sign*2*j) for j in range(3)]
            if ds[0] is None:
                assert all(v is None for v in ds); continue
            assert ds[0]>cut and ds[1]>ds[0] and ds[2]-2*ds[1]+ds[0]>=0
    out={}
    for n in range(-bound,bound+1):
        deg=degree(data,u,n)
        if deg is None or deg>cut:continue
        v=[u[j]+n*h[j] for j in range(len(h))]
        mins=[seed_degree(v[j+1]//2,-v[j]//2) for j in range(0,len(v),2)]
        weight=-2*n*data['chi_A']+n*omega(h,u)
        # Track each factor's minimum so intermediate truncations cannot lose
        # terms that a later negative degree would lower into the window.
        term={0:(1,0)}; remain=sum(mins)
        for b,j in enumerate(range(0,len(v),2)):
            remain-=mins[b]
            limit=cut-weight-remain
            term=mul(term,seed(v[j+1]//2,-v[j]//2,cut-deg+mins[b]),limit)
        out=add(out,scale(term,weight,ipow(-n*data['chi_phase']),cut))
    return out

def rational_first(a,b,cut):
    m=b-a
    if m<=0:return {}
    exponent=4*b*(m-1)
    return scale(invpoch(m-1,cut-exponent),exponent,ipow(2*b),cut)

checks=[]; results={}
for name,data in DATA.items():
    h=data['h']
    for row in data['rows']:
        v=P(row)
        n=next((v[j]//h[j] for j in range(len(h)) if h[j]),0)
        assert v==tuple(n*x for x in h)
        p=sum(row[1::3])
        assert (2-p-n*data['chi_phase'])%4==0
        assert 2-p==n*data['chi_A']
        checks.append('native row character '+name)
    unit=average(data,(0,)*len(h)); detected=average(data,data['u'])
    assert not unit
    expected=seed(0,0,64)
    if name=='solid':expected=scale(expected,g=(-1,0))
    else:expected=mul(expected,expected,64)
    assert detected==expected
    results[name]={'unit':unit,'detector_input':data['u'],
                   'detector_q_coefficients':[detected.get(8*j,(0,0))[0] for j in range(9)]}
    # Original primitive edge relation on even and odd root inputs.
    for u in ((0,)*len(h),data['u'],tuple(-x for x in h),
              tuple(data['u'][j]-h[j] for j in range(len(h)))):
        shift=omega(h,u)
        v=tuple(u[j]+h[j] for j in range(len(h)))
        left=scale(average(data,v,80),shift,cut=64)
        right=scale(average(data,u,80),2*data['chi_A'],ipow(data['chi_phase']),64)
        assert left==right
        checks.append('primitive edge '+name)
    # Correct-sided tetrahedron relation W(u) f_1.
    for u in ((0,)*len(h),data['u'],tuple(-x for x in h)):
        v=list(u);v[0]-=2
        w=list(u);w[1]+=2
        lhs=add(scale(average(data,tuple(v),80),-2*u[1],cut=64),
                scale(average(data,tuple(w),80),-2*u[0],cut=64))
        assert lhs==average(data,u)
        checks.append('tetrahedron '+name)

for a in range(-2,3):
    for b in range(-2,4):
        actual=average(DATA['solid'],(2*a,2*b,0,0),48)
        expected=mul(rational_first(a,b,80),seed(0,0,80),48)
        assert actual==expected,(a,b,actual,expected)
        checks.append('solid closed evaluation')

# Direct integer unit bridge and native sign for a finite range of actual rows.
for data in DATA.values():
    row=data['rows'][1]
    for n in range(-8,9):
        R=tuple(n*v for v in row); lam=P(R)
        if any(x%2 for x in lam):continue
        Y=tuple(v for j in range(0,len(lam),2) for v in (-lam[j]//2,0,-lam[j+1]//2))
        low=[min(Y[j:j+3]) for j in range(0,len(Y),3)]
        p=sum(R[1::3]); assert p%2==0
        chi=n-p//2+sum(low)
        assert (2*n+p-2*sum(low)+2*chi)%4==0
        assert -2*n+p-2*sum(low)==-2*chi
        checks.append('unit phase bridge')

print(json.dumps({'status':'passed','checks':len(checks),'x_cutoff':64,
 'results':results,'scope':'two cited actual incidence data; exact series controls with tail bounds; no formal or independent topological verification'},indent=2))