Re: [sage-support] Quickly compute signs of eigenvalues

51 views
Skip to first unread message

William Stein

unread,
May 20, 2013, 8:37:55 PM5/20/13
to sage-s...@googlegroups.com
On Mon, May 20, 2013 at 2:19 PM, Theo Belaire
<tyr.god....@gmail.com> wrote:
> I have a large computation where I need to compute the number of positive eigenvalues of a matrix.
> I am currently computing all the eigenvalues then counting how many are positive, but I see when profiling that "{method 'roots' of 'sage.rings.polynomial.polynomial_element.Polynomial' objects}" is taking up most of the time.
> Any way to just get sign?
> Or to use estimates to speed up most of the computation, and I can recompute any that have more than 1 positive eigenvalues exactly. I'd really like to avoid any false negatives though.
>

Can you somehow give an example of the sort of polynomials you're considering?

What's the base field? RR, RealField(100), RDF, CDF, etc.?

For example:

R.<x> = RR['x']
f = x^50 + 7*x + 5
%timeit f.roots()
5 loops, best of 3: 97.6 ms per loop

versus

R.<x> = RDF['x']
f = x^50 + 7*x + 5
%timeit f.roots()
125 loops, best of 3: 5.74 ms per loop

> --
> You received this message because you are subscribed to the Google Groups "sage-support" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-support...@googlegroups.com.
> To post to this group, send email to sage-s...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sage-support?hl=en.
> For more options, visit https://groups.google.com/groups/opt_out.
>
>



--
William Stein
Professor of Mathematics
University of Washington
http://wstein.org

Theo Belaire

unread,
May 21, 2013, 8:39:48 AM5/21/13
to sage-s...@googlegroups.com
I haven't explicitly set the ring it's working over, but all the entries of the matrix are integral.

Jason Grout

unread,
May 21, 2013, 8:53:14 AM5/21/13
to sage-s...@googlegroups.com
On 5/21/13 7:39 AM, Theo Belaire wrote:
> I haven't explicitly set the ring it's working over, but all the entries
> of the matrix are integral.

Are your eigenvalues guaranteed to be a little away from 0, so you could
use numeric approximations? If so, computing the eigenvalues over the
RDF field, which uses some numeric approximation algorithms, may be much
faster than computing roots of integral characteristic polynomials.

Thanks,

Jason


Theo Belaire

unread,
May 21, 2013, 3:29:05 PM5/21/13
to sage-s...@googlegroups.com
I have a fixed number of 0 eigenvalues, and no real guarantee that the rest will be more than δ away for any δ.
I wouldn't mind only recomputing close ones exactly, if I can bound the error so I know I won't miss any.  I know that having more than one positive eigenvalues will be very rare, if it ever happens, but I really don't want to miss any of those cases.

Michael Orlitzky

unread,
May 21, 2013, 6:24:50 PM5/21/13
to sage-s...@googlegroups.com
On 05/20/2013 05:19 PM, Theo Belaire wrote:
> I have a large computation where I need to compute the number of
> positive eigenvalues of a matrix. I am currently computing all the
> eigenvalues then counting how many are positive, but I see when
> profiling that "{method 'roots' of
> 'sage.rings.polynomial.polynomial_element.Polynomial' objects}" is
> taking up most of the time. Any way to just get sign? Or to use
> estimates to speed up most of the computation, and I can recompute
> any that have more than 1 positive eigenvalues exactly. I'd really
> like to avoid any false negatives though.
>

Singular (included in sage) has a library that can quickly count the
number of real roots in an interval:

http://mate.dm.uba.ar/~alidick/etobis/rrc.pdf

I don't see an easy way to interface with it, but you could probably
come up with something based on other code that uses singular.

This seems to do what you want:

sturmha(p, 0, maxabs(p))

Jeroen Demeyer

unread,
May 22, 2013, 1:10:40 AM5/22/13
to sage-s...@googlegroups.com
On 05/22/2013 12:24 AM, Michael Orlitzky wrote:
> Singular (included in sage) has a library that can quickly count the
> number of real roots in an interval:
>
> http://mate.dm.uba.ar/~alidick/etobis/rrc.pdf
>
> I don't see an easy way to interface with it, but you could probably
> come up with something based on other code that uses singular.
>
> This seems to do what you want:
>
> sturmha(p, 0, maxabs(p))
>

Same for PARI, where the function is called polsturm().

gp> ?polsturm()
polsturm(pol,{a},{b}): number of real roots of the squarefree polynomial
pol in the interval ]a,b] (which are respectively taken to be -oo or +oo
when omitted).

Theo Belaire

unread,
May 22, 2013, 12:01:24 PM5/22/13
to sage-s...@googlegroups.com
This has cut my computation time to a hundredth of that it was, thanks.
Reply all
Reply to author
Forward
0 new messages