sagemath / sage

Main repository of SageMath
https://www.sagemath.org
Other
1.43k stars 479 forks source link

"look up" a face in the face lattice of a polyhedron #29683

Closed kliem closed 3 years ago

kliem commented 4 years ago

We implement two methods that look up a face in the face lattice of a polyhedron:

This allows an easy answer for ​https://ask.sagemath.org/question/34485/what-is-the-most-efficient-way-to-look-up-a-face-in-the-face-lattice-of-a-polyhedron/#50965

Depends on #31834

CC: @jplab @LaisRast @videlec @yuan-zhou

Component: geometry

Keywords: polyhedron, faces, meet, join

Author: Jonathan Kliem

Branch/Commit: 8190917

Reviewer: Yuan Zhou, Matthias Koeppe

Issue created by migration from https://trac.sagemath.org/ticket/29683

kliem commented 4 years ago

Branch: pubic/29683

kliem commented 4 years ago

Changed branch from pubic/29683 to public/29683

kliem commented 4 years ago

Commit: 521f9e0

kliem commented 4 years ago

Last 10 new commits:

295039adocumentation
d36da4acoverage and small improvement
2d0f0d9method `reset` for the face iterator
6d99a4ctypo
5711f6bmethod `ignore_subsets`
f6633bdmethod current and fix for reset
e8c17c3join_of_Vrep and meet_of_facets for face iterator
b371a73expose in combinatorial_polyhedron
954b5c8raise index error for index error
521f9e0expose in Polyhedron_base
kliem commented 4 years ago
comment:3

There are failing tests. I would wait until the dependices are taken care of to make things work again.

kliem commented 4 years ago

New commits:

19d2742join_of_Vrep and meet_of_facets for face iterator
c5c17ffaccount for changes in the newest develop
7a4b950expose in combinatorial_polyhedron
f108b06raise index error for index error; fix doctests
a69ba9fexpose in Polyhedron_base
a6f53affix doctests
kliem commented 4 years ago

Changed commit from 521f9e0 to a6f53af

kliem commented 4 years ago

Changed branch from public/29683 to public/29683-reb

kliem commented 4 years ago
comment:6

I want to redesign some of the setup before making everything more complicated.

More precisely, creating strucutures face_struct and faces_list_struct that take care of the details. Then this overly long argument list of get_next_level will just reduce to three arguments or so. This would also make future changes easier including the transition to bitsets.pxi.

kliem commented 3 years ago

Changed commit from a6f53af to 7e0a25e

kliem commented 3 years ago

New commits:

2a14433join_of_Vrep and meet_of_facets for face iterator
530d85cexpose in combinatorial_polyhedron
73a8be1raise index error for index error; fix doctests
c7ce411expose in Polyhedron_base
7e0a25efix doctests
kliem commented 3 years ago

Changed dependencies from #29681 to none

kliem commented 3 years ago

Changed branch from public/29683-reb to public/29683-reb2

mkoeppe commented 3 years ago
comment:8

Some quick comments:

kliem commented 3 years ago
comment:9

I discovered a bug by accident:

sage: P = polytopes.permutahedron(3, backend='cdd')                                                                                                                                 
sage: [x.ambient_Hrepresentation() for x in P.facets()]                                                                                                                             
[(An inequality (1, 1, 0) x - 3 >= 0, An inequality (1, 0, 1) x - 3 >= 0),
 (An inequality (1, 1, 0) x - 3 >= 0, An inequality (1, 0, 0) x - 1 >= 0),
 (An inequality (1, 1, 0) x - 3 >= 0, An inequality (0, 1, 1) x - 3 >= 0),
 (An inequality (1, 1, 0) x - 3 >= 0, An inequality (0, 1, 0) x - 1 >= 0),
 (An inequality (1, 1, 0) x - 3 >= 0, An inequality (0, 0, 1) x - 1 >= 0),
 (An inequality (1, 1, 0) x - 3 >= 0, An equation (1, 1, 1) x - 6 == 0)]

The problem is that the backend cdd doesn't put equations in a stable place, which my code assumed. It depends whether you initialize from Vrep or from Hrep I guess.

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

69740e7improved documentation
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from 7e0a25e to 69740e7

mkoeppe commented 3 years ago
comment:11

Setting new milestone based on a cursory review of ticket status, priority, and last modification date.

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from 69740e7 to 84d05ce

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

84d05cefix a variable name and add a doctest
kliem commented 3 years ago

Changed branch from public/29683-reb2 to public/29683-reb3

kliem commented 3 years ago

Changed commit from 84d05ce to deca8e0

kliem commented 3 years ago
comment:14

rebased


New commits:

901b789join_of_Vrep and meet_of_facets for face iterator
2a5a782expose in combinatorial_polyhedron
6cecfd7raise index error for index error; fix doctests
9a66e32expose in Polyhedron_base
977c088fix doctests
15e2c58improved documentation
deca8e0fix a variable name and add a doctest
kliem commented 3 years ago
comment:15

From #30469 comment:25

  • :meth:~sage.geometry.polyhedron.base.Polyhedron_base.meet_of_facets | largest face contained specified facets` --> contained in the specified facets
  • Polyhedron_base.join_of_Vrep(self, *Vrepresentatives): ... in case of... --> "in the case of" * 2
  • FaceIterator_base(SageObject).join_of_Vrep(self, *indices): .. NOTE:: ... may not te well defined. --> "be"
  • Why is it meet_of_facets, but not meet_of_Hrep? If rays/lines are allowed on the primal side, why doesn't it make sense to consider equations on the dual side?
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from deca8e0 to ac9ff1e

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

ac9ff1etypos
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

e9769abanother typo
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from ac9ff1e to e9769ab

kliem commented 3 years ago
comment:18

Replying to @kliem:

From #30469 comment:25

  • :meth:~sage.geometry.polyhedron.base.Polyhedron_base.meet_of_facets | largest face contained specified facets` --> contained in the specified facets
  • Polyhedron_base.join_of_Vrep(self, *Vrepresentatives): ... in case of... --> "in the case of" * 2
  • FaceIterator_base(SageObject).join_of_Vrep(self, *indices): .. NOTE:: ... may not te well defined. --> "be"

All fixed.

  • Why is it meet_of_facets, but not meet_of_Hrep? If rays/lines are allowed on the primal side, why doesn't it make sense to consider equations on the dual side?

Lines are basically ignored as are equations.

However, rays do make a difference. This is how I came up with it 12 months ago. join_of_Vrep is just better than join_of_vertices_and_rays, isn't it?

In the end, I don't know what is best.

yuan-zhou commented 3 years ago
comment:19

I'm very new to the topic. I just started learning about combinatorial polyhedra today, and haven't read all the code changes on this ticket yet. But I'm wondering if it makes more sense to start by faces of the (possibly unbounded) polyhedron as input. Let

P = Polyhedron(vertices=[[1,0], [0,1]], rays=[[1,1]])
sage: V = P.Vrepresentation()
sage: V                                                                        
(A vertex at (0, 1), A vertex at (1, 0), A ray in the direction (1, 1))

Since "A ray in the direction (1, 1)" is not a face of P, V[2] shouldn't be taken as input for the look-up. Instead, the input could be something like [(1, 2), 0], where each element (= a single index or a tuple of indices) in the list must describe an actual face of the polyhedron.

For example, we can define the functions join_of_faces and meet_of_faces with input

Just my two cents.

kliem commented 3 years ago
comment:20

As for the combinatorial polytope: It's just a coatom(ist)ic and atom(ist)ic lattice with some extra properties. Hence the name join_of and meet_of. As the lattice is coatom(ist)ic and atom(ist)ic the concept of meet_of_faces and join_of_faces is redundant. To compute the meet of faces it suffices to compute the meet of the coatoms/facets. As it is only one depth-first lookup the computation time of doing this is usually neglectable.

The unbounded polyhedron (no lines/linear subspaces assumed) is a bit different. It is basically a polytope with marked far face. The face lattice is given by [0, 1] \ [0, F]. Where [0,1] denotes the lattice of the bounded polytope and F is the far face.

The rays are vertices of that bounded polytope and thus it makes sense to include them in the notion of a join. Computing the join is always well defined with the exception for the case where you only include elements of the far face, e.g. rays.

However, I don't have a strong opinion on a name scheme. In the combinatorial polyhedron module it should be join and meet I think, but in polyhedron/base.py it might as well be look_up_face or find_face.

yuan-zhou commented 3 years ago
comment:21

Given two faces (as polyhedra), if I'm interested in finding the largest face (as polyhedron) that is contained in both faces, is it better to call intersection or to apply join_of_facets?

kliem commented 3 years ago
comment:22

If you want the largest face contained in both of them, you need to do the intersection, which translates to the meet operation in the face lattice.

I would propose not defining the meet of the faces but instead the meet of the coatoms/facets. If you give me the faces as polyhedron, I need to first determine how your facet-description of the facets corresponds to the facets of the original polyhedron.

Meet should be much faster, because all the information is already encoded in the combinatorial polyhedron. It is just a simple depth first search in a lattice, doing a couple of intersections and subset checks of bitsets.

But it also depends on the kind of output you want.

kliem commented 3 years ago
comment:23

Actually what you describe is really easy to do, but neither implemented nor intended for in this ticket:

In your case you know (at least in the bounded case) that the intersection of both polytopes is exactly the convex hull of the intersection of the vertices. That is of course by far not true in the general case.

To determine the equations is really easy. It's just computing the kernel of a matrix. The possibly redundant inequalities are given by the union of the inequalities. So the only thing left to do (at least with backend field) is to figure out, which of the inequalities is redundant.

But this is exactly the functionality that face_as_combinatorial_polyhedron tries to implement as well.

yuan-zhou commented 3 years ago
comment:24

Thanks!

Another possibly irrelevant question: for an empty face, does its ambient_H_indices or ambient_Hrepresentation mean anything?

P =  Polyhedron(rays=[(1,0),(0,1)])
F = P.faces(-1)[0]
F.ambient_H_indices()  #  returns  (0, 1)
F.ambient_Hrepresentation() #  (An inequality (1, 0) x + 0 >= 0, An inequality (0, 1) x + 0 >= 0)
P.faces(0)[0].ambient_Hrepresentation() # same output as above
kliem commented 3 years ago
comment:25

Stupid corner cases.

Yes, usually it ambient_H_indices or ambient_Hrepresentation does mean something for the empty face as well. But not in this case. There is at least two different definitions for faces:

Usually they agree. But not in some case. In your example, there are two inequalities that are not redundant. They define those two unbounded 1-dimension faces, the 2-dimensional face is the empty intersection and the vertex. The empty face cannot be represented by those inequalities.

Different story for the second definition.

If you have at least two vertices, then things are fine. Sometimes cones are so much nicer...

That is why the CombinatorialPolyhedron itself just ignores the empty face and the full-dimensional face for the face iterator.

yuan-zhou commented 3 years ago
comment:26

I'm reading base.py and find the following corner case:

sage: P = Polyhedron(vertices=[[1,0,0], [0,1,0]], rays=[[1,1,0]])               
sage: Q = Polyhedron(vertices=[[0,0],[0,1],[1,1],[1,0]])
sage: a, b, c = Q.Hrepresentation()[:3]
sage: P.meet_of_facets(b, c)
sage: P.meet_of_facets(a, b)

That is, join_of_Vrep and meet_of_facets probably need to check if v and facet given as PolyhedronFace or V/H-representation actually correspond to self.

Some typos:

yuan-zhou commented 3 years ago
comment:27

In /src/sage/geometry/polyhedron/combinatorial_polyhedron/base.pyx, do the methods join_of_Vrep and meet_of_facets of CombinatorialPolyhedron have input indices corresponding to Vrepresentation/Hrepresentation or .vertices()/.facets()?

I guess it is the latter, but the documentation doesn't say. (And the name .._Vrepr somehow suggests the former?)

I also got a SIGSEGV when trying the following

sage: Q = Polyhedron(lines=[[1,1]])                                             
sage: CQ = CombinatorialPolyhedron(Q)
sage: CQ.join_of_Vrep(0)
kliem commented 3 years ago
comment:28

Replying to @yuan-zhou:

In /src/sage/geometry/polyhedron/combinatorial_polyhedron/base.pyx, do the methods join_of_Vrep and meet_of_facets of CombinatorialPolyhedron have input indices corresponding to Vrepresentation/Hrepresentation or .vertices()/.facets()?

Same thing.

I guess it is the latter, but the documentation doesn't say. (And the name .._Vrepr somehow suggests the former?)

I also got a SIGSEGV when trying the following

sage: Q = Polyhedron(lines=[[1,1]])                                             
sage: CQ = CombinatorialPolyhedron(Q)
sage: CQ.join_of_Vrep(0)

Bad.

kliem commented 3 years ago
comment:29

Replying to @yuan-zhou:

I'm reading base.py and find the following corner case:

sage: P = Polyhedron(vertices=[[1,0,0], [0,1,0]], rays=[[1,1,0]])               
sage: Q = Polyhedron(vertices=[[0,0],[0,1],[1,1],[1,0]])
sage: a, b, c = Q.Hrepresentation()[:3]
sage: P.meet_of_facets(b, c)
sage: P.meet_of_facets(a, b)

That is, join_of_Vrep and meet_of_facets probably need to check if v and facet given as PolyhedronFace or V/H-representation actually correspond to self.

Some typos:

  • Line 6960: empty line but ::

This is common practice in sage. I don't know if it is good practice. It just means that this is another example explained by the previous comment.

  • Line 7025: maybe you meant that the index cannot correspond to an equation in the Hrepresentation.
  • Line 7030: ValueError: 0 is not a facet --> An equation (1, 1, 1, 1, 1) x - 15 == 0 is not a facet?

Everything else will be fixed with the next commit.

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

bec98dfimproved documentation and check for elements to be of the correct polyhedron
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from e9769ab to bec98df

yuan-zhou commented 3 years ago
comment:31

Replying to @kliem:

Replying to @yuan-zhou:

In /src/sage/geometry/polyhedron/combinatorial_polyhedron/base.pyx, do the methods join_of_Vrep and meet_of_facets of CombinatorialPolyhedron have input indices corresponding to Vrepresentation/Hrepresentation or .vertices()/.facets()?

Same thing.

They are different in some cases, such as

sage: CQ.Hrepresentation()                                                      
(An equation (1, -1) x + 0 == 0,)
sage: CQ.facets()                                                               
()
kliem commented 3 years ago
comment:32

Ok. Forgot that.

yuan-zhou commented 3 years ago
comment:33

This is certainly detail, but it bothers me that Polyhedron.join_of_Vrep can take a line as Vrep (it simply ignores the line) whereas Polyhedron.meet_of_facets would raise an error for equation or index corresponding to an equation. It is not symmetric. Also, I would give facet.dim() + 1 == self.dim() a separate error message.

kliem commented 3 years ago
comment:34

Replying to @yuan-zhou:

Replying to @kliem:

Replying to @yuan-zhou:

In /src/sage/geometry/polyhedron/combinatorial_polyhedron/base.pyx, do the methods join_of_Vrep and meet_of_facets of CombinatorialPolyhedron have input indices corresponding to Vrepresentation/Hrepresentation or .vertices()/.facets()?

Same thing.

They are different in some cases, such as

sage: CQ.Hrepresentation()                                                      
(An equation (1, -1) x + 0 == 0,)
sage: CQ.facets()                                                               
()

The whole equation thing is so messed up, I don't even know where to start.

It was never a good decision to permit those in the first place.

Sometimes it's first the equations and then the inequalities, sometimes it's the other way around.

Maybe have same consistently last throughout the module? What do you think.

I have the feeling that there are some really bad design decisions that are haunting me.

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from bec98df to bfb08f5

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

bfb08f5typo
kliem commented 3 years ago
comment:36

I mean the ambient_H_indices of faces are not matching up the Hrepresentation. That is one thing that is annoying.

kliem commented 3 years ago
comment:37

Replying to @yuan-zhou:

This is certainly detail, but it bothers me that Polyhedron.join_of_Vrep can take a line as Vrep (it simply ignores the line) whereas Polyhedron.meet_of_facets would raise an error for equation or index corresponding to an equation. It is not symmetric. Also, I would give facet.dim() + 1 == self.dim() a separate error message.

As stated above (at least I think I said this). It don't have a strong opinion on this. So you are proposing to just ignore equations? Or raise an error in both cases?