15

I have polynomials of nontrivial degree (4+) and need to robustly and efficiently determine whether or not they have a root in the interval [0,T]. The precise location or number of roots don't concern me, I just need to know if there is at least one.

Right now I'm using interval arithmetic as a quick check to see if I can prove that no roots can exist. If I can't, I'm using Jenkins-Traub to solve for all of the polynomial roots. This is obviously inefficient since it's checking for all real roots and finding their exact positions, information I don't end up needing.

Is there a standard algorithm I should be using? If not, are there any other efficient checks I could do before doing a full Jenkins-Traub solve for all roots?

For example, one optimization I could do is to check if my polynomial f(t) has the same sign at 0 and T. If not, there is obviously a root in the interval. If so, I can solve for the roots of f'(t) and evaluate f at all roots of f' in the interval [0,T]. f(t) has no root in that interval if and only if all of these evaluations have the same sign as f(0) and f(T). This reduces the degree of the polynomial I have to root-find by one. Not a huge optimization, but perhaps better than nothing.

moinudin
  • 134,091
  • 45
  • 190
  • 216
user168715
  • 5,469
  • 1
  • 31
  • 42
  • 3
    I'd recommend http://math.stackexchange.com – Jim Brissom Dec 22 '10 at 00:59
  • I'm trying to imagine a recursive solution where you recurse on the derivative until you can solve it, and then use that result to determine the roots of the integral. But ultimately you will still need to find *actual* roots rather than just approximations; estimations are not very useful with high degree polynomials, which can fluctuate wildly. What degree polynomials are you talking about - only 4, or does it go much higher? – Kirk Broadhurst Dec 22 '10 at 05:07

5 Answers5

17

Sturm's theorem lets you calculate the number of real roots in the range (a, b). Given the number of roots, you know if there is at least one. From the bottom half of page 4 of this paper:

Let f(x) be a real polynomial. Denote it by f0(x) and its derivative f′(x) by f1(x). Proceed as in Euclid's algorithm to find

f0(x) = q1(x) · f1(x) − f2(x),
f1(x) = q2(x) · f2(x) − f3(x),
.
.
.
fk−2(x) = qk−1(x) · fk−1(x) − fk,

where fk is a constant, and for 1 ≤ i ≤ k, fi(x) is of degree lower than that of fi−1(x). The signs of the remainders are negated from those in the Euclid algorithm.

Note that the last non-vanishing remainder fk (or fk−1 when fk = 0) is a greatest common divisor of f(x) and f′(x). The sequence f0, f1,. . ., fk (or fk−1 when fk = 0) is called a Sturm sequence for the polynomial f.

Theorem 1 (Sturm's Theorem) The number of distinct real zeros of a polynomial f(x) with real coefficients in (a, b) is equal to the excess of the number of changes of sign in the sequence f0(a), ..., fk−1(a), fk over the number of changes of sign in the sequence f0(b), ..., fk−1(b), fk.

Ryan Dougherty
  • 528
  • 6
  • 21
moinudin
  • 134,091
  • 45
  • 190
  • 216
3

You could certainly do binary search on your interval arithmetic. Start with [0,T] and substitute it into your polynomial. If the result interval does not contain 0, you're done. If it does, divide the interval in 2 and recurse on each half. This scheme will find the approximate location of each root pretty quickly.

If you eventually get 4 separate intervals with a root, you know you are done. Otherwise, I think you need to get to intervals [x,y] where f'([x,y]) does not contain zero, meaning that the function is monotonically increasing or decreasing and hence contains at most one zero. Double roots might present a problem, I'd have to think more about that.

Edit: if you suspect a multiple root, find roots of f' using the same procedure.

Keith Randall
  • 22,985
  • 2
  • 35
  • 54
  • 1
    I've read this a few times and it's stumping me: "Start with [0,T] and substitute it into your polynomial. If the result interval does not contain 0, you're done." I can imagine plenty of polynomials with roots between 0 & T, but where p(0) and p(T) are both positive. Am I misunderstanding you? – Kirk Broadhurst Dec 22 '10 at 04:58
  • I'm talking about interval arithmetic as in http://en.m.wikipedia.org/wiki/Interval_arithmetic?wasRedirected=true. This doesn't evaluate at just 0 and T, but a (conservative approximation of) all points in between. – Keith Randall Dec 22 '10 at 06:15
1

Use Descartes rule of signs to glean some information. Just count the number of sign changes in the coefficients. This gives you an upper bound on the number of positive real roots. Consider the polynomial P.

P = 131.1 - 73.1*x + 52.425*x^2 - 62.875*x^3 - 69.225*x^4 + 11.225*x^5 + 9.45*x^6 + x^7

In fact, I've constructed P to have a simple list of roots. They are...

{-6, -4.75, -2, 1, 2.3, -i, +i}

Can we determine if there is a root in the interval [0,3]? Note that there is no sign change in the value of P at the endpoints.

P(0) = 131.1
P(3) = 4882.5

How many sign changes are there in the coefficients of P? There are 4 sign changes, so there may be as many as 4 positive roots.

But, now substitute x+3 for x into P. Thus

Q(x) = P(x+3) = ...
  4882.5 + 14494.75*x + 15363.9*x^2 + 8054.675*x^3 + 2319.9*x^4 + 370.325*x^5 + 30.45*x^6 + x^7

See that Q(x) has NO sign changes in the coefficients. All of the coefficients are positive values. Therefore there can be no roots larger than 3.

So there MAY be either 2 or 4 roots in the interval [0,3].

At least this tells you whether to bother looking at all. Of course, if the function has opposite signs on each end of the interval, we know there are an odd number of roots in that interval.

  • "So there MAY be either 2 or 4 roots in the interval [0,3]." - Or there may be zero roots in [0,3]. If there were sign changes in Q(x) then it only means that there may be roots > 3 but not that such roots necessarily exist. – Kirk Broadhurst Dec 23 '10 at 01:05
  • This is restating some of what the OP said in his "same sign at 0 and T" statement. Good when it shows that a root exists, but no good for determining when a root doesn't exist. – Kirk Broadhurst Dec 23 '10 at 01:08
  • The point is that the OP can at least some times exclude a root from existing. So Some of the time he need do no search. This is far better than looking for a root, but not being sure that he just did not search carefully enough. And since the cost is small, why not do so? –  Dec 23 '10 at 04:14
  • 1
    @Kirk - if the rule of signs says that no root exists, then there is no need to search. So it IS useful as a tool. –  Dec 23 '10 at 04:16
1

It's not that efficient, but is quite reliable. You can construct the polynomial's Companion Matrix (A sparse matrix whose eigenvalues are the polynomial's roots).

There are efficient eigenvalue algorithms that can find eigenvalues in a given interval. One of them is the inverse iteration (Can find eigenvalues closest to some input value. Just give the middle point of the interval as the above value).

Alex Shtoff
  • 2,520
  • 1
  • 25
  • 53
0

If the value f(0)*f(t)<=0 then you are guaranteed to have a root. Otherwise you can start splitting the domain into two parts (bisection) and check the values in the ends until you are confident there is no root in that segment.

if f(0)*f(t)>0 you either have no, two, four, .. roots. Your limit is the polynomial order. if f(0)*f(t)<0 you may have one, three, five, .. roots.

John Alexiou
  • 28,472
  • 11
  • 77
  • 133