Re: [eigen] Polynomial solver, eigenvalues of companion matrix and balancing

[Gael Guennebaud <gael.guennebaud@xxxxxxxxx>]

just to say that I've tried your unit test using octave's eigen solver
I get the same behavior... sometime there is a very large error
showing up...

gael.


On Mon, Jan 18, 2010 at 7:39 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx> wrote:
> 2010/1/18 Manuel Yguel <manuel.yguel@xxxxxxxxx>:
>> On Mon, Jan 18, 2010 at 6:51 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx> wrote:
>>> 2010/1/18 Manuel Yguel <manuel.yguel@xxxxxxxxx>:
>>>> Hello,
>>>> I have written a class to compute the roots of a real polynomial.
>>>> I build the companion matrix and compute its eigenvalues.
>>>> You can see the code in the attached file (not a patch yet ... see the
>>>> following).
>>>> I have some problems with that method:
>>>>
>>>> A] When testing: I encounter almost systematic failure for polynomials
>>>> with deg greater than 7 (see test file attached),
>>>> therefore I wonder:
>>>>
>>>> 1) Do I use the right eigensolver ?
>>>
>>> If you only want to support real polynomials, then yes.
>>>
>>>> 2) Does the eigensolver balance the input matrix ?
>>>
>>> Not as far as I can see.
>> thanks, that what I guessed too.
>>>
>>>> I have written some code to balance a companion matrix explicitly (and
>>>> I am testing this stuff at the moment) but before investing more time
>>>> in that direction I need to know if doing the balancing is redundant.
>>>
>>> Ask Gael to be sure, but I don't think that we have this right now.
>>>
>>>> B] A somehow related question is: the eigensolver make a copy of the
>>>> matrix however, the companion matrix is sparse and for high degree
>>>> polynomial this copy could be avoided.
>>>> I have thought about writing the companion matrix class as a
>>>> cwiseNullaryOp. Do I follow the right thread?
>>>
>>> In all these decomposition algorithms, the work matrix is initialized
>>> at the start by copying the input matrix into it. This is how these
>>> algorithms work. If you want to preserve sparsity, you need a
>>> completely different algorithm: this one won't preserve sparsity at
>>> all. But I am not sure that companion matrices offer real
>>> opportunities to take advantage of sparsity here.
>> Sorry I was not clear enough, so let me explain it again:
>> a companion matrix (balanced) needs 2d-1 coefficients for a polynomial
>> of degree d.
>> There is a big memory gain, by not building an entire matrix then
>> copying it, but just initializing the matrix inside the eigensolver
>> algorithm by the 2d-1 coefficients.
>> After, it is up to the eigensolver to do whatever is needed with its matrix.
>> I just do not want to build an entire dxd matrix that will just be
>> copied afterward.
>
> Yes, this is what Gael is pointing out when he says that the companion
> matrix is already hessenberg. We need to provide an API for
> eigenvalues of a Hessenberg matrix. In the same vein, right now i'm
> adding SVD of a bidiagonal matrix.
>
> Benoit
>
>>
>> Manuel
>>>
>>> Benoit
>>>
>>>> P.S. I am also writing a solver with a Bézier bissection, it is
>>>> claimed to perform faster.
>>>
>>> Interesting!
>>>
>>>
>>>
>>
>>
>>
>
>
>