CoCalc Shared Filescode / non-primitive.sagews
Author: Maarten Derickx


load("compute_lambda.sage")

Compiling ./modular_symbol_map.pyx...
E = EllipticCurve("37a1")
M = E.modular_symbol_space(sign=1)
f = ModularSymbolMap(M)
gn = E.modular_symbol_numerical()
g = E.modular_symbol()
h = ModularSymbolMap(g)
inf_zero = M.rational_period_mapping()([oo,0])[0]

1
a=compute_dist(E,3,100000)
h = plot_histogram(a)

1 0 There are 4783 primes to use up to 100000 Starting... X X X X X X X X X X
a1=compute_dist1(E,3,20000)
h1 = plot_histogram(a1,color='red')

There are 1123 primes to use up to 20000 Starting... X X X X X X X X X X
Error in lines 2-2 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> TypeError: plot_histogram() got an unexpected keyword argument 'color'
h = histogram(a,bins=200,color = 'red');h1 = histogram([i/sqrt(3.5) for i in a1[:len(a)]],bins=200)

h+h1;h1+h

h.plot?

File: /projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/plot/graphics.py
Signature : h.plot(self)
Docstring :
Draw a 2D plot of this graphics object, which just returns this
object since this is already a 2D graphics object.

EXAMPLES:

sage: S = circle((0,0), 2)
sage: S.plot() is S
True

It does not accept any argument (https://trac.sagemath.org/19539):

sage: S.plot(1)
Traceback (most recent call last):
...
TypeError: plot() takes exactly 1 argument (2 given)
sage: S.plot(hey="hou")
Traceback (most recent call last):
...
TypeError: plot() got an unexpected keyword argument 'hey'

a=compute_dist(E,5,10000)
plot_histogram(a)

There are 306 primes to use up to 10000 Starting... X X X X X X X X X X
a=compute_dist(E,7,10000)
plot_histogram(a)

There are 203 primes to use up to 10000 Starting... X X X
Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> File "", line 32, in compute_dist File "", line 21, in alphas File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/schemes/elliptic_curves/ell_modular_symbols.py", line 349, in __call__ if r != oo: File "sage/structure/element.pyx", line 1125, in sage.structure.element.Element.__richcmp__ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10384) return coercion_model.richcmp(self, other, op) File "sage/structure/coerce.pyx", line 1868, in sage.structure.coerce.CoercionModel_cache_maps.richcmp (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/coerce.c:20376) return PyObject_RichCompare(x, y, op) File "sage/structure/element.pyx", line 1123, in sage.structure.element.Element.__richcmp__ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10360) return (<Element>self)._richcmp_(other, op) File "sage/structure/element.pyx", line 1127, in sage.structure.element.Element._richcmp_ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10462) cpdef _richcmp_(left, right, int op): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/modular/cusps.py", line 444, in _richcmp_ return richcmp(s, o, op) File "stringsource", line 67, in cfunc.to_py.__Pyx_CFunc_object____object____object____int___to_py.wrap (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/sage_object.c:18617) File "sage/structure/element.pyx", line 1125, in sage.structure.element.Element.__richcmp__ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10384) return coercion_model.richcmp(self, other, op) File "sage/structure/coerce.pyx", line 1860, in sage.structure.coerce.CoercionModel_cache_maps.richcmp (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/coerce.c:20259) x, y = self.canonical_coercion(x, y) File "sage/structure/coerce.pyx", line 1167, in sage.structure.coerce.CoercionModel_cache_maps.canonical_coercion (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/coerce.c:10890) x_elt = (<Map>x_map)._call_(x) File "sage/structure/coerce_maps.pyx", line 105, in sage.structure.coerce_maps.DefaultConvertMap_unique._call_ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/coerce_maps.c:4757) return C._element_constructor(x) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/rings/infinity.py", line 1181, in _element_constructor_ return FiniteNumber(self, cmp(x, 0)) File "src/cysignals/signals.pyx", line 252, in cysignals.signals.python_check_interrupt (build/src/cysignals/signals.c:2854) File "src/cysignals/signals.pyx", line 97, in cysignals.signals.sig_raise_exception (build/src/cysignals/signals.c:1303) KeyboardInterrupt
a=compute_dist(E,11,10000)
plot_histogram(a)

There are 945 primes to use up to 100000 Starting... X X X X X X X X X X
a=compute_dist(E,13,10000)
plot_histogram(a)

There are 1008 primes to use up to 130000 Starting... X X X
Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> File "", line 32, in compute_dist File "", line 21, in alphas File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/schemes/elliptic_curves/ell_modular_symbols.py", line 349, in __call__ if r != oo: File "sage/structure/element.pyx", line 1125, in sage.structure.element.Element.__richcmp__ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10384) return coercion_model.richcmp(self, other, op) File "sage/structure/coerce.pyx", line 1868, in sage.structure.coerce.CoercionModel_cache_maps.richcmp (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/coerce.c:20376) return PyObject_RichCompare(x, y, op) File "sage/structure/element.pyx", line 1123, in sage.structure.element.Element.__richcmp__ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10360) return (<Element>self)._richcmp_(other, op) File "sage/structure/element.pyx", line 1127, in sage.structure.element.Element._richcmp_ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10462) cpdef _richcmp_(left, right, int op): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/modular/cusps.py", line 444, in _richcmp_ return richcmp(s, o, op) File "stringsource", line 67, in cfunc.to_py.__Pyx_CFunc_object____object____object____int___to_py.wrap (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/sage_object.c:18617) File "sage/structure/element.pyx", line 1125, in sage.structure.element.Element.__richcmp__ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10384) return coercion_model.richcmp(self, other, op) File "sage/structure/coerce.pyx", line 1868, in sage.structure.coerce.CoercionModel_cache_maps.richcmp (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/coerce.c:20376) return PyObject_RichCompare(x, y, op) File "sage/structure/element.pyx", line 1123, in sage.structure.element.Element.__richcmp__ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10360) return (<Element>self)._richcmp_(other, op) File "sage/structure/element.pyx", line 1158, in sage.structure.element.Element._richcmp_ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10559) c = left._cmp_(right) File "sage/structure/element.pyx", line 1174, in sage.structure.element.Element._cmp_ (/projects/sage/sage-7.6/src/build/cythonized/sage/structure/element.c:10948) return left_cmp(right) File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/sage/rings/infinity.py", line 1272, in __cmp__ def __cmp__(self, other): File "src/cysignals/signals.pyx", line 252, in cysignals.signals.python_check_interrupt (build/src/cysignals/signals.c:2854) File "src/cysignals/signals.pyx", line 97, in cysignals.signals.sig_raise_exception (build/src/cysignals/signals.c:1303) KeyboardInterrupt
def compute_dist1(E, d, m_max):
M = E.modular_symbol_space(sign=1)
assert d%2 == 1
N = M.level()
f = E.modular_symbol()

# much more ms, since this code is massively faster...
ms = [m for m in prime_range(3,m_max+1) if \
gcd(m, N) == 1 and euler_phi(m) % d == 0]

print 'There are %s primes to use up to %s'%(len(ms), m_max)

def alphas(m, d):
assert d%2 == 1
R = Integers(m)
gen = R(primitive_root(m))
n = euler_phi(m)//d
b = gen
h = gen^d
denom = float(sqrt(euler_phi(m)*log(m)))
alphas = []
for i in range(0, d):
s = 0
for j in range(n):
period = f((b^i * h^j).lift()/ m)
s += period
alphas.append(s / denom)
return [alphas[i]-alphas[j] for i in range(d) for j in range(d) if i!=j]

t0 = walltime()
print "Starting..."
data = []
progress = len(ms) // 10 + 1
for i, m in enumerate(ms):
if i % progress == 0:
print 'X',; sys.stdout.flush()
data += alphas(m, d)
return stats.TimeSeries(data)

%time E.modular_symbol_numerical()(1/2)

-0.800000000000001 CPU time: 0.04 s, Wall time: 0.04 s
f._eval1(3,25)[0]/ZZ(f.denom)+inf_zero;g(3,25)

-3/2 1/5
randint(0,N)

2287724
N=15
for i in srange(N):
print "test"
f._eval1(i,N)[0]/25/f.denom,g(i/N)+g(-i/N)

test (0, 2/5) test (1/10, 7/5) test (1/10, 7/5) test (1/5, 12/5) test (1/10, 7/5) test (-1/10, -3/5) test (-3/10, -13/5) test (-2/5, -18/5) test (-2/5, -18/5) test (-3/10, -13/5) test (-1/10, -3/5) test (1/10, 7/5) test (1/5, 12/5) test (1/10, 7/5) test (1/10, 7/5)
g(3/25)

-3/10
gn(3/7)

-1.80000000000001
%time f._eval1(1,2)

[-10] CPU time: 0.00 s, Wall time: 0.00 s
f.denom*f.C;


[ -2] [ 2] [ 0] [ 10] [ 5] [ -5] [-10] [-10] [ -5] [ 5] [ 10] [ 0]
f(1/4)

Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.6/local/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 995, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> TypeError: '_projects_68c8b2b8_03ba_44d4_a0d1_5d771c8cb465_code_modular_symbol_map_pyx._projects_68c8b2b8_03ba_44d4_a0d1_5d771c8cb465_code_modular_symbol_map_pyx_0.ModularSymbolMap' object is not callable
ModularSymbolMap??

   File:
Source:
cdef class ModularSymbolMap:
cdef long d, N
cdef public long denom
cdef long* X  # coefficients of linear map from P^1 to Q^d.
cdef public object C
cdef public P1List P1

def __cinit__(self):
self.X = NULL

def __repr__(self):
return "Modular symbols map for modular symbols factor of dimension %s and level %s"%(self.d, self.N)

def __init__(self, A):
"""
EXAMPLES::

We illustrate a few ways to construct the modular symbol map for 389a.

sage: from psage.modform.rational.modular_symbol_map import ModularSymbolMap
sage: A = ModularSymbols(389,sign=1).cuspidal_subspace().new_subspace().decomposition()[0]
sage: f = ModularSymbolMap(A)
sage: f._eval1(-3,7)
[-2]
sage: f.denom
2

sage: E = EllipticCurve('389a'); g = E.modular_symbol()
sage: h = ModularSymbolMap(g)
sage: h._eval1(-3,7)
[-2]
sage: h.denom
1
sage: g(-3/7)
-2

sage: f.denom*f.C == h.denom*h.C
True
"""
import sage.schemes.elliptic_curves.ell_modular_symbols as ell_modular_symbols
from sage.modular.modsym.space import ModularSymbolsSpace

if isinstance(A, ell_modular_symbols.ModularSymbol):
C = self._init_using_ell_modular_symbol(A)

elif isinstance(A, ModularSymbolsSpace):
C = self._init_using_modsym_space(A)

else:

raise NotImplementedError, "Creating modular symbols from object of type %s not implemented"%type(A)

# useful for debugging only -- otherwise is a waste of memory/space
self.C = C

X, self.denom = C._clear_denom()
# Now store in a C data structure the entries of X (as long's)
self.X = <long*>sage_malloc(sizeof(long*)*X.nrows()*X.ncols())
cdef Py_ssize_t i, j, n
n = 0
for a in X.list():
self.X[n] = a
n += 1

def _init_using_ell_modular_symbol(self, f):
# Input f is an elliptic curve modular symbol map
assert f.sign() != 0
self.d = 1  # the dimension

E = f.elliptic_curve()
self.N = E.conductor()
self.P1 = P1List(self.N)

# Make a matrix whose rows are the images of the Manin symbols
# corresponding to the elements of P^1 under f.
n = len(self.P1)
C = matrix(QQ, n, 1)
for i in range(n):
# If the ith element of P1 is (u,v) for some u,v modulo N.
# We need to turn this into something we can evaluate f on.
# 1. Lift to a 2x2 SL2Z matrix M=(a,b;c,d) with c=u, d=v (mod N).
# 2. The modular symbol is M{0,oo} = {b/d,a/c}.
# 3. So {b/d, a/c} = {b/d,oo} + {oo,a/c} = -{oo,b/d} + {oo,a/c} = f(a/c)-f(b/d).
# 4. Thus x |--> f(a/c)-f(b/d).
a,b,c,d = self.P1.lift_to_sl2z(i)  # output are Python ints, so careful!
C[i,0] = (f(Integer(a)/c) if c else 0) - (f(Integer(b)/d) if d else 0)
return C

def _init_using_modsym_space(self, A):
# Very slow generic setup code.  This is "O(1)" in that we
# care mainly about evaluation time being fast, at least in
# this code.
if A.sign() == 0:
raise ValueError, "A must have sign +1 or -1"
self.d = A.dimension()

self.N = A.level()
self.P1 = P1List(self.N)

# The key data we need from the modular symbols space is the
# map that assigns to an element of P1 the corresponding
# element of ZZ^n.  That's it.  We forget everything else.
M = A.ambient_module()
W = matrix([M.manin_symbol(x).element() for x in self.P1])
B = A.dual_free_module().basis_matrix().transpose()

# Return matrix whose rows are the images of the Manin symbols
# corresponding to the elements of P^1 under the modular
# symbol map.
return W*B

def __dealloc__(self):
if self.X:
sage_free(self.X)

cdef int evaluate(self, long v[MAX_DEG], long a, long b) except -1:
cdef long q[MAX_CONTFRAC]
cdef int i, j, k, n, sign=1
cdef long* x

# initialize answer vector to 0
for i in range(self.d):
v[i] = 0

# compute continued fraction
n = contfrac_q(q, a, b)

# compute corresponding modular symbols, mapping over...
for i in range(1,n):
j = self.P1.index((sign*q[i])%self.N, q[i-1]%self.N)
# map over, adding a row of the matrix self.X
# to the answer vector v.
x = self.X + j*self.d
for k in range(self.d):
v[k] += x[k]
# change sign, so q[i] is multiplied by (-1)^(i-1)
sign *= -1

def dimension(self):
return self.d

def _eval0(self, a, b):
cdef long v[MAX_DEG]
self.evaluate(v, a, b)

def _eval1(self, a, b):
"""
EXAMPLE::

sage: from psage.modform.rational.modular_symbol_map import ModularSymbolMap
sage: A = ModularSymbols(188,sign=1).cuspidal_subspace().new_subspace().decomposition()[-1]
sage: f = ModularSymbolMap(A)
sage: f._eval1(-3,7)
[-3, 0]

"""
cdef long v[MAX_DEG]
self.evaluate(v, a, b)
cdef int i
return [v[i] for i in range(self.d)]