Utilities for Gross-Zagier L-series¶
- sage.modular.modform.l_series_gross_zagier_coeffs.bqf_theta_series(Q, bound, var=None)¶
Return the theta series associated to a positive definite quadratic form.
For a given form \(f = ax^2 + bxy + cy^2\) this is the sum
\[\sum_{(x,y) \in \ZZ^2} q^{f(x,y)} = \sum_{n=-\infty}^{\infty} r(n)q^n\]where \(r(n)\) give the number of way \(n\) is represented by \(f\).
INPUT:
Q– a positive definite quadratic formbound– how many terms to computevar– (optional) the variable in which to express this power series
OUTPUT: a power series in
var, or list of ints ifvaris unspecifiedEXAMPLES:
sage: from sage.modular.modform.l_series_gross_zagier_coeffs import bqf_theta_series sage: bqf_theta_series([2,1,5], 10) [1, 0, 2, 0, 0, 2, 2, 0, 4, 0, 0] sage: Q = BinaryQF([41,1,1]) sage: bqf_theta_series(Q, 50, ZZ[['q']].gen()) 1 + 2*q + 2*q^4 + 2*q^9 + 2*q^16 + 2*q^25 + 2*q^36 + 4*q^41 + 4*q^43 + 4*q^47 + 2*q^49 + O(q^51)
>>> from sage.all import * >>> from sage.modular.modform.l_series_gross_zagier_coeffs import bqf_theta_series >>> bqf_theta_series([Integer(2),Integer(1),Integer(5)], Integer(10)) [1, 0, 2, 0, 0, 2, 2, 0, 4, 0, 0] >>> Q = BinaryQF([Integer(41),Integer(1),Integer(1)]) >>> bqf_theta_series(Q, Integer(50), ZZ[['q']].gen()) 1 + 2*q + 2*q^4 + 2*q^9 + 2*q^16 + 2*q^25 + 2*q^36 + 4*q^41 + 4*q^43 + 4*q^47 + 2*q^49 + O(q^51)
- sage.modular.modform.l_series_gross_zagier_coeffs.gross_zagier_L_series(an_list, Q, N, u, var=None)¶
Compute the coefficients of the Gross-Zagier \(L\)-series.
INPUT:
an_list– list of coefficients of the \(L\)-series of an elliptic curveQ– a positive definite quadratic formN– conductor of the elliptic curveu– number of roots of unity in the field associated withQvar– (optional) the variable in which to express this power series
OUTPUT: a power series in
var, or list of ints ifvaris unspecifiedThe number of terms is the length of the input
an_list.EXAMPLES:
sage: from sage.modular.modform.l_series_gross_zagier_coeffs import gross_zagier_L_series sage: E = EllipticCurve('37a') sage: N = 37 sage: an = E.anlist(60); len(an) 61 sage: K.<a> = QuadraticField(-40) sage: A = K.class_group().gen(0) sage: Q = A.ideal().quadratic_form().reduced_form() sage: Q2 = (A**2).ideal().quadratic_form().reduced_form() sage: u = K.zeta_order() sage: t = PowerSeriesRing(ZZ,'t').gen() sage: LA = gross_zagier_L_series(an,Q,N,u,t); LA -2*t^2 - 2*t^5 - 2*t^7 - 4*t^13 - 6*t^18 - 4*t^20 + 20*t^22 + 4*t^23 - 4*t^28 + 8*t^32 - 2*t^37 - 6*t^45 - 18*t^47 + 2*t^50 - 8*t^52 + 2*t^53 + 20*t^55 + O(t^61) sage: len(gross_zagier_L_series(an,Q,N,u)) 61 sage: LA + gross_zagier_L_series(an,Q2,N,u,t) t - 2*t^2 + 2*t^4 - 2*t^5 - 2*t^7 + 3*t^9 + 4*t^10 - 10*t^11 - 4*t^13 + 4*t^14 - 4*t^16 - 6*t^18 - 4*t^20 + 20*t^22 + 4*t^23 - t^25 + 8*t^26 - 4*t^28 + 8*t^32 + 4*t^35 + 6*t^36 - 2*t^37 - 18*t^41 - 20*t^44 - 6*t^45 - 8*t^46 - 18*t^47 - 11*t^49 + 2*t^50 - 8*t^52 + 2*t^53 + 20*t^55 + 16*t^59 + O(t^61)
>>> from sage.all import * >>> from sage.modular.modform.l_series_gross_zagier_coeffs import gross_zagier_L_series >>> E = EllipticCurve('37a') >>> N = Integer(37) >>> an = E.anlist(Integer(60)); len(an) 61 >>> K = QuadraticField(-Integer(40), names=('a',)); (a,) = K._first_ngens(1) >>> A = K.class_group().gen(Integer(0)) >>> Q = A.ideal().quadratic_form().reduced_form() >>> Q2 = (A**Integer(2)).ideal().quadratic_form().reduced_form() >>> u = K.zeta_order() >>> t = PowerSeriesRing(ZZ,'t').gen() >>> LA = gross_zagier_L_series(an,Q,N,u,t); LA -2*t^2 - 2*t^5 - 2*t^7 - 4*t^13 - 6*t^18 - 4*t^20 + 20*t^22 + 4*t^23 - 4*t^28 + 8*t^32 - 2*t^37 - 6*t^45 - 18*t^47 + 2*t^50 - 8*t^52 + 2*t^53 + 20*t^55 + O(t^61) >>> len(gross_zagier_L_series(an,Q,N,u)) 61 >>> LA + gross_zagier_L_series(an,Q2,N,u,t) t - 2*t^2 + 2*t^4 - 2*t^5 - 2*t^7 + 3*t^9 + 4*t^10 - 10*t^11 - 4*t^13 + 4*t^14 - 4*t^16 - 6*t^18 - 4*t^20 + 20*t^22 + 4*t^23 - t^25 + 8*t^26 - 4*t^28 + 8*t^32 + 4*t^35 + 6*t^36 - 2*t^37 - 18*t^41 - 20*t^44 - 6*t^45 - 8*t^46 - 18*t^47 - 11*t^49 + 2*t^50 - 8*t^52 + 2*t^53 + 20*t^55 + 16*t^59 + O(t^61)
- sage.modular.modform.l_series_gross_zagier_coeffs.to_series(L, var)¶
Create a power series element out of a list
Lin the variablevar.EXAMPLES:
sage: from sage.modular.modform.l_series_gross_zagier_coeffs import to_series sage: to_series([1,10,100], 't') 1 + 10*t + 100*t^2 + O(t^3) sage: to_series([0..5], CDF[['z']].0) 0.0 + 1.0*z + 2.0*z^2 + 3.0*z^3 + 4.0*z^4 + 5.0*z^5 + O(z^6)
>>> from sage.all import * >>> from sage.modular.modform.l_series_gross_zagier_coeffs import to_series >>> to_series([Integer(1),Integer(10),Integer(100)], 't') 1 + 10*t + 100*t^2 + O(t^3) >>> to_series((ellipsis_range(Integer(0),Ellipsis,Integer(5))), CDF[['z']].gen(0)) 0.0 + 1.0*z + 2.0*z^2 + 3.0*z^3 + 4.0*z^4 + 5.0*z^5 + O(z^6)