Incorrect roots found
I have been testing the library on some polynomials with integer coefficients. On many of them incorrect roots are found. One of the examples:
Coefficients: [1,0,-1,0,1,1,1,1,1,1,1,0,-1,0,1]
Roots found:
1.1910 + 0.6054i
1.1910 - 0.6054i
-1.2766 + 0.3311i
-1.2766 - 0.3311i
-0.9114 + 0.8715i
-0.9114 - 0.8715i
-0.2781 + 1.1229i
-0.2781 - 1.1229i
0.4173 + 0.9860i
0.4173 - 0.9860i
0.6739 + 0.3241i
0.6739 - 0.3241i
0.1838 + 0.2115i
0.1838 - 0.2115i
True roots:
1.2051 + 0.5795i
1.2051 - 0.5795i
-0.5464 + 0.8375i
-0.5464 - 0.8375i
-0.9793 + 0.2023i
-0.9793 - 0.2023i
-0.8508 + 0.5254i
-0.8508 - 0.5254i
0.5333 + 0.8459i
0.5333 - 0.8459i
-0.0358 + 0.9994i
-0.0358 - 0.9994i
0.6739 + 0.3241i
0.6739 - 0.3241i
Roots are sorted by magnitude.
I cannot reproduce the issue. This is the result I get:
(0.673943,0.3241)
(0.673943,-0.3241)
(-0.85084,0.525425)
(-0.85084,-0.525425)
(-0.979328,0.202277)
(-0.979328,-0.202277)
(0.533286,0.845935)
(0.533286,-0.845935)
(-0.0357928,0.999359)
(-0.0357928,-0.999359)
(-0.546374,0.837541)
(-0.546374,-0.837541)
(1.20511,0.579537)
(1.20511,-0.579537)
I did few other tests, accuracy varies greatly depending on parameters (tolerances, number of shifts, etc.). What parameters do you use? What are residuals for each root? I use defaults (Intel C++ on Windows).
I work with reciprocal polynomials with integer coefficients, most of them have ill-conditioned roots (especially cyclotomic polys).
So far I was not able to get consistent accuracy with RpolyPlusPlus (one set of parameters work for one polynomial but not for the other). I compare it with DHSEQR from MKL/LAPACK, which is very stable (but slow).
One more thing, it is more numerically robust to use std::hypot to compute magnitude of a complex number. It takes care of the overflow/underflow and other accuracy losses.
I used the default parameters. Did you enable any particular optimization flags when compiling using Intel C++?