
def genber_numerical(n, chi, prec, primes):
    if not chi.is_primitive():
        raise ValueError, "chi (=%s) must be primitive"%chi
    CC = ComplexField(prec)
    K = genber_K(n, chi, CC)
    a = CC.pi()**n
    g = chi.gauss_sum_numerical(prec)
    L = genber_L(n, chi**(-ZZ(1)), CC, primes)
    d = genber_d(n, chi, prec)
    print "d = ", d
    print "K = ", K
    print "L = ", L
    print "a = ", a
    print "g = ", g
    num = d*K*L/(a*g)
    return num

def genber_d(n, chi, prec):
    f = chi.conductor()
    if f == ZZ(1):
        raise NotImplementedError
    F = f.factor()
    assert len(F) >= ZZ(1)
    if len(F) > ZZ(1):
        return ZZ(1)
    elif f == ZZ(4):
        return ZZ(2)
    P = F[0]
    if P[0] == ZZ(2) and P[1] > ZZ(2):
        return ZZ(1)
    elif P[1] == ZZ(1):
        return n*P[0]
    elif P[1] > ZZ(1):
        R = chi.base_ring()
        phi = R.complex_embedding(prec)
        return phi(ZZ(1) - chi(ZZ(1)+p))


def genber_K(n, chi, CC):
    f = chi.modulus()
    i = CC.gen(0)
    z = ZZ(2)*factorial(n)*(f/(ZZ(2)*i))**n
    if n % ZZ(2) == ZZ(0):
        z *= -ZZ(1)
    return z

def genber_L(n, chi, CC, max_prime):
    R = chi.base_ring()
    m = chi.modulus()
    phi = R.complex_embedding(CC.prec())
    n = phi(n)
    z = ZZ(1)
    zero = ZZ(2)**(-CC.prec())
    z_prev = z
    for p in prime_range(max_prime):
        h = z * phi(chi(p)) / (CC(p)**n)
        z = z - h
        print p, abs(z - z_prev)
        z_prev = z
    return ZZ(1)/z

def test1(prec, primes):
    e = DirichletGroup(ZZ(13)).gen(0)
    print "group: ", e.parent()
    print "chi = ", e
    B = genber_numerical(ZZ(389), e, prec, primes)
    print "B3 = ", B
    return e, B
    
    
