Dokchitser’s \(L\)-functions calculator#
Todo
add more examples from SAGE_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).
AUTHORS:
Tim Dokchitser (2002): original PARI code and algorithm (and the documentation below is based on Dokchitser’s docs).
William Stein (2006-03-08): Sage interface
- class sage.lfunctions.dokchitser.Dokchitser(*args, **kwargs)[source]#
Bases:
SageObject
Dokchitser’s \(L\)-functions Calculator
Create a Dokchitser \(L\)-series with
Dokchitser(conductor, gammaV, weight, eps, poles, residues, init, prec)
where
conductor
– integer, the conductorgammaV
– 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/\(\QQ\)eps
– complex number; sign in functional equationpoles
– (default:[]
) list of points where \(L^*(s)\) has (simple) poles; only poles with \(Re(s)>weight/2\) should be includedresidues
– 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 sage: L(2) 1.64493406684823 sage: L(2, 1.1) 1.64493406684823 sage: L.derivative(2) -0.937548254315844 sage: h = RR('0.0000000000001') sage: (zeta(2+h) - zeta(2.))/h -0.937028232783632 sage: L.taylor_series(2, k=5) 1.64493406684823 - 0.937548254315844*z + 0.994640117149451*z^2 - 1.00002430047384*z^3 + 1.00006193307...*z^4 + O(z^5)
>>> from sage.all import * >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0)], weight=Integer(1), eps=Integer(1), poles=[Integer(1)], residues=[-Integer(1)], init='1') >>> L Dokchitser L-series of conductor 1 and weight 1 >>> L(Integer(1)) Traceback (most recent call last): ... ArithmeticError >>> L(Integer(2)) 1.64493406684823 >>> L(Integer(2), RealNumber('1.1')) 1.64493406684823 >>> L.derivative(Integer(2)) -0.937548254315844 >>> h = RR('0.0000000000001') >>> (zeta(Integer(2)+h) - zeta(RealNumber('2.')))/h -0.937028232783632 >>> L.taylor_series(Integer(2), k=Integer(5)) 1.64493406684823 - 0.937548254315844*z + 0.994640117149451*z^2 - 1.00002430047384*z^3 + 1.00006193307...*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(algorithm='gp'); L Dokchitser L-function associated to Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field sage: L(1) 0.000000000000000 sage: L.derivative(1) 0.305999773834052 sage: L.derivative(1,2) 0.373095594536324 sage: L.num_coeffs() 48 sage: L.taylor_series(1,4) 0.000000000000000 + 0.305999773834052*z + 0.186547797268162*z^2 - 0.136791463097188*z^3 + O(z^4) sage: L.check_functional_equation() # abs tol 1e-19 6.04442711160669e-18
>>> from sage.all import * >>> E = EllipticCurve('37a') >>> L = E.lseries().dokchitser(algorithm='gp'); L Dokchitser L-function associated to Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field >>> L(Integer(1)) 0.000000000000000 >>> L.derivative(Integer(1)) 0.305999773834052 >>> L.derivative(Integer(1),Integer(2)) 0.373095594536324 >>> L.num_coeffs() 48 >>> L.taylor_series(Integer(1),Integer(4)) 0.000000000000000 + 0.305999773834052*z + 0.186547797268162*z^2 - 0.136791463097188*z^3 + O(z^4) >>> L.check_functional_equation() # abs tol 1e-19 6.04442711160669e-18
RANK 2 ELLIPTIC CURVE:
We compute the leading coefficient and Taylor expansion of the \(L\)-series of a rank \(2\) elliptic curve.
sage: E = EllipticCurve('389a') sage: L = E.lseries().dokchitser(algorithm='gp') sage: L.num_coeffs() 156 sage: L.derivative(1,E.rank()) 1.51863300057685 sage: L.taylor_series(1,4) -1.27685190980159e-23 + (7.23588070754027e-24)*z + 0.759316500288427*z^2 - 0.430302337583362*z^3 + O(z^4) # 32-bit -2.72911738151096e-23 + (1.54658247036311e-23)*z + 0.759316500288427*z^2 - 0.430302337583362*z^3 + O(z^4) # 64-bit
>>> from sage.all import * >>> E = EllipticCurve('389a') >>> L = E.lseries().dokchitser(algorithm='gp') >>> L.num_coeffs() 156 >>> L.derivative(Integer(1),E.rank()) 1.51863300057685 >>> L.taylor_series(Integer(1),Integer(4)) -1.27685190980159e-23 + (7.23588070754027e-24)*z + 0.759316500288427*z^2 - 0.430302337583362*z^3 + O(z^4) # 32-bit -2.72911738151096e-23 + (1.54658247036311e-23)*z + 0.759316500288427*z^2 - 0.430302337583362*z^3 + O(z^4) # 64-bit
NUMBER FIELD:
We compute with the Dedekind zeta function of a number field.
sage: x = var('x') sage: K = NumberField(x**4 - x**2 - 1,'a') sage: L = K.zeta_function(algorithm='gp') sage: L.conductor 400 sage: L.num_coeffs() 264 sage: L(2) 1.10398438736918 sage: L.taylor_series(2,3) 1.10398438736918 - 0.215822638498759*z + 0.279836437522536*z^2 + O(z^3)
>>> from sage.all import * >>> x = var('x') >>> K = NumberField(x**Integer(4) - x**Integer(2) - Integer(1),'a') >>> L = K.zeta_function(algorithm='gp') >>> L.conductor 400 >>> L.num_coeffs() 264 >>> L(Integer(2)) 1.10398438736918 >>> L.taylor_series(Integer(2),Integer(3)) 1.10398438736918 - 0.215822638498759*z + 0.279836437522536*z^2 + O(z^3)
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)
>>> from sage.all import * >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0),Integer(1)], weight=Integer(12), eps=Integer(1)) >>> 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))' >>> 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
>>> from sage.all import * >>> L.num_coeffs() 12 >>> L.set_coeff_growth('2*n^(11/2)') >>> L.num_coeffs() 11
Now we’re ready to evaluate, etc.
sage: L(1) 0.0374412812685155 sage: L(1, 1.1) 0.0374412812685155 sage: L.taylor_series(1,3) 0.0374412812685155 + 0.0709221123619322*z + 0.0380744761270520*z^2 + O(z^3)
>>> from sage.all import * >>> L(Integer(1)) 0.0374412812685155 >>> L(Integer(1), RealNumber('1.1')) 0.0374412812685155 >>> L.taylor_series(Integer(1),Integer(3)) 0.0374412812685155 + 0.0709221123619322*z + 0.0380744761270520*z^2 + O(z^3)
- check_functional_equation(T=1.2)[source]#
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).if what this function returns does not look like 0 at all, probably the functional equation is wrong (i.e. some of the parameters gammaV, conductor etc., or the coefficients are wrong),
if checkfeq(T) is to be used, more coefficients have to be generated (approximately T times more), e.g. call cflength(1.3), initLdata(“a(k)”,1.3), checkfeq(1.3)
T=1 always (!) returns 0, so T has to be away from 1
default value \(T=1.2\) seems to give a reasonable balance
if you don’t have to verify the functional equation or the L-values, call num_coeffs(1) and initLdata(“a(k)”,1), you need slightly less coefficients.
EXAMPLES:
sage: L = Dokchitser(conductor=1, gammaV=[0], weight=1, eps=1, poles=[1], residues=[-1], init='1') sage: L.check_functional_equation() # abs tol 1e-19 -2.71050543121376e-20
>>> from sage.all import * >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0)], weight=Integer(1), eps=Integer(1), poles=[Integer(1)], residues=[-Integer(1)], init='1') >>> L.check_functional_equation() # abs tol 1e-19 -2.71050543121376e-20
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.73967861488124
>>> from sage.all import * >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0)], weight=Integer(1), eps=-Integer(11), poles=[Integer(1)], residues=[-Integer(1)], init='1') >>> L.check_functional_equation() -9.73967861488124
- derivative(s, k=1)[source]#
Return the \(k\)-th derivative of the \(L\)-series at \(s\).
Warning
If \(k\) is greater than the order of vanishing of \(L\) at \(s\) you may get nonsense.
EXAMPLES:
sage: E = EllipticCurve('389a') sage: L = E.lseries().dokchitser(algorithm='gp') sage: L.derivative(1,E.rank()) 1.51863300057685
>>> from sage.all import * >>> E = EllipticCurve('389a') >>> L = E.lseries().dokchitser(algorithm='gp') >>> L.derivative(Integer(1),E.rank()) 1.51863300057685
- gp()[source]#
Return the gp interpreter that is used to implement this Dokchitser L-function.
EXAMPLES:
sage: E = EllipticCurve('11a') sage: L = E.lseries().dokchitser(algorithm='gp') sage: L(2) 0.546048036215014 sage: L.gp() PARI/GP interpreter
>>> from sage.all import * >>> E = EllipticCurve('11a') >>> L = E.lseries().dokchitser(algorithm='gp') >>> L(Integer(2)) 0.546048036215014 >>> L.gp() PARI/GP interpreter
- init_coeffs(v, cutoff=1, w=None, pari_precode='', max_imaginary_part=0, max_asymp_coeffs=40)[source]#
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 initLdatamax_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).
EXAMPLES:
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)
>>> from sage.all import * >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0),Integer(1)], weight=Integer(12), eps=Integer(1)) >>> 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))' >>> L.init_coeffs('tau(k)', pari_precode=pari_precode)
Evaluate the resulting L-function at a point, and compare with the answer that one gets “by definition” (of L-function attached to a modular form):
sage: L(14) 0.998583063162746 sage: a = delta_qexp(1000) sage: sum(a[n]/float(n)^14 for n in reversed(range(1,1000))) 0.9985830631627461
>>> from sage.all import * >>> L(Integer(14)) 0.998583063162746 >>> a = delta_qexp(Integer(1000)) >>> sum(a[n]/float(n)**Integer(14) for n in reversed(range(Integer(1),Integer(1000)))) 0.9985830631627461
Illustrate that one can give a list of complex numbers for v (see Issue #10937):
sage: L2 = Dokchitser(conductor=1, gammaV=[0,1], weight=12, eps=1) sage: L2.init_coeffs(list(delta_qexp(1000))[1:]) sage: L2(14) 0.998583063162746
>>> from sage.all import * >>> L2 = Dokchitser(conductor=Integer(1), gammaV=[Integer(0),Integer(1)], weight=Integer(12), eps=Integer(1)) >>> L2.init_coeffs(list(delta_qexp(Integer(1000)))[Integer(1):]) >>> L2(Integer(14)) 0.998583063162746
- num_coeffs(T=1)[source]#
Return number of coefficients \(a_n\) that are needed in order to perform most relevant \(L\)-function computations to the desired precision.
EXAMPLES:
sage: E = EllipticCurve('11a') sage: L = E.lseries().dokchitser(algorithm='gp') sage: L.num_coeffs() 26 sage: E = EllipticCurve('5077a') sage: L = E.lseries().dokchitser(algorithm='gp') 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
>>> from sage.all import * >>> E = EllipticCurve('11a') >>> L = E.lseries().dokchitser(algorithm='gp') >>> L.num_coeffs() 26 >>> E = EllipticCurve('5077a') >>> L = E.lseries().dokchitser(algorithm='gp') >>> L.num_coeffs() 568 >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0)], weight=Integer(1), eps=Integer(1), poles=[Integer(1)], residues=[-Integer(1)], init='1') >>> L.num_coeffs() 4
Verify that
num_coeffs
works with non-real spectral parameters, e.g. for the L-function of the level 10 Maass form with eigenvalue 2.7341055592527126:sage: ev = 2.7341055592527126 sage: L = Dokchitser(conductor=10, gammaV=[ev*i, -ev*i],weight=2,eps=1) sage: L.num_coeffs() 26
>>> from sage.all import * >>> ev = RealNumber('2.7341055592527126') >>> L = Dokchitser(conductor=Integer(10), gammaV=[ev*i, -ev*i],weight=Integer(2),eps=Integer(1)) >>> L.num_coeffs() 26
- set_coeff_growth(coefgrow)[source]#
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.
EXAMPLES:
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) sage: L.set_coeff_growth('2*n^(11/2)') sage: L(1) 0.0374412812685155
>>> from sage.all import * >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0),Integer(1)], weight=Integer(12), eps=Integer(1)) >>> 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))' >>> L.init_coeffs('tau(k)', pari_precode=pari_precode) >>> L.set_coeff_growth('2*n^(11/2)') >>> L(Integer(1)) 0.0374412812685155
- taylor_series(a=0, k=6, var='z')[source]#
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 viewvar
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 expandk
– integer (default: 6), series is \(O(``var``^k)\)var
– string (default: ‘z’), variable of power series
EXAMPLES:
sage: L = Dokchitser(conductor=1, gammaV=[0], weight=1, eps=1, poles=[1], residues=[-1], init='1') sage: L.taylor_series(2, 3) 1.64493406684823 - 0.937548254315844*z + 0.994640117149451*z^2 + O(z^3) sage: E = EllipticCurve('37a') sage: L = E.lseries().dokchitser(algorithm='gp') sage: L.taylor_series(1) 0.000000000000000 + 0.305999773834052*z + 0.186547797268162*z^2 - 0.136791463097188*z^3 + 0.0161066468496401*z^4 + 0.0185955175398802*z^5 + O(z^6)
>>> from sage.all import * >>> L = Dokchitser(conductor=Integer(1), gammaV=[Integer(0)], weight=Integer(1), eps=Integer(1), poles=[Integer(1)], residues=[-Integer(1)], init='1') >>> L.taylor_series(Integer(2), Integer(3)) 1.64493406684823 - 0.937548254315844*z + 0.994640117149451*z^2 + O(z^3) >>> E = EllipticCurve('37a') >>> L = E.lseries().dokchitser(algorithm='gp') >>> L.taylor_series(Integer(1)) 0.000000000000000 + 0.305999773834052*z + 0.186547797268162*z^2 - 0.136791463097188*z^3 + 0.0161066468496401*z^4 + 0.0185955175398802*z^5 + O(z^6)
We compute a Taylor series where each coefficient is to high precision.
sage: E = EllipticCurve('389a') sage: L = E.lseries().dokchitser(200, algorithm='gp') sage: L.taylor_series(1,3) ...e-82 + (...e-82)*z + 0.75931650028842677023019260789472201907809751649492435158581*z^2 + O(z^3)
>>> from sage.all import * >>> E = EllipticCurve('389a') >>> L = E.lseries().dokchitser(Integer(200), algorithm='gp') >>> L.taylor_series(Integer(1),Integer(3)) ...e-82 + (...e-82)*z + 0.75931650028842677023019260789472201907809751649492435158581*z^2 + O(z^3)
Check that Issue #25402 is fixed:
sage: L = EllipticCurve("24a1").modular_form().lseries() sage: L.taylor_series(-1, 3) 0.000000000000000 - 0.702565506265199*z + 0.638929001045535*z^2 + O(z^3)
>>> from sage.all import * >>> L = EllipticCurve("24a1").modular_form().lseries() >>> L.taylor_series(-Integer(1), Integer(3)) 0.000000000000000 - 0.702565506265199*z + 0.638929001045535*z^2 + O(z^3)
Check that Issue #25965 is fixed:
sage: L2 = EllipticCurve("37a1").modular_form().lseries(); L2 L-series associated to the cusp form q - 2*q^2 - 3*q^3 + 2*q^4 - 2*q^5 + O(q^6) sage: L2.taylor_series(0,4) 0.000000000000000 - 0.357620466127498*z + 0.273373112603865*z^2 + 0.303362857047671*z^3 + O(z^4) sage: L2.taylor_series(0,1) O(z^1) sage: L2(0) 0.000000000000000
>>> from sage.all import * >>> L2 = EllipticCurve("37a1").modular_form().lseries(); L2 L-series associated to the cusp form q - 2*q^2 - 3*q^3 + 2*q^4 - 2*q^5 + O(q^6) >>> L2.taylor_series(Integer(0),Integer(4)) 0.000000000000000 - 0.357620466127498*z + 0.273373112603865*z^2 + 0.303362857047671*z^3 + O(z^4) >>> L2.taylor_series(Integer(0),Integer(1)) O(z^1) >>> L2(Integer(0)) 0.000000000000000