Context Navigation

← Previous Ticket
Next Ticket →

#17895 closed defect (fixed)

Computing all roots is faster than computing a single one

Reported by:	gagern	Owned by:
Priority:	major	Milestone:	sage-8.9
Component:	number fields	Keywords:	qqbar polynomial_root roots performance
Cc:	sbrandhorst	Merged in:
Authors:	Vincent Delecroix	Reviewers:	Simon Brandhorst
Report Upstream:	N/A	Work issues:
Branch:	d5105b6 (Commits, GitHub, GitLab)	Commit:	d5105b625060d4b2a9cfe9cee4316d8b7c96ea46
Dependencies:	#28191	Stopgaps:

Status badges

Description (last modified by vdelecroix)

The following example comes from comment:27:ticket:16964 via comment:2:ticket:17886.

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: ival = CIF((0.77, 0.78), (-0.08, -0.07))
sage: %time z1, = [r for r in p4.roots(QQbar, False) if r in ival]
CPU times: user 1.43 s, sys: 195 ms, total: 1.62 s
Wall time: 1.47 s
sage: %time z2 = QQbar.polynomial_root(p4, ival)
CPU times: user 1min 5s, sys: 212 ms, total: 1min 5s
Wall time: 1min 5s

The computation for z1 works reasonably well and completes in under 2 seconds, but the one for z2 takes over a minute. Which is definitely wrong, since computing all roots and then choosing the right one should be more work, not less than just computing a single one!

The reason for this is the time spent in the square-free decomposition (called from sage/rings/polynomial/complex_roots.py) which behave differently depending whether the polynomial is defined with coefficients in QQ or AA

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 somehow independent problem about gcd is tracked at #28199).

With the branch attached the timings become

sage: sage: %time z1, = [r for r in p4.roots(QQbar, False) if r in ival]
CPU times: user 411 ms, sys: 152 µs, total: 411 ms
Wall time: 410 ms
sage: sage: %time z2 = QQbar.polynomial_root(p4, ival)
CPU times: user 294 ms, sys: 1 µs, total: 294 ms
Wall time: 296 ms

Change History (20)

comment:1 Changed 7 years ago by pernici

If the interval is close enough to the root, QQbar.polynomial_root is fast

sage: z1
0.7745442035157198? - 0.07841573936877289?*I
sage: ival = CIF((0.77454420351571, 0.77454420351572), (-0.078415739368773, -0.078415739368772))
sage: %time z2 = QQbar.polynomial_root(p4, ival); z2._more_precision()
Wall time: 26.7 ms
sage: z2
0.7745442035157198? - 0.07841573936877289?*I

Last edited 3 years ago by vdelecroix (previous) (diff)

comment:2 Changed 3 years ago by vdelecroix

For the record, here is what profiling shows

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
     6320   23.769    0.004   23.769    0.004 {built-in method _operator.sub}
     6481   17.266    0.003   17.266    0.003 {built-in method _operator.mul}
        2    0.225    0.112    0.231    0.116 {method 'roots' of 'sage.rings.polynomial.polynomial_element.Polynomial' objects}
    16711    0.188    0.000    0.188    0.000 sets_cat.py:974(_element_constructor_from_element_class)
       89    0.069    0.001    0.070    0.001 {sage.rings.polynomial.refine_root.refine_root}
    15053    0.037    0.000    0.272    0.000 qqbar.py:3074(__init__)
    15287    0.031    0.000    0.227    0.000 qqbar.py:5723(_interval_fast)
        1    0.028    0.028   41.429   41.429 fields.py:218(_gcd_univariate_polynomial)
    14889    0.024    0.000   41.096    0.003 qqbar.py:7736(an_binop_rational)
    15052    0.024    0.000    0.036    0.000 qqbar.py:5653(__init__)
    47328    0.022    0.000    0.022    0.000 {built-in method builtins.isinstance}
    15052    0.017    0.000    0.297    0.000 qqbar.py:4712(__init__)
     6481    0.014    0.000   17.451    0.003 qqbar.py:3240(_mul_)
     6320    0.012    0.000   23.946    0.004 qqbar.py:3284(_sub_)
    15449    0.009    0.000    0.009    0.000 {built-in method sage.rings.real_mpfi.RealIntervalField}
    15287    0.006    0.000    0.008    0.000 qqbar.py:4730(_ensure_real)
     2088    0.003    0.000    0.022    0.000 qqbar.py:3271(_add_)
    15450    0.003    0.000    0.003    0.000 {sage.rings.complex_interval.is_ComplexIntervalFieldElement}
  448/160    0.002    0.000    0.004    0.000 complex_field.py:389(_element_constructor_)
        4    0.002    0.000    0.024    0.006 {method 'derivative' of 'sage.rings.polynomial.polynomial_element.Polynomial' objects}
     2088    0.002    0.000    0.002    0.000 {built-in method _operator.add}
      484    0.001    0.000    0.005    0.000 qqbar.py:4023(interval)
1322/1000    0.001    0.000    0.005    0.000 complex_interval_field.py:427(__call__)
      484    0.001    0.000    0.003    0.000 qqbar.py:3999(interval_diameter)
      216    0.001    0.000    0.003    0.000 {method 'n' of 'sage.structure.element.Element' objects}
        1    0.001    0.001   41.766   41.766 qqbar.py:6487(_complex_refine_interval)
        4    0.001    0.000    0.005    0.001 {built-in method builtins.sorted}
       81    0.001    0.000    0.001    0.000 qqbar.py:5809(invert)
  628/324    0.001    0.000    0.005    0.000 complex_field.py:352(__call__)
      561    0.001    0.000    0.001    0.000 {built-in method builtins.hasattr}
      240    0.001    0.000    0.004    0.000 misc.py:5558(_key_complex_for_display)
      144    0.001    0.000    0.001    0.000 {method '__copy__' of 'sage.categories.map.Map' objects}
      719    0.000    0.000    0.000    0.000 {method 'diameter' of 'sage.rings.real_mpfi.RealIntervalFieldElement' objects}
       15    0.000    0.000    0.007    0.000 polynomial_element_generic.py:1052(__init__)
  161/155    0.000    0.000    0.001    0.000 homset.py:84(Hom)
      312    0.000    0.000    0.001    0.000 copy.py:66(copy)
        2    0.000    0.000    0.086    0.043 complex_roots.py:47(interval_roots)
      235    0.000    0.000    0.001    0.000 qqbar.py:3905(_more_precision)
        1    0.000    0.000   41.757   41.757 complex_roots.py:147(complex_roots)
      162    0.000    0.000    0.001    0.000 qqbar.py:5244(real_number)
      161    0.000    0.000    0.001    0.000 qqbar.py:3970(interval_fast)
     1067    0.000    0.000    0.000    0.000 complex_field.py:289(_real_field)
       11    0.000    0.000    0.000    0.000 homset.py:582(__init__)
       81    0.000    0.000    0.004    0.000 qqbar.py:3262(__invert__)
      235    0.000    0.000    0.001    0.000 qqbar.py:4751(_more_precision)
        2    0.000    0.000    0.000    0.000 polynomial_ring.py:234(__init__)
      469    0.000    0.000    0.001    0.000 <frozen importlib._bootstrap>:1009(_handle_fromlist)
       85    0.000    0.000    0.000    0.000 qqbar.py:3386(__bool__)
        2    0.000    0.000    0.000    0.000 complex_field.py:193(__init__)
       82    0.000    0.000    0.001    0.000 qqbar.py:720(_element_constructor_)
      671    0.000    0.000    0.000    0.000 complex_interval_field.py:266(prec)
      160    0.000    0.000    0.000    0.000 {method 'real' of 'sage.rings.complex_number.ComplexNumber' objects}
        1    0.000    0.000   41.766   41.766 qqbar.py:6320(refine_interval)
      218    0.000    0.000    0.001    0.000 complex_field.py:91(ComplexField)
       80    0.000    0.000    0.003    0.000 fields.py:767(inverse_of_unit)
        1    0.000    0.000    0.001    0.001 {method 'change_ring' of 'sage.rings.polynomial.polynomial_element.Polynomial' objects}
      216    0.000    0.000    0.000    0.000 {method 'abs' of 'sage.structure.element.RingElement' objects}
        6    0.000    0.000    0.001    0.000 polynomial_ring_constructor.py:52(PolynomialRing)
        1    0.000    0.000   41.436   41.436 fields.py:441(_squarefree_decomposition_univariate_polynomial)
      160    0.000    0.000    0.000    0.000 {method 'imag' of 'sage.rings.complex_number.ComplexNumber' objects}
        1    0.000    0.000    0.001    0.001 complex_roots.py:85(intervals_disjoint)
       17    0.000    0.000    0.000    0.000 dynamic_class.py:127(dynamic_class)
      160    0.000    0.000    0.000    0.000 {method 'real' of 'sage.rings.complex_interval.ComplexIntervalFieldElement' objects}
        1    0.000    0.000    0.000    0.000 {method 'monic' of 'sage.rings.polynomial.polynomial_element.Polynomial' objects}
        9    0.000    0.000    0.008    0.001 polynomial_ring.py:320(_element_constructor_)
       80    0.000    0.000    0.000    0.000 complex_roots.py:118(row_disjoint)
        6    0.000    0.000    0.000    0.000 polynomial_ring_constructor.py:677(_single_variate)
        1    0.000    0.000    0.000    0.000 dynamic_class.py:346(dynamic_class_internal)
      6/5    0.000    0.000    0.000    0.000 {method 'has_coerce_map_from' of 'sage.structure.parent.Parent' objects}
      496    0.000    0.000    0.000    0.000 complex_interval_field.py:255(is_exact)
        2    0.000    0.000    0.000    0.000 category.py:499(__repr_object_names)
       91    0.000    0.000    0.000    0.000 weakref.py:341(__init__)
        6    0.000    0.000    0.000    0.000 homset.py:40(RingHomset)
      468    0.000    0.000    0.000    0.000 {method 'precision' of 'sage.rings.real_mpfi.RealIntervalFieldElement' objects}
        8    0.000    0.000    0.000    0.000 rings.py:322(_Hom_)
        1    0.000    0.000   41.757   41.757 qqbar.py:6060(complex_roots)
      160    0.000    0.000    0.000    0.000 {method 'imag' of 'sage.rings.complex_interval.ComplexIntervalFieldElement' objects}
       36    0.000    0.000    0.000    0.000 complex_roots.py:112(column_disjoint)
       85    0.000    0.000    0.000    0.000 {method 'contains_zero' of 'sage.rings.real_mpfi.RealIntervalFieldElement' objects}
        3    0.000    0.000    0.000    0.000 {method '_coerce_map_via' of 'sage.structure.parent.Parent' objects}
       97    0.000    0.000    0.000    0.000 weakref.py:336(__new__)
        5    0.000    0.000    0.001    0.000 rings.py:845(__getitem__)
       80    0.000    0.000    0.001    0.000 misc.py:5634(<lambda>)
      152    0.000    0.000    0.000    0.000 {built-in method builtins.issubclass}
        8    0.000    0.000    0.000    0.000 pushout.py:817(__init__)
       97    0.000    0.000    0.000    0.000 {built-in method __new__ of type object at 0x7fc5f0f21a00}
        6    0.000    0.000    0.000    0.000 {built-in method sage.structure.category_object.normalize_names}
      333    0.000    0.000    0.000    0.000 {method 'append' of 'list' objects}
      155    0.000    0.000    0.000    0.000 homset.py:1174(domain)
      318    0.000    0.000    0.000    0.000 {method 'get' of 'dict' objects}
        5    0.000    0.000    0.000    0.000 polynomial_ring.py:661(_coerce_map_from_)
      306    0.000    0.000    0.000    0.000 complex_field.py:43(late_import)
        1    0.000    0.000   41.768   41.768 {built-in method builtins.exec}
        1    0.000    0.000    0.000    0.000 qqbar.py:6545(<listcomp>)
      144    0.000    0.000    0.000    0.000 {method 'precision' of 'sage.rings.real_mpfr.RealNumber' objects}
        5    0.000    0.000    0.000    0.000 {method 'coerce_map_from' of 'sage.structure.parent.Parent' objects}
        2    0.000    0.000    0.000    0.000 polynomial_singular_interface.py:355(can_convert_to_singular)
        4    0.000    0.000    0.000    0.000 {method '_populate_coercion_lists_' of 'sage.structure.parent.Parent' objects}
      162    0.000    0.000    0.000    0.000 {method 'precision' of 'sage.rings.real_mpfr.RealField_class' objects}
        3    0.000    0.000    0.000    0.000 complex_field.py:431(_coerce_map_from_)
        1    0.000    0.000    0.000    0.000 sequence.py:80(Sequence)
        8    0.000    0.000    0.000    0.000 polynomial_ring.py:628(construction)
        2    0.000    0.000    0.000    0.000 polynomial_ring.py:1650(__init__)
        3    0.000    0.000    0.000    0.000 rational_field.py:349(_coerce_map_from_)
       16    0.000    0.000    0.000    0.000 <frozen importlib._bootstrap>:416(parent)
        1    0.000    0.000   41.768   41.768 <string>:1(<module>)
        1    0.000    0.000    0.000    0.000 qqbar.py:6533(<listcomp>)
      146    0.000    0.000    0.000    0.000 {built-in method builtins.getattr}
      200    0.000    0.000    0.000    0.000 {built-in method builtins.len}
        2    0.000    0.000    0.000    0.000 polynomial_ring.py:1826(__init__)
        4    0.000    0.000    0.000    0.000 polynomial_ring_constructor.py:672(_save_in_cache)
      155    0.000    0.000    0.000    0.000 homset.py:1189(codomain)
      149    0.000    0.000    0.000    0.000 copy.py:111(_copy_immutable)
        2    0.000    0.000    0.000    0.000 category.py:525(_repr_object_names)
        1    0.000    0.000   41.757   41.757 qqbar.py:6060(complex_roots)
        1    0.000    0.000   41.768   41.768 qqbar.py:6171(__init__)
        1    0.000    0.000   41.436   41.436 {method 'squarefree_decomposition' of 'sage.rings.polynomial.polynomial_element.Polynomial' objects}
        1    0.000    0.000   41.768   41.768 qqbar.py:1361(polynomial_root)
        1    0.000    0.000   41.757   41.757 qqbar.py:6637(_complex_isolate_interval)

comment:3 Changed 3 years ago by vdelecroix

Description modified (diff)

Last edited 3 years ago by vdelecroix (previous) (diff)

comment:4 Changed 3 years ago by vdelecroix

Authors set to Vincent Delecroix
Branch set to u/vdelecroix/17895
Commit set to 2984ecb978aca6e1f3d3926f3e4e9afbe954e3f7
Milestone changed from sage-6.6 to sage-8.9
Status changed from new to needs_review

With the change in the branch, timings go back to normal! Though, I have some strange doctest failures in qqbar.py. Before trying to understand them I would like to see if I did break anything else. So I am just waiting for a more complete patchbot report.

New commits:

3893ba3	`rewrite root refinement in qqbar`
2984ecb	`allow rational polys in AlgebraicPolynomialTracker`

comment:5 Changed 3 years ago by vdelecroix

Dependencies set to #28191

comment:6 Changed 3 years ago by git

Commit changed from 2984ecb978aca6e1f3d3926f3e4e9afbe954e3f7 to 659b7a6ce931a998df29fca979772f604a78f677

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

242a41a	`allow rational polys in AlgebraicPolynomialTracker`
659b7a6	`fix doctests in qqbar.py`

comment:7 Changed 3 years ago by vdelecroix

Description modified (diff)

comment:8 Changed 3 years ago by vdelecroix

Cc sbrandhorst added

comment:9 Changed 3 years ago by vdelecroix

Description modified (diff)

comment:10 Changed 3 years ago by git

Commit changed from 659b7a6ce931a998df29fca979772f604a78f677 to 49f917a3b74b3784e5ced901404b3778c94ecf92

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

49f917a also handle ZZ

comment:11 Changed 3 years ago by vdelecroix

ready for review :-)

comment:12 follow-up: ↓ 14 Changed 3 years ago by sbrandhorst

seems to look good. One question:

Is there a reason to hardcode QQy in several places?

comment:13 Changed 3 years ago by sbrandhorst

Reviewers set to Simon Brandhorst
Status changed from needs_review to needs_info

comment:14 in reply to: ↑ 12 Changed 3 years ago by vdelecroix

Replying to sbrandhorst:

seems to look good. One question:

Is there a reason to hardcode QQy in several places?

Mainly for consistency and the fact that you don't want to multiply parents. This is an internal class, so the name of the variable does not matter from the user perspective. Before, it used to be canoncalized with a y (after the call to complex roots). I kept this convention.

comment:15 Changed 3 years ago by sbrandhorst

You could add a doctest to show the speedup. Then I am happy.

comment:16 Changed 3 years ago by git

Commit changed from 49f917a3b74b3784e5ced901404b3778c94ecf92 to d5105b625060d4b2a9cfe9cee4316d8b7c96ea46

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

d5105b6 add a doctest for speedup

comment:17 Changed 3 years ago by vdelecroix

Status changed from needs_info to needs_review

comment:18 Changed 3 years ago by sbrandhorst

Status changed from needs_review to positive_review

comment:19 Changed 3 years ago by vdelecroix

Thank you.

comment:20 Changed 3 years ago by vbraun

Branch changed from u/vdelecroix/17895 to d5105b625060d4b2a9cfe9cee4316d8b7c96ea46
Resolution set to fixed
Status changed from positive_review to closed

Note: See TracTickets for help on using tickets.

Download in other formats: