back to table · edit · history · where entries came from · files · download
9358 bytes, as of the version from 2026-11-10 00:44 (current). Recorded here, not run.
"""Independent checks for Bloch elliptic dilogarithm -- numberdb.org/T700.
sage -python checks.py # pre-submission numerical evidence
sage -python checks.py --stored # repeat on values read from the draft
Needs generate.py beside it, SageMath, numberdb and mpmath. The independent
Fourier calculation uses mpmath, no polylogarithm, and sums the geometric
q-series analytically before summing the Fourier index. A finite Hurwitz-zeta
sum supplies the unit-circle term. This is an agreement check, not the source
of the table's rigorous bounds. The generator's Arb balls certify the digits.
The Mahler check uses a different integral, with an exact nonzero integration
control. The lattice sum only settles the sign at the square-lattice 4-torsion
point; its polynomial tail does not verify a hundred digits. All-entry Fourier
checks and symmetry checks provide the full-range tests.
"""
import os
import sys
import json
import math
import hashlib
from pathlib import Path
import numberdb.sage as numberdb
from sage.rings.real_arb import RealBallField
from sage.rings.complex_arb import ComplexBallField
from sage.rings.rational_field import QQ
from numberdb._write import to_text
from numberdb._generate import _stored_entries
import mpmath as mp
from generate import Table, DISCRIMINANTS, elliptic_d
MANIFEST = 'agents/runs/20261110T001313Z-build-skills/manifest.json'
SKILL_SHA = '886390324dd1971e381a43cf3b46b822ab0bc67d62f31d8481e29ec79556746b'
def fourier(delta, N, r, s, dps=150):
with mp.workdps(dps):
L=mp.pi*mp.sqrt(-delta)
a=mp.mpf(r)/N
theta=2*mp.pi*(mp.mpf(s)-mp.mpf((delta%2)*r)/2)/N
rho=mp.exp(-L)
q=((-1)**(delta%2))*rho
x=mp.exp(-L*a+1j*theta)
if r == 0:
# A complete residue-class grouping of the absolutely convergent
# Clausen series, not a library dilogarithm evaluation.
total=mp.fsum(mp.sin(2*mp.pi*s*j/N)*mp.zeta(2,mp.mpf(j)/N)
for j in range(1,N+1))/(N*N)
rate=L
else:
total=mp.mpf(0)
rate=L*min(a,1-a)
K=int(mp.ceil((dps+15)*mp.log(10)/rate))+3
X=mp.mpc(1);Q=mp.mpf(1)
for k in range(1,K+1):
X*=x;Q*=q
if r==0:
term=(X-1/X)*(Q/(k*k*(1-Q))+L*Q/(k*(1-Q)**2))
else:
term=(X*((1/(mp.mpf(k)**2)+L*a/k)/(1-Q)+L*Q/(k*(1-Q)**2))
-(1/X)*((1/(mp.mpf(k)**2)-L*a/k)*Q/(1-Q)+L*Q/(k*(1-Q)**2)))
total+=mp.im(term)
return +total
def as_mp_ball(value):
lo,hi=value.endpoints()
# Sage's endpoint conversion is directed; retain all displayed digits.
return mp.mpf(str(lo)),mp.mpf(str(hi))
def parse_written(text):
text=str(text).rstrip('?')
from decimal import Decimal
d=Decimal(text)
# A plain stored decimal means plus or minus one last-place unit.
u=mp.power(10,d.as_tuple().exponent)
m=mp.mpf(text)
return m-u,m+u,u
def mahler_check():
C=ComplexBallField(460);R=RealBallField(460)
control=C.integral(lambda z,analytic:z*z,0,1)
assert control.real().is_finite() and control.imag().is_finite() and (control.real()-R(QQ(1)/3)).contains_zero()
assert control.real().accuracy()>400 and control.imag().contains_zero()
q=(-2*R.pi()).exp()
def theta_phi(sign):
v=R(1)+2*sum((sign*q)**(n*n) for n in range(1,16))
tail=2*q**256/(1-q**33)
return v.add_error(tail)
k=4*(theta_phi(1)/theta_phi(-1))**2
def integrand(z,analytic):
w=(C(k)+2*z.cos())/2
return (w+(w-1).sqrt(analytic=analytic)*(w+1).sqrt(analytic=analytic)).log(analytic=analytic)
integral=C.integral(integrand,0,2*C.pi(),abs_tol=R(10)**(-120),rel_tol=R(10)**(-120))
m=(integral/(2*C.pi())).real()
target=4*elliptic_d(-4,4,0,1,110)/R.pi()
assert m.is_finite() and m.accuracy()>360 and integral.imag().contains_zero()
assert (m-target).contains_zero(), (m,target)
# Float lattice sum, with an explicit analytic remainder bound for tau=i.
# Outside max(|m|,|n|)<=B, each shell has 8k points and |lambda|>=k,
# so sum |lambda|^-3 <= 8/B. Multiplication by 1/pi gives 8/(pi*B).
B=240
terms=[]
for mindex in range(-B,B+1):
sine=(0,1,0,-1)[mindex%4]
if not sine:continue
for nindex in range(-B,B+1):
terms.append(-mindex*sine/(mindex*mindex+nindex*nindex)**2)
K=math.fsum(terms)
signed=-K/math.pi
remainder=8/(math.pi*B)+1e-10
v=float(elliptic_d(-4,4,0,1,100).mid())
assert abs(signed-v)<remainder and abs(-signed-v)>remainder
return {'control_integral_x_squared':str(control),'theta_mahler_parameter':str(k),
'mahler_integral':str(m),'four_over_pi_times_D':str(target),
'mahler_accuracy_bits':int(m.accuracy()),'lattice_box':B,
'negative_sign_lattice_value':signed,'absolute_tail_bound_with_roundoff_margin':remainder,
'positive_sign_rejected':True,
'sign_resolution':'The common printed phase gives D=-Im(tau)^2 Re(K_2,1)/pi, agreeing with Brunault. Samart prints the opposite sign.'}
def main():
if os.environ.get('NUMBERDB_KEY_FROM_STDIN')=='1':
os.environ['NUMBERDB_API_KEY']=sys.stdin.read().strip()
mp.mp.dps=160
g=Table();params=list(g.enumerate())
assert len(params)==362
assert set(DISCRIMINANTS)=={-3,-12,-27,-4,-16,-7,-28,-8,-11,-19,-43,-67,-163}
# Independent full-group enumeration, quotient as frozensets, then filter.
expected=set()
for d in DISCRIMINANTS:
for N in (3,4,5,6):
orbits={frozenset(((r,s),((-r)%N,(-s)%N)))
for r in range(N) for s in range(N) if math.gcd(N,r,s)==1}
for orbit in orbits:
r,s=min(orbit)
if (2*s-(d%2)*r)%N:expected.add((d,N,r,s))
actual={tuple(p[k] for k in g.parameters) for p in params}
assert actual==expected
stored=None
if '--stored' in sys.argv:
stored=_stored_entries(numberdb.table('T700'))
assert set(stored)=={','.join(str(p[k]) for k in g.parameters) for p in params}
rows=[];least=10000;max_relative=mp.mpf(0);max_units=mp.mpf(0)
for p in params:
key=tuple(p[k] for k in g.parameters)
v,diag=elliptic_d(*key,100,diagnostics=True)
assert v.is_finite() and not v.contains_zero() and v.accuracy()>360
v2=elliptic_d(*key,130)
assert v2.is_finite() and v2.accuracy()>450 and (v-v2).contains_zero()
ind=fourier(*key)
lo,hi=as_mp_ball(v)
midpoint=(lo+hi)/2
error=abs(ind-midpoint)/abs(ind)
assert error < mp.mpf('1e-115'), (key, str(error))
# The independently computed point lies in the 100-digit interval.
text=stored[','.join(map(str,key))] if stored is not None else to_text(v,100)
wl,wh,ulp=parse_written(text)
assert wl <= ind <= wh, (key,text,mp.nstr(ind,110))
assert wl<=lo and hi<=wh
units=abs(ind-(wl+wh)/2)/ulp
# Negative control: changing the last place by four units must fail.
assert not wl+4*ulp <= ind <= wh+4*ulp
delta,N,r,s=key
opposite=elliptic_d(delta,N,(-r)%N,(-s)%N,100)
assert opposite.is_finite() and (v+opposite).contains_zero()
# Complex conjugation of x means s -> b*r-s.
conjugate=elliptic_d(delta,N,r,((delta%2)*r-s)%N,100)
assert conjugate.is_finite() and (v+conjugate).contains_zero()
least=min(least,diag['accuracy_bits']);max_relative=max(max_relative,error)
max_units=max(max_units,units)
rows.append({'params':p,'number':text,**diag,
'fourier_relative_difference':mp.nstr(error,8)})
zero_count=0
for delta in DISCRIMINANTS:
for N in range(1,7):
for r in range(N):
for s in range(N):
if math.gcd(N,r,s)!=1 or (r,s)>((-r)%N,(-s)%N):continue
if N<=2 or (2*s-(delta%2)*r)%N==0:
v=elliptic_d(delta,N,r,s,100)
lo,hi=as_mp_ball(v)
assert v.is_finite() and v.contains_zero() and hi-lo<mp.mpf('1e-115')
assert abs(fourier(delta,N,r,s))<mp.mpf('1e-130')
zero_count+=1
evidence={'table':'T700','mode':'stored' if stored is not None else 'pre-fill',
'skill_manifest':MANIFEST,'skill_sha256_local':SKILL_SHA,'skill_sha256_live':SKILL_SHA,
'source_sha256':{p:hashlib.sha256(Path(__file__).with_name(p).read_bytes()).hexdigest()
for p in ['generate.py','checks.py']},
'versions':{'sage':__import__('sage.version',fromlist=['version']).version,
'numberdb':numberdb.__version__,'mpmath':mp.__version__},
'entries':len(rows),'min_accuracy_bits':least,
'max_fourier_relative_difference':mp.nstr(max_relative,12),
'max_written_last_place_error':mp.nstr(max_units,12),
'shifted_last_place_rejection_controls':len(rows),
'oddness_checks':len(rows),'conjugation_checks':len(rows),
'excluded_zero_checks':zero_count,'mahler_and_lattice':mahler_check(),'rows':rows}
print(json.dumps(evidence,indent=2))
if __name__=='__main__':main()