arb_fmpz_poly_complex_roots: explicitly guarantee that roots are isolated (was true in practice, but could theoretically fail when the deflation hack is used)

This commit is contained in:
fredrik 2021-05-19 12:17:32 +02:00
parent 9c240a1e69
commit 2a1172bf3e
3 changed files with 56 additions and 4 deletions

View file

@ -12,7 +12,7 @@
#include "arb_fmpz_poly.h"
#include "flint/profiler.h"
int check_accuracy(acb_ptr vec, slong len, slong prec)
static int check_accuracy(acb_ptr vec, slong len, slong prec)
{
slong i;
@ -25,6 +25,28 @@ int check_accuracy(acb_ptr vec, slong len, slong prec)
return 1;
}
static int check_isolation(acb_srcptr roots, slong len)
{
slong i, j;
for (i = 0; i < len; i++)
{
if (arf_sgn(arb_midref(acb_imagref(roots + i))) >= 0)
{
for (j = i + 1; j < len; j++)
{
if (arf_sgn(arb_midref(acb_imagref(roots + j))) >= 0)
{
if (acb_overlaps(roots + i, roots + j))
return 0;
}
}
}
}
return 1;
}
void
arb_fmpz_poly_complex_roots(acb_ptr roots, const fmpz_poly_t poly, int flags, slong target_prec)
{
@ -149,6 +171,15 @@ arb_fmpz_poly_complex_roots(acb_ptr roots, const fmpz_poly_t poly, int flags, sl
arb_zero(acb_imagref(roots + i));
}
if (!check_isolation(roots, deg))
{
/* extremely unlikely */
if (flags & ARB_FMPZ_POLY_ROOTS_VERBOSE)
flint_printf("isolation failure!\n");
continue;
}
if (flags & ARB_FMPZ_POLY_ROOTS_VERBOSE)
flint_printf("done!\n");

View file

@ -20,7 +20,7 @@ check_roots(const fmpz_poly_t poly, acb_srcptr roots, slong prec)
arb_poly_t rpoly;
arb_t lead;
slong i, num_real, num_upper, deg;
slong i, j, num_real, num_upper, deg;
deg = fmpz_poly_degree(poly);
@ -64,6 +64,23 @@ check_roots(const fmpz_poly_t poly, acb_srcptr roots, slong prec)
flint_abort();
}
for (i = 0; i < deg; i++)
{
for (j = i + 1; j < deg; j++)
{
if (acb_overlaps(roots + i, roots + j))
{
flint_printf("FAIL! (isolation)\n");
flint_printf("deg = %wd, num_real = %wd, num_upper = %wd\n\n", deg, num_real, num_upper);
for (i = 0; i < deg; i++)
{
acb_printn(roots + i, 30, 0);
flint_printf("\n");
}
}
}
}
_arb_vec_clear(real, num_real);
_acb_vec_clear(upper, num_upper);
arb_poly_clear(rpoly);

View file

@ -66,7 +66,10 @@ Polynomial roots
.. function:: void arb_fmpz_poly_complex_roots(acb_ptr roots, const fmpz_poly_t poly, int flags, slong prec)
Writes to *roots* all the real and complex roots of the polynomial *poly*,
computed to *prec* accurate bits.
computed to at least *prec* accurate bits.
The root enclosures are guaranteed to be disjoint, so that
all roots are isolated.
The real roots are written first in ascending order (with
the imaginary parts set exactly to zero). The following
nonreal roots are written in arbitrary order, but with conjugate pairs
@ -94,7 +97,8 @@ Polynomial roots
All roots are refined to a relative accuracy of at least *prec* bits.
The output values will generally have higher actual precision,
depending on the precision used internally by the algorithm.
depending on the precision needed for isolation and the
precision used internally by the algorithm.
This implementation should be adequate for general use, but it is not
currently competitive with state-of-the-art isolation