CoCalc Shared Filescode / eulerprod.pyOpen in CoCalc with one click!
Authors: Jennifer Balakrishnan, Alex Best, Maarten Derickx, Nicholas Triantafillou, Jan Vonk
1
"""
2
General L-series
3
4
AUTHOR:
5
- William Stein
6
7
TODO:
8
- Symmetric powers (and modular degree -- see trac 9758)
9
- Triple product L-functions: Gross-Kudla, Zhang, etc -- see the code in triple_prod/triple.py
10
- Support L-calc L-function
11
- Make it so we use exactly one GP session for *all* of the Dokchitser L-functions
12
- Tensor products
13
- Genus 2 curves, via smalljac and genus2reduction
14
- Fast L-series of elliptic curves over number fields (not just sqrt(5)), via smalljac
15
- Inverse of number_of_coefficients function.
16
"""
17
18
#################################################################################
19
#
20
# (c) Copyright 2011 William Stein
21
#
22
# This file is part of PSAGE
23
#
24
# PSAGE is free software: you can redistribute it and/or modify
25
# it under the terms of the GNU General Public License as published by
26
# the Free Software Foundation, either version 3 of the License, or
27
# (at your option) any later version.
28
#
29
# PSAGE is distributed in the hope that it will be useful,
30
# but WITHOUT ANY WARRANTY; without even the implied warranty of
31
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
32
# GNU General Public License for more details.
33
#
34
# You should have received a copy of the GNU General Public License
35
# along with this program. If not, see <http://www.gnu.org/licenses/>.
36
#
37
#################################################################################
38
39
import copy, math, types
40
41
from sage.all import prime_range, cached_method, sqrt, SR, vector, gcd
42
from sage.rings.all import ZZ, Integer, QQ, O, ComplexField, CDF, infinity as oo
43
from sage.arith.all import primes
44
from sage.rings.rational_field import is_RationalField
45
from sage.schemes.elliptic_curves.ell_generic import is_EllipticCurve
46
47
from psage.ellcurve.lseries.helper import extend_multiplicatively_generic
48
#import pyximport; pyximport.install()
49
#from helper import extend_multiplicatively_generic
50
51
from sage.misc.all import prod
52
from sage.modular.abvar.abvar import is_ModularAbelianVariety
53
from sage.modular.dirichlet import is_DirichletCharacter
54
from sage.lfunctions.dokchitser import Dokchitser
55
from sage.modular.modsym.space import is_ModularSymbolsSpace
56
from sage.modular.abvar.abvar import is_ModularAbelianVariety
57
from sage.rings.number_field.number_field_base import is_NumberField
58
import sage.modular.modform.element
59
from sage.modular.all import Newform
60
from sage.structure.factorization import Factorization
61
from sage.misc.mrange import cartesian_product_iterator
62
63
I = sqrt(-1)
64
65
def prec(s):
66
"""
67
Return precision of s, if it has a precision attribute. Otherwise
68
return 53. This is useful internally in this module.
69
70
EXAMPLES::
71
72
sage: from psage.lseries.eulerprod import prec
73
sage: prec(ComplexField(100)(1))
74
100
75
sage: prec(RealField(125)(1))
76
125
77
sage: prec(1/3)
78
53
79
"""
80
if hasattr(s, 'prec'):
81
return s.prec()
82
return 53
83
84
def norm(a):
85
"""
86
Return the norm of a, for a in either a number field or QQ.
87
88
This is a function used internally in this module, mainly because
89
elements of QQ and ZZ have no norm method.
90
91
EXAMPLES::
92
93
sage: from psage.lseries.eulerprod import norm
94
sage: K.<a> = NumberField(x^2-x-1)
95
sage: norm(a+5)
96
29
97
sage: (a+5).norm()
98
29
99
sage: norm(17)
100
17
101
"""
102
try:
103
return a.norm()
104
except AttributeError:
105
return a
106
107
def tiny(prec):
108
"""
109
Return a number that we consider tiny to the given precision prec
110
in bits. This is used in various places as "zero" to the given
111
precision for various checks, e.g., of correctness of the
112
functional equation.
113
"""
114
return max(1e-8, 1.0/2**(prec-1))
115
116
def prime_below(P):
117
"""
118
Return the prime in ZZ below the prime P (or element of QQ).
119
120
EXAMPLES::
121
122
sage: from psage.lseries.eulerprod import prime_below
123
sage: K.<a> = NumberField(x^2-x-1)
124
sage: prime_below(K.prime_above(11))
125
11
126
sage: prime_below(K.prime_above(5))
127
5
128
sage: prime_below(K.prime_above(3))
129
3
130
sage: prime_below(7)
131
7
132
"""
133
try:
134
return P.smallest_integer()
135
except AttributeError:
136
return ZZ(P)
137
138
139
class LSeriesDerivative(object):
140
"""
141
The formal derivative of an L-series.
142
143
EXAMPLES::
144
145
sage: from psage.lseries.eulerprod import LSeries
146
sage: L = LSeries('delta')
147
sage: L.derivative()
148
First derivative of L-function associated to Ramanujan's Delta (a weight 12 cusp form)
149
sage: L.derivative()(11/2)
150
0.125386233743526
151
152
We directly create an instance of the class (users shouldn't need to do this)::
153
154
sage: from psage.lseries.eulerprod import LSeriesDerivative
155
sage: Ld = LSeriesDerivative(L, 2); Ld
156
Second derivative of L-function associated to Ramanujan's Delta (a weight 12 cusp form)
157
sage: type(Ld)
158
<class 'psage.lseries.eulerprod.LSeriesDerivative'>
159
"""
160
def __init__(self, lseries, k):
161
"""
162
INPUT:
163
- lseries -- any LSeries object (derives from LseriesAbstract)
164
- k -- positive integer
165
"""
166
k = ZZ(k)
167
if k <= 0:
168
raise ValueError, "k must be a positive integer"
169
self._lseries = lseries
170
self._k = k
171
172
def __cmp__(self, right):
173
return cmp((self._lseries,self._k), (right._lseries, right._k))
174
175
def __call__(self, s):
176
"""
177
Return the value of this derivative at s, which must coerce to a
178
complex number. The computation is to the same precision as s,
179
or to 53 bits of precision if s is exact.
180
181
As usual, if s has large imaginary part, then there could be
182
substantial precision loss (a warning is printed in that case).
183
184
EXAMPLES::
185
186
sage: from psage.lseries.eulerprod import LSeries
187
sage: L = LSeries(EllipticCurve('389a'))
188
sage: L1 = L.derivative(1); L1(1)
189
-1.94715429754927e-20
190
sage: L1(2)
191
0.436337613850735
192
sage: L1(I)
193
-19.8890471908356 + 31.2633280771869*I
194
sage: L2 = L.derivative(2); L2(1)
195
1.51863300057685
196
sage: L2(I)
197
134.536162459604 - 62.6542402272310*I
198
"""
199
return self._lseries._function(prec(s)).derivative(s, self._k)
200
201
def __repr__(self):
202
"""
203
String representation of this derivative of an L-series.
204
205
EXAMPLES::
206
207
sage: from psage.lseries.eulerprod import LSeries
208
sage: L = LSeries('delta')
209
sage: L.derivative(1).__repr__()
210
"First derivative of L-function associated to Ramanujan's Delta (a weight 12 cusp form)"
211
sage: L.derivative(2).__repr__()
212
"Second derivative of L-function associated to Ramanujan's Delta (a weight 12 cusp form)"
213
sage: L.derivative(3).__repr__()
214
"Third derivative of L-function associated to Ramanujan's Delta (a weight 12 cusp form)"
215
sage: L.derivative(4).__repr__()
216
"4-th derivative of L-function associated to Ramanujan's Delta (a weight 12 cusp form)"
217
sage: L.derivative(2011).__repr__()
218
"2011-th derivative of L-function associated to Ramanujan's Delta (a weight 12 cusp form)"
219
"""
220
k = self._k
221
if k == 1:
222
kth = 'First'
223
elif k == 2:
224
kth = 'Second'
225
elif k == 3:
226
kth = 'Third'
227
else:
228
kth = '%s-th'%k
229
return "%s derivative of %s"%(kth, self._lseries)
230
231
def derivative(self, k=1):
232
"""
233
Return the k-th derivative of this derivative object.
234
235
EXAMPLES::
236
237
sage: from psage.lseries.eulerprod import LSeries
238
sage: f = Newforms(43,2,names='a')[1]; f
239
q + a1*q^2 - a1*q^3 + (-a1 + 2)*q^5 + O(q^6)
240
sage: L = LSeries(f); L1 = L.derivative()
241
sage: L1(1)
242
0.331674007376949
243
sage: L(1)
244
0.620539857407845
245
sage: L = LSeries(f); L1 = L.derivative(); L1
246
First derivative of L-series of a degree 2 newform of level 43 and weight 2
247
sage: L1.derivative()
248
Second derivative of L-series of a degree 2 newform of level 43 and weight 2
249
sage: L1.derivative(3)
250
4-th derivative of L-series of a degree 2 newform of level 43 and weight 2
251
252
"""
253
if k == 0:
254
return self
255
return LSeriesDerivative(self._lseries, self._k + k)
256
257
class LSeriesParentClass(object):
258
def __contains__(self, x):
259
return isinstance(x, (LSeriesAbstract, LSeriesProduct))
260
261
def __call__(self, x):
262
"""
263
EXAMPLES::
264
265
sage: from psage.lseries.eulerprod import LSeries
266
sage: L = LSeries('zeta')
267
sage: P = L.parent(); P
268
All L-series objects and their products
269
sage: P(L) is L
270
True
271
272
sage: P(L^3)
273
(Riemann Zeta function viewed as an L-series)^3
274
sage: (L^3).parent()
275
All L-series objects and their products
276
277
You can make the L-series attached to an object by coercing
278
that object into the parent::
279
280
sage: P(EllipticCurve('11a'))
281
L-series of Elliptic Curve defined by y^2 + y = x^3 - x^2 - 10*x - 20 over Rational Field
282
283
We also allow coercing in 1, because this is very useful for
284
formal factorizations and products::
285
286
sage: P(1)
287
1
288
289
Other numbers do not coerce in, of course::
290
291
sage: P(2)
292
Traceback (most recent call last):
293
...
294
TypeError
295
"""
296
if isinstance(x, LSeriesAbstract):
297
return x
298
elif isinstance(x, LSeriesProduct):
299
return x
300
elif x == 1:
301
return x
302
else:
303
try:
304
return LSeries(x)
305
except NotImplementedError:
306
raise TypeError
307
308
def __repr__(self):
309
"""
310
Return string representation of this parent object.
311
312
sage: from psage.lseries.eulerprod import LSeriesParent
313
sage: LSeriesParent.__repr__()
314
'All L-series objects and their products'
315
"""
316
return "All L-series objects and their products"
317
318
def __cmp__(self, right):
319
"""
320
Returns equality if right is an instance of
321
LSeriesParentClass; otherwise compare memory locations.
322
323
EXAMPLES::
324
325
sage: from psage.lseries.eulerprod import LSeriesParent, LSeriesParentClass
326
sage: cmp(LSeriesParent, LSeriesParentClass())
327
0
328
sage: cmp(LSeriesParent, 7) != 0
329
True
330
sage: cmp(7, LSeriesParent) == -cmp(LSeriesParent, 7)
331
True
332
"""
333
if isinstance(right, LSeriesParentClass):
334
return 0
335
return cmp(id(self), id(right))
336
337
338
LSeriesParent = LSeriesParentClass()
339
340
class LSeriesAbstract(object):
341
r"""
342
L-series defined by an Euler product.
343
344
The parameters that define the 'shape' of the L-series are:
345
346
conductor, hodge_numbers, weight, epsilon, poles, base_field
347
348
Let gamma(s) = prod( Gamma((s+h)/2) for h in hodge_numbers ).
349
Denote this L-series by L(s), and let `L^*(s) = A^s \gamma(s) L(s)`, where
350
`A = \sqrt(N)/\pi^{d/2}`, where d = len(hodge_numbers) and N = conductor.
351
Then the functional equation is
352
353
Lstar(s) = epsilon * Lstar(weight - s).
354
355
To actually use this class we create a derived class that in
356
addition implements a method _local_factor(P), that takes as input
357
a prime ideal P of K=base_field, and returns a polynomial, which
358
is typically the reversed characteristic polynomial of Frobenius
359
at P of Gal(Kbar/K) acting on the maximal unramified quotient of
360
some Galois representation. This class automatically computes the
361
Dirichlet series coefficients `a_n` from the local factors of the
362
`L`-function.
363
364
The derived class may optionally -- and purely as an optimization
365
-- define a method self._precompute_local_factors(bound,
366
prec=None), which is typically called before
367
368
[_local_factor(P) for P with norm(P) < bound]
369
370
is called in the course of various computations.
371
"""
372
def __init__(self,
373
conductor,
374
hodge_numbers,
375
weight,
376
epsilon,
377
poles,
378
residues,
379
base_field,
380
is_selfdual=True,
381
prec=53):
382
"""
383
INPUT:
384
- ``conductor`` -- list or number (in a subset of the
385
positive real numbers); if the conductor is a list, then
386
each conductor is tried in order (to the precision prec
387
below) until we find one that works.
388
- ``hodge_numbers`` -- list of numbers (in a subring of the complex numbers)
389
- ``weight`` -- number (in a subset of the positive real numbers)
390
- ``epsilon`` -- number (in a subring of the complex numbers)
391
- ``poles`` -- list of numbers (in subring of complex numbers); poles of the *completed* L-function
392
- ``residues`` -- list of residues at each pole given in poles or string "automatic"
393
- ``base_field`` -- QQ or a number field; local L-factors
394
correspond to nonzero prime ideals of this field.
395
- ``is_selfdual`` -- bool (default: True)
396
- ``prec`` -- integer (default: 53); precision to use when trying to figure
397
out parameters using the functional equation
398
399
400
EXAMPLES::
401
402
sage: from psage.lseries.eulerprod import LSeriesAbstract
403
sage: L = LSeriesAbstract(conductor=1, hodge_numbers=[0], weight=1, epsilon=1, poles=[1], residues=[-1], base_field=QQ)
404
sage: type(L)
405
<class 'psage.lseries.eulerprod.LSeriesAbstract'>
406
sage: L._conductor
407
1
408
sage: L._hodge_numbers
409
[0]
410
sage: L._weight
411
1
412
sage: L._epsilon
413
1
414
sage: L._poles
415
[1]
416
sage: L._residues
417
[-1]
418
sage: L._base_field
419
Rational Field
420
sage: L
421
Euler Product L-series with conductor 1, Hodge numbers [0], weight 1, epsilon 1, poles [1], residues [-1] over Rational Field
422
"""
423
self._anlist = {None:[]}
424
425
(self._conductor, self._hodge_numbers, self._weight, self._epsilon,
426
self._poles, self._residues, self._base_field, self._is_selfdual) = (
427
conductor, hodge_numbers, weight, epsilon, poles, residues,
428
base_field, is_selfdual)
429
430
# the following parameters allow for specifying a list of
431
# possibilities:
432
#
433
# conductor -- list of possible conductors
434
# hodge_numbers -- give a list of lists
435
# weight -- list of possible weights
436
# poles (residues must be the default 'automatic')
437
# epsilon -- list of possible epsilon's.
438
439
# 1. Figure out for which parameters we have multiple options
440
# 2. Run through them until checking each until one is found
441
# that works.
442
v = []
443
if isinstance(conductor, list):
444
v.append('_conductor')
445
if len(hodge_numbers)>0 and isinstance(hodge_numbers[0], list):
446
v.append('_hodge_numbers')
447
if isinstance(weight, list):
448
v.append('_weight')
449
if len(poles) > 0 and isinstance(poles[0], list):
450
v.append('_poles')
451
if isinstance(epsilon, list):
452
v.append('_epsilon')
453
454
if len(v) > 0:
455
found_params = False
456
for X in cartesian_product_iterator([getattr(self, attr) for attr in v]):
457
kwds = dict((v[i],X[i]) for i in range(len(v)))
458
if self._is_valid_parameters(prec=prec, save=True, **kwds):
459
found_params = True
460
break
461
if not found_params:
462
raise RuntimeError, "no choice of values for %s works"%(', '.join(v))
463
464
def _is_valid_parameters(self, prec=53, save=True, **kwds):
465
valid = False
466
try:
467
old = [(k, getattr(self, k)) for k in kwds.keys()]
468
for k,v in kwds.iteritems():
469
setattr(self, k, v)
470
self._function.clear_cache()
471
self._function(prec=prec)
472
try:
473
self._function(prec=prec)
474
valid = True
475
except RuntimeError:
476
pass
477
finally:
478
if not save:
479
for k, v in old:
480
setattr(self, k, v)
481
return valid
482
483
def __cmp__(self, right):
484
if self is right:
485
return 0
486
for a in ['degree', 'weight', 'conductor', 'epsilon', 'base_field']:
487
c = cmp(getattr(self, a)(), getattr(right, a)())
488
if c: return c
489
c = cmp(type(self), type(right))
490
return self._cmp(right)
491
492
def _cmp(self, right):
493
raise TypeError
494
495
def parent(self):
496
"""
497
Return parent of this L-series, which is the collection of all
498
L-series and their products.
499
500
EXAMPLES::
501
502
sage: from psage.lseries.eulerprod import LSeries
503
sage: L = LSeries('delta'); P = L.parent(); P
504
All L-series objects and their products
505
sage: L in P
506
True
507
"""
508
return LSeriesParent
509
510
def __pow__(self, n):
511
"""
512
Return the n-th power of this L-series, where n can be any integer.
513
514
EXAMPLES::
515
516
sage: from psage.lseries.eulerprod import LSeries
517
sage: L = LSeries('delta');
518
sage: L3 = L^3; L3
519
(L-function associated to Ramanujan's Delta (a weight 12 cusp form))^3
520
sage: L3(1)
521
0.0000524870430366548
522
sage: L(1)^3
523
0.0000524870430366548
524
sage: M = L^(-3); M(1)
525
19052.3211471761
526
sage: L(1)^(-3)
527
19052.3211471761
528
529
Higher precision::
530
531
sage: M = L^(-3); M(RealField(100)(1))
532
19052.321147176093380952680193
533
sage: L(RealField(100)(1))^(-3)
534
19052.321147176093380952680193
535
536
Special case -- 0th power -- is not allowed::
537
538
sage: L^0
539
Traceback (most recent call last):
540
...
541
ValueError: product must be nonempty
542
543
"""
544
return LSeriesProduct([(self, ZZ(n))])
545
546
def __mul__(self, right):
547
"""
548
Multiply two L-series, or an L-series times a formal product of L-series.
549
550
EXAMPLES::
551
552
sage: from psage.lseries.eulerprod import LSeries
553
sage: d = LSeries('delta'); z = LSeries('zeta')
554
sage: d * z
555
(Riemann Zeta function viewed as an L-series) * (L-function associated to Ramanujan's Delta (a weight 12 cusp form))
556
sage: d * (d * z)
557
(Riemann Zeta function viewed as an L-series) * (L-function associated to Ramanujan's Delta (a weight 12 cusp form))^2
558
"""
559
if isinstance(right, LSeriesAbstract):
560
return LSeriesProduct([(self, 1), (right, 1)])
561
elif isinstance(right, LSeriesProduct):
562
return right * self
563
raise TypeError
564
565
def __div__(self, right):
566
"""
567
Divide two L-series or formal L-series products.
568
569
EXAMPLES::
570
571
sage: from psage.lseries.eulerprod import LSeries
572
sage: d = LSeries('delta'); z = LSeries('zeta')
573
sage: d / z
574
(Riemann Zeta function viewed as an L-series)^-1 * (L-function associated to Ramanujan's Delta (a weight 12 cusp form))
575
sage: d / (z^3)
576
(Riemann Zeta function viewed as an L-series)^-3 * (L-function associated to Ramanujan's Delta (a weight 12 cusp form))
577
"""
578
if isinstance(right, LSeriesAbstract):
579
return LSeriesProduct([(self, 1), (right, -1)])
580
elif isinstance(right, LSeriesProduct):
581
return LSeriesProduct(Factorization([(self, 1)]) / right._factorization)
582
raise TypeError
583
584
def conductor(self):
585
"""
586
Return the conductor of this L-series.
587
588
EXAMPLES::
589
590
sage: from psage.lseries.eulerprod import LSeries
591
sage: LSeries('zeta').conductor()
592
1
593
sage: LSeries('delta').conductor()
594
1
595
sage: LSeries(EllipticCurve('11a')).conductor()
596
11
597
sage: LSeries(Newforms(33)[0]).conductor()
598
33
599
sage: LSeries(DirichletGroup(37).0).conductor()
600
37
601
sage: LSeries(kronecker_character(7)).conductor()
602
28
603
sage: kronecker_character(7).conductor()
604
28
605
sage: L = LSeries(EllipticCurve('11a3').base_extend(QQ[sqrt(2)]), prec=5)
606
sage: L.conductor().factor()
607
2^6 * 11^2
608
"""
609
return self._conductor
610
611
def hodge_numbers(self):
612
"""
613
Return the Hodge numbers of this L-series. These define the local Gamma factors.
614
615
EXAMPLES::
616
617
sage: from psage.lseries.eulerprod import LSeries
618
sage: LSeries('zeta').hodge_numbers()
619
[0]
620
sage: LSeries('delta').hodge_numbers()
621
[0, 1]
622
sage: LSeries(EllipticCurve('11a')).hodge_numbers()
623
[0, 1]
624
sage: LSeries(Newforms(43,names='a')[1]).hodge_numbers()
625
[0, 1]
626
sage: LSeries(DirichletGroup(37).0).hodge_numbers()
627
[1]
628
sage: LSeries(EllipticCurve(QQ[sqrt(-1)],[1,2]), prec=5).hodge_numbers() # long time
629
[0, 0, 1, 1]
630
"""
631
return list(self._hodge_numbers)
632
633
def weight(self):
634
"""
635
Return the weight of this L-series.
636
637
EXAMPLES::
638
639
sage: from psage.lseries.eulerprod import LSeries
640
sage: LSeries('zeta').weight()
641
1
642
sage: LSeries('delta').weight()
643
12
644
sage: LSeries(EllipticCurve('389a')).weight()
645
2
646
sage: LSeries(Newforms(43,names='a')[1]).weight()
647
2
648
sage: LSeries(Newforms(6,4)[0]).weight()
649
4
650
sage: LSeries(DirichletGroup(37).0).weight()
651
1
652
sage: L = LSeries(EllipticCurve('11a3').base_extend(QQ[sqrt(2)]),prec=5); L.weight()
653
2
654
"""
655
return self._weight
656
657
def poles(self):
658
"""
659
Poles of the *completed* L-function with the extra Gamma
660
factors included.
661
662
WARNING: These are not just the poles of self.
663
664
OUTPUT:
665
- list of numbers
666
667
EXAMPLES::
668
669
sage: from psage.lseries.eulerprod import LSeries
670
sage: LSeries('zeta').poles()
671
[0, 1]
672
sage: LSeries('delta').poles()
673
[]
674
"""
675
return list(self._poles)
676
677
def residues(self, prec=None):
678
"""
679
Residues of the *completed* L-function at each pole.
680
681
EXAMPLES::
682
683
sage: from psage.lseries.eulerprod import LSeries
684
sage: LSeries('zeta').residues()
685
[9, 8]
686
sage: LSeries('delta').residues()
687
[]
688
sage: v = LSeries('zeta').residues()
689
sage: v.append(10)
690
sage: LSeries('zeta').residues()
691
[9, 8]
692
693
The residues of the Dedekind Zeta function of a field are dynamically computed::
694
695
sage: K.<a> = NumberField(x^2 + 1)
696
sage: L = LSeries(K); L
697
Dedekind Zeta function of Number Field in a with defining polynomial x^2 + 1
698
699
If you just call residues you get back that they are automatically computed::
700
701
sage: L.residues()
702
'automatic'
703
704
But if you call with a specific precision, they are computed using that precision::
705
706
sage: L.residues(prec=53)
707
[-0.886226925452758]
708
sage: L.residues(prec=200)
709
[-0.88622692545275801364908374167057259139877472806119356410690]
710
"""
711
if self._residues == 'automatic':
712
if prec is None:
713
return self._residues
714
else:
715
C = ComplexField(prec)
716
return [C(a) for a in self._function(prec=prec).gp()('Lresidues')]
717
else:
718
return list(self._residues)
719
720
def base_field(self):
721
"""
722
EXAMPLES::
723
724
sage: from psage.lseries.eulerprod import LSeries
725
sage: LSeries('zeta').base_field()
726
Rational Field
727
sage: L = LSeries(EllipticCurve('11a3').base_extend(QQ[sqrt(2)]), prec=5); L.base_field()
728
Number Field in sqrt2 with defining polynomial x^2 - 2
729
"""
730
return self._base_field
731
732
def epsilon(self, prec=None):
733
"""
734
EXAMPLES::
735
736
sage: from psage.lseries.eulerprod import LSeries
737
sage: LSeries('zeta').epsilon()
738
1
739
sage: LSeries('delta').epsilon()
740
1
741
sage: LSeries(EllipticCurve('389a')).epsilon()
742
1
743
sage: LSeries(EllipticCurve('37a')).epsilon()
744
-1
745
sage: LSeries(Newforms(389,names='a')[1]).epsilon()
746
-1
747
sage: LSeries(Newforms(6,4)[0]).epsilon()
748
1
749
sage: LSeries(DirichletGroup(7).0).epsilon()
750
1/7*I*((((((e^(2/21*I*pi) + 1)*e^(2/21*I*pi) + 1)*e^(1/21*I*pi) - 1)*e^(1/21*I*pi) - 1)*e^(2/21*I*pi) - 1)*e^(1/21*I*pi) - 1)*sqrt(7)*e^(1/21*I*pi)
751
sage: LSeries(DirichletGroup(7).0).epsilon(prec=53)
752
0.386513572759156 + 0.922283718859307*I
753
754
In this example, the epsilon factor is computed when the curve
755
is created. The prec parameter determines the floating point precision
756
used in computing the epsilon factor::
757
758
sage: L = LSeries(EllipticCurve('11a3').base_extend(QQ[sqrt(2)]), prec=5); L.epsilon()
759
-1
760
761
Here is extra confirmation that the rank is really odd over the quadratic field::
762
763
sage: EllipticCurve('11a').quadratic_twist(2).rank()
764
1
765
766
We can compute with the L-series too::
767
768
sage: L(RealField(5)(2))
769
0.53
770
771
For L-functions of newforms with nontrivial character, the
772
epsilon factor is harder to find (we don't have a good
773
algorithm implemented to find it) and might not even be 1 or
774
-1, so it is set to 'solve'. In this case, the functional
775
equation is used to determine the solution.::
776
777
sage: f = Newforms(kronecker_character_upside_down(7),3)[0]; f
778
q - 3*q^2 + 5*q^4 + O(q^6)
779
sage: L = LSeries(f)
780
sage: L.epsilon()
781
'solve'
782
sage: L(3/2)
783
0.332981771482934
784
sage: L.epsilon()
785
1
786
787
Here is an example with nontrivial character::
788
789
sage: f = Newforms(DirichletGroup(7).0, 5, names='a')[0]; f
790
q + a0*q^2 + ((zeta6 - 2)*a0 - zeta6 - 1)*q^3 + (-4*zeta6*a0 + 2*zeta6 - 2)*q^4 + ((4*zeta6 - 2)*a0 + 9*zeta6 - 18)*q^5 + O(q^6)
791
sage: L = LSeries(f)
792
793
First trying to evaluate with the default (53 bits) of
794
precision fails, since for some reason (that I do not
795
understand) the program is unable to find a valid epsilon
796
factor::
797
798
sage: L(0)
799
Traceback (most recent call last):
800
...
801
RuntimeError: unable to determine epsilon from functional equation working to precision 53, since we get epsilon=0.806362085925390 - 0.00491051026156292*I, which is not sufficiently close to 1
802
803
However, when we evaluate to 100 bits of precision it works::
804
805
sage: L(RealField(100)(0))
806
0
807
808
The epsilon factor is *not* known to infinite precision::
809
810
sage: L.epsilon()
811
'solve'
812
813
But it is now known to 100 bits of precision, and here it is::
814
815
sage: L.epsilon(100)
816
0.42563106101692403875896879406 - 0.90489678963824790765479396740*I
817
818
When we try to compute to higher precision, again Sage solves for the epsilon factor
819
numerically::
820
821
sage: L(RealField(150)(1))
822
0.26128389551787271923496480408992971337929665 - 0.29870133769674001421149135036267324347896657*I
823
824
And now it is known to 150 bits of precision. Notice that
825
this is consistent with the value found above, and has
826
absolute value (very close to) 1.
827
828
sage: L.epsilon(150)
829
0.42563106101692403875896879406038776338921622 - 0.90489678963824790765479396740501409301704122*I
830
sage: abs(L.epsilon(150))
831
1.0000000000000000000000000000000000000000000
832
"""
833
if self._epsilon == 'solve' or (
834
hasattr(self._epsilon, 'prec') and (prec is None or self._epsilon.prec() < prec)):
835
return 'solve'
836
if prec is not None:
837
C = ComplexField(prec)
838
if isinstance(self._epsilon, list):
839
return [C(x) for x in self._epsilon]
840
return C(self._epsilon)
841
return self._epsilon
842
843
def is_selfdual(self):
844
"""
845
Return True if this L-series is self dual; otherwise, return False.
846
847
EXAMPLES::
848
849
Many L-series are self dual::
850
851
sage: from psage.lseries.eulerprod import LSeries
852
sage: LSeries('zeta').is_selfdual()
853
True
854
sage: LSeries('delta').is_selfdual()
855
True
856
sage: LSeries(Newforms(6,4)[0]).is_selfdual()
857
True
858
859
Nonquadratic characters have non-self dual L-series::
860
861
sage: LSeries(DirichletGroup(7).0).is_selfdual()
862
False
863
sage: LSeries(kronecker_character(7)).is_selfdual()
864
True
865
866
Newforms with non-quadratic characters also have non-self dual L-seris::
867
868
sage: L = LSeries(Newforms(DirichletGroup(7).0, 5, names='a')[0]); L.is_selfdual()
869
False
870
"""
871
return self._is_selfdual
872
873
@cached_method
874
def degree(self):
875
"""
876
Return the degree of this L-function, which is by definition the number of Gamma
877
factors (e.g., the number of Hodge numbers) divided by the degree of the base field.
878
879
EXAMPLES::
880
881
sage: from psage.lseries.eulerprod import LSeries
882
sage: LSeries('zeta').degree()
883
1
884
sage: LSeries(DirichletGroup(5).0).degree()
885
1
886
sage: LSeries(EllipticCurve('11a')).degree()
887
2
888
889
The L-series attached to this modular symbols space of dimension 2 is a product
890
of 2 degree 2 L-series, hence has degree 4::
891
892
sage: M = ModularSymbols(43,2,sign=1).cuspidal_subspace()[1]; M.dimension()
893
2
894
sage: L = LSeries(M); L
895
L-series attached to Modular Symbols subspace of dimension 2 of Modular Symbols space of dimension 4 for Gamma_0(43) of weight 2 with sign 1 over Rational Field
896
sage: L.factor()
897
(L-series of a degree 2 newform of level 43 and weight 2) * (L-series of a degree 2 newform of level 43 and weight 2)
898
sage: L.degree()
899
4
900
901
sage: x = var('x'); K.<a> = NumberField(x^2-x-1); LSeries(EllipticCurve([0,-a,a,0,0])).degree()
902
2
903
"""
904
n = len(self.hodge_numbers())
905
d = self.base_field().degree()
906
assert n % d == 0, "degree of base field must divide the number of Hodge numbers"
907
return n//d
908
909
def twist(self, chi, conductor=None, epsilon=None, prec=53):
910
r"""
911
Return the quadratic twist of this L-series by the character chi, which
912
must be a character of self.base_field(). Thus chi should take as input
913
prime ideals (or primes) of the ring of integers of the base field, and
914
output something that can be coerced to the complex numbers.
915
916
INPUT:
917
- `\chi` -- 1-dimensional character
918
919
EXAMPLES::
920
921
sage: from psage.lseries.eulerprod import LSeries
922
sage: E = EllipticCurve('11a')
923
sage: L = LSeries(E); L
924
L-series of Elliptic Curve defined by y^2 + y = x^3 - x^2 - 10*x - 20 over Rational Field
925
sage: L3 = L.twist(DirichletGroup(3).0); L3
926
Twist of L-series of Elliptic Curve defined by y^2 + y = x^3 - x^2 - 10*x - 20 over Rational Field by Dirichlet character modulo 3 of conductor 3 mapping 2 |--> -1
927
sage: L3._chi
928
Dirichlet character modulo 3 of conductor 3 mapping 2 |--> -1
929
sage: L3._L
930
L-series of Elliptic Curve defined by y^2 + y = x^3 - x^2 - 10*x - 20 over Rational Field
931
sage: L3(1)
932
1.68449633297548
933
sage: F = E.quadratic_twist(-3)
934
sage: L3.conductor()
935
99
936
sage: F.conductor()
937
99
938
sage: F.lseries()(1)
939
1.68449633297548
940
sage: L3.anlist(20)
941
[0, 1, 2, 0, 2, -1, 0, -2, 0, 0, -2, -1, 0, 4, -4, 0, -4, 2, 0, 0, -2]
942
sage: F.anlist(20)
943
[0, 1, 2, 0, 2, -1, 0, -2, 0, 0, -2, -1, 0, 4, -4, 0, -4, 2, 0, 0, -2]
944
sage: L3.anlist(1000) == F.anlist(1000)
945
True
946
sage: L3.local_factor(11)
947
T + 1
948
949
A higher degree twist::
950
951
sage: L = LSeries(EllipticCurve('11a'))
952
sage: L5 = L.twist(DirichletGroup(5).0); L5
953
Twist of L-series of Elliptic Curve defined by y^2 + y = x^3 - x^2 - 10*x - 20 over Rational Field by Dirichlet character modulo 5 of conductor 5 mapping 2 |--> zeta4
954
sage: L5(1)
955
1.28009593569230 - 0.681843202124309*I
956
sage: L5.epsilon()
957
'solve'
958
sage: L5.epsilon(53)
959
0.989989082587826 + 0.141143956147310*I
960
sage: L5.conductor()
961
275
962
sage: L5.taylor_series(center=1, degree=3)
963
1.28009593569230 - 0.681843202124309*I + (-0.536450338806282 + 0.166075270978779*I)*z + (0.123743053129226 + 0.320802890011298*I)*z^2 + O(z^3)
964
965
966
WARNING!! Twisting is not implemented in full generality when
967
the conductors are not coprime. One case where we run into
968
trouble is when twisting lowers the level of a newform. Below
969
we take the form of level 11 and weight 2, twist it by the
970
character chi of conductor 3 to get a form of level 99. Then
971
we take the L-series of the level 99 form, and twist that by
972
chi, which should be the L-series attached to the form of
973
level 11. Unfortunately, our code for working out the local
974
L-factors doesn't succeed in this case, hence the local factor
975
is wrong, so the functional equation is not satisfied.
976
977
sage: f = Newform('11a'); f
978
q - 2*q^2 - q^3 + 2*q^4 + q^5 + O(q^6)
979
sage: L = LSeries(f)
980
sage: chi = DirichletGroup(3).0
981
sage: Lc = L.twist(chi); Lc
982
Twist of L-series of a degree 1 newform of level 11 and weight 2 by Dirichlet character modulo 3 of conductor 3 mapping 2 |--> -1
983
sage: Lc.anlist(20)
984
[0, 1, 2, 0, 2, -1, 0, -2, 0, 0, -2, -1, 0, 4, -4, 0, -4, 2, 0, 0, -2]
985
sage: g = Newform('99d'); g
986
q + 2*q^2 + 2*q^4 - q^5 + O(q^6)
987
sage: list(g.qexp(20))
988
[0, 1, 2, 0, 2, -1, 0, -2, 0, 0, -2, -1, 0, 4, -4, 0, -4, 2]
989
sage: Lt = Lc.twist(chi, conductor=11)
990
Traceback (most recent call last):
991
...
992
RuntimeError: no choice of values for _epsilon works
993
sage: Lt = Lc.twist(chi,conductor=11, epsilon=1)
994
sage: Lt(1)
995
Traceback (most recent call last):
996
...
997
RuntimeError: invalid L-series parameters: functional equation not satisfied
998
999
This is because the local factor is wrong::
1000
1001
sage: Lt.local_factor(3)
1002
1
1003
sage: L.local_factor(3)
1004
3*T^2 + T + 1
1005
1006
1007
"""
1008
return LSeriesTwist(self, chi=chi, conductor=conductor, epsilon=epsilon, prec=prec)
1009
1010
@cached_method
1011
def local_factor(self, P, prec=None):
1012
"""
1013
Return the local factor of the L-function at the prime P of
1014
self._base_field. The result is cached.
1015
1016
INPUT:
1017
- a prime P of the ring of integers of the base_field
1018
- prec -- None or positive integer (bits of precision)
1019
1020
OUTPUT:
1021
- a polynomial, e.g., something like "1-a*T+p*T^2".
1022
1023
EXAMPLES:
1024
1025
You must overload this in the derived class::
1026
1027
sage: from psage.lseries.eulerprod import LSeriesAbstract
1028
sage: L = LSeriesAbstract(conductor=1, hodge_numbers=[0], weight=1, epsilon=1, poles=[1], residues=[-1], base_field=QQ)
1029
sage: L.local_factor(2)
1030
Traceback (most recent call last):
1031
...
1032
NotImplementedError: must be implemented in the derived class
1033
"""
1034
return self._local_factor(P, prec=prec)
1035
1036
def _local_factor(self, P, prec=None):
1037
"""
1038
Compute local factor at prime P. This must be overwritten in the derived class.
1039
1040
EXAMPLES::
1041
1042
sage: from psage.lseries.eulerprod import LSeriesAbstract
1043
sage: L = LSeriesAbstract(conductor=1, hodge_numbers=[0], weight=1, epsilon=1, poles=[1], residues=[-1], base_field=QQ)
1044
sage: L.local_factor(2)
1045
Traceback (most recent call last):
1046
...
1047
NotImplementedError: must be implemented in the derived class
1048
"""
1049
raise NotImplementedError, "must be implemented in the derived class"
1050
1051
def __repr__(self):
1052
"""
1053
EXAMPLES::
1054
1055
sage: from psage.lseries.eulerprod import LSeriesAbstract
1056
sage: L = LSeriesAbstract(conductor=1, hodge_numbers=[0], weight=1, epsilon=1, poles=[1], residues=[-1], base_field=QQ)
1057
sage: L.__repr__()
1058
'Euler Product L-series with conductor 1, Hodge numbers [0], weight 1, epsilon 1, poles [1], residues [-1] over Rational Field'
1059
"""
1060
return "Euler Product L-series with conductor %s, Hodge numbers %s, weight %s, epsilon %s, poles %s, residues %s over %s"%(
1061
self._conductor, self._hodge_numbers, self._weight, self.epsilon(), self._poles, self._residues, self._base_field)
1062
1063
def _precompute_local_factors(self, bound, prec):
1064
"""
1065
Derived classes may use this as a 'hint' that _local_factors
1066
will soon get called for primes of norm less than the bound.
1067
1068
In the base class, this is a no-op, and it is not necessary to
1069
overload this class.
1070
1071
INPUT:
1072
- ``bound`` -- integer
1073
- ``prec`` -- integer
1074
1075
EXAMPLES::
1076
1077
sage: from psage.lseries.eulerprod import LSeriesAbstract
1078
sage: L = LSeriesAbstract(conductor=1, hodge_numbers=[0], weight=1, epsilon=1, poles=[1], residues=[-1], base_field=QQ)
1079
sage: L._precompute_local_factors(100, 53)
1080
"""
1081
pass
1082
1083
def _primes_above(self, p):
1084
"""
1085
Return the primes of the ring of integers of the base field above the integer p.
1086
1087
INPUT:
1088
- p -- prime integer (no type checking necessarily done)
1089
1090
EXAMPLES::
1091
1092
sage: from psage.lseries.eulerprod import LSeriesAbstract
1093
sage: L = LSeriesAbstract(conductor=1, hodge_numbers=[0], weight=1, epsilon=1, poles=[1], residues=[-1], base_field=QQ)
1094
sage: L._primes_above(3)
1095
[3]
1096
sage: L = LSeriesAbstract(conductor=1, hodge_numbers=[0], weight=1, epsilon=1, poles=[1], residues=[-1], base_field=QQ[sqrt(-1)])
1097
sage: L._primes_above(5)
1098
[Fractional ideal (I + 2), Fractional ideal (-I + 2)]
1099
sage: L._primes_above(3)
1100
[Fractional ideal (3)]
1101
"""
1102
K = self._base_field
1103
if is_RationalField(K):
1104
return [p]
1105
else:
1106
return K.primes_above(p)
1107
1108
def zeros(self, n):
1109
"""
1110
Return the imaginary parts of the first n nontrivial zeros of
1111
this L-function on the critical line in the upper half plane,
1112
as 32-bit real numbers.
1113
1114
INPUT:
1115
- n -- nonnegative integer
1116
1117
EXAMPLES::
1118
1119
"""
1120
return self._lcalc().zeros(n)
1121
1122
def _lcalc(self):
1123
"""
1124
Return Rubinstein Lcalc object attached to this L-series. This
1125
is useful both for evaluating the L-series, especially when
1126
the imaginary part is large, and for computing the zeros in
1127
the critical strip.
1128
1129
EXAMPLES::
1130
"""
1131
# this will require using Rubinstein's L-calc
1132
raise NotImplementedError
1133
1134
def anlist(self, bound, prec=None):
1135
"""
1136
Return list `v` of Dirichlet series coefficients `a_n`for `n` up to and including bound,
1137
where `v[n] = a_n`. If at least one of the coefficients is ambiguous, e.g., the
1138
local_factor method returns a list of possibilities, then this function instead
1139
returns a generator over all possible `a_n` lists.
1140
1141
In particular, we include v[0]=0 as a convenient place holder
1142
to avoid having `v[n-1] = a_n`.
1143
1144
INPUT:
1145
- ``bound`` -- nonnegative integer
1146
1147
EXAMPLES::
1148
1149
sage: from psage.lseries.eulerprod import LSeries; L = LSeries('zeta')
1150
sage: L.anlist(30)
1151
[0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
1152
sage: from psage.lseries.eulerprod import LSeries; L = LSeries(EllipticCurve('11a'))
1153
sage: L.anlist(30)
1154
[0, 1, -2, -1, 2, 1, 2, -2, 0, -2, -2, 1, -2, 4, 4, -1, -4, -2, 4, 0, 2, 2, -2, -1, 0, -4, -8, 5, -4, 0, 2]
1155
sage: K.<a> = NumberField(x^2-x-1); L = LSeries(EllipticCurve([0,-a,a,0,0]))
1156
sage: L.anlist(30)
1157
[0, 1, 0, 0, -2, -1, 0, 0, 0, -4, 0, 3, 0, 0, 0, 0, 0, 0, 0, 5, 2, 0, 0, 0, 0, -4, 0, 0, 0, 11, 0]
1158
"""
1159
# First check if we know anlist to infinite bit precision up to given bound:
1160
if len(self._anlist[None]) > bound:
1161
if prec is None:
1162
# request it to infinite precision
1163
return self._anlist[None][:bound+1]
1164
else:
1165
# finite precision request
1166
C = ComplexField(prec)
1167
return [C(a) for a in self._anlist[None]]
1168
1169
if prec is not None:
1170
# check if numerically computed already to at least this precision
1171
t = [z for z in self._anlist.iteritems() if z[0] >= prec and len(z[1]) > bound]
1172
if len(t) > 0:
1173
C = ComplexField(prec)
1174
return [C(a) for a in t[0][1][:bound+1]]
1175
1176
self._precompute_local_factors(bound+1, prec=prec)
1177
compute_anlist_multiple = False
1178
LF = []
1179
for p in prime_range(bound+1):
1180
lf = []
1181
for P in self._primes_above(p):
1182
if norm(P) <= bound:
1183
F = self._local_factor(P, prec)
1184
if isinstance(F, list):
1185
compute_anlist_multiple = True
1186
lf.append(F)
1187
LF.append((p, lf))
1188
1189
coefficients = self._compute_anlist(LF, bound, prec)
1190
if not compute_anlist_multiple:
1191
coefficients = list(coefficients)[0]
1192
# save answer in cache
1193
self._anlist[prec] = coefficients
1194
return coefficients
1195
1196
def _compute_anlist(self, LF, bound, prec):
1197
"""
1198
Iterator over possible anlists, given LF, bound, and prec.
1199
1200
INPUT:
1201
- ``LF`` -- list of pairs (p, [local factors (or lists of them) at primes over p])
1202
- ``bound`` -- positive integer
1203
- ``prec`` -- positive integer (bits of precision)
1204
"""
1205
K = self._base_field
1206
coefficients = [0,1] + [0]*(bound-1)
1207
1208
for i, (p, v) in enumerate(LF):
1209
if len(v) > 0:
1210
# Check for the possibility of several different
1211
# choices of Euler factor at a given prime. If this happens,
1212
# we switch gears.
1213
some_list = False
1214
for j, lf in enumerate(v):
1215
if isinstance(lf, list):
1216
some_list = True
1217
for f in list(lf):
1218
LF0 = copy.deepcopy(LF)
1219
LF0[i][1][j] = f
1220
for z in self._compute_anlist(LF0, bound, prec):
1221
yield z
1222
if some_list:
1223
return
1224
# Not several factors -- just compute the a_{p^r} up to the required bound:
1225
f = prod(v)
1226
accuracy_p = int(math.floor(math.log(bound)/math.log(p))) + 1
1227
T = f.parent().gen()
1228
series_p = (f + O(T**accuracy_p))**(-1)
1229
for j in range(1, accuracy_p):
1230
coefficients[p**j] = series_p[j]
1231
1232
# fill in non-prime power coefficients
1233
extend_multiplicatively_generic(coefficients)
1234
yield list(coefficients)
1235
1236
def _symbolic_(self, R, bound=10, prec=None):
1237
"""
1238
Convert self into the symbolic ring as a truncated Dirichleter series, including
1239
terms up to `n^s` where n=bound.
1240
1241
EXAMPLES::
1242
1243
sage: from psage.lseries.eulerprod import LSeries
1244
sage: L = LSeries(EllipticCurve('37a'))
1245
sage: SR(L)
1246
-2/2^s - 3/3^s + 2/4^s - 2/5^s + 6/6^s - 1/7^s + 6/9^s + 4/10^s + 1
1247
sage: L._symbolic_(SR, 20)
1248
-2/2^s - 3/3^s + 2/4^s - 2/5^s + 6/6^s - 1/7^s + 6/9^s + 4/10^s - 5/11^s - 6/12^s - 2/13^s + 2/14^s + 6/15^s - 4/16^s - 12/18^s - 4/20^s + 1
1249
"""
1250
s = R.var('s')
1251
a = self.anlist(bound, prec)
1252
return sum(a[n]/n**s for n in range(1,bound+1))
1253
1254
def __call__(self, s):
1255
"""
1256
Return value of this L-function at s. If s is a real or
1257
complex number to prec bits of precision, then the result is
1258
also computed to prec bits of precision. If s has infinite or
1259
unknown precision, then the L-value is computed to 53 bits of
1260
precision.
1261
1262
EXAMPLES::
1263
1264
sage: from psage.lseries.eulerprod import LSeries
1265
sage: L = LSeries(EllipticCurve('37a'))
1266
sage: L(1)
1267
0
1268
sage: L(2)
1269
0.381575408260711
1270
sage: z = L(RealField(100)(2)); z
1271
0.38157540826071121129371040958
1272
sage: z.prec()
1273
100
1274
sage: L(RealField(150)(2))
1275
0.38157540826071121129371040958008663667709753
1276
1277
WARNING: There will be precision loss (with a warning) if the imaginary
1278
part of the input is large::
1279
1280
sage: L = LSeries('zeta')
1281
sage: L(1/2 + 40*I)
1282
verbose -1 (...: dokchitser.py, __call__) Warning: Loss of 14 decimal digits due to cancellation
1283
0.793046013671137 - 1.04127377821427*I
1284
sage: L(ComplexField(200)(1/2 + 40*I))
1285
verbose -1 (...: dokchitser.py, __call__) Warning: Loss of 14 decimal digits due to cancellation
1286
0.79304495256192867196489258889793696080572220439302833315881 - 1.0412746146510650200518905953910554313275550685861559488384*I
1287
1288
An example with a pole::
1289
1290
sage: L = LSeries('zeta')
1291
sage: L(2) # good
1292
1.64493406684823
1293
sage: L(1) # a pole!
1294
Traceback (most recent call last):
1295
...
1296
ZeroDivisionError: pole at 1
1297
"""
1298
if s in self._poles:
1299
raise ZeroDivisionError, "pole at %s"%s
1300
return self._function(prec(s))(s)
1301
1302
def derivative(self, k=1):
1303
"""
1304
Return the k-th derivative of self.
1305
1306
INPUT:
1307
- k -- (default: 1) nonnegative integer
1308
1309
EXAMPLES::
1310
1311
sage: from psage.lseries.eulerprod import LSeries
1312
sage: L = LSeries('zeta')
1313
sage: L.derivative()
1314
First derivative of Riemann Zeta function viewed as an L-series
1315
1316
We numerically approximate the derivative at two points and compare
1317
with evaluating the derivative object::
1318
1319
sage: eps=1e-10; (L(2+eps) - L(2))/eps
1320
-0.937547817159157
1321
sage: Lp = L.derivative(); Lp(2)
1322
-0.937548254315844
1323
sage: eps=1e-10; (L(2+I+eps) - L(2+I))/eps
1324
0.0624900131640516 + 0.489033813444451*I
1325
sage: Lp(2+I)
1326
0.0624900021906470 + 0.489033591679899*I
1327
1328
Higher derivatives::
1329
1330
sage: L.derivative(2)
1331
Second derivative of Riemann Zeta function viewed as an L-series
1332
sage: L.derivative(2)(2)
1333
1.98928023429890
1334
sage: L.derivative(3)
1335
Third derivative of Riemann Zeta function viewed as an L-series
1336
sage: L.derivative(4)
1337
4-th derivative of Riemann Zeta function viewed as an L-series
1338
sage: L.derivative(5)
1339
5-th derivative of Riemann Zeta function viewed as an L-series
1340
1341
Derivative of derivative::
1342
1343
sage: L.derivative().derivative()
1344
Second derivative of Riemann Zeta function viewed as an L-series
1345
1346
Using the derivative function in Sage works::
1347
1348
sage: derivative(L)
1349
First derivative of Riemann Zeta function viewed as an L-series
1350
sage: derivative(L,2)
1351
Second derivative of Riemann Zeta function viewed as an L-series
1352
"""
1353
if k == 0:
1354
return self
1355
return LSeriesDerivative(self, k)
1356
1357
def taylor_series(self, center=None, degree=6, variable='z', prec=53):
1358
"""
1359
Return the Taylor series expansion of self about the given
1360
center point to the given degree in the specified variable
1361
numerically computed to the precision prec in bits. If the
1362
center is not specified it defaults to weight / 2.
1363
1364
INPUT:
1365
- ``center`` -- None or number that coerces to the complex numbers
1366
- ``degree`` -- integer
1367
- ``variable`` -- string or symbolic variable
1368
- ``prec`` -- positive integer (floating point bits of precision)
1369
1370
EXAMPLES:::
1371
1372
sage: from psage.lseries.eulerprod import LSeries; L = LSeries('zeta')
1373
sage: L.taylor_series()
1374
-1.46035450880959 - 3.92264613920915*z - 8.00417850696433*z^2 - 16.0005515408865*z^3 - 31.9998883216853*z^4 - 64.0000050055172*z^5 + O(z^6)
1375
sage: RealField(53)(zeta(1/2))
1376
-1.46035450880959
1377
sage: L.taylor_series(center=2, degree=4, variable='t', prec=30)
1378
1.6449341 - 0.93754825*t + 0.99464012*t^2 - 1.0000243*t^3 + O(t^4)
1379
sage: RealField(30)(zeta(2))
1380
1.6449341
1381
1382
"""
1383
if center is None:
1384
center = ComplexField(prec)(self._weight)/2
1385
return self._function(prec).taylor_series(center, degree, variable)
1386
1387
def analytic_rank(self, tiny=1e-8, prec=53):
1388
center = ComplexField(prec)(self._weight) / 2
1389
degree = 4
1390
while True:
1391
f = self.taylor_series(center, degree, prec=prec)
1392
i = 0
1393
while i < degree:
1394
if abs(f[i]) > tiny:
1395
return i
1396
i += 1
1397
degree += 2
1398
1399
@cached_method
1400
def _function(self, prec=53, T=1.2):
1401
"""
1402
Return Dokchitser object that allows for computation of this
1403
L-series computed to enough terms so that the functional
1404
equation checks out with the given value of T and precision.
1405
1406
This is used behind the scenes for evaluation and computation
1407
of Taylor series.
1408
"""
1409
eps = self.epsilon(prec)
1410
return self._dokchitser(prec, eps, T=T)
1411
1412
def _dokchitser_unitialized(self, prec, epsilon):
1413
# Create the Dokchitser object
1414
if epsilon == 'solve':
1415
eps = 'X'
1416
else:
1417
eps = epsilon
1418
return Dokchitser(conductor = self.conductor(), gammaV = self.hodge_numbers(), weight = self.weight(),
1419
eps = eps, poles = self.poles(), residues = self.residues(),
1420
prec = prec)
1421
1422
def number_of_coefficients(self, prec=53, T=1.2):
1423
"""
1424
Return the number of Dirichlet series coefficients that will
1425
be needed in order to evaluate this L-series (near the real
1426
line) to prec bits of precision and have the functional
1427
equation test pass with the given value of T.
1428
1429
INPUT:
1430
- prec -- integer
1431
1432
EXAMPLES::
1433
1434
sage: from psage.lseries.eulerprod import LSeries
1435
sage: L = LSeries(DirichletGroup(5).0)
1436
sage: L.number_of_coefficients(20)
1437
8
1438
sage: L.number_of_coefficients()
1439
11
1440
sage: L.number_of_coefficients(1000)
1441
43
1442
"""
1443
return self._dokchitser_unitialized(prec, self.epsilon(prec)).num_coeffs(T)
1444
1445
def _dokchitser(self, prec, epsilon, T=1.2):
1446
L = self._dokchitser_unitialized(prec, epsilon)
1447
1448
# Find out how many coefficients of the Dirichlet series are needed
1449
# to compute to the requested precision.
1450
n = L.num_coeffs(T=T)
1451
#if n >= 500: # TODO: for debugging only -- remove later
1452
# print "num coeffs =", n
1453
1454
# Compute the Dirichlet series coefficients
1455
X = self.anlist(n, prec)
1456
if isinstance(X, types.GeneratorType):
1457
# Several possible coefficients -- we try them until finding on that works.
1458
coeff_lists = X
1459
else:
1460
# Only one to try
1461
coeff_lists = [X]
1462
1463
tiny0 = tiny(prec)
1464
for coeffs in coeff_lists:
1465
# Define a string that when evaluated in PARI defines a function
1466
# a(k), which returns the Dirichlet coefficient a_k.
1467
s = 'v=[%s]; a(k)=v[k];'%','.join([str(z) if isinstance(z, (int,long,Integer)) else z._pari_init_() for z in coeffs[1:]])
1468
1469
# Tell the L-series / PARI about the coefficients.
1470
1471
if self.is_selfdual():
1472
L.init_coeffs('a(k)', pari_precode = s)
1473
else:
1474
# Have to directly call gp_eval, since case of functional equation having two different
1475
# (conjugate) L-functions isn't supported in Dokchitser class (yet).
1476
L._Dokchitser__init = True
1477
L._gp_eval(s)
1478
L._gp_eval('initLdata("a(k)",1,"conj(a(k))")')
1479
1480
if epsilon == 'solve':
1481
cmd = "sgneq = Vec(checkfeq()); sgn = -sgneq[2]/sgneq[1]; sgn"
1482
epsilon = ComplexField(prec)(L._gp_eval(cmd))
1483
if abs(abs(epsilon)-1) > tiny0:
1484
raise RuntimeError, "unable to determine epsilon from functional equation working to precision %s, since we get epsilon=%s, which is not sufficiently close to 1"%(prec, epsilon)
1485
# 1, -1 are two common special cases, where it is clear what the
1486
# infinite precision version is.
1487
if epsilon == 1:
1488
self._epsilon = 1
1489
elif epsilon == -1:
1490
self._epsilon = -1
1491
else:
1492
self._epsilon = epsilon
1493
1494
fe = L.check_functional_equation()
1495
if abs(fe) <= tiny0:
1496
# one worked!
1497
self._anlist[prec] = coeffs
1498
return L
1499
else:
1500
pass
1501
1502
# They all failed.
1503
raise RuntimeError, "invalid L-series parameters: functional equation not satisfied"
1504
1505
def check_functional_equation(self, T, prec=53):
1506
return self._function(prec=prec).check_functional_equation(T)
1507
1508
class LSeriesProductEvaluator(object):
1509
def __init__(self, factorization, prec):
1510
self._factorization = factorization
1511
self._prec = prec
1512
1513
def __call__(self, s):
1514
try:
1515
v = self._functions
1516
except AttributeError:
1517
self._functions = [(L._function(self._prec),e) for L,e in self._factorization]
1518
v = self._functions
1519
return prod(f(s)**e for f,e in v)
1520
1521
class LSeriesProduct(object):
1522
"""
1523
A formal product of L-series.
1524
"""
1525
def __init__(self, F):
1526
"""
1527
INPUT:
1528
- `F` -- list of pairs (L,e) where L is an L-function and e is a nonzero integer.
1529
"""
1530
if not isinstance(F, Factorization):
1531
F = Factorization(F)
1532
F.sort(cmp)
1533
if len(F) == 0:
1534
raise ValueError, "product must be nonempty"
1535
self._factorization = F
1536
1537
def __cmp__(self, right):
1538
# TODO: make work even if right not a product
1539
return cmp(self._factorization, right._factorization)
1540
1541
def is_selfdual(self):
1542
"""
1543
Return True if every factor of self is self dual; otherwise, return False.
1544
1545
EXAMPLES::
1546
1547
sage: from psage.lseries.eulerprod import LSeries
1548
sage: L1 = LSeries('zeta'); L1.is_selfdual()
1549
True
1550
sage: L2 = LSeries(DirichletGroup(9).0); L2.is_selfdual()
1551
False
1552
sage: (L1*L1).is_selfdual()
1553
True
1554
sage: (L1*L2).is_selfdual()
1555
False
1556
sage: (L2*L2).is_selfdual()
1557
False
1558
"""
1559
is_selfdual = set(L.is_selfdual() for L,e in self._factorization)
1560
if len(is_selfdual) > 1:
1561
return False
1562
else:
1563
return list(is_selfdual)[0]
1564
1565
def conductor(self):
1566
"""
1567
Return the conductor of this product, which we define to be
1568
the product with multiplicities of the conductors of the
1569
factors of self.
1570
1571
EXAMPLES::
1572
1573
sage: from psage.lseries.eulerprod import LSeries
1574
sage: J = J0(33); L = LSeries(J)
1575
sage: L.conductor()
1576
3993
1577
sage: L.conductor().factor()
1578
3 * 11^3
1579
sage: J.decomposition()
1580
[
1581
Simple abelian subvariety 11a(1,33) of dimension 1 of J0(33),
1582
Simple abelian subvariety 11a(3,33) of dimension 1 of J0(33),
1583
Simple abelian subvariety 33a(1,33) of dimension 1 of J0(33)
1584
]
1585
1586
Of course, the conductor as we have defined it need not be an integer::
1587
1588
sage: L = LSeries(J[0])/LSeries(J[2])^2; L
1589
(L-series of a degree 1 newform of level 11 and weight 2) * (L-series of a degree 1 newform of level 33 and weight 2)^-2
1590
sage: L.conductor()
1591
1/99
1592
"""
1593
return prod(L.conductor()**e for L, e in self._factorization)
1594
1595
def __pow__(self, n):
1596
"""
1597
Return the n-th power of this formal L-series product, where n can be any integer.
1598
1599
EXAMPLES::
1600
1601
sage: from psage.lseries.eulerprod import LSeries
1602
sage: L = LSeries('zeta') * LSeries(DirichletGroup(9).0)
1603
sage: L(2)
1604
1.80817853715812 + 0.369298119218816*I
1605
sage: L = LSeries('zeta') * LSeries(DirichletGroup(9).0)
1606
sage: L3 = L^3; L3
1607
(Riemann Zeta function viewed as an L-series)^3 * (L-series attached to Dirichlet character modulo 9 of conductor 9 mapping 2 |--> zeta6)^3
1608
sage: L3(2)
1609
5.17205298762567 + 3.57190597873829*I
1610
sage: CC((zeta(2)*LSeries(DirichletGroup(9).0)(2))^3)
1611
5.17205298762567 + 3.57190597873829*I
1612
sage: L^(-1999)
1613
(Riemann Zeta function viewed as an L-series)^-1999 * (L-series attached to Dirichlet character modulo 9 of conductor 9 mapping 2 |--> zeta6)^-1999
1614
sage: (L^(-1999)) (2)
1615
8.90248311986228e-533 - 6.20123089437732e-533*I
1616
"""
1617
return LSeriesProduct(self._factorization**ZZ(n))
1618
1619
def __mul__(self, right):
1620
if isinstance(right, LSeriesAbstract):
1621
return LSeriesProduct(self._factorization * Factorization([(right,1)]))
1622
elif isinstance(right, LSeriesProduct):
1623
return LSeriesProduct(self._factorization * right._factorization)
1624
else:
1625
raise TypeError
1626
1627
def __div__(self, right):
1628
if isinstance(right, LSeriesAbstract):
1629
return LSeriesProduct(self._factorization * Factorization([(right,-1)]))
1630
elif isinstance(right, LSeriesProduct):
1631
return LSeriesProduct(self._factorization / right._factorization)
1632
raise TypeError
1633
1634
def parent(self):
1635
return LSeriesParent
1636
1637
def factor(self):
1638
return self._factorization
1639
1640
def hodge_numbers(self):
1641
"""
1642
Return the Hodge numbers of this product of L-series.
1643
"""
1644
1645
def degree(self):
1646
return sum(e*L.degree() for L,e in self._factorization)
1647
1648
def local_factor(self, P, prec=None):
1649
return prod(L.local_factor(P,prec)**e for L,e in self._factorization)
1650
1651
def __getitem__(self, *args, **kwds):
1652
return self._factorization.__getitem__(*args, **kwds)
1653
1654
def __repr__(self):
1655
return self.factor().__repr__()
1656
1657
def __call__(self, s):
1658
return self._function(prec(s))(s)
1659
1660
def derivative(self, k=1):
1661
raise NotImplementedError
1662
1663
def taylor_series(self, center=None, degree=6, variable='z', prec=53):
1664
"""
1665
EXAMPLE::
1666
1667
sage: from psage.lseries.eulerprod import LSeries
1668
sage: L1 = LSeries('zeta'); L2 = LSeries('delta')
1669
sage: f = L1 * L2; f
1670
(Riemann Zeta function viewed as an L-series) * (L-function associated to Ramanujan's Delta (a weight 12 cusp form))
1671
sage: f.taylor_series(center=2, degree=4, variable='w', prec=30)
1672
0.24077647 + 0.10066485*w + 0.061553731*w^2 - 0.041923238*w^3 + O(w^4)
1673
sage: L1.taylor_series(2, 4, 'w', 30) * L2.taylor_series(2, 4, 'w', 30)
1674
0.24077647 + 0.10066485*w + 0.061553731*w^2 - 0.041923238*w^3 + O(w^4)
1675
sage: f = L1 / L2; f
1676
(Riemann Zeta function viewed as an L-series) * (L-function associated to Ramanujan's Delta (a weight 12 cusp form))^-1
1677
sage: f.taylor_series(center=2, degree=4, variable='w', prec=30)
1678
11.237843 - 17.508629*w + 21.688182*w^2 - 24.044641*w^3 + O(w^4)
1679
sage: L1.taylor_series(2, 4, 'w', 30) / L2.taylor_series(2, 4, 'w', 30)
1680
11.237843 - 17.508629*w + 21.688182*w^2 - 24.044641*w^3 + O(w^4)
1681
1682
"""
1683
return prod(L.taylor_series(center, degree, variable, prec)**e for L,e in self._factorization)
1684
1685
def analytic_rank(self, prec=53):
1686
"""
1687
Return sum of the order of vanishing counted with
1688
multiplicities of each factors at their center point.
1689
1690
WARNING: The analytic rank is computed numerically, so is
1691
definitely not provably correct.
1692
1693
EXAMPLES::
1694
1695
We compute the analytic rank of the non-simple 32-dimensional modular abelian variety `J_0(389)`.
1696
1697
sage: from psage.lseries.eulerprod import LSeries
1698
sage: M = ModularSymbols(389,sign=1).cuspidal_subspace()
1699
sage: L = LSeries(M); L
1700
L-series attached to Modular Symbols subspace of dimension 32 of Modular Symbols space of dimension 33 for Gamma_0(389) of weight 2 with sign 1 over Rational Field
1701
1702
We first attempt computation of the analytic rank with the default of 53 bits precision::
1703
1704
sage: L.analytic_rank()
1705
Traceback (most recent call last):
1706
...
1707
RuntimeError: invalid L-series parameters: functional equation not satisfied
1708
1709
The above failed because trying to compute one of the degree
1710
20 newforms resulting in some internal error when double
1711
checking the functional equation. So we try with slightly more precision::
1712
1713
sage: L.analytic_rank(70)
1714
13
1715
1716
This works, since the factors have dimensions 1,2,3,6,20, and
1717
the one of degree 1 has rank 2, the ones of degree 2,3,6 have
1718
rank 2,3,6, respectively, and the one of degree 20 has rank 0::
1719
1720
sage: 2*1 + 2 + 3 + 6
1721
13
1722
"""
1723
return sum(e*L.analytic_rank(prec=prec) for L,e in self._factorization)
1724
1725
def weight(self):
1726
"""
1727
Return the weight of this L-series, which is the sum of the weights
1728
of the factors counted with multiplicity.
1729
"""
1730
return sum(e*L.weight() for L,e in self._factorization)
1731
1732
@cached_method
1733
def _function(self, prec=53):
1734
"""
1735
Return Dokchitser object that allows for computation of this
1736
L-series. This is used behind the scenes for evaluation and
1737
computation of Taylor series.
1738
"""
1739
return LSeriesProductEvaluator(self._factorization, prec)
1740
1741
1742
class LSeriesZeta(LSeriesAbstract):
1743
"""
1744
EXAMPLES::
1745
1746
sage: from psage.lseries.eulerprod import LSeries
1747
sage: L = LSeries('zeta'); L
1748
Riemann Zeta function viewed as an L-series
1749
sage: L(2)
1750
1.64493406684823
1751
sage: L(2.000000000000000000000000000000000000000000000000000)
1752
1.64493406684822643647241516664602518921894990120680
1753
sage: zeta(2.000000000000000000000000000000000000000000000000000)
1754
1.64493406684822643647241516664602518921894990120680
1755
sage: L.local_factor(3)
1756
-T + 1
1757
sage: L.local_factor(5)
1758
-T + 1
1759
sage: L.anlist(30)
1760
[0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
1761
"""
1762
def __init__(self):
1763
LSeriesAbstract.__init__(self, conductor=1, hodge_numbers=[0], weight=1, epsilon=1,
1764
poles=[0,1], residues=[9,8], base_field=QQ)
1765
T = ZZ['T'].gen()
1766
self._lf = 1 - T
1767
1768
def _cmp(self, right):
1769
return 0
1770
1771
def _local_factor(self, P, prec):
1772
return self._lf
1773
1774
def __repr__(self):
1775
return "Riemann Zeta function viewed as an L-series"
1776
1777
class LSeriesDelta(LSeriesAbstract):
1778
"""
1779
EXAMPLES::
1780
1781
sage: from psage.lseries.eulerprod import LSeriesDelta; L = LSeriesDelta()
1782
sage: L.anlist(10)
1783
[0, 1, -24, 252, -1472, 4830, -6048, -16744, 84480, -113643, -115920]
1784
sage: list(delta_qexp(11))
1785
[0, 1, -24, 252, -1472, 4830, -6048, -16744, 84480, -113643, -115920]
1786
sage: L.anlist(10^4) == list(delta_qexp(10^4+1))
1787
True
1788
"""
1789
def __init__(self):
1790
LSeriesAbstract.__init__(self, conductor=1, hodge_numbers=[0,1], weight=12, epsilon=1,
1791
poles=[], residues=[], base_field=QQ)
1792
self._T = ZZ['T'].gen()
1793
self._lf = {}
1794
1795
def _local_factor(self, P, prec):
1796
"""
1797
EXAMPLES::
1798
1799
sage: from psage.lseries.eulerprod import LSeriesDelta; L = LSeriesDelta()
1800
sage: L.local_factor(2)
1801
2048*T^2 + 24*T + 1
1802
sage: L._local_factor(11, None) # really this is called
1803
285311670611*T^2 - 534612*T + 1
1804
1805
The representation is reducible modulo 691::
1806
1807
sage: L.local_factor(2).factor_mod(691)
1808
(666) * (T + 387) * (T + 690)
1809
sage: L.local_factor(3).factor_mod(691)
1810
(251) * (T + 234) * (T + 690)
1811
sage: L.local_factor(11).factor_mod(691)
1812
(468) * (T + 471) * (T + 690)
1813
1814
... because of the 691 here::
1815
1816
sage: bernoulli(12)
1817
-691/2730
1818
"""
1819
try:
1820
return self._lf[P]
1821
except KeyError:
1822
pass
1823
self._precompute_local_factors(P+1)
1824
return self._lf[P]
1825
1826
def _precompute_local_factors(self, bound, prec=None):
1827
"""
1828
Precompute local factors up to the given bound.
1829
1830
EXAMPLES::
1831
1832
sage: from psage.lseries.eulerprod import LSeriesDelta; L = LSeriesDelta()
1833
sage: L._lf
1834
{}
1835
sage: L._precompute_local_factors(10)
1836
sage: L._lf
1837
{2: 2048*T^2 + 24*T + 1, 3: 177147*T^2 - 252*T + 1, 5: 48828125*T^2 - 4830*T + 1, 7: 1977326743*T^2 + 16744*T + 1}
1838
"""
1839
from sage.modular.all import delta_qexp
1840
T = self._T
1841
T2 = T**2
1842
f = delta_qexp(bound)
1843
for p in prime_range(bound):
1844
if not self._lf.has_key(p):
1845
self._lf[p] = 1 - f[p]*T + (p**11)*T2
1846
1847
def __repr__(self):
1848
return "L-function associated to Ramanujan's Delta (a weight 12 cusp form)"
1849
1850
def _cmp(self, right):
1851
return 0
1852
1853
1854
class LSeriesEllipticCurve(LSeriesAbstract):
1855
def __init__(self, E, prec=53):
1856
"""
1857
EXAMPLES::
1858
sage: from psage.lseries.eulerprod import LSeriesEllipticCurve
1859
sage: L = LSeriesEllipticCurve(EllipticCurve('389a'))
1860
sage: L(2)
1861
0.360092863578881
1862
"""
1863
E = E.global_minimal_model()
1864
self._E = E
1865
K = E.base_field()
1866
d = K.degree()
1867
self._lf = {}
1868
self._T = ZZ['T'].gen()
1869
self._N = E.conductor()
1870
if gcd(norm(self._N), K.discriminant())!=1:
1871
raise NotImplementedError, "Computation of conductor only implemented when the discrimenant of the base field and the norm of the conductor are coprime"
1872
LSeriesAbstract.__init__(self, conductor = norm(self._N) * K.discriminant()**2,
1873
hodge_numbers = [0]*d+[1]*d, weight = 2, epsilon = 1,#hardcoded to be 1 here
1874
poles = [], residues=[], base_field = K, prec=prec)
1875
1876
def elliptic_curve(self):
1877
return self._E
1878
1879
def _cmp(self, right):
1880
return cmp(self.elliptic_curve(), right.elliptic_curve())
1881
1882
def __repr__(self):
1883
return "L-series of %s"%self.elliptic_curve()
1884
1885
def _lf0(self, P):
1886
R = ZZ['T']
1887
T = R.gen()
1888
q = norm(P)
1889
p = prime_below(P)
1890
f = ZZ(q).ord(p)
1891
if P.divides(self._N):
1892
a = self._E.local_data(P).bad_reduction_type()
1893
return 1 - a*(T**f)
1894
else:
1895
a = q + 1 - self._E.reduction(P).count_points()
1896
return 1 - a*(T**f) + q*(T**(2*f))
1897
1898
def _local_factor(self, P, prec):
1899
if not self._lf.has_key(P):
1900
self._lf[P]=self._lf0(P)
1901
return self._lf[P]
1902
1903
1904
def _precompute_local_factors(self, bound, prec=None):
1905
for p in primes(bound):
1906
for q in self._primes_above(p):
1907
if q.norm()>=bound:
1908
continue
1909
if not self._lf.has_key(q):
1910
self._lf[q] = self._lf0(q)
1911
1912
1913
1914
1915
class LSeriesEllipticCurveQQ(LSeriesEllipticCurve):
1916
def __init__(self, E):
1917
E = E.global_minimal_model()
1918
self._E = E
1919
K = E.base_field()
1920
self._N = E.conductor()
1921
self._lf = {}
1922
self._T = ZZ['T'].gen()
1923
LSeriesAbstract.__init__(self, conductor = self._N,
1924
hodge_numbers = [0,1], weight = 2, epsilon = E.root_number(),
1925
poles = [], residues=[], base_field = QQ)
1926
1927
def _lf0(self, p):
1928
a = self._E.ap(p)
1929
T = self._T
1930
if self._N%p == 0:
1931
if self._N%(p*p) == 0:
1932
return T.parent()(1)
1933
else:
1934
return 1 - a*T
1935
else:
1936
return 1 - a*T + p*T*T
1937
1938
def _precompute_local_factors(self, bound, prec=None):
1939
for p in primes(bound):
1940
if not self._lf.has_key(p):
1941
self._lf[p] = self._lf0(p)
1942
1943
1944
1945
1946
class LSeriesEllipticCurveSqrt5(LSeriesEllipticCurve):
1947
def _precompute_local_factors(self, bound, prec=None):
1948
E = self._E
1949
1950
# Compute all of the prime ideals of the ring of integers up to the given bound
1951
from psage.number_fields.sqrt5.prime import primes_of_bounded_norm, Prime
1952
primes = primes_of_bounded_norm(bound)
1953
1954
# Compute the traces of Frobenius: this is supposed to be the hard part
1955
from psage.ellcurve.lseries.aplist_sqrt5 import aplist
1956
v = aplist(E, bound)
1957
1958
# Compute information about the primes of bad reduction, in
1959
# particular the integers i such that primes[i] is a prime of bad
1960
# reduction.
1961
bad_primes = set([Prime(a.prime()) for a in E.local_data()])
1962
1963
# Compute the local factors of the L-series.
1964
P = ZZ['T']
1965
T = P.gen()
1966
# Table of powers of T, so we don't have to compute T^4 (say) thousands of times.
1967
Tp = [T**i for i in range(5)]
1968
1969
# For each prime, we write down the local factor.
1970
if not hasattr(self, '_lf'):
1971
self._lf = {}
1972
for i, P in enumerate(primes):
1973
inertial_deg = 2 if P.is_inert() else 1
1974
a_p = v[i]
1975
if P in bad_primes:
1976
# bad reduction
1977
f = 1 - a_p*Tp[inertial_deg]
1978
else:
1979
# good reduction
1980
q = P.norm()
1981
f = 1 - a_p*Tp[inertial_deg] + q*Tp[2*inertial_deg]
1982
self._lf[P] = f
1983
1984
def _local_factor(self, P, prec):
1985
from psage.number_fields.sqrt5.prime import primes_of_bounded_norm, Prime
1986
if not isinstance(P, Prime):
1987
P = Prime(P)
1988
if self._lf.has_key(P):
1989
return self._lf[P]
1990
else:
1991
return LSeriesEllipticCurve._local_factor(self, P.sage_ideal())
1992
1993
def _primes_above(self, p):
1994
"""
1995
Return the primes above p. This function returns a special
1996
optimized prime of the ring of integers of Q(sqrt(5)).
1997
"""
1998
from prime import primes_above
1999
return primes_above(p)
2000
2001
class LSeriesDedekindZeta(LSeriesAbstract):
2002
"""
2003
EXAMPLES::
2004
2005
sage: from psage.lseries.eulerprod import LSeries
2006
sage: K.<a> = NumberField(x^3 - 2)
2007
sage: L = LSeries(K); L
2008
Dedekind Zeta function of Number Field in a with defining polynomial x^3 - 2
2009
sage: L(2)
2010
1.60266326190044
2011
sage: L.residues()
2012
'automatic'
2013
sage: L.residues(prec=53)
2014
[-4.77632833933856]
2015
sage: L.residues(prec=100)
2016
[-4.7763283393385594030639875094]
2017
"""
2018
def __init__(self, K):
2019
if not K.is_absolute():
2020
K = K.absolute_field(names='a')
2021
self._K = K
2022
d = K.degree()
2023
sigma = K.signature()[1]
2024
LSeriesAbstract.__init__(self,
2025
conductor = abs(K.discriminant()),
2026
hodge_numbers = [0]*(d-sigma) + [1]*sigma,
2027
weight = 1,
2028
epsilon = 1,
2029
poles = [1],
2030
residues = 'automatic',
2031
base_field = K,
2032
is_selfdual = True)
2033
self._T = ZZ['T'].gen()
2034
2035
def _cmp(self, right):
2036
return cmp(self.number_field(), right.number_field())
2037
2038
def number_field(self):
2039
return self._K
2040
2041
def __repr__(self):
2042
return "Dedekind Zeta function of %s"%self._K
2043
2044
def _local_factor(self, P, prec):
2045
T = self._T
2046
return 1 - T**P.residue_class_degree()
2047
2048
2049
class LSeriesDirichletCharacter(LSeriesAbstract):
2050
"""
2051
EXAMPLES::
2052
2053
sage: from psage.lseries.eulerprod import LSeries; L = LSeries(DirichletGroup(5).0)
2054
sage: L(3)
2055
0.988191681624057 + 0.0891051883457395*I
2056
"""
2057
def __init__(self, chi):
2058
if not chi.is_primitive():
2059
raise NotImplementedError, "chi must be primitive"
2060
if chi.is_trivial():
2061
raise NotImplementedError, "chi must be nontrivial"
2062
if chi.base_ring().characteristic() != 0:
2063
raise ValueError, "base ring must have characteristic 0"
2064
self._chi = chi
2065
LSeriesAbstract.__init__(self, conductor = chi.conductor(),
2066
hodge_numbers = [1] if chi.is_odd() else [0],
2067
weight = 1,
2068
epsilon = None,
2069
poles = [], # since primitive
2070
residues = [], # since primitive
2071
base_field = QQ,
2072
is_selfdual = chi.order() <= 2)
2073
self._T = ZZ['T'].gen()
2074
2075
def _cmp(self, right):
2076
return cmp(self.character(), right.character())
2077
2078
def __repr__(self):
2079
"""
2080
EXAMPLES::
2081
2082
sage: from psage.lseries.eulerprod import LSeries; L = LSeries(DirichletGroup(3).0)
2083
sage: L.__repr__()
2084
'L-series attached to Dirichlet character modulo 3 of conductor 3 mapping 2 |--> -1'
2085
"""
2086
return "L-series attached to %s"%self._chi
2087
2088
def character(self):
2089
return self._chi
2090
2091
def epsilon(self, prec=None):
2092
chi = self._chi
2093
if prec is None:
2094
# answer in symbolic ring
2095
return (sqrt(-1) * SR(chi.gauss_sum())) / chi.modulus().sqrt()
2096
else:
2097
C = ComplexField(prec)
2098
x = C(chi.modulus()).sqrt() / chi.gauss_sum_numerical(prec=prec)
2099
if chi.is_odd():
2100
x *= C.gen()
2101
return 1/x
2102
2103
def _local_factor(self, P, prec):
2104
a = self._chi(P)
2105
if prec is not None:
2106
a = ComplexField(prec)(a)
2107
return 1 - a*self._T
2108
2109
class LSeriesModularSymbolsAbstract(LSeriesAbstract):
2110
def _cmp(self, right):
2111
return cmp(self.modular_symbols(), right.modular_symbols())
2112
2113
def modular_symbols(self):
2114
return self._M
2115
2116
def __repr__(self):
2117
"""
2118
EXAMPLES::
2119
2120
sage: from psage.lseries.eulerprod import LSeries
2121
sage: f = Newforms(43,2,names='a')[1]; f
2122
q + a1*q^2 - a1*q^3 + (-a1 + 2)*q^5 + O(q^6)
2123
sage: LSeries(f).__repr__()
2124
'L-series of a degree 2 newform of level 43 and weight 2'
2125
"""
2126
return "L-series of a degree %s newform of level %s and weight %s"%(self._M.dimension(), self._M.level(), self._M.weight())
2127
2128
def _precompute_local_factors(self, bound, prec):
2129
primes = [p for p in prime_range(bound) if not self._lf.has_key(p) or self._lf[p][0] < prec]
2130
self._do_precompute(primes, prec)
2131
2132
def _do_precompute(self, primes, prec):
2133
E, v = self._M.compact_system_of_eigenvalues(primes)
2134
if prec == 53:
2135
C = CDF
2136
elif prec is None or prec==oo:
2137
if v.base_ring() == QQ:
2138
C = QQ
2139
else:
2140
C = CDF
2141
else:
2142
C = ComplexField(prec)
2143
2144
phi = v.base_ring().embeddings(C)[self._conjugate]
2145
v = vector(C, [phi(a) for a in v])
2146
aplist = E.change_ring(C) * v
2147
T = C['T'].gen(); T2 = T**2
2148
chi = self._M.character()
2149
k = self.weight()
2150
for i in range(len(primes)):
2151
p = primes[i]
2152
s = chi(p)
2153
if s != 0: s *= p**(k-1)
2154
F = 1 - aplist[i]*T + s*T2
2155
self._lf[p] = (prec, F)
2156
2157
def _local_factor(self, P, prec):
2158
# TODO: ugly -- get rid of all "prec=None" in whole program -- always use oo.
2159
if prec is None: prec = oo
2160
if self._lf.has_key(P) and self._lf[P][0] >= prec:
2161
return self._lf[P][1]
2162
else:
2163
self._do_precompute([P],prec)
2164
return self._lf[P][1]
2165
2166
class LSeriesModularSymbolsNewformGamma0(LSeriesModularSymbolsAbstract):
2167
def _cmp(self, right):
2168
return cmp((self._M, self._conjugate), (right._M, right._conjugate))
2169
2170
def __init__(self, M, conjugate=0, check=True, epsilon=None):
2171
"""
2172
INPUT:
2173
- M -- a simple, new, cuspidal modular symbols space with
2174
sign 1
2175
- conjugate -- (default: 0), integer between 0 and dim(M)-1
2176
- check -- (default: True), if True, checks that M is
2177
simple, new, cuspidal, which can take a very long time,
2178
depending on how M was defined
2179
- epsilon -- (default: None), if not None, should be the sign
2180
in the functional equation, which is -1 or 1. If this is
2181
None, then epsilon is computed by computing the sign of
2182
the main Atkin-Lehner operator on M. If you have a faster
2183
way to determine epsilon, use it.
2184
2185
EXAMPLES::
2186
2187
sage: from psage.lseries.eulerprod import LSeriesModularSymbolsNewformGamma0
2188
sage: M = ModularSymbols(43,sign=1).cuspidal_subspace()[1]; M
2189
Modular Symbols subspace of dimension 2 of Modular Symbols space of dimension 4 for Gamma_0(43) of weight 2 with sign 1 over Rational Field
2190
sage: L0 = LSeriesModularSymbolsNewformGamma0(M,0,check=False,epsilon=1); L0
2191
L-series of a degree 2 newform of level 43 and weight 2
2192
sage: L0.taylor_series()
2193
0.620539857407845 + 0.331674007376949*z - 0.226392184536589*z^2 + 0.0960519649929789*z^3 - 0.00451826124421802*z^4 - 0.0203363026933833*z^5 + O(z^6)
2194
sage: L1 = LSeriesModularSymbolsNewformGamma0(M,1,check=False,epsilon=1); L1
2195
L-series of a degree 2 newform of level 43 and weight 2
2196
sage: L1(1)
2197
0.921328017272472
2198
sage: L1.taylor_series()
2199
0.921328017272472 + 0.492443075089339*z - 0.391019352704047*z^2 + 0.113271812405127*z^3 + 0.0213067052584679*z^4 - 0.0344198080536274*z^5 + O(z^6)
2200
"""
2201
if M.dimension() == 0:
2202
raise ValueError, "modular symbols space must positive dimension"
2203
chi = M.character()
2204
if chi is None or not chi.is_trivial():
2205
raise ValueError, "modular symbols space must have trivial character"
2206
self._M = M
2207
N = M.level()
2208
2209
if check:
2210
if not M.is_simple():
2211
raise ValueError, "modular symbols space must be simple"
2212
if not M.is_new():
2213
raise ValueError, "modular symbols space must be new"
2214
if not M.is_cuspidal():
2215
raise ValueError, "modular symbols space must be cuspidal"
2216
2217
k = M.weight()
2218
if epsilon is None:
2219
w = M.atkin_lehner_operator(N).matrix()
2220
if w not in [-1, 1]:
2221
raise ValueError, "modular symbols space must have constant Atkin-Lehner operator"
2222
epsilon = (-1)**(k/2) * w[0,0]
2223
2224
conjugate = ZZ(conjugate)
2225
if conjugate < 0 or conjugate >= M.dimension():
2226
raise ValueError, "conjugate must a nonnegative integer less than the dimension"
2227
self._conjugate = conjugate
2228
self._lf = {}
2229
2230
LSeriesAbstract.__init__(self,
2231
conductor = N,
2232
hodge_numbers = [0,1],
2233
weight = k,
2234
epsilon = epsilon,
2235
poles = [], # since primitive
2236
residues = [], # since primitive
2237
base_field = QQ)
2238
2239
def _is_valid_modsym_space(M):
2240
if not is_ModularSymbolsSpace(M):
2241
raise TypeError, "must be a modular symbols space"
2242
if M.dimension() == 0:
2243
raise ValueError, "modular symbols space must positive dimension"
2244
if M.character() is None:
2245
raise ValueError, "modular symbols space must have associated character"
2246
if not M.is_simple():
2247
raise ValueError, "modular symbols space must be simple"
2248
if not M.is_new():
2249
raise ValueError, "modular symbols space must be new"
2250
if not M.is_cuspidal():
2251
raise ValueError, "modular symbols space must be cuspidal"
2252
2253
2254
class LSeriesModularSymbolsNewformCharacter(LSeriesModularSymbolsAbstract):
2255
def _cmp(self, right):
2256
return cmp((self._M, self._conjugate), (right._M, right._conjugate))
2257
2258
def __init__(self, M, conjugate=0):
2259
_is_valid_modsym_space(M)
2260
chi = M.character()
2261
self._M = M
2262
N = M.level()
2263
2264
# See Remark 5.0.2 in [Diamond-Im] which says: "Let f be a newform of level N
2265
# which is a common eigenform under all the Hecke operators T_p. Then
2266
# w_N(f) = c*fbar, where fbar = sum bar(a_n) q^n and c is a scalar. The functional
2267
# equation may be rewritten as Lambda(s, f) = c * i^k * Lambda(k-s, fbar).
2268
# That said, this seems hard to compute, so we just solve using the
2269
# functional equation.
2270
epsilon = 'solve'
2271
2272
k = M.weight()
2273
conjugate = ZZ(conjugate)
2274
if conjugate < 0 or conjugate >= M.dimension():
2275
raise ValueError, "conjugate must a nonnegative integer less than the dimension"
2276
self._conjugate = conjugate
2277
2278
2279
LSeriesAbstract.__init__(self,
2280
conductor = N,
2281
hodge_numbers = [0,1],
2282
weight = k,
2283
epsilon = epsilon,
2284
poles = [], # since primitive
2285
residues = [], # since primitive
2286
base_field = QQ,
2287
is_selfdual = chi.order() <= 2)
2288
self._lf = {}
2289
2290
class LSeriesModularEllipticCurveSqrt5(LSeriesAbstract):
2291
def __init__(self, E, **kwds):
2292
self._E = E
2293
self._N = E.conductor()
2294
self._T = ZZ['T'].gen()
2295
LSeriesAbstract.__init__(self,
2296
conductor = norm(self._N) * 25,
2297
hodge_numbers = [0,0,1,1],
2298
weight = 2,
2299
epsilon = [1,-1],
2300
poles = [], residues = [],
2301
base_field = E.base_field(),
2302
is_selfdual = True,
2303
**kwds)
2304
2305
def elliptic_curve(self):
2306
return self._E
2307
2308
def _cmp(self, right):
2309
return cmp(self.elliptic_curve(), right.elliptic_curve())
2310
2311
def __repr__(self):
2312
return "L-series attached to %s"%self._E
2313
2314
def _local_factor(self, P, prec):
2315
T = self._T
2316
v = self._N.valuation(P)
2317
if v >= 2:
2318
# additive reduction -- local factor is 1
2319
return T.parent()(1)
2320
elif v >= 1:
2321
# multiplicative reduction -- we have no algorithm yet to compute this
2322
# local factor, so we let the functional equation do the work.
2323
d = P.residue_class_degree()
2324
return [1 - T**d, 1 + T**d]
2325
else:
2326
# good reduction -- use Hecke operator
2327
a = self._E.ap(P)
2328
d = P.residue_class_degree()
2329
q = P.norm()
2330
return 1 - a*T**d + q*T**(2*d)
2331
2332
def _new_modsym_space_with_multiplicity(M):
2333
"""
2334
Returns a simple new modular symbols space N and an integer d such
2335
that M is isomorphic to `N^d` as a module over the anemic Hecke
2336
algebra.
2337
2338
INPUT:
2339
- M -- a sign=1 modular simple space for the full Hecke
2340
algebra (including primes dividing the level) that can't be
2341
decomposed further by the Hecke operators. None of the
2342
conditions on M are explicitly checked.
2343
2344
OUTPUT:
2345
- N -- a simple new modular symbols space
2346
- d -- a positive integer
2347
"""
2348
if M.is_new():
2349
return [(M,1)]
2350
raise NotImplementedError
2351
2352
def LSeriesModularSymbolsNewform(M, i=0):
2353
chi = M.character()
2354
if chi is None:
2355
raise NotImplementedError
2356
elif chi.is_trivial():
2357
return LSeriesModularSymbolsNewformGamma0(M, i)
2358
else:
2359
return LSeriesModularSymbolsNewformCharacter(M, i)
2360
2361
class LSeriesModularSymbolsMotive(LSeriesProduct):
2362
"""
2363
The product of L-series attached to the modular symbols space M.
2364
"""
2365
def __init__(self, M):
2366
self._M = M
2367
if not is_ModularSymbolsSpace(M):
2368
raise TypeError, "X must be a modular symbols space or have a modular symbols method"
2369
self._M = M
2370
D = M.decomposition()
2371
for A in D:
2372
_is_valid_modsym_space(A)
2373
F = []
2374
for A in D:
2375
for X in _new_modsym_space_with_multiplicity(A):
2376
N, d = X
2377
chi = N.character()
2378
for i in range(N.dimension()):
2379
F.append( (LSeriesModularSymbolsNewform(N,i), d))
2380
LSeriesProduct.__init__(self, F)
2381
2382
def modular_symbols(self):
2383
return self._M
2384
2385
def __repr__(self):
2386
return "L-series attached to %s"%self._M
2387
2388
class LSeriesModularAbelianVariety(LSeriesProduct):
2389
"""
2390
The product of L-series attached to the modular abelian variety A.
2391
2392
EXAMPLES::
2393
2394
sage: from psage.lseries.eulerprod import LSeries
2395
sage: L = LSeries(J0(54)); L
2396
L-series attached to Abelian variety J0(54) of dimension 4
2397
sage: L.factor()
2398
(L-series of a degree 1 newform of level 27 and weight 2)^2 * (L-series of a degree 1 newform of level 54 and weight 2) * (L-series of a degree 1 newform of level 54 and weight 2)
2399
sage: L(1)
2400
0.250717238804658
2401
sage: L.taylor_series(prec=20)
2402
0.25072 + 0.59559*z + 0.15099*z^2 - 0.35984*z^3 + 0.056934*z^4 + 0.17184*z^5 + O(z^6)
2403
2404
Independent check of L(1)::
2405
2406
sage: prod(EllipticCurve(lbl).lseries()(1) for lbl in ['54a', '54b', '27a', '27a'])
2407
0.250717238804658
2408
2409
Different check that totally avoids using Dokchitser::
2410
2411
sage: prod(EllipticCurve(lbl).lseries().at1()[0] for lbl in ['54a', '54b', '27a', '27a'])
2412
0.250848605530185
2413
"""
2414
def __init__(self, A):
2415
self._A = A
2416
D = A.decomposition()
2417
F = None
2418
for A in D:
2419
# TODO: This is ugly, but I don't know a cleaner way to do it yet.
2420
# Could be painfully inefficient in general.
2421
f = Newform(A.newform_label(), names='a')
2422
M = f.modular_symbols(sign=1)
2423
d = ZZ(A.dimension() / M.dimension())
2424
L = LSeriesModularSymbolsMotive(M)**d
2425
if F is None:
2426
F = L
2427
else:
2428
F *= L
2429
if F is None:
2430
raise ValueError, "abelian variety must have positive dimension"
2431
LSeriesProduct.__init__(self, F.factor())
2432
2433
def abelian_variety(self):
2434
return self._A
2435
2436
def __repr__(self):
2437
return "L-series attached to %s"%self._A
2438
2439
2440
#import psage.modform.hilbert.sqrt5.hmf
2441
2442
class LSeriesTwist(LSeriesAbstract):
2443
"""
2444
Twist of an L-series by a character.
2445
"""
2446
def __init__(self, L, chi, conductor=None, epsilon=None, prec=53):
2447
"""
2448
INPUT:
2449
- `L` -- an L-series
2450
- ``chi`` -- a character of the base field of L
2451
- ``conductor`` -- None, or a list of conductors to try
2452
- ``prec`` -- precision to use when trying conductors, if
2453
conductor is a list
2454
"""
2455
self._L = L
2456
self._chi = chi
2457
2458
if not chi.is_primitive():
2459
raise ValueError, "character must be primitive"
2460
2461
A = ZZ(L.conductor())
2462
B = chi.conductor()
2463
if conductor is None:
2464
if A.gcd(B) != 1:
2465
# Make a list of all possible conductors, and let the
2466
# functional equation figure it out.
2467
smallest = ZZ(A)
2468
while smallest.gcd(B) != 1:
2469
smallest = smallest // smallest.gcd(B)
2470
biggest = A * (B**L.degree())
2471
assert biggest % smallest == 0
2472
#
2473
# TODO: improve this using the theorem stated
2474
# on page 1 of http://wstein.org/papers/padictwist/
2475
#
2476
conductor = [smallest*d for d in divisors(biggest//smallest)]
2477
else:
2478
conductor = A * (B**L.degree())
2479
hodge_numbers = L.hodge_numbers()
2480
weight = L.weight()
2481
if epsilon is None:
2482
if L.epsilon() != 'solve':
2483
if chi.order() <= 2:
2484
if A.gcd(B) == 1:
2485
epsilon = L.epsilon() * chi(-A)
2486
else:
2487
epsilon = [L.epsilon(), -L.epsilon()]
2488
else:
2489
epsilon = 'solve'
2490
else:
2491
epsilon = 'solve'
2492
is_selfdual = L.is_selfdual()
2493
poles = [] # ??? TODO -- no clue here.
2494
residues = 'automatic'
2495
base_field = L.base_field()
2496
LSeriesAbstract.__init__(self, conductor, hodge_numbers, weight, epsilon,
2497
poles, residues, base_field, is_selfdual, prec)
2498
2499
def _local_factor(self, P, prec):
2500
L0 = self._L.local_factor(P, prec)
2501
chi = self._chi
2502
T = L0.parent().gen()
2503
c = chi(P)
2504
if prec is not None:
2505
c = ComplexField(prec)(c)
2506
return L0(c*T)
2507
2508
def __repr__(self):
2509
return "Twist of %s by %s"%(self._L, self._chi)
2510
2511
def _cmp(self, right):
2512
return cmp((self._L, self._chi), (right._L, right._chi))
2513
2514
def untwisted_lseries(self):
2515
return self._L
2516
2517
def twist_character(self):
2518
return self._chi
2519
2520
##############
2521
# TODO: General tensor products and symmetric powers: see
2522
# http://magma.maths.usyd.edu.au/magma/handbook/text/1392#15272
2523
##############
2524
2525
2526
def LSeries(X, *args, **kwds):
2527
"""
2528
Return the L-series of X, where X can be any of the following:
2529
2530
- elliptic curve over a number field (including QQ)
2531
- Dirichlet character
2532
- cuspidal newform
2533
- new cuspidal modular symbols space -- need not be simple
2534
- string: 'zeta' (Riemann Zeta function), 'delta'
2535
- modular elliptic curve attached to Hilbert modular forms space
2536
2537
For convenience, if L is returned, then L._X is set to X.
2538
2539
EXAMPLES::
2540
2541
The Dedekind Zeta function of a number field::
2542
2543
sage: from psage.lseries.eulerprod import LSeries
2544
sage: K.<a> = NumberField(x^2 + 1)
2545
sage: L = LSeries(K); L
2546
Dedekind Zeta function of Number Field in a with defining polynomial x^2 + 1
2547
sage: L(2)
2548
1.50670300992299
2549
2550
sage: K.zeta_coefficients(100) == L.anlist(100)[1:]
2551
True
2552
2553
sage: L = LSeries(ModularSymbols(43, weight=2,sign=1).cuspidal_subspace().decomposition()[1])
2554
sage: L(1)
2555
0.571720756464112
2556
sage: L.factor()[0][0](1)
2557
0.620539857407845
2558
sage: L.factor()[1][0](1)
2559
0.921328017272472
2560
sage: L._X
2561
Modular Symbols subspace of dimension 2 of Modular Symbols space of dimension 4 for Gamma_0(43) of weight 2 with sign 1 over Rational Field
2562
sage: L = LSeries(ModularSymbols(DirichletGroup(13).0^2, weight=2,sign=1).cuspidal_subspace())
2563
sage: L(1)
2564
0.298115272465799 - 0.0402203326076733*I
2565
2566
sage: from psage.modform.hilbert.sqrt5.hmf import F, HilbertModularForms
2567
sage: D = HilbertModularForms(-13*F.0+5).elliptic_curve_factors()
2568
sage: D
2569
[
2570
Isogeny class of elliptic curves over QQ(sqrt(5)) attached to form number 0 in Hilbert modular forms of dimension 4, level -13*a+5 (of norm 209=11*19) over QQ(sqrt(5)),
2571
Isogeny class of elliptic curves over QQ(sqrt(5)) attached to form number 1 in Hilbert modular forms of dimension 4, level -13*a+5 (of norm 209=11*19) over QQ(sqrt(5)),
2572
Isogeny class of elliptic curves over QQ(sqrt(5)) attached to form number 2 in Hilbert modular forms of dimension 4, level -13*a+5 (of norm 209=11*19) over QQ(sqrt(5))
2573
]
2574
sage: L = LSeries(D[0])
2575
sage: L(RealField(10)(1))
2576
0
2577
sage: L.epsilon()
2578
-1
2579
2580
The L-series of a modular abelian variety with both new and old parts::
2581
2582
sage: L = LSeries(J0(33)); L
2583
L-series attached to Abelian variety J0(33) of dimension 3
2584
sage: L.factor()
2585
(L-series of a degree 1 newform of level 11 and weight 2)^2 * (L-series of a degree 1 newform of level 33 and weight 2)
2586
sage: L.local_factor(2, prec=oo)
2587
8*T^6 + 12*T^5 + 12*T^4 + 8*T^3 + 6*T^2 + 3*T + 1
2588
sage: L(1)
2589
0.0481553138900504
2590
2591
We check the above computation of L(1) via independent methods (and implementations)::
2592
2593
sage: prod(EllipticCurve(lbl).lseries().at1()[0] for lbl in ['11a', '11a', '33a'])
2594
0.0481135342926321
2595
sage: prod(EllipticCurve(lbl).lseries()(1) for lbl in ['11a', '11a', '33a'])
2596
0.0481553138900504
2597
2598
A nonsimple new modular symbols space of level 43::
2599
2600
sage: L = LSeries(ModularSymbols(43,sign=1).cuspidal_subspace())
2601
sage: L
2602
L-series attached to Modular Symbols subspace of dimension 3 of Modular Symbols space of dimension 4 for Gamma_0(43) of weight 2 with sign 1 over Rational Field
2603
sage: L(1)
2604
0
2605
sage: L.taylor_series()
2606
0.196399786632435*z + 0.314922741074845*z^2 - 0.0797083673829092*z^3 - 0.161630566287135*z^4 + 0.123939472976207*z^5 + O(z^6)
2607
sage: L.factor()
2608
(L-series of a degree 1 newform of level 43 and weight 2) * (L-series of a degree 2 newform of level 43 and weight 2) * (L-series of a degree 2 newform of level 43 and weight 2)
2609
sage: L.analytic_rank()
2610
1
2611
sage: D = ModularSymbols(43,sign=1).cuspidal_subspace().decomposition()
2612
sage: L0 = LSeries(D[0]); L1 = LSeries(D[1])
2613
sage: L0.taylor_series() * L1.taylor_series()
2614
0.196399786632435*z + 0.314922741074845*z^2 - 0.0797083673829091*z^3 - 0.161630566287135*z^4 + 0.123939472976207*z^5 + O(z^6)
2615
sage: L0.factor()
2616
L-series of a degree 1 newform of level 43 and weight 2
2617
sage: L1.factor()
2618
(L-series of a degree 2 newform of level 43 and weight 2) * (L-series of a degree 2 newform of level 43 and weight 2)
2619
2620
"""
2621
L = _lseries(X, *args, **kwds)
2622
L._X = X
2623
return L
2624
2625
def _lseries(X, *args, **kwds):
2626
"""
2627
Helper function used by LSeries function.
2628
"""
2629
if is_EllipticCurve(X):
2630
K = X.base_ring()
2631
if is_RationalField(K):
2632
return LSeriesEllipticCurveQQ(X, *args, **kwds)
2633
elif list(K.defining_polynomial()) == [-1,-1,1]:
2634
return LSeriesEllipticCurveSqrt5(X, *args, **kwds)
2635
else:
2636
return LSeriesEllipticCurve(X)
2637
2638
if is_DirichletCharacter(X):
2639
if X.is_trivial() and X.is_primitive():
2640
return LSeriesZeta(*args, **kwds)
2641
else:
2642
return LSeriesDirichletCharacter(X, *args, **kwds)
2643
2644
if is_NumberField(X):
2645
return LSeriesDedekindZeta(X, *args, **kwds)
2646
2647
if isinstance(X, sage.modular.modform.element.Newform):
2648
return LSeriesModularSymbolsNewform(X.modular_symbols(sign=1), *args, **kwds)
2649
2650
if is_ModularSymbolsSpace(X):
2651
if X.sign() != 1:
2652
raise NotImplementedError
2653
return LSeriesModularSymbolsMotive(X, *args, **kwds)
2654
2655
if isinstance(X, str):
2656
y = X.lower()
2657
if y == 'zeta':
2658
return LSeriesZeta(*args, **kwds)
2659
elif y == 'delta':
2660
return LSeriesDelta(*args, **kwds)
2661
else:
2662
raise ValueError, 'unknown L-series "%s"'%y
2663
2664
if isinstance(X, psage.modform.hilbert.sqrt5.hmf.EllipticCurveFactor):
2665
return LSeriesModularEllipticCurveSqrt5(X, *args, **kwds)
2666
2667
if is_ModularAbelianVariety(X):
2668
return LSeriesModularAbelianVariety(X, *args, **kwds)
2669
2670
raise NotImplementedError
2671
2672