solve_cubic_equation_real.sf

#!/usr/bin/ruby

# Find a real solution to a cubic equation, using reduction to a depressed cubic, followed by the Cardano formula.

# Dividing ax^3 + bx^2 + cx + d = 0 by `a` and substituting `t - b/(3a)` for x we get the equation:
#   t^3 + pt + q = 0

# This allows us to use the Cardano formula to solve for `t`, which gives us:
#   x = t - b/(3a)

# Example (with x = 79443853):
#    15 x^3 - 22 x^2 + 8 x - 7520940423059310542039581 = 0

# See also:
#   https://en.wikipedia.org/wiki/Cubic_function

func solve_cubic_equation (a, b, c, d) {

    var p = (3*a*c - b*b)/(3*a*a)
    var q = (2*b**3 - 9*a*b*c + 27*a*a*d)/(27 * a**3)

    var t = (cbrt(-(q/2) + sqrt((q**2 / 4) + (p**3 / 27))) +
             cbrt(-(q/2) - sqrt((q**2 / 4) + (p**3 / 27))))

    var x = (t - b/(3*a))

    return x
}

say solve_cubic_equation(15, -22, 8, -7520940423059310542039581).round    #=> 79443853