cesaro.sagews

theta = var('theta')
def f(C):
    T = SR(0)
    for q in prime_powers(C+1):
        if q > 1:
            p, n = factor(q)[0]
            T += 2 * p^(-n/2) * log(p) * cos(n*log(p)*theta)
    return T

h = f(5)

plot(h^2, 5, 40)

k = h^2

T = var('T')
assume(T>0)
g(T) = 1/T * k.integrate(theta, 0, T)

plot(g, 0, 40, color='red') + plot(k,0,40)