back to table · edit · history · where entries came from · files · download
7069 bytes, as of the version from 2026-11-02 14:37 (current). Recorded here, not run.
"""Independent checks for T679. Run: sage -python checks.py.
No hypergeometric library is used by the alternative computations. Smooth
values use Taylor continuation of the differential equation. Conifold values
use a finite coefficient sum with a Bernoulli-polynomial gamma-ratio tail,
summed using Hurwitz zeta, independently of mpmath's Richardson/Shanks method.
A rank-two elliptic control, a telescoping tail control, and a modular L-value
provide checks outside the 14 rank-four systems. Heuristic checks, not proofs
of error bounds. The serialized JSON lines are the numerical evidence.
"""
import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from sage.libs.pari import pari
from sage.rings.real_mpfr import RealField
import mpmath as mp
import math, json, time
from generate import HolomorphicPeriods, components
# This explicit pair list is independent of the generator's cyclotomic parsing.
PAIRS=[('2.2.2.2',2,2),('2.2.3',2,3),('2.2.4',2,4),('2.2.6',2,6),
('3.3',3,3),('3.4',3,4),('3.6',3,6),('4.4',4,4),('4.6',4,6),('6.6',6,6),
('5',5,QQ(5)/2),('8',8,QQ(8)/3),('10',10,QQ(10)/3),('12',12,QQ(12)/5)]
def alpha_of(a,b):
r=mp.mpf(str(QQ(1)/a)); s=mp.mpf(str(QQ(1)/b))
return [r,1-r,s,1-s]
def product(xs):
ans=mp.mpf(1)
for x in xs: ans*=x
return ans
def falling(n,k):
return math.prod(range(n-k+1,n+1))
def ode_polynomials(alpha):
# theta^r - z prod(theta+alpha); coefficient of ordinary derivative j.
r=len(alpha)
symmetric=[mp.mpf(1)]
for a in alpha:
new=[mp.mpf(0)]*(len(symmetric)+1)
for j,c in enumerate(symmetric): new[j]+=a*c; new[j+1]+=c
symmetric=new
stirling=[[0]*(r+1) for _ in range(r+1)]; stirling[0][0]=1
for n in range(1,r+1):
for j in range(1,n+1): stirling[n][j]=stirling[n-1][j-1]+j*stirling[n-1][j]
polynomials=[]
for j in range(r+1):
p=[mp.mpf(0)]*(j+2)
p[j]=stirling[r][j]
p[j+1]=-sum(symmetric[k]*stirling[k][j] for k in range(j,r+1))
polynomials.append(p)
return polynomials
def initial_jet(alpha,x,terms):
r=len(alpha); jet=[mp.mpf(0)]*r; c=mp.mpf(1)
for n in range(terms):
for j in range(min(r,n+1)):
jet[j]+=c*math.comb(n,j)*x**(n-j)
c*=product(a+n for a in alpha)/(mp.mpf(n+1)**r)
return jet
def step(polynomials,jet,x,h,terms):
r=len(jet)
shifted=[]
for p in polynomials:
shifted.append([sum(p[k]*math.comb(k,l)*x**(k-l) for k in range(l,len(p)))
for l in range(len(p))])
coef=list(jet)
for n in range(terms-r):
rhs=mp.mpc(0)
for j,p in enumerate(shifted):
for l,v in enumerate(p):
k=n-l+j
if n>=l and (j,l)!=(r,0):
rhs+=v*falling(k,j)*coef[k]
coef.append(-rhs/(shifted[r][0]*falling(n+r,r)))
return [sum(coef[n]*math.comb(n,j)*h**(n-j) for n in range(j,len(coef)))
for j in range(r)]
def continuation(alpha,z,working=150):
# Fixed rational paths avoid 0 and 1; each step is under half the nearest
# singularity distance. This handles repeated exponents without limits.
terms=4*working+80
if z<0:
x=mp.mpf('-0.25'); targets=[mp.mpf('-0.375'),mp.mpf('-0.5'),mp.mpf('-0.75'),mp.mpf(-1)]
else:
x=mp.mpf('0.25')
if z==mp.mpf('0.5'): targets=[mp.mpf('0.375'),mp.mpf('0.5')]
else: targets=[mp.mpc('0.25','0.125'),mp.mpc('0.25','0.25'),mp.mpc('0.25','0.375'),
mp.mpc('0.25','0.5'),mp.mpc('0.5','0.5'),mp.mpc('0.75','0.5'),
mp.mpc('1','0.5'),mp.mpc('1.25','0.5'),mp.mpc('1.5','0.5'),
mp.mpc('1.75','0.5'),mp.mpc('2','0.5'),mp.mpc('2','0.25'),mp.mpf(2)]
p=ode_polynomials(alpha); jet=initial_jet(alpha,x,terms)
for target in targets:
jet=step(p,jet,x,target-x,terms); x=target
return jet[0]
def conifold(alpha,lower=None,N=1400,K=130):
if lower is None: lower=[mp.mpf(1)]*(len(alpha)-1)
beta=list(lower)+[mp.mpf(1)]
exponent=sum(beta)-sum(alpha)
logcoef=[mp.mpf(0)]
for k in range(1,K+1):
v=sum(mp.bernpoly(k+1,a) for a in alpha)-sum(mp.bernpoly(k+1,b) for b in beta)
logcoef.append((-1)**(k+1)*v/(mp.mpf(k)*(k+1)))
coeff=[mp.mpf(1)]
for n in range(1,K+1):
coeff.append(sum(k*logcoef[k]*coeff[n-k] for k in range(1,n+1))/n)
c=mp.mpf(1); head=mp.mpf(0)
for n in range(N):
head+=c
c*=product(a+n for a in alpha)/product(b+n for b in beta)
C=product(mp.gamma(b) for b in beta)/product(mp.gamma(a) for a in alpha)
tail=C*sum(coeff[k]*mp.zeta(exponent+k,N) for k in range(K+1))
return head+tail
def close(actual,expected,digits,label):
err=abs(actual-expected)/max(abs(expected),mp.mpf('0.01'))
if not mp.isfinite(actual) or not mp.isfinite(expected) or err>=mp.mpf(10)**(-digits):
raise AssertionError((label,mp.nstr(err,15)))
return float(-mp.log10(err)) if err else 999.0
def main():
mp.mp.dps=165
start=time.monotonic()
from sage.env import SAGE_VERSION
print(json.dumps({'versions':{'Sage':SAGE_VERSION,'mpmath':mp.__version__,'PARI':str(pari('version()'))},'skill_manifest':'agents/runs/20261102T134642Z-build-skills/manifest.json'}),flush=True)
# Known telescoping sum: 2F1(1,1;3;1)=2.
control=conifold([mp.mpf(1)]*2,[mp.mpf(3)])
print(json.dumps({'control':'telescoping conifold sum','agreement_digits':close(control,mp.mpf(2),125,'tail control')}),flush=True)
# Rank-two control on all smooth paths, including upper-bank continuation.
for z in [mp.mpf('0.5'),mp.mpf(-1),mp.mpf(2)]:
y=continuation([mp.mpf('0.5')]*2,z,150)
ref=2*mp.ellipk(z)/mp.pi
if z>1: ref=mp.conj(ref)
print(json.dumps({'control':'elliptic continuation','z':str(z),'agreement_digits':close(y,ref,125,'elliptic')}),flush=True)
# PARI constructs the weight-four level-eight newform, not a hypergeometric.
pari.set_real_precision(165)
L=pari('lfun(concat(lfunmf(mfinit([8,4],0)))[1],2)')
modular=16*mp.mpf(str(L))/mp.pi**2
y=conifold([mp.mpf('0.5')]*4)
print(json.dumps({'control':'16 L(eta(2tau)^4 eta(4tau)^4,2)/pi^2','value':mp.nstr(modular,140),'agreement_digits':close(y,modular,125,'modular')}),flush=True)
for template,a,b in PAIRS:
alpha=alpha_of(a,b)
for z in ['1','1/2','-1','2']:
x=mp.mpf(z)
y=conifold(alpha) if x==1 else continuation(alpha,x,150)
c=components(template,z,190); target=mp.mpc(*c)
agreement=close(y,target,120,(template,z))
# Record both independent values for later comparison with database.
print(json.dumps({'template':template,'z':z,'independent_re':mp.nstr(y.real,140),'independent_im':mp.nstr(y.imag,140),'generator_re':c[0],'generator_im':c[1],'agreement_digits':agreement}),flush=True)
print(json.dumps({'passed':True,'entries_checked':56,'seconds':time.monotonic()-start}),flush=True)
if __name__=='__main__': main()