from pathlib import Path
import re, math, random

HERE=Path(__file__).resolve().parent
PAPER=HERE/"CSM_RH_Paper_64_Selberg_Symmetry_Pseudo_Prime_Boundary_Models_and_Ordinary_Arithmetic_Residual_v0.1_2026-09-09.md"

def F(x,d,gamma,eps):
    return x + eps*(x**(1-d))*math.cos(gamma*math.log(x))

def greedy_model(N,d=0.2,gamma=3.7,eps=0.02,n0=100):
    S=0.0
    selected=[]
    Fbase=F(n0-1,d,gamma,eps)
    errs=[]
    for n in range(n0,N+1):
        q=math.log(n)
        delta=F(n,d,gamma,eps)-F(n-1,d,gamma,eps)
        eprev=S-(F(n-1,d,gamma,eps)-Fbase)
        trial=eprev-delta
        if trial >= -q/2:
            b=0.0
        else:
            b=q
            selected.append(n)
        S += b
        e=S-(F(n,d,gamma,eps)-Fbase)
        errs.append((n,e,q))
    return selected,errs,S,Fbase

def check_greedy():
    selected,errs,S,Fbase=greedy_model(200000)
    worst=0.0
    for n,e,q in errs[100:]:
        ratio=abs(e)/(q+1)
        worst=max(worst,ratio)
        assert abs(e) <= q/2+3.0,(n,e,q)
    return worst,len(selected)

def check_density():
    selected,errs,S,Fbase=greedy_model(500000)
    N=500000
    expected=N/math.log(N)
    ratio=len(selected)/expected
    # finite-N rough sanity only
    assert 0.5<ratio<1.5,ratio
    return ratio

def check_boundary_derivative():
    d=0.2
    gamma=3.7
    eps=0.02
    for x in [1e3,1e5,1e8]:
        pert=eps*x**(-d)*((1-d)*math.cos(gamma*math.log(x))-gamma*math.sin(gamma*math.log(x)))
        assert abs(pert)<1
    return True

def check_BDH_scale():
    for d in [0.05,0.2,0.45]:
        for A in [2,10,100]:
            # compare logs: x^(2)/log^A x vs x^(2-2d) log x
            # ratio = x^(2d)/log^(A+1)x -> infinity.
            # sanity at extremely large log x chosen algebraically.
            L=max(1000,(A+1)*100/d)
            logr=2*d*L-(A+1)*math.log(L)
            assert logr>0
    return True

def check_selberg_exponent():
    # x^rho * x^(1-rho) = x at real exponent level.
    random.seed(64)
    for _ in range(1000):
        beta=random.uniform(0.01,0.99)
        assert abs(beta+(1-beta)-1)<1e-14
    return True

def check_source():
    s=PAPER.read_text(encoding="utf-8")
    forbidden=[
        r"(?<!\\)\\\(",
        r"(?<!\\)\\\)",
        r"(?<!\\)\\\[",
        r"(?<!\\)\\\]",
    ]
    for pat in forbidden:
        assert re.search(pat,s) is None,pat
    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("greedy",check_greedy())
    print("density_ratio",check_density())
    print("boundary_derivative",check_boundary_derivative())
    print("BDH_scale",check_BDH_scale())
    print("Selberg_exponent",check_selberg_exponent())
    print("source_delimiters",check_source())
    print("PASS")