Trac #28506: Direct sum of polyhedron is broken, so is minkowski diff…

…erence and face truncation

The following returns a `TypeError`:

{{{
sage: s2 = polytopes.simplex(2)
sage: s3 = polytopes.simplex(3)
sage: t = s2.direct_sum(s3)
}}}

H-Polyhedra might have non-integral vertices and it is hard to tell from
the H-Representation, whether this this is the case.

We fix `direct_sum`, `minkowski_difference` and `face_truncation` to try
the quotient ring, if the given base ring does not work.

URL: https://trac.sagemath.org/28506
Reported by: jipilab
Ticket author(s): Jonathan Kliem
Reviewer(s): Jean-Philippe Labbé

src/doc/en/thematic_tutorials/geometry/tips.rst

@@ -27,7 +27,7 @@ You can obtain different operations using natural symbols:
     sage: Cube * Octahedron   # Cartesian product
     A 6-dimensional polyhedron in QQ^6 defined as the convex hull of 48 vertices
     sage: Cube - Polyhedron(vertices=[[-1,0,0],[1,0,0]])  # Minkowski difference
-    A 2-dimensional polyhedron in ZZ^3 defined as the convex hull of 4 vertices
+    A 2-dimensional polyhedron in QQ^3 defined as the convex hull of 4 vertices
 
 .. end of output
 
@@ -64,19 +64,19 @@ the latex presentation, there is a method for that!
 
     sage: Nice_repr = TCube.Hrepresentation_str()
     sage: print(Nice_repr)
-                  -2*x0 >= -1 
-                  -2*x1 >= -1 
-    -6*x0 + 6*x1 - 6*x2 >= -7 
-                   2*x0 >= -1 
-                   2*x1 >= -1 
-     6*x0 - 6*x1 - 6*x2 >= -7 
-     6*x0 + 6*x1 - 6*x2 >= -7 
-     6*x0 - 6*x1 + 6*x2 >= -7 
-     6*x0 + 6*x1 + 6*x2 >= -7 
-                   2*x2 >= -1 
-                  -2*x2 >= -1 
-    -6*x0 + 6*x1 + 6*x2 >= -7 
-    -6*x0 - 6*x1 + 6*x2 >= -7 
+                  -2*x0 >= -1
+                  -2*x1 >= -1
+    -6*x0 + 6*x1 - 6*x2 >= -7
+                   2*x0 >= -1
+                   2*x1 >= -1
+     6*x0 - 6*x1 - 6*x2 >= -7
+     6*x0 + 6*x1 - 6*x2 >= -7
+     6*x0 - 6*x1 + 6*x2 >= -7
+     6*x0 + 6*x1 + 6*x2 >= -7
+                   2*x2 >= -1
+                  -2*x2 >= -1
+    -6*x0 + 6*x1 + 6*x2 >= -7
+    -6*x0 - 6*x1 + 6*x2 >= -7
     -6*x0 - 6*x1 - 6*x2 >= -7
 
     sage: print(TCube.Hrepresentation_str(latex=True))
@@ -94,7 +94,7 @@ the latex presentation, there is a method for that!
                               -2 \, x_{2} & \geq & -1 \\
     -6 \, x_{0} + 6 \, x_{1} + 6 \, x_{2} & \geq & -7 \\
     -6 \, x_{0} - 6 \, x_{1} + 6 \, x_{2} & \geq & -7 \\
-    -6 \, x_{0} - 6 \, x_{1} - 6 \, x_{2} & \geq & -7 
+    -6 \, x_{0} - 6 \, x_{1} - 6 \, x_{2} & \geq & -7
     \end{array}
 
     sage: Latex_repr = LatexExpr(TCube.Hrepresentation_str(latex=True))

src/sage/geometry/polyhedron/base.py

@@ -3735,7 +3735,7 @@ def minkowski_difference(self, other):
             sage: four_cube = polytopes.hypercube(4)
             sage: four_simplex = Polyhedron(vertices = [[0, 0, 0, 1], [0, 0, 1, 0], [0, 1, 0, 0], [1, 0, 0, 0]])
             sage: four_cube - four_simplex
-            A 4-dimensional polyhedron in ZZ^4 defined as the convex hull of 16 vertices
+            A 4-dimensional polyhedron in QQ^4 defined as the convex hull of 16 vertices
             sage: four_cube.minkowski_difference(four_simplex) == four_cube - four_simplex
             True
 
@@ -3762,6 +3762,13 @@ def minkowski_difference(self, other):
             True
             sage: (X-Y)+Y == X
             True
+
+        Testing that :trac:`28506` is fixed::
+
+            sage: Q = Polyhedron([[1,0],[0,1]])
+            sage: S = Polyhedron([[0,0],[1,2]])
+            sage: S.minkowski_difference(Q)
+            A 1-dimensional polyhedron in QQ^2 defined as the convex hull of 2 vertices
         """
         if other.is_empty():
             return self.parent().universe()   # empty intersection = everything
@@ -3780,7 +3787,9 @@ def minkowski_difference(self, other):
             ieq = list(ieq)
             ieq[0] += min(values)   # shift constant term
             new_ieqs.append(ieq)
-        P = self.parent()
+
+        # Some vertices might need fractions.
+        P = self.parent().change_ring(self.base_ring().fraction_field())
         return P.element_class(P, None, [new_ieqs, new_eqns])
 
     def __sub__(self, other):
@@ -4125,17 +4134,27 @@ def direct_sum(self, other):
             :meth:`join`
             :meth:`subdirect_sum`
 
-        TESTS::
+        TESTS:
+
+        Check that the backend is preserved::
 
             sage: P = polytopes.simplex(backend='cdd')
             sage: Q = polytopes.simplex(backend='ppl')
             sage: P.direct_sum(Q).backend()
             'cdd'
             sage: Q.direct_sum(P).backend()
             'ppl'
+
+        Check that :trac:`28506` is fixed::
+
+            sage: s2 = polytopes.simplex(2)
+            sage: s3 = polytopes.simplex(3)
+            sage: s2.direct_sum(s3)
+            A 5-dimensional polyhedron in QQ^7 defined as the convex hull of 7 vertices
         """
         try:
-            new_ring = self.parent()._coerce_base_ring(other)
+            # Some vertices might need fractions.
+            new_ring = self.parent()._coerce_base_ring(other).fraction_field()
         except TypeError:
             raise TypeError("no common canonical parent for objects with parents: " + str(self.parent())
                      + " and " + str(other.parent()))
@@ -4557,12 +4576,22 @@ def face_truncation(self, face, linear_coefficients=None, cut_frac=None):
              sage: face_trunc.face_lattice().is_isomorphic(Cube.face_lattice())
              True
 
-        TESTS::
+        TESTS:
+
+        Testing that the backend is preserved::
 
             sage: Cube = polytopes.cube(backend='field')
             sage: face_trunc = Cube.face_truncation(Cube.faces(2)[0])
             sage: face_trunc.backend()
             'field'
+
+        Testing that :trac:`28506` is fixed::
+
+            sage: P = polytopes.twenty_four_cell()
+            sage: P = P.dilation(6)
+            sage: P = P.change_ring(ZZ)
+            sage: P.face_truncation(P.faces(2)[0], cut_frac=1)
+            A 4-dimensional polyhedron in QQ^4 defined as the convex hull of 27 vertices
         """
         if cut_frac is None:
             cut_frac = ZZ.one() / 3
@@ -4599,7 +4628,8 @@ def face_truncation(self, face, linear_coefficients=None, cut_frac=None):
         new_ieqs = self.inequalities_list() + [ineq_vector]
         new_eqns = self.equations_list()
 
-        parent = self.parent().base_extend(cut_frac)
+        # Some vertices might need fractions.
+        parent = self.parent().base_extend(cut_frac/1)
         return parent.element_class(parent, None, [new_ieqs, new_eqns])
 
     def stack(self, face, position=None):

src/sage/geometry/polyhedron/base_ZZ.py

@@ -764,6 +764,7 @@ def is_known_summand(poly):
             if is_known_summand(X):
                 continue
             Y = self - X
+            Y = Y.change_ring(ZZ)  # Minkowski difference returns QQ-polyhedron
             if X+Y != self:
                 continue
             decompositions.append((X, Y))