I'm trying to get solution of cubic equations analytically in R, not numerically.
I looked up on the internet and get the formula for cubic roots and wrote the following code: The link is: http://www.math.vanderbilt.edu/~schectex/courses/cubic/
cub <- function(a,b,c,d) {
p <- -b/3/a
q <- p^3 + (b*c-3*a*(d))/(6*a^2)
r <- c/3/a
x <- (q+(q^2+(r-p^2)^3)^0.5)^(1/3)+(q-(q^2+(r-p^2)^3)^0.5)^(1/3)+p
x
}
However this function doesn't work in most cases and I guess it's because of the power of negative numbers inside the formula, for example I noticed R cannot get the real root of (-8)^(1/3) which is -2. But Im not sure how I could fix my code so that it can be used to solve for exact cubic solutions in general.
Thanks.
I'd use
polyroot()
. See here for more details.Here is a function to compute all the analytical solutions: 'cubsol' . Any comments would be most welcome. One question - at the moment the code searches rather inefficiently for which the real solution is amongst the three complex ones produced by ... s2 = cuberoot(q-s0^0.5); xtemp[1:3] <- s1+ s2 +p; Is there a more efficient way of knowing which one it would be before calculating it?
Try this:
If you only want the real root then here is another option:
Or you may be interested in the Ryacas package or the polynom package for other options.