#!/usr/bin/env python3 """Exact, standard-library-only certificates for the Gaussian norm unit. Run with Python 3.9+. Results are printed as JSON; no input/output files are changed. All numerical checks use integers or Fraction; drawings do not certify equality. Verification gates remain active in normal and optimized Python. """ from fractions import Fraction from math import gcd,isqrt import json def verify(condition,message): """Keep mathematical verification active under python -O.""" if not condition:raise AssertionError(message) def add(z,w):return (z[0]+w[0],z[1]+w[1]) def sub(z,w):return (z[0]-w[0],z[1]-w[1]) def mul(z,w):return (z[0]*w[0]-z[1]*w[1],z[0]*w[1]+z[1]*w[0]) def norm(z):return z[0]*z[0]+z[1]*z[1] def conj(z):return z[0],-z[1] def power(z,k): out=(1,0) for _ in range(k):out=mul(out,z) return out def divrem(z,w): if w==(0,0):raise ZeroDivisionError('zero Gaussian divisor') u=mul(z,conj(w));d=norm(w) def nearest(v):return (2*v+d)//(2*d) # ties toward +infinity, including negatives q=(nearest(u[0]),nearest(u[1]));r=sub(z,mul(q,w)) verify(norm(r)*2<=norm(w), 'divrem certificate or invariant failed: norm(r)*2<=norm(w)') return q,r def exact_div(z,w): q,r=divrem(z,w) if r!=(0,0):raise ValueError('not divisible in Z[i]') return q def xgcd(z,w): if z==w==(0,0):raise ValueError('both Gaussian inputs zero') a,b=z,w;u,v=(1,0),(0,0);s,t=(0,0),(1,0);trace=[] while b!=(0,0): q,r=divrem(a,b);trace.append({'a':a,'b':b,'q':q,'r':r}) a,b=b,r;u,v=v,sub(u,mul(q,v));s,t=t,sub(s,mul(q,t)) verify(add(mul(u,z),mul(s,w))==a, 'xgcd certificate or invariant failed: add(mul(u,z),mul(s,w))==a') exact_div(z,a);exact_div(w,a) return a,u,s,trace def factor(n): if n<1:raise ValueError('positive factorization input required') out={};p=2 while p*p<=n: while n%p==0:out[p]=out.get(p,0)+1;n//=p p=3 if p==2 else p+2 if n>1:out[n]=out.get(n,0)+1 return out def r2_formula(n): f=factor(n) if any(p%4==3 and e%2 for p,e in f.items()):return 0 v=4 for p,e in f.items(): if p%4==1:v*=e+1 return v def primitive_formula(n): f=factor(n) if f.get(2,0)>=2 or any(p%4==3 for p in f):return 0 return 4*2**sum(p%4==1 for p in f) def all_representations(n): if n<0:return set() out=set() for x in range(-isqrt(n),isqrt(n)+1): q=n-x*x;y=isqrt(q) if y*y==q:out.add((x,y));out.add((x,-y)) return out def primitive_representations(n):return {z for z in all_representations(n) if gcd(*z)==1} def positive_unordered(reps):return {tuple(sorted(map(abs,z))) for z in reps if 0 not in z} def chi4(n):return 0 if n%2==0 else (1 if n%4==1 else -1) def divisors(n): out=[] for d in range(1,isqrt(n)+1): if n%d==0: out.append(d) if d*d!=n:out.append(n//d) return out def cornacchia(m,t): if m<=1:raise ValueError('m must exceed 1') if (t*t+1)%m:raise ValueError('invalid square root of -1') t%=m;t=min(t,m-t) a,b,u,v=m,t,0,1;j=1 trace=[{'j':0,'r':m,'s':0},{'j':1,'r':t,'s':1}] verify(gcd(m,t)==1, 'cornacchia certificate or invariant failed: gcd(m,t)==1') while True: verify(a*v+b*u==m, 'cornacchia certificate or invariant failed: a*v+b*u==m') verify((b-(-1)**(j-1)*t*v)%m==0, 'cornacchia certificate or invariant failed: (b-(-1)**(j-1)*t*v)%m==0') if b*b0, 'cornacchia certificate or invariant failed: b>0') verify(a*a>=m and 0=m and 0n>0 and gcd(m,n)==1 and (m-n)%2==1):raise ValueError('invalid primitive parameters') return m*m-n*n,2*m*n,m*m+n*n def recover(a,b,c): if not(a>0 and b>0 and c>0 and a*a+b*b==c*c and gcd(gcd(a,b),c)==1):raise ValueError('not primitive triple') if b%2:a,b=b,a m2,n2=(c+a)//2,(c-a)//2;m,n=isqrt(m2),isqrt(n2) verify(m*m==m2 and n*n==n2 and triple(m,n)==(a,b,c), 'recover certificate or invariant failed: m*m==m2 and n*n==n2 and triple(m,n)==(a,b,c)') return m,n def main(): checks={};z=(5,2);w=(5,-2) verify(mul(z,w)==(29,0), 'main certificate or invariant failed: mul(z,w)==(29,0)') verify(add(mul(z,(-2,2)),mul((3,0),w))==(1,0), 'main certificate or invariant failed: add(mul(z,(-2,2)),mul((3,0),w))==(1,0)') verify(mul((16,0),(3,1))==add((1,0),mul(w,(7,6))), 'main certificate or invariant failed: mul((16,0),(3,1))==add((1,0),mul(w,(7,6)))') delta,u,v,tr=xgcd((29,0),(12,1)) verify(delta==(5,-2) and u==(1,0) and v==(-2,0), 'main certificate or invariant failed: delta==(5,-2) and u==(1,0) and v==(-2,0)') checks['gaussian_gcd_29_12_plus_i']={'gcd':delta,'bezout':[u,v],'trace':tr} grid=[(a,b) for a in range(-3,4) for b in range(-3,4)];gc=0 for a in grid: for b in grid: if a==b==(0,0):continue d,u,v,_=xgcd(a,b);gc+=1 verify(add(mul(u,a),mul(v,b))==d, 'main certificate or invariant failed: add(mul(u,a),mul(v,b))==d') verify(norm(mul(a,b))==norm(a)*norm(b), 'main certificate or invariant failed: norm(mul(a,b))==norm(a)*norm(b)') checks['gaussian_grid_gcd_cases']=gc # Kernel equivalence is tested independently by coordinate congruence and exact division. kernel=0 for a in range(-30,31): for b in range(-30,31): in_kernel=(a-12*b)%29==0 divisible=(5*a-2*b)%29==0 and (2*a+5*b)%29==0 verify(in_kernel==divisible, 'main certificate or invariant failed: in_kernel==divisible');kernel+=1 checks['quotient_kernel_cases']=kernel count_n=3000 for n in range(1,count_n+1): reps=all_representations(n);prim={z for z in reps if gcd(*z)==1} verify(len(reps)==r2_formula(n), 'main certificate or invariant failed: len(reps)==r2_formula(n)') verify(len(prim)==primitive_formula(n), 'main certificate or invariant failed: len(prim)==primitive_formula(n)') verify(len(reps)==4*sum(chi4(d) for d in divisors(n)), 'main certificate or invariant failed: len(reps)==4*sum(chi4(d) for d in divisors(n))') A=int(isqrt(n)**2==n);D=int(n%2==0 and isqrt(n//2)**2==n//2) verify(len(positive_unordered(reps))*8==len(reps)-4*A+4*D, 'main certificate or invariant failed: len(positive_unordered(reps))*8==len(reps)-4*A+4*D') verify(len({tuple(sorted(map(abs,z))) for z in reps})*8==len(reps)+4*A+4*D, 'main certificate or invariant failed: len({tuple(sorted(map(abs,z))) for z in reps})*8==len(reps)+4*A+4*D') checks['count_formula_n_range']=[1,count_n] r325=all_representations(325);p325=primitive_representations(325) verify(positive_unordered(r325)=={(1,18),(6,17),(10,15)}, 'main certificate or invariant failed: positive_unordered(r325)=={(1,18),(6,17),(10,15)}') verify(positive_unordered(p325)=={(1,18),(6,17)}, 'main certificate or invariant failed: positive_unordered(p325)=={(1,18),(6,17)}') allocations=set() for h in range(3): for k in range(2): z=mul(mul(power((2,1),h),power((2,-1),2-h)),mul(power((3,2),k),power((3,-2),1-k))) for u in [(1,0),(-1,0),(0,1),(0,-1)]:allocations.add(mul(u,z)) verify(allocations==r325, 'main certificate or invariant failed: allocations==r325') checks['325']={'r2':len(r325),'primitive_r2':len(p325),'positive_unordered':sorted(positive_unordered(r325)),'primitive_pairs':sorted(positive_unordered(p325)),'all_signed_ordered':sorted(r325)} corn_cases=0 for m in range(2,2001): produced=set() for t in range(1,m//2+1): if (t*t+1)%m:continue z,_=cornacchia(m,t);produced.add(tuple(sorted(z)));corn_cases+=1 target=positive_unordered(primitive_representations(m)) verify(produced==target, (m,produced,target)) checks['cornacchia_all_root_cases']=corn_cases checks['cornacchia_moduli_range']=[2,2000] checks['cornacchia_325']={str(t):{'answer':cornacchia(325,t)[0],'trace':cornacchia(325,t)[1]} for t in [18,57]} checks['roots_minus_one_mod_325']=[t for t in range(325) if (t*t+1)%325==0] verify(checks['roots_minus_one_mod_325']==[18,57,268,307], "main certificate or invariant failed: checks['roots_minus_one_mod_325']==[18,57,268,307]") verify(pow(2,7,29)==12 and cornacchia(29,12)[0]==(5,2), 'main certificate or invariant failed: pow(2,7,29)==12 and cornacchia(29,12)[0]==(5,2)') # Independent brute-force triples vs parameter-generated triples, with odd leg first. triples=set() for c in range(1,501): for a in range(1,c): b=isqrt(c*c-a*a) if b*b!=c*c-a*a or gcd(a,b)!=1:continue odd,even=(a,b) if a%2 else (b,a) triples.add((odd,even,c));recover(odd,even,c) generated=set() for m in range(2,isqrt(500)+1): for n in range(1,m): if m*m+n*n<=500 and gcd(m,n)==1 and (m-n)%2:generated.add(triple(m,n)) verify(generated==triples, 'main certificate or invariant failed: generated==triples') fixed=sorted(t for t in triples if t[2]==325) verify(fixed==[(253,204,325),(323,36,325)], 'main certificate or invariant failed: fixed==[(253,204,325),(323,36,325)]') checks['primitive_triples_c_le_500']=len(triples) checks['primitive_triples_c_325']=fixed rational=0 for num in range(-20,21): for den in range(1,21): t=Fraction(num,den);X=(1-t*t)/(1+t*t);Y=2*t/(1+t*t) verify(X*X+Y*Y==1 and X!=-1 and Y/(1+X)==t, 'main certificate or invariant failed: X*X+Y*Y==1 and X!=-1 and Y/(1+X)==t');rational+=1 checks['rational_circle_exact_cases']=rational # Mutated input/certificates must actually be rejected, not silently tolerated. rejected=0 for f in [lambda:cornacchia(325,56),lambda:cornacchia(1,0),lambda:xgcd((0,0),(0,0)),lambda:exact_div((5,2),(5,-2)),lambda:triple(3,1),lambda:triple(4,2),lambda:recover(9,12,15)]: try:f() except (ValueError,ZeroDivisionError):rejected+=1 else:raise AssertionError('invalid input unexpectedly accepted') verify(40*40+5*5!=325 and 17*17+5*5!=325, 'main certificate or invariant failed: 40*40+5*5!=325 and 17*17+5*5!=325') verify(1+5*2*2==21 and not any(x*x+5*y*y==7 for x in range(3) for y in range(2)), 'main certificate or invariant failed: 1+5*2*2==21 and not any(x*x+5*y*y==7 for x in range(3) for y in range(2))') checks['rejected_invalid_inputs']=rejected checks['mutated_stop_and_coefficient_rejected']=2 checks['general_d_counterexample']={'d':5,'m':7,'root':3,'terminal_r':1,'terminal_s':2,'weighted_norm':21} return {'status':'PASS','arithmetic':'integer and Fraction only','checks':checks} if __name__=='__main__':print(json.dumps(main(),ensure_ascii=False,indent=2))