from pathlib import Path
import re, math
import mpmath as mp

HERE=Path(__file__).resolve().parent
PAPER=HERE/"CSM_RH_Paper_84_Prime_Mertens_Cross_Spectrum_Derivative_Novelty_and_Horizontal_Exponent_Neutrality_v0.1_2026-09-09.md"
mp.mp.dps=70

def mobius(n):
    x=n; cnt=0; p=2
    while p*p<=x:
        if x%p==0:
            x//=p; cnt+=1
            if x%p==0:return 0
        p+=1
    if x>1:cnt+=1
    return -1 if cnt%2 else 1

def von_mangoldt(n):
    if n<2:return mp.mpf('0')
    x=n; p=2
    while p*p<=x:
        if x%p==0:
            while x%p==0:x//=p
            return mp.log(p) if x==1 else mp.mpf('0')
        p+=1
    return mp.log(n)

def divisors(n):
    out=[]
    for d in range(1,int(math.isqrt(n))+1):
        if n%d==0:
            out.append(d)
            if d*d!=n:out.append(n//d)
    return out

def check_convolution():
    vals=[]
    for n in range(1,300):
        lhs=sum(von_mangoldt(d)*mobius(n//d) for d in divisors(n))
        rhs=-mobius(n)*mp.log(n) if n>1 else mp.mpf('0')
        assert abs(lhs-rhs)<mp.mpf('1e-50'),(n,lhs,rhs)
        if n in [2,6,30,60,210]:
            vals.append((n,float(lhs),float(rhs)))
    return vals

def check_simple_zero_residues():
    rho=mp.zetazero(1)
    gamma=mp.gamma(rho)
    zp=mp.diff(mp.zeta,rho)
    vals=[]
    for p in [5,7,9]:
        eps=mp.mpf(10)**(-p)
        s=rho+eps
        prime=eps*mp.gamma(s)*(-mp.diff(mp.zeta,s)/mp.zeta(s))
        mert=eps*mp.gamma(s)/mp.zeta(s)
        vals.append((p,complex(prime),complex(mert)))
    assert abs(vals[-1][1]+complex(gamma))<1e-6
    assert abs(vals[-1][2]-complex(gamma/zp))<1e-6
    return vals

def check_gradient():
    s=mp.mpf('1.37')+mp.mpf('2.4')*1j
    z=mp.zeta(s)
    zp=mp.diff(mp.zeta,s)
    F=1/z
    Fp=-zp/(z*z)
    cross=(-zp/z)*(abs(F)**2)
    assert abs(cross-Fp*mp.conj(F))<mp.mpf('1e-55')
    h=mp.mpf('1e-6')
    def energy_sigma(sig):
        return abs(1/mp.zeta(sig+mp.im(s)*1j))**2
    ds=(energy_sigma(mp.re(s)+h)-energy_sigma(mp.re(s)-h))/(2*h)
    assert abs(mp.re(cross)-ds/2)<mp.mpf('1e-7')
    def energy_t(tt):
        return abs(1/mp.zeta(mp.re(s)+tt*1j))**2
    dt=(energy_t(mp.im(s)+h)-energy_t(mp.im(s)-h))/(2*h)
    assert abs(mp.im(cross)+dt/2)<mp.mpf('1e-7')
    return complex(cross),float(ds),float(dt)

def check_exponent():
    for beta in [0.51,0.63,0.8,0.95]:
        s=2*(1-beta)
        kappa=2*(1-beta)
        assert abs(s-kappa)<1e-14
        energy=4*(1-beta)
        assert abs(energy/2-kappa)<1e-14
    return True

def check_multiple_leading_factor():
    # algebraic coefficient m!/((m-1)! zeta^(m)) = m/zeta^(m)
    for m in [2,3,4,5]:
        assert math.factorial(m)/math.factorial(m-1)==m
    return True

def check_source():
    s=PAPER.read_text(encoding="utf-8")
    for pat in [r"(?<!\\)\\\(",r"(?<!\\)\\\)",r"(?<!\\)\\\[",r"(?<!\\)\\\]"]:
        assert re.search(pat,s) is None
    assert s.count("$$")%2==0
    tmp=re.sub(r"\$\$.*?\$\$","",s,flags=re.S)
    assert len(re.findall(r"(?<!\\)\$",tmp))%2==0
    return True

if __name__=="__main__":
    print("convolution",check_convolution())
    print("simple_zero_residues",check_simple_zero_residues())
    print("gradient",check_gradient())
    print("exponent",check_exponent())
    print("multiple_factor",check_multiple_leading_factor())
    print("source_delimiters",check_source())
    print("PASS")