# Isogeny computation for larger numbers doesn't finish

I have the following code segment in Sage:

from sage.all import *

proof.arithmetic(False)

f = 1
lA = 2
lB = 3
eA = 372
eB = 239

p = f*lA**eA*lB**eB-1
assert p.is_prime()

Fp = GF(p)
R.<x> = PolynomialRing(Fp)
# The quadratic extension via x^2 + 1 since p = 3 mod 4
Fp2.<j> = Fp.extension(x^2 + 1)

E_0 = EllipticCurve(Fp2, [1,0])
assert E_0.is_supersingular()

l = [lA, lB]
e = [eA, eB]

def starting_node(E_i, count):
if count == 96:
return E_i
kers = E_i(0).division_points(3)
P = E_i(0)
while P == E_i(0):
P = kers[randint(0, 8)]
psi = E_i.isogeny(P)
return starting_node(psi.codomain(), count+1)

def walk(E_i, j,lvl, prev_j ,max_lvl):
children = {}
if lvl < max_lvl:
for P in E_i(0).division_points(2):
if P == E_i(0):
continue
psi = E_i.isogeny(P)
E_child = psi.codomain()
j_child = str(E_child.j_invariant())
lvl_child = lvl+1
if j_child in children:
children[j_child] += 1
else:
children[j_child] = 1
if (j_child != prev_j) or (children[j_child] > 1):
walk(E_child, j_child, lvl_child,j, max_lvl)

idx = 0
E = starting_node(E_0, 0)
print "start: ", str(E.j_invariant())
walk(E, str(E.j_invariant()), 0, "", e[idx])


The above code seems to execute fine when I use smaller values, but still it executes very slowly. Whereas, when I use the above numbers it does not finish with the execution. There already seems to be a similar problem mentioned here, but the solution there seems to be just for 3-isogenies, but what about the general case, for different isogenies and using kernels? Also Magma seems to handle isogeny computations quite efficiently, why is Sage having such a hard time? Any ideas how can I modify the above code so that it actually finishes the execution in a "reasonable" time?

edit retag close merge delete

Sort by » oldest newest most voted

There was a related question here:

41236

So i will use the ad-hoc routine pointIsogeny, the curve E, and point $R$ defined in there for the start. If i understand the situation, then the line E = starting_node(E_0, 0) wants to build (randomly) the following:

p = 2^372*3^239 - 1
F = GF( p, proof=0 )
R.<x> = PolynomialRing(F)
# The quadratic extension via x^2 + 1 since p = 3 mod 4
Fp2.<j> = Fp.extension(x^2 + 1)

E0 = EllipticCurve(Fp2, [1,0])
R = E0.point( ( 9007518108646169717637829256143902727256908604612852170262845383236308734546752870948023665818612994373405439104311563180515971827416888758364379345147971116263603311381076594736408857657724917306603115510356333363208849059629*j + 8143875544876203102731118758998042519548033956324586056599014675159430178641351639084698635604334996400279245810044652728374901773305503117205094200107841651156165349992200320758569566012680521517849444975619314122642286738078, 5933067621852699133314119054291511797259450704514751366342623924502189539964363940282453138760679452407770908256363554867263728524097776132882938143129922521585626385240622607283870971683009886348379499590584807528366215593257*j + 3243060684426863808401390460086135176972427334603406644358264254553792355536601776050361287848102972929373427788393162068323805718475829130551068182582167101080584984966674872491237953303465307684436948645319932524248022850554 ) )

R_list, degrees_list = [], []
S = 3^(239-96)*R
while S != E0(0):
R_list.append( S )
degrees_list.append(3)
S *= 3

X = ( E0, R_list, [], degrees_list )
count = 0
while X[1]:
count += 1
print count
X = pointIsogeny( *X )

E = X[0]
E


So we have a starting reproducible elliptic curve E

Elliptic Curve defined by y^2 = x^3 + (7365844174710734349703979824267397716070848379207150714853939305480598177974648401456642446481782170109092245224349222954361643535826935062105068360613554750604758205242084022642407354099811489418896293593385968974015370653730*j+211392337950744742300479627443770128587762706472203430775009140876045287163506520387404594567561228376780127226866671410674341193885329801042247857950441952859225536775626695894662031453971315468618629095593065757279048175980)*x + (191888858730120017972798004127909028669268717193436606467893596092138461980127494443633833240727656111845218913366368409973705945852937546861824400761190208020898528637753121592516776112625907992364462314541787272525903553709*j+577381321884404637808391205760891583039563709656261523188365644721839518395090059951637547964468822095873175263496035533602755866120293764901657313023000531067504187145886064396313750664978511489728895045390110744492542111986) over Finite Field in j of size 10354717741769305252977768237866805321427389645549071170116189679054678940682478846502882896561066713624553211618840202385203911976522554393044160468771151816976706840078913334358399730952774926980235086850991501872665651576831^2


for the further walk process. (That i maybe do not understand in this second.) The walking part may be given by the following:

MAXLEVEL = 372

def walk(E, jString, level, previous_jString ):
children_DIC = {}    # dictionary jString-value -> number of occurences
if level >= MAXLEVEL:
return

for P in E(0).division_points(2):
if P == E(0):
continue
# psi = E.isogeny(P)
# EE  = psi.codomain()
EE = pointIsogeny( E, [P], [], [2] )[0]
jj = str( EE.j_invariant() )    # ?a string (here and elsewhere in similar places)
level_child = level + 1
if jj in children_DIC:
children_DIC[jj] += 1
else:
children_DIC[jj]  = 1

print "Level = %s children_DIC[%s] = %s\n" % ( level, jj, children_DIC[jj] )
if (jj != previous_jString) or (children_DIC[jj] > 1):
walk(EE, jj, level+1, jString)

print "start: ", str(E.j_invariant())
walk(E, str(E.j_invariant()), 0, "")


It may be that the children dictionary is / should be a global constant. (Then it should be defined outside, and declared as global inside the walk.) I inserted the print in order to see something.

At any rate, the 2-torsion isogeny problem should be solved so far.

more

@dan_fulea Thank you for the answer. But apparently, this code does not finish for the above primes. I believe the problem is not the isogeny computation, but how I walk in the graph recursively. Do you have any idea how to handle this more efficiently?

( 2018-03-02 16:15:21 -0500 )edit

In your starting_node function you are constructing points of order 3 but it seems that the function E_i(0).division_points(3) returns points whose order attribute is not set. That should be considered a bug, and probably my fault, certainly fixable. You can get around it rather easily if instead of constructing the isogenies as you do you call the method E_i.isogenies_prime_degree(3) which returns (quickly) a list of all 3-isogenies from E_i. This will include 3-isogenies whose kernel points are not in GF(p) (only their x-coordinates will certainly be), so if that matters you'll need to add a little more.

In the code as is Sage is computing the order of the point P before constructing the 3-isogeny, with no clue as to what that order might be or any bound on it, which may well cause the computation of the cardinality of the curve and its factorization.

more

I tried the above code, and I think the main bottleneck is inside the walk function, in the for loop when it tries to compute E_i(0).division_points(2), it doesn't seem to finish. From what I understood, you ask to change the starting_node function, with something like this:

isogs = E_i.isogenies_prime_degree(3)
psi = isogs[randint(0, len(isogs)-1)]
return starting_node(psi.codomain(), count+1)


Am I right? Although, this will work for the starting_node function, I believe the problem still remains inside the walk function. Maybe one can do the same thing there too, and remove the for loop and use this:

isogs = E_i.isogenies_prime_degree(2)
psi = isogs[randint(0, len(isogs)-1)]


But I don't know whether that will give the same results that OP asks for.

( 2018-02-26 11:59:16 -0500 )edit