# A Brief Introduction to Polytopes in Sage#

*Author: sarah-marie belcastro <smbelcas@toroidalsnark.net>*

If you already know some convex geometry *a la* Grünbaum or
Brøndsted, then you may have itched to get your hands dirty with some
polytope calculations. Here is a mini-guide to doing just that.

## Basics#

First, let’s define a polytope as the convex hull of a set of points, i.e. given \(S\) we compute \(P={\rm conv}(S)\):

```
sage: P1 = Polyhedron(vertices = [[-5,2], [4,4], [3,0], [1,0], [2,-4], [-3,-1], [-5,-3]])
sage: P1
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 4 vertices
```

```
>>> from sage.all import *
>>> P1 = Polyhedron(vertices = [[-Integer(5),Integer(2)], [Integer(4),Integer(4)], [Integer(3),Integer(0)], [Integer(1),Integer(0)], [Integer(2),-Integer(4)], [-Integer(3),-Integer(1)], [-Integer(5),-Integer(3)]])
>>> P1
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 4 vertices
```

Notice that Sage tells you the dimension of the polytope as well as the dimension of the ambient space.

Of course, you want to know what this object looks like:

```
sage: P1.plot()
Graphics object consisting of 6 graphics primitives
```

```
>>> from sage.all import *
>>> P1.plot()
Graphics object consisting of 6 graphics primitives
```

Even in only 2 dimensions, it’s a pain to figure out what the supporting hyperplanes are. Luckily Sage will take care of that for us.

```
sage: for q in P1.Hrepresentation():
....: print(q)
An inequality (-4, 1) x + 12 >= 0
An inequality (1, 7) x + 26 >= 0
An inequality (1, 0) x + 5 >= 0
An inequality (2, -9) x + 28 >= 0
```

```
>>> from sage.all import *
>>> for q in P1.Hrepresentation():
... print(q)
An inequality (-4, 1) x + 12 >= 0
An inequality (1, 7) x + 26 >= 0
An inequality (1, 0) x + 5 >= 0
An inequality (2, -9) x + 28 >= 0
```

That notation is not immediately parseable, because seriously, those do not look like equations of lines (or of halfspaces, which is really what they are).

`(-4, 1) x + 12 >= 0`

really means \((-4, 1)\cdot\vec{x} + 12 \geq 0\).

So… if you want to define a polytope via inequalities, you have to translate each inequality into a vector. For example, \((-4, 1)\cdot\vec{x} + 12 \geq 0\) becomes (12, -4, 1).

```
sage: altP1 = Polyhedron(ieqs=[(12, -4, 1), (26, 1, 7),(5,1,0), (28, 2, -9)])
sage: altP1.plot()
Graphics object consisting of 6 graphics primitives
```

```
>>> from sage.all import *
>>> altP1 = Polyhedron(ieqs=[(Integer(12), -Integer(4), Integer(1)), (Integer(26), Integer(1), Integer(7)),(Integer(5),Integer(1),Integer(0)), (Integer(28), Integer(2), -Integer(9))])
>>> altP1.plot()
Graphics object consisting of 6 graphics primitives
```

Other information you might want to pull out of Sage about a polytope is the vertex list, which can be done in two ways:

```
sage: for q in P1.Vrepresentation():
....: print(q)
A vertex at (-5, -3)
A vertex at (-5, 2)
A vertex at (4, 4)
A vertex at (2, -4)
```

```
>>> from sage.all import *
>>> for q in P1.Vrepresentation():
... print(q)
A vertex at (-5, -3)
A vertex at (-5, 2)
A vertex at (4, 4)
A vertex at (2, -4)
```

```
sage: P1.vertices()
(A vertex at (-5, -3), A vertex at (-5, 2), A vertex at (4, 4), A vertex at (2, -4))
```

```
>>> from sage.all import *
>>> P1.vertices()
(A vertex at (-5, -3), A vertex at (-5, 2), A vertex at (4, 4), A vertex at (2, -4))
```

## Polar duals#

Surely you want to compute the polar dual:

```
sage: P1dual = P1.polar()
sage: P1dual
A 2-dimensional polyhedron in QQ^2 defined as the convex hull of 4 vertices
```

```
>>> from sage.all import *
>>> P1dual = P1.polar()
>>> P1dual
A 2-dimensional polyhedron in QQ^2 defined as the convex hull of 4 vertices
```

Check it out---we started with an integer-lattice polytope and dualized to a rational-lattice polytope. Let’s look at that.

```
sage: P1dual.plot()
Graphics object consisting of 6 graphics primitives
```

```
>>> from sage.all import *
>>> P1dual.plot()
Graphics object consisting of 6 graphics primitives
```

```
sage: P1.plot() + P1dual.plot()
Graphics object consisting of 12 graphics primitives
```

```
>>> from sage.all import *
>>> P1.plot() + P1dual.plot()
Graphics object consisting of 12 graphics primitives
```

Oh, yeah, unless the polytope is unit-sphere-sized, the dual will be a very different size. Let’s rescale.

```
sage: ((1/4)*P1).plot() + (4*P1dual).plot()
Graphics object consisting of 12 graphics primitives
```

```
>>> from sage.all import *
>>> ((Integer(1)/Integer(4))*P1).plot() + (Integer(4)*P1dual).plot()
Graphics object consisting of 12 graphics primitives
```

If you think that looks a little bit shady, you’re correct. Here is an example that makes the issue a bit clearer.

```
sage: P2 = Polyhedron(vertices = [[-5,0], [-1,1], [-2,0], [1,0], [-2,-1], [-3,-1], [-5,-1]])
sage: P2
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 5 vertices
sage: P2dual = P2.polar(); P2dual
A 2-dimensional polyhedron in QQ^2 defined as the convex hull of 5 vertices
sage: P2.plot() + P2dual.plot()
Graphics object consisting of 14 graphics primitives
```

```
>>> from sage.all import *
>>> P2 = Polyhedron(vertices = [[-Integer(5),Integer(0)], [-Integer(1),Integer(1)], [-Integer(2),Integer(0)], [Integer(1),Integer(0)], [-Integer(2),-Integer(1)], [-Integer(3),-Integer(1)], [-Integer(5),-Integer(1)]])
>>> P2
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 5 vertices
>>> P2dual = P2.polar(); P2dual
A 2-dimensional polyhedron in QQ^2 defined as the convex hull of 5 vertices
>>> P2.plot() + P2dual.plot()
Graphics object consisting of 14 graphics primitives
```

That is clearly not computing what we think of as the polar dual. But look at this…

```
sage: P2.plot() + (-1*P2dual).plot()
Graphics object consisting of 14 graphics primitives
```

```
>>> from sage.all import *
>>> P2.plot() + (-Integer(1)*P2dual).plot()
Graphics object consisting of 14 graphics primitives
```

Here is what’s going on.

If a polytope `P`

is in \(\ZZ\), then…

…the dual is inverted in some way, which is vertically for polygons.

…the dual is taken of P itself.

…if the origin is not in P, then an error is returned.

However, if a polytope is *not* in \(\ZZ\), for example if it’s in \(\QQ\) or
`RDF`

, then…

(1’) …the dual is not inverted.

(2’) …the dual is taken of P-translated-so-barycenter-is-at-origin.

Keep all of this in mind as you take polar duals.

## Polytope Constructions#

Minkowski sums! Now with two syntaxes!

```
sage: P1+P2
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 8 vertices
```

```
>>> from sage.all import *
>>> P1+P2
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 8 vertices
```

```
sage: P1.minkowski_sum(P2)
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 8 vertices
```

```
>>> from sage.all import *
>>> P1.minkowski_sum(P2)
A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 8 vertices
```

Okay, fine. We should have some 3-dimensional examples, at least. (Note that in order to display polytopes effectively you’ll need visualization software such as Javaview and Jmol installed.)

```
sage: P3 = Polyhedron(vertices=[(0,0,0), (0,0,1/2), (0,1/2,0), (1/2,0,0), (3/4,1/5,3/2)]); P3
A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 5 vertices
sage: P4 = Polyhedron(vertices=[(-1,1,0),(1,1,0),(-1,0,1), (1,0,1),(0,-1,1),(0,1,1)]); P4
A 3-dimensional polyhedron in ZZ^3 defined as the convex hull of 6 vertices
sage: P3.plot() + P4.plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> P3 = Polyhedron(vertices=[(Integer(0),Integer(0),Integer(0)), (Integer(0),Integer(0),Integer(1)/Integer(2)), (Integer(0),Integer(1)/Integer(2),Integer(0)), (Integer(1)/Integer(2),Integer(0),Integer(0)), (Integer(3)/Integer(4),Integer(1)/Integer(5),Integer(3)/Integer(2))]); P3
A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 5 vertices
>>> P4 = Polyhedron(vertices=[(-Integer(1),Integer(1),Integer(0)),(Integer(1),Integer(1),Integer(0)),(-Integer(1),Integer(0),Integer(1)), (Integer(1),Integer(0),Integer(1)),(Integer(0),-Integer(1),Integer(1)),(Integer(0),Integer(1),Integer(1))]); P4
A 3-dimensional polyhedron in ZZ^3 defined as the convex hull of 6 vertices
>>> P3.plot() + P4.plot()
Graphics3d Object
```

```
sage: (P3+P4).plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> (P3+P4).plot()
Graphics3d Object
```

We can also find the intersection of two polytopes… and this too has two syntaxes!

```
sage: int12 = P1.intersection(P2*.5); int12.plot()
Graphics object consisting of 7 graphics primitives
```

```
>>> from sage.all import *
>>> int12 = P1.intersection(P2*RealNumber('.5')); int12.plot()
Graphics object consisting of 7 graphics primitives
```

```
sage: int34 = P3 & P4; int34.plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> int34 = P3 & P4; int34.plot()
Graphics3d Object
```

Should one wish to translate, one can.

```
sage: transP2 = P2.translation([2,1])
sage: P2.plot() + transP2.plot()
Graphics object consisting of 14 graphics primitives
```

```
>>> from sage.all import *
>>> transP2 = P2.translation([Integer(2),Integer(1)])
>>> P2.plot() + transP2.plot()
Graphics object consisting of 14 graphics primitives
```

Then of course we can take prisms, pyramids, and bipyramids of polytopes…

```
sage: P2.prism().plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> P2.prism().plot()
Graphics3d Object
```

```
sage: P1.pyramid().plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> P1.pyramid().plot()
Graphics3d Object
```

```
sage: P2dual.bipyramid().plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> P2dual.bipyramid().plot()
Graphics3d Object
```

Okay, fine. Yes, Sage has some kinds of polytopes built in.
If you type `polytopes.`

and then press `TAB`

after the period, you’ll get a
list of pre-built polytopes.

```
sage: P5 = polytopes.hypercube(5)
sage: P6 = polytopes.cross_polytope(3)
sage: P7 = polytopes.simplex(7)
```

```
>>> from sage.all import *
>>> P5 = polytopes.hypercube(Integer(5))
>>> P6 = polytopes.cross_polytope(Integer(3))
>>> P7 = polytopes.simplex(Integer(7))
```

Let’s look at a 4-dimensional polytope.

```
sage: P8 = polytopes.hypercube(4)
sage: P8.plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> P8 = polytopes.hypercube(Integer(4))
>>> P8.plot()
Graphics3d Object
```

We can see it from a different perspective:

```
sage: P8.schlegel_projection(position=1/2).plot()
Graphics3d Object
```

```
>>> from sage.all import *
>>> P8.schlegel_projection(position=Integer(1)/Integer(2)).plot()
Graphics3d Object
```

## Queries to polytopes#

Once you’ve constructed some polytope, you can ask Sage questions about it.

```
sage: P1.contains([1,0])
True
```

```
>>> from sage.all import *
>>> P1.contains([Integer(1),Integer(0)])
True
```

```
sage: P1.interior_contains([3,0])
False
```

```
>>> from sage.all import *
>>> P1.interior_contains([Integer(3),Integer(0)])
False
```

```
sage: P3.contains([1,0,0])
False
```

```
>>> from sage.all import *
>>> P3.contains([Integer(1),Integer(0),Integer(0)])
False
```

Face information can be useful.

```
sage: int34.f_vector()
(1, 8, 12, 6, 1)
```

```
>>> from sage.all import *
>>> int34.f_vector()
(1, 8, 12, 6, 1)
```

Well, geometric information might be *more* helpful…
Here we are told which of the vertices form each 2-face:

```
sage: [f.ambient_V_indices() for f in int34.faces(2)]
[(2, 6, 7), (0, 1, 3, 5), (1, 3, 4), (0, 5, 6, 7), (0, 1, 2, 4, 6), (2, 3, 4, 5, 7)]
```

```
>>> from sage.all import *
>>> [f.ambient_V_indices() for f in int34.faces(Integer(2))]
[(2, 6, 7), (0, 1, 3, 5), (1, 3, 4), (0, 5, 6, 7), (0, 1, 2, 4, 6), (2, 3, 4, 5, 7)]
```

Yeah, that isn’t so useful as it is. Let’s figure out the vertex and hyperplane representations of the first face in the list.

```
sage: first2faceofint34 = int34.faces(2)[0]
sage: first2faceofint34.ambient_Hrepresentation(); first2faceofint34.vertices()
(An inequality (0, 0, -1) x + 1 >= 0,)
(A vertex at (2/3, 2/15, 1), A vertex at (3/8, 1/10, 1), A vertex at (1/2, 3/10, 1))
```

```
>>> from sage.all import *
>>> first2faceofint34 = int34.faces(Integer(2))[Integer(0)]
>>> first2faceofint34.ambient_Hrepresentation(); first2faceofint34.vertices()
(An inequality (0, 0, -1) x + 1 >= 0,)
(A vertex at (2/3, 2/15, 1), A vertex at (3/8, 1/10, 1), A vertex at (1/2, 3/10, 1))
```

If you want more… Base class for polyhedra: Miscellaneous methods is the first place you want to go.