22.3 Dokchitser's L-functions Calculator

Module: sage.lfunctions.dokchitser

Author Log:

TODO: - add more examples from data/extcode/pari/dokchitser that illustrate use with Eisenstein series, number fields, etc. - plug this code into number fields and modular forms code (elliptic curves are done).

Module-level Functions

reduce_load_dokchitser( D)

Class: Dokchitser

class Dokchitser
Dokchitser's $ L$ -functions Calculator

Create a Dokchitser $ L$ -series with

Dokchitser(conductor, gammaV, weight, eps, poles, residues, init, prec)

where

conductor - integer, the conductor gammaV - list of Gamma-factor parameters, e.g. [0] for Riemann zeta, [0,1] for ell.curves, (see examples). weight - positive real number, usually an integer e.g. 1 for Riemann zeta, 2 for $ H^1$ of curves/ $ \mathbf{Q}$ eps - complex number; sign in functional equation poles - (default: []) list of points where $ L^*(s)$ has (simple) poles; only poles with Re(s)>weight/2 should be included residues - vector of residues of $ L^*(s)$ in those poles or set residues='automatic' (default value) prec - integer (default: 53) number of *bits* of precision

RIEMANN ZETA FUNCTION:

We compute with the Riemann Zeta function.

sage: L = Dokchitser(conductor=1, gammaV=[0], weight=1, eps=1, poles=[1], residues=[-1], init='1')
sage: L
Dokchitser L-series of conductor 1 and weight 1
sage: L(1)
Traceback (most recent call last):
...
ArithmeticError:   ###   user error: L*(s) has a pole at
s=1.000000000000000000
sage: L(2)
1.6449340668482264
sage: L(2, 1.1)
1.6449340668482264
sage: L.derivative(2)
-0.93754825431584377
sage: h = RR('0.0000000000001')
sage: (zeta(2+h) - zeta(2))/h
-0.93702823278363212
sage: L.taylor_series(2, k=5)
1.6449340668482264 + -0.93754825431584377*z + 0.99464011714945055*z^2 +
-1.0000243004738407*z^3 + 1.0000619330723526*z^4 + O(z^5)

RANK 1 ELLIPTIC CURVE:

We compute with the $ L$ -series of a rank $ 1$ curve.

sage: E = EllipticCurve('37a')
sage: L = E.Lseries_dokchitser(); L
Dokchitser L-function associated to Elliptic Curve defined by y^2 + y = x^3
- x over Rational Field
sage: L(1)
0
sage: L.derivative(1)
0.30599977383405230
sage: L.derivative(1,2)
0.37309559453632390
sage: L.num_coeffs()
48
sage: L.taylor_series(1,4)
0.30599977383405230*z + 0.18654779726816195*z^2 + -0.13679146309718768*z^3
+ O(z^4)
sage: L.check_functional_equation()
0.0000000000000000061121897480000001            # 32-bit
0.0000000000000000060444271116099998            # 64-bit

RANK 2 ELLIPTIC CURVE:

We compute the leading coefficient and Taylor expansion of the $ L$ -series of a rank $ 2$ curve.

sage: E = EllipticCurve('389a')
sage: L = E.Lseries_dokchitser()
sage: L.num_coeffs ()
156
sage: L.derivative(1,E.rank())
1.5186330005768536
sage: L.taylor_series(1,4)
-0.000000000000000000000012815814569193140 +
0.0000000000000000000000072626829063558658*z + 0.75931650028842679*z^2 +
-0.43030233758336200*z^3 + O(z^4)  # 32-bit
-0.000000000000000000000026912956656300001 +
0.000000000000000000000015251490196900001*z + 0.75931650028842679*z^2 +
-0.43030233758336200*z^3 + O(z^4)   # 64-bit

RAMANUJAN DELTA L-FUNCTION:

The coefficients are given by Ramanujan's tau function:

sage: L = Dokchitser(conductor=1, gammaV=[0,1], weight=12, eps=1)
sage: pari_precode = 'tau(n)=(5*sigma(n,3)+7*sigma(n,5))*n/12 - 35*sum(k=1,n-1,(6*k-4*(n-k))*sigma(k,3)*sigma(n-k,5))'
sage: L.init_coeffs('tau(k)', pari_precode=pari_precode)

We redefine the default bound on the coefficients: Deligne's estimate on tau(n) is better than the default coefgrow(n)= $ (4n)^{11/2}$ (by a factor 1024), so re-defining coefgrow() improves efficiency (slightly faster).

sage: L.num_coeffs()
12
sage: L.set_coeff_growth('2*n^(11/2)')
sage: L.num_coeffs()
11

Now we're ready to evaluate, etc.

sage: L(1)
0.037441281268515543
sage: L(1, 1.1)
0.037441281268515543
sage: L.taylor_series(1,3)
0.037441281268515543 + 0.070922112361932230*z + 0.038074476127051969*z^2 +
O(z^3)
Dokchitser( self, conductor, gammaV, weight, eps, [poles=None], [residues=53], [prec=automatic], [init=[]])

Functions: check_functional_equation,$  $ derivative,$  $ gp,$  $ init_coeffs,$  $ num_coeffs,$  $ set_coeff_growth,$  $ taylor_series

check_functional_equation( self, [T=1.2])

Verifies how well numerically the functional equation is satisfied, and also determines the residues if self.poles != [] and residues='automatic'.

More specifically: for $ T>1$ (default 1.2), self.check_functional_equation(T) should ideally return 0 (to the current precision).

sage: L = Dokchitser(conductor=1, gammaV=[0], weight=1, eps=1, poles=[1], residues=[-1], init='1')
sage: L.check_functional_equation ()
-0.000000000000000000027105054319999997      # 32-bit
-0.000000000000000000027105054312100001      # 64-bit

If we choose the sign in functional equation for the $ \zeta$ function incorrectly, the functional equation doesn't check out.

sage: L = Dokchitser(conductor=1, gammaV=[0], weight=1, eps=-11, poles=[1], residues=[-1], init='1')
sage: L.check_functional_equation ()
-9.7396786148812371

derivative( self, s, [k=1])

Return the k-th derivative of the L-series at s.

gp( self)

Return the gp interpreter that is used to implement this Dokchitser L-function.

sage: E = EllipticCurve('11a')
sage: L = E.Lseries_dokchitser()
sage: L(2)
0.54604803621501352
sage: L.gp()
GP/PARI interpreter

init_coeffs( self, v, [cutoff=40], [w=0], [pari_precode=], [max_imaginary_part=None], [max_asymp_coeffs=1])

Set the coefficients $ a_n$ of the $ L$ -series. If $ L(s)$ is not equal to its dual, pass the coefficients of the dual as the second optional argument.

INPUT:
    v -- list of complex numbers or string (pari function of k)
    cutoff -- real number >= 1 (default: 1)
    w -- list of complex numbers or string (pari function of k)
    pari_precode -- some code to execute in pari before calling initLdata
    max_imaginary_part -- (default: 0): redefine if you want to compute
                          L(s) for s having large imaginary part,
    max_asymp_coeffs --   (default: 40): at most this many terms are
                          generated in asymptotic series for phi(t)
                          and G(s,t).

num_coeffs( self, [T=1])

Return number of coefficients $ a_n$ that are needed in order to perform most relevant $ L$ -function computations to the desired precision.

sage: E = EllipticCurve('11a')
sage: L = E.Lseries_dokchitser()
sage: L.num_coeffs()
26
sage: E = EllipticCurve('5077a')
sage: L = E.Lseries_dokchitser()
sage: L.num_coeffs()
568
sage: L = Dokchitser(conductor=1, gammaV=[0], weight=1, eps=1, poles=[1], residues=[-1], init='1')
sage: L.num_coeffs()
4

set_coeff_growth( self, coefgrow)

You might have to redefine the coefficient growth function if the $ a_n$ of the $ L$ -series are not given by the following PARI function:

            coefgrow(n) = if(length(Lpoles),   
                              1.5*n^(vecmax(real(Lpoles))-1),  
                              sqrt(4*n)^(weight-1));

INPUT:
    coefgrow -- string that evaluates to a PARI function of n
                that defines a coefgrow function.

taylor_series( self, [a=z], [k=6], [var=0])

Return the first k terms of the Taylor series expansion of the $ L$ -series about $ a$ .

This is returned as a series in var, where you should view var as equal to $ s-a$ . Thus this function returns the formal power series whose coefficients are $ L^{(n)}(a)/n!$ .

INPUT:
    a -- complex number (default: 0); point about which to expand
    k -- integer (default: 6), series is O(var^k)
    var -- string (default: 'z'), variable of power series

sage: L = Dokchitser(conductor=1, gammaV=[0], weight=1, eps=1, poles=[1], residues=[-1], init='1')
sage: L.taylor_series(2, 3)
1.6449340668482264 + -0.93754825431584377*z + 0.99464011714945055*z^2 +
O(z^3)
sage: E = EllipticCurve('37a')
sage: L = E.Lseries_dokchitser()
sage: L.taylor_series(1)
0.30599977383405230*z + 0.18654779726816195*z^2 + -0.13679146309718768*z^3
+ 0.016106646849640055*z^4 + 0.018595517539880219*z^5 + O(z^6)

Special Functions: __call__,$  $ __reduce__,$  $ _Dokchitser__check_init,$  $ _Dokchitser__to_CC,$  $ _repr_

__call__( self, s, [c=None])

INPUT:
    s -- complex number
    c -- (optional); cutoff which should be a real number > 1
         that is chosen close to 1.

sage: E = EllipticCurve('5077a')
sage: L = E.Lseries_dokchitser(100)
sage: L(1)
0
sage: L(1+I)
-1.3085436607849493358323930438319 + 0.81298000036784359634835412129371*I

See About this document... for information on suggesting changes.