8438: implement general cubic solution

The approach below cannot be used without special handling
to make sure that a valid pair of roots for _S and T are chosen
such that _S*T = -Q. [Thanks to Kalevi Suominen for noting why
the formulas below did not work as written for all cases].

** eqns (54)-(56) from http://mathworld.wolfram.com/CubicFormula.html

f = Poly([1, -S(3)/2 - 7*I/2, -2, -3], x)
_, a2, a1, a0 = f.monic().all_coeffs()
coeff = I*sqrt(3)/2
Q = (3*a1 - a2**2)/9
R = (9*a2*a1 - 27*a0 - 2*a2**3)/54
D = Q**3 + R**2
_S = root(R + sqrt(D), 3)
T = root(R - sqrt(D), 3)
B = _S + T
A = _S - T
r = -a2/3 - B/2
bad = [-a2/3 + B, r + coeff*A, r - coeff*A]

Although Vieta's approach seems to work, it produces largers expressions
and thereore takes longer:

def vieta(eq, x):
 # https://en.wikipedia.org/wiki/Cubic_function, Vieta's substitution
 f = Poly(eq, x).monic()
 _,a,b,c = f.all_coeffs()
 d = f.as_expr().subs(x, x-b/3/a) if a else f.as_expr()
 q = d.subs(x,0)
 p = d.subs(x,1)-q-1
 w3 = solve(x**2 + q*x - p**3/27, x)
 w = solve(x**3 - w3[0], x)
 t = [w - p/3/w for w in w]
 return [t - b/3/a for t in t] if a else t

If interested, compare these solutions with those produced
by solve:

vieta(Poly([1, y, -2, -3], x), x)
vieta(Poly([1, 4, -2, -3], x), x)
vieta(Poly([1, I, -2, -3], x), x)