Trac #14222: Coercion/conversion from Cone and other objects to Polyh…

…edron

Mixed intersections between Cones and Polyhedron returning a Polyhedron:

   {{{
      sage: cone = Cone([(1,0), (-1,0), (0,1)])
      sage: p = polytopes.hypercube(2)
      sage: p & cone
      A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 3
vertices
   }}}

In a similar direction, we extend the `Polyhedron` constructor so it can
take various types as objects as input, for example:
{{{
        sage: H.<x,y> = HyperplaneArrangements(QQ)
        sage: h = x + y - 1; h
        Hyperplane x + y - 1
        sage: Polyhedron(h, base_ring=ZZ)
        A 1-dimensional polyhedron in ZZ^2 defined as the convex hull of
1 vertex and 1 line
        sage: Polyhedron(h)
}}}

URL: https://trac.sagemath.org/14222
Reported by: nthiery
Ticket author(s): Matthias Koeppe
Reviewer(s): Jonathan Kliem

src/sage/geometry/cone.py

@@ -3017,6 +3017,13 @@ def intersection(self, other):
             N( 2, 5)
             in 2-d lattice N
 
+        The intersection can also be expressed using the operator ``&``::
+
+            sage: (cone1 & cone2).rays()
+            N(-1, 3),
+            N( 2, 5)
+            in 2-d lattice N
+
         It is OK to intersect cones living in sublattices of the same ambient
         lattice::
 
@@ -3042,7 +3049,18 @@ def intersection(self, other):
             N(3, 1),
             N(0, 1)
             in 2-d lattice N
+
+        An intersection with a polyhedron returns a polyhedron::
+
+            sage: cone = Cone([(1,0), (-1,0), (0,1)])
+            sage: p = polytopes.hypercube(2)
+            sage: cone & p
+            A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 4 vertices
+            sage: sorted(_.vertices_list())
+            [[-1, 0], [-1, 1], [1, 0], [1, 1]]
         """
+        if not isinstance(other, ConvexRationalPolyhedralCone):
+            return self.polyhedron().intersection(other)
         if self._ambient is other._ambient:
             # Cones of the same ambient cone or fan intersect nicely/quickly.
             # Can we maybe even return an element of the cone lattice?..
@@ -3058,6 +3076,8 @@ def intersection(self, other):
         p.add_constraints(other._PPL_cone().constraints())
         return _Cone_from_PPL(p, self.lattice().intersection(other.lattice()))
 
+    __and__ = intersection
+
     def is_equivalent(self, other):
         r"""
         Check if ``self`` is "mathematically" the same as ``other``.
@@ -3495,7 +3515,7 @@ def plot(self, **options):
         result += tp.plot_walls(walls)
         return result
 
-    def polyhedron(self):
+    def polyhedron(self, **kwds):
         r"""
         Return the polyhedron associated to ``self``.
 
@@ -3523,7 +3543,7 @@ def polyhedron(self):
             A 0-dimensional polyhedron in ZZ^2 defined as the convex hull
             of 1 vertex
         """
-        return Polyhedron(rays=self.rays(), vertices=[self.lattice()(0)])
+        return Polyhedron(rays=self.rays(), vertices=[self.lattice()(0)], **kwds)
 
     def an_affine_basis(self):
         r"""

src/sage/geometry/hyperplane_arrangement/hyperplane.py

@@ -282,7 +282,7 @@ def __contains__(self, q):
         return self.A() * q + self._const == 0
 
     @cached_method
-    def polyhedron(self):
+    def polyhedron(self, **kwds):
         """
         Return the hyperplane as a polyhedron.
 
@@ -304,8 +304,10 @@ def polyhedron(self):
              A vertex at (0, 0, 4/3))
         """
         from sage.geometry.polyhedron.constructor import Polyhedron
-        R = self.parent().base_ring()
-        return Polyhedron(eqns=[self.coefficients()], base_ring=R)
+        R = kwds.pop('base_ring', None)
+        if R is None:
+            R = self.parent().base_ring()
+        return Polyhedron(eqns=[self.coefficients()], base_ring=R, **kwds)
 
     @cached_method
     def linear_part(self):

src/sage/geometry/lattice_polytope.py

@@ -3649,7 +3649,7 @@ def plot3d(self,
                 pplot += text3d(i+self.nvertices(), bc+index_shift*(p-bc), rgbcolor=pindex_color)
         return pplot
 
-    def polyhedron(self):
+    def polyhedron(self, **kwds):
         r"""
         Return the Polyhedron object determined by this polytope's vertices.
 
@@ -3660,7 +3660,7 @@ def polyhedron(self):
             A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 4 vertices
         """
         from sage.geometry.polyhedron.constructor import Polyhedron
-        return Polyhedron(vertices=[list(v) for v in self._vertices])
+        return Polyhedron(vertices=[list(v) for v in self._vertices], **kwds)
 
     def show3d(self):
         """

src/sage/geometry/polyhedron/constructor.py

@@ -291,6 +291,8 @@
 #                  https://www.gnu.org/licenses/
 ########################################################################
 
+import sage.geometry.abc
+
 from sage.rings.integer_ring import ZZ
 from sage.rings.real_double import RDF
 from sage.rings.real_mpfr import RR
@@ -313,11 +315,16 @@ def Polyhedron(vertices=None, rays=None, lines=None,
 
     INPUT:
 
-    - ``vertices`` -- list of points. Each point can be specified as
+    - ``vertices`` -- iterable of points. Each point can be specified as
       any iterable container of ``base_ring`` elements. If ``rays`` or
       ``lines`` are specified but no ``vertices``, the origin is
       taken to be the single vertex.
 
+      Instead of vertices, the first argument can also be an object
+      that can be converted to a :func:`Polyhedron` via an :meth:`as_polyhedron`
+      or :meth:`polyhedron` method. In this case, the following 5 arguments
+      cannot be provided.
+
     - ``rays`` -- list of rays. Each ray can be specified as any
       iterable container of ``base_ring`` elements.
 
@@ -332,17 +339,17 @@ def Polyhedron(vertices=None, rays=None, lines=None,
       any iterable container of ``base_ring`` elements. An entry equal to
       ``[-1,7,3,4]`` represents the equality `7x_1+3x_2+4x_3= 1`.
 
+    - ``ambient_dim`` -- integer. The ambient space dimension. Usually
+      can be figured out automatically from the H/Vrepresentation
+      dimensions.
+
     - ``base_ring`` -- a sub-field of the reals implemented in
       Sage. The field over which the polyhedron will be defined. For
       ``QQ`` and algebraic extensions, exact arithmetic will be
       used. For ``RDF``, floating point numbers will be used. Floating
       point arithmetic is faster but might give the wrong result for
       degenerate input.
 
-    - ``ambient_dim`` -- integer. The ambient space dimension. Usually
-      can be figured out automatically from the H/Vrepresentation
-      dimensions.
-
     - ``backend`` -- string or ``None`` (default). The backend to use. Valid choices are
 
       * ``'cdd'``: use cdd
@@ -465,17 +472,45 @@ def Polyhedron(vertices=None, rays=None, lines=None,
         ...
         ValueError: invalid base ring
 
-    Create a mutable polyhedron::
-
-        sage: P = Polyhedron(vertices=[[0, 1], [1, 0]], mutable=True)
-        sage: P.is_mutable()
-        True
-        sage: hasattr(P, "_Vrepresentation")
-        False
-        sage: P.Vrepresentation()
-        (A vertex at (0, 1), A vertex at (1, 0))
-        sage: hasattr(P, "_Vrepresentation")
-        True
+    Converting from a given polyhedron::
+
+        sage: cb = polytopes.cube(); cb
+        A 3-dimensional polyhedron in ZZ^3 defined as the convex hull of 8 vertices
+        sage: Polyhedron(cb, base_ring=QQ)
+        A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 8 vertices
+
+    Converting from other objects to a polyhedron::
+
+        sage: quadrant = Cone([(1,0), (0,1)])
+        sage: Polyhedron(quadrant)
+        A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 1 vertex and 2 rays
+        sage: Polyhedron(quadrant, base_ring=QQ)
+        A 2-dimensional polyhedron in QQ^2 defined as the convex hull of 1 vertex and 2 rays
+
+        sage: o = lattice_polytope.cross_polytope(2)
+        sage: Polyhedron(o)
+        A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 4 vertices
+        sage: Polyhedron(o, base_ring=QQ)
+        A 2-dimensional polyhedron in QQ^2 defined as the convex hull of 4 vertices
+
+        sage: p = MixedIntegerLinearProgram(solver='PPL')
+        sage: x, y = p['x'], p['y']
+        sage: p.add_constraint(x <= 1)
+        sage: p.add_constraint(x >= -1)
+        sage: p.add_constraint(y <= 1)
+        sage: p.add_constraint(y >= -1)
+        sage: Polyhedron(o)
+        A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 4 vertices
+        sage: Polyhedron(o, base_ring=QQ)
+        A 2-dimensional polyhedron in QQ^2 defined as the convex hull of 4 vertices
+
+        sage: H.<x,y> = HyperplaneArrangements(QQ)
+        sage: h = x + y - 1; h
+        Hyperplane x + y - 1
+        sage: Polyhedron(h, base_ring=ZZ)
+        A 1-dimensional polyhedron in ZZ^2 defined as the convex hull of 1 vertex and 1 line
+        sage: Polyhedron(h)
+        A 1-dimensional polyhedron in QQ^2 defined as the convex hull of 1 vertex and 1 line
 
     .. NOTE::
 
@@ -486,7 +521,6 @@ def Polyhedron(vertices=None, rays=None, lines=None,
         input data - the results can depend upon the tolerance
         setting of cdd.
 
-
     TESTS:
 
     Check that giving ``float`` input gets converted to ``RDF`` (see :trac:`22605`)::
@@ -569,10 +603,53 @@ def Polyhedron(vertices=None, rays=None, lines=None,
         sage: Polyhedron(ambient_dim=2, vertices=[], rays=[], lines=[], base_ring=QQ)
         The empty polyhedron in QQ^2
 
+    Create a mutable polyhedron::
+
+        sage: P = Polyhedron(vertices=[[0, 1], [1, 0]], mutable=True)
+        sage: P.is_mutable()
+        True
+        sage: hasattr(P, "_Vrepresentation")
+        False
+        sage: P.Vrepresentation()
+        (A vertex at (0, 1), A vertex at (1, 0))
+        sage: hasattr(P, "_Vrepresentation")
+        True
+
     .. SEEALSO::
 
         :mod:`Library of polytopes <sage.geometry.polyhedron.library>`
     """
+    # Special handling for first argument, for coercion-like uses
+    constructor = None
+    first_arg = vertices
+    if isinstance(first_arg, sage.geometry.abc.Polyhedron):
+        constructor = first_arg.change_ring
+    try:
+        # PolyhedronFace.as_polyhedron (it also has a "polyhedron" method with a different purpose)
+        constructor = first_arg.as_polyhedron
+    except AttributeError:
+        try:
+            # ConvexRationalPolyhedralCone, LatticePolytopeClass, MixedIntegerLinearProgram, Hyperplane
+            constructor = first_arg.polyhedron
+        except AttributeError:
+            pass
+    if constructor:
+        if not all(x is None for x in (rays, lines, ieqs, eqns, ambient_dim)):
+            raise ValueError('if a polyhedron is given, cannot provide H- and V-representations objects')
+        # Only pass non-default arguments
+        kwds = {}
+        if base_ring is not None:
+            kwds['base_ring'] = base_ring
+        if verbose is not False:
+            kwds['verbose'] = verbose
+        if backend is not None:
+            kwds['backend'] = backend
+        if minimize is not True:
+            kwds['minimize'] = minimize
+        if mutable is not False:
+            kwds['mutable'] = mutable
+        return constructor(**kwds)
+
     got_Vrep = not ((vertices is None) and (rays is None) and (lines is None))
     got_Hrep = not ((ieqs is None) and (eqns is None))
 

src/sage/geometry/polyhedron/face.py

@@ -750,7 +750,7 @@ def is_compact(self):
                        for V in self.ambient_Vrepresentation())
 
     @cached_method
-    def as_polyhedron(self):
+    def as_polyhedron(self, **kwds):
         """
         Return the face as an independent polyhedron.
 
@@ -774,7 +774,12 @@ def as_polyhedron(self):
         P = self._polyhedron
         parent = P.parent()
         Vrep = (self.vertices(), self.rays(), self.lines())
-        return P.__class__(parent, Vrep, None)
+        result = P.__class__(parent, Vrep, None)
+        if any(kwds.get(kwd) is not None
+               for kwd in ('base_ring', 'backend')):
+            from .constructor import Polyhedron
+            return Polyhedron(result, **kwds)
+        return result
 
     def _some_elements_(self):
         r"""