The Quartic spreadsheet presented here previously (most recently here) uses algebraic techniques to provide an “exact” solution to polynomial equations of up to fouth order. For polynomials of higher order there is no general algebraic solution, and numerical methods are required. Also under some circumstances the “exact” solution, when evaluated with floating point numbers of limited precision, produces results that are very far from exact. In these circumstances a well written iterative numerical procedure can give much better results.
One such procedure is the Jenkins-Traub algorithm, which apart from the advantage of one author having an excellent family name, is also said to be both fast and stable. The Wikipedia article says that it “has been described as practically a standard in black-box polynomial root-finders”.
I have added a VBA user defined function (UDF) implementation of the Jenkins-Traub method to the Quartic spreadsheet, which has been re-named Polynomial.xlsb, and is available for free download, including open-source code.
The VBA code is a translation of Fortran77 code at the NetLib library. Similar code in Fortran90 format can be found at Alan Miller’s Fortran Software. Both of these sites are invaluable sources of open source scientific Fortran code. The chosen routine provides all real and complex solutions to any polynomial with real coefficients. Fortran routines are also available for polynomials with complex coefficients, but I have not yet translated these.
The screen-shots below show output from the Jenkins-Traub UDF (RPolyJT) for a quartic polynomial, followed by the results from the Quartic UDF for the same input. Note that the Jenkins-Traub result is close to machine precision, but in this case the Quartic output has errors in the third decimal place.
Output from RPolyJT for higher order polynomials is shown in the screenshots below, where the function has been used to solve the equation : 2.x^60 – 1 = 0. The roots are 60 equally spaced points in the complex plane, with a modulus of 0.5^(1/60). Note that the input is the range A26:A86; the cells containing zero coefficients could have been left blank, with the same results. The top and bottom of the input and output results are shown in the two screenshots below:
Great! I started using your code right away. I have a nonlinear system of equation which can be reduced to a quartic equation in a single unknown, and both Quartic and PolyRJT gave accurate and fast results. Since the overall computational time is similar, in my case, and you warn that PolyRJT can be more accurate, I’m using it instead of Quartic, but both worked really fine. I’m extremely grateful, and booked your website for future reference.
LikeLike
Pingback: Ping Back Ping Backs | Newton Excel Bach, not (just) an Excel Blog
Pingback: Daily Download 19: Solving polynomials | Newton Excel Bach, not (just) an Excel Blog
Pingback: what Microsoft think VBA is good for … | Newton Excel Bach, not (just) an Excel Blog
Great function! Solved my problem. Many many many thanks
LikeLike
Pingback: Solving Quadratic, Cubic, Quartic and higher order equations; examples | Newton Excel Bach, not (just) an Excel Blog
Pingback: Solving polynomials – update | Newton Excel Bach, not (just) an Excel Blog
I’ve been working on the couch problem. (what is the widest rectangle which will turn a corner between a 3 foot hallway and a 4 foot hallway. A classic problem. it requires solving a 6th degree equation. Normal Excel won’t do it, but your program worked like a charm. The complete solution has further steps, and now I can do the whole package. Thank you very much.
LikeLike
Tom – edited as requested (I don’t think visitors can edit).
Glad you found it useful!
(But now I’m going to have to try to work out how to solve it with a 6th degree polynomial)
LikeLike
@Tom – Interesting problem… Do you mean the sofa is 3 x W where the maximum value of W is somewhere between 3 and 4?
If so, probably overlooked something but from the back of an envelope i get the following formula for W:
=7*t/2+1/(2*t)-3*(t-1)/(t+1)
where
t = tan(x/2)
and x is the angle of rotation. The first order condition for the maximum is a quartic function of t with a real root which can also be calculated using Excel functions:=1/(1+IRR({-1,-2,-6,14,7}))
LikeLike
Anyone wondering how IRR solves polynomials (as I did when I first saw it) might like a look at:
https://newtonexcelbach.wordpress.com/2011/12/01/linest-npv-irr-and-solving-polynomials/
LikeLike
Pingback: Year 10 Report | Newton Excel Bach, not (just) an Excel Blog