Constructions
The standard way to define a polyhedron is by either giving a $V$-representation or an $H$-representation. But polyhedra may also be constructed through other means: by name, via operations on other polyhedra, or from other objects in OSCAR.
$H$- and $V$-representations
Intersecting halfspaces: $H$-representation
Polyhedron
— MethodPolyhedron{T}(A::AnyVecOrMat, b) where T<:scalar_types
The (convex) polyhedron defined by
\[P(A,b) = \{ x | Ax ≤ b \}.\]
see Def. 3.35 and Section 4.1. of Michael Joswig, Thorsten Theobald (2013)
Examples
The following lines define the square $[0,1]^2 \subset \mathbb{R}^2$:
julia> A = [1 0; 0 1; -1 0 ; 0 -1];
julia> b = [1, 1, 0, 0];
julia> Polyhedron(A,b)
A polyhedron in ambient dimension 2
Polyhedron
— MethodPolyhedron{T}(I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing) where T<:scalar_types
The (convex) polyhedron obtained intersecting the halfspaces I
(inequalities) and the hyperplanes E
(equations).
Examples
The following lines define the square $[0,1]^2 \subset \mathbb{R}^2$:
julia> A = [1 0; 0 1; -1 0 ; 0 -1];
julia> b = [1, 1, 0, 0];
julia> Polyhedron((A,b))
A polyhedron in ambient dimension 2
As an example for a polyhedron constructed from both inequalities and equations, we construct the polytope $[0,1]\times\{0\}\subset\mathbb{R}^2$
julia> P = Polyhedron(([-1 0; 1 0], [0,1]), ([0 1], [0]))
A polyhedron in ambient dimension 2
julia> is_feasible(P)
true
julia> dim(P)
1
julia> vertices(P)
2-element SubObjectIterator{PointVector{fmpq}}:
[1, 0]
[0, 0]
The complete $H$-representation can be retrieved using facets
and affine_hull
:
julia> P = Polyhedron(([-1 0; 1 0], [0,1]), ([0 1], [0]))
A polyhedron in ambient dimension 2
julia> facets(P)
2-element SubObjectIterator{AffineHalfspace{fmpq}} over the Halfspaces of R^2 described by:
-x₁ ≦ 0
x₁ ≦ 1
julia> affine_hull(P)
1-element SubObjectIterator{AffineHyperplane{fmpq}} over the Hyperplanes of R^2 described by:
x₂ = 0
julia> Q0 = Polyhedron(facets(P))
A polyhedron in ambient dimension 2
julia> P == Q0
false
julia> Q1 = Polyhedron(facets(P), affine_hull(P))
A polyhedron in ambient dimension 2
julia> P == Q1
true
Computing convex hulls: $V$-representation
convex_hull
— Methodconvex_hull([::Type{T} = fmpq,] V [, R [, L]]; non_redundant::Bool = false)
Construct the convex hull of the vertices V
, rays R
, and lineality L
. If R
or L
are omitted, then they are assumed to be zero.
Arguments
V::AbstractCollection[PointVector]
: Points whose convex hull is to be computed.R::AbstractCollection[RayVector]
: Rays completing the set of points.L::AbstractCollection[RayVector]
: Generators of the Lineality space.
If an argument is given as a matrix, its content has to be encoded row-wise.
R
can be given as an empty matrix or as nothing
if the user wants to compute the convex hull only from V
and L
.
If it is known that V
and R
only contain extremal points and that the description of the lineality space is complete, set non_redundant = true
to avoid unnecessary redundancy checks.
See Def. 2.11 and Def. 3.1 of Michael Joswig, Thorsten Theobald (2013).
Examples
The following lines define the square $[0,1]^2 \subset \mathbb{R}^2$:
julia> Square = convex_hull([0 0; 0 1; 1 0; 1 1])
A polyhedron in ambient dimension 2
To construct the positive orthant, rays have to be passed:
julia> V = [0 0];
julia> R = [1 0; 0 1];
julia> PO = convex_hull(V, R)
A polyhedron in ambient dimension 2
The closed-upper half plane can be constructed by passing rays and a lineality space:
julia> V = [0 0];
julia> R = [0 1];
julia> L = [1 0];
julia> UH = convex_hull(V, R, L)
A polyhedron in ambient dimension 2
To obtain the x-axis in $\mathbb{R}^2$:
julia> V = [0 0];
julia> R = nothing;
julia> L = [1 0];
julia> XA = convex_hull(V, R, L)
A polyhedron in ambient dimension 2
This is a standard triangle, defined via a (redundant) $V$-representation and its unique minimal $H$-representation:
julia> T = convex_hull([ 0 0 ; 1 0 ; 0 1; 0 1/2 ])
A polyhedron in ambient dimension 2
julia> halfspace_matrix_pair(facets(T))
(A = [-1 0; 0 -1; 1 1], b = Polymake.RationalAllocated[0, 0, 1])
The complete $V$-representation can be retrieved using vertices
, rays
and lineality_space
:
julia> P = convex_hull([0 0], [1 0], [0 1])
A polyhedron in ambient dimension 2
julia> Q0 = convex_hull(vertices(P))
A polyhedron in ambient dimension 2
julia> P == Q0
false
julia> Q1 = convex_hull(vertices(P), rays(P))
A polyhedron in ambient dimension 2
julia> P == Q1
false
julia> Q0 == Q1
false
julia> Q2 = convex_hull(vertices(P), rays(P), lineality_space(P))
A polyhedron in ambient dimension 2
julia> P == Q2
true
Named polyhedra
archimedean_solid
— Functionarchimedean_solid(s)
Construct an Archimedean solid with the name given by String s
from the list below. The polytopes are realized with floating point numbers and thus not exact; Vertex-facet-incidences are correct in all cases.
Arguments
s::String
: The name of the desired Archimedean solid. Possible values:- "truncated_tetrahedron" : Truncated tetrahedron. Regular polytope with four triangular and four hexagonal facets.
- "cuboctahedron" : Cuboctahedron. Regular polytope with eight triangular and six square facets.
- "truncated_cube" : Truncated cube. Regular polytope with eight triangular and six octagonal facets.
- "truncated_octahedron" : Truncated Octahedron. Regular polytope with six square and eight hexagonal facets.
- "rhombicuboctahedron" : Rhombicuboctahedron. Regular polytope with eight triangular and 18 square facets.
- "truncated_cuboctahedron" : Truncated Cuboctahedron. Regular polytope with 12 square, eight hexagonal and six octagonal facets.
- "snub_cube" : Snub Cube. Regular polytope with 32 triangular and six square facets. The vertices are realized as floating point numbers. This is a chiral polytope.
- "icosidodecahedron" : Icosidodecahedon. Regular polytope with 20 triangular and 12 pentagonal facets.
- "truncated_dodecahedron" : Truncated Dodecahedron. Regular polytope with 20 triangular and 12 decagonal facets.
- "truncated_icosahedron" : Truncated Icosahedron. Regular polytope with 12 pentagonal and 20 hexagonal facets.
- "rhombicosidodecahedron" : Rhombicosidodecahedron. Regular polytope with 20 triangular, 30 square and 12 pentagonal facets.
- "truncated_icosidodecahedron" : Truncated Icosidodecahedron. Regular polytope with 30 square, 20 hexagonal and 12 decagonal facets.
- "snub_dodecahedron" : Snub Dodecahedron. Regular polytope with 80 triangular and 12 pentagonal facets. The vertices are realized as floating point numbers. This is a chiral polytope.
Examples
julia> T = archimedean_solid("cuboctahedron")
A polyhedron in ambient dimension 3
julia> sum([nvertices(F) for F in faces(T, 2)] .== 3)
8
julia> sum([nvertices(F) for F in faces(T, 2)] .== 4)
6
julia> nfacets(T)
14
birkhoff
— Functionbirkhoff(n::Integer, even::Bool = false)
Construct the Birkhoff polytope of dimension $n^2$.
This is the polytope of $n \times n$ stochastic matrices (encoded as row vectors of length $n^2$), i.e., the matrices with non-negative real entries whose row and column entries sum up to one. Its vertices are the permutation matrices.
Use even = true
to get the vertices only for the even permutation matrices.
Examples
julia> b = birkhoff(3)
A polyhedron in ambient dimension 9
julia> vertices(b)
6-element SubObjectIterator{PointVector{fmpq}}:
[1, 0, 0, 0, 1, 0, 0, 0, 1]
[0, 1, 0, 1, 0, 0, 0, 0, 1]
[0, 0, 1, 1, 0, 0, 0, 1, 0]
[1, 0, 0, 0, 0, 1, 0, 1, 0]
[0, 1, 0, 0, 0, 1, 1, 0, 0]
[0, 0, 1, 0, 1, 0, 1, 0, 0]
catalan_solid
— Functioncatalan_solid(s::String)
Construct a Catalan solid with the name s
from the list below. The polytopes are realized with floating point coordinates and thus are not exact. However, vertex-facet-incidences are correct in all cases.
Arguments
s::String
: The name of the desired Archimedean solid. Possible values:- "triakis_tetrahedron" : Triakis Tetrahedron. Dual polytope to the Truncated Tetrahedron, made of 12 isosceles triangular facets.
- "triakis_octahedron" : Triakis Octahedron. Dual polytope to the Truncated Cube, made of 24 isosceles triangular facets.
- "rhombic_dodecahedron" : Rhombic dodecahedron. Dual polytope to the cuboctahedron, made of 12 rhombic facets.
- "tetrakis_hexahedron" : Tetrakis hexahedron. Dual polytope to the truncated octahedron, made of 24 isosceles triangluar facets.
- "disdyakis_dodecahedron" : Disdyakis dodecahedron. Dual polytope to the truncated cuboctahedron, made of 48 scalene triangular facets.
- "pentagonal_icositetrahedron" : Pentagonal Icositetrahedron. Dual polytope to the snub cube, made of 24 irregular pentagonal facets. The vertices are realized as floating point numbers.
- "pentagonal_hexecontahedron" : Pentagonal Hexecontahedron. Dual polytope to the snub dodecahedron, made of 60 irregular pentagonal facets. The vertices are realized as floating point numbers.
- "rhombic_triacontahedron" : Rhombic triacontahedron. Dual polytope to the icosidodecahedron, made of 30 rhombic facets.
- "triakis_icosahedron" : Triakis icosahedron. Dual polytope to the icosidodecahedron, made of 30 rhombic facets.
- "deltoidal_icositetrahedron" : Deltoidal Icositetrahedron. Dual polytope to the rhombicubaoctahedron, made of 24 kite facets.
- "pentakis_dodecahedron" : Pentakis dodecahedron. Dual polytope to the truncated icosahedron, made of 60 isosceles triangular facets.
- "deltoidal_hexecontahedron" : Deltoidal hexecontahedron. Dual polytope to the rhombicosidodecahedron, made of 60 kite facets.
- "disdyakis_triacontahedron" : Disdyakis triacontahedron. Dual polytope to the truncated icosidodecahedron, made of 120 scalene triangular facets.
Examples
julia> T = catalan_solid("triakis_tetrahedron");
julia> count(F -> nvertices(F) == 3, faces(T, 2))
12
julia> nfacets(T)
12
cross
— Functioncross([::Type{T} = fmpq,] d::Int [,n::Rational])
Construct a $d$-dimensional cross polytope around origin with vertices located at $\pm e_i$ for each unit vector $e_i$ of $R^d$, scaled by $n$.
Examples
Here we print the facets of a non-scaled and a scaled 3-dimensional cross polytope:
julia> C = cross(3)
A polyhedron in ambient dimension 3
julia> facets(C)
8-element SubObjectIterator{AffineHalfspace{fmpq}} over the Halfspaces of R^3 described by:
x₁ + x₂ + x₃ ≦ 1
-x₁ + x₂ + x₃ ≦ 1
x₁ - x₂ + x₃ ≦ 1
-x₁ - x₂ + x₃ ≦ 1
x₁ + x₂ - x₃ ≦ 1
-x₁ + x₂ - x₃ ≦ 1
x₁ - x₂ - x₃ ≦ 1
-x₁ - x₂ - x₃ ≦ 1
julia> D = cross(3, 2)
A polyhedron in ambient dimension 3
julia> facets(D)
8-element SubObjectIterator{AffineHalfspace{fmpq}} over the Halfspaces of R^3 described by:
x₁ + x₂ + x₃ ≦ 2
-x₁ + x₂ + x₃ ≦ 2
x₁ - x₂ + x₃ ≦ 2
-x₁ - x₂ + x₃ ≦ 2
x₁ + x₂ - x₃ ≦ 2
-x₁ + x₂ - x₃ ≦ 2
x₁ - x₂ - x₃ ≦ 2
-x₁ - x₂ - x₃ ≦ 2
cube
— Functioncube([::Type{T} = fmpq,] d::Int , [l::Rational = -1, u::Rational = 1])
Construct the $[l,u]$-cube in dimension $d$.
Examples
In this example the 5-dimensional unit cube is constructed to ask for one of its properties:
julia> C = cube(5,0,1);
julia> normalized_volume(C)
120
cyclic_polytope
— Functioncyclic_polytope(d::Int, n::Int)
Construct the cyclic polytope that is the convex hull of $n$ points on the moment curve in dimension $d$.
Examples
julia> cp = cyclic_polytope(3, 20)
A polyhedron in ambient dimension 3
julia> nvertices(cp)
20
del_pezzo_polytope
— Functiondel_pezzo_polytope(d::Int)
Produce the d-dimensional del Pezzo polytope, which is the convex hull of the cross polytope together with the all-ones and minus all-ones vector.
julia> DP = del_pezzo_polytope(4)
A polyhedron in ambient dimension 4
julia> f_vector(DP)
4-element Vector{fmpz}:
10
40
60
30
fano_simplex
— Functionfano_simplex(d::Int)
Construct a lattice simplex such that the origin is the unique interior lattice point. The normal toric variety associated with its face fan is smooth.
julia> S = fano_simplex(3)
A polyhedron in ambient dimension 3
julia> X = NormalToricVariety(face_fan(S))
A normal toric variety
julia> is_smooth(X)
true
fractional_cut_polytope
— Functionfractional_cut_polytope(G::Graph{Undirected})
Construct the fractional cut polytope of the graph $G$.
Examples
julia> G = complete_graph(4);
julia> fractional_cut_polytope(G)
A polyhedron in ambient dimension 6
fractional_matching_polytope
— Functionfractional_matching_polytope(G::Graph{Undirected})
Construct the fractional matching polytope of the graph $G$.
Examples
julia> G = complete_graph(4);
julia> fractional_matching_polytope(G)
A polyhedron in ambient dimension 6
gelfand_tsetlin
— Functiongelfand_tsetlin(lambda::AbstractVector)
Construct the Gelfand Tsetlin polytope indexed by a weakly decreasing vector lambda
.
julia> P = gelfand_tsetlin([5,3,2])
A polyhedron in ambient dimension 6
julia> is_fulldimensional(P)
false
julia> p = project_full(P)
A polyhedron in ambient dimension 3
julia> is_fulldimensional(p)
true
julia> volume(p)
3
simplex
— Functionsimplex([::Type{T} = fmpq,] d::Int [,n::Rational])
Construct the simplex which is the convex hull of the standard basis vectors along with the origin in $\mathbb{R}^d$, scaled by $n$.
Examples
Here we take a look at the facets of the 7-simplex and a scaled 7-simplex:
julia> s = simplex(7)
A polyhedron in ambient dimension 7
julia> facets(s)
8-element SubObjectIterator{AffineHalfspace{fmpq}} over the Halfspaces of R^7 described by:
-x₁ ≦ 0
-x₂ ≦ 0
-x₃ ≦ 0
-x₄ ≦ 0
-x₅ ≦ 0
-x₆ ≦ 0
-x₇ ≦ 0
x₁ + x₂ + x₃ + x₄ + x₅ + x₆ + x₇ ≦ 1
julia> t = simplex(7, 5)
A polyhedron in ambient dimension 7
julia> facets(t)
8-element SubObjectIterator{AffineHalfspace{fmpq}} over the Halfspaces of R^7 described by:
-x₁ ≦ 0
-x₂ ≦ 0
-x₃ ≦ 0
-x₄ ≦ 0
-x₅ ≦ 0
-x₆ ≦ 0
-x₇ ≦ 0
x₁ + x₂ + x₃ + x₄ + x₅ + x₆ + x₇ ≦ 5
Operations on polyhedra
Polyhedra can be produced through operations on other polyhedra. For example, they can be added using Minkowski addition or scaled; each of which results in a new polyhedron.
+
— Method+(P::Polyhedron, Q::Polyhedron)
Return the Minkowski sum $P + Q = \{ x+y\ |\ x∈P, y∈Q\}$ of P
and Q
(see also minkowski_sum
).
Examples
The Minkowski sum of a square and the 2-dimensional cross-polytope is an octagon:
julia> P = cube(2);
julia> Q = cross(2);
julia> M = minkowski_sum(P, Q)
A polyhedron in ambient dimension 2
julia> nvertices(M)
8
*
— Method*(k::Int, Q::Polyhedron)
Return the scaled polyhedron $kQ = \{ kx\ |\ x∈Q\}$.
Note that k*Q = Q*k
.
Examples
Scaling an $n$-dimensional bounded polyhedron by the factor $k$ results in the volume being scaled by $k^n$. This example confirms the statement for the 6-dimensional cube and $k = 2$.
julia> C = cube(6);
julia> SC = 2*C
A polyhedron in ambient dimension 6
julia> volume(SC)//volume(C)
64
*
— Method*(P::Polyhedron, Q::Polyhedron)
Return the Cartesian product of P
and Q
(see also product
).
Examples
The Cartesian product of a triangle and a line segment is a triangular prism.
julia> T=simplex(2)
A polyhedron in ambient dimension 2
julia> S=cube(1)
A polyhedron in ambient dimension 1
julia> length(vertices(T*S))
6
bipyramid
— Functionbipyramid(P::Polyhedron, z::Number = 1, z_prime::Number = -z)
Make a bipyramid over a pointed polyhedron P
.
The bipyramid is the convex hull of the input polyhedron P
and two apexes (v
, z
), (v
, z_prime
) on both sides of the affine span of P
. For bounded polyhedra, the projections of the apexes v
to the affine span of P
is the vertex barycenter of P
.
Examples
julia> c = cube(2)
A polyhedron in ambient dimension 2
julia> vertices(bipyramid(c,2))
6-element SubObjectIterator{PointVector{fmpq}}:
[-1, -1, 0]
[1, -1, 0]
[-1, 1, 0]
[1, 1, 0]
[0, 0, 2]
[0, 0, -2]
intersect
— Methodintersect(P::Polyhedron, Q::Polyhedron)
Return the intersection $P \cap Q$ of P
and Q
.
Examples
The positive orthant of the plane is the intersection of the two halfspaces with $x≥0$ and $y≥0$ respectively.
julia> UH1 = convex_hull([0 0],[1 0],[0 1]);
julia> UH2 = convex_hull([0 0],[0 1],[1 0]);
julia> PO = intersect(UH1, UH2)
A polyhedron in ambient dimension 2
julia> rays(PO)
2-element SubObjectIterator{RayVector{fmpq}}:
[1, 0]
[0, 1]
pyramid
— Functionpyramid(P::Polyhedron, z::Number = 1)
Make a pyramid over the given polyhedron P
.
The pyramid is the convex hull of the input polyhedron P
and a point v
outside the affine span of P
. For bounded polyhedra, the projection of v
to the affine span of P
coincides with the vertex barycenter of P
. The scalar z
is the distance between the vertex barycenter and v
.
Examples
julia> c = cube(2)
A polyhedron in ambient dimension 2
julia> vertices(pyramid(c,5))
5-element SubObjectIterator{PointVector{fmpq}}:
[-1, -1, 0]
[1, -1, 0]
[-1, 1, 0]
[1, 1, 0]
[0, 0, 5]
The convex hull of two polytopes can be computed via convex_hull
.
convex_hull
— Methodconvex_hull(P::Polyhedron, Q::Polyhedron)
Return the convex_hull of P
and Q
.
Examples
The convex hull of the following two line segments in $R^3$ is a tetrahedron.
julia> L₁ = convex_hull([-1 0 0; 1 0 0])
A polyhedron in ambient dimension 3
julia> L₂ = convex_hull([0 -1 0; 0 1 0])
A polyhedron in ambient dimension 3
julia> T=convex_hull(L₁,L₂);
julia> f_vector(T)
2-element Vector{fmpz}:
4
4
Polyhedra from other mathematical objects
orbit_polytope
— Functionorbit_polytope(V::AbstractCollection[PointVector], G::PermGroup)
Construct the convex hull of the orbit of one or several points (given row-wise in V
) under the action of G
.
Examples
This will construct the $3$-dimensional permutahedron:
julia> V = [1 2 3];
julia> G = symmetric_group(3);
julia> P = orbit_polytope(V, G)
A polyhedron in ambient dimension 3
julia> vertices(P)
6-element SubObjectIterator{PointVector{fmpq}}:
[1, 2, 3]
[1, 3, 2]
[2, 1, 3]
[2, 3, 1]
[3, 1, 2]
[3, 2, 1]
newton_polytope
— Functionnewton_polytope(poly::Polynomial)
Compute the Newton polytope of the multivariate polynomial poly
.
Examples
julia> S, (x, y) = PolynomialRing(ZZ, ["x", "y"])
(Multivariate Polynomial Ring in x, y over Integer Ring, fmpz_mpoly[x, y])
julia> f = x^3*y + 3x*y^2 + 1
x^3*y + 3*x*y^2 + 1
julia> NP = newton_polytope(f)
A polyhedron in ambient dimension 2
julia> vertices(NP)
3-element SubObjectIterator{PointVector{fmpq}}:
[3, 1]
[1, 2]
[0, 0]