"""profile_sums.py -- the profile sums  Sigma* = sum_k |W(k/(r1 s1 E))| / k,
   computed through the collapsed cotangent formula with a tail bound from
   |W(xi)| <= 2(1+log rs)(r+s)/(rs |xi|).  Empirically Sigma* <= 2 across
   structures, approx  min(1, s1 E / g) * O(1).
"""
from wlaw import W_cot
from math import gcd, log

def sigma_star(r, s, E, K=400):
    g = gcd(r, s); r1, s1 = r // g, s // g
    den = r1 * s1 * E
    tot = 0.0
    for k in range(1, K + 1):
        if (r * s * k) % den: continue
        tot += abs(W_cot(r, s, k, den)) / k
    L = 1 + log(r * s)
    tail = 2 * L * (r + s) * den / (r * s) * sum(1.0 / k**2 for k in range(K + 1, 10 * K))
    return tot, tail

if __name__ == "__main__":
    cases = [(10,15,1),(6,15,1),(21,35,1),(14,21,1),(194,291,1),(194,291,97),
             (115,145,1),(115,145,5),(101,103,1),(303,505,1),(26,39,13),
             (30,105,1),(30,105,15),(229,233,1),(458,687,1)]
    print(" r     s     g    E   Sigma*(head)  tail-bound")
    for r, s, E in cases:
        t, tl = sigma_star(r, s, E)
        print(f"{r:5d} {s:5d} {gcd(r,s):4d} {E:4d}   {t:9.4f}   {tl:8.3f}")