repair-program.py

back to table · edit · history · where entries came from · files · download

815 bytes, as of the version from 2026-11-02 17:54. Recorded here, not run.

import numberdb.sage
from sage.rings.integer_ring import ZZ
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField
from sage.rings.complex_arb import ComplexBallField


def moment(n, k, bits=512):
    n, k = ZZ(n), ZZ(k)
    if n < 2 or k not in (2, 3):
        raise ValueError('Require integer n >= 2 and k = 2 or 3.')
    R, C = RealBallField(bits), ComplexBallField(bits)
    pi = R.pi()
    if k == 2:
        return pi**2 * QQ(n - 1) / (6*n)
    total = R(0)
    for r in range(1, n):
        if 2*r == n:
            continue  # cot(pi/2) is exactly zero.
        a = QQ(r) / n
        theta = pi * R(a)
        total += theta.cos()/theta.sin() * C(2).zeta(C(a)).real()
    return (QQ(9*n - 6 - 3*n**3)/(2*n**2) * R(3).zeta()
            + 3*pi*total/(2*n**2))


moment(6, 3)