Opened 5 years ago

Closed 5 years ago

# Suggested improvement to computing Ihara zeta functions

Reported by: Owned by: kimball minor sage-6.9 graph theory Ihara zeta function chapoton, ncohen Kimball Martin, Frédéric Chapoton Nathann Cohen, Frédéric Chapoton N/A bdd724c (Commits) bdd724c11d70f09ddabdc6d0e7161a27e6c0d8cc

### Description

I needed to compute a lot of Ihara zeta functions awhile ago, and the built-in function ihara_zeta_function_inverse() was far too slow for my needs. I suppose the current function uses something like the Bass determinant formula, though I haven't checked the code. My approach is the following:

1. Given a graph G, "prune" it, i.e., successively delete all vertices of degree < 2 until one is left with a minimum degree (>=) 2 graph or the null graph (no vertices). Call the pruned graph H. Then G and H will have the same zeta function.
1. Compute the (Hashimoto) edge matrix T of H.
1. Return the reverse (in the sense of Sage's polynomial reverse function) of the characteristic polynomial of T. This is the reciprocal of the Ihara zeta function of H, and therefore G, by the Hashimoto determinant formula.

I did various personal testing and my code seems much faster than the built-in function almost all of the time. I would like to suggest some form of this for a future release of Sage.

Here are some sample runs from Sage 6.3 (my function being ihara_inverse()).

A sparse graph:

```sage: G = graphs.CycleGraph(50)
sage: time(G.ihara_zeta_function_inverse())
CPU times: user 6.93 s, sys: 2.91 ms, total: 6.93 s
Wall time: 6.94 s
t^100 - 2*t^50 + 1
sage: time(ihara_inverse(G))
CPU times: user 30.6 ms, sys: 3 ms, total: 33.6 ms
Wall time: 31.7 ms
t^100 - 2*t^50 + 1
```

Some non-sparse graphs:

```sage: G = graphs.CompleteBipartiteGraph(5,5)
sage: time(G.ihara_zeta_function_inverse())
CPU times: user 947 ms, sys: 1.55 ms, total: 948 ms
Wall time: 958 ms
-1048576*t^50 + 14745600*t^48 - 95027200*t^46 + 369868800*t^44 - 962560000*t^42 + 1744135680*t^40 - 2201382400*t^38 + 1831972800*t^36 - 787721200*t^34 - 154949375*t^32 + 428125040*t^30 - 204331400*t^28 - 24665200*t^26 + 60598300*t^24 - 13690000*t^22 - 7728440*t^20 + 3527600*t^18 + 564550*t^16 - 433200*t^14 - 35000*t^12 + 31920*t^10 + 3100*t^8 - 1200*t^6 - 200*t^4 + 1
sage: time(ihara_inverse(G))
CPU times: user 13.1 ms, sys: 3.26 ms, total: 16.4 ms
Wall time: 38.8 ms
-1048576*t^50 + 14745600*t^48 - 95027200*t^46 + 369868800*t^44 - 962560000*t^42 + 1744135680*t^40 - 2201382400*t^38 + 1831972800*t^36 - 787721200*t^34 - 154949375*t^32 + 428125040*t^30 - 204331400*t^28 - 24665200*t^26 + 60598300*t^24 - 13690000*t^22 - 7728440*t^20 + 3527600*t^18 + 564550*t^16 - 433200*t^14 - 35000*t^12 + 31920*t^10 + 3100*t^8 - 1200*t^6 - 200*t^4 + 1
```
```sage: G = graphs.CompleteBipartiteGraph(6,6)
sage: time(G.ihara_zeta_function_inverse())
CPU times: user 5.25 s, sys: 5.38 ms, total: 5.25 s
Wall time: 5.27 s
244140625*t^72 - 5625000000*t^70 + 61699218750*t^68 - 428250000000*t^66 + 2108422265625*t^64 - 7820550000000*t^62 + 22648537650000*t^60 - 52343499600000*t^58 + 97765790872500*t^56 - 148342725328000*t^54 + 182431141590600*t^52 - 179637252149376*t^50 + 137507735376900*t^48 - 76129085692800*t^46 + 23737562794800*t^44 + 3520304265600*t^42 - 8655860024370*t^40 + 4403201990400*t^38 - 494274924300*t^36 - 606594787200*t^34 + 322933258350*t^32 - 16849733760*t^30 - 39070587600*t^28 + 11389507200*t^26 + 1979295300*t^24 - 1398038400*t^22 + 1778760*t^20 + 102646400*t^18 - 5820300*t^16 - 5443200*t^14 + 267600*t^12 + 213120*t^10 + 2025*t^8 - 4800*t^6 - 450*t^4 + 1
sage: time(ihara_inverse(G))
CPU times: user 20.9 ms, sys: 4.92 ms, total: 25.8 ms
Wall time: 22.2 ms
244140625*t^72 - 5625000000*t^70 + 61699218750*t^68 - 428250000000*t^66 + 2108422265625*t^64 - 7820550000000*t^62 + 22648537650000*t^60 - 52343499600000*t^58 + 97765790872500*t^56 - 148342725328000*t^54 + 182431141590600*t^52 - 179637252149376*t^50 + 137507735376900*t^48 - 76129085692800*t^46 + 23737562794800*t^44 + 3520304265600*t^42 - 8655860024370*t^40 + 4403201990400*t^38 - 494274924300*t^36 - 606594787200*t^34 + 322933258350*t^32 - 16849733760*t^30 - 39070587600*t^28 + 11389507200*t^26 + 1979295300*t^24 - 1398038400*t^22 + 1778760*t^20 + 102646400*t^18 - 5820300*t^16 - 5443200*t^14 + 267600*t^12 + 213120*t^10 + 2025*t^8 - 4800*t^6 - 450*t^4 + 1
```
```sage: G = graphs.CompleteBipartiteGraph(7,6)
sage: time(G.ihara_zeta_function_inverse())
CPU times: user 11 s, sys: 4.44 ms, total: 11.1 s
Wall time: 11.1 s
-3645000000*t^84 + 102060000000*t^82 - 1373472450000*t^80 + 11821912200000*t^78 - 73058147977500*t^76 + 344926224417600*t^74 - 1292328995339375*t^72 + 3939476291089776*t^70 - 9936410853590550*t^68 + 20970989057296320*t^66 - 37290794329153575*t^64 + 56044177460026560*t^62 - 71135923125029280*t^60 + 75847829145762240*t^58 - 67114948712186700*t^56 + 48106249389322880*t^54 - 26524408167721800*t^52 + 9751404304200384*t^50 - 849449754020700*t^48 - 1689661308484800*t^46 + 1245020652426600*t^44 - 369917663198400*t^42 - 35867256926130*t^40 + 75387616466400*t^38 - 24087116417700*t^36 - 1812737707200*t^34 + 3444792879150*t^32 - 682752530880*t^30 - 191262470400*t^28 + 95440161600*t^26 + 997964100*t^24 - 6891091200*t^22 + 551436984*t^20 + 347144000*t^18 - 39257820*t^16 - 13910400*t^14 + 1241940*t^12 + 443520*t^10 - 6615*t^8 - 8400*t^6 - 630*t^4 + 1
sage: time(ihara_inverse(G))
CPU times: user 27.4 ms, sys: 6.38 ms, total: 33.7 ms
Wall time: 29 ms
-3645000000*t^84 + 102060000000*t^82 - 1373472450000*t^80 + 11821912200000*t^78 - 73058147977500*t^76 + 344926224417600*t^74 - 1292328995339375*t^72 + 3939476291089776*t^70 - 9936410853590550*t^68 + 20970989057296320*t^66 - 37290794329153575*t^64 + 56044177460026560*t^62 - 71135923125029280*t^60 + 75847829145762240*t^58 - 67114948712186700*t^56 + 48106249389322880*t^54 - 26524408167721800*t^52 + 9751404304200384*t^50 - 849449754020700*t^48 - 1689661308484800*t^46 + 1245020652426600*t^44 - 369917663198400*t^42 - 35867256926130*t^40 + 75387616466400*t^38 - 24087116417700*t^36 - 1812737707200*t^34 + 3444792879150*t^32 - 682752530880*t^30 - 191262470400*t^28 + 95440161600*t^26 + 997964100*t^24 - 6891091200*t^22 + 551436984*t^20 + 347144000*t^18 - 39257820*t^16 - 13910400*t^14 + 1241940*t^12 + 443520*t^10 - 6615*t^8 - 8400*t^6 - 630*t^4 + 1
```

My code is attached.

### comment:1 follow-up: ↓ 2 Changed 5 years ago by chapoton

Thanks a lot!

Indeed, this seems to be much faster. Do you want me to take on, or do you want to go on doing (part of) the job yourself ?

First step would be to rewrite your code as a single function, by inclusion of your two auxiliary functions.

Next step would be to replace the code of the current ihara_zeta_function_inverse by this single function code of yours. For this you need to have a source install of sage, and to have set up git as explained in the develop manual.

I am ready to help, but this nevertheless will require some effort from you.

### comment:2 in reply to: ↑ 1 Changed 5 years ago by kimball

I would be happy if you started taking over, and then let me know what you need from me.

### comment:3 Changed 5 years ago by ncohen

Just a note:

If you want to remove vertices of degree `<2` from your graph until none is left, you can do this:

```g.subgraph(vertices=g.cores(k=2))
```

Nathann

### comment:4 Changed 5 years ago by chapoton

• Branch set to u/chapoton/19118

New commits:

 ​83ad065 `trac #19118 new and faster algo for Ihara zeta function of graphs`

### comment:5 Changed 5 years ago by git

• Commit changed from 83ad065e0d3d5e211d91705dfe5bc5217a5d2938 to 200e8f2460981464b4ff9011758aca4d0d77115f

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

 ​200e8f2 `trac #19118 simplify core construction`

### comment:6 Changed 5 years ago by git

• Commit changed from 200e8f2460981464b4ff9011758aca4d0d77115f to 5aa4d7ab554983e1e830b7eda7b467177e04d27c

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

 ​5aa4d7a `trac #19118 fixing mistakes`

### comment:7 Changed 5 years ago by ncohen

Also, I am not totally sure that I understand the definition of the Hashimoto matrix. To me it sounds like `g.line_graph().adjacency_matrix()` (i.e. the adjacency matrix of the line graph), but I am not totally sure `O_o`

It sounds a bit weird to me that this matrix seems to depend on whether an edge `(i,j)` is returned as `i,j` or as `j,i`.

Nathann

### comment:8 Changed 5 years ago by chapoton

Yes, maybe the author of this code can provide us with explanations and references ?

I guess that one needs to pick an arbitrary orientation of the graph first. Then compute this matrix, and the resulting charpoly will not depend on the orientation choice.

### comment:9 Changed 5 years ago by git

• Commit changed from 5aa4d7ab554983e1e830b7eda7b467177e04d27c to 2372bee1f8d3cc9ae56e11fa6b93489726852d4a

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

 ​55b1244 `Merge branch 'u/chapoton/19118' into 6.9.b6` ​2372bee `trac #19118 final clean-up, removal of old algorithm`

### comment:10 Changed 5 years ago by chapoton

• Authors set to Kimball Martin, Frédéric Chapoton
• Status changed from new to needs_review

This looks good to me, and is really much more efficient. Nathann, do you approve ?

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

• Commit changed from 2372bee1f8d3cc9ae56e11fa6b93489726852d4a to bdd724c11d70f09ddabdc6d0e7161a27e6c0d8cc

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

 ​9fbd65d `Merge branch 'u/chapoton/19118' into 6.9.b7` ​bdd724c `trac #19118 refiewer's comment : elif instead of if`

### comment:12 Changed 5 years ago by ncohen

I have nothing to add to the algorithn. If you are all okay with the mathematical side, then it can go.

Nathann

### comment:13 Changed 5 years ago by chapoton

• Reviewers set to Nathann Cohen, Frédéric Chapoton
• Status changed from needs_review to positive_review

ok, then let it be.

Thanks Nathann

### comment:14 Changed 5 years ago by vbraun

• Branch changed from u/chapoton/19118 to bdd724c11d70f09ddabdc6d0e7161a27e6c0d8cc
• Resolution set to fixed
• Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.