Opened 13 months ago
Closed 12 months ago
#29837 closed enhancement (fixed)
Improvement for incidence matrix of polyhedra
Reported by:  ghkliem  Owned by:  

Priority:  major  Milestone:  sage9.2 
Component:  geometry  Keywords:  polyhedra, incidence matrix 
Cc:  jipilab, ghLaisRast  Merged in:  
Authors:  Jonathan Kliem  Reviewers:  Travis Scrimshaw 
Report Upstream:  N/A  Work issues:  
Branch:  955f12b (Commits, GitHub, GitLab)  Commit:  955f12baa3ecf69aac9e7c6fb7773993f3eb9e80 
Dependencies:  #29838, #29839  Stopgaps: 
Description (last modified by )
This ticket collects some improvements for obtaining incidence matrices:
 The method
self._is_zero
is called many times and we can save plenty of time by obtaining it beforehand.  We can obtain the slack matrix by multiplying the Vrepresentation matrix with the homogenous Vrepresentation. From this we can obtain the incidence matrix (strangely this last step takes up most of the time).
Unfortunately, matrix multiplication for algebraic fields isn't very efficient. So we need to keep the old method in this case.
Before this ticket:
sage: P = polytopes.permutahedron(7) sage: %time _ = P.incidence_matrix() CPU times: user 751 ms, sys: 7.95 ms, total: 759 ms Wall time: 758 ms sage: P = polytopes.associahedron(['A', 9]) sage: %time _ = P.incidence_matrix() CPU times: user 1.94 s, sys: 4.03 ms, total: 1.94 s Wall time: 1.94 s sage: P = polytopes.hypercube(14) sage: %time _ = P.incidence_matrix() CPU times: user 670 ms, sys: 0 ns, total: 670 ms Wall time: 669 ms sage: K.<sqrt2> = QuadraticField(2) sage: P = polytopes.permutahedron(7).dilation(sqrt2) sage: %time _ = P.incidence_matrix() CPU times: user 3.82 s, sys: 0 ns, total: 3.82 s Wall time: 3.82 s
With this ticket:
sage: P = polytopes.permutahedron(7) sage: %time _ = P.incidence_matrix() CPU times: user 64.4 ms, sys: 7 µs, total: 64.4 ms Wall time: 63.6 ms sage: P = polytopes.associahedron(['A', 9]) sage: %time _ = P.incidence_matrix() CPU times: user 325 ms, sys: 12.1 ms, total: 337 ms Wall time: 335 ms sage: P = polytopes.hypercube(14) sage: %time _ = P.incidence_matrix() CPU times: user 139 ms, sys: 12 ms, total: 151 ms Wall time: 149 ms sage: K.<sqrt2> = QuadraticField(2) sage: P = polytopes.permutahedron(7).dilation(sqrt2) sage: %time _ = P.incidence_matrix() CPU times: user 3.69 s, sys: 3.69 ms, total: 3.69 s Wall time: 3.69 s
Change History (32)
comment:1 Changed 13 months ago by
 Branch set to public/29837
 Commit set to 6dbe5d24132b62c38a1e77c09a471510d9ab1776
 Status changed from new to needs_review
 Type changed from PLEASE CHANGE to enhancement
comment:2 Changed 13 months ago by
 Dependencies set to #29838
 Status changed from needs_review to needs_work
This should be rebased on #29838 and another ticket that speeds up getting the incidence matrix from the slack matrix.
comment:3 Changed 13 months ago by
 Dependencies changed from #29838 to #29838, #29839
comment:4 Changed 13 months ago by
 Dependencies changed from #29838, #29839 to #29838, #29839, #29840
comment:5 Changed 13 months ago by
 Description modified (diff)
comment:6 Changed 13 months ago by
 Branch changed from public/29837 to public/29837reb
 Commit changed from 6dbe5d24132b62c38a1e77c09a471510d9ab1776 to 89aeab39856360ffe80a8e62a0e4886aa3d8d622
 Description modified (diff)
 Status changed from needs_work to needs_review
New commits:
43e6dfb  implement _zero_matrix

05eba1e  fix bug

7b70830  rename

47a19c6  implement slack matrix

54755d7  Merge branch 'public/29838' of git://trac.sagemath.org/sage into public/29837reb

a4576d6  set up incidence/adjacency matrix with GF(2)

30c4715  Merge branch 'public/29840' of git://trac.sagemath.org/sage into public/29837reb

fcd62d2  determine is_zero beforehand

13b8583  small fix

89aeab3  obtain incidence matrix from slack matrix

comment:7 followup: ↓ 8 Changed 13 months ago by
I don't think changing the interface so that GF(2) matrices are returned is a good idea. How about leaving it as ZZ by default, and using a keyword argument for the base ring so it can be overridden.
comment:8 in reply to: ↑ 7 Changed 13 months ago by
Replying to mkoeppe:
I don't think changing the interface so that GF(2) matrices are returned is a good idea. How about leaving it as ZZ by default, and using a keyword argument for the base ring so it can be overridden.
I don't think the change to GF(2)
was needed to make things better. So I can rewrite some stuff and remove #29839, that would also do it, I think. sage.matrix.matrix_integer_dense
is probably just as fast or faster than the other one.
I really only thought that this is the more natural thing for a 0/1 matrix. If the more natural thing is base ring integer or GF(2)
really performs awful, I can just revert that decision.
comment:9 Changed 13 months ago by
See my comment in #29840
comment:10 Changed 13 months ago by
 Dependencies changed from #29838, #29839, #29840 to #29838, #29839
 Status changed from needs_review to needs_work
Changing the base ring of incidence matrices was rejected in #29340.
comment:11 followup: ↓ 12 Changed 13 months ago by
I wish this was implemented using numpy matrices (knows as arrays), which provide faster linear algebra and have builtin methods to manipulate zero patterns of their matrices.
comment:12 in reply to: ↑ 11 Changed 13 months ago by
I don't know how fast conversion back and forth works. It is true, that numpy is a bit faster, but that doesn't help me if I cannot convert things.
E.g.
sage: P = polytopes.associahedron(['A', 9]) sage: %time a = np.array(P.Vrepresentation()) CPU times: user 383 ms, sys: 0 ns, total: 383 ms Wall time: 382 ms
This is already longer than it takes to obtain the entire incidence matrix. Of course this is most likely not an ideal way to convert, but I don't know any better way at the moment.
And conversion back is also slow.
Replying to dimpase:
I wish this was implemented using numpy matrices (knows as arrays), which provide faster linear algebra and have builtin methods to manipulate zero patterns of their matrices.
comment:13 Changed 13 months ago by
 Branch changed from public/29837reb to public/29837reb2
 Commit changed from 89aeab39856360ffe80a8e62a0e4886aa3d8d622 to 3d28d09c5c74cdc8cb79f3b17ef50b1bff8b828d
 Status changed from needs_work to needs_review
comment:14 Changed 13 months ago by
 Commit changed from 3d28d09c5c74cdc8cb79f3b17ef50b1bff8b828d to 7ef6a111bf70637a01b6789819b87dde647bdc4c
Branch pushed to git repo; I updated commit sha1. New commits:
7ef6a11  use slack matrix to obtain incidence matrix quicker

comment:15 Changed 13 months ago by
 Commit changed from 7ef6a111bf70637a01b6789819b87dde647bdc4c to 759d1f20427b7a213789e1c114e427824f0884c2
comment:16 Changed 12 months ago by
 Branch changed from public/29837reb2 to public/29837reb3
 Commit changed from 759d1f20427b7a213789e1c114e427824f0884c2 to 764a99b9a8c476f1d73e120f90b68652ff230eae
Last 10 new commits:
9462376  Adding optimized get_is_zero_unsafe() to rational sparse matrices.

94febde  Custom get_is_zero_unsafe() for modn sparse matrices.

bf0352f  See if we can check trivial zeros first for symbolic matrices.

2a08abe  small mistakes

333c4ba  wrong doctest

82c6b4a  docstring

211f4f1  Merge branch 'public/29839' of git://trac.sagemath.org/sage into public/29837reb3

4c5433e  determine is_zero beforehand

27a1b90  use slack matrix to obtain incidence matrix quicker

764a99b  account for change in #29839

comment:17 Changed 12 months ago by
comment:18 Changed 12 months ago by
 Description modified (diff)
comment:19 Changed 12 months ago by
You might be able to get a little more speed out by doing
return x == 0 +return (not x)
comment:20 Changed 12 months ago by
 Reviewers set to Travis Scrimshaw
However, whether you make that change or not, you can set a positive review on my behalf.
comment:21 Changed 12 months ago by
About three percent for P = polytopes.permutahedron(7).dilation(sqrt2)
.
Thank you.
comment:22 Changed 12 months ago by
 Commit changed from 764a99b9a8c476f1d73e120f90b68652ff230eae to f70e1c5d22c28afd9d3145ce43f26a4eb8b98504
Branch pushed to git repo; I updated commit sha1. New commits:
f70e1c5  improvement suggested by reviewer

comment:23 Changed 12 months ago by
 Status changed from needs_review to positive_review
comment:24 Changed 12 months ago by
I am getting this with this ticket in sageongentoo
sage t long warnlong 101.0 /usr/lib/python3.7/sitepackages/sage/geometry/polyhedron/backend_cdd.py ********************************************************************** File "/usr/lib/python3.7/sitepackages/sage/geometry/polyhedron/backend_cdd.py", line 508, in sage.geometry.polyhedron.backend_cdd.Polyhedron_RDF_cdd.__init__ Failed example: TestSuite(p).run() Expected nothing Got: Failure in _test_an_affine_basis: Traceback (most recent call last): File "/usr/lib/python3.7/sitepackages/sage/misc/sage_unittest.py", line 296, in run test_method(tester=tester) File "/usr/lib/python3.7/sitepackages/sage/geometry/polyhedron/base.py", line 2122, in _test_an_affine_basis b = self.an_affine_basis() File "/usr/lib/python3.7/sitepackages/sage/geometry/polyhedron/base.py", line 2085, in an_affine_basis chain = self.a_maximal_chain()[1:] # we exclude the empty face File "/usr/lib/python3.7/sitepackages/sage/geometry/polyhedron/base.py", line 3098, in a_maximal_chain comb_chain = self.combinatorial_polyhedron().a_maximal_chain() File "sage/misc/cachefunc.pyx", line 2310, in sage.misc.cachefunc.CachedMethodCallerNoArgs.__call__ (/dev/shm/portage/scimathematics/sage9999/work/sage9999/srcpython3_7/build/cythonized/sage/misc/cachefunc.c:12875) self.cache = f(self._instance) File "/usr/lib/python3.7/sitepackages/sage/geometry/polyhedron/base.py", line 3362, in combinatorial_polyhedron return CombinatorialPolyhedron(self) File "sage/geometry/polyhedron/combinatorial_polyhedron/base.pyx", line 367, in sage.geometry.polyhedron.combinatorial_polyhedron.base.CombinatorialPolyhedron.__init__ (/dev/shm/portage/scimathematics/sage9999/work/sage9999/srcpython3_7/build/cythonized/sage/geometry/polyhedron/combinatorial_polyhedron/base.c:6786) data = data.incidence_matrix() File "sage/misc/cachefunc.pyx", line 2310, in sage.misc.cachefunc.CachedMethodCallerNoArgs.__call__ (/dev/shm/portage/scimathematics/sage9999/work/sage9999/srcpython3_7/build/cythonized/sage/misc/cachefunc.c:12875) self.cache = f(self._instance) File "/usr/lib/python3.7/sitepackages/sage/geometry/polyhedron/base.py", line 2749, in incidence_matrix is_zero == self._is_zero UnboundLocalError: local variable 'is_zero' referenced before assignment  The following tests failed: _test_an_affine_basis **********************************************************************
seems circular.
comment:25 Changed 12 months ago by
 Status changed from positive_review to needs_work
is_zero == self._is_zero +is_zero = self._is_zero
comment:26 Changed 12 months ago by
The problem is that usually the cdd backend figures out the incidence matrix itself for inexact polyhedra. So you don't notice, because incidence_matrix
is never run (except for the empty polyhedron, hence the error). I should add a doctest for this.
comment:27 Changed 12 months ago by
 Branch changed from public/29837reb3 to public/29837reb4
 Commit changed from f70e1c5d22c28afd9d3145ce43f26a4eb8b98504 to da9ec61583cd63e246b8bbd326c3958f4679d12d
 Status changed from needs_work to needs_review
comment:28 followup: ↓ 30 Changed 12 months ago by
Trivial tweak, but the formatting is better:
Test that this method works for inexact base ring  (`cdd` sets the cache already):: + (``cdd`` sets the cache already)::
comment:29 Changed 12 months ago by
 Commit changed from da9ec61583cd63e246b8bbd326c3958f4679d12d to 955f12baa3ecf69aac9e7c6fb7773993f3eb9e80
Branch pushed to git repo; I updated commit sha1. New commits:
955f12b  sphinx syntax

comment:30 in reply to: ↑ 28 Changed 12 months ago by
Done.
Replying to tscrim:
Trivial tweak, but the formatting is better:
Test that this method works for inexact base ring  (`cdd` sets the cache already):: + (``cdd`` sets the cache already)::
comment:32 Changed 12 months ago by
 Branch changed from public/29837reb4 to 955f12baa3ecf69aac9e7c6fb7773993f3eb9e80
 Resolution set to fixed
 Status changed from positive_review to closed
New commits:
determine is_zero beforehand
small fix
use matrix multiplication to obtain incidence matrix