Opened 4 years ago

Closed 4 years ago

# infinite loop in spectral_radius for trees

Reported by: Owned by: dcoudert major sage-8.4 graph theory vdelecroix, maurimo Vincent Delecroix David Coudert N/A 7c8240d 7c8240d7c4c034c33dab4de89b3ee8673f56d0f0

### Description

The following instruction will never end

```sage: g = graphs.CompleteBipartiteGraph(1,3)
```

Adding a print in the while loop to get `e_max`, `e_min`, `c_prec`, `e_max-e_min`, `e_max * c_prec`, we get

```sage: g.spectral_radius()
3.0 1.0 1e-10 2.0 3e-10
3.0 1.0 1e-10 2.0 3e-10
3.0 1.0 1e-10 2.0 3e-10
...
```

and its the same for every tree. It's working well if the graph has a single edge, but for the path graph of order 3, we already get the infinite loop.

### comment:1 Changed 4 years ago by dcoudert

I don't know how to fix that.

### comment:2 Changed 4 years ago by vdelecroix

Indeed, this is very bad for any bipartite graph because the two eigenvalues are `lambda1` (the spectral radius) and `-lambda1`. Something special needs to be done in that case. One simple option is to consider the square of the adjacency matrix (reduced to one of the bi partition). I will try to provide a fix.

### comment:3 Changed 4 years ago by vdelecroix

Though, in the case of trees as requested on ask I am sure there are much more clever things to do.

### comment:4 Changed 4 years ago by vdelecroix

More precisely, here is the problem. Iterating the incidence matrix gets you into a cycle of length two

```sage: g = graphs.CompleteBipartiteGraph(1,3)
sage: v = vector([1,2,3,4])
sage: v = A*v; v/= v.norm(); v
(0.9819805060619657, 0.1091089451179962, 0.1091089451179962, 0.1091089451179962)
sage: v = A*v; v/= v.norm(); v
(0.18898223650461365, 0.566946709513841, 0.566946709513841, 0.566946709513841)
sage: v = A*v; v/= v.norm(); v
(0.9819805060619656, 0.10910894511799618, 0.10910894511799618, 0.10910894511799618)
sage: v = A*v; v/= v.norm(); v
(0.18898223650461363, 0.5669467095138409, 0.5669467095138409, 0.5669467095138409)
```

But with `A^2` everything is fine

```sage: w = A**2 * v
sage: [w[i] / v[i] for i in range(4)]
[3.0, 3.0, 3.0, 3.0]
```

(and the eigenvalue is indeed sqrt(3) for the initial graph).

### comment:5 follow-up: ↓ 7 Changed 4 years ago by vdelecroix

@dcoudert: do you know why `is_bipartite` is not there on directed graph? The notion of bipartiteness does not care about orientation at all.

### comment:6 Changed 4 years ago by vdelecroix

• Authors set to Vincent Delecroix
• Branch set to u/vdelecroix/26148
• Commit set to cb000cd3bd5b0ca05cbec61b9d3ce6f7a8e335ca
• Status changed from new to needs_review

New commits:

 ​cb000cd `26148: fix spectral radius for bipartite graphs`

### comment:7 in reply to: ↑ 5 Changed 4 years ago by dcoudert

@dcoudert: do you know why `is_bipartite` is not there on directed graph? The notion of bipartiteness does not care about orientation at all.

It's certainly because noone asked for this feature yet?

We can move the method to `generic_graph.py`. The part that tests bipartiteness is the same. However, we must update the method for finding a cycle when the digraph is not bipartite as shortest path methods care about orientation.

### comment:8 follow-up: ↓ 9 Changed 4 years ago by dcoudert

The test `if not G` treats the empty graph only, so it's not appropriate for graph without edges. The correct test is `if not G.size()`.

The current behavior is that a graph without edges is considered bipartite. I agree with that. For the empty graph, it's always ambiguous in definitions...

Concerning the method you propose for bipartite graphs in fact you use a method similar to `project_left` in `BipartiteGraph`, except that the projection don't accept loops/multiedges.

### comment:9 in reply to: ↑ 8 Changed 4 years ago by vdelecroix

The test `if not G` treats the empty graph only, so it's not appropriate for graph without edges. The correct test is `if not G.size()`.

The current behavior is that a graph without edges is considered bipartite. I agree with that. For the empty graph, it's always ambiguous in definitions...

Concerning the method you propose for bipartite graphs in fact you use a method similar to `project_left` in `BipartiteGraph`, except that the projection don't accept loops/multiedges.

Indeed. Note that `project_left` also produces a graph and not a digraph. Why on earth is this called `project_left`?

### comment:10 Changed 4 years ago by dcoudert

Well, there is also method `project_right` ;)

### comment:11 Changed 4 years ago by git

• Commit changed from cb000cd3bd5b0ca05cbec61b9d3ce6f7a8e335ca to c4aeaca9d4b88c1ced0beb1571d4e1abff339eda

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

 ​5544dda `26148: move is_bipartite to generic_graph.py` ​c4aeaca `26148: use is_bipartite when g is oriented`

### comment:12 Changed 4 years ago by vdelecroix

I removed the call to shortest path.

### comment:13 Changed 4 years ago by dcoudert

• Reviewers set to David Coudert
• Why not writing directly `return sqrt(e_min), sqrt(e_max)` ?
• `while s is not None:` -> `while s:`
• `u_to_root` -> `v_to_root` (you start from v)
• before the loop `while u_to_root and w_to_root and u_to_root[-1] == w_to_root[-1]` you could add a comment like `# Search for the nearest common ancestor of v and w`. No obligation.
• You must add special case for the graph with a single vertex and no edge (don't know the expected result).
```sage: Graph(1).spectral_radius()
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-8-c8d7e35649a7> in <module>()

985
--> 986         e_min, e_max = spectral_radius(H, prec)
987         e_min = sqrt(e_min)
988         e_max = sqrt(e_max)

957     """
958     if not G:
--> 959         raise ValueError("empty graph")
960     if G.is_directed():
961         if not G.is_strongly_connected():

ValueError: empty graph
```
Note that `Graph([(0, 0)], loops=True).spectral_radius()` returns `(1.0, 1.0)`

### comment:14 Changed 4 years ago by git

• Commit changed from c4aeaca9d4b88c1ced0beb1571d4e1abff339eda to 7c8240d7c4c034c33dab4de89b3ee8673f56d0f0

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

 ​8150104 `26148: renaming and comments in is_bipartite` ​7c8240d `26148: fix one vertex case for spectral_radius`

### comment:15 follow-up: ↓ 16 Changed 4 years ago by vdelecroix

I did implement everything excepted

• `while s is not None:` -> `while s:`

Nope, `s = 0` or `s = False` may be valid vertices.

### comment:16 in reply to: ↑ 15 Changed 4 years ago by dcoudert

Nope, `s = 0` or `s = False` may be valid vertices.

Right :P

### comment:17 Changed 4 years ago by dcoudert

• Status changed from needs_review to positive_review

LGTM.

### comment:18 Changed 4 years ago by vbraun

• Branch changed from u/vdelecroix/26148 to 7c8240d7c4c034c33dab4de89b3ee8673f56d0f0
• Resolution set to fixed
• Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.