Trac #33611: Improve edges/ridges for simple/simplicial polytopes

In #30040 we have improved the face iterator for simple/simplicial
polytopes, but have not applied this to our heuristic, whether we should
find edges/ridges in dual/non-dual mode.

This is noticeable for edges of the permutahedron (and theoretically
also for other edges of simple polytopes) and does not make a difference
otherwise.

Before:

{{{
sage: C = polytopes.hypercube(12).combinatorial_polyhedron()
sage: %time _ = C.edges()
CPU times: user 86.1 ms, sys: 9 µs, total: 86.1 ms
Wall time: 86 ms
sage: %time _ = C.ridges()
CPU times: user 617 µs, sys: 3 µs, total: 620 µs
Wall time: 622 µs
sage: C = polytopes.cyclic_polytope(6, 20).combinatorial_polyhedron()
sage: %time _ = C.edges()
CPU times: user 114 µs, sys: 0 ns, total: 114 µs
Wall time: 117 µs
sage: %time _ = C.ridges()
CPU times: user 1.09 ms, sys: 5 µs, total: 1.09 ms
Wall time: 1.1 ms
sage: C = polytopes.permutahedron(7).combinatorial_polyhedron()
sage: %time _ = C.edges()
CPU times: user 16.1 s, sys: 0 ns, total: 16.1 s
Wall time: 16.1 s
sage: %time _ = C.ridges()
CPU times: user 2.35 ms, sys: 0 ns, total: 2.35 ms
Wall time: 2.36 ms
}}}

After:

{{{
sage: C = polytopes.hypercube(12).combinatorial_polyhedron()
sage: %time _ = C.edges()
CPU times: user 39.5 ms, sys: 0 ns, total: 39.5 ms
Wall time: 39.4 ms
sage:  %time _ = C.ridges()
CPU times: user 725 µs, sys: 60 µs, total: 785 µs
Wall time: 791 µs
sage: C = polytopes.cyclic_polytope(6, 20).combinatorial_polyhedron()
sage: %time _ = C.edges()
CPU times: user 154 µs, sys: 12 µs, total: 166 µs
Wall time: 168 µs
sage: %time _ = C.ridges()
CPU times: user 1.16 ms, sys: 95 µs, total: 1.25 ms
Wall time: 1.26 ms
sage: C = polytopes.permutahedron(7).combinatorial_polyhedron()
sage: %time _ = C.edges()
CPU times: user 63.2 ms, sys: 0 ns, total: 63.2 ms
Wall time: 63.2 ms
sage: %time _ = C.ridges()
CPU times: user 4.56 ms, sys: 0 ns, total: 4.56 ms
Wall time: 4.57 ms
}}}

The heuristic isn't perfect, but if it fails, you can now select a
better algorithm:

{{{
sage: P = polytopes.permutahedron(7)
sage: P1 = P.stack(next(P.face_generator(1)))
sage: %time P1.vertex_adjacency_matrix()
CPU times: user 20.6 s, sys: 7.96 ms, total: 20.6 s
Wall time: 20.6 s
5041 x 5041 dense matrix over Integer Ring (use the '.str()' method to
see the entries)
sage: P1 = P.stack(next(P.face_generator(1)))
sage: %time P1.vertex_adjacency_matrix(algorithm='primal')
CPU times: user 206 ms, sys: 16 ms, total: 222 ms
Wall time: 222 ms
5041 x 5041 dense matrix over Integer Ring (use the '.str()' method to
see the entries)
}}}

URL: https://trac.sagemath.org/33611
Reported by: gh-kliem
Ticket author(s): Jonathan Kliem
Reviewer(s): Travis Scrimshaw, Yuan Zhou

src/sage/geometry/polyhedron/base3.py

@@ -781,8 +781,8 @@ def facets(self):
             return ()
         return self.faces(self.dimension()-1)
 
-    @cached_method(do_pickle=True)
-    def f_vector(self, num_threads=None, parallelization_depth=None):
+    @cached_method(do_pickle=True, key=lambda self, x, y, z: None)
+    def f_vector(self, num_threads=None, parallelization_depth=None, algorithm=None):
         r"""
         Return the f-vector.
 
@@ -794,6 +794,12 @@ def f_vector(self, num_threads=None, parallelization_depth=None):
         - ``parallelization_depth`` -- integer (optional); specify
           how deep in the lattice the parallelization is done
 
+        - ``algorithm`` -- string (optional);
+          specify whether the face generator starts with facets or vertices:
+          * ``'primal'`` -- start with the facets
+          * ``'dual'`` -- start with the vertices
+          * ``None`` -- choose automatically
+
         OUTPUT:
 
         Return a vector whose `i`-th entry is the number of
@@ -849,7 +855,7 @@ def f_vector(self, num_threads=None, parallelization_depth=None):
             sage: Q.f_vector.is_in_cache()
             True
         """
-        return self.combinatorial_polyhedron().f_vector(num_threads, parallelization_depth)
+        return self.combinatorial_polyhedron().f_vector(num_threads, parallelization_depth, algorithm=algorithm)
 
     def bounded_edges(self):
         """
@@ -879,10 +885,18 @@ def bounded_edges(self):
                 yield (obj[i], obj[j])
 
     @cached_method
-    def vertex_adjacency_matrix(self):
+    def vertex_adjacency_matrix(self, algorithm=None):
         """
         Return the binary matrix of vertex adjacencies.
 
+        INPUT:
+
+        - ``algorithm`` -- string (optional);
+          specify whether the face generator starts with facets or vertices:
+          * ``'primal'`` -- start with the facets
+          * ``'dual'`` -- start with the vertices
+          * ``None`` -- choose automatically
+
         EXAMPLES::
 
             sage: polytopes.simplex(4).vertex_adjacency_matrix()
@@ -1003,15 +1017,23 @@ def vertex_adjacency_matrix(self):
                 sage: P.adjacency_matrix().is_immutable()
                 True
         """
-        return self.combinatorial_polyhedron().vertex_adjacency_matrix()
+        return self.combinatorial_polyhedron().vertex_adjacency_matrix(algorithm=algorithm)
 
     adjacency_matrix = vertex_adjacency_matrix
 
     @cached_method
-    def facet_adjacency_matrix(self):
+    def facet_adjacency_matrix(self, algorithm=None):
         """
         Return the adjacency matrix for the facets.
 
+        INPUT:
+
+        - ``algorithm`` -- string (optional);
+          specify whether the face generator starts with facets or vertices:
+          * ``'primal'`` -- start with the facets
+          * ``'dual'`` -- start with the vertices
+          * ``None`` -- choose automatically
+
         EXAMPLES::
 
             sage: s4 = polytopes.simplex(4, project=True)
@@ -1053,7 +1075,7 @@ def facet_adjacency_matrix(self):
             [1 0 1]
             [1 1 0]
         """
-        return self.combinatorial_polyhedron().facet_adjacency_matrix()
+        return self.combinatorial_polyhedron().facet_adjacency_matrix(algorithm=algorithm)
 
     def a_maximal_chain(self):
         r"""

src/sage/geometry/polyhedron/base4.py

@@ -139,11 +139,23 @@ def vertex_facet_graph(self, labels=True):
         """
         return self.combinatorial_polyhedron().vertex_facet_graph(names=labels)
 
-    def vertex_graph(self):
+    def vertex_graph(self, **kwds):
         """
         Return a graph in which the vertices correspond to vertices
         of the polyhedron, and edges to edges.
 
+        INPUT:
+
+        - ``names`` -- boolean (default: ``True``); if ``False``,
+          then the nodes of the graph are labeld by the
+          indices of the Vrepresentation
+
+        - ``algorithm`` -- string (optional);
+          specify whether the face generator starts with facets or vertices:
+          * ``'primal'`` -- start with the facets
+          * ``'dual'`` -- start with the vertices
+          * ``None`` -- choose automatically
+
         ..NOTE::
 
             The graph of a polyhedron with lines has no vertices,
@@ -184,7 +196,7 @@ def vertex_graph(self):
             sage: polytopes.simplex(1).graph().edges()
             [(A vertex at (0, 1), A vertex at (1, 0), None)]
         """
-        return self.combinatorial_polyhedron().vertex_graph()
+        return self.combinatorial_polyhedron().vertex_graph(**kwds)
 
     graph = vertex_graph
 

src/sage/geometry/polyhedron/combinatorial_polyhedron/base.pxd

@@ -51,13 +51,14 @@ cdef class CombinatorialPolyhedron(SageObject):
     cdef ListOfFaces bitrep_facets(self)
     cdef ListOfFaces bitrep_Vrep(self)
     cdef tuple far_face_tuple(self)
+    cdef int _algorithm_to_dual(self, algorithm) except -2
 
     # Methods to obtain a different combinatorial polyhedron.
     cpdef CombinatorialPolyhedron dual(self)
     cpdef CombinatorialPolyhedron pyramid(self, new_vertex=*, new_facet=*)
 
     cdef FaceIterator _face_iter(self, bint dual, int dimension)
-    cdef int _compute_f_vector(self, size_t num_threads, size_t parallelization_depth) except -1
+    cdef int _compute_f_vector(self, size_t num_threads, size_t parallelization_depth, int dual) except -1
 
     cdef inline int _compute_edges(self, dual) except -1:
         return self._compute_edges_or_ridges(dual, True)