Opened 3 years ago

Closed 11 months ago

#23946 closed enhancement (fixed)

Exhaust over Weil polynomials

Reported by: kedlaya Owned by:
Priority: minor Milestone: sage-9.1
Component: number theory Keywords: Weil polynomials, sd91
Cc: roed Merged in:
Authors: Kiran Kedlaya Reviewers: David Roe
Report Upstream: N/A Work issues:
Branch: e02ae3c (Commits) Commit: e02ae3cfc9bd4da34ee8f3b3b9a9b94e7666fa6f
Dependencies: #24016 Stopgaps:

Description

I have reasonably stable code (mostly in C, depending on FLINT) for exhausting over Weil polynomials:

https://github.com/kedlaya/root-unitary

The goal of this ticket is to incorporate this code into Sage in some fashion.

Change History (50)

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

  • Dependencies set to #24016

Status update: I've pushed many recent changes to the branch repackage_for_sage of the Github repo. In particular, I have been working to make the call syntax much more transparent; for instance, to compute the Weil polynomials corresponding to abelian varieties of dimension g over the finite field of order q, one can now use standard Python iterator syntax:

sage: wp = WeilPolynomials(2*g, q)
sage: for i in wp:
...

Some to-do items: check that all my old test scripts still work; write more tests; add more documentation; add more input sanitization.

In the meantime, I'd appreciate some feedback about where this might belong in Sage.

comment:2 in reply to: ↑ 1 Changed 3 years ago by jdemeyer

Replying to kedlaya:

In the meantime, I'd appreciate some feedback about where this might belong in Sage.

First of all, your github repo contains both C code and .sage files. Are the .sage files considered to be part of your package or are they really only use cases/examples/tests/documentation?

There are two possibilities, each with its advantages and disadvantages:

(A) Make it an external (presumably optional) package and interface it from Sage.

(B) Make it actually a part of Sage itself.

Advantages of (A) are that you keep control of the package: you can easily change the package at will (with (B) every change would need to go through the Sage Trac) and you don't need to care about Sage coding standards. Also important: with (A) your package could be used without Sage. With (B), there might be more initial work needed to get it into Sage but then it might be less work to maintain.

comment:3 Changed 3 years ago by kedlaya

I just merged a pull request in from another branch, which has the effect of significantly reconfiguring the file structure. In particular, there are no longer any .sage files in the main directory; the ones that exist are all in subdirectories and are indeed examples/tests.

The main code base now consists of two .c files (and their .h headers) and one .pyx file. I still don't have any sort of independent build structure set up; moreover, the .pyx file is more than just a simple wrapper, it is where the parallelism is currently implemented. So at the moment there is really no way to use the code without Sage, and I don't particularly intend to work on making that possible. (Exception: if someone wants to use the code from some other platform like Julia, I'd be willing to cooperate with that.)

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

I think that the code is stable enough that we should just add it to Sage. Maybe create a folder sage/rings/polynomial/weil/ and put them in there?

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

Replying to roed:

I think that the code is stable enough that we should just add it to Sage. Maybe create a folder sage/rings/polynomial/weil/ and put them in there?

Do you mean put everything in there? Or put the Cython in sage/rings/polynomials/weil_polynomials.pyx and the C files in sage/rings/polynomial/weil?

Also, if someone can point me to some guidance about a module to the Sage library (e.g., how to format __init__.py and all.py and the like to make the new code discoverable), that would be helpful.

comment:6 in reply to: ↑ 5 Changed 3 years ago by roed

Replying to kedlaya:

Replying to roed:

I think that the code is stable enough that we should just add it to Sage. Maybe create a folder sage/rings/polynomial/weil/ and put them in there?

Do you mean put everything in there? Or put the Cython in sage/rings/polynomials/weil_polynomials.pyx and the C files in sage/rings/polynomial/weil?

I meant putting everything in there. We should also talk about the best way for a user to access these methods; perhaps a method weil_polynomials on ZZ[x] and QQ[x] that returns your iterator?

Also, if someone can point me to some guidance about a module to the Sage library (e.g., how to format __init__.py and all.py and the like to make the new code discoverable), that would be helpful.

See http://doc.sagemath.org/html/en/developer/coding_basics.html#files-and-directory-structure for some details. You'll also need to add lines to src/module_list.py for each extension module.

comment:7 Changed 3 years ago by kedlaya

  • Branch set to u/kedlaya/exhaust_over_weil_polynomials

comment:8 Changed 3 years ago by kedlaya

  • Commit set to d12b0ea0aa53cbf558e6ac873414778dd12796a4

Pushed some partial progress, but this doesn't compile yet.


New commits:

da3d483Add directory src/sage/rings/polynomial/weil/
d12b0eaUpdate module list, all.py

comment:9 Changed 12 months ago by kedlaya

New commits:

da3d483Add directory src/sage/rings/polynomial/weil/
d12b0eaUpdate module list, all.py

comment:10 Changed 12 months ago by kedlaya

  • Branch u/kedlaya/exhaust_over_weil_polynomials deleted
  • Commit d12b0ea0aa53cbf558e6ac873414778dd12796a4 deleted

comment:11 Changed 12 months ago by kedlaya

  • Branch set to u/kedlaya/weilpoly_search

comment:12 Changed 12 months ago by kedlaya

  • Commit set to 2a62e66aa58ccd7e82740078ce2c71d62a10ea52

OK, here is a new branch which seems to compile. Still a lot of polishing to be done, but it should be possible to at least play with the code now.


New commits:

2bfaf04Initial commit
d7e1e2aSuccesful compilation
2a62e66Bug fix from upstream (incorrect handling of leading coeffs)

comment:13 Changed 12 months ago by chapoton

does not work for me:

sage: from sage.rings.polynomial.weil.weil_polynomials import *
---------------------------------------------------------------------------
ImportError                               Traceback (most recent call last)
<ipython-input-1-5b4f8184ba57> in <module>()
----> 1 from sage.rings.polynomial.weil.weil_polynomials import *

ImportError: /home/chapoton/sage3/local/lib/python3.7/site-packages/sage/rings/polynomial/weil/weil_polynomials.cpython-37m-x86_64-linux-gnu.so: undefined symbol: GOMP_loop_nonmonotonic_dynamic_start

and there are missing empty lines in the doc after lines ending with ::

comment:14 Changed 12 months ago by kedlaya

  • Milestone changed from sage-8.1 to sage-9.1

That is the code looking for OpenMP (which is optional, since we don't ship it with Sage). Does it help if you take out the compiler directive in line 2 of sage/rings/polynomial/weil/weil_polynomials.pyx?

The docstrings need a fair bit of work to get to 100% doctest coverage, let alone proper formatting. I can of course work on that independently.

comment:15 Changed 12 months ago by chapoton

yes, removing the distutils directive make it work

comment:16 Changed 12 months ago by git

  • Commit changed from 2a62e66aa58ccd7e82740078ce2c71d62a10ea52 to 1438c238cc5ab5aca7b0f558ff0261f7112ca9ec

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

cd2f537Eliminate compiler warnings about return from a void function
1438c23Suppress OpenMP, add doctests

comment:17 Changed 12 months ago by kedlaya

  • Cc roed added

In addition to suppressing OpenMP, I also added doctests to get to 100% coverage.

Cc'ing roed, who is somewhat familiar to the code.

comment:18 Changed 12 months ago by chapoton

please remove the "next" method.

and in the test of __next__, use next(it) to call the iterator.

comment:19 Changed 12 months ago by git

  • Commit changed from 1438c238cc5ab5aca7b0f558ff0261f7112ca9ec to 7516836408314c97f176d62d677b7d2d924e861b

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

7516836Removed "next" method

comment:20 Changed 12 months ago by git

  • Commit changed from 7516836408314c97f176d62d677b7d2d924e861b to 1a6feb9edc7b3621738f1dac42ea0dbeea7282c3

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

1a6feb9Merge edits from Edgar Costa

comment:21 Changed 11 months ago by chapoton

needs review ?

comment:22 Changed 11 months ago by kedlaya

Standing by for some feedback from roed first.

comment:23 Changed 11 months ago by roed

  • Branch changed from u/kedlaya/weilpoly_search to u/roed/weilpoly_search

comment:24 Changed 11 months ago by roed

  • Commit changed from 1a6feb9edc7b3621738f1dac42ea0dbeea7282c3 to c1be1fbff4f97dcd107be3e92c6e3b676afdbf00

About to board a flight, and I want to make a few more changes, but can you take a look at the documentation I added and make sure it looks right? It would also be nice to see some examples illustrating the different inputs available to the iterator.


New commits:

c1be1fbSome Weil polynomial changes

comment:25 Changed 11 months ago by git

  • Commit changed from c1be1fbff4f97dcd107be3e92c6e3b676afdbf00 to d88c654e80d164044cfb3b946f2475ef6dd1ece2

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

d88c654Fix some raise/return problems, add weil_polynomials method for polynomial rings, add some documentation

comment:26 follow-up: Changed 11 months ago by roed

Alright, I made the other changes I wanted to make (most notably adding a method to polynomial rings). Here are some documentation questions that I'd like to see addressed, but I think this ticket is basically ready to go.

  • I don't understand how specifying congruence constraints on leading coefficients works. I tried the following:
    sage: R.<T> = ZZ[]
    sage: R.weil_polynomials(4,2,lead=((1,0),(0,2)))
    [T^4 + 4, T^4 - T^2 + 4, T^4 - 2*T^2 + 4, T^4 - 3*T^2 + 4, T^4 - 4*T^2 + 4]
    

I would expect this to give polynomials where the coefficient of T^3 was even, but I'm not sure what I'm getting. Is it ignoring the 2 and giving me those where the coefficient of T^3 is exactly 0? Kiran, can you give me an idea of how this functionality is supposed to work?

  • I don't quite understand what node_limit does. Raises a RuntimeError if the number of terminal nodes ever exceeds the given value?
  • Similarly, I guessed that squarefree omits polynomials that are not squarefree.

Overall, can you check/correct the changes I made? Also, if there's anything else useful to be said about how to recompile with OpenMP support, feel free to augment the documentation that I added.

comment:27 Changed 11 months ago by kedlaya

  • Branch changed from u/roed/weilpoly_search to u/kedlaya/weilpoly_search

comment:28 Changed 11 months ago by git

  • Commit changed from d88c654e80d164044cfb3b946f2475ef6dd1ece2 to 0aaec041320d09285a7d95171147d32f4304e850

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

0aaec04Updated some docstrings

comment:29 Changed 11 months ago by git

  • Commit changed from 0aaec041320d09285a7d95171147d32f4304e850 to 6f48bf4e43e3afd232f719a760fb877c97951853

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

6f48bf4Fix squarefree search, add doctest

comment:30 in reply to: ↑ 26 Changed 11 months ago by kedlaya

  • Authors set to Kiran Kedlaya
  • Status changed from new to needs_review

Replying to roed:

Alright, I made the other changes I wanted to make (most notably adding a method to polynomial rings). Here are some documentation questions that I'd like to see addressed, but I think this ticket is basically ready to go.

I made a few tweaks in the latest commits. I really like having a method for polynomial rings; I'm guessing most users will only use the code this way. The iterator is really only needed for power use cases.

Since there is already an is_weil_polynomial method for polynomial elements, I edited the docstring for that method to point to weil_polynomials and vice versa.

  • I don't understand how specifying congruence constraints on leading coefficients works. I tried the following:
    sage: R.<T> = ZZ[]
    sage: R.weil_polynomials(4,2,lead=((1,0),(0,2)))
    [T^4 + 4, T^4 - T^2 + 4, T^4 - 2*T^2 + 4, T^4 - 3*T^2 + 4, T^4 - 4*T^2 + 4]
    

I would expect this to give polynomials where the coefficient of T^3 was even, but I'm not sure what I'm getting. Is it ignoring the 2 and giving me those where the coefficient of T^3 is exactly 0? Kiran, can you give me an idea of how this functionality is supposed to work?

You just caught a bug in the C code which I have now fixed (here and upstream). I added a doctest documenting how this is supposed to work.

  • I don't quite understand what node_limit does. Raises a RuntimeError if the number of terminal nodes ever exceeds the given value?

Correct. It is there as a failsafe in case you don't want your code to run forever.

  • Similarly, I guessed that squarefree omits polynomials that are not squarefree.

Correct, except that there was a (long-known) bug in that too. I fixed that and added a doctest.

Overall, can you check/correct the changes I made? Also, if there's anything else useful to be said about how to recompile with OpenMP support, feel free to augment the documentation that I added.

One issue is that on my system, I have to put the distutils lines at the top of the code when uncommenting them, which I think runs contrary to Sage conventions. I added a comment about this.

Since it seems we are basically into the review phase anyway, let me go ahead and set state accordingly. I would be particularly open to suggestions for additional doctests; there are probably many use cases I haven't thought to test yet, and I'd like to make sure not to break those in future versions of the code.

comment:31 Changed 11 months ago by roed

  • Branch changed from u/kedlaya/weilpoly_search to u/roed/weilpoly_search

comment:32 Changed 11 months ago by git

  • Commit changed from 6f48bf4e43e3afd232f719a760fb877c97951853 to d8a813f7df1d503ef41e16c28beb2f77907fe444

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

d8a813fAdd some weil polynomial tests

comment:33 Changed 11 months ago by roed

  • Reviewers set to David Roe

I've added some more tests. I think I'm happy to give this positive review if all tests are passing. I'm also happy to try to keep brainstorming some more doctests if you'd like.

comment:34 Changed 11 months ago by kedlaya

  • Branch changed from u/roed/weilpoly_search to u/kedlaya/weilpoly_search

comment:35 Changed 11 months ago by git

  • Commit changed from d8a813f7df1d503ef41e16c28beb2f77907fe444 to dbb60955f2e17b0b73dbd0353f4f8bfed223916b

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

dbb6095Fix compiler warnings

comment:36 Changed 11 months ago by git

  • Commit changed from dbb60955f2e17b0b73dbd0353f4f8bfed223916b to e0b9ff9bcd8d3b5ef3f1c6b40fd6e6f6a1c175a9

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

e0b9ff9Add ref to github repo for Weil polynomials

comment:37 Changed 11 months ago by kedlaya

A few last-minute changes:

  • There was a trivial bug in polynomial_ring.py about checking the base ring (which patchbot also noticed); I fixed that and added some doctests to confirm.
  • I fixed some compiler warnings on the backend.
  • I added a link to the github repo in the main package docstring.

That's probably enough doctests for the moment; I'll inevitably have to add more in response to later bugfixes. Let's see what patchbot has to say next.

comment:38 Changed 11 months ago by roed

  • Milestone changed from sage-9.1 to sage-9.0
  • Status changed from needs_review to positive_review

Patchbot is green; I think we're good to go! I'll optimistically set the milestone to 9.0. :-)

comment:39 follow-up: Changed 11 months ago by vbraun

On Python 2 there are various errors of the form

**********************************************************************
File "src/sage/rings/polynomial/weil/weil_polynomials.pyx", line 168, in sage.rings.polynomial.weil.weil_polynomials.dfs_manager.node_count
Failed example:
    _ = next(it)
Exception raised:
    Traceback (most recent call last):
      File "/var/lib/buildbot/slave/sage2_git/build/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 681, in _run
        self.compile_and_execute(example, compiler, test.globs)
      File "/var/lib/buildbot/slave/sage2_git/build/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 1123, in compile_and_execute
        exec(compiled, globs)
      File "<doctest sage.rings.polynomial.weil.weil_polynomials.dfs_manager.node_count[3]>", line 1, in <module>
        _ = next(it)
    TypeError: instance has no next() method
**********************************************************************

comment:40 Changed 11 months ago by vbraun

  • Status changed from positive_review to needs_work

comment:41 Changed 11 months ago by git

  • Commit changed from e0b9ff9bcd8d3b5ef3f1c6b40fd6e6f6a1c175a9 to 7f44a49a42de8291e091ca2f4c516111282888b0

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

7f44a49Reinstate next() method for Py2 compatibility

comment:42 in reply to: ↑ 39 Changed 11 months ago by kedlaya

  • Status changed from needs_work to needs_review

Replying to vbraun:

On Python 2 there are various errors of the form

**********************************************************************
File "src/sage/rings/polynomial/weil/weil_polynomials.pyx", line 168, in sage.rings.polynomial.weil.weil_polynomials.dfs_manager.node_count
Failed example:
    _ = next(it)
Exception raised:
    Traceback (most recent call last):
      File "/var/lib/buildbot/slave/sage2_git/build/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 681, in _run
        self.compile_and_execute(example, compiler, test.globs)
      File "/var/lib/buildbot/slave/sage2_git/build/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 1123, in compile_and_execute
        exec(compiled, globs)
      File "<doctest sage.rings.polynomial.weil.weil_polynomials.dfs_manager.node_count[3]>", line 1, in <module>
        _ = next(it)
    TypeError: instance has no next() method
**********************************************************************

I was asked to remove that method (see comment:18), but anyway here it is again.

Last edited 11 months ago by kedlaya (previous) (diff)

comment:43 Changed 11 months ago by git

  • Commit changed from 7f44a49a42de8291e091ca2f4c516111282888b0 to 6a3e93d0c200fa3d381c7d6b6433134439e04b33

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

6a3e93dFix typos in doctest

comment:44 Changed 11 months ago by kedlaya

Looking at the patchbot logs, it seems that at least one of the bots failed to build because of nested function declarations in the C code. I thought we were always compiling the Sage library with gcc, which supports nested functions even though they are not in the ANSI standard...?

comment:45 Changed 11 months ago by vbraun

on osx we now use clang by default

comment:46 Changed 11 months ago by git

  • Commit changed from 6a3e93d0c200fa3d381c7d6b6433134439e04b33 to e02ae3cfc9bd4da34ee8f3b3b9a9b94e7666fa6f

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

e02ae3cUn-nest functions in C code for non-GNU compilers

comment:47 Changed 11 months ago by kedlaya

Ah, it seems clang does not (yet) support nested functions as in GNU C. So I redid that bit of the C code; let's see what the patchbots have to say now.

comment:48 Changed 11 months ago by kedlaya

Some of the patchbots are set up strangely, but I'm not seeing any systemic issues with this branch.

comment:49 Changed 11 months ago by roed

  • Milestone changed from sage-9.0 to sage-9.1
  • Status changed from needs_review to positive_review

I agree. We have successful tests run on Darwin, so I'm setting this back to positive review.

comment:50 Changed 11 months ago by vbraun

  • Branch changed from u/kedlaya/weilpoly_search to e02ae3cfc9bd4da34ee8f3b3b9a9b94e7666fa6f
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.