Îlu algorithm: asymptotically faster elliptic-curve isogenies

[yyyyx4]

Sage currently computes an isogeny of prime degree ℓ using ϴ(ℓ) base-field operations. We can achieve the same task in Õ(√ℓ) operations using the √élu algorithm; see https://velusqrt.isogeny.org.

This patch adds a fairly reasonable implementation of the algorithm to Sage. It works for odd-degree isogenies over all finite fields and can outperform the current EllipticCurveIsogeny class starting from degrees around 100, depending on the concrete performance of the underlying base-field and polynomial arithmetic.

The new code is marked as experimental and isn't used anywhere by default. After some real-world testing and careful benchmarking, we can (and should) let Sage automatically decide which implementation to use for a given input.

CC: @defeo @JohnCremona @categorie @sagetrac-bensmith @sagetrac-tanja @tscrim @videlec @slel @kwankyu

Component: elliptic curves

Author: Lorenz Panny

Branch: a0f4c4a

Reviewer: John Cremona, Travis Scrimshaw, Kwankyu Lee

Issue created by migration from https://trac.sagemath.org/ticket/34303

[yyyyx4]

Author: Lorenz Panny

[yyyyx4]

Commit: 4b24a77

[yyyyx4]

Branch: public/velusqrt_algorithm_for_elliptic_curve_isogenies

[JohnCremona]

comment:2

I have a question: although I know that the application of greatest interest involves curves over finite fields, I assume that the faster algorithm also works in characteristic 0, so that isogenies between curves defined over number fields will also be affected by this? I hope so.

[yyyyx4]

comment:3

That's a good question. The Îlu speedup works best if you have an actual rational kernel point as input, which requires big field extensions for large isogeny degrees over number fields. There is a variant of the algorithm for rational kernels made up of irrational points, but the speedup is smaller in that case.

Do you have an application in mind? Which field and isogeny degree, and how is the input given? It shouldn't be hard to make this work for number fields, but I'm not sure if it helps at all, which is why I restricted the current implementation to finite fields.

[yyyyx4]

comment:4

I was hoping we could get this into the upcoming release so that other people can play around with it. This is opt-in and marked as experimental, so a minimum viable review basically consists of checking that it doesn't do any harm. Could someone please spare a few minutes to have a look?

[JohnCremona]

comment:5

Replying to @yyyyx4:

That's a good question. The Îlu speedup works best if you have an actual rational kernel point as input, which requires big field extensions for large isogeny degrees over number fields. There is a variant of the algorithm for rational kernels made up of irrational points, but the speedup is smaller in that case.

Do you have an application in mind? Which field and isogeny degree, and how is the input given? It shouldn't be hard to make this work for number fields, but I'm not sure if it helps at all, which is why I restricted the current implementation to finite fields.

The isogenies I compute are mostly over number fields and of small prime degree ell, where "small" mostly means <50. And for most of these ell, I have implemented special code for them (which works over finte fields too, except characteristic ell). As for the input, it's just the curve and ell and I need all the ell-isogenous curves; there are easy necessary conditions for an ell-isogeny to exist so I would rarely try to compute one when it does not exist.

[JohnCremona]

comment:6

Replying to @yyyyx4:

I was hoping we could get this into the upcoming release so that other people can play around with it. This is opt-in and marked as experimental, so a minimum viable review basically consists of checking that it doesn't do any harm. Could someone please spare a few minutes to have a look?

I had a first look, and will come back to it tomorrow.

[JohnCremona]

Reviewer: John Cremona

[JohnCremona]

Work Issues: document parameters for all methods and functions

[JohnCremona]

comment:7

Here are some comments made as I read the code, so not in order of importance. Point 4 about documenting parameters is the most important -- otherwise docstrings, examples ans tests are very good.

1. Instead of

    sage: sum(iso(psiQ) == phiQ for iso in isos)
    1

I would suggest

    sage: any(iso(psiQ) == phiQ for iso in isos)
    True

2. There's a comment which says that the code works on "short Weierstrass curves", so I assume that menas that the characteristic should not be 2 (and perhaps also not 3). This should probably be stated at the top of the file hom_velusqrt.py. I haven't yet checked whether this condition is checked in the code.

3. I see the comment about this implementation only being over finite fields. I think that is fine, given that isogenies of high degree over number fields are rarely encountered (since the field would have to have high degree over QQ), as I said in an earlier comment.

4. Methods' docstrings should have an INPUT block documenting all parameters (apart from self), for example the _init_ method, of ProductTree and more (not itemised here!).

I see that the code here is not called from elsewhere -- in particular, not called by default by isogeny methods of elliptic curve classes -- so is certainly harmless to include.

[sagetrac-git]

Changed commit from 4b24a77 to 986462f

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

b0b8cfb	clarify that curves without a short Weierstraß model only exist in characteristics 2 and 3
9f9a77c	improve documentation for Îlu code
8bdcd45	move comparison by evaluation up from hom_composite.py to hom.py
1e80198	implement more flexible comparison for EllipticCurveHom children including Îlu
ab4b4bf	ensure correct base rings (example from failing random doctest)
986462f	add glue code in EllipticCurve_field.isogeny() to use the Îlu implementation

[yyyyx4]

comment:9

Thank you for the comments! I added INPUT: blocks to (I believe) all functions that need them.

I also made comparison work using the same code as for composite isogenies (by reducing to evaluation), and added a (purely opt-in) call path from the .isogeny() method for elliptic curves so that it's easy to use the new algorithm without typing long import statements.

About your other remarks:

1. sum(...) == 1 is a stronger statement than any(...), which is why I wrote it this way: The former checks unique existence, the latter just existence.

2. The implementation does work for some curves in characteristic 3, but not all of them. In the code, it simply tries converting the curve to short Weierstraß form and throws an exception if that fails. This is slightly more general than checking the characteristic.

3. Indeed, the situation over number fields mentioned in #34303 comment:5 does not seem like it would benefit from this. Should it turn out that there's something we can do after all, this can of course be taken care of at a later time.

[yyyyx4]

Changed work issues from document parameters for all methods and functions to none

[tscrim]

comment:10

Replying to @JohnCremona:

4. Methods' docstrings should have an INPUT block documenting all parameters (apart from self), for example the _init_ method, of ProductTree and more (not itemised here!).

One comment on this: To make it more publicly visible, it can be useful to have the inputs for __init__ to be at the class-level documentation and not (repeated) at the method-level since it does not appear in the docstrings. However, this is mostly a suggestion when you expect the user to directly call the class to build the object (rather than, e.g., an element of a parent which is through the parent's _element_constructor_).

One additional comment to John's, the INPUT: items should always be in code format, e.g.:

     INPUT:
 
-    - `n` -- odd integer `\geq 5`
+    - ``n`` -- odd integer `\geq 5`

In _choose_IJK, it would also be good to specify that this needs to be an :class:`Integer` otherwise the method will break. In fact, since you are just wanting the integer part, you could do

(n-1).sqrtrem()[0] // 2

As a general programming note, I prefer to avoid similar variable names, such as b and b_. It is too easy to make a mistake with them.

[sagetrac-git]

Changed commit from 986462f to 368a18e

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

156c820	prevent crash in _random_example_for_testing() when group is trivial
3d63a50	move INPUT blocks from .__init__() docstrings to class docstrings
c874dc3	explicitly convert given n to Integer in FastEllipticPolynomial
e499ea6	clarify that _choose_IJK needs an Integer; use .isqrt(); rename b_ to c
9c7b344	change all INPUT parameters to code formatting
368a18e	ensure correct base rings (for real this time); cf. ab4b4bf817b

[yyyyx4]

comment:12

Thanks! All done.

[JohnCremona]

comment:13

Thanks for making all these changes, and answering my comments. I am happy with this, so if Travis is too he can add his name to the reviewers list and set it to "positive review".

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

5fe7ef4	add "short_weierstrass" to list of valid inputs for compute_model()
a0f4c4a	remove experimental warning

[yyyyx4]

comment:20

Thanks, done! I've also created #34409, which (inspired by your comments here) suggests to remove the experimental warning for composite elliptic-curve isogenies.

[kwankyu]

comment:21

Replying to @tscrim:

You have 2 (2)'s leading to a non-bijective map. ;)

Fixed. Thanks. It was awful.

[kwankyu]

Changed reviewer from John Cremona, Travis Scrimshaw to John Cremona, Travis Scrimshaw, Kwankyu Lee

[kwankyu]

comment:22

Replying to @yyyyx4:

Thanks, done! I've also created #34409, which (inspired by your comments here) suggests to remove the experimental warning for composite elliptic-curve isogenies.

Thank you. It was perhaps your being modest. But Sage should always deliver reliable code, and I think we should avoid including experimental code in Sage.

[yyyyx4]

comment:23

Thanks a lot everyone!

(Indeed, I'm reasonably convinced that the code is correct. I added the experimental warning not because I feared things were broken, but because it'd make it easier to tweak the API later. Presumably though, this won't be necessary anyway.)

[vbraun]

Changed branch from public/velusqrt_algorithm_for_elliptic_curve_isogenies to a0f4c4a

[antonio-rojas]

Changed commit from a0f4c4a to none

[antonio-rojas]

comment:25

Got a test failure with a specific random seed

sage -t --long --random-seed=66297844331759291946775140613791774160 /usr/lib/python3.10/site-packages/sage/schemes/elliptic_curves/hom_velusqrt.py
**********************************************************************
File "/usr/lib/python3.10/site-packages/sage/schemes/elliptic_curves/hom_velusqrt.py", line 829, in sage.schemes.elliptic_curves.hom_velusqrt.EllipticCurveHom_velusqrt
Failed example:
    any(map(check, psi.codomain().isomorphisms(phi.codomain())))
Expected:
    True
Got:
    False
**********************************************************************

[fchapoton]

comment:26

and the linter now requires a blank line before nested definitions

[yyyyx4]

comment:27

Thanks both of you. See #34467 and #34466.