Opened 2 years ago

Closed 7 months ago

Last modified 7 months ago

#23222 closed defect (fixed)

Add test to EllipticCurve.isogeny() to check for validity of kernel polynomials

Reported by: defeo Owned by:
Priority: major Milestone: sage-8.7
Component: elliptic curves Keywords: isogenies, vélu's formulas
Cc: cremona Merged in:
Authors: John Cremona Reviewers: Luca De Feo, Frédéric Chapoton
Report Upstream: N/A Work issues:
Branch: 28e369f (Commits) Commit:
Dependencies: Stopgaps:

Description (last modified by defeo)

In this example

sage: K = GF(47^2)
sage: E = EllipticCurve([0, K.gen()])
sage: f = E.division_polynomial(7).factor()[4][0]
sage: E.is_supersingular()
True
sage: E.isogeny(f).codomain().is_supersingular()
False
sage: phi = E.isogeny(E.division_polynomial(7).factor()[8][0])
sage: phi(E.random_point())
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
...
TypeError: Coordinates [16*z2 + 30, 3*z2 + 22, 1] do not define a point on Elliptic Curve defined by y^2 = x^3 + (3*z2+36)*x over Finite Field in z2 of size 47^2

The factors of the seventh division polynomial are not detected as bad kernel polynomials, and a wrong isogeny is computed instead.

Tests to detect these cases are in .isogenies_prime_degree_general(). They should be moved to an independent function and called before applying Vélu's formulas.

Change History (31)

comment:1 follow-up: Changed 2 years ago by cremona

Remember that the Trace of Frobenius is an integer, not an element of the field, so seeing 94 which is 0 mod 47 is ok. Secondly, you have omitted the isogeny whose kernel is the product of the 3 linear factors of the 7-division polynomial. Thirdly, if you want to see the 7-isogenous curves it is better to do

E.isogenies_prime_degree(7)

since it is non-trivial to combine the factors of the division polynomial correctly. In this case there are only two 7-isogenies. Not that E.isogeny(f) does *not* actually check that f is a valid kernel polynomial even when check=True is give (the default). It should! When it was written we did not have code to do that but now we do. It is slightly buried in the function isogenies_prime_degree_general().

comment:2 in reply to: ↑ 1 ; follow-up: Changed 2 years ago by defeo

Replying to cremona:

Remember that the Trace of Frobenius is an integer, not an element of the field, so seeing 94 which is 0 mod 47 is ok.

I agree with that, but doesn't Sato-Tate say that two curves are isogenous over a finite field iff they have the same trace? As I said, j=0 and j=36 are isogenous over GF(472), but the exact curves in the example are not. And, indeed:

sage: phi = E.isogeny(E.division_polynomial(7).factor()[8][0])
sage: phi(E.random_point())
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
...
TypeError: Coordinates [16*z2 + 30, 3*z2 + 22, 1] do not define a point on Elliptic Curve defined by y^2 = x^3 + (3*z2+36)*x over Finite Field in z2 of size 47^2

Secondly, you have omitted the isogeny whose kernel is the product of the 3 linear factors of the 7-division polynomial.

Yes. I just wanted to highlight the bug. One of the buggy factors would have been enough.

Thirdly, if you want to see the 7-isogenous curves it is better to do

E.isogenies_prime_degree(7)

since it is non-trivial to combine the factors of the division polynomial correctly.

Absolutely. However, I found the bug while doing a more general search.

Not that E.isogeny(f) does *not* actually check that f is a valid kernel polynomial even when check=True is give (the default). It should! When it was written we did not have code to do that but now we do. It is slightly buried in the function isogenies_prime_degree_general().

Well, that's not completely correct:

sage: E.isogeny(E.base_field().polynomial_ring().random_element(3))
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
...
ValueError: The polynomial does not define a finite subgroup of the elliptic curve.

So there are some checks, but they are buggy, at least for j=0,1728, I suspect.

I do not have time to dig this further right now. Hopefully in a couple of weeks I can come up with a fix.

comment:3 Changed 2 years ago by cremona

The j-invariant does not determine the trace since it only determines the curve up to twist; in your case there are more than 2 twists as well.

Yes, some checking is done for validity of kernel polys, but that is misleading since bad input passes those tests as they are not complete.

We'll come back to that.

comment:4 Changed 17 months ago by cremona

Can we close this as invalid / won't fix? I just re-read the comments and think that I explained why none of the output which had been thought to be suspect was incorrect.

comment:5 in reply to: ↑ 2 Changed 17 months ago by defeo

I still don't get it. Admittedly, E.isogeny(f) does not test that f is a valid kernel polynomial. But, why is f not a valid kernel polynomial in this example:

sage: f = E.division_polynomial(7).factor()[4][0]
sage: E.is_supersingular()
True
sage: E.isogeny(f).codomain().is_supersingular()
False

Also, the TypeError in comment 2 is not nice and would deserve fixing.

comment:6 Changed 17 months ago by cremona

The kernel poly of a 7-isogeny is a degree 3 factor of the 7-division poly but not every degree 3 factor is a kernel poly! In this case the 3 roots of your degree 3 factor, which all lie in GF(476), are the x-coordinates of points in 3 different subgroups of order 7.

In the function isogenies_prime_degree_general() in sage.schemes.elliptic_curves.isogeny_small_degree, this is implemented properly: factor the l-division poly (here l is an odd prime); keep only factors of degree dividing (l-1)/2; for each such factor work out whether or not it is an irreducible factor of a kernel poly.

It would be possible to extract from the code of that function a test for whether a polynomial of degree (l-1)/2 actually is a kernel poly. Feel free to make a ticket for that! The relevent theory is in Kimi Tsukazaki's thesis (cited in that file).

comment:7 Changed 17 months ago by defeo

Oh, right!

Can't we just edit this ticket so that it points to the problem, instead of closing it and opening a new one?

comment:8 Changed 17 months ago by cremona

Sure, go ahead. Once we have a stand-alone function test_kernel_poly(E,f,ell) it is clear where to put a call to it: in place of lines 2450-2454 where now it only tests if the polynomial divides the division poly. I should have done this ages ago.

comment:9 Changed 17 months ago by defeo

  • Description modified (diff)
  • Keywords vélu's formulas added; supersingular removed
  • Summary changed from Incorrect output from Vélu's formulas on j=0? to Add test to EllipticCurve.isogeny() to check for validity of kernel polynomials
  • Type changed from defect to enhancement

comment:10 Changed 16 months ago by cremona

  • Type changed from enhancement to defect

Working on it...

comment:11 Changed 16 months ago by cremona

  • Authors set to John Cremona
  • Branch set to u/cremona/23222
  • Commit set to 0c120af3a0cb8323954bd5994984bc433adadfe8
  • Milestone changed from sage-8.0 to sage-8.2
  • Status changed from new to needs_review

New commits:

0c120af#23222: check validity of isogeny kernel polynomials

comment:12 Changed 16 months ago by git

  • Commit changed from 0c120af3a0cb8323954bd5994984bc433adadfe8 to f7518f7da457d7dc25a32a0ba9e92ab58a189b45

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

f7518f7#23222: fix docstring formatting error

comment:13 follow-up: Changed 9 months ago by chapoton

  • Status changed from needs_review to needs_work

red banch, needs rebasing

comment:14 in reply to: ↑ 13 Changed 9 months ago by cremona

Replying to chapoton:

red branch, needs rebasing

I will rebase it, but it would be good to have a commitment from someone to review it, since the actual work on this was done over 6 months ago.

comment:15 Changed 9 months ago by git

  • Commit changed from f7518f7da457d7dc25a32a0ba9e92ab58a189b45 to 0683290a5a04626d57eed1e28910bfc9346a6f34

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

0683290Merge branch 'develop' into 23222

comment:16 Changed 9 months ago by cremona

  • Status changed from needs_work to needs_review

comment:17 Changed 9 months ago by defeo

It's been in my todo list for 6 months indeed. But I won't have time before Christmas.

comment:18 Changed 7 months ago by defeo

I have one doctest failure:

sage -t --long --warn-long 83.4 src/sage/schemes/elliptic_curves/ell_field.py                                                          
**********************************************************************                                                                 
File "src/sage/schemes/elliptic_curves/ell_field.py", line 1023, in sage.schemes.elliptic_curves.ell_field.EllipticCurve_field.isogeni$
s_prime_degree
Failed example:
    E.isogenies_prime_degree(5)
Exception raised:
    Traceback (most recent call last):
      File "/data/dfl/sage/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 671, in _run
        self.compile_and_execute(example, compiler, test.globs)
      File "/data/dfl/sage/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 1086, in compile_and_execute
        exec(compiled, globs)
      File "<doctest sage.schemes.elliptic_curves.ell_field.EllipticCurve_field.isogenies_prime_degree[28]>", line 1, in <module>
        E.isogenies_prime_degree(Integer(5))
      File "/data/dfl/sage/local/lib/python2.7/site-packages/sage/schemes/elliptic_curves/ell_field.py", line 1091, in isogenies_prime$
degree
        return sum([isogenies_prime_degree(self, d) for d in L], [])
      File "/data/dfl/sage/local/lib/python2.7/site-packages/sage/schemes/elliptic_curves/isogeny_small_degree.py", line 2328, in isog$
nies_prime_degree
        return isogenies_prime_degree_general(E,l, minimal_models=minimal_models)
      File "/data/dfl/sage/local/lib/python2.7/site-packages/sage/schemes/elliptic_curves/isogeny_small_degree.py", line 2197, in isog$
nies_prime_degree_general
        return [E.isogeny(k) for k in ker]
      File "/data/dfl/sage/local/lib/python2.7/site-packages/sage/schemes/elliptic_curves/ell_field.py", line 880, in isogeny
        raise RuntimeError("Unable to construct isogeny: %s" % e)
    RuntimeError: Unable to construct isogeny: 'sage.rings.fraction_field_FpT.FpTElement' object has no attribute 'divides'
**********************************************************************

...investigating.

comment:19 Changed 7 months ago by defeo

The doctest failure is due to #27064. No issue to be addressed here.

comment:20 Changed 7 months ago by defeo

In is_kernel_polynomial, there is this test:

psi_m = E.division_polynomial(m)
if not f.divides(psi_m):
    return False

isn't it (much) more efficient to do

if E.division_polynomial(m, x=f.parent().quo(f).gen()) != 0:
    return False
Last edited 7 months ago by defeo (previous) (diff)

comment:21 Changed 7 months ago by defeo

The doc still says

   * "check" (default: True) does some partial checks that the

        input is valid (e.g., that the kernel polynomial provided is
        valid); however, invalid input can in some cases still pass,
        since that the points define a group is not checked.

what are the invalid inputs that still pass?

comment:22 Changed 7 months ago by cremona

  • Milestone changed from sage-8.2 to sage-8.7
  • Status changed from needs_review to needs_work

I will correct the docstring and also avoid using f.divides as suggested.

comment:23 Changed 7 months ago by git

  • Commit changed from 0683290a5a04626d57eed1e28910bfc9346a6f34 to 72a61bbd024b6b089b71aa8a3405f07e6fb55ae2

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

f123620Merge branch 'develop' of trac.sagemath.org:sage into 23222
0c30e53Merge branch 'develop' into 23222
7e7593f#23222 small change to avoid #27064, and update docstring
72a61bbMerge remote-tracking branch 'trac/u/cremona/23222' into 23222

comment:24 Changed 7 months ago by cremona

  • Reviewers set to defeo, chapoton
  • Status changed from needs_work to needs_review

comment:25 Changed 7 months ago by defeo

  • Branch changed from u/cremona/23222 to u/defeo/23222

comment:26 Changed 7 months ago by git

  • Commit changed from 72a61bbd024b6b089b71aa8a3405f07e6fb55ae2 to 28e369f252fa8339ba379c5619096016afdc4493

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

f123620Merge branch 'develop' of trac.sagemath.org:sage into 23222
0c30e53Merge branch 'develop' into 23222
7e7593f#23222 small change to avoid #27064, and update docstring
72a61bbMerge remote-tracking branch 'trac/u/cremona/23222' into 23222
28e369fMerge branch 'u/cremona/23222' of git://trac.sagemath.org/sage into t/23222/23222

comment:27 Changed 7 months ago by defeo

  • Reviewers changed from defeo, chapoton to Luca De Feo, Fréderic Chapoton
  • Status changed from needs_review to positive_review

Doctests pass. I just fixed a minor problem in the docstring (r""" instead of """). Thanks, and sorry for the long time.

comment:28 Changed 7 months ago by cremona

Thanks!

comment:29 Changed 7 months ago by chapoton

  • Reviewers changed from Luca De Feo, Fréderic Chapoton to Luca De Feo, Frédéric Chapoton

comment:30 Changed 7 months ago by vbraun

  • Branch changed from u/defeo/23222 to 28e369f252fa8339ba379c5619096016afdc4493
  • Resolution set to fixed
  • Status changed from positive_review to closed

comment:31 Changed 7 months ago by chapoton

  • Commit 28e369f252fa8339ba379c5619096016afdc4493 deleted

one py3-incompatible import was added here ; fix at #27107

Note: See TracTickets for help on using tickets.