Module complex_roots

We are given a polynomial and its derivative (with complex
interval coefficients), an estimated root, and a complex interval
field to use in computations.  We use interval arithmetic to
refine the root and prove that we have in fact isolated a unique
root.

If we succeed, we return the isolated root; if we fail, we return
None.

EXAMPLES:
    sage: from sage.rings.polynomial.complex_roots import *
    sage: x = polygen(ZZ)
    sage: p = x^9 - 1
    sage: ip = CIF['x'](p); ip
    [1.0000000000000000 .. 1.0000000000000000]*x^9 + [-1.0000000000000000 .. -1.0000000000000000]
    sage: ipd = CIF['x'](p.derivative()); ipd
    [9.0000000000000000 .. 9.0000000000000000]*x^8
    sage: irt = CIF(CC(cos(2*pi/9), sin(2*pi/9))); irt
    [0.76604444311897801 .. 0.76604444311897802] + [0.64278760968653925 .. 0.64278760968653926]*I
    sage: ip(irt)
    [-1.3322676295501879e-15 .. 6.6613381477509393e-16] + [-1.2212453270876722e-15 .. 6.6613381477509393e-16]*I
    sage: ipd(irt)
    [6.8943999880707931 .. 6.8943999880708056] - [5.7850884871788474 .. 5.7850884871788600]*I
    sage: refine_root(ip, ipd, irt, CIF)
    [0.76604444311897779 .. 0.76604444311897824] + [0.64278760968653914 .. 0.64278760968653948]*I

We are given a squarefree polynomial p, a list of estimated roots,
and a precision.

We attempt to verify that the estimated roots are in fact distinct
roots of the polynomial, using interval arithmetic of precision prec.
If we succeed, we return a list of intervals bounding the roots; if we
fail, we return None.

EXAMPLES:
    sage: x = polygen(ZZ)
    sage: p = x^3 - 1
    sage: rts = [CC.zeta(3)^i for i in range(0, 3)]
    sage: from sage.rings.polynomial.complex_roots import interval_roots
    sage: interval_roots(p, rts, 53)
    [[1.0000000000000000 .. 1.0000000000000000], [-0.50000000000000012 .. -0.49999999999999983] + [0.86602540378443848 .. 0.86602540378443882]*I, [-0.50000000000000012 .. -0.49999999999999988] - [0.86602540378443848 .. 0.86602540378443882]*I]
    sage: interval_roots(p, rts, 200)
    [[1.0000000000000000000000000000000000000000000000000000000000000 .. 1.0000000000000000000000000000000000000000000000000000000000000], [-0.50000000000000000000000000000000000000000000000000000000000063 .. -0.49999999999999999999999999999999999999999999999999999999999968] + [0.86602540378443864676372317075293618347140262690519031402790313 .. 0.86602540378443864676372317075293618347140262690519031402790377]*I, [-0.50000000000000000000000000000000000000000000000000000000000063 .. -0.49999999999999999999999999999999999999999999999999999999999968] - [0.86602540378443864676372317075293618347140262690519031402790313 .. 0.86602540378443864676372317075293618347140262690519031402790377]*I]

Given a list of complex intervals, check whether they are pairwise
disjoint.

EXAMPLES:
    sage: from sage.rings.polynomial.complex_roots import intervals_disjoint
    sage: a = CIF(RIF(0, 3), 0)
    sage: b = CIF(0, RIF(1, 3))
    sage: c = CIF(RIF(1, 2), RIF(1, 2))
    sage: d = CIF(RIF(2, 3), RIF(2, 3))
    sage: intervals_disjoint([a,b,c,d])
    False
    sage: d2 = CIF(RIF(2, 3), RIF(2.001, 3))
    sage: intervals_disjoint([a,b,c,d2])
    True

Compute the complex roots of a given polynomial with exact
coefficients (integer, rational, Gaussian rational, and algebraic
coefficients are supported).  Returns a list of pairs of a root
and its multiplicity.

Roots are returned as a ComplexIntervalFieldElement; each interval
includes exactly one root, and the intervals are disjoint.

By default, the algorithm will do a squarefree decomposition
to get squarefree polynomials.  The skip_squarefree parameter
lets you skip this step.  (If this step is skipped, and the polynomial
has a repeated root, then the algorithm will loop forever!)

You can specify retval='interval' (the default) to get roots as
complex intervals.  The other options are retval='algebraic' to
get elements of QQbar, or retval='algebraic_real' to get only
the real roots, and to get them as elements of AA.

EXAMPLES:
    sage: from sage.rings.polynomial.complex_roots import complex_roots
    sage: x = polygen(ZZ)
    sage: complex_roots(x^5 - x - 1)
    [([1.1673039782614185 .. 1.16730397826141...], 1), ([0.18123244446987518 .. 0.18123244446987558] + [1.0839541013177103 .. 1.0839541013177110]*I, 1), ([0.181232444469875... .. 0.1812324444698755...] - [1.083954101317710... .. 1.0839541013177110]*I, 1), ([-0.76488443360058489 .. -0.76488443360058455] + [0.35247154603172609 .. 0.3524715460317264...]*I, 1), ([-0.76488443360058489 .. -0.76488443360058455] - [0.35247154603172609 .. 0.35247154603172643]*I, 1)]

    sage: K.<im> = NumberField(x^2 + 1)
    sage: eps = 1/2^100
    sage: x = polygen(K)
    sage: p = (x-1)*(x-1-eps)*(x-1+eps)*(x-1-eps*im)*(x-1+eps*im)

This polynomial actually has all-real coefficients, and is very, very
close to (x-1)^5:
    sage: [RR(QQ(a)) for a in list(p - (x-1)^5)]
    [3.87259191484932e-121, -3.87259191484932e-121]
    sage: rts = complex_roots(p)
    sage: [ComplexIntervalField(10)(rt[0] - 1) for rt in rts]
    [0, [7.8886e-31 .. 7.8887e-31]*I, [7.8886e-31 .. 7.8887e-31], [-7.8887e-31 .. -7.8886e-31]*I, [-7.8887e-31 .. -7.8886e-31]]

We can get roots either as intervals, or as elements of QQbar or AA.
    sage: p = (x^2 + x - 1)
    sage: p = p * p(x*im)
    sage: p
    (-1)*x^4 + (im - 1)*x^3 + im*x^2 + (-im - 1)*x + 1

Two of the roots have a zero real component; two have a zero
imaginary component.  These zero components will be found slightly
inaccurately, and the exact values returned are very sensitive to
the (non-portable) results of numpy.  So we post-process the roots
for printing, to get predictable doctest results.
    sage: def tiny(x):
    ...       return x.contains_zero() and x.absolute_diameter() <  1e-14
    sage: def smash(x):
    ...       x = CIF(x[0]) # discard multiplicity
    ...       if tiny(x.imag()): return x.real()
    ...       if tiny(x.real()): return CIF(0, x.imag())
    sage: rts = complex_roots(p); type(rts[0][0]), sorted(map(smash, rts))
    (<type 'sage.rings.complex_interval.ComplexIntervalFieldElement'>, [[-1.6180339887498952 .. -1.618033988749894...], [-0.61803398874989502 .. -0.61803398874989468]*I, [1.618033988749894... .. 1.6180339887498952]*I, [0.61803398874989468 .. 0.61803398874989502]])
    sage: rts = complex_roots(p, retval='algebraic'); type(rts[0][0]), sorted(map(smash, rts))
    (<class 'sage.rings.qqbar.AlgebraicNumber'>, [[-1.6180339887498950 .. -1.6180339887498946], [-0.61803398874989491 .. -0.61803398874989479]*I, [1.6180339887498946 .. 1.6180339887498950]*I, [0.61803398874989479 .. 0.61803398874989491]])
    sage: rts = complex_roots(p, retval='algebraic_real'); type(rts[0][0]), rts
    (<class 'sage.rings.qqbar.AlgebraicReal'>, [([-1.6180339887498950 .. -1.6180339887498946], 1), ([0.61803398874989479 .. 0.61803398874989491], 1)])