Finding more than one root with root finding algorithms

Question

What is the best way to find the roots of an equation with more than one root. I understand that no one method which can solve every equation, and that you have to use more than one, but I can't find a root finding algorithm that can solve for more than one root in even the simplest instance.

For example: y = x^2

Although a root solving algorithm to solve a basic equation like this is helpful, it would need to be something I could adapt to solve an equation with more than two roots.

One more thing to note is that the equations wouldn't be your typical polynomials, but could be something such as ln(x^2) + x - 15 = 0

What is a root finding algorithm that could solve for this, or how could you edit a root finding algorithm such as the Bisection/Newton/Brent method to solve this problem, (Assuming I'm correct in that Newton and Brent's method can only solve for one root).

Answer 1

As everyone has said, it is impossible to provide a general algorithm to find all, or some of the roots (where some is greater than one.) And how many roots is some? You cannot find all of the roots in general, since many functions will have infinitely many roots.

Even methods like Newton do not always converge to a solution. I tend to like a good, fairly stable method that will converge to a solution under reasonable circumstances, such as a bracket where the function is known to change sign. You can find such a code that has good convergence behavior on single roots, yet still is protected to behave basically like a bisection scheme when the function is less well behaved.

So, given a decent root finding scheme, you can try simple things like deflation. Thus, consider a simple function, like a first kind Bessel function. I'll do all my examples using MATLAB, but any tool that has a stable well written rootfinder like fzero in MATLAB will suffice.

ezplot(@(x) besselj(0,x),[0,10])
grid on

f0 = @(x) besselj(0,x);
xroots(1) = fzero(f0,1)
xroots =
    2.4048

From the plot, we can see there is a second root around 5 or 6.

Now, deflate f0 for that root, creating a new function based on f0, but one that lacks a root at xroots(1).

f1 = @(x) f0(x)./(x-xroots(1));
ezplot(f1,[0,10])
grid on

Note that in this curve, the root of f0 at xroots(1) has now been zapped away, as if it did not exist. Can we find a second root?

xroots(2) = fzero(f1,2)
xroots =
    2.4048    5.5201

We can go on of course, but at some point this scheme will fail due to numerical issues. And that failure won't take too terribly long either.

A better scheme might be to (for 1-d problems) use a bracketing scheme, conjoined with a sampling methodology. I say better because it does not require modifying the initial function to deflate the roots. (For 2-d or higher, things get far more hairy of course.)

xlist = (0:1:10)';
flist = f0(xlist);

[xlist,flist]
ans =
         0    1.0000
    1.0000    0.7652
    2.0000    0.2239
    3.0000   -0.2601
    4.0000   -0.3971
    5.0000   -0.1776
    6.0000    0.1506
    7.0000    0.3001
    8.0000    0.1717
    9.0000   -0.0903
   10.0000   -0.2459

As you can see, the function has sign changes in the intervals [2,3], [5,6], and [8,9]. A rootfinder that can search in a bracket will do here.

fzero(f0,[2,3])
ans =
    2.4048

fzero(f0,[5,6])
ans =
    5.5201

fzero(f0,[8,9])
ans =
    8.6537

Just look for sign changes, then throw the known bracket into a root finder. This will give as many solutions as you can find brackets.

Be advised, there are serious problems with the above scheme. It will completely fail to find the double root of a simple function like f(x)=x^2, because no sign change exists. And if you choose too coarse of a sampling, then you may have an interval with TWO roots in it, but you won't see a sign change at the endpoints.

For example, consider the function f(x) = x^2-x, which has single roots at 0 and at 1. But if you sample that function at -1 and 2, you will find that it is positive at both points. There is no sign change, but there are two roots.

Again, NO method can be made perfect. You can ALWAYS devise a function that will cause any such numerical method to fail.

Answer 2

Bottom line: You need to isolate the roots yourself.

Details depend on the algorithm: If you're using bisection or Brent's method, you need to come up with a set of intervals each containing a unique root. If using the Newton's method, you need to come up with a set of starting estimates (since it converges to a root given a starting point, and with different starting points it may or may not converge to different roots).

Answer 3

I'd say that there's no general method to find all the roots of a general equation. However, one can try and devise methodologies once sufficient conditions have been specified. Even simple quadratic equations ax² + bx + c = 0 aren't completely trivial, because the existence of real roots depends on the sign of b²-4ac, which isn't immediately obvious. So there are lots of techniques to apply, e.g Newton-Raphson, but no general method for the general case, especially for equations like ln(x²)+x-15 = 0.