How to solve the polynomial eigenvalue in python?

Question

In my python code, I would like to solve the polynomial eigenvalue problem:

A0 + lambda*A1 + lambda^2*A2 + lambda^3*A3 + .... = 0

where An are dense matrices, and lambda is a constant. In matlab it is possible to solve this problem using the polyeig function. It seems that there is no equivalent functionality in scipy. So far the only way I can think to do it is to form the corresponding companion matrix. This creates an equivalent linear eigenvalue problem which can be given to existing scipy solvers, however it is much larger and I believe it can be quite ill-conditioned.

Can anyone suggest an existing, open source or freely-available library which can solve this problem? I'd be quite happy with a fortran library which could be linked via f2py or C/C++ library to link via cython.

Edit: For anyone interested in solving nonlinear eigenvalue problems in python, the code I wrote myself to solve this problem can be found here. Note that I deal with the more general case of a nonlinear eigenvalue problem (in the sense that it has a nonlinear dependence on lambda). To understand the method, please read the paper mentioned in the code comments.

Answer 1

This discussion points to a general method for turning a polynomial eigenvalue problem into a generalized eigenvalue problem, which can later be solved using scipy's linear algebra functions. Hope this helps!