Opened 7 years ago
Closed 6 years ago
#18916 closed enhancement (fixed)
Use Kedlaya algorithm to count points on hyperelliptic curves
Reported by:  jpflori  Owned by:  

Priority:  major  Milestone:  sage7.3 
Component:  number fields  Keywords:  hyperelliptic curves, matrix of Frobenius 
Cc:  defeo, katestange  Merged in:  
Authors:  Kiran Kedlaya  Reviewers:  JeanPierre Flori, Frédéric Chapoton 
Report Upstream:  N/A  Work issues:  
Branch:  e79317e (Commits, GitHub, GitLab)  Commit:  e79317e72f2e50cb94e00e2a91fd643da644875e 
Dependencies:  #20219, #20581  Stopgaps: 
Description (last modified by )
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 7 years ago by
comment:2 Changed 7 years ago by
+1 to use PARI/GP
comment:3 followup: ↓ 5 Changed 6 years ago by
 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 6 years ago by
Addendum: naive point counting is already implemented directly in Sage, see #11980.
comment:5 in reply to: ↑ 3 ; followup: ↓ 6 Changed 6 years ago by
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 ; followup: ↓ 7 Changed 6 years ago by
Replying to jpflori:
I confirm that
hyperellpadicfrobenius
only works for prime fields whereas thenf...
one handles extensions. Though both of them are not markedstatic
and usegerepile
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 6 years ago by
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 followup: ↓ 9 Changed 6 years ago by
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 ; followup: ↓ 14 Changed 6 years ago by
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 insrc/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 6 years ago by
If you do make any changes to the PARI interface, please base them on #20219.
comment:11 Changed 6 years ago by
 Dependencies set to #20219
comment:12 Changed 6 years ago by
 Branch set to u/kedlaya/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves
comment:13 Changed 6 years ago by
 Commit set to 6a4e932eef76137b6856c3b84f2fe82de807cc3a
Branch pushed to git repo; I updated commit sha1. Last 10 new commits:
c0ed97a  Stop using deprecated PARI factoring features

edc5ce2  Merge branch 't/20205/get_rid_of_factorint_withproof_sage_in_pari_interface' into HEAD

5fb408d  Replace pari_catch_sig_on by sig_on

d5c934c  Deprecate PARI nth_prime and prime_list

d7d2d7d  Remove redundant functions from pari_instance.pyx

130fc90  Remove redundant functions from gen.pyx

2299743  Use SR.symbol() instead var()

6acc79a  Add back order() as deprecated alias of multiplicative_order()

0f2386a  Fix basic t_LIST functionality

6a4e932  Merge 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 ; followup: ↓ 15 Changed 6 years ago by
Replying to jpflori:
I mentioned above that
nfhyperellpadicfrobenius
is not currently exposed. However, this function does appear insrc/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 6 years ago by
Replying to kedlaya:
Replying to jpflori:
I mentioned above that
nfhyperellpadicfrobenius
is not currently exposed. However, this function does appear insrc/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 whyhyperellpadicfrobenius
(which does not have its own explicitly defined wrapper as far as I can find) shows up this way whilenfhyperellpadicfrobenius
does not. I had thoughtparidecl.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 autogenerated 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 6 years ago by
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 ColemanGross height pairing. But I think the padic sigma function on an elliptic curve does constitute an example where the Frobenius matrix is strictly more useful than its charpoly.
comment:17 Changed 6 years ago by
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 followup: ↓ 20 Changed 6 years ago by
 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 toplevel: hyperellcharpoly(x^3 *** ^ *** hyperellcharpoly: bug in PARI/GP (Segmentation Fault), please report.
Off to go report this upstream...
comment:19 Changed 6 years ago by
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 6 years ago by
 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 toplevel: 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.ubordeaux1.fr> Date: Fri Mar 11 12:38:30 2016 +0100 hyperellcharpoly: fix typo in ZXX_to_FpXC
comment:21 Changed 6 years ago by
 Commit changed from 6a4e932eef76137b6856c3b84f2fe82de807cc3a to e510540635d871837fd5de0ae47f3fdbf51f07e6
Branch pushed to git repo; I updated commit sha1. New commits:
e510540  Implement wrapper for PARI hyperellcharpoly

comment:22 Changed 6 years ago by
 Commit changed from e510540635d871837fd5de0ae47f3fdbf51f07e6 to 95ede4bc4354e8717bfb38143d8b5ebb168bf878
Branch pushed to git repo; I updated commit sha1. New commits:
95ede4b  Correctly handle frobenius_matrix() for an even degree polynomial

comment:23 Changed 6 years ago by
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 6 years ago by
 Commit changed from 95ede4bc4354e8717bfb38143d8b5ebb168bf878 to 4c92503eeaed4d9cca5016ee5fec2c6082eb65d8
Branch pushed to git repo; I updated commit sha1. New commits:
4c92503  Fix count_points for evendegree hyperelliptics

comment:26 Changed 6 years ago by
 Commit changed from 4c92503eeaed4d9cca5016ee5fec2c6082eb65d8 to 6ca4b8c808ba9a53a17ce061a96ecbabea58d595
Branch pushed to git repo; I updated commit sha1. New commits:
6ca4b8c  Add doctest regarding #20391

comment:27 Changed 6 years ago by
... and a doctest to indicate that #20391 is addressed.
comment:28 followup: ↓ 29 Changed 6 years ago by
@Kiran: is it in state of review or do you plan further changes?
comment:29 in reply to: ↑ 28 Changed 6 years ago by
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/b8cc019c120444b1bea9eb81c119388e/sage/local/lib/python2.7/sitepackages/sage/doctest/forker.py", line 496, in _run self.compile_and_execute(example, compiler, test.globs) File "/projects/b8cc019c120444b1bea9eb81c119388e/sage/local/lib/python2.7/sitepackages/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/b8cc019c120444b1bea9eb81c119388e/sage/local/lib/python2.7/sitepackages/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/b8cc019c120444b1bea9eb81c119388e/sage/src/build/cythonized/sage/misc/cac hefunc.c:12462) self.cache = f(self._instance) File "/projects/b8cc019c120444b1bea9eb81c119388e/sage/local/lib/python2.7/sitepackages/sage/schemes/hyperelliptic_curves/hyperelliptic_finite_field.py", line 551, in frobenius_polyno mial return self.frobenius_polynomial_pari() File "/projects/b8cc019c120444b1bea9eb81c119388e/sage/local/lib/python2.7/sitepackages/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/b8cc019c120444b1bea9eb81c119388e/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 6 years ago by
Created #20581 for the PARI/GP upgrade.
comment:31 Changed 6 years ago by
 Dependencies changed from #20219 to #20219, #20581
comment:32 Changed 6 years ago by
 Milestone changed from sage6.8 to sage7.3
comment:33 followup: ↓ 34 Changed 6 years ago by
comment:34 in reply to: ↑ 33 Changed 6 years ago by
 Status changed from new to needs_review
comment:35 Changed 6 years ago by
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 6 years ago by
@Frédéric: Do you plan on pushing your changes as a reviewer's action?
comment:37 Changed 6 years ago by
not today, so please go on
comment:38 Changed 6 years ago by
 Commit changed from 6ca4b8c808ba9a53a17ce061a96ecbabea58d595 to 5a046fe8c581a3e29894ba50706efc82e2af7462
Branch pushed to git repo; I updated commit sha1. New commits:
f6e4aa1  Implement wrapper for PARI hyperellcharpoly

834b5c5  Correctly handle frobenius_matrix() for an even degree polynomial

2690f0d  Fix count_points for evendegree hyperelliptics

c2bd95d  Add doctest regarding #20391

625e636  Merge 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

1f90f95  Upgrade PARI to latest master

186c4b0  Merge 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

5a046fe  Correct doctest formatting

comment:39 Changed 6 years ago by
 Commit changed from 5a046fe8c581a3e29894ba50706efc82e2af7462 to cf5ad9834607c2295df2c592b173a184b7091fdb
Branch pushed to git repo; I updated commit sha1. New commits:
cf5ad98  Merge branch 'develop' into t/18916/use_kedlaya_algorithm_to_count_points_on_hyperelliptic_curves

comment:40 Changed 6 years ago by
 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 (post7.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 6 years ago by
 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:
c814597  trac #18916 reviewer's commit (details)

comment:42 Changed 6 years ago by
 Description modified (diff)
 Reviewers set to JeanPierre Flori, Frédéric Chapoton
 Status changed from needs_review to positive_review
The changes look good to me.
comment:43 Changed 6 years ago by
 Status changed from positive_review to needs_work
[sagelib7.4.beta0] [ 8/311] gcc fnostrictaliasing I/Users/buildslavesage/slave/sage_git/build/local/var/tmp/sage/build/python22.7.10.p2/include DNDEBUG g fwrapv O3 Wall Wnounused I/Users/buildslavesage/slave/sage_git/build/local/lib/python2.7/sitepackages/cysignals I/Users/buildslavesage/slave/sage_git/build/local/include I/Users/buildslavesage/slave/sage_git/build/local/include/python2.7 I/Users/buildslavesage/slave/sage_git/build/local/lib/python2.7/sitepackages/numpy/core/include I/Users/buildslavesage/slave/sage_git/build/src I/Users/buildslavesage/slave/sage_git/build/src/sage/ext I/Users/buildslavesage/slave/sage_git/build/src/build/cythonized I/Users/buildslavesage/slave/sage_git/build/src/build/cythonized/sage/ext I/Users/buildslavesage/slave/sage_git/build/local/include/python2.7 c /Users/buildslavesage/slave/sage_git/build/src/build/cythonized/sage/libs/pari/handle_error.c o build/temp.macosx10.9x86_642.7/Users/buildslavesage/slave/sage_git/build/src/build/cythonized/sage/libs/pari/handle_error.o fnostrictaliasing [sagelib7.4.beta0] In file included from /Users/buildslavesage/slave/sage_git/build/local/include/c++/4.9.3/backward/strstream:51:0, [sagelib7.4.beta0] from /Users/buildslavesage/slave/sage_git/build/local/include/libLfunction/L.h:34, [sagelib7.4.beta0] from /Users/buildslavesage/slave/sage_git/build/src/sage/libs/lcalc/lcalc_sage.h:1, [sagelib7.4.beta0] from /Users/buildslavesage/slave/sage_git/build/src/build/cythonized/sage/libs/lcalc/lcalc_Lfunction.cpp:338: [sagelib7.4.beta0] /Users/buildslavesage/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 nondeprecated interface with equivalent functionality instead. For a listing of replacement headers and interfaces, consult the file backward_warning.h. To disable this warning use Wnodeprecated. [Wcpp] [sagelib7.4.beta0] #warning \ [sagelib7.4.beta0] ^ [sagelib7.4.beta0] In file included from /Users/buildslavesage/slave/sage_git/build/src/build/cythonized/sage/libs/pari/gen.c:341:0: [sagelib7.4.beta0] /Users/buildslavesage/slave/sage_git/build/local/include/pari/paripriv.h:16:48: error: unknown type name 'ulong'
on 32bit
comment:44 Changed 6 years ago by
 Commit changed from c814597edef29a03cb1ae787fb484d2e84dcc461 to e79317e72f2e50cb94e00e2a91fd643da644875e
comment:45 Changed 6 years ago by
 Status changed from needs_work to needs_review
Between pulling from the PARIrelated 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 6 years ago by
 Status changed from needs_review to positive_review
I can confirm it now builds on my setup with 7.4.beta0.
comment:47 Changed 6 years ago by
 Branch changed from public/18916 to e79317e72f2e50cb94e00e2a91fd643da644875e
 Resolution set to fixed
 Status changed from positive_review to closed
The magic line reads something like:
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)