Opened 7 years ago

Last modified 3 years ago

#17895 closed defect

Computing all roots is faster than computing a single one — at Version 3

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: Reviewers:
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: 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

Change History (3)

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)
Note: See TracTickets for help on using tickets.