Implementation of quartic equation solver doesn't work

Question

What's going wrong with my implementation of a quartic equation solver?
Here is my code on GitHub. I followed this: http://www.1728.com/quartic2.htm
In fact the real implementation starts at line 271, where I create the monic poly.

If I try it with a polynomial with 4 real roots it works fine (for example with 3x^4 + 6x^3 - 123x^2 - 126x + 1,080), otherwise gives wrong roots.

Thanks,
rubik

P.S. I called the function __quartic because it is still in development

Answer 1

If it isn't a problem with your cubic solver (as phkahler suggested), then I'm betting that it is a case of int division in these lines:

f = float(c - (3*b**2 / 8))
g = float(d + (b**3/ 8) - (b*c / 2))
h = e - (3*b**4 / 256) + (b**2 * c / 16) - (b*d / 4)

If you are using a version of Python that is less than 3.0 and not importing division from __future__, then you can have major problems here. For instance, if b = 3, then (3*b**2/8) is equal to 3 by int division rather than the correct float division value of 3.375. Probably the best way to correct this is to convert a,b,c,d and e to floats before you start doing calculations. Your conversion to float after all of the calculations for f and g are really unnecessary after you've made this correction and they aren't really helping you much now.

For more information about this division issue and its history see:

http://www.python.org/dev/peps/pep-0238/

http://python-history.blogspot.com/2009/03/problem-with-integer-division.html

Answer 2

Does your cubic root solver work? You should have it reduced to only returning a single real root for use within the quartic solver. Run a number of test cases on the cubic solver first, including equations with both 1 and 3 real roots. Only after this is verified should you try to debug the quartic solver.

Answer 3

Perhaps I can help out. I wrote a program and the procedure for solving quartic equations.

I'm guessing you are not writing the program in JavaScript. I'd be willing to help out (if I can).