I understand the logic behind the polynomial rootfinding algorithms that use companion matrices. I checked the numpy implementation and it's pretty straightforward but I don't understand this comment:
rotated companion matrix reduces error
Why is this true mathematically? I don't even understand why would this give the same roots as the original companion matrix.
Thanks in advance for the help!
From binomial theorems we know that the polynomial division $$ q(a,b)=\frac{p(b)-p(a)}{b-a} $$ is without remainder. So if $a$ is a root and $b=x$, then $$ p(x)=(x-a)q(a,x) $$ is a factorization into a linear factor and a polynomial of one degree less. One can write this equation also equivalently as $$ xq(a,x)\equiv aq(a,x)\pmod{p(x)} $$ or $$ (x\cdot q(a,x))\bmod p(x) = a\cdot q(a,x) $$ This is a linear system for the coefficients of $q_a(x)=q(a,x)$, especially a system for an eigenvalue-eigenvector pair and the basis for the root finders based on the companion matrix.
For some mix of tradition and convenience the usual basis is the monomial basis, ordered for increasing degrees. This give a system matrix with a filled last row and entries directly above the diagonal.
Most eigenvalue finders apply a variant of the QR algorithm. This gets simplified if the matrix is in Hessenberg form (upper triangular including one sub-diagonal), as that pattern is preserved by the algorithm.
If the (monomial) basis is ordered for falling degree, the companion matrix is automatically in Hessenberg form. This may reduce the number of operations and thus increase the accuracy, as each numerical matrix transformation moves the actually represented (characteristic) polynomial a little more away from the original polynomial.
But generally the roots found by
numpy.rootscan be substantially improved by one or two Newton iterations. While the eigenvalue algorithm finds all roots of a perturbed polynomial, there is a minimal risk that such a Newton step will map roots in a cluster to the same root of the original polynomial. But multiple roots and and small clusters are always a challenge for root-finding algorithms.