How to find the first root of a cubic equation in matlab

Question

Question: Given a cubic equation where the coefficients are real but can take any extreme conditions (e.g very large or very small number), write a program in Matlab that finds all the roots of this equation. You can't use the built-in functions roots and fzero

Here's my approach: By the Fundamental Theorem of Algebra, f(x) must have at least one real root. If I can find this root, I can factor it out and find the other two roots by using the quadratic formula. But I stuck badly on finding the first real root. I'm thinking about the bisection method, but how can I even find the interval that contains the first root. Can anyone help me at this step?

Thanks!

Answer 1

Let's assume the coeff of $x^3$ is positive; the negative case is similar.

The derivative is a quadratic, whose roots $r_1 < r_2$ you can find. (If it has no roots, see below).

To the right of $x = r_2$, the function is increasing, so if $f(r^2) > 0$, then there's no root there, and if $f(r_2) < 0$, then there's at most one root there. A similar argument applies to $f(r_1)$. So we have three cases:

1. $f(r_1) < 0, f(r_2) > 0$: there's a root between them.

2. $f(r_1)$ and $f(r_2)$ are both positive: there's a root to the left of $r_1$.

3. $f(r_1)$ and $f(r_2)$ are both negative: there's a root to the right of $r_2$.

Here's what you can do in case 3: look at $f(r_2 + 1)$, $f(r_2 + 2)$, $f(r_2 + 4)$, ... $f(r_2 + 2^n)$. Eventually one of these will be positive. Then between $r^2$ and this point, there's a root, which you can find by the intermediate value theorem.

A similar thing works for case 2. Case 1 is easy.

And what if the roots of the derivative are both complex? Then the slope is everywhere positive. Look at the points $x =\pm 1$, $x = \pm 2$, ..., $x = \pm 2^n$, increasing $n$ until the signs differ. Now use bisection to find the root.

Answer 2

Normalize so that the leading coefficient is $1$. By the triangle inequality:

$$|x^3+ax^2+bx+c| \geq |x|^3-|a||x|^2-|b||x|+|c|.$$

Now let $R=\max \{ 1,|a|+|b|+|c| \}$. Then the inequality above implies that all roots must lie in $[-R,R]$. Thus we can find one of the real roots (or the only one) to within an accuracy of $\varepsilon$ with essentially $\log_2(R/\varepsilon)$ bisection steps with the initial interval $[-R,R]$.

Now this might be too slow, if $R$ is big or $\varepsilon$ is small. Indeed, this simple method is not "secure", in the sense that I can cook up a problem which is within the domain of your problem statement and which will make this method take an arbitrarily long time to run. Yet somehow these problems with huge $R$ should actually be very easy.

And indeed, better methods exist which can handle large $R$. As usual it is helpful to know something analytical about the problem in advance. In this case we can find the extrema explicitly using the quadratic formula. If none exist, then the cubic is increasing, which means we can apply Newton's method starting anywhere we want. It is impossible for only one extremum to exist in a cubic (why?) If two extrema $r_1 < r_2$ exist, then you can compute the value of the cubic at the two extrema. If the signs differ, then the root is in between them, so you can run Newton's method starting at $\frac{r_1+r_2}{2}$. If they are both negative, then the root is to the right of them, so you can run Newton's method starting at $\frac{r_2+R}{2}$. If they are both positive, then the root is to the left of them, so you can run Newton's method starting at $\frac{r_1-R}{2}$.

As an example, with $a=-169,b=45,c=-29$, we have two extrema, one around $0$ and one at about $110$. Both are negative and $R$ is about $250$, so we start Newton's method at about $180$ and get a root of about $169$ in four steps. Bisection would have required 40-50 steps.

Answer 3

This may help: If $x^3 +A x^2+B x+C=0$ then $$|x|<1+\max (|A|,|B|,|C|).$$ And if $C\ne 0$ then $$|x|>(1+\max (|B/C|,|A/C|,|1/C|))^{-1}.$$ Sometimes this can be sharpened by letting $x=k y$ with $k>0$ because we have $$|x|<k(1+\max (|A/k|,|B/k^2|,|C/k^3|)=k+\max (|A|,|B/k|,|C/k^2|),$$ with a corresponding lower bound for |x| in terms of $A,B,C,k$ from the 2nd inequality. These inequalities are verifiable by entirely elementary algebra.