Opened 7 years ago

Closed 6 years ago

# permanental_minor_vector, matching polynomial

Reported by: Owned by: pernici major sage-6.4 linear algebra ncohen, jsp Mario Pernici Vincent Delecroix N/A 901d36c (Commits) 901d36c702cfab6c13e1d14083cef7e8ee931c37

Added function `permanental_minor_polynomial`, which computes the polynomial of the sums of the permanental minors in an efficient way, see arXiv:1406.5337.

As a benchmark, consider dancing.sage. On my computer `dance(10)` took 3 hours using sage-6.2; it takes 0.08 seconds with this patch, using the Butera-Pernici algorithm, 13 seconds using the Godsil algorithm.

In the case of banded matrices, the speedup is greater, since for fixed bandwidth the new algorithm is polynomial; here is a permanent which cannot be computed with sage-6.2

```sage: from sage.matrix.matrix_misc import permanental_minor_polynomial
sage: n, w = 100, 5
sage: m = matrix([[i*j + 1 if abs(i-j) <= w else 0
....:              for i in range(n)] for j in range(n)])
sage: time p = permanental_minor_polynomial(m, permanent_only=True)
Wall time: 773 ms
```

This can be used to compute the matching polynomial of bipartite graphs much faster than matching_polynomial for certain graphs. For example, using the attached file to compute the reduced adjacency matrix of a 2d grid (n1, n2):

```sage: from grid import sq_mat
sage: from sage.matrix.matrix_misc import permanental_minor_polynomial
sage: m = matrix(sq_mat(6, 6))
sage: time a = permanental_minor_polynomial(m)
Wall time: 10.6 ms

sage: g = BipartiteGraph(m)
sage: time p = g.matching_polynomial()
Wall time: 46.9 s
```

For n1,n2=10,10 permanental_minor_vector(m) takes 0.3s For n1,n2=50,6 it takes 0.5s

### comment:1 Changed 7 years ago by pernici

• Authors set to Mario Pernici

### comment:2 Changed 7 years ago by pernici

• Branch set to u/pernici/ticket/16603
• Created changed from 07/01/14 18:24:37 to 07/01/14 18:24:37
• Modified changed from 07/01/14 18:29:25 to 07/01/14 18:29:25

### comment:3 Changed 7 years ago by git

• Commit set to b0d714572f15dd915b1caabbb11f5f96f4020c85

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

 ​b0d7145 `fixed bug in `permanental_minor_vector` in case of vanishing permanent`

### comment:4 Changed 7 years ago by pernici

• Description modified (diff)

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

• Description modified (diff)

### comment:6 Changed 7 years ago by pernici

• Description modified (diff)
• Summary changed from permanental_minor_vector to permanental_minor_vector, matching polynomial

### comment:7 Changed 6 years ago by vbraun_spam

• Milestone changed from sage-6.3 to sage-6.4

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

• Branch changed from u/pernici/ticket/16603 to public/ticket/16603
• Commit changed from b0d714572f15dd915b1caabbb11f5f96f4020c85 to f65acdedf20554ed8cbd88fde232da4229d96f85

I have had a quick look and made a review patch. You need to document the auxiliary function.

New commits:

 ​55d328f `Merge branch 'u/pernici/ticket/16603' in 6.4.rc0` ​f65acde `trac #16603 first review commit`

### comment:9 Changed 6 years ago by pernici

• Branch changed from public/ticket/16603 to u/pernici/ticket/16603
• Modified changed from 10/31/14 19:35:58 to 10/31/14 19:35:58

### comment:10 Changed 6 years ago by pernici

• Commit changed from f65acdedf20554ed8cbd88fde232da4229d96f85 to 947ddd16da731019f25a951e43bf551ad997a700

I made a commit in which I documented the auxiliary function, but by mistake I committed it to u/pernici/ticket/16603 instead of public/ticket/16603

New commits:

 ​947ddd1 `Added missing documentation in ``_prm_mul``;`

### comment:11 follow-up: ↓ 12 Changed 6 years ago by chapoton

Do you want this to be reviewed ? If yes, please turn it into "needs review".

### comment:12 in reply to: ↑ 11 Changed 6 years ago by pernici

Do you want this to be reviewed ? If yes, please turn it into "needs review".

How do I turn it into "needs review" ?

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

• Status changed from new to needs_review

Click on "Modify ticket", then on the button "needs review". I did it for you.

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

• Status changed from needs_review to needs_work

Hello,

This code is very cool! You can have a look at my branch `u/vdelecroix/16603` where I provide two commits that I describe below.

• I edited a bit the code, I hope it is cleaner (and actually also a bit faster)
• Please, could you fill the section `ALGORITHM` of `permanental_minor_vector`. It would be nice to have here a rough description of the algorithm.
• I edited a lot the documentation. First of all the file `matrix_misc.py` was not part of the matrix documentation! After my commits it is. If you compile the documentation (with either `make doc` or `sage -docbuild reference html`) then you will see the result in the section `Matrices and Spaces of Matrices` -> `Miscellaneous matrix functions`. I changed the name `_prm_mul` to `prm_mul` in order to make it appear in the documentation.
• Currently the situtation for benchmarks is unfair:
• the code of Jaap Spies to compute the permanent and the permanent minors can be optimized a lot!
• your code is in Python while Jaap's one is in Cython
• As you mentioned this algorithm can be used in several contexts. It would be nice to have some methods (of graphs and matrices) that call it. But if you prefer, we can do it in a further ticket.

Vincent

### comment:15 follow-up: ↓ 16 Changed 6 years ago by pernici

Hello,

how can I access your branch u/vdelecroix/16603 ?

### comment:16 in reply to: ↑ 15 Changed 6 years ago by vdelecroix

Hello,

how can I access your branch u/vdelecroix/16603 ?

How do you do with your branch ` u/pernici/ticket/16603`? All users have write access to `public/*` and `u/username/*` and read access everywhere. Just do in the Sage root:

```\$ git fetch trac u/vdelecroix/16603
\$ git checkout -b 16603_vincent_branch FETCH_HEAD
```

The first command download my two commits. The special temporary branch `FETCH_HEAD` is a pointer to the last commit. The second command creates a new branch called `16603_vincent_branch` and move to it.

You can then merge these commits (if you like them) in your branch with

```\$ git checkout my_original_branch
\$ git merge 16603_vincent_branch
```

Is that clear enough?

Vincent

### comment:17 follow-up: ↓ 18 Changed 6 years ago by pernici

Previously I used ./sage -dev checkout --ticket 16603 Following your instructions I got your branch; thanks.

### comment:18 in reply to: ↑ 17 Changed 6 years ago by vdelecroix

Previously I used ./sage -dev checkout --ticket 16603 Following your instructions I got your branch; thanks.

I see. The `sage -dev` commands are a bit outdated (and give you much less freedom). You might better use `git` or the enhanced `git trac` that Volker Braun wrote. Have a look at the documentation.

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

• Description modified (diff)

I rewrote the description (formatting the code).

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

• Commit changed from 947ddd16da731019f25a951e43bf551ad997a700 to f511623cce7bde7f022fa212a505526b3b62f0db

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

 ​a6fd53b `trac #16603: speedup + doc` ​f611f58 `trac #16603: add matrix_misc.py to the doc` ​f511623 `Added documentation to `permanental_minor_vector``

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

Hello,

I will push several commits in a minute on my branch `u/vdelecroix/16603`. Again you are free to use/discuss/reject them. I discuss the most important changes below.

Change in your original code and documentation:

• I changed your initial function into `permanental_minor_polynomial`. It looks more natural as it is what you compute... I moved the part of the code that makes it a list in `matrix.rook_vector`
• I added a keyword `var` (as it is the case for `matrix.charpoly` for example) and change the argument name from `m` to `A` (since A is used in the docstring)
• I optimized a bit the method `prm_mul` (it was possible to remove the loop computing `exp`) and introduced the list `vars_to_do` inplace of the set `done_vars`.
• Sometimes it might happen that some coefficient of the partial product become zero. I added a test in `prm_mul`. If you think this is not worth it, I can revert that change.
• Again, I edited a lot the documentation of `permanental_minor_polynomial`. Please reread carefully, I tried to make it more understandable but I might have made many mistakes.
• Do you know an adaptation of your method to compute only a specific permanental minor? (I guess that to compute the `k`-th permanental minor it is just enough to ignore the terms `t^l^` with `l > k`)

Methods of matrices:

• does your algorithm have a name? I used Butera-Pernici (see the item just below, it has some importance).
• To allow the computation using your algorith, I made changes in `matrix.permanental_minor` and `matrix.permanent`. After my commits you can do
```sage: m = matrix(ZZ,4,4,range(16))
sage: m.permanent(algorithm="Ryser")
22432
sage: m.permanent(algorithm="ButeraPernici")
22432
```
it will be easier to make benchmarks (but I repeat that right now it is not fair for both, see also #17457).
• Right now, in all methods except `matrix.rook_vector`, the default choice for algorithm is Ryser (what was before). Where do you think we should put yours in first?
• why `matrix.rook_vector` is only defined for {0,1} matrices... it is well defined for any matrix, isn't it?
• possibly we can add a method `matrix.permanental_minor_polynomial` but it will overlap a lot with rook vector. I am not sure what to do as these methods are very specific to matrix combinatorics. Perhaps you have a better overview than me on that.

Despite this big list of remarks, I am sure this is not far for being ready for integration! This is a really nice piece of work.

Vincent

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

• Status changed from needs_work to needs_info

changes pushed on `u/vdelecroix/16603`.

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

• Commit changed from f511623cce7bde7f022fa212a505526b3b62f0db to 19350a4574f4caa2e0b5ad08f66e51930b647a42

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

 ​97f2062 `trac #16603: documentation + argument modifs` ​63e37cb `trac #16603: various optimizations` ​8a228d8 `trac #16603: modify methods of matrices` ​c76177e `Added `prec` parameter to `permanental_minor_polynomial`;` ​6f2c1c7 `used default "ButeraPernici" algorithm in `rook_vector`` ​19350a4 `Added "Godsil" algorithm to `rook_vector`.`

### comment:24 Changed 6 years ago by pernici

Hello Vincent,

Do you know an adaptation of your method to compute only a specific permanental minor? (I guess that to compute the k-th permanental minor it is just enough to ignore the terms tl with l > k)

I adapted the algorithm to compute fastly k-th permanental minor for small `k`; testing on `A = matrix(n, k, range(n)` using this algorithm `A.permanent_minor(k)` is faster for `2 <= k <= n-2

does your algorithm have a name? I used Butera-Pernici (see the item just below, it has some importance).

No, it did not have a name.

Right now, in all methods except matrix.rook_vector, the default choice for algorithm is Ryser (what was before). Where do you think we should put yours in first?

In the example above, for n x n matrices in `permanental_minor` the Ryser algorithm is generally faster for `k=n` (permanent), slightly faster for `k=1,n-1`; maybe one could make a default option which uses the Ryser algorithm in those cases, the new algorithm in the other cases.

why matrix.rook_vector is only defined for {0,1} matrices... it is well defined for any matrix, isn't it?

Usually `rook vector` is defined only on {0,1} matrices; I added the "Godsil" algorithm which is faster than the "Ryser" algorithm.

### comment:25 Changed 6 years ago by pernici

• Description modified (diff)

### comment:26 Changed 6 years ago by pernici

• Description modified (diff)

### comment:27 Changed 6 years ago by pernici

• Description modified (diff)

### comment:28 Changed 6 years ago by pernici

• Status changed from needs_info to needs_review

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

• Branch changed from u/pernici/ticket/16603 to u/vdelecroix/16603
• Commit changed from 19350a4574f4caa2e0b5ad08f66e51930b647a42 to 3e7297b35f94986b8e0f37f30a3b3d6ee3d3b699
• Reviewers set to Vincent Delecroix

Hello,

I really think that the documentation is half of a software: it provides specification of functions and help the user to orient himself/herself. Next time you implement a function in Sage pay attention that your new functionality is referenced from other part of the code and that it also make references to relevant functions. There is a lot of possibilities for organization of the documentation (see the section "Docstring" of the developer guide).

I was not happy with the current behavior and documentation of `rook_vector` and `permanental_minor` (these methods should work for rectangular matrices whatever their shape). I guess that this is my last commit for that ticket. Just have a look and if you agree we can set the ticket to positive review.

Next step would be optimisation of Ryser algorithm: #17457.

Vincent

New commits:

 ​3e7297b `trac #16603: documentation commit`

### comment:30 Changed 6 years ago by pernici

• Branch changed from u/vdelecroix/16603 to u/pernici/ticket/16603
• Modified changed from 12/18/14 23:23:23 to 12/18/14 23:23:23

### comment:31 Changed 6 years ago by pernici

• Commit changed from 3e7297b35f94986b8e0f37f30a3b3d6ee3d3b699 to e4359bb3d2444375eed077c72274e4f6a4ed87bc

Hello, I added the matching_polynomial method to the BipartiteGraph? class, so one can optionally use the "rook" algorithm. Now rook_vector is much faster when the matrix is made mostly of ones. I fixed a couple of bugs and added documentation.

A simple optimization I did not include here is the one in case of matrices which can be put in block form permuting rows and columns. I did not include it because it does not seem to be a commonly occurring case. The same optimization can be used in the computation of determinants of large matrices of this kind. Do you think that it is worth opening a ticket on it, or is it too obvious and useless?

Sooner or later I will open a ticket for a faster algorithm for the matching polynomial of graphs (non bipartite).

New commits:

 ​dd4636a `fixed bug in ``permanental_minor_poly``;` ​a60cc74 `added argument `complement`, eliminated argument `check` in ``rook_vector``;` ​e4359bb `added ``matching_polynomial`` method to ``BipartiteGraph`` class`

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

• Branch changed from u/pernici/ticket/16603 to u/vdelecroix/16603
• Commit changed from e4359bb3d2444375eed077c72274e4f6a4ed87bc to e03f93ae3016870383b50bf2a3abea4b5234e51b

Hello,

First of all, there is already a lot of material on this ticket. It would be good to stop here and open tickets for new features (please put me on cc on them).

I push a commit as I do not like the `complement` keyword that sounds like you want the rook vector of the complement matrix. I replaced it with `use_complement`. I know that it is different from the specification for `BipartiteGraphs` but it fits with `IndependentSets`

```sage: from sage.graphs.independent_sets import IndependentSets
sage: I = IndependentSets(graphs.GridGraph((5,5)))
sage: I.cardinality()
55447
sage: I = IndependentSets(graphs.GridGraph((5,5)), complement=True)
sage: I.cardinality()
66
```

Nevertheless, I kept the keyword `complement` to do the following as it is very useful

```sage: identity_matrix(6).rook_vector(complement=True)
[1, 30, 315, 1420, 2715, 1854, 265]
```

Vincent

New commits:

 ​e03f93a `trac #16603: change keywords for `rook_vector` + doc`

### comment:33 Changed 6 years ago by pernici

• Branch changed from u/vdelecroix/16603 to u/pernici/ticket/16603
• Modified changed from 01/09/15 09:22:05 to 01/09/15 09:22:05

### comment:34 Changed 6 years ago by pernici

• Commit changed from e03f93ae3016870383b50bf2a3abea4b5234e51b to 127418fc436cad3561c298a8f5c685c1e4b4fb81

Hello, I made one more small change. Thank you for reviewing and improving the code and documentation.

New commits:

 ​127418f `fixed a bug in ``matching_polynomial``; added test.`

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

• Status changed from needs_review to positive_review

Hello,

Let us stop there with positive review! If you want to go any further open a new ticket (and please put me in cc).

Vincent

### comment:36 Changed 6 years ago by vbraun

• Status changed from positive_review to needs_work

Doctest failure in src/sage/misc/sagedoc.py

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

• Branch changed from u/pernici/ticket/16603 to u/vdelecroix/16603
• Commit changed from 127418fc436cad3561c298a8f5c685c1e4b4fb81 to 901d36c702cfab6c13e1d14083cef7e8ee931c37
• Status changed from needs_work to needs_review

Funny! Done...

New commits:

 ​901d36c `fix doctest in misc/sagedoc.py`

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

• Status changed from needs_review to positive_review

### comment:39 Changed 6 years ago by vbraun

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