CoCalc Shared Filescode / non-primitive.sagewsOpen in CoCalc with one click!
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)]