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