Opened 3 years ago

Last modified 3 months ago

#28199 new defect

gcd can be very slow in AA[x]

Reported by: vdelecroix Owned by:
Priority: major Milestone: sage-9.7
Component: algebra Keywords:
Cc: bruno, mmezzarobba Merged in:
Authors: Reviewers:
Report Upstream: N/A Work issues:
Branch: u/bruno/gcd_via_UFD (Commits, GitHub, GitLab) Commit: 62db08ca0b85c535f3d6ef5bacac19ecefa8908b
Dependencies: Stopgaps:

Status badges

Description

Here is a x50 slower example

sage: x,y = polygens(QQ,"x,y")
sage: p1 = x^5 + 6*x^4 - 42*x^3 - 142*x^2 + 467*x + 422
sage: p2 = p1(x=(x-1)^2)
sage: p3 = p2(x=x*y).resultant(p2,x).univariate_polynomial()
sage: p4, = [f[0] for f in p3.factor() if f[0].degree() == 80]
sage: %time _ = p4.squarefree_decomposition()
CPU times: user 807 µs, sys: 0 ns, total: 807 µs
Wall time: 883 µs
sage: %time _ = p4.change_ring(AA).squarefree_decomposition()
CPU times: user 40.1 s, sys: 3.21 ms, total: 40.1 s
Wall time: 40.1 s

This problem originally appeared in #17895 (where a better workaround has been found).

Change History (14)

comment:1 Changed 3 years ago by bruno

  • Branch set to u/bruno/gcd_via_UFD

comment:2 Changed 3 years ago by git

  • Commit set to 62db08ca0b85c535f3d6ef5bacac19ecefa8908b

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

62db08c28199: Try to use subresultants algorithm

comment:3 Changed 3 years ago by bruno

I've made a very quick try: Instead of using the plain old Euclidean algorithm, one can use the so-called subresultant algorithm. This algorithm is already provided for UFDs so I tried to branch this algorithms for fields also. The resulting timings are still bad, but slightly less so:

sage: x,y = polygens(QQ,"x,y")
sage: p1 = x^5 + 6*x^4 - 42*x^3 - 142*x^2 + 467*x + 422
sage: p2 = p1(x=(x-1)^2)
sage: p3 = p2(x=x*y).resultant(p2,x).univariate_polynomial()
sage: p4, = [f[0] for f in p3.factor() if f[0].degree() == 80]
sage: %time _ = p4.squarefree_decomposition()
CPU times: user 1.89 ms, sys: 0 ns, total: 1.89 ms
Wall time: 1.91 ms
sage: %time _ = p4.change_ring(AA).squarefree_decomposition()
CPU times: user 3.14 s, sys: 15.3 ms, total: 3.16 s
Wall time: 3.16 s

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

You were quick :-) This is great: 10x improvement! There is no hope to beat the C-implementation from flint anyway.

Do you know why the subresultant is not the default method over fields? Could it be always better?

This reminds me the problem for computing determinant of matrices with rational coefficients. Going via Gauss elimination (which is a sort of "default" for fields) will create a terrible coefficient blowup.

comment:5 in reply to: ↑ 4 Changed 3 years ago by bruno

Replying to vdelecroix:

You were quick :-) This is great: 10x improvement! There is no hope to beat the C-implementation from flint anyway.

Of course, but I would hope for less than a second...

Do you know why the subresultant is not the default method over fields? Could it be always better?

That's why my first try. The problem is that it requires a _floordiv_ and not all fields have it. For instance:

sage: R.<x> = RR[]
sage: p = (x^2-1)^2
sage: p.is_squarefree()
Traceback (most recent call last):
...
TypeError: unsupported operand parent(s) for //: 'Real Field with 53 bits of precision' and 'Real Field with 53 bits of precision'

or

sage: S.<y> = GF(5)[]
sage: Q = S.quotient_ring(y^2+y+1)
sage: T.<z> = Q[]
sage: gcd(y*z+1, z + y)
Traceback (most recent call last):
...
TypeError: unsupported operand parent(s) for //: 'Univariate Quotient Polynomial Ring in ybar over Finite Field of size 5 with modulus y^2 + y + 1' and 'Univariate Quotient Polynomial Ring in ybar over Finite Field of size 5 with modulus y^2 + y + 1'

One possibility to check is whether adding a _floordiv_ method that simply calls self._div_(other) in the Fields category is sufficient and yields better timings.

comment:6 follow-up: Changed 3 years ago by bruno

Actually I checked: It is not faster in general (you can for instance test in the ring T of my example of the previous comment). This is reasonable: The only cases where it is faster is when there is coefficient growth. So there is a difficulty to decide when to use each algorithm. Maybe I should add a test as: if not self.base_ring().is_finite(). Note that in that case, RR would use subresultants while it is slightly slower. (But not that computing (non-approximate) gcds for polynomials over RR is mainly nonsense because of approximation issues.)

comment:7 in reply to: ↑ 6 Changed 3 years ago by vdelecroix

I think that exact algorithm should raise an error with non-approximate base ring.

comment:8 Changed 2 years ago by embray

  • Milestone changed from sage-8.9 to sage-9.1

Ticket retargeted after milestone closed

comment:9 Changed 2 years ago by mkoeppe

  • Milestone changed from sage-9.1 to sage-9.2

Batch modifying tickets that will likely not be ready for 9.1, based on a review of the ticket title, branch/review status, and last modification date.

comment:10 Changed 22 months ago by mkoeppe

  • Milestone changed from sage-9.2 to sage-9.3

comment:11 Changed 17 months ago by mkoeppe

  • Milestone changed from sage-9.3 to sage-9.4

Setting new milestone based on a cursory review of ticket status, priority, and last modification date.

comment:12 Changed 11 months ago by mkoeppe

  • Milestone changed from sage-9.4 to sage-9.5

comment:13 Changed 6 months ago by mkoeppe

  • Milestone changed from sage-9.5 to sage-9.6

comment:14 Changed 3 months ago by mkoeppe

  • Milestone changed from sage-9.6 to sage-9.7
Note: See TracTickets for help on using tickets.