Opened 6 years ago
Closed 6 years ago
#17695 closed enhancement (fixed)
Slightly fasten is_cyclotomic()
Reported by:  bruno  Owned by:  

Priority:  minor  Milestone:  sage6.5 
Component:  number theory  Keywords:  cyclotomic polynomials 
Cc:  chapoton  Merged in:  
Authors:  Bruno Grenet  Reviewers:  Frédéric Chapoton 
Report Upstream:  N/A  Work issues:  
Branch:  a2a6f02 (Commits, GitHub, GitLab)  Commit:  a2a6f026794fb12e1d557d8ca6738c45d5ed937f 
Dependencies:  Stopgaps: 
Description (last modified by )
In src/sage/rings/polynomial/polynomial_element.pyx
, the method is_cyclotomic()
begins by checking the irreducibility of the input polynomial, and then checks the leading and constant coefficients. I simply propose to invert those tests to fasten a bit the test for polynomials that are not monic or whose constant coefficient is not 1
.
Change History (31)
comment:1 Changed 6 years ago by
 Branch set to u/bruno/slightly_fasten_is_cyclotomic__
comment:2 Changed 6 years ago by
 Commit set to f79eeaaa630b1bd966a25a2f7b0f7a5f1dda3378
 Description modified (diff)
 Keywords cyclotomic polynomials added
 Status changed from new to needs_review
comment:3 followup: ↓ 5 Changed 6 years ago by
comment:4 Changed 6 years ago by
 Cc chapoton added
comment:5 in reply to: ↑ 3 ; followup: ↓ 6 Changed 6 years ago by
Replying to vdelecroix:
Hello Bruno,
Actually, I discussed this with Karim Belabas in Bordeaux few weeks ago. We should either call Pari or stick to their implementation. It is much faster:
sage: p = cyclotomic_polynomial(150) sage: pp = p._gp_() sage: timeit("pp.poliscyclo()") 625 loops, best of 3: 697 µs per loop sage: timeit("p.is_cyclotomic()") 125 loops, best of 3: 2.78 ms per loopMoreover they do have the function
poliscycloprod
that checks if all roots are roots of unity (which is even faster).Vincent
It is not surprising indeed that Pari is faster for this. I don't understand what you mean by
either call Pari or stick to their implementation.
Do you mean there are two ways to use Pari's implementation?
comment:6 in reply to: ↑ 5 ; followup: ↓ 9 Changed 6 years ago by
Replying to bruno:
Replying to vdelecroix: It is not surprising indeed that Pari is faster for this. I don't understand what you mean by
either call Pari or stick to their implementation.
Do you mean there are two ways to use Pari's implementation?
Nope. Either we should copy their implementation within Sage (I do not think this is justified) or just call Pari.
And by the way there are two ways of using pari. Either with gp
(this is a pexpect interface with a lot of latency in communication) and with libpari
(which is very efficient in C). But sadly there is not yet poliscyclo
and poliscycloprod
in the libpari interface (i.e. in the file libs/pari/gen.pyx
). I guess that we should wait for #17631 and not modify gen.pyx
directly.
Vincent
comment:7 followup: ↓ 8 Changed 6 years ago by
A side remark:
I would like to have a certificate (something like : this is the cyclotomic polynomial Phi_20)
Is this allowed by the pari function ?
comment:8 in reply to: ↑ 7 Changed 6 years ago by
Replying to chapoton:
A side remark:
I would like to have a certificate (something like : this is the cyclotomic polynomial Phi_20)
Is this allowed by the pari function ?
Of course ;) From pari documentation
poliscyclo(f): returns 0 if f is not a cyclotomic polynomial, and n > 0 if f = Phi_n, the nth cyclotomic polynomial.
and you can already use it
sage: R.<x> = ZZ[] sage: (x^2+x+1)._gp_().poliscyclo() 3 sage: (x^21)._gp_().poliscyclo() 0
Vincent
comment:9 in reply to: ↑ 6 ; followup: ↓ 10 Changed 6 years ago by
Replying to vdelecroix:
Replying to bruno:
Replying to vdelecroix: It is not surprising indeed that Pari is faster for this. I don't understand what you mean by
either call Pari or stick to their implementation.
Do you mean there are two ways to use Pari's implementation?
Nope. Either we should copy their implementation within Sage (I do not think this is justified) or just call Pari.
And by the way there are two ways of using pari. Either with
gp
(this is a pexpect interface with a lot of latency in communication) and withlibpari
(which is very efficient in C). But sadly there is not yetpoliscyclo
andpoliscycloprod
in the libpari interface (i.e. in the filelibs/pari/gen.pyx
). I guess that we should wait for #17631 and not modifygen.pyx
directly.Vincent
For concreteness, what do you propose to do now? If I understand correctly, you propose to wait for #17631 and then make is_cyclotomic()
directly call the pari functions, and not introduce the tiny change I made? I am fine with that, just wanted to be sure!
comment:10 in reply to: ↑ 9 ; followup: ↓ 11 Changed 6 years ago by
Replying to bruno:
For concreteness, what do you propose to do now? If I understand correctly, you propose to wait for #17631 and then make
is_cyclotomic()
directly call the pari functions, and not introduce the tiny change I made? I am fine with that, just wanted to be sure!
For now, we can already do the thing I proposed with the gp interface (and return the certificate instead of True
/False
). And we can also add a TODO
section in the documentation saying that we should move to the pari interface as soon as #17631 is ready.
Vincent
comment:11 in reply to: ↑ 10 ; followup: ↓ 12 Changed 6 years ago by
Replying to vdelecroix:
Replying to bruno:
For concreteness, what do you propose to do now? If I understand correctly, you propose to wait for #17631 and then make
is_cyclotomic()
directly call the pari functions, and not introduce the tiny change I made? I am fine with that, just wanted to be sure!For now, we can already do the thing I proposed with the gp interface (and return the certificate instead of
True
/False
). And we can also add aTODO
section in the documentation saying that we should move to the pari interface as soon as #17631 is ready.
I am not sure to agree since it actually makes the function slower:
sage: p = cyclotomic_polynomial(150) sage: %timeit p.is_cyclotomic() 100 loops, best of 3: 2.75 ms per loop sage: %timeit p._gp_().poliscyclo() 100 loops, best of 3: 3.96 ms per loop
As an aside, isn't it a problem for a method named is_cyclotomic()
to return an integer rather than a boolean? Of course, the certificate is interesting to have, and one can still test with is_cyclotomic()
as before. It is just that I expect a method is_...()
to return a boolean. Maybe I am wrong...
comment:12 in reply to: ↑ 11 ; followup: ↓ 13 Changed 6 years ago by
Replying to bruno:
Replying to vdelecroix:
Replying to bruno:
For concreteness, what do you propose to do now? If I understand correctly, you propose to wait for #17631 and then make
is_cyclotomic()
directly call the pari functions, and not introduce the tiny change I made? I am fine with that, just wanted to be sure!For now, we can already do the thing I proposed with the gp interface (and return the certificate instead of
True
/False
). And we can also add aTODO
section in the documentation saying that we should move to the pari interface as soon as #17631 is ready.I am not sure to agree since it actually makes the function slower:
sage: p = cyclotomic_polynomial(150) sage: %timeit p.is_cyclotomic() 100 loops, best of 3: 2.75 ms per loop sage: %timeit p._gp_().poliscyclo() 100 loops, best of 3: 3.96 ms per loop
Right. Most of the time seems to be lost in the conversion
sage: l = [cyclotomic_polynomial(i) for i in range(2,50)] sage: timeit("for p in l: _ = p.is_cyclotomic()") 5 loops, best of 3: 51.7 ms per loop sage: timeit("for p in l: _ = p._gp_().poliscyclo()") 5 loops, best of 3: 303 ms per loop sage: lgp = [p._gp_() for p in l] sage: timeit("for p in lgp: _ = p.poliscyclo()") 25 loops, best of 3: 34.8 ms per loop
As an aside, isn't it a problem for a method named
is_cyclotomic()
to return an integer rather than a boolean? Of course, the certificate is interesting to have, and one can still test withis_cyclotomic()
as before. It is just that I expect a methodis_...()
to return a boolean. Maybe I am wrong...
You are right. The best option would be:
def is_cyclotomic(self, certificate=False): ... return ans if certificate else bool(ans)
A compromise can be:
 add the keyword
certificate
(if it isTrue
we useself._gp_().polisyclco()
)  do your initial suggested change
 add a
TODO
to the documentation to switch from gp to pari (and open a ticket for it). We will then consider the fact of making it the default.
Vincent
comment:13 in reply to: ↑ 12 Changed 6 years ago by
Replying to vdelecroix:
A compromise can be:
 add the keyword
certificate
(if it isTrue
we useself._gp_().polisyclco()
) do your initial suggested change
 add a
TODO
to the documentation to switch from gp to pari (and open a ticket for it). We will then consider the fact of making it the default.
I like this proposition, I'll do it.
comment:14 Changed 6 years ago by
 Commit changed from f79eeaaa630b1bd966a25a2f7b0f7a5f1dda3378 to 053a30cc574fb122ec8527c53849419ccf871e68
Branch pushed to git repo; I updated commit sha1. New commits:
053a30c  Add certificate keyword

comment:15 Changed 6 years ago by
 Status changed from needs_review to needs_work
Just typographical comment. The TODO
block must be indented. In other words
.. TODO:: As soon as #XYZ is finished we should do blih and blah.
Have a look to Docstring section of the developer manual.
It is _gp_()
and not __gp__()
(i.e. one underscore).
I would move the ._gp_().poliscyclo()
below the tests for characteristic.
Vincent
comment:16 Changed 6 years ago by
I pushed too fast, and noticed my mistake on __gp__()
. Actually, there are other ones. I'm correcting this.
Thanks anyway for the feedback, and especially for the documentation!
comment:17 Changed 6 years ago by
 Commit changed from 053a30cc574fb122ec8527c53849419ccf871e68 to c7b5f6861c5bbabadd9a9e68cdc8c6d7a5673c50
Branch pushed to git repo; I updated commit sha1. New commits:
c7b5f68  Correcting mistakes + docstring formatting

comment:18 Changed 6 years ago by
 Status changed from needs_work to needs_review
comment:19 followup: ↓ 21 Changed 6 years ago by
One more thing for the doc: for trac tickets you should use :trac:`17730`
. That will create automatically a link in the compiled documentation.
More seriously, we definitely do not want to test irreducibility at the beginning. You can possibly add this test in between step 2 and 3. That would prevent some useless factorization. But if you do so, I want to see some benchmarks first.
Vincent
comment:20 Changed 6 years ago by
 Status changed from needs_review to needs_work
comment:21 in reply to: ↑ 19 ; followup: ↓ 22 Changed 6 years ago by
Replying to vdelecroix:
One more thing for the doc: for trac tickets you should use
:trac:`17730`
. That will create automatically a link in the compiled documentation.
OK.
More seriously, we definitely do not want to test irreducibility at the beginning. You can possibly add this test in between step 2 and 3. That would prevent some useless factorization. But if you do so, I want to see some benchmarks first.
I do not understand your remark: irreducibility was tested as the first step before, and moving this test later is precisely the purpose of my ticket.
comment:22 in reply to: ↑ 21 Changed 6 years ago by
Replying to bruno:
More seriously, we definitely do not want to test irreducibility at the beginning. You can possibly add this test in between step 2 and 3. That would prevent some useless factorization. But if you do so, I want to see some benchmarks first.
I do not understand your remark: irreducibility was tested as the first step before, and moving this test later is precisely the purpose of my ticket.
Right ;) I will try a bit more optimization on my side.
comment:23 Changed 6 years ago by
 Commit changed from c7b5f6861c5bbabadd9a9e68cdc8c6d7a5673c50 to 1b01fad201f72a254845e763cf8b072d28cde8b8
Branch pushed to git repo; I updated commit sha1. New commits:
1b01fad  improved formatting (trac ticket)

comment:24 Changed 6 years ago by
 Status changed from needs_work to needs_review
comment:25 followup: ↓ 27 Changed 6 years ago by
 Branch changed from u/bruno/slightly_fasten_is_cyclotomic__ to u/vdelecroix/17695
 Commit changed from 1b01fad201f72a254845e763cf8b072d28cde8b8 to 90ac73a342710459415a32442d418751e3560be6
Actually, the factor at the end is overkill as we only want a square root. I did have a look at the method .is_square()
and it was dirty enough to be slow. So I just make it simpler/clearer. The speedup is very small and concerns mostly positive cases.
I am happy with your commits. Just review mine.
Vincent
New commits:
90ac73a  trac #17695: simplify factor step + faster .sqrt()

comment:26 Changed 6 years ago by
 Branch changed from u/vdelecroix/17695 to u/bruno/17695
comment:27 in reply to: ↑ 25 Changed 6 years ago by
 Commit changed from 90ac73a342710459415a32442d418751e3560be6 to 85d8ffcdd69d911dcc042e4dd3ab8374f90e4d56
Replying to vdelecroix:
Actually, the factor at the end is overkill as we only want a square root. I did have a look at the method
.is_square()
and it was dirty enough to be slow. So I just make it simpler/clearer. The speedup is very small and concerns mostly positive cases.I am happy with your commits. Just review mine.
I am fine with your commit too, but for one typo that I've corrected in the docstring of is_square()
:
def is_square(self, root=False): """ Returns whether or not polynomial is square. If the optional  argument ``root`` is set to ``True, then also returns the square root + argument ``root`` is set to ``True``, then also returns the square root (or ``None``, if the polynomial is not square).
I do not know whether it is acceptable for me to set "positive review". I let you do so if you consider everything's OK. Thanks for the patient feedback!
New commits:
85d8ffc  Review commit

comment:28 Changed 6 years ago by
 Branch changed from u/bruno/17695 to public/ticket/17695
 Commit changed from 85d8ffcdd69d911dcc042e4dd3ab8374f90e4d56 to a2a6f026794fb12e1d557d8ca6738c45d5ed937f
 Status changed from needs_review to positive_review
Looks good to me. I have just made a few typographical changes, so I allow myself to put this into positive review.
New commits:
a2a6f02  trac #17695 some doc and code formatting details

comment:29 Changed 6 years ago by
Ho great! Nice catch!
Vincent
comment:30 Changed 6 years ago by
 Reviewers set to Frédéric Chapoton
comment:31 Changed 6 years ago by
 Branch changed from public/ticket/17695 to a2a6f026794fb12e1d557d8ca6738c45d5ed937f
 Resolution set to fixed
 Status changed from positive_review to closed
Hello Bruno,
Actually, I discussed this with Karim Belabas in Bordeaux few weeks ago. We should either call Pari or stick to their implementation. It is much faster:
Moreover they do have the function
poliscycloprod
that checks if all roots are roots of unity (which is even faster).Vincent