Dear all, I am happy to announce the first version of CmplxRoots.jl <https://github.com/giordano/CmplxRoots.jl>, a fast root finder of real and complex polynomials.
This is a Julia implementation of the algorithm described in - J. Skowron & A. Gould, 2012, "General Complex Polynomial Root Solver and Its Further Optimization for Binary Microlenses", arXiv:1203.1034 <http://arxiv.org/abs/1203.1034> I first wrote the package as a wrapper around the Fortran implementation available at http://www.astrouw.edu.pl/~jskowron/cmplx_roots_sg/ but then I rewrote it entirely in Julia, finding that it's even faster than ccall'ing functions in the shared object. CmplxRoots.jl is a registered Julia package, so you can install it with the package manager: Pkg.add("CmplxRoots") Two functions are exposed to the users: roots(polynomial[, roots, polish=true]) roots5(polynomial[, roots]) The first function can be used to solve polynomials of any order, the second one is specialized and optimized for (and works only for) polynomials of fifth-order, but I found it to be a bit slower than roots(). The only mandatory argument is the list of coefficients of the polynomial, in ascending order. So, if you want to find all the roots of polynomial (x - i)*(x - 2)*(x - 3*i)*(x - 4)*(x - 5*i) = x^5 - (9*i + 6)*x^4 + (54*i - 15)*x^3 + (138 - 57*i)*x^2 - (184 + 90*i)*x + 120*i you can use julia> roots([120im, -(184 + 90im), (138 - 57im), (54im - 15), -(6 + 9im), 1 ]) 5-element Array{Complex{Float64},1}: -1.55431e-15+5.0im 4.0+1.77636e-16im 1.55028e-15+3.0im -1.04196e-16+1.0im 2.0-2.00595e-16im julia> roots5([120im, -(184 + 90im), (138 - 57im), (54im - 15), -(6 + 9im), 1]) 5-element Array{Complex{Float64},1}: 5.92119e-16+5.0im 4.0-1.4803e-16im 2.0+1.88202e-16im -1.88738e-15+3.0im 2.10942e-15+1.0im Optional roots argument is the array of starting values that can be used to find the final roots (it won't be changed in-place). Note that even if the roots are all real, the output array is always of type Complex{Float64}. Also Polynomials.jl package has a root finding function, but according to some random tests performed with @benchmark (from Benchmarks.jl package) it's roughly 20 times slower than roots() function provided by CmplxRoots.jl. Of course do not expect results with machine precision (you can have catastrophic cancellation even for simple quadratic polynomials), but some random tests indicate that CmplxRoots.jl is slightly more precise than Polynomials.jl (just compute abs2() of polynomials at the root). The package CmplxRoots.jl is licensed under the GNU Lesser General Public License version 3 or any later version, as well as under a "customary scientific license", which implies that if this code was important in the scientific process or for the results of your scientific work, you are asked for the appropriate citation of the paper Skowron & Gould 2012 ( http://arxiv.org/abs/1203.1034). Happy root-finding, Mosè
