Automatic differentiation using dual numbers

[Albert Chan]

An example that c = (1 - 1/2 * (f/f') * (f''/f')) can turn negative.

f = x^4 - 30*x^3 +63*x^2 - 12*x + 33      // real roots ≈ [2.31090, 27.7432]
f' = 4*x^3 - 90*x^2 + 126*x - 12              // real roots ≈ [0.102744, 1.38994, 21.0073]
f'' = 12*x^2 -180*x + 126                         // real roots ≈ [0.736125, 14.2639]

m = f'^2 - 1/2*f*f''                                    // real roots ≈ [-0.34128, 0.484147]
c = m / f'^2

c = ±∞ when guess = f' roots (Newton's correction also ±∞)
c = 0 when guess = m roots (Halley's method hit by divide-by-zero)
c < 0 when guess within m roots (iterate opposite direction of Newton's)

Comment:

Above example is easier to solve for roots of f(1/x), then flip the results.

g = 33*x^4 - 12*x^3 + 63*x^2 - 30*x + 1
g' = 132*x^3 - 36*x^2 + 126*x - 30
g'' = 396*x^2 - 72*x + 126

g'', with negative quadratic discriminant, have no root
g'' > 0 ⇒ g is concave up.