Sharedwww / talks / bernoulli / genber.sageOpen in CoCalc
Author: William A. Stein
1
2
def genber_numerical(n, chi, prec, primes):
3
if not chi.is_primitive():
4
raise ValueError, "chi (=%s) must be primitive"%chi
5
CC = ComplexField(prec)
6
K = genber_K(n, chi, CC)
7
a = CC.pi()^n
8
g = chi.gauss_sum_numerical(prec)
9
L = genber_L(n, chi^(-1), CC, primes)
10
d = genber_d(n, chi, prec)
11
print "d = ", d
12
print "K = ", K
13
print "L = ", L
14
print "a = ", a
15
print "g = ", g
16
num = d*K*L/(a*g)
17
return num
18
19
def genber_d(n, chi, prec):
20
f = chi.conductor()
21
if f == 1:
22
raise NotImplementedError
23
F = f.factor()
24
assert len(F) >= 1
25
if len(F) > 1:
26
return 1
27
elif f == 4:
28
return 2
29
P = F[0]
30
if P[0] == 2 and P[1] > 2:
31
return 1
32
elif P[1] == 1:
33
return n*P[0]
34
elif P[1] > 1:
35
R = chi.base_ring()
36
phi = R.complex_embedding(prec)
37
return phi(1 - chi(1+p))
38
39
40
def genber_K(n, chi, CC):
41
f = chi.modulus()
42
i = CC.0
43
z = 2*factorial(n)*(f/(2*i))^n
44
if n % 2 == 0:
45
z *= -1
46
return z
47
48
def genber_L(n, chi, CC, max_prime):
49
R = chi.base_ring()
50
m = chi.modulus()
51
phi = R.complex_embedding(CC.prec())
52
n = phi(n)
53
z = 1
54
zero = 2^(-CC.prec())
55
z_prev = z
56
for p in prime_range(max_prime):
57
h = z * phi(chi(p)) / (CC(p)^n)
58
z = z - h
59
print p, abs(z - z_prev)
60
z_prev = z
61
return 1/z
62
63
def test1(prec, primes):
64
e = DirichletGroup(13).0
65
print "group: ", e.parent()
66
print "chi = ", e
67
B = genber_numerical(389, e, prec, primes)
68
print "B3 = ", B
69
return e, B
70
71
72