CoCalc Shared Filescode / maarten-scratch.sagews
Authors: Alex Best, Maarten Derickx
%load modular_symbol_map.pyx

Compiling ./modular_symbol_map.pyx...
sage_server.MAX_OUTPUT = 1000000
for mf in ["431a", "431b", "433a"]:
fp=open(mf,"w")
M =  EllipticCurve(mf+"1").modular_symbol_space(sign=1)
MSinfz = M.rational_period_mapping()([oo,0])[0]
d = M.dimension()
assert d == 1
f = ModularSymbolMap(M)
X,MSden = f.C._clear_denom()
X=[a[0] for a in X]
N = M.level()
P1 = f.P1
MStab = [[X[P1.index(a,b)] for b in Integers(N)] for a in Integers(N)]
fp.write("int MStab_%s[] = {\n"%mf)
for i in range(N):
fp.write("%d,"%MStab[i][0])
for j in range(1,N):
fp.write("%d,"%MStab[i][j])
fp.write("\n")
fp.write("};\n\n")
fp.write("int MSden_%s = %d;\n"%(mf,MSden))
fp.write("int MSinfz_%s = %d;\n"%(mf,MSinfz))
fp.write("\n")
fp.close()

def contfrac_q(a, b):
"""
Given coprime integers a, b, compute the sequence q_i of
denominators in the continued fraction expansion of a/b
"""
qi = [1]
if a == 0:
return qi
if b == 0:
raise ZeroDivisionError

# one iteration isn't needed, since we are only computing the q_i.
a,b = b, a%b
i=1
while b:
if i >= 2:
qi.append(qi[-1]*(a//b) + qi[-2])
else:
qi.append(a//b)
a,b = b, a%b
i += 1
return qi

def eval(a,b,tab):
q = contfrac_q(a,b)
v = 0
sign = 1
for i in range(1,len(q)):
v += tab[(sign*q[i])%N][q[i-1]%N]
sign *= -1
return v

def alphas(m,d,N,tab):
assert is_prime(m) and d%2 == 1
R = Integers(m)
Npow = R(N)^((d-1)//2)
gen = R(primitive_root(m))
n = (m-1)//d
b = gen
h = gen^d
return [sum([eval((Npow*b^i*h^j).lift(),m,tab) for j in range(n)]) for i in range(1,(d-1)//2+1)]

def normalize_alphas(alist,m,MSden,MSinfz):
# should be normalized sqrt(euler-phi(m)/d * C_E + D_E), but computing C_E and D_E is a pain, and d is fixed
return [(a/MSden+((m-1)//d)*MSinfz)/float(sqrt((m-1)*log(m))) for a in alist]

MOMENTS=10

def moments(alist):
M=[0.0 for i in range(MOMENTS+1)]
for a in alist:
apow = 1.0
for i in range(MOMENTS+1):
M[i]+=apow
apow *= a
return [M[i]/len(alist) for i in range(MOMENTS+1)]


d = 3
ms = [m for m in prime_range(3,100) if gcd(m,N)==1 and m % d == 1]
a = []
for m in ms:
alist = alphas(m,d,N,MStab)
print m, alist, normalize_alphas(alist,m,MSden,MSinfz)
a += normalize_alphas(alist,m,MSden,MSinfz)
print len(a), moments(a)


7 [2] [0.5853194810641169] 13 [-2] [-0.3604956057297047] 19 [2] [0.27472138285709213] 31 [4] [0.39409418584560885] 37 [-8] [-0.7016651856311937] 43 [2] [0.1591263605933982] 61 [4] [0.25469316510042145] 67 [-2] [-0.12005794395973912] 73 [-2] [-0.11379203050314328] 79 [6] [0.32500546873918773] 97 [20] [0.9543607031040101] 11 [1.00000000000000, 0.150119089225459, 0.211787029485880, 0.0737345846554369, 0.113856308766378, 0.0636422712804702, 0.0839059400262076, 0.0601951443657269, 0.0692500441872164, 0.0567079792107614, 0.0600532480700175]
%time
a=compute_dist(EllipticCurve("37a1"),29,1000000)
plot_histogram(a)
print ""
print "%d:"%len(a),
M = moments(a)
for i in range(1,MOMENTS+1):
print "%.6f "%M[i],


2 1 There are 2816 primes to use up to 1000000 Starting... X X X X X X X X X X X
39424: 0.000072 0.636646 0.000728 1.645757 -0.010270 9.505303 1.751655 111.743050 109.795762 2240.739724 CPU time: 2827.92 s, Wall time: 2828.03 s
capture?


File: /projects/sage/sage-7.5/local/lib/python2.7/site-packages/smc_sagews/sage_salvus.py
Signature : capture(self, code=None, stdout=None, stderr=None, append=False, echo=False)
Docstring :
Capture or ignore the output from evaluating the given code.
(SALVUS only).

Use capture as a block decorator by placing either %capture or
%capture(optional args) at the beginning of a cell or at the
beginning of a line.  If you use just plain %capture then stdout
and stderr are completely ignored.  If you use %capture(args) you
can redirect or echo stdout and stderr to variables or files.  For
example if you start a cell with this line:

%capture(stdout='output', stderr=open('error','w'), append=True, echo=True)

then stdout is appended (because append=True) to the global
variable output, stderr is written to the file 'error', and the
output is still displayed in the output portion of the cell
(echo=True).

INPUT:

* stdout -- string (or object with write method) to send stdout
output to (string=name of variable)

* stderr -- string (or object with write method) to send stderr
output to (string=name of variable)

* append -- (default: False) if stdout/stderr are a string,
append to corresponding variable

* echo   -- (default: False) if True, also echo stdout/stderr to
the output cell.



File: /projects/sage/sage-7.5/local/lib/python2.7/site-packages/smc_sagews/sage_server.py
Docstring :
SageMathCloud functionality.

An instance of this object is created each time you execute a cell.
It has various methods for sending different types of output
messages, links to files, etc.   Type 'help(smc)' for more details.

OUTPUT LIMITATIONS -- There is an absolute limit on the number of
messages output for a given cell, and also the size of the output
message for each cell.  You can access or change those limits
dynamically in a worksheet as follows by viewing or changing any of
the following variables:

sage_server.MAX_STDOUT_SIZE       # max length of each stdout output message
sage_server.MAX_STDERR_SIZE       # max length of each stderr output message
sage_server.MAX_MD_SIZE           # max length of each md (markdown) output message
sage_server.MAX_HTML_SIZE         # max length of each html output message
sage_server.MAX_TEX_SIZE          # max length of tex output message
sage_server.MAX_OUTPUT_MESSAGES   # max number of messages output for a cell.

And:

sage_server.MAX_OUTPUT            # max total character output for a single cell; computation
# terminated/truncated if sum of above exceeds this.

"%s/431a.txt"%os.environ['HOME']

'/projects/68c8b2b8-03ba-44d4-a0d1-5d771c8cb465/431a.txt'
fp=open("foo.txt","w")
fp.write("hi")
fp.close()

︠d47c1e60-002b-477f-b5c9-42a39c7b7711︠