"""Independent computations and stored-value checks for numberdb.org/T527. Run alongside generate.py and the unmodified campaign check.py: sage -python checks.py sage -python checks.py --stored # NUMBERDB_API_KEY for the private draft The polynomial resolvent computation and factorial-series quotient share no coefficient code with generate.py. See derivation.md for the ODE derivation. The script prints JSON evidence; redirect stdout to preserve it. No writes. """ import os import sys import json from decimal import Decimal from fractions import Fraction from math import factorial, comb import numberdb.sage as numberdb from numberdb._write import to_text from sage.rings.rational_field import QQ from sage.rings.real_arb import RealBallField from sage.version import version as sage_version import generate from check import exactness, measure, agrees_with, names_its_rings MANIFEST = 'agents/runs/20261026T204019Z-build-skills/manifest.json' SKILL_SHA = 'b4d51ded0d52cfe3dcc0859de0845c3879b0c2363273ba7e809a7e06ded40b18' def derivative(a): return [j*a[j] for j in range(1,len(a))] or [QQ(0)] def add(a,b): out=[QQ(0)]*max(len(a),len(b)) for j,x in enumerate(a):out[j]+=x for j,x in enumerate(b):out[j]+=x while len(out)>1 and out[-1]==0:out.pop() return out def scaled(a,c): return [x*c for x in a] def resolvent_q(count): # v_n = P_n + exp(-x)R_n; (d²-1)v_n = x v_(n-1), v_n(0)=0. p,r=[QQ(1)],[-QQ(1)] out=[QQ(1)] for n in range(1,count+1): a=[QQ(0)]+p pnew=[QQ(0)] while any(a): pnew=add(pnew,scaled(a,-1)) a=derivative(derivative(a)) a=[QQ(0)]+r sprime=[QQ(0)] power=QQ(2) while any(a): sprime=add(sprime,scaled(a,-QQ(1)/power)) a=derivative(a);power*=2 rnew=[-pnew[0]]+[sprime[j]/QQ(j+1) for j in range(len(sprime))] # Verify the polynomial identities, independently of the recurrence extraction. assert add(derivative(derivative(pnew)),scaled(pnew,-1))==add([QQ(0)]+p,[QQ(0)]) assert add(derivative(derivative(rnew)),scaled(derivative(rnew),-2))==add([QQ(0)]+r,[QQ(0)]) assert pnew[0]+rnew[0]==0 slope=derivative(pnew)[0]+derivative(rnew)[0]-rnew[0] out.append((-1)**n*slope) p,r=pnew,rnew return out def quotient_q(count): alpha=[QQ(factorial(6*n))/QQ(576**n*factorial(3*n)*factorial(2*n)) for n in range(count+1)] beta=[QQ(1)] q=[QQ(1)] for n in range(1,count+1): beta.append(alpha[n]+QQ(6*n-3)/QQ(4)*beta[-1]) q.append(beta[n]-sum((alpha[j]*q[n-j] for j in range(1,n+1)),QQ(0))) return q def direct_moment(n,q,field): w=sum((q[j]*q[n-j] for j in range(n+1)),QQ(0)) return field(2)**(-field(n)/2)*field(factorial(n))*field(w)/field(QQ(3*n+2)/QQ(2)).gamma() def written_interval(text): if '.' not in text and 'e' not in text.lower(): v=QQ(text);return v,v d=Decimal(text);f=Fraction(d);center=QQ(f.numerator)/QQ(f.denominator) unit=QQ(10)**int(d.as_tuple().exponent) return center-unit,center+unit def check_inside(text,ball): assert ball.is_finite() and ball.accuracy()>=700 lo,hi=written_interval(text) field=ball.parent() if lo==hi: assert lo in ball else: assert ball>field(lo) and ball=750 # Published values: Majumdar–Comtet (74); Janson Table 5 adds order 4. controls={0:QQ(1),1:QQ(1),2:QQ(17)/24,3:QQ(123)/140,4:QQ(2963)/3840} assert all(generate.rational_factor(n)==x for n,x in controls.items()) k,gq=generate.coefficients(count) assert not agrees_with(dict(enumerate(gq)),lambda n:q[n]) assert not agrees_with(dict(enumerate(gq)),lambda n:quotient_q(count)[n]) assert not exactness(dict(enumerate(k))) assert not exactness(dict(enumerate(gq))) assert not names_its_rings(generate.__file__) gen=generate.DoubleMeanderMoments() values={int(p['n']):gen.value(p,100) for p in gen.enumerate()} assert not exactness(values) written={n:to_text(v.get('number') if isinstance(v,dict) else v,100) for n,v in values.items()} fresh={n:direct_moment(n,q,high) for n in range(count+1)} for n in range(1,count): assert fresh[n-1]*fresh[n+1]-fresh[n]**2>0 # Probabilistic arcsine splitting, with independently computed meander moments. me=[high.pi().sqrt()*high(2)**(-high(n)/2)*factorial(n)*q[n] /high(QQ(3*n+1)/2).gamma() for n in range(count+1)] for n in range(count+1): conv=high(0) for j in range(n+1): a,b=high(QQ(3*j+1)/2),high(QQ(3*(n-j)+1)/2) conv+=comb(n,j)*me[j]*me[n-j]*a.gamma()*b.gamma()/(a+b).gamma()/high.pi() assert conv.is_finite() and conv.accuracy()>700 assert conv.overlaps(fresh[n]) rows=[] for n,text in written.items(): check_inside(text,fresh[n]) v=generate.moment(n) rows.append({'n':n,'exact_rational_factor':str(generate.rational_factor(n)), 'factor_multiplies_sqrt_2_over_pi':bool(n%2),'written':text, 'relative_accuracy_bits':int(v.accuracy()) if n%2 else None, 'independent_relative_accuracy_bits':int(fresh[n].accuracy()), 'written_interval_contains_independent_enclosure':True}) result={'table':'T527','skill_manifest':MANIFEST,'skill_sha256':SKILL_SHA, 'sage_version':sage_version,'working_bits':numberdb.bits(100,losing=generate.WORKING_GUARD), 'independent_bits':768,'minimum_odd_relative_accuracy_bits':min(x['relative_accuracy_bits'] for x in rows if x['n']%2), 'checks':{'resolvent_polynomial_identities':'all orders 1..12', 'coefficient_agreement_resolvent':'all orders 0..12', 'coefficient_agreement_factorial_series_quotient':'all orders 0..12', 'gamma_normalization_controls':True,'reflection_principle_mean':True, 'published_low_moments':'orders 0..4','arcsine_meander_convolution':'orders 0..12', 'strict_log_convexity':'orders 1..11, with M0=1', 'exactness_and_named_imports':True}, 'serialized_measurement':measure(written),'rows':rows} if stored: report=gen.verify(sample=None) assert report.ok,str(report) doc=numberdb.table('T527') entries=doc['Numbers'] if isinstance(entries,list): texts={int(x['params']['n']):str(x['number']) for x in entries} else: texts={int(n):str(v['number'] if isinstance(v,dict) else v) for n,v in entries.items()} assert list(texts)==list(range(1,count+1)) intervals={0:high(1)} for n,text in texts.items(): assert text==written[n] check_inside(text,fresh[n]) lo,hi=written_interval(text) intervals[n]=high((lo+hi)/2, (hi-lo)/2) # Reconstruct W_n from the database value, testing every displayed formula. wread=intervals[n]*high(2)**(high(n)/2)*high(QQ(3*n+2)/2).gamma()/factorial(n) wexact=sum((q[j]*q[n-j] for j in range(n+1)),QQ(0)) assert wread.is_finite() and wread.accuracy()>300 assert 0 in (wread-high(wexact)) for n in range(1,count): assert intervals[n-1]*intervals[n+1]-intervals[n]**2>0 result['stored_checks']={'generator_verify_all':True,'count':len(texts), 'independent_enclosures_contained':True, 'reconstructed_W_matches_independent_convolution':True, 'strict_log_convexity':True,'enumeration_matches':True} return result if __name__=='__main__': if os.environ.get('NUMBERDB_KEY_FROM_STDIN')=='1': key=sys.stdin.read().strip() if key.startswith('NUMBERDB_API_KEY='):key=key.split('=',1)[1].strip().strip('\"\'') os.environ['NUMBERDB_API_KEY']=key use_stored='--stored' in sys.argv or os.environ.get('NUMBERDB_RESTATING')=='1' print(json.dumps(run(use_stored),indent=2))