sagemath / sage

Main repository of SageMath
https://www.sagemath.org
Other
1.47k stars 487 forks source link

Spectral radius of graphs #20228

Closed videlec closed 8 years ago

videlec commented 8 years ago

There are no current way to compute the spectral radius of a graph (the Perron-Frobenius eigenvalue of its adjacency matrix). We provide a method spectral_radius that provides an interval of floating approximations.

It is much faster than previously available methods (where the exact methods are unusable because of characteristic polynomial computation)

sage: G = digraphs.RandomDirectedGNM(10,40)
sage: %timeit max(G.adjacency_matrix().charpoly().roots(AA, multiplicities=False))
100 loops, best of 3: 6.23 ms per loop
sage: %timeit G.spectral_radius()
1000 loops, best of 3: 178 µs per loop

And for a larger graph

sage: G = digraphs.RandomDirectedGNM(400, 6000)
sage: %timeit max(G.adjacency_matrix().charpoly().roots(AA, multiplicities=False))
1 loop, best of 3: 4.63 s per loop
sage: %timeit G.spectral_radius()
10 loops, best of 3: 13 ms per loop

The new method spectral_radius can be used for much larger graph as soon as the spectral gap is large enough to ensure the convergence.

Depends on #20210

Component: graph theory

Keywords: days71

Author: Vincent Delecroix

Branch/Commit: 89f2c10

Reviewer: Maurizio Monge

Issue created by migration from https://trac.sagemath.org/ticket/20228

videlec commented 8 years ago

New commits:

dd2dd58Trac 20228: spectral radius
videlec commented 8 years ago

Commit: dd2dd58

videlec commented 8 years ago

Branch: u/vdelecroix/20228

videlec commented 8 years ago

Description changed:

--- 
+++ 
@@ -18,4 +18,4 @@
 sage: %timeit G.spectral_radius()
 10 loops, best of 3: 13 ms per loop

-The new method spectral_radius can be used for much larger graph as soon as the spectral graph is large enough to ensure the convergence. +The new method spectral_radius can be used for much larger graph as soon as the spectral gap is large enough to ensure the convergence.

234549f8-15c0-4ce6-b5ca-e861340d063a commented 8 years ago

Reviewer: Maurizio Monge

234549f8-15c0-4ce6-b5ca-e861340d063a commented 8 years ago
comment:3

Is there a reason for not using Sage's interval arithmetic? http://doc.sagemath.org/html/en/reference/rings_numerical/sage/rings/real_mpfi.html Sage's interval library is already taking care of lowlevel floating point rounding problems, and could remove complexity from the function. Apart from this everything looks fine.

234549f8-15c0-4ce6-b5ca-e861340d063a commented 8 years ago

Changed keywords from none to days71

videlec commented 8 years ago
comment:4

Hello Maurizio,

Thanks for having a look!

Interval arithmetic does not use machine floating point but mpfr. It would be much slower than the version I wrote. This would only be a constant time slower which on my computer seems to be a factor x8. The intended application of this method is for huge sparse graphs (let say with ~1000000 vertices) and a x8 speedup is not negligible. Similarly, I could have used ball arithmetic arb which is also shipped with Sage (and should be faster than interval arithmetic).

Moreover, it is a certain amount of work to write matrix/vector code using either mpir or arb. Using mpir or arb would also decrease readability since multiplying numbers with mpfr is done via mpfr_mul(x, a, b, MPFR_RNDN) which is much less readable than x = a*b. Though the approach using either mpir or arb would have the advantage of making available higher precision. This would indeed be interesting for small graphs.

On the other hand, the function I made is written in C. It is not reasonable to use the Cython wrapper we have in Sage as allocating/deallocating million of Python classes is expensive.

The only thing that can easily be done would be to return an element in RIF instead of a tuple of Python float.

What I wrote is not the "best and fastest algorithm" to get the Perron-Frobenius eigenvalue. For example, it suffers from slow convergence in absence of spectral gap. But it was very helpful to me and might as well be helpful for others. I consider than having an arbitrary precision version using mpfr, mpir or arb would be indeed good for another ticket. It will need some care because you want to increase slowly the precision (it makes no sense to use the full precision on the first iterations). Moreover, it is not immediate to determine what is the precision needed for the computation to obtain a given precision for the answer.

I need to fix a mistake (the condition in the loop should be (e_max - e_min) > e_max * c_prec instead of e_max * (e_max - e_min) > c_prec). I would also like to move the sage_free inside the finally block (to avoid unfreed memory). I am waiting for your comments before submitting the commit.

234549f8-15c0-4ce6-b5ca-e861340d063a commented 8 years ago
comment:5

Hi Vincent, Sure, yours is a really good reason. When you fix the mistakes you are talking about you have my approval for the patch.

videlec commented 8 years ago
comment:6

Replying to @sagetrac-maurimo:

Hi Vincent, Sure, yours is a really good reason. When you fix the mistakes you are talking about you have my approval for the patch.

Nope. The mistakes have to be fixed. It was about the rest. For example, I am not sure whether it makes more sense to return a pair of float or an element of RIF (the standard field for interval arithmetic).

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 8 years ago

Changed commit from dd2dd58 to 64b2d44

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 8 years ago

Branch pushed to git repo; I updated commit sha1. New commits:

64b2d44Trac 20228: cleaning
234549f8-15c0-4ce6-b5ca-e861340d063a commented 8 years ago
comment:9

Well, In my opinion a pair of float is fine, if you are not using intervals internally we can avoid an unnecessary dependency, a user can easily create an interval from the float pair if necessary. Waiting for your last word, and I will accept your patch.

videlec commented 8 years ago
comment:10

I have nothing to add. Thanks.

videlec commented 8 years ago
comment:11

There is a conflict with #20210. I need to rebase this ticket over it.

videlec commented 8 years ago

Dependencies: #20210

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 8 years ago

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

03458eaUpgrade cysignals package
dce67fcMove memory functions to cysignals
4bb8337Rename sage_malloc -> sig_malloc and friends
1d89b0fTrac 20228: spectral radius
bca45c5Trac 20228: cleaning
89f2c10Trac 20228: memory allocation: sage_ -> sig_
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 8 years ago

Changed commit from 64b2d44 to 89f2c10

vbraun commented 8 years ago

Changed branch from u/vdelecroix/20228 to 89f2c10