Is there reliable polyhedra library?

Question

I tried to use qhalf, but found bug in its output for 8 dimensions. Parma polyhedra library seemed to be more credible, but it gives wrong results even in 6 dimensions. Input central hyperplanes for chamber:

[-1, 1, 1, 1, 1, -1], 
[-1, -1, 1, 1, 1, 1], 
[1, -1, -1, 1, 1, 1], 
[1, 1, -1, -1, 1, 1], 
[1, 1, 1, -1, -1, 1], 
[1, 1, 1, 1, -1, -1]

Facets are the same, and polyhedron.vrep() is equal to

[[ 0  1 -1  1 -1  1 -1]
 [ 1  0  0  0  0  0  0]
 [ 0  1  0  1  0  1 -1]
 [ 0  1  0  0  1  0  0]
 [ 0  1  0  1  0  0  0]
 [ 0  1  0  1 -1  1  0]
 [ 0  0  0  1  0  0  1]
 [ 0  0  0  1  0  1  0]
 [ 0  1  0  1 -1  2 -1]
 [ 0  0  0  0  1  0  1]
 [ 0  1  0  0  0  1  0]]

which is clearly wrong, because vertices should be symmetric relative to cyclic change of coordinate axes.
What library can you suggest which is more reliable?

Answer 1

No, there is no bug in PPL (Parma Polyhedral Library) here. Here is why.

Let $A$ be the first of your matrices. Note that for $\ell=(1, -1, 1, -1, 1, -1)$ one has $A\ell=0$, and in fact $\ell$ spans the kernel of $A$. Thus your polyhedron is the Minkowski sum of the line spanned by $\ell$ and pointed cone $C$ given by the 12 inequalities, 6 of the form $\langle a_i,x\rangle\geq 0$ coming from the rows $a_i$ of $A$, and 6 of the form $x_j\geq 0$, for $1\leq j\leq 6$.

To enter your data in a typical polyhedral solver in "H-representation", you'd have to add the 1st column of 0s to $A$- as the rows are interpreted as inequalities, with the 1st entry of each row being the right hand side of the inequality. That is, as follows:

A=[[0, -1, 1, 1, 1, 1, -1],
   [0, -1, -1, 1, 1, 1, 1],
   [0, 1, -1, -1, 1, 1, 1],
   [0, 1, 1, -1, -1, 1, 1],
   [0, 1, 1, 1, -1, -1, 1],
   [0, 1, 1, 1, 1, -1, -1]]

Now, I'll give the input for Sagemath system, where you can in fact choose from a number of polyhedral backends (the default is PPL, one you accuse of being buggy.)

sage: p=Polyhedron(ieqs=A); p
A 6-dimensional polyhedron in QQ^6 defined as the convex hull of 1 vertex, 9 rays, 1 line
sage: p.lines()
(A line in the direction (1, -1, 1, -1, 1, -1),)
sage: p.vertices()
(A vertex at (0, 0, 0, 0, 0, 0),)
sage: p.rays()
(A ray in the direction (1, 0, 1, 0, 1, -1),
 A ray in the direction (1, 0, 0, 1, 0, 0),
 A ray in the direction (1, 0, 1, 0, 0, 0),
 A ray in the direction (1, 0, 1, -1, 1, 0),
 A ray in the direction (0, 0, 1, 0, 0, 1),
 A ray in the direction (0, 0, 1, 0, 1, 0),
 A ray in the direction (1, 0, 1, -1, 2, -1),
 A ray in the direction (0, 0, 0, 1, 0, 1),
 A ray in the direction (1, 0, 0, 0, 1, 0))

At this point one might say: "oh, the output is wrong, as the cyclic symmetry is lost!". But in fact it is OK, as each ray is only defined up to a linear combination with $\ell$; and indeed you can see that the 1st ray $(1, 0, 1, 0, 1, -1)$ becomes $(0,1,0,1,0,0)$ if one subtracts $\ell$, similarly the 4th becomes $(0,1,0,0,0,1)$ and the 7th becomes $(0,1,0,0,1,0)$. That is, one can represent all the 9 rays as 0-1 vectors with just two 1s, and the cyclic symmetry is still there; there is no bug here, it's just the way the output is encoded is a bit counter-intuitive. In fact, one can check that this is all correct, by adding the inequalities specifying nonnegativity of the variables to $A$, computing the corresponding pointed cone $C$, and the Minkowski sum of $C$ and $\ell$, as follows:

sage: B=matrix(A).stack(matrix.zero(1,6).stack(matrix.identity(6)).T)
sage: pp=Polyhedron(ieqs=B); pp; pp.rays()
A 6-dimensional polyhedron in QQ^6 defined as the convex hull of 1 vertex and 9 rays
(A ray in the direction (0, 0, 0, 1, 0, 1),
 A ray in the direction (0, 0, 1, 0, 1, 0),
 A ray in the direction (1, 0, 0, 0, 1, 0),
 A ray in the direction (0, 1, 0, 0, 1, 0),
 A ray in the direction (0, 0, 1, 0, 0, 1),
 A ray in the direction (0, 1, 0, 0, 0, 1),
 A ray in the direction (1, 0, 1, 0, 0, 0),
 A ray in the direction (0, 1, 0, 1, 0, 0),
 A ray in the direction (1, 0, 0, 1, 0, 0))
sage: pp.Minkowski_sum(Polyhedron(lines=p.lines())).inequalities()
(An inequality (-1, -1, 1, 1, 1, 1) x + 0 >= 0,
 An inequality (1, 1, -1, -1, 1, 1) x + 0 >= 0,
 An inequality (1, -1, -1, 1, 1, 1) x + 0 >= 0,
 An inequality (-1, 1, 1, 1, 1, -1) x + 0 >= 0,
 An inequality (1, 1, 1, 1, -1, -1) x + 0 >= 0,
 An inequality (1, 1, 1, -1, -1, 1) x + 0 >= 0)

And the latter is exactly the original 6 inequalities.