CoCalc Shared Filescode / twistedlseries-maarten.sagewsOpen in CoCalc with one click!
Authors: Jennifer Balakrishnan, Alex Best, Maarten Derickx, Maarten Derickx, Michael Rubinstein, Jan Vonk
from psage.eulerprod import LSeriesEllipticCurve, LSeriesTwist
def rep(alpha, primes): vals = [] for p in primes: vals += [p.ideallog(alpha)[0]] return vals def _iter__OK_mod_p(OK,P): M = P.pari_hnf().sage().diagonal() from itertools import product reps = product(*map(range,M)) B = OK.basis() return [sum(i*j for i,j in zip(r,B)) for r in reps] def totally_positive_generator(p): K = p.ring() units = K.unit_group().gens() x = p.gens_reduced() assert len(x)==1, "p not pricipally generated!" x = x[0] M = Matrix(GF(2),[[phi(g)<0 for phi in K.real_embeddings()] for g in units]) v = vector(GF(2),[phi(x)<0 for phi in K.real_embeddings()]) x *= prod(g for g,e in zip(units,M.solve_left(v)) if e) assert x.is_totally_positive(), "Totally positive generator does not exist!" return x class NFChar(): def __init__(self, factored_conductor, rootlist, order, zeta = None, prec = None): """ Currently we allow only tamely ramified characters. """ assert all(p.is_prime() for p in factored_conductor), "factored_conductor contains an element that is not prime" assert len(set(factored_conductor)) == len(factored_conductor), "factored_conductor contains duplicate elements" self._factored_conductor = factored_conductor self._primes_conductor = [p.smallest_integer() for p in factored_conductor] self._conductor = prod(factored_conductor) self._m = totally_positive_generator(self._conductor) self._rootlist = rootlist assert gcd(rootlist+[d])==1, "the order of the character is not %s"%(order) self._order = order self._field = factored_conductor[0].number_field() self._ring = self._field.maximal_order() self._degree = self._field.absolute_degree() if prec == None: self._CC = CC else: self._CC = ComplexField(prec) if zeta is None: self._zeta = self._CC.zeta(d) else: #assert zeta**d==1 self._zeta = zeta assert self._field.class_number()==1 for g in self._field.unit_group().gens(): assert self._call_nf_elt(self._field(g)) == 1, " character is not 1 on units" def __call__(self, p): x = p.gens_reduced() assert len(x)==1, "Do not know how to evaluate on a non principal ideal!" return self._call_nf_elt(p.gens_reduced()[0]) @cached_method def _iter_rep_mod_m(self): return _iter__OK_mod_p(self._ring,self._conductor) def gauss_sum(self): assert self._ring.discriminant()==5, "only implemented for sqrt5" tau = 0 sqrt5 = (self._field*5).factor()[0][0].gens_reduced()[0] print sqrt5 for alpha in self._iter_rep_mod_m(): #print tau a = alpha/self._m#*sqrt5 a = a.trace() den = a.denominator() num = a.numerator() tau += self._call_nf_elt(alpha)*self._CC.zeta(den)^(num%den) return tau def _call_nf_elt(self,alpha): fact = self._factored_conductor rootlist = self._rootlist for p,q in zip(self._primes_conductor,fact): if p.divides(alpha.norm()) and alpha in q: return 0 #return m.ideallog(p.gens_reduced()[0]) reps = rep(alpha, fact) #print "#rootlist %s, #reps %s"%(len(rootlist),len(reps)) return self._zeta**(sum([rootlist[i]*reps[i] for i in range(len(rootlist))]) % self._order) def is_primitive(self): return all(i%self._order != 0 for i in self._rootlist) def conductor(self): return prod(self._factored_conductor).norm() def order(self): return self._order
R.<x> = QQ[] K.<phi>=NumberField(x^2 -x-1) E = EllipticCurve(K, [1, phi + 1, phi, phi, 0]) m = K*13 fact = [prime for (prime, i) in m.factor()] d = 3
foundlist = len(fact)*[1] for t in Integers(d)^len(m.factor()): try: if 0 in list(t): continue rootlist = list(t.lift()) chi = NFChar(fact,rootlist,d,prec=1000) foundlist = rootlist print rootlist x = chi.gauss_sum(); print print x.abs(),x except: pass
[1] -2*phi + 1 13.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 11.3407561726636152699368162869855117983731253191903113812827727652124965276757522368493039894057371978570715000611418213033343214441382821742420971004022136131129716550720474800350052564067483123669351089814797540077428600188606608171977226854509076353420967694562397343116132495426591104271695619234 + 6.35509633539823757082843897576929635164123227418046373191723445222979581773290744537857331917542382091652147604151185486267441580032280296429497772415175454074114827575353571701097120106343964546518212091823462667259056273414807931764641831386718012747116711234367738464348211389670063607384427700001*I [2] -2*phi + 1 13.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 11.3407561726636152699368162869855117983731253191903113812827727652124965276757522368493039894057371978570715000611418213033343214441382821742420971004022136131129716550720474800350052564067483123669351089814797540077428600188606608171977226854509076353420967694562397343116132495426591104271695619234 - 6.35509633539823757082843897576929635164123227418046373191723445222979581773290744537857331917542382091652147604151185486267441580032280296429497772415175454074114827575353571701097120106343964546518212091823462667259056273414807931764641831386718012747116711234367738464348211389670063607384427700001*I
x
11.3407561726636152699368162869855117983731253191903113812827727652124965276757522368493039894057371978570715000611418213033343214441382821742420971004022136131129716550720474800350052564067483123669351089814797540077428600188606608171977226854509076353420967694562397343116132495426591104271695619234 - 6.35509633539823757082843897576929635164123227418046373191723445222979581773290744537857331917542382091652147604151185486267441580032280296429497772415175454074114827575353571701097120106343964546518212091823462667259056273414807931764641831386718012747116711234367738464348211389670063607384427700001*I
CC1 = ComplexField(1000) I = CC1.0 omega = prod([E.period_lattice(phi).basis(prec=1000)[0] for phi in K.embeddings(RR)])
#magma LValues for cubic twist of conductor 13 prec e-74 l1,l2 =[ 1.5457956238192788683557199097437729233665818071610358735968195553567568888501269808125460910379551385036304017986579868697753000481659891333945500875045540441514902685477548634003391025711156928736075 - 0.86622796175514932234718383126453984183577137642394022731205644771173084410073936907261536775497232496651092083801009064197224548531397645425498120806474732224747486052910601468664988353184945233634762*I, 1.5457956238192788683557199097437729233665818071610358735968195553567568888501269808125460910379551385036304017986579868697753000481659891333945500875045540441514902685477548634003391025711156928736077 + 0.86622796175514932234718383126453984183577137642394022731205644771173084410073936907261536775497232496651092083801009064197224548531397645425498120806474732224747486052910601468664988353184945233634753*I ]
13.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
#mamgma LValues for cubic twist of conductor 47 prec e74 l1,l2 = [ 0.24505801487087417016068854563803448172521943510270968636805221995783706098457669244070126558322132754114006237893747704077409789994033154666761763371421995661540622180222286491609840322732421598650475 + 0.42445293274705574412757500887055189552420282612221945972766975611605748380421821057087580584770905794678215958395810249645747226385011392829719555699234058808465069312539939557361013672581656837204405*I, 0.24505801487087417016068854563803448172521943510270968636805221995783706098457669244070126558322132754114006237893747704077409789994033154666761763371421995661540622180222286491609840322732421598650477 - 0.42445293274705574412757500887055189552420282612221945972766975611605748380421821057087580584770905794678215958395810249645747226385011392829719555699234058808465069312539939557361013672581656837204386*I ]
omega = prod([E.period_lattice(phi).basis(prec=1000)[0] for phi in K.embeddings(RR)])
c1 = l1*47*CC1(sqrt(5))/omega c2 = l2*47.conjugate() *CC1(sqrt(5))/omega
c1,c2
(0.9999999995653503062529576591124754896370504363434188055463447835913000245126430422276394574958995879729683943441266440810124106473200871203105837568300491479678574046284886618745443521709166986551736 + 1.732050807586195067941157851961213937336650671437163840656006336307774130379298105089823768284308618777381632838540174758579615723491693438187234643639901697333509474185442000836983390248527316165499*I, 0.9999999995653503062529576591124754896370504363434188055463447835913000245126430422276394574958995879729683943441266440810124106473200871203105837568300491479678574046284886618745443521709166986551737 - 1.732050807586195067941157851961213937336650671437163840656006336307774130379298105089823768284308618777381632838540174758579615723491693438187234643639901697333509474185442000836983390248527316165498*I)
x.conjugate()
11.3407561726636152699368162869855117983731253191903113812827727652124965276757522368493039894057371978570715000611418213033343214441382821742420971004022136131129716550720474800350052564067483123669351089814797540077428600188606608171977226854509076353420967694562397343116132495426591104271695619234 - 6.35509633539823757082843897576929635164123227418046373191723445222979581773290744537857331917542382091652147604151185486267441580032280296429497772415175454074114827575353571701097120106343964546518212091823462667259056273414807931764641831386718012747116711234367738464348211389670063607384427700001*I
1/c1;1/c2;
0.4999999999999999999999999999999999999999999999999999999999999999999999999934956231565666574343088924407419896235212233272735000147214157716184640643318561169074097207246090575533634257256677689681113 - 9.789279086055920531037253731441653548566460090198335763674624714507240225895177402992985655384451524612420671959350614999606032086253312249764159225105263460065631461269354850272811495630508557477052e-75*I 0.2610221927089224350365631712415115144995162664082569178958147777508118762544302932184219095265161745286366069433103338206566987549141885123237098148801496534879660519721023329528268952448575047226335 - 0.4264591597250857028365957388617469633260920057454412677493119493117616729542769917171234683549142379885358787137148975693919254844183663218964443829974242399178036667969166857387601436367253872148488*I
c1, c2;
(0.9999999995653503062529576591124754896370504363434188055463447835913000245126430422276394574958995879729683943441266440810124106473200871203105837568300491479678574046284886618745443521709166986551736 + 1.732050807586195067941157851961213937336650671437163840656006336307774130379298105089823768284308618777381632838540174758579615723491693438187234643639901697333509474185442000836983390248527316165499*I, 0.9999999995653503062529576591124754896370504363434188055463447835913000245126430422276394574958995879729683943441266440810124106473200871203105837568300491479678574046284886618745443521709166986551737 - 1.732050807586195067941157851961213937336650671437163840656006336307774130379298105089823768284308618777381632838540174758579615723491693438187234643639901697333509474185442000836983390248527316165498*I)
c1*c2-4;c1+c2-2
-8.093088569877350958144983871293747909312945882759218380628346650543790566425014413834834802396732744152051493872087143862134447449917591500152645793415664931919724196162966488372064625571333503992959e-10 + 9.144941236411822604101211875669864313850200804989490725343365823468150650983299713086503530538041448730974292603275057833515443317510437438730834983953523435630821104602005456475703659443420968491774e-199*I -8.692993874940846817750490207258991273131623889073104328173999509747139155447210850082008240540632113117467118379751787053598257593788324863399017040642851907430226762509112956581666026896526704936875e-10 + 7.707879042118536194885307152350314207388026392776856468503694051208869834400209758172910118596349221073249760908474691602534445081901654412644560915046541181460263502450261741886664512959454816300209e-199*I
for d in range(1,13): algdep(c1,d,known_digits=10).factor()
1 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4 x^2 - 2*x + 4
for d in range(1,13): algdep(c2,d,known_digits=60).factor()
x - 2 x - 2 x - 2 x - 2 x - 2 x - 2 x - 2 x - 2 x - 2 x - 2 x - 2 x - 2
R.<x> = QQ[] K.<phi>=NumberField(x^2 -x-1) E = EllipticCurve(K, [1, phi + 1, phi, phi, 0]) m = K*47 fact = [prime for (prime, i) in m.factor()] d = 3
foundlist = len(fact)*[1] for t in Integers(d)^len(m.factor()): try: if 0 in list(t): continue rootlist = list(t.lift()) chi = NFChar(fact,rootlist,d) foundlist = rootlist print rootlist x = chi.gauss_sum(); print print x.abs(),x except: pass
[1] -2*phi + 1 47.0000000000000 47.0000000000000 + 1.30007116183606e-13*I [2] -2*phi + 1 46.9999999999999 46.9999999999999 + 3.00870439673417e-14*I
newchi = NFChar(fact, foundlist, d) newchi(K.prime_above(2))
1.00000000000000
x = newchi.gauss_sum(); print x.norm(),x
169.000000000000 11.3407561726636 - 6.35509633539824*I
p = K.prime_above(2) K = p.ring() units = K.unit_group().gens() x = p.gens_reduced() print x assert len(x)==1, "p not pricipally generated!" x = x[0] M = Matrix(GF(2),[[phi(g)<0 for phi in K.real_embeddings()] for g in units]) v = vector(GF(2),[phi(x)<0 for phi in K.real_embeddings()]) print v x *= prod(g for g,e in zip(units,M.solve_left(v)) if e) print x assert x.is_totally_positive(), "Totally positive generator does not exist!"
Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.5/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> NameError: name 'K' is not defined
K
Error in lines 1-1 Traceback (most recent call last): File "/projects/sage/sage-7.5/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> NameError: name 'K' is not defined
M.solve_left(v)
(0, 1)
v;x
(1, 0) 2*phi
totally_positive_generator??
File: /projects/68c8b2b8-03ba-44d4-a0d1-5d771c8cb465/code Unable to read source code (source code not available)
totally_positive_generator(K.prime_above(2))
2 (0, 0) (0, 0) 2
L = LSeriesEllipticCurve(E) print "done1" LEchi = LSeriesTwist(L, newchi, epsilon='solve', conductor=ZZ(5^2*E.conductor().norm()*(newchi.conductor()^2))) print E.conductor().norm() print newchi.conductor() LEchi print "done2"
1.00000000000000 done1 31 169 Twist of L-series of Elliptic Curve defined by y^2 + x*y + phi*y = x^3 + (phi+1)*x^2 + phi*x over Number Field in phi with defining polynomial x^2 - x - 1 by <__builtin__.NFChar instance at 0x7fd03418dd40> done2
LEchi.anlist(100) L.conductor() L.epsilon() L.number_of_coefficients(2) L.degree() L.hodge_numbers() L(3/2) LEchi.number_of_coefficients(5) LEchi.degree() LEchi.hodge_numbers() LEchi(ComplexField(5)(3/2))
[0, 1, -0.000000000000000, -0.000000000000000, -3.00000000000000, 2.00000000000000, 0.000000000000000, -0.000000000000000, -0.000000000000000, 2.00000000000000, -0.000000000000000, 0.000000000000000, 0.000000000000000, 0, 0.000000000000000, -0.000000000000000, 5.00000000000000, 0, -0.000000000000000, 0.000000000000000, -6.00000000000000, 0.000000000000000, 0.000000000000000, 0, 0.000000000000000, -1.00000000000000, 0.000000000000000, -0.000000000000000, 0.000000000000000, -4.00000000000000, 0.000000000000000, -7.00000000000000, -0.000000000000000, 0.000000000000000, 0.000000000000000, -0.000000000000000, -6.00000000000000, 0, 0.000000000000000, 0.000000000000000, -0.000000000000000, 12.0000000000000, 0.000000000000000, 0, 0.000000000000000, 4.00000000000000, 0.000000000000000, 0, -0.000000000000000, 2.00000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, -8.00000000000000, 0.000000000000000, 4.00000000000000, 0.000000000000000, -0.000000000000000, -3.00000000000000, 0.000000000000000, 0.000000000000000, 0, 0.000000000000000, 0.000000000000000, 0.000000000000000, 8.00000000000000, -0.000000000000000, 0, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 16.0000000000000, 10.0000000000000, -5.00000000000000, -0.000000000000000, 0, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, -4.00000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0, -0.000000000000000, 0.000000000000000, 3.00000000000000] 775 1 534 2 [0, 0, 1, 1] 0.607704196335035 99729 2 [0, 0, 1, 1]
#better def is_prime(n): """returns True if n is a prime false otherwise""" if n<2: return False for i in range(2,floor(sqrt(n))+1): if n%i == 0: return False return True def primes_in_range(a,b): """returns the list of all primes in the range(a,b)""" return [i for i in range(a,b) if is_prime(i)] print primes_in_range(2,20)
[2, 3, 5, 7, 11, 13, 17, 19]
L = LSeriesEllipticCurve(E) import cProfile t = cputime() cProfile.runctx("L.anlist(100000)",None, locals(), filename="timing_data/profile_test", sort='cumulative') print cputime(t)
478.148 478.148
LEchi = LSeriesTwist(L, newchi, epsilon='solve', conductor=ZZ(5^2*E.conductor().norm()*(newchi.conductor()^2))) t = cputime() cProfile.runctx("LEchi.anlist(100000)",None, locals(), filename="timing_data/profile_test1", sort='cumulative') print cputime(t)
26.52
print cputime(t)
L = LSeriesEllipticCurve(E) L.anlist(100) print "L done" LEchi = LSeriesTwist(L, newchi, epsilon='solve', conductor=ZZ(5^2*E.conductor().norm()*(newchi.conductor()^2))) LEchi.anlist(100)
%sh python gprof2dot/gprof2dot.py -f pstats timing_data/profile_test | dot -Tpng -o timing_data/profile_LSeriesEllipticCurve_anlist100000_cachetest_opt.png 2> /dev/null
%sh python gprof2dot/gprof2dot.py -f pstats timing_data/profile_test1 | dot -Tpng -o timing_data/profile_LSeriesTwist_anlist100000_cachtest_opt.png 2> /dev/null
L
L-series of Elliptic Curve defined by y^2 + x*y + phi*y = x^3 + (phi+1)*x^2 + phi*x over Number Field in phi with defining polynomial x^2 - x - 1
96c00861-087e-464f-b834-00be02de9d9f