# eigen-decomposition of a special companion matrix

I have a special type of companion matrix, where the "special" part is that each element in the matrix are matrices. For instance, the diagonal with "1":s is instead a diagonal with identity matrices, and so on... My question is if anyone knows of a solid algorithm for finding the eigen-decomposition of such a matrix? If you don't know the exact algorithm, perhaps you know of some work that has been done in the area?

If anyone is interested, I can clarify what my eigendecomposition will be used; AND this poses another problem, which is really the core problem I am facing.

I want to simulate a companion matrix (p x p), with each element being a matrix (n x n). The only restriction I have is that all eigenvalues must satisfy |λ|<1 . Ultimately, I would like to simulate the eigenvalues from (−1,1) , and then produce a companion matrix from that. I realize there might be (most likely are) multiple solutions, but that does not matter. I would settle for any solution. :)

My idea of solving this is to look at the general eigendecomposition of a companion matrix, simulate the eigenvalues, possible simulate the eigenvectors, and then reproduce the companion matrix.

Am I making any sense with this? I hope you understand the issue I have, and if not, please feel free to ask. As I said above, this area is fairly new to me...

-

You can get some savings with respect to the naive $O(n^3d^3)$ by using the semiseparable structure, but I don't think you can get anything faster than $O(n^3d^2)$ (stable) or $O(n^3d\log d)$ (maybe unstable). Algorithms for semiseparable matrices aren't exactly easy to start working with unless you work in numerical linear algebra, so unless $d$ is large, this isn't typically worth the trouble.
Inside Matlab there are some computations that use this problem (see polyeig), but as far as I know they run unoptimized QR.