Article 5XVHN Numerically finding roots of a cubic

Numerically finding roots of a cubic

by
John
from John D. Cook on (#5XVHN)

The analog of the quadratic formula for cubic equations is cumbersome. A lot of people naturally say Forget all that. If I need to find the roots of a cubic, I'll just use a numerical method like Newton's method." Sounds good.

Where to start?

But how do you know where to look for the roots? You have to have some idea where a root is before you can use a numerical method.

For the bisection method, you have to know an interval where the polynomial changes sign. For Newton's method, if you don't start close enough, you diverge. In fact, you can make complex fractals by plotting the points where Newton's method converges and diverges.

Cornelius Lanczos gives an elegant explanation of how to solve cubic equations at the beginning of his book Applied Analysis. His book was published in 1956, before the days of ubiquitous computers, and so he couldn't say Plot the function and scroll around until you can see roughly where the zeros are."

Lanczos does a sequence of change of variables that reduce the problem to finding a unique root in the unit interval.

Setup

Let

f(x) = ax^3 + bx^2 + cx + d

where the coefficients are all real numbers. We want to find all real values of x where f(x) = 0.

The main task is to find a root of f(x). Once we've found a root r, we can divide by (x - r) to reduce the problem to a quadratic.

Reduce to standard form

By a change of variables we can reduce the problem to finding a root of

x^3 + bx^2 + cx - 1 = 0.

How? First of all, you can always make the leading coefficient 1 by dividing by a if necessary. Second, you can make the constant term negative by replacing x by a new variable equal to -x. Third, you can make the constant term -1 by replacing x with d-1/3x.

We've divided by a and d above. How do we know they're not zero? If they were zero, that would be good news because it would make things easier.

We assume a is not zero because otherwise we don't have a cubic equation. And we assume d is not zero because otherwise we can factor out an x and we're left with a quadratic equation.

Bracketing a root

After our changes of variable,

f(x) = x^3 + bx^2 + cx - 1.

So f(0) = -1 and f(x) goes to infinity as x goes to infinity, and this means f(x) must have at least one positive root.

Evaluate f(1). If f(1) = 0 then you're very lucky because you've found your root.

If f(1) > 0 then there is a root in (0, 1) since the polynomial changes signs over the interval.

If f(-1) < 0 then the function has a root in (1,) . In that case, replace x with 1/x and then the new polynomial has a root somewhere in (0, 1).

Uniqueness

If f(x) equals zero somewhere in the interval (0, 1), could it be zero at two or three places in the interval?

It cannot have just two roots in the interval because the sign is different at the two ends of the interval.

It also cannot have three roots in the interval. We know from Vieta's formula that the product of the roots of f equals 1. It cannot have three roots in the interval (0, 1) because the product of three positive numbers less than 1 is less than 1.

If f(x) had a zero somewhere in (1,) before replacing x with 1/x, an analogous argument applies. It can't have two zeros because it changes signs over the interval, and it can't have three zeros in the interval because if three numbers are each greater than 1, their product is greater than 1.

In either case, possibly after replacing x with 1/x, we know f has a unique root in the interval (0, 1), and now we can apply a numerical method.

More general functions

This post looked at a very special function, a cubic polynomial with real coefficients, and showed that after a change of variables it has a unique root in an interval. The next post shows how to determine how many zeros a more general, complex-valued function has in a region of the complex plane.

Related postsThe post Numerically finding roots of a cubic first appeared on John D. Cook.
External Content
Source RSS or Atom Feed
Feed Location http://feeds.feedburner.com/TheEndeavour?format=xml
Feed Title John D. Cook
Feed Link https://www.johndcook.com/blog
Reply 0 comments