#  Copyright (c) 1997-2024
#  Ewgenij Gawrilow, Michael Joswig, and the polymake team
#  Technische Universität Berlin, Germany
#  https://polymake.org
#
#  This program is free software; you can redistribute it and/or modify it
#  under the terms of the GNU General Public License as published by the
#  Free Software Foundation; either version 2, or (at your option) any
#  later version: http://www.gnu.org/licenses/gpl.txt.
#
#  This program is distributed in the hope that it will be useful,
#  but WITHOUT ANY WARRANTY; without even the implied warranty of
#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#  GNU General Public License for more details.
#-------------------------------------------------------------------------------


# @topic objects/Polytope/specializations/Polytope<Rational>
# A rational polyhedron realized in Q^d

# @topic objects/Polytope/specializations/Polytope<Float>
# A pointed polyhedron with float coordinates realized in R<sup>d</sup>.
#
# It mainly exists for visualization.
#
# Convex hull and related algorithms use floating-point arithmetics.
# Due to numerical errors inherent to this kind of computations, the resulting
# combinatorial description can be arbitrarily far away from the truth, or even
# not correspond to any valid polytope.  You have been warned.
#
# None of the standard construction clients produces objects of this type.
# If you want to get one, create it with the explicit constructor or [[convert_to]].

object Polytope {

file_suffix poly

# @category Input property
# Points such that the polyhedron is their convex hull.
# Redundancies are allowed.
# The vector (x<sub>0</sub>, x<sub>1</sub>, ... x<sub>d</sub>) represents a point in d-space given in homogeneous coordinates.
# Affine points are identified by x<sub>0</sub> > 0.
# Points with x<sub>0</sub> = 0 can be interpreted as rays.
#
# polymake automatically normalizes each coordinate vector, dividing them by the first non-zero element.
# The clients and rule subroutines can always assume that x<sub>0</sub> is either 0 or 1.
# All vectors in this section must be non-zero.
# Dual to [[INEQUALITIES]].
#
# Input section only.  Ask for [[VERTICES]] if you want to compute a V-representation from an H-representation.
# @example Given some (homogeneous) points in 3-space we first construct a matrix containing them. Assume we don't know wether these are all
# vertices of their convex hull or not. To safely produce a polytope from these points, we set the input to the matrix representing them.
# In the following the points under consideration are the vertices of the 3-simplex together with their barycenter, which will be no vertex:
# > $M = new Matrix([[1,0,0,0],[1,1,0,0],[1,0,1,0],[1,0,0,1],[1,1/4,1/4,1/4]]);
# > $p = new Polytope(POINTS=>$M);
# > print $p->VERTICES;
# | 1 0 0 0
# | 1 1 0 0
# | 1 0 1 0
# | 1 0 0 1
property POINTS = override INPUT_RAYS;


# @category Geometry
# Vertices of the polyhedron. No redundancies are allowed.
# All vectors in this section must be non-zero.
# The coordinates are normalized the same way as [[POINTS]]. Dual to [[FACETS]].
# This section is empty if and only if the polytope is empty.
# The property [[VERTICES]] appears only in conjunction with the property [[LINEALITY_SPACE]].
# The specification of the property [[VERTICES]] requires the specification of [[LINEALITY_SPACE]], and vice versa.
# @example To print the vertices (in homogeneous coordinates) of the standard 2-simplex, i.e. a right-angled isoceles triangle, type this:
# > print simplex(2)->VERTICES;
# | (3) (0 1)
# | 1 1 0
# | 1 0 1
# @example If we know some points to be vertices of their convex hull, we can store them as rows in a Matrix and construct a new polytope with it.
# The following produces a 3-dimensioanl pyramid over the standard 2-simplex with the specified vertices:
# > $M = new Matrix([[1,0,0,0],[1,1,0,0],[1,0,1,0],[1,0,0,3]]);
# > $p = new Polytope(VERTICES=>$M);
# @example The following adds a (square) pyramid to one facet of a 3-cube. We do this by extracting the vertices of the cube via the built-in
# method and then attach the apex of the pyramid to the matrix.
# > $v = new Vector([1,0,0,3/2]);
# > $M = cube(3)->VERTICES / $v;
# > $p = new Polytope(VERTICES=>$M);
property VERTICES = override RAYS;

# @topic //properties//INEQUALITIES
# Inequalities that describe half-spaces such that the polyhedron is their intersection.
# Redundancies are allowed.  Dual to [[POINTS]].
#
# A vector (A<sub>0</sub>, A<sub>1</sub>, ..., A<sub>d</sub>) defines the
# (closed affine) half-space of points (1, x<sub>1</sub>, ..., x<sub>d</sub>) such that
# A<sub>0</sub> + A<sub>1</sub> x<sub>1</sub> + ... + A<sub>d</sub> x<sub>d</sub> >= 0.
#
# Input section only.  Ask for [[FACETS]] and [[AFFINE_HULL]] if you want to compute an H-representation
# from a V-representation.


# @topic //properties//EQUATIONS
# Equations that hold for all points of the polyhedron.
#
# A vector (A<sub>0</sub>, A<sub>1</sub>, ..., A<sub>d</sub>) describes the hyperplane
# of all points (1, x<sub>1</sub>, ..., x<sub>d</sub>) such that A<sub>0</sub> + A<sub>1</sub> x<sub>1</sub> + ... + A<sub>d</sub> x<sub>d</sub> = 0.
# All vectors in this section must be non-zero.
#
# Input section only.  Ask for [[AFFINE_HULL]] if you want to see an irredundant description of the affine span.

# @category Geometry
# Dual basis of the affine hull of the polyhedron.
# The property [[AFFINE_HULL]] appears only in conjunction with the property [[FACETS]].
# The specification of the property [[FACETS]] requires the specification of [[AFFINE_HULL]], and vice versa.
property AFFINE_HULL = override LINEAR_SPAN {
   sub canonical { &orthogonalize_affine_subspace; }
}

method canonical_ineq {
   add_extra_polytope_ineq($_[1]);
}


# @category Geometry
# The i-th row is the normal vector of a hyperplane separating the i-th vertex from the others.
# This property is a by-product of redundant point elimination algorithm.
# All vectors in this section must be non-zero.
# @example [require bundled:cdd] This prints a matrix in which each row represents a normal vector of a hyperplane seperating one vertex of a centered square
# with side length 2 from the other ones. The first and the last hyperplanes as well as the second and third hyperplanes are the same
# up to orientation.
# > print cube(2)->VERTEX_NORMALS;
# | 0 1/2 1/2
# | 0 -1/2 1/2
# | 0 1/2 -1/2
# | 0 -1/2 -1/2
property VERTEX_NORMALS = override RAY_SEPARATORS;

# @category Geometry
# Some point belonging to the polyhedron.
# @example This stores a (homogeneous) point belonging to the 3-cube as a vector and prints its coordinates:
# > $v = cube(3)->VALID_POINT;
# > print $v;
# | 1 -1 -1 -1
property VALID_POINT : Vector<Scalar> {
   sub canonical { &canonicalize_rays; }
}

# @category Geometry
# True if the polyhedron is not empty.
property FEASIBLE : Bool;

# @topic //properties//CONE_DIM
# One more than the dimension of the affine hull of the polyhedron
#   = one more than the dimension of the polyhedron.
#   = dimension of the homogenization of the polyhedron
# If the polytope is given purely combinatorially, this is the dimension of a minimal embedding space
# @example This prints the cone dimension of a 3-cube. Since the dimension of its affine closure is 3, the result is 4.
# > print cube(3)->CONE_DIM;
# | 4

# @topic //properties//CONE_AMBIENT_DIM
# One more than the dimension of the space in which the polyhedron lives.
#   = dimension of the space in which the homogenization of the polyhedron lives

# @topic //properties//POINTED
# True if the polyhedron does not contain an affine line.
# @example A square does not contain an affine line and is therefore pointed. Removing one facet does not change this, although
# it is no longer bounded.  After removing two opposing facets, it contains infinitely many affine lines orthogonal to the
# removed facets.
# > $p = cube(2);
# > print $p->POINTED;
# | true
# > print facet_to_infinity($p,0)->POINTED;
# | true
# > print new Polytope(INEQUALITIES=>$p->FACETS->minor([0,1],All))->POINTED;
# | false

# @category Geometry
# A vertex of a pointed polyhedron.
# @example This prints the first vertex of the 3-cube (corresponding to the first row in the vertex matrix).
# > print cube(3)->ONE_VERTEX;
# | 1 -1 -1 -1
property ONE_VERTEX = override ONE_RAY;


# @category Geometry
# True if and only if [[LINEALITY_SPACE]] trivial and [[FAR_FACE]] is trivial.
# @example A pyramid over a square is bounded. Removing the base square yields an unbounded pointed polyhedron
# (the vertices with first entry equal to zero correspond to rays).
# > $p = pyramid(cube(2));
# > print $p->BOUNDED;
# | true
# > $q = facet_to_infinity($p,4);
# > print $q->BOUNDED;
# | false
property BOUNDED : Bool;


# @category Geometry
# True if (1, 0, 0, ...) is in the relative interior.
# If full-dimensional then polar to [[BOUNDED]].
# @example The cube [0,1]^3 is not centered, since the origin is on the boundary. By a small translation we can make it centered:
# > $p = cube(3,0,0);
# > print $p->CENTERED;
# | false
# > $t = new Vector([-1/2,-1/2,-1/2]);
# > print translate($p,$t)->CENTERED;
# | true
property CENTERED : Bool;


# @category Geometry
# True if (1, 0, 0, ...) is contained (possibly in the boundary).
# @example The cube [0,1]^3 is only weakly centered, since the origin is on the boundary.
# > $p = cube(3,0,0);
# > print $p->WEAKLY_CENTERED;
# | true
# > print $p->CENTERED;
# | false
property WEAKLY_CENTERED : Bool;


# @category Geometry
# The following is defined for [[CENTERED]] polytopes only:
# A facet is special if the cone over that facet with the origin as the apex contains the [[VERTEX_BARYCENTER]].
# Motivated by Obro's work on Fano polytopes.
property SPECIAL_FACETS : Set;

# @category Geometry
# True if P = -P.
# @example A centered 3-cube is centrally symmetric. By stacking a single facet (5), this property is lost. We can
# recover it by stacking the opposing facet (4) as well.
# > $p = cube(3);
# > print $p->CENTRALLY_SYMMETRIC;
# | true
# > print stack($p,5)->CENTRALLY_SYMMETRIC;
# | false
# > print stack($p,new Set<Int>(4,5))->CENTRALLY_SYMMETRIC;
# | true
property CENTRALLY_SYMMETRIC : Bool;

# @category Geometry
# The permutation induced by the central symmetry, if present.
property CS_PERMUTATION : Array<Int>;

# @category Geometry
# Number of [[POINTS]].
property N_POINTS = override N_INPUT_RAYS;

# @category Combinatorics
# Number of [[VERTICES]].
# @example The following prints the number of vertices of a 3-dimensional cube:
# > print cube(3)->N_VERTICES;
# | 8
# @example The following prints the number of vertices of the convex hull of 10 specific points lying in the unit square [0,1]^2:
# > print rand_box(2,10,1,seed=>4583572)->N_VERTICES;
# | 4
property N_VERTICES = override N_RAYS;

# @category Combinatorics
# Number of pairs of incident vertices and facets.
property N_VERTEX_FACET_INC = override N_RAY_FACET_INC;

# @category Combinatorics
# Vertex-facet incidence matrix, with rows corresponding to facets and columns
# to vertices. Vertices and facets are numbered from 0 to [[N_VERTICES]]-1 rsp.
# [[N_FACETS]]-1, according to their order in [[VERTICES]] rsp. [[FACETS]].
# 
# This property is at the core of all combinatorial properties.  It has the following semantics:
# (1) The combinatorics of an unbounded and pointed polyhedron is defined to be the combinatorics
#     of the projective closure.
# (2) The combiantorics of an unbounded polyhedron which is not pointed is defined to be the
#     combinatorics of the quotient modulo the lineality space.
# Therefore: [[VERTICES_IN_FACETS]] and each other property which is grouped under "Combinatorics"
# always refers to some polytope.
# @example [prefer cdd] [require bundled:cdd]
# The following prints the vertex-facet incidence matrix of a 5-gon by listing all facets as a set of contained vertices
# in a cyclic order (each line corresponds to an edge):
# > print n_gon(5)->VERTICES_IN_FACETS;
# | {1 2}
# | {2 3}
# | {3 4}
# | {0 4}
# | {0 1}
# @example [prefer cdd] [require bundled:cdd] The following prints the Vertex_facet incidence matrix of the standard 3-simplex together with the facet numbers:
# > print rows_numbered(simplex(3)->VERTICES_IN_FACETS);
# | 0:1 2 3
# | 1:0 2 3
# | 2:0 1 3
# | 3:0 1 2
property VERTICES_IN_FACETS = override RAYS_IN_FACETS;

# @category Geometry
property VERTICES_IN_RIDGES = override RAYS_IN_RIDGES;

# @category Geometry
# Similar to [[VERTICES_IN_FACETS]], but with columns corresponding to [[POINTS]] instead of [[VERTICES]].
# This property is a byproduct of convex hull computation algorithms.
# It is discarded as soon as [[VERTICES_IN_FACETS]] is computed.
property POINTS_IN_FACETS = override INPUT_RAYS_IN_FACETS;


# @category Combinatorics
# transposed [[VERTICES_IN_FACETS]]
# Notice that this is a temporary property; it will not be stored in any file.
property FACETS_THRU_VERTICES = override FACETS_THRU_RAYS;

# @category Geometry
# similar to [[FACETS_THRU_VERTICES]], but with [[POINTS]] instead of [[VERTICES]]
# Notice that this is a temporary property; it will not be stored in any file.
property FACETS_THRU_POINTS = override FACETS_THRU_INPUT_RAYS;

# @category Geometry
# Similar to [[VERTICES_IN_FACETS]], but with rows corresponding to [[INEQUALITIES]] instead of [[FACETS]].
# This property is a byproduct of convex hull computation algorithms.
# It is discarded as soon as [[VERTICES_IN_FACETS]] is computed.
property VERTICES_IN_INEQUALITIES = override RAYS_IN_INEQUALITIES;

# @category Geometry
# transposed [[VERTICES_IN_INEQUALITIES]]
property INEQUALITIES_THRU_VERTICES = override INEQUALITIES_THRU_RAYS;

# @category Combinatorics
# Number of incident facets for each vertex.
# @example [require bundled:cdd] The following prints the number of incident facets for each vertex of the elongated pentagonal pyramid (Johnson solid 9)
# > print johnson_solid(9)->VERTEX_SIZES;
# | 5 4 4 4 4 4 3 3 3 3 3
property VERTEX_SIZES = override RAY_SIZES;

# @category Combinatorics
# Measures the deviation of the cone from being simple in terms of the [[GRAPH]].
# @example The excess vertex degree of an egyptian pyramid is one.
# > print pyramid(cube(2))->EXCESS_VERTEX_DEGREE;
# | 1
property EXCESS_VERTEX_DEGREE = override EXCESS_RAY_DEGREE;

# @topic //properties//HASSE_DIAGRAM
# Number of incident vertices for each facet.

# @category Unbounded polyhedra
# Indices of vertices that are rays.
property FAR_FACE : Set;

rule FAR_FACE : VertexPerm.FAR_FACE, VertexPerm.PERMUTATION {
   $this->FAR_FACE=permuted($this->VertexPerm->FAR_FACE, $this->VertexPerm->PERMUTATION);
}
weight 1.10;

# @category Combinatorics
# Vertex-edge graph.

property GRAPH {

   # @category Geometry
   # Difference of the vertices for each edge (only defined up to signs).
   property EDGE_DIRECTIONS : EdgeMap<Undirected, Vector<Scalar>> : construct(ADJACENCY);

   # @category Geometry
   # Squared Euclidean length of each edge
   property SQUARED_EDGE_LENGTHS : EdgeMap<Undirected, Scalar> : construct(ADJACENCY);

}

property DUAL_GRAPH {

   # @category Geometry
   # Dihedral angles (in radians) between the two facets corresponding to
   # each edge of the dual graph, i.e. the ridges of the polytope.
   property DIHEDRAL_ANGLES : EdgeMap<Undirected, Float> : construct(ADJACENCY);

}


# @category Combinatorics
# Lists for each occurring size (= number of incident vertices or edges) of a 2-face how many there are.
# @example This prints the number of facets spanned by 3,4 or 5 vertices a truncated 3-dimensional cube has.
# > $p = truncation(cube(3),5);
# > print $p->TWO_FACE_SIZES;
# | {(3 1) (4 3) (5 3)}
property TWO_FACE_SIZES : Map<Int,Int>;

# @category Combinatorics
# Lists for each occurring size (= number of incident facets or ridges) of a subridge how many there are.
property SUBRIDGE_SIZES : Map<Int,Int>;


# @topic //properties//F_VECTOR
# f<sub>k</sub> is the number of k-faces.
# @example This prints the f-vector of a 3-dimensional cube. The first entry represents the vertices.
# > print cube(3)->F_VECTOR;
# | 8 12 6
# @example This prints the f-vector of the 3-dimensional cross-polytope. Since the cube and the cross polytope
# of equal dimension are dual, their f-vectors are the same up to reversion.
# > print cross(3)->F_VECTOR;
# | 6 12 8
# @example After truncating the first standard basis vector of the 3-dimensional cross-polytope the f-vector changes.
# Only segments of the incident edges of the cut off vertex remain and the intersection of these with the new hyperplane
# generate four new vertices. These also constitute four new edges and a new facet.
# > print truncation(cross(3),4)->F_VECTOR;
# | 9 16 9


# @topic //properties//F2_VECTOR
# f<sub>ik</sub> is the number of incident pairs of i-faces and k-faces; the main diagonal contains the [[F_VECTOR]].
# @example The following prints the f2-vector of a 3-dimensional cube:
# > print cube(3)->F2_VECTOR;
# | 8 24 24
# | 24 12 24
# | 24 24 6


# @topic //properties//SIMPLICIAL
# True if the polytope is simplicial.
# @example A polytope with random vertices uniformly distributed on the unit sphere is simplicial. The following checks
# this property and prints the result for 8 points in dimension 3:
# > print rand_sphere(3,8)->SIMPLICIAL;
# | true

# @topic //properties//SIMPLE
# True if the polytope is simple. Dual to [[SIMPLICIAL]].
# @example This determines if a 3-dimensional cube is simple or not:
# > print cube(3)->SIMPLE;
# | true

# @topic //properties//SELF_DUAL
# True if the polytope is self-dual.
# @example The following checks if the centered square with side length 2 is self dual:
# > print cube(2)->SELF_DUAL;
# | true
# @example The elongated square pyramid (Johnson solid 8) is dual to itself, since the apex of the square pyramid attachted to the cube
# and the opposing square of the cube swap roles. The following checks this property and prints the result:
# > print johnson_solid(8)->SELF_DUAL;
# | true


# @category Combinatorics
# Maximal dimension in which all faces are simplices.
# @example The 3-dimensional cross-polytope is simplicial, i.e. its simplicity is 2. After truncating an arbitrary vertex
# the simplicity is reduced to 1.
# > print cross(3)->SIMPLICIALITY;
# | 2
# > print truncation(cross(3),4)->SIMPLICIALITY;
# | 1
property SIMPLICIALITY : Int;

# @category Combinatorics
# Maximal dimension in which all dual faces are simplices.
# @example This checks the 3-dimensional cube for simplicity. Since the cube is dual to the cross-polytope of equal dimension and all its faces are simplices,
# the result is 2.
# > print cube(3)->SIMPLICITY;
# | 2
property SIMPLICITY : Int;

# @category Combinatorics
# Maximal dimension in which all faces are simple polytopes.
# This checks the 3-dimensional cube for face simplicity. Since the cube is dual to the cross-polytope of equal dimension and it is simplicial,
# the result is 3.
# > print cube(3)->SIMPLICITY;
# | 3
property FACE_SIMPLICITY : Int;

# @category Combinatorics
# True if all facets are cubes.
# @example A k-dimensional cube has k-1-dimensional cubes as facets and is therefore cubical. The following checks if this holds for the
# 3-dimensional case:
# > print cube(3)->CUBICAL;
# | true
# @example This checks if a zonotope generated by 4 random points on the 3-dimensional sphere is cubical, which is always the case.
# > print zonotope(rand_sphere(3,4)->VERTICES)->CUBICAL;
# | true
property CUBICAL : Bool;

# @category Combinatorics
# Maximal dimension in which all facets are cubes.
# @example We will modify the 3-dimensional cube in two different ways. While stacking some facets (in this case facets 4 and 5) preserves the cubicality up to
# dimension 2, truncating an arbitrary vertex reduces the cubicality to 1.
# > print stack(cube(3),[4,5])->CUBICALITY;
# | 2
# > print truncation(cube(3),5)->CUBICALITY;
# | 1
property CUBICALITY : Int;

# @category Combinatorics
# Dual to [[CUBICAL]].
# @example Since the cross-polytope is dual to a cube of same dimension, it is cocubical. The following checks this for the 3-dimensional case:
# > print cross(3)->COCUBICAL;
# | true
property COCUBICAL : Bool;

# @category Combinatorics
# Dual to [[CUBICALITY]].
# @example After stacking a facet of the 3-dimensional cube, its cubicality is lowered to 2. Hence its dual polytope has cocubicality 2 as well. The
# following produces such a stacked cube and asks for its cocubicality after polarization:
# > $p = stack(cube(3),5);
# > print polarize($p)->COCUBICALITY;
# | 2
property COCUBICALITY : Int;

# @category Combinatorics
# True if the polytope is neighborly.
# @example This checks the 4-dimensional cyclic polytope with 6 points on the moment curve for neighborliness, i.e. if it is ⌊dim/2⌋ neighborly:
# > print cyclic(4,6)->NEIGHBORLY;
# | true
property NEIGHBORLY : Bool;

# @category Combinatorics
# Maximal dimension in which all facets are neighborly.
# @example [nocompare] This determines that the full dimensional polytope given by 10 specific vertices on the 8-dimensional sphere is 3-neighborly, i.e.
# all 3-dimensional faces are tetrahedra. Hence the polytope is not neighborly.
# > print rand_sphere(8,10,seed=>8866463)->NEIGHBORLINESS;
# | 3
property NEIGHBORLINESS : Int;

# @category Combinatorics
# Dual to [[NEIGHBORLY]].
# @example Since cyclic polytopes generated by vertices on the moment curve are neighborly, their dual polytopes are balanced. The following checks this
# for the 4-dimensional case by centering the cyclic polytope and then polarizing it:
# > $p = cyclic(4,6);
# > $q = polarize(center($p));
# > print $q->BALANCED;
# | true
property BALANCED : Bool;

# @category Combinatorics
# Maximal dimension in which all facets are balanced.
# @example [nocompare] The following full dimensional polytope given by 10 specific vertices on the 8-dimensional sphere is 3-neighborly. Hence the dual polytope is
# 3-balanced, where we first center and then polarize it.
# > $p = rand_sphere(8,10,seed=>8866463);
# > $q = polarize(center($p));
# > print $q->BALANCE;
# | 3
property BALANCE : Int;

# @category Combinatorics
# Parameter describing the shape of the face-lattice of a 4-polytope.
property FATNESS : Float;

# @category Combinatorics
# Parameter describing the shape of the face-lattice of a 4-polytope.
property COMPLEXITY : Float;

# @category Geometry
# The center of gravity of the vertices of a bounded polytope.
# @example This prints the vertex barycenter of the standard 3-simplex:
# > print simplex(3)->VERTEX_BARYCENTER;
# | 1 1/4 1/4 1/4
property VERTEX_BARYCENTER : Vector<Scalar>;

# @category Geometry
# The minimal angle between any two vertices (seen from the [[VERTEX_BARYCENTER]]).
property MINIMAL_VERTEX_ANGLE : Float;

# @category Geometry
# The rows of this matrix contain a configuration of affine points in homogeneous cooordinates.
# The zonotope is obtained as the Minkowski sum of all rows, normalized to x_0 = 1.
# Thus, if the input matrix has n columns, the ambient affine dimension of the resulting zonotope is n-1.
property ZONOTOPE_INPUT_POINTS : Matrix<Scalar>;

# @category Geometry
# is the zonotope calculated from ZONOTOPE_INPUT_POINTS or ZONOTOPE_INPUT_VECTORS to be centered at the origin?
# The zonotope is always calculated as the Minkowski sum of all segments conv {x,v}, where
# * v ranges over the ZONOTOPE_INPUT_POINTS or ZONOTOPE_INPUT_VECTORS, and
# * x = -v if CENTERED_ZONOTOPE = 1,
# * x = 0  if CENTERED_ZONOTOPE = 0.
# Input section only.
property CENTERED_ZONOTOPE : Bool;

# @category Geometry
# The minimal Ball which contains a the Polytope. The ball is 
# given by its center c and the square of it radius r.
# @example
# > $P = new Polytope(POINTS=>[[1,0],[1,4]]);
# > print $P->MINIMAL_BALL;
# | 4 <1 2>
property MINIMAL_BALL : Pair<Scalar, Vector<Scalar>>;

rule MINIMAL_BALL : POINTS | VERTICES | INEQUALITIES | FACETS{
   $this->MINIMAL_BALL = minimal_ball($this);
}



} # /Polytope

# @category Geometry
# a lattice that is displaced from the origin, i.e.,
# a set of the form x + L, where x is a non-zero vector and L a (linear) lattice
# @tparam Scalar coordinate type
declare object AffineLattice<Scalar=Rational> {

# an affine point, ie first coordinate non-zero
property ORIGIN : Vector<Scalar>;

# rows are vectors, ie first coordinate zero
property BASIS : Matrix<Scalar>;


property SQUARED_DETERMINANT : Scalar;


rule SQUARED_DETERMINANT : BASIS {
   $this->SQUARED_DETERMINANT = det($this->BASIS * transpose($this->BASIS));
}

} # /AffineLattice


object Polytope {

# @category Geometry
# An affine lattice L such that P + L tiles the affine span of P
property TILING_LATTICE : AffineLattice<Scalar>;


# @topic category properties/Triangulation and volume
# Everything in this group is defined for [[BOUNDED]] polytopes only.


# @category Triangulation and volume
# Volume of the polytope.
# @example The following prints the volume of the centered 3-dimensional cube with side length 2:
# > print cube(3)->VOLUME;
# | 8
property VOLUME : Scalar;

# @category Triangulation and volume
# Mahler volume (or volume product) of the polytope.
# Defined as the volume of the polytope and the volume of its polar (for [[BOUNDED]], [[CENTERED]] and [[FULL_DIM]] polytopes only).
# Often studied for centrally symmetric convex bodies, where the regular cubes are conjectured to be the global minimiers.
# @example The following prints the Mahler volume of the centered 2-cube:
# > print cube(2)->MAHLER_VOLUME;
# | 8
property MAHLER_VOLUME : Scalar;

# @category Triangulation and volume
# Array of the squared relative //k//-dimensional volumes of the simplices in
# a triangulation of a //d//-dimensional polytope.
property SQUARED_RELATIVE_VOLUMES : Array<Scalar>;

# @category Geometry
# Centroid (center of mass) of the polytope.
property CENTROID : Vector<Scalar>;

property TRIANGULATION {

   # @category Triangulation and volume
   # GKZ-vector
   # 	See Chapter 7 in Gelfand, Kapranov, and Zelevinsky:
   # 	Discriminants, Resultants and Multidimensional Determinants, Birkhäuser 1994
   property GKZ_VECTOR : Vector<Scalar>;

}


# @category Combinatorics
# Minimal non-faces of a [[SIMPLICIAL]] polytope.
property MINIMAL_NON_FACES : IncidenceMatrix;

rule MINIMAL_NON_FACES : VertexPerm.MINIMAL_NON_FACES, VertexPerm.PERMUTATION {
   $this->MINIMAL_NON_FACES=permuted_cols($this->VertexPerm->MINIMAL_NON_FACES, $this->VertexPerm->PERMUTATION);
}
weight 1.20;


# @category Visualization
# Reordered [[VERTICES_IN_FACETS]] for 2d and 3d-polytopes.
# Vertices are listed in the order of their appearance
# when traversing the facet border counterclockwise seen from outside of the polytope.
#
# For a 2d-polytope (which is a closed polygon), lists all vertices in the border traversing order.
property VIF_CYCLIC_NORMAL = override RIF_CYCLIC_NORMAL;

# @category Visualization
# Reordered transposed [[VERTICES_IN_FACETS]]. Dual to [[VIF_CYCLIC_NORMAL]].
property FTV_CYCLIC_NORMAL = override FTR_CYCLIC_NORMAL;

# @category Visualization
# Reordered [[GRAPH]]. Dual to [[NEIGHBOR_FACETS_CYCLIC_NORMAL]].
property NEIGHBOR_VERTICES_CYCLIC_NORMAL = override NEIGHBOR_RAYS_CYCLIC_NORMAL;


# @category Visualization
# Unique names assigned to the [[VERTICES]].
# If specified, they are shown by visualization tools instead of vertex indices.
#
# For a polytope build from scratch, you should create this property by yourself,
# either manually in a text editor, or with a client program.
#
# If you build a polytope with a construction function
# taking some other input polytope(s), the labels are created the labels automatically
# except if you call the function with a //no_labels// option. The exact format of the
# abels is dependent on the construction, and is described in the corresponding help topic.
# @option Bool no_labels
property VERTEX_LABELS = override RAY_LABELS;

# @category Visualization
# Unique names assigned to the [[POINTS]], analogous to [[VERTEX_LABELS]].
property POINT_LABELS = override INPUT_RAY_LABELS;

# @category Geometry
# Find the vertices by given labels.
# @param String label ... vertex labels
# @return Set<Int> vertex indices
user_method labeled_vertices {
   my $this=shift;
   if (defined (my $labels=$this->lookup("VERTEX_LABELS"))) {
      new Set([ find_labels($labels, @_) ]);
   } else {
      die "No VERTEX_LABELS stored in this complex";
   }
}

# @topic //properties//FACET_LABELS
# Unique names assigned to the [[FACETS]], analogous to [[VERTEX_LABELS]].

# @topic //properties//INEQUALITY_LABELS
# Unique names assigned to the [[INEQUALITIES]], analogous to [[VERTEX_LABELS]].

# @category Geometry
# The splits of the polytope, i.e., hyperplanes cutting the polytope in
# two parts such that we have a regular subdivision.
property SPLITS : Matrix<Scalar>;

# @category Geometry
# Two [[SPLITS]] are compatible if the defining hyperplanes do not intersect in the
# interior of the polytope.  This defines a graph.
property SPLIT_COMPATIBILITY_GRAPH : Graph;

# @topic category properties/Matroid properties
# Properties which belong to the corresponding (oriented) matroid

# @category Matroid properties
# Chirotope corresponding to the [[VERTICES]]. TOPCOM format.
# @depends topcom
property CHIROTOPE : Text;

rule CHIROTOPE : VERTICES {
   $this->CHIROTOPE = chirotope($this->VERTICES);
}
precondition : FULL_DIM && BOUNDED;
weight 6.11;

# @category Matroid properties
# Circuits in [[VECTORS]]
# @depends topcom
property CIRCUITS : Set<Pair<Set<Int>,Set<Int>>>;

# @category Matroid properties
# Cocircuits in [[VECTORS]]
# @depends topcom
property COCIRCUITS : Set<Pair<Set<Int>,Set<Int>>>;

# @category Geometry
# The slack matrix of the polytope. The (i,j)-th entry is the value of the j-th
# facet on the i-th vertex.
#
# See
# > João Gouveia, Antonio Macchia, Rekha R. Thomas, Amy Wiebe:
# > The Slack Realization Space of a Polytope
# > (https://arxiv.org/abs/1708.04739)
property SLACK_MATRIX : Matrix<Scalar>;

rule SLACK_MATRIX : VERTICES, FACETS {
   $this->SLACK_MATRIX = $this->VERTICES * transpose($this->FACETS);
}

# @category Geometry
# Returns the inner description of a Polytope:
# [V,L] where V are the vertices and L is the lineality space
# @return Array<Matrix<Scalar>>
user_method INNER_DESCRIPTION : VERTICES, LINEALITY_SPACE {
    my $this=shift;
    new Array<Matrix>([$this->VERTICES, $this->LINEALITY_SPACE]);
}

# @category Geometry
# Returns the outer description of a Polytope:
# [F,A] where F are the facets and A is the affine hull
# @return Array<Matrix<Scalar>>
user_method OUTER_DESCRIPTION : FACETS, AFFINE_HULL {
    my $this=shift;
    new Array<Matrix>([$this->FACETS, $this->AFFINE_HULL]);
}

# @category Geometry
# checks if a given Ball B(c,r) is contained in a Polytope
# @param Vector<Scalar> c the center of the ball
# @param Scalar r the radius of the ball
# @return Bool
user_method contains_ball {
    my ($self, $c, $r) = @_;
    return polytope_contains_ball($c, $r, $self);
}

# @category Geometry
# check if a given Ball B(c,r) contains a Polytope
# @param Vector<Scalar> c the center of the ball
# @param Scalar r the radius of the ball
# @return Bool
user_method contained_in_ball {
    my  ($self, $c, $r) = @_;
    return polytope_contained_in_ball($self, $c, $r);
}

}

# @category Coordinate conversions
# provide a Polytope object with desired coordinate type
# @tparam Coord target coordinate type
# @param Polytope P source object
# @return Polytope<Coord> //P// if it already has the requested type, a new object otherwise
# @example
# > print cube(2)->type->full_name;
# | Polytope<Rational>
# > $pf = convert_to<Float>(cube(2));
# > print $pf->type->full_name;
# | Polytope<Float>

user_function convert_to<Coord>(Polytope) {
   my $target_type=typeof Coord;
   if ($_[0]->type->params->[0] == $target_type) {
      $_[0]
   } else {
      new Polytope<Coord>($_[0]);
   }
}

# @category Coordinate conversions
# Dehomogenize the [[VERTICES|vertex coordinates]] and convert them to Float
# @param Polytope P source object
# @return Matrix<Float> the dehomogenized vertices converted to Float
# @example
# > print cube(2,1/2)->VERTICES;
# | 1 -1/2 -1/2
# | 1 1/2 -1/2
# | 1 -1/2 1/2
# | 1 1/2 1/2
# > print affine_float_coords(cube(2,1/2));
# | -0.5 -0.5
# | 0.5 -0.5
# | -0.5 0.5
# | 0.5 0.5
user_function affine_float_coords(Polytope) {
   my $P=shift;
   if ($P->FAR_FACE->size()!=0) {
      croak("polytope must be bounded");
   }
   dehomogenize(convert_to<Float>($P->VERTICES));
}

# @category Combinatorics
# Provide a basic combinatorial description (unless there exists one already).
# Only the options //extended// and //matchthenet// (for 3-polytopes) can trigger any computation.
# @param Polytope P source object
# @return String the description
# @option Bool extended
# @option Bool matchthenet
# @option String language
# @example
# > print describe(cube(3));
# | cube of dimension 3
# > $P = new Polytope(POINTS=>cube(3)->VERTICES); print describe($P);
# | polytope with POINTS, CONE_AMBIENT_DIM
# > $P = new Polytope(POINTS=>cube(3)->VERTICES); print describe($P, extended=>1, language=>"de");
# | einfaches 3-Polytop im 3-Raum mit f-Vektor 8 12 6
# > $P = new Polytope(POINTS=>cube(3)->VERTICES); print describe($P, matchthenet=>1);
# | { "en": "A simple polytope with 8 vertices, 12 edges and 6 facets.", "de": "Ein einfaches Polytop mit 8 Ecken, 12 Kanten und 6 Seiten."}
user_function describe(Polytope; { extended => 0, matchthenet => 0, language => "en" }) {
   my ($P, $options) = @_;
   my ($simplicial, $simple, $cubical, $polytope, $unbounded, $polyhedron, $pointed, $with, $dimensional, $linealityspace, $in, $space, $fvector, $vertices, $edges, $facets, $and)
     = ("simplicial", "simple", "cubical", "polytope", "unbounded", "polyhedron", "pointed", "with", "dimensional", "lineality space", "in", "space", "f-vector", "vertices", "edges", "facets", "and");
   if ($options->{language} eq "de") {
       ($simplicial, $simple, $cubical, $polytope, $unbounded, $polyhedron, $pointed, $with, $dimensional, $linealityspace, $in, $space, $fvector, $vertices, $edges, $facets)
           = ("simpliziales", "einfaches", "kubisches", "Polytop", "unbeschränktes", "Polyeder", "spitzes", "mit", "dimensionalem", "Linealitätsraum", "im", "Raum", "f-Vektor", "Ecken", "Kanten", "Facetten", "und");
   }
   if ($P->description() eq "") {
      my $desc = "";
      if ($options->{extended}) {
         my $ad = $P->AMBIENT_DIM;
         my $d = $P->DIM;
         my $fv = $P->F_VECTOR;
         if ($P->SIMPLICIAL) {
            $desc = "$simplicial ";
         } elsif ($P->SIMPLE) {
            $desc = "$simple ";
         } elsif ($P->CUBICAL) {
            $desc = "$cubical ";
         }
         if ($P->BOUNDED) {
            $desc .= "$d-$polytope";
         } else {
            $desc .= "$unbounded $d-$polyhedron";
            my $ld = $P->LINEALITY_DIM;
            if ($ld==0) {
               $desc = "$pointed " . $desc;
            } else {
               $desc .= ", $with $ld-$dimensional $linealityspace,";
            }
         }
         $desc .= " $in $ad-$space $with $fvector $fv";
      } elsif ($options->{matchthenet}) {
         my $desc_en = "";
         my $desc_de = "";
         my $ad = $P->AMBIENT_DIM;
         my $d = $P->DIM;
         my ($n,$e,$f) = @{ $P->F_VECTOR };
         if ($P->SIMPLICIAL) {
            $desc_en = "simplicial ";
            $desc_de = "simpliziales ";
         } elsif ($P->SIMPLE) {
            $desc_en = "simple ";
            $desc_de = "einfaches ";
         } elsif ($P->CUBICAL) {
            $desc_en = "simple ";
            $desc_de = "kubisches ";
         }
         die "option 'matchthenet' only valid for 3-polytopes" unless $ad==3 and $d==3;
         $desc_en .= "polytope with $n vertices, $e edges and $f facets";
         $desc_de .= "Polytop mit $n Ecken, $e Kanten und $f Seiten";
         $desc = "{ \"en\": \"A $desc_en.\", \"de\": \"Ein $desc_de.\"}";
         } else {
         $desc = "$polytope $with " . (join ", ", $P->list_properties());
      }
      $P->description() = $desc;
   }
   return $P->description();
}

# Local Variables:
# mode: perl
# cperl-indent-level:3
# indent-tabs-mode:nil
# End: