"""Independent Arb enclosure of C4 and the C3(delta) hand-off.
Different prime cutoffs and a direct coefficient recurrence; imports no author code.
"""
import os
for k in ['OPENBLAS_NUM_THREADS','OMP_NUM_THREADS','MKL_NUM_THREADS','NUMEXPR_NUM_THREADS']:os.environ[k]='1'
from pathlib import Path
import math,json
from flint import arb,ctx
ctx.prec=256;ctx.threads=1
OUT=Path(__file__).resolve().parent/'resume1'
def primes(n):
    s=bytearray(b'\x01')*(n+1);s[:2]=b'\x00\x00'
    for p in range(2,math.isqrt(n)+1):
        if s[p]:s[p*p::p]=b'\x00'*len(s[p*p::p])
    return [i for i in range(13,n+1) if s[i]]
BIG,SMALL,SPLIT=200000,3000,20
pr=primes(BIG); ps=[p for p in pr if p<=SMALL]
def bounds(x): return [x.lower().str(28),x.upper().str(28)]
def point_up(x):return arb(x.upper())
def coefficients(LX):
    km=int(float(LX)/math.log(13)); lx=arb(LX)
    assert km*arb(13).log()<lx<(km+1)*arb(13).log()
    S=point_up(lx.log()+arb('0.26149721284765')+1/(10*lx**2)+arb(4)/(15*lx**3)-sum((arb(1)/p for p in [2,3,5,7,11]),arb(0)))
    P=[arb(0),S]; A=[arb(0),arb(0)]
    lp={p:arb(p).log() for p in pr}
    for k in range(2,km+1):
        cutoff,plist=(BIG,pr) if k<=SPLIT else (SMALL,ps)
        vals=[arb(p)**(-k) for p in plist]
        finite=sum(vals,arb(0)); logsum=sum((lp[p]*v for p,v in zip(plist,vals)),arb(0))
        cut=arb(cutoff)
        P.append(point_up(finite+cut**(1-k)/(k-1)))
        A.append(point_up(k*(logsum+cut**(1-k)*(cut.log()/(k-1)+1/arb(k-1)**2))))
    R=[arb(1)]
    for n in range(1,km+1):
        R.append(point_up(sum((P[k]*R[n-k] for k in range(1,n+1)),arb(0))/n))
    return km,R,A,S
ct=1+arb(1)/36260
results={}
for LX,Li,K in [('2736.86','1800',18),('1800',None,18)]:
    km,R,A,S=coefficients(LX)
    if Li:
        li=arb(Li); d=d2=arb('.02');logY=arb(8000000000).log()
        assert li<=arb(LX) and d+d2<1 and d*li>=1 and d2*li>=2 and li>logY
        rsum=sum(R[K:km],arb(0))
        windows=sum((ct*(ct*(1+logY)*R[k-2]+sum((A[a]*R[k-1-a] for a in range(2,k)),arb(0)))/(li-logY) for k in range(K+1,km+1)),arb(0))
        small=li*(-d*li).exp();sq=2*arb('1.5').zeta()*li*(-d2*li/2).exp()
        main=(1+arb(480)/(840*logY))*rsum/(1-d-d2)
        win=arb(480)/135*windows/(1-d-d2)
        rest=arb(480)/135*(small+sq)+win
        lhs=main+rest
        rhs=(1-arb(480)/(840*(arb('4e18')-2).log()))*(1-arb('1e-18'))-arb(480)/135*arb('1e-18')
        record={k:v.str(28) for k,v in [('S',S),('main',main),('window',win),('small',small),('squarefull',sq),('rest',rest),('lhs_majorant',lhs),('rhs',rhs)]}
        assert lhs<arb('.4068037') and lhs<rhs
    else:
        li=arb('4e18').log(); d=arb('.2')
        total=sum((ct*R[k-1]+sum((A[a]*R[k-a] for a in range(2,k+1)),arb(0)) for k in range(K+1,km+1)),arb(0))
        lhs=li*(-d*li).exp()+total/(1-d)
        rhs=135*(arb(1)/480-1/(840*(arb('4e18')-2).log()))*(1-arb('1e-18'))-arb('1e-18')
        record={'S':S.str(28),'lhs_majorant':lhs.str(28),'rhs':rhs.str(28)}
        assert lhs<arb('.209801') and lhs<rhs
    record['kmax']=km
    results[LX]=record;print(LX,json.dumps(record),flush=True)
(OUT/'audit_c4.json').write_text(json.dumps(results,indent=2)+'\n')