## `numeric::polyroots` -- numerical roots of a univariate polynomial

#### Introduction

`numeric::polyroots``(p, ..)` returns numerical approximations of all real and complex roots of the univariate polynomial `p`.

#### Call(s)

`numeric::polyroots(p <, FixedPrecision> <, SquareFree> <, Factor> )`

#### Parameters

`p` | - | a univariate polynomial expression or a univariate polynomial of domain type `DOM_POLY`.

#### Options

`FixedPrecision` | - | launches a quick numerical search with fixed precision
`SquareFree` | - | a squarefree factorization of `p` is computed before the numerical search starts
`Factor` | - | a factorization of `p` is computed before the numerical search starts

#### Returns

a list of numerical roots.

#### Side Effects

The function is sensitive to the environment variable `DIGITS`, which determines the numerical working precision.

#### Related Functions

`RootOf`, `numeric::fsolve`, `numeric::polysysroots`, `numeric::realroot`, `numeric::realroots`, `polylib::realroots`, `solve`

#### Details

- The coefficients may be real or complex numbers. Also symbolic coefficients are accepted if they can be converted to floats.
- The trivial polynomial p=0 results in an error message. The empty list is returned for constant polynomials p<>0.
- Multiple roots are listed according to their multiplicities, i.e., the length of the root list coincides with the degree of p.
- The root list is sorted by `numeric::sort`.
- Up to roundoff effects, the numerical roots should be accurate to `DIGITS` significant digits, unless the option FixedPrecision is used.
- All floating point entries in `p` are internally approximated by rational numbers: `numeric::polyroots``(p)` computes the roots of `numeric::rationalize``(p, Minimize)`.
- For polynomial expressions in factored form, the numerical search is applied to each factor separately.
- With `setuserinfo(numeric::polyroots, 2)`, detailed information on the internal search is provided.
- It is recommended to use `polylib::realroots` if `p` is a real polynomial and only real roots are of interest.

#### Option: FixedPrecision

- This option provides the fastest way to obtain approximations of the roots by a numerical search with a fixed internal precision of 2 `DIGITS` decimal places.
- Note that badly isolated roots or multiple roots will usually not be approximated to `DIGITS` decimals when using this option. The problem of finding such roots is numerically ill-conditioned, i.e., such roots cannot be found to full precision with fixed precision arithmetic. Typically, a q-fold root will be approximated only to about 2*DIGITS/q decimal places. Cf. example 3.
- Without this option, `numeric::polyroots` internally increases the working precision until all roots are found to `DIGITS` decimal places.

#### Option: SquareFree

- With this option, a symbolic square free factorization is computed via `polylib::sqrfree(p)`. The numerical root finding algorithm is then applied to each square free factor.
- This option is recommended, when `p` is known to have multiple roots. Such roots force `numeric::polyroots` to increase the working precision internally increasing the costs of the numerical search. A square free factorization reduces the multiplicity of each root to one, speeding up the final numerical search. Cf. example 3.
- For polynomials with real rational coefficients, a square free factorization is always used, i.e., this option does not have any effect for such polynomials. For all other types of coefficients, a square free factorization may be costly and must be requested by this option.
- Multiple roots of `p` can be successfully dealt with by this option. However, for badly separated distinct roots the square free factorization will not improve the performance of the numerical search.

#### Option: Factor

- With this option, a symbolic factorization of `p` via `factor` is computed. The numerical root finding algorithm is then applied to each factor.
- This option is useful, when `p` can be successfully factorized (e.g., when `p` has multiple roots). The numerical search on the factors is much more efficient than the search on the original polynomial. On the other hand, symbolic factorization of `p` may be costly.

#### Example 1

Both polynomial expressions as well as `DOM_POLY` objects may be used to specify the polynomial:

    >> numeric::polyroots(x^3 - 3*x - sqrt(2))

                    [-1.414213562, -0.5176380902, 1.931851653]

    >> numeric::polyroots(PI*z^4 + I*z + 0.1)

          [- 0.5936837297 - 0.3729248918 I, 0.1003181767 I,

             0.6455316068 I, 0.5936837297 - 0.3729248918 I]

    >> numeric::polyroots(poly(x^5 - x^2, [x]))

          [- 0.5 - 0.8660254038 I, - 0.5 + 0.8660254038 I, 0.0, 0.0, 1.0]

#### Example 2

The following polynomial has exact coefficients:

    >> p := poly((x - 1)*(x - PI)^3, [x]):

    >> numeric::polyroots(p)

                   [1.0, 3.141592654, 3.141592654, 3.141592654]

Note that roundoff errors in the coefficients of `p` have a dramatic effect on multiple roots:

    >> p := poly((x - 1.0)*(x - float(PI))^3, [x]):

    >> numeric::polyroots(p)

          [0.9999999995, 3.140177788 - 0.00244957836 I,

             3.140177788 + 0.00244957836 I, 3.144422386]

These are the roots of the rationalized polynomial

    >> numeric::rationalize(p, Minimize)

                4                 3                2
          poly(x  - 461286/44249 x  + 176627/4525 x  - 2515405/41498 x +

             1837649/59267, [x])

    >> delete p:

#### Example 3

The multiple root `I/3` of the following polynomial can only be computed with restricted precision by fixed precision arithmetic:

    >> p := poly((x^2 - 6*x +8)*(x - I/3)^5, [x]):

    >> numeric::polyroots(p, FixedPrecision)

          [- 0.006109919675 + 0.3322082825 I,

             - 0.002938026857 + 0.3387257645 I,

             - 0.0007618302219 + 0.3271290976 I,

             0.004162114315 + 0.3378408591 I,

             0.005647662459 + 0.330762663 I, 2.0, 4.0]

Without the option FixedPrecision, the working precision is increased internally to compute better approximations:

    >> setuserinfo(numeric::polyroots, 1): numeric::polyroots(p)

          Info: computing roots of factor poly(x^7 - .. , [x])
          Info: increasing working precision to DIGITS=20
          Info: increasing working precision to DIGITS=40
          Info: increasing working precision to DIGITS=80
          Info: increasing working precision to DIGITS=160
          Info: accepting last approximation

          [0.3333333333 I, 0.3333333333 I, 0.3333333333 I,

             0.3333333333 I, 0.3333333333 I, 2.0, 4.0]

A square free factorization reduces the multiplicity of the root `I/3`:

    >> numeric::polyroots(p, SquareFree)

          Info: starting squarefree factorization
          Info: computing roots of factor poly(3*x - I, [x])
          Info: increasing working precision to DIGITS=20
          Info: computing roots of factor poly(x^2 - 6*x + 8, [x])
          Info: increasing working precision to DIGITS=20

          [0.3333333333 I, 0.3333333333 I, 0.3333333333 I,

             0.3333333333 I, 0.3333333333 I, 2.0, 4.0]

    >> setuserinfo(numeric::polyroots, 0): delete p:

#### Example 4

The following polynomial has badly separated roots. `numeric::polyroots` does not manage to separate them properly:

    >> p := poly(_mult((x - 1 - i/10^9) $ i=0..5), [x]):

    >> numeric::polyroots(p)

          [0.9999999987, 1.000000003, 1.000000003, 1.000000003,

             1.000000003, 1.000000003]

One can preprocess the polynomial by a symbolic factorization:

    >> numeric::polyroots(p, Factor)

          [1.0, 1.000000001, 1.000000002, 1.000000003, 1.000000004,

             1.000000005]

Alternatively, one can increase the working precision to separate the roots:

    >> DIGITS := 20: numeric::polyroots(p)

          [1.0, 1.000000001, 1.000000002, 1.000000003, 1.000000004,

             1.000000005]

    >> delete p, DIGITS:

#### Background

- The numerical root finding algorithm implemented by `numeric::polyroots` is Laguerre's method: W.H. Press, B.P. Flannery, S.A. Teukolsky and W.T. Vetterling: Numerical Recipes in C, Cambridge University Press, 1988.

#### Changes

- new function

Do you have questions or comments?
Copyright © SciFace Software GmbH & Co. KG 2000