Opened 2 years ago

Closed 2 years ago

Last modified 2 years ago

#27534 closed enhancement (fixed)

Implement Lawrence extension for polytopes

Reported by: gh-LaisRast Owned by:
Priority: major Milestone: sage-8.8
Component: geometry Keywords: polytopes, lawrence extension
Cc: jipilab Merged in:
Authors: Laith Rastanawi Reviewers: Jonathan Kliem
Report Upstream: N/A Work issues:
Branch: bc1b9cd (Commits, GitHub, GitLab) Commit:
Dependencies: Stopgaps:

Status badges

Description (last modified by gh-LaisRast)

Adding the following methods to base.py:

  • lawrence_extension(self, v): Return the Lawrence extension of self on the point v.
  • lawrence_polytope(self): Return the Lawrence polytope of self.
  • is_lawrence_polytope(self): Return true if self is a Lawrence polytope.

for definitions, see Section 6.6 of Lectures on Polytopes, Günter M. Ziegler, [Zie2007]

Change History (22)

comment:1 Changed 2 years ago by gh-LaisRast

  • Branch set to public/27534
  • Commit set to a7849b3ecea4f36f7233940ee6d18a0af943013c

New commits:

a7849b3adding Lawrence polytope

comment:2 Changed 2 years ago by gh-LaisRast

  • Status changed from new to needs_review

comment:3 Changed 2 years ago by gh-kliem

  • Reviewers set to Jonathan Kliem
  • Status changed from needs_review to needs_work

A few small comments.

  1. Single quotes (or whatever the name of those things is) is for latex, double for code-style. There is at least one place (using self), where this is inconsistent.
  2. The Lawrence extension of P, what is P? How about The Lawrence extension of a polytope P?
  3. A comment I just received myself: "These bullet points for INPUT: do not, by general Sage convention, end in a period/full-stop."
  4. It seems to me that the requirement of being full dimensional is not necessary. Instead one can require 0 to be contained, if its not full dimensional. If a polytope P satisfies "ax = 0" for all x in P, then the lawrence extension will just satisfie "(a,0)x = 0". I.e. if all hyperplanes containing P are linear, then the Lawrence extension will be contained in the "same" hyperplanes. I propose something like:
    -        if not self.is_full_dimensional():
    -            raise NotImplementedError("`self` must be full dimensional")
    +        if not self.is_full_dimensional():
    +            if not self.contains([0]*self.ambient_dim()):
    +                raise NotImplementedError("``self`` must be full-dimensional or contain the origin, try ``self.translation``")
    

Then, one can add an example that explains how to proceed in case of a not full-dimensional polytope:

sage: P = polytopes.permutahedron(4)
sage: Q = lawrence_polytope(P)
Traceback (most recent call last):
...
NotImplementedError: ``self`` must be full-dimensional or contain the origin, try ``self.translation``
sage: T = P.translation(-vector(P.vertices()[0]))
sage: Q = lawrence_polytope(T)
sage: Q
A 26-dimensional polyhedron in ZZ^28 defined as the convex hull of 47 vertices

Maybe add another example like this to the Lawrence extension (maybe Birkhoff_polytope(3)).

  1. I find too many commands in one line hard to read:
    -            sage: P = polytopes.cube(); Q = P.lawrence_polytope(); Q.is_lawrence_polytope()
    -            True
    +            sage: P = polytopes.cube()
    +            sage: Q = P.lawrence_polytope()
    +            sage: Q.is_lawrence_polytope()
    +            True
    
  2. I would use True instead of true.
  3. You never use d = self.dim().

comment:4 Changed 2 years ago by gh-LaisRast

Thanks for your comments.

The functions need to be rewritten entirely: for instance, my function lawrence_polytope does produce a Lawrence polytope, but not the Lawrence polytope of self according to the usual definitions (there are two different definitions for it).

I will come back to this ticket later and fix everything.

comment:5 follow-ups: Changed 2 years ago by gh-kliem

One more remark.

is_lawrence_polytope is very slow in cases, where the polytope does not have few vertices. Try checking this for permutahedron(6).

Instead of computing the Gale-transform, one could just grap the information from the facets.

  1. Create a list of vertices V.
  2. Iterate through all facets F containing exactly n-1 vertices.
  3. Remove the vertex not contained in F from V (this vertex "is" the origin in the Gale transform).
  4. Iterate through all facets F containing exactly n-2 vertices.
  5. Let v_1,v_2 not be in F. If v_1 and v_2 in V, remove them.

If at the end of this V is empty, the Gale diagram is centrally symmetric. If not, then it is not.

comment:6 in reply to: ↑ 5 ; follow-ups: Changed 2 years ago by jipilab

Replying to gh-kliem:

One more remark.

is_lawrence_polytope is very slow in cases, where the polytope does not have few vertices. Try checking this for permutahedron(6).

It is hard to have an algorithm that gets faster as the input gets larger!

In spite of that, yes, it might not be optimized, but I would not go into intricate implementation of the function unless there is yet a known use case.

One suggestion:

Perles' example of a non-rational polytope uses lawrence extension to be constructed. It would be nice to add the example to the library of polytopes (either using lawrence extension or just hard-coding the vertices. Or even better, hard-coding the vertices and in the test of the function for the polytope, construct the polytope using lawrence extension and check that it is really combinatorially isomorphic.

comment:7 Changed 2 years ago by jipilab

The Title of the ticket should be more precise about what the ticket provides.

In order to make this ticket more accessible, it would be nice to have related keywords: polytopes and lawrence extension.

comment:8 in reply to: ↑ 6 Changed 2 years ago by gh-kliem

Replying to jipilab:

Replying to gh-kliem:

One more remark.

is_lawrence_polytope is very slow in cases, where the polytope does not have few vertices. Try checking this for permutahedron(6).

It is hard to have an algorithm that gets faster as the input gets larger!

Sure. But its also nice to not recalculate something that is already known.

In spite of that, yes, it might not be optimized, but I would not go into intricate implementation of the function unless there is yet a known use case.

Whether a polytope is a Lawrence polytope can simply be decided from the incidence matrix. It is simply the question, whether the vertices V can be labeled v1,...,vn,w1,...,w2m such that V\{vi} and V\{w2j,w2j+1} is a facet.

Moreover, the current version has issues. We discovered that gale_transform() is not normalized. At the moment the answer depends on the realization of the polytope. Also, the current code will return True for the bipyramid over a triangle. To my understanding this is not a Lawrence polytope as the Gale diagram is not symmetric with respect to multiplicity.

comment:9 Changed 2 years ago by gh-kliem

Is a pyramid over a Lawrence polytope a Lawrence polytope? I think the definition should allow for an arbitrary number of points at the origin of the Gale diagram.

polymake's construction of a Lawrence polytope of a polytope with a vertex at the origin will even have an odd number of points. This construction is according to Bernd Sturmfels' definition, I think.

Let V be the vertex matrix of a polytope. The Lawrence polytope will have the following vertices according to Günter Ziegler (Lectures on Polytopes):

            \begin{pmatrix}
             V      &   V    \\
             I_n    &   2*I_n
            \end{pmatrix},

According to Bernd Sturmfels (https://arxiv.org/pdf/math/0202104.pdf):

            \begin{pmatrix}
             V      &   0    \\
             I_n    &   I_n
            \end{pmatrix},

Both construction yield centrally symmetric Gale diagrams, but the second construction might not have an odd number of points, as pointed out.

Only the first construction will "add" for each point x in the Gale diagram a point -x.

Btw, the first construction will not have gale_transform be centrally symmetric, before normalizing. E.g. when applying it for the cube, those are the coordinates of the Gale transform:

sage: P = polytopes.cube()
sage: I_n = matrix.identity(8)
sage: V = P.vertices_matrix()
sage: lambda_V = block_matrix([[V,V],[I_n, 2*I_n]])
sage: Polyhedron(lambda_V.transpose()).gale_transform()

[(2, 0, 0, 0),
 (-1, 0, 0, 0),
 (0, 2, 0, 0),
 (0, -1, 0, 0),
 (0, 0, 2, 0),
 (0, 0, -1, 0),
 (-2, -2, -2, 0),
 (1, 1, 1, 0),
 (0, 0, 0, 2),
 (0, 0, 0, -1),
 (-2, -2, 0, -2),
 (1, 1, 0, 1),
 (-2, 0, -2, -2),
 (1, 0, 1, 1),
 (4, 2, 2, 2),
 (-2, -1, -1, -1)]

comment:10 Changed 2 years ago by git

  • Commit changed from a7849b3ecea4f36f7233940ee6d18a0af943013c to 69cd2c3f76444244edc5fa2d9d38001ac695ba36

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

69cd2c3improve is_lawrence and fix lawrence_polytope

comment:11 in reply to: ↑ 5 Changed 2 years ago by gh-LaisRast

  • Status changed from needs_work to needs_review

Replying to gh-kliem:

One more remark.

is_lawrence_polytope is very slow in cases, where the polytope does not have few vertices. Try checking this for permutahedron(6).

Instead of computing the Gale-transform, one could just grap the information from the facets.

  1. Create a list of vertices V.
  2. Iterate through all facets F containing exactly n-1 vertices.
  3. Remove the vertex not contained in F from V (this vertex "is" the origin in the Gale transform).
  4. Iterate through all facets F containing exactly n-2 vertices.
  5. Let v_1,v_2 not be in F. If v_1 and v_2 in V, remove them.

If at the end of this V is empty, the Gale diagram is centrally symmetric. If not, then it is not.

This is really helpful and fast. I implemented this in is_lawrence_polytope.

Both construction yield centrally symmetric Gale diagrams, but the second construction might not have an odd number of points, as pointed out.

I decided to only put the construction given by Ziegler. It works in the case where the polytope is not full-dimensional.

comment:12 in reply to: ↑ 6 Changed 2 years ago by gh-LaisRast

Replying to jipilab:

One suggestion:

Perles' example of a non-rational polytope uses lawrence extension to be constructed. It would be nice to add the example to the library of polytopes (either using lawrence extension or just hard-coding the vertices. Or even better, hard-coding the vertices and in the test of the function for the polytope, construct the polytope using lawrence extension and check that it is really combinatorially isomorphic.

This would be a very nice example. I think this needs its own ticket "adding Perles's example to polytopes library". After that we can add it to the test of the function lawrence_extension .

comment:13 Changed 2 years ago by gh-LaisRast

  • Description modified (diff)
  • Keywords polytopes lawrence extension added
  • Summary changed from Lawrence Extension to Implement Lawrence extension for polytopes

comment:14 Changed 2 years ago by git

  • Commit changed from 69cd2c3f76444244edc5fa2d9d38001ac695ba36 to a975b55312cfd8227a88b934901005dbb6b7ef95

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

a975b55change the name of the variable non_vertices to facet_non_vertices

comment:15 Changed 2 years ago by gh-LaisRast

Sorry. I know this is on "needs_review". I only changed the name of a variable from non_vertices to facet_non_vertices.

comment:16 Changed 2 years ago by gh-kliem

+src/sage/geometry/polyhedron/base.py:509: undefined name 'x'
+src/sage/geometry/polyhedron/base.py:4030: 'sage.matrix.constructor.block_matrix' imported but unused

Maybe one could add x = None, before line 503. Then we don't get that pyflakes error anymore.

Delete line 4030.

comment:17 Changed 2 years ago by gh-kliem

Tiny things:

-        I_n= matrix.identity(n)
-        lambda_V = block_matrix([[V, I_n],[V, 2*I_n]])
+        I_n = matrix.identity(n)
+        lambda_V = block_matrix([[V, I_n], [V, 2*I_n]])
+            For more information, see [BaSt1990]_.
-
+        """
-

I think, there are in general no blank lines after the end of the docstring """.

Also there is a double blank line.

sage: Q = polytopes.regular_polygon(4).pyramid()

Instead of Q, I would maybe use something more descriptive as egyptian or egyptian_py. Then one does not have to read all of polytopes.regular_polygon(4).pyramid() to understand the construction.

This is all that I can come up with.

Last edited 2 years ago by gh-kliem (previous) (diff)

comment:18 Changed 2 years ago by git

  • Commit changed from a975b55312cfd8227a88b934901005dbb6b7ef95 to e5c6cad97afb517d6ccfee370f81869c6a492314

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

e5c6cadfixing pyflakes undefined variable error/cleaning the code

comment:19 Changed 2 years ago by git

  • Commit changed from e5c6cad97afb517d6ccfee370f81869c6a492314 to bc1b9cd71307b02eee6ff17e33aa3fe05089230e

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

bc1b9cdchange variable name from egyption_py to egyption_pyramid

comment:20 Changed 2 years ago by gh-kliem

  • Status changed from needs_review to positive_review

Looks fine to me.

comment:21 Changed 2 years ago by vbraun

  • Branch changed from public/27534 to bc1b9cd71307b02eee6ff17e33aa3fe05089230e
  • Resolution set to fixed
  • Status changed from positive_review to closed

comment:22 Changed 2 years ago by jipilab

  • Commit bc1b9cd71307b02eee6ff17e33aa3fe05089230e deleted

Here is one things that could have been changed in this ticket:

-     raise NotImplementedError("self must be a polytope")
+     raise NotImplementedError("Only polytopes (compact polyhedra) are allowed.")

like in bounding_box or

-     raise NotImplementedError("self must be a polytope")
+     raise NotImplementedError("The polytope has to be compact.")

like in barycentric_subdivision. Usually self is replaced by something more informative about the actual object and making the requirement precise (in this case 'compactness'). Since polytope is not a universally defined word, mentioning compactness makes thing completely unambiguous.

This could be changed some other time when practical.

Note: See TracTickets for help on using tickets.