Opened 5 years ago

Closed 4 years ago

#18916 closed enhancement (fixed)

Use Kedlaya algorithm to count points on hyperelliptic curves

Reported by: jpflori Owned by:
Priority: major Milestone: sage-7.3
Component: number fields Keywords: hyperelliptic curves, matrix of Frobenius
Cc: defeo, katestange Merged in:
Authors: Kiran Kedlaya Reviewers: Jean-Pierre Flori, Frédéric Chapoton
Report Upstream: N/A Work issues:
Branch: e79317e (Commits) Commit: e79317e72f2e50cb94e00e2a91fd643da644875e
Dependencies: #20219, #20581 Stopgaps:

Description (last modified by jpflori)

There is an implementation in monsky_washnitzer.py and one in PARI/GP. Let's use the PARI/GP one.

Change History (47)

comment:1 Changed 5 years ago by jpflori

The magic line reads something like:

monsky_washnitzer.matrix_of_frobenius_hyperelliptic(self.hyperelliptic_polynomials()[0], p, prec)[0]

The defining polynomial needs to be monicified which is also required by hypellfrob, so it makes sense to refactor the code a little bit. And it only works for prime fields.

There is also a PARI/GP implem without these restriction (only p != 2 I'd say)

comment:2 Changed 5 years ago by jdemeyer

+1 to use PARI/GP

comment:3 follow-up: Changed 4 years ago by kedlaya

  • Keywords hyperelliptic curves matrix of Frobenius added

From a quick inspection of the PARI/GP implementation:

The command hyperellpadicfrobenius appears to only work for prime fields. There is an internal command nfhyperellpadicfrobenius which seems to handle nonprime fields, which is called by hyperellcharpoly to get the characteristic polynomial of Frobenius; but it does not appear to be exposed.

Also, the Frobenius matrix commands have a similar restriction on p as in Harvey's code. The charpoly command seems to be doing a naive point count for smaller p.

In any case, it would be worth updating HyperellipticCurve?.frobenius_polynomial and related methods to include calls to PARI as appropriate.

comment:4 Changed 4 years ago by kedlaya

Addendum: naive point counting is already implemented directly in Sage, see #11980.

comment:5 in reply to: ↑ 3 ; follow-up: Changed 4 years ago by jpflori

To make things clear, I'm looking at PARI git version, though the version in Sage is quite recent and should be similar as far as hyperell.c is concerned.

Replying to kedlaya:

From a quick inspection of the PARI/GP implementation:

The command hyperellpadicfrobenius appears to only work for prime fields. There is an internal command nfhyperellpadicfrobenius which seems to handle nonprime fields, which is called by hyperellcharpoly to get the characteristic polynomial of Frobenius; but it does not appear to be exposed.

I confirm that hyperellpadicfrobenius only works for prime fields whereas the nf... one handles extensions. Though both of them are not marked static and use gerepile magic to leave a clean stack so I guess they are possible and safe to use from outside the PARI library.

Also, the Frobenius matrix commands have a similar restriction on p as in Harvey's code. The charpoly command seems to be doing a naive point count for smaller p.

Do you mean those in the PARI library? I don't see at first glance such restriction. There is some naive point counting involved in hyperellcharpoly but that is when the curve is of (very) low genus and the characteristic is small. Did I miss something else?

For p == 2 it indeed seems it only works for F_2 using naive point counting.

In any case, it would be worth updating HyperellipticCurve?.frobenius_polynomial and related methods to include calls to PARI as appropriate.

comment:6 in reply to: ↑ 5 ; follow-up: Changed 4 years ago by kedlaya

Replying to jpflori:

I confirm that hyperellpadicfrobenius only works for prime fields whereas the nf... one handles extensions. Though both of them are not marked static and use gerepile magic to leave a clean stack so I guess they are possible and safe to use from outside the PARI library.

Bill Allombert confirmed this. The nf... command is not (currently) available in GP, but we can call it from the library.

Also, the Frobenius matrix commands have a similar restriction on p as in Harvey's code. The charpoly command seems to be doing a naive point count for smaller p.

Do you mean those in the PARI library? I don't see at first glance such restriction. There is some naive point counting involved in hyperellcharpoly but that is when the curve is of (very) low genus and the characteristic is small. Did I miss something else?

Oh, maybe that is just an optimization for cases where naive point counting is actually faster?

comment:7 in reply to: ↑ 6 Changed 4 years ago by jpflori

Replying to kedlaya:

Did I miss something else?

Oh, maybe that is just an optimization for cases where naive point counting is actually faster?

I guess so, it just looks like hand optimized random values.

comment:8 follow-up: Changed 4 years ago by kedlaya

I think I need some assistance with the PARI library interface here.

I mentioned above that nfhyperellpadicfrobenius is not currently exposed. However, this function does appear in src/sage/libs/pari/paridecl.pxd. Does something else need to be done in order to make it visible as a method of a PARI object?

comment:9 in reply to: ↑ 8 ; follow-up: Changed 4 years ago by jpflori

Replying to kedlaya:

I think I need some assistance with the PARI library interface here.

I mentioned above that nfhyperellpadicfrobenius is not currently exposed. However, this function does appear in src/sage/libs/pari/paridecl.pxd. Does something else need to be done in order to make it visible as a method of a PARI object?

The methods of PARI objecgt are in gen.pxd/pyx. You should add a method there calling the C function. Then in the hyperelliptic curve class, something like

pari(H.polynomial()).frobenius()

could be called.

Maybe a better option is to add a method directly in Sage's hyperelliptic curve over finite field class.

I'll have some free time in a few hours, I can have a look as well.

comment:10 Changed 4 years ago by jdemeyer

If you do make any changes to the PARI interface, please base them on #20219.

comment:11 Changed 4 years ago by kedlaya

  • Dependencies set to #20219

comment:12 Changed 4 years ago by kedlaya

  • Branch set to u/kedlaya/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves

comment:13 Changed 4 years ago by git

  • Commit set to 6a4e932eef76137b6856c3b84f2fe82de807cc3a

Branch pushed to git repo; I updated commit sha1. Last 10 new commits:

c0ed97aStop using deprecated PARI factoring features
edc5ce2Merge branch 't/20205/get_rid_of_factorint_withproof_sage_in_pari_interface' into HEAD
5fb408dReplace pari_catch_sig_on by sig_on
d5c934cDeprecate PARI nth_prime and prime_list
d7d2d7dRemove redundant functions from pari_instance.pyx
130fc90Remove redundant functions from gen.pyx
2299743Use SR.symbol() instead var()
6acc79aAdd back order() as deprecated alias of multiplicative_order()
0f2386aFix basic t_LIST functionality
6a4e932Merge branch 'u/jdemeyer/remove_redundant_functions_from_gen_pyx' of git://trac.sagemath.org/sage into t/18916/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves

comment:14 in reply to: ↑ 9 ; follow-up: Changed 4 years ago by kedlaya

Replying to jpflori:

I mentioned above that nfhyperellpadicfrobenius is not currently exposed. However, this function does appear in src/sage/libs/pari/paridecl.pxd. Does something else need to be done in order to make it visible as a method of a PARI object?

The methods of PARI objecgt are in gen.pxd/pyx. You should add a method there calling the C function.

Maybe I wasn't clear here. I think there is some mechanism by which most PARI functions appear "automatically" as methods of a PARI object; most of them are not explicitly declared in src/sage/libs/pari/gen.pyx. I was wondering why hyperellpadicfrobenius (which does not have its own explicitly defined wrapper as far as I can find) shows up this way while nfhyperellpadicfrobenius does not. I had thought paridecl.pxd was the master source for the autogenerated methods, but it includes both of these functions and yet only one makes it through the rest of the process.

comment:15 in reply to: ↑ 14 Changed 4 years ago by jpflori

Replying to kedlaya:

Replying to jpflori:

I mentioned above that nfhyperellpadicfrobenius is not currently exposed. However, this function does appear in src/sage/libs/pari/paridecl.pxd. Does something else need to be done in order to make it visible as a method of a PARI object?

The methods of PARI objecgt are in gen.pxd/pyx. You should add a method there calling the C function.

Maybe I wasn't clear here. I think there is some mechanism by which most PARI functions appear "automatically" as methods of a PARI object; most of them are not explicitly declared in src/sage/libs/pari/gen.pyx. I was wondering why hyperellpadicfrobenius (which does not have its own explicitly defined wrapper as far as I can find) shows up this way while nfhyperellpadicfrobenius does not. I had thought paridecl.pxd was the master source for the autogenerated methods, but it includes both of these functions and yet only one makes it through the rest of the process.

Oh I see. I'd say it's exactly because the nf... function is not in GP. IIRC functions in GP are described in pari.desc which is also used for autogenerating our wrappers. This puts Cython methods in the gen_auto class in some pxi auto generated files.

The paridecl.pxd file is still handwritten/copied from the PARI library headers and is not involved during the generation of the wrappers (though it is needed to help Cython at compilation time I'd say). Part of it could surely be auto-generated from pari.desc (or even all of it from the C headers though that would be more involved).

Therefore in the end we have Cython access to the nf... but no method wrapping it.

See #17631 and #17860 and https://github.com/sagemath/sage/blob/master/src/sage_setup/autogen/pari/generator.py

Jeroen, please correct me if I'm wrong.

comment:16 Changed 4 years ago by kedlaya

Bill Allombert is wondering what use case(s) we have in mind for the Frobenius matrix (rather than just its characteristic polynomial). Keep in mind that they don't (currently) provide Coleman integration, which rules out the Coleman-Gross height pairing. But I think the p-adic sigma function on an elliptic curve does constitute an example where the Frobenius matrix is strictly more useful than its charpoly.

comment:17 Changed 4 years ago by kedlaya

Sorry, I seem to have dragged this ticket somewhat off topic. I created the new ticket #20309 to discuss Frobenius matrices.

The original topic of this ticket was point counting (and by extension zeta functions), and for that by far the easiest solution would be to use PARI's hyperellcharpoly. Better would be to intelligently combine this with Harvey's code for Frobenius matrices when it applies (it's very fast but not as general), but maybe that can be a separate ticket later.

comment:18 follow-up: Changed 4 years ago by kedlaya

  • Report Upstream changed from N/A to Reported upstream. No feedback yet.

I'm having trouble getting the PARI code to work in characteristic 2, but even in odd characteristic there seem to be some issues:

pari: T = ffinit(3, 2);
pari: y = ffgen(T, 'y);
pari: hyperellcharpoly(x^3 + y*x)
***   at top-level: hyperellcharpoly(x^3
  ***                 ^--------------------
  *** hyperellcharpoly: bug in PARI/GP (Segmentation Fault), please report.

Off to go report this upstream...

comment:19 Changed 4 years ago by jpflori

IIRC PARI supposes that the equation is y² = f(x) and the function eats f(x), in particular in char 2 this is a singular curve and so says PARI.

comment:20 in reply to: ↑ 18 Changed 4 years ago by kedlaya

  • Report Upstream changed from Reported upstream. No feedback yet. to Fixed upstream, but not in a stable release.

Replying to kedlaya:

I'm having trouble getting the PARI code to work in characteristic 2, but even in odd characteristic there seem to be some issues:

pari: T = ffinit(3, 2);
pari: y = ffgen(T, 'y);
pari: hyperellcharpoly(x^3 + y*x)
***   at top-level: hyperellcharpoly(x^3
  ***                 ^--------------------
  *** hyperellcharpoly: bug in PARI/GP (Segmentation Fault), please report.

Off to go report this upstream...

Bill Allombert reports that this particular issue (the one not in characteristic 2) is fixed upstream:

commit 31253df93771c72e44c19c245dcf514822a3bdb4
Author: Bill Allombert <Bill.Allombert@math.u-bordeaux1.fr>
Date:   Fri Mar 11 12:38:30 2016 +0100

    hyperellcharpoly: fix typo in ZXX_to_FpXC

comment:21 Changed 4 years ago by git

  • Commit changed from 6a4e932eef76137b6856c3b84f2fe82de807cc3a to e510540635d871837fd5de0ae47f3fdbf51f07e6

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

e510540Implement wrapper for PARI hyperellcharpoly

comment:22 Changed 4 years ago by git

  • Commit changed from e510540635d871837fd5de0ae47f3fdbf51f07e6 to 95ede4bc4354e8717bfb38143d8b5ebb168bf878

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

95ede4bCorrectly handle frobenius_matrix() for an even degree polynomial

comment:23 Changed 4 years ago by kedlaya

The previous commits call to PARI in all cases where the characteristic is odd and David Harvey's code is not applicable (either the field is not prime, or the characteristic is too small, or h != 0, or f is of even degree). The naive code is still used in characteristic 2.

There is one failing doctest due to the PARI issue described above. Presumably this will be fixed once we pull from upstream.

comment:24 Changed 4 years ago by git

  • Commit changed from 95ede4bc4354e8717bfb38143d8b5ebb168bf878 to 4c92503eeaed4d9cca5016ee5fec2c6082eb65d8

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

4c92503Fix count_points for even-degree hyperelliptics

comment:25 Changed 4 years ago by kedlaya

  • Authors set to Kiran Kedlaya

The last commit addresses #20391.

comment:26 Changed 4 years ago by git

  • Commit changed from 4c92503eeaed4d9cca5016ee5fec2c6082eb65d8 to 6ca4b8c808ba9a53a17ce061a96ecbabea58d595

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

6ca4b8cAdd doctest regarding #20391

comment:27 Changed 4 years ago by kedlaya

... and a doctest to indicate that #20391 is addressed.

comment:28 follow-up: Changed 4 years ago by jpflori

@Kiran: is it in state of review or do you plan further changes?

comment:29 in reply to: ↑ 28 Changed 4 years ago by kedlaya

Replying to jpflori:

@Kiran: is it in state of review or do you plan further changes?

Even after rebasing to 7.2.rc1, there is still one failing doctest:

File "src/sage/schemes/hyperelliptic_curves/hyperelliptic_finite_field.py", line 1348, in sage.schemes.hyperelliptic_curves.hyperelliptic_finite_field.HyperellipticCurve_finite_field.zeta_funct
ion
Failed example:
    H.zeta_function()
Exception raised:
    Traceback (most recent call last):
      File "/projects/b8cc019c-1204-44b1-bea9-eb81c119388e/sage/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 496, in _run
        self.compile_and_execute(example, compiler, test.globs)
      File "/projects/b8cc019c-1204-44b1-bea9-eb81c119388e/sage/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 858, in compile_and_execute
        exec(compiled, globs)
      File "<doctest sage.schemes.hyperelliptic_curves.hyperelliptic_finite_field.HyperellipticCurve_finite_field.zeta_function[8]>", line 1, in <module>
        H.zeta_function()
      File "/projects/b8cc019c-1204-44b1-bea9-eb81c119388e/sage/local/lib/python2.7/site-packages/sage/schemes/hyperelliptic_curves/hyperelliptic_finite_field.py", line 1364, in zeta_function
        P = self.frobenius_polynomial()
      File "sage/misc/cachefunc.pyx", line 2235, in sage.misc.cachefunc.CachedMethodCallerNoArgs.__call__ (/projects/b8cc019c-1204-44b1-bea9-eb81c119388e/sage/src/build/cythonized/sage/misc/cac
hefunc.c:12462)
        self.cache = f(self._instance)
      File "/projects/b8cc019c-1204-44b1-bea9-eb81c119388e/sage/local/lib/python2.7/site-packages/sage/schemes/hyperelliptic_curves/hyperelliptic_finite_field.py", line 551, in frobenius_polyno
mial
        return self.frobenius_polynomial_pari()
      File "/projects/b8cc019c-1204-44b1-bea9-eb81c119388e/sage/local/lib/python2.7/site-packages/sage/schemes/hyperelliptic_curves/hyperelliptic_finite_field.py", line 473, in frobenius_polyno
mial_pari
        return ZZ['x'](pari([f, h]).hyperellcharpoly())
      File "sage/libs/pari/auto_gen.pxi", line 8650, in sage.libs.pari.gen.gen_auto.hyperellcharpoly (/projects/b8cc019c-1204-44b1-bea9-eb81c119388e/sage/src/build/cythonized/sage/libs/pari/gen
.c:47546)
        sig_on()
      File "src/cysignals/signals.pyx", line 108, in cysignals.signals.sig_raise_exception (build/src/cysignals/signals.c:1448)
    SignalError: Segmentation fault

As noted above, this is due to a bug in PARI which has been fixed upstream. To resolve this, we need to pull from upstream, which is probably worth a separate ticket.

comment:30 Changed 4 years ago by jpflori

Created #20581 for the PARI/GP upgrade.

comment:31 Changed 4 years ago by kedlaya

  • Dependencies changed from #20219 to #20219, #20581

comment:32 Changed 4 years ago by kedlaya

  • Milestone changed from sage-6.8 to sage-7.3

comment:33 follow-up: Changed 4 years ago by slelievre

Now that #20581 is in, is this good to go? Or does it need the new Pari update at #21005?

comment:34 in reply to: ↑ 33 Changed 4 years ago by kedlaya

  • Status changed from new to needs_review

Replying to slelievre:

Now that #20581 is in, is this good to go? Or does it need the new Pari update at #21005?

It should be good to go with just #20581, but I haven't yet had a chance to recompile to see if the broken doctest is now fixed. But let's go ahead and make that part of the review...

comment:35 Changed 4 years ago by chapoton

missing blank line after

Over prime fields of odd characteristic, `f` may have even degree::

missing blank line after

This example shows that ticket #20391 is resolved:

which also should instead be

This example shows that :trac:`20391` is resolved::

comment:36 Changed 4 years ago by jpflori

@Frédéric: Do you plan on pushing your changes as a reviewer's action?

comment:37 Changed 4 years ago by chapoton

not today, so please go on

comment:38 Changed 4 years ago by git

  • Commit changed from 6ca4b8c808ba9a53a17ce061a96ecbabea58d595 to 5a046fe8c581a3e29894ba50706efc82e2af7462

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

f6e4aa1Implement wrapper for PARI hyperellcharpoly
834b5c5Correctly handle frobenius_matrix() for an even degree polynomial
2690f0dFix count_points for even-degree hyperelliptics
c2bd95dAdd doctest regarding #20391
625e636Merge branch 'u/kedlaya/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves' of git://trac.sagemath.org/sage into t/18916/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves
1f90f95Upgrade PARI to latest master
186c4b0Merge branch 'u/jdemeyer/upgrade_pari_to_latest_master' of git://trac.sagemath.org/sage into t/18916/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves
5a046feCorrect doctest formatting

comment:39 Changed 4 years ago by git

  • Commit changed from 5a046fe8c581a3e29894ba50706efc82e2af7462 to cf5ad9834607c2295df2c592b173a184b7091fdb

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

cf5ad98Merge branch 'develop' into t/18916/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves

comment:40 Changed 4 years ago by kedlaya

  • Report Upstream changed from Fixed upstream, but not in a stable release. to N/A

I fixed the formatting issues that were reported. Separately, I updated to the current (post-7.3) develop, which includes a PARI update, and can confirm that the broken doctest is now fixed.

I would set this back to needs_review except that it was never unset...

comment:41 Changed 4 years ago by chapoton

  • Branch changed from u/kedlaya/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves to public/18916
  • Commit changed from cf5ad9834607c2295df2c592b173a184b7091fdb to c814597edef29a03cb1ae787fb484d2e84dcc461

Here is a public branch with just a small commit on top of your branch,

with just some enhanced details in code and doc formatting.

This is not a full review, sorry.


New commits:

c814597trac #18916 reviewer's commit (details)

comment:42 Changed 4 years ago by jpflori

  • Description modified (diff)
  • Reviewers set to Jean-Pierre Flori, Frédéric Chapoton
  • Status changed from needs_review to positive_review

The changes look good to me.

comment:43 Changed 4 years ago by vbraun

  • Status changed from positive_review to needs_work
[sagelib-7.4.beta0] [  8/311] gcc -fno-strict-aliasing -I/Users/buildslave-sage/slave/sage_git/build/local/var/tmp/sage/build/python2-2.7.10.p2/include -DNDEBUG -g -fwrapv -O3 -Wall -Wno-unused -I/Users/buildslave-sage/slave/sage_git/build/local/lib/python2.7/site-packages/cysignals -I/Users/buildslave-sage/slave/sage_git/build/local/include -I/Users/buildslave-sage/slave/sage_git/build/local/include/python2.7 -I/Users/buildslave-sage/slave/sage_git/build/local/lib/python2.7/site-packages/numpy/core/include -I/Users/buildslave-sage/slave/sage_git/build/src -I/Users/buildslave-sage/slave/sage_git/build/src/sage/ext -I/Users/buildslave-sage/slave/sage_git/build/src/build/cythonized -I/Users/buildslave-sage/slave/sage_git/build/src/build/cythonized/sage/ext -I/Users/buildslave-sage/slave/sage_git/build/local/include/python2.7 -c /Users/buildslave-sage/slave/sage_git/build/src/build/cythonized/sage/libs/pari/handle_error.c -o build/temp.macosx-10.9-x86_64-2.7/Users/buildslave-sage/slave/sage_git/build/src/build/cythonized/sage/libs/pari/handle_error.o -fno-strict-aliasing
[sagelib-7.4.beta0] In file included from /Users/buildslave-sage/slave/sage_git/build/local/include/c++/4.9.3/backward/strstream:51:0,
[sagelib-7.4.beta0]                  from /Users/buildslave-sage/slave/sage_git/build/local/include/libLfunction/L.h:34,
[sagelib-7.4.beta0]                  from /Users/buildslave-sage/slave/sage_git/build/src/sage/libs/lcalc/lcalc_sage.h:1,
[sagelib-7.4.beta0]                  from /Users/buildslave-sage/slave/sage_git/build/src/build/cythonized/sage/libs/lcalc/lcalc_Lfunction.cpp:338:
[sagelib-7.4.beta0] /Users/buildslave-sage/slave/sage_git/build/local/include/c++/4.9.3/backward/backward_warning.h:32:2: warning: #warning This file includes at least one deprecated or antiquated header which may be removed without further notice at a future date. Please use a non-deprecated interface with equivalent functionality instead. For a listing of replacement headers and interfaces, consult the file backward_warning.h. To disable this warning use -Wno-deprecated. [-Wcpp]
[sagelib-7.4.beta0]  #warning \
[sagelib-7.4.beta0]   ^
[sagelib-7.4.beta0] In file included from /Users/buildslave-sage/slave/sage_git/build/src/build/cythonized/sage/libs/pari/gen.c:341:0:
[sagelib-7.4.beta0] /Users/buildslave-sage/slave/sage_git/build/local/include/pari/paripriv.h:16:48: error: unknown type name 'ulong'

on 32-bit

comment:44 Changed 4 years ago by git

  • Commit changed from c814597edef29a03cb1ae787fb484d2e84dcc461 to e79317e72f2e50cb94e00e2a91fd643da644875e

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

ed5c343Remove spurious changes to PARI headers
e79317eMerge branch 'public/18916' of git://trac.sagemath.org/sage into t/18916/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves

comment:45 Changed 4 years ago by kedlaya

  • Status changed from needs_work to needs_review

Between pulling from the PARI-related tickets and updating develop, a couple of unwanted changes to the PARI headers seem to have crept in, causing the build warning. I just removed these; here we go again.

comment:46 Changed 4 years ago by jpflori

  • Status changed from needs_review to positive_review

I can confirm it now builds on my setup with 7.4.beta0.

comment:47 Changed 4 years ago by vbraun

  • Branch changed from public/18916 to e79317e72f2e50cb94e00e2a91fd643da644875e
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.