Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
r""" Lattice and reflexive polytopes
This module provides tools for work with lattice and reflexive polytopes. A *convex polytope* is the convex hull of finitely many points in `\RR^n`. The dimension `n` of a polytope is the smallest `n` such that the polytope can be embedded in `\RR^n`.
A *lattice polytope* is a polytope whose vertices all have integer coordinates.
If `L` is a lattice polytope, the dual polytope of `L` is
.. MATH::
\{y \in \ZZ^n : x\cdot y \geq -1 \text{ all } x \in L\}
A *reflexive polytope* is a lattice polytope, such that its polar is also a lattice polytope, i.e. it is bounded and has vertices with integer coordinates.
This Sage module uses Package for Analyzing Lattice Polytopes (PALP), which is a program written in C by Maximilian Kreuzer and Harald Skarke, which is freely available under the GNU license terms at http://hep.itp.tuwien.ac.at/~kreuzer/CY/. Moreover, PALP is included standard with Sage.
PALP is described in the paper :arxiv:`math.SC/0204356`. Its distribution also contains the application nef.x, which was created by Erwin Riegler and computes nef-partitions and Hodge data for toric complete intersections.
ACKNOWLEDGMENT: polytope.py module written by William Stein was used as an example of organizing an interface between an external program and Sage. William Stein also helped Andrey Novoseltsev with debugging and tuning of this module.
Robert Bradshaw helped Andrey Novoseltsev to realize plot3d function.
.. note::
IMPORTANT: PALP requires some parameters to be determined during compilation time, i.e., the maximum dimension of polytopes, the maximum number of points, etc. These limitations may lead to errors during calls to different functions of these module. Currently, a ValueError exception will be raised if the output of poly.x or nef.x is empty or contains the exclamation mark. The error message will contain the exact command that caused an error, the description and vertices of the polytope, and the obtained output.
Data obtained from PALP and some other data is cached and most returned values are immutable. In particular, you cannot change the vertices of the polytope or their order after creation of the polytope.
If you are going to work with large sets of data, take a look at ``all_*`` functions in this module. They precompute different data for sequences of polynomials with a few runs of external programs. This can significantly affect the time of future computations. You can also use dump/load, but not all data will be stored (currently only faces and the number of their internal and boundary points are stored, in addition to polytope vertices and its polar).
AUTHORS:
- Andrey Novoseltsev (2007-01-11): initial version
- Andrey Novoseltsev (2007-01-15): ``all_*`` functions
- Andrey Novoseltsev (2008-04-01): second version, including:
- dual nef-partitions and necessary convex_hull and minkowski_sum
- built-in sequences of 2- and 3-dimensional reflexive polytopes
- plot3d, skeleton_show
- Andrey Novoseltsev (2009-08-26): dropped maximal dimension requirement
- Andrey Novoseltsev (2010-12-15): new version of nef-partitions
- Andrey Novoseltsev (2013-09-30): switch to PointCollection.
- Maximilian Kreuzer and Harald Skarke: authors of PALP (which was also used to obtain the list of 3-dimensional reflexive polytopes)
- Erwin Riegler: the author of nef.x
"""
#***************************************************************************** # Copyright (C) 2007-2013 Andrey Novoseltsev <novoselt@gmail.com> # Copyright (C) 2007-2013 William Stein <wstein@gmail.com> # # Distributed under the terms of the GNU General Public License (GPL) # as published by the Free Software Foundation; either version 2 of # the License, or (at your option) any later version. # http://www.gnu.org/licenses/ #***************************************************************************** from __future__ import print_function, absolute_import from six.moves import range
from sage.arith.all import gcd from sage.combinat.posets.posets import FinitePoset from sage.env import POLYTOPE_DATA_DIR from sage.geometry.cone import _ambient_space_point, integral_length from sage.geometry.hasse_diagram import Hasse_diagram_from_incidences from sage.geometry.point_collection import PointCollection,\ is_PointCollection, read_palp_point_collection from sage.geometry.toric_lattice import ToricLattice, is_ToricLattice from sage.graphs.graph import DiGraph, Graph from sage.groups.perm_gps.permgroup_element import PermutationGroupElement from sage.libs.ppl import (C_Polyhedron, Generator_System, Linear_Expression, point as PPL_point) from sage.matrix.constructor import matrix from sage.structure.element import is_Matrix from sage.misc.all import cached_method, flatten, tmp_filename from sage.misc.superseded import deprecated_function_alias from sage.modules.all import vector from sage.numerical.mip import MixedIntegerLinearProgram from sage.plot.plot3d.index_face_set import IndexFaceSet from sage.plot.plot3d.all import line3d, point3d from sage.plot.plot3d.shapes2 import text3d from sage.rings.all import Integer, ZZ, QQ from sage.sets.set import Set_generic from sage.structure.all import Sequence from sage.structure.sage_object import SageObject
from copy import copy import collections from six.moves import copyreg import os import subprocess import warnings from six import StringIO from functools import reduce
class SetOfAllLatticePolytopesClass(Set_generic): def _repr_(self): r""" Return a string representation.
TESTS::
sage: lattice_polytope.SetOfAllLatticePolytopesClass()._repr_() 'Set of all Lattice Polytopes' """
def __call__(self, x): r""" TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: lattice_polytope.SetOfAllLatticePolytopesClass().__call__(o) 3-d reflexive polytope in 3-d lattice M """ raise TypeError
SetOfAllLatticePolytopes = SetOfAllLatticePolytopesClass()
def LatticePolytope(data, compute_vertices=True, n=0, lattice=None): r""" Construct a lattice polytope.
INPUT:
- ``data`` -- points spanning the lattice polytope, specified as one of:
* a :class:`point collection <sage.geometry.point_collection.PointCollection>` (this is the preferred input and it is the quickest and the most memory efficient one);
* an iterable of iterables (for example, a list of vectors) defining the point coordinates;
* a file with matrix data, opened for reading, or
* a filename of such a file, see :func:`~sage.geometry.point_collection.read_palp_point_collection` for the file format;
- ``compute_vertices`` -- boolean (default: ``True``). If ``True``, the convex hull of the given points will be computed for determining vertices. Otherwise, the given points must be vertices;
- ``n`` -- an integer (default: 0) if ``data`` is a name of a file, that contains data blocks for several polytopes, the ``n``-th block will be used;
- ``lattice`` -- the ambient lattice of the polytope. If not given, a suitable lattice will be determined automatically, most likely the :class:`toric lattice <sage.geometry.toric_lattice.ToricLatticeFactory>` `M` of the appropriate dimension.
OUTPUT:
- a :class:`lattice polytope <LatticePolytopeClass>`.
EXAMPLES::
sage: points = [(1,0,0), (0,1,0), (0,0,1), (-1,0,0), (0,-1,0), (0,0,-1)] sage: p = LatticePolytope(points) sage: p 3-d reflexive polytope in 3-d lattice M sage: p.vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M
We draw a pretty picture of the polytope in 3-dimensional space::
sage: p.plot3d().show()
Now we add an extra point, which is in the interior of the polytope...
::
sage: points.append((0,0,0)) sage: p = LatticePolytope(points) sage: p.nvertices() 6
You can suppress vertex computation for speed but this can lead to mistakes::
sage: p = LatticePolytope(points, compute_vertices=False) ... sage: p.nvertices() 7
Given points must be in the lattice::
sage: LatticePolytope([[1/2], [3/2]]) Traceback (most recent call last): ... ValueError: points [[1/2], [3/2]] are not in 1-d lattice M!
But it is OK to create polytopes of non-maximal dimension::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (0,0,0), ....: (-1,0,0), (0,-1,0), (0,0,0), (0,0,0)]) sage: p 2-d lattice polytope in 3-d lattice M sage: p.vertices() M(-1, 0, 0), M( 0, -1, 0), M( 1, 0, 0), M( 0, 1, 0) in 3-d lattice M
An empty lattice polytope can be considered as well::
sage: p = LatticePolytope([], lattice=ToricLattice(3).dual()); p -1-d lattice polytope in 3-d lattice M sage: p.lattice_dim() 3 sage: p.npoints() 0 sage: p.nfacets() 0 sage: p.points() Empty collection in 3-d lattice M sage: p.faces() ((-1-d lattice polytope in 3-d lattice M,),) """ (lattice is None or lattice is data.module())): f = open(data) skip_palp_matrix(f, n) data = read_palp_point_collection(data) f.close() data = read_palp_point_collection(data) try: data = list(data) except TypeError: raise TypeError("cannot construct a polytope from\n%s" % data) raise ValueError("lattice must be given explicitly for " "empty polytopes!") except TypeError: raise TypeError("cannot construct a polytope from\n%s" % data)
copyreg.constructor(LatticePolytope) # "safe for unpickling"
def ReflexivePolytope(dim, n): r""" Return the `n`-th 2- or 3-dimensional reflexive polytope.
.. NOTE::
#. Numeration starts with zero: `0 \leq n \leq 15` for `{\rm dim} = 2` and `0 \leq n \leq 4318` for `{\rm dim} = 3`.
#. During the first call, all reflexive polytopes of requested dimension are loaded and cached for future use, so the first call for 3-dimensional polytopes can take several seconds, but all consecutive calls are fast.
#. Equivalent to ``ReflexivePolytopes(dim)[n]`` but checks bounds first.
EXAMPLES:
The 3rd 2-dimensional polytope is "the diamond"::
sage: ReflexivePolytope(2, 3) 2-d reflexive polytope #3 in 2-d lattice M sage: lattice_polytope.ReflexivePolytope(2,3).vertices() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M
There are 16 reflexive polygons and numeration starts with 0::
sage: ReflexivePolytope(2,16) Traceback (most recent call last): ... ValueError: there are only 16 reflexive polygons!
It is not possible to load a 4-dimensional polytope in this way::
sage: ReflexivePolytope(4,16) Traceback (most recent call last): ... NotImplementedError: only 2- and 3-dimensional reflexive polytopes are available! """ raise ValueError("there are only 4319 reflexive 3-polytopes!") else:
# Sequences of reflexive polytopes _rp = [None]*4
def ReflexivePolytopes(dim): r""" Return the sequence of all 2- or 3-dimensional reflexive polytopes.
.. NOTE::
During the first call the database is loaded and cached for future use, so repetitive calls will return the same object in memory.
:param dim: dimension of required reflexive polytopes :type dim: 2 or 3 :rtype: list of lattice polytopes
EXAMPLES:
There are 16 reflexive polygons::
sage: len(ReflexivePolytopes(2)) 16
It is not possible to load 4-dimensional polytopes in this way::
sage: ReflexivePolytopes(4) Traceback (most recent call last): ... NotImplementedError: only 2- and 3-dimensional reflexive polytopes are available! """ global _rp os.path.join(POLYTOPE_DATA_DIR, "reflexive_polytopes_%dd" % dim)) # Data files have normal form of reflexive polytopes # Prevents dimension computation later
def is_LatticePolytope(x): r""" Check if ``x`` is a lattice polytope.
INPUT:
- ``x`` -- anything.
OUTPUT:
- ``True`` if ``x`` is a :class:`lattice polytope <LatticePolytopeClass>`, ``False`` otherwise.
EXAMPLES::
sage: from sage.geometry.lattice_polytope import is_LatticePolytope sage: is_LatticePolytope(1) False sage: p = LatticePolytope([(1,0), (0,1), (-1,-1)]) sage: p 2-d reflexive polytope #0 in 2-d lattice M sage: is_LatticePolytope(p) True """
class LatticePolytopeClass(SageObject, collections.Hashable): r""" Create a lattice polytope.
.. WARNING::
This class does not perform any checks of correctness of input nor does it convert input into the standard representation. Use :func:`LatticePolytope` to construct lattice polytopes.
Lattice polytopes are immutable, but they cache most of the returned values.
INPUT:
The input can be either:
- ``points`` -- :class:`~sage.geometry.point_collection.PointCollection`;
- ``compute_vertices`` -- boolean.
or (these parameters must be given as keywords):
- ``ambient`` -- ambient structure, this polytope *must be a face of* ``ambient``;
- ``ambient_vertex_indices`` -- increasing list or tuple of integers, indices of vertices of ``ambient`` generating this polytope;
- ``ambient_facet_indices`` -- increasing list or tuple of integers, indices of facets of ``ambient`` generating this polytope.
OUTPUT:
- lattice polytope.
.. NOTE::
Every polytope has an ambient structure. If it was not specified, it is this polytope itself. """
def __init__(self, points=None, compute_vertices=None, ambient=None, ambient_vertex_indices=None, ambient_facet_indices=None): r""" Construct a lattice polytope.
See :func:`LatticePolytope` for documentation.
TESTS::
sage: LatticePolytope([(1,2,3), (4,5,6)]) # indirect test 1-d lattice polytope in 3-d lattice M """ [PPL_point(Linear_Expression(p, 0)) for p in points])) else:
def __contains__(self, point): r""" Check if ``point`` is contained in ``self``.
See :meth:`_contains` (which is called by this function) for documentation.
TESTS::
sage: p = lattice_polytope.cross_polytope(2) sage: (1,0) in p True sage: [1,0] in p True sage: (-2,0) in p False """
def __eq__(self, other): r""" Compare ``self`` with ``other``.
INPUT:
- ``other`` -- anything.
OUTPUT:
- ``True`` if ``other`` is a :class:`lattice polytope <LatticePolytopeClass>` equal to ``self``, ``False`` otherwise.
.. NOTE::
Two lattice polytopes are equal if they have the same vertices listed in the same order.
TESTS::
sage: p1 = LatticePolytope([(1,0), (0,1), (-1,-1)]) sage: p2 = LatticePolytope([(1,0), (0,1), (-1,-1)]) sage: p3 = LatticePolytope([(0,1), (1,0), (-1,-1)]) sage: p1 == p1 True sage: p1 == p2 True sage: p1 is p2 False sage: p1 == p3 False sage: p1 == 0 False """ and self._vertices == other._vertices)
@cached_method def __hash__(self): r""" Return the hash of ``self``.
OUTPUT:
- an integer.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: hash(o) == hash(o) True """ # FIXME: take into account other things that may be preset?..
def __ne__(self, other): r""" Compare ``self`` with ``other``.
INPUT:
- ``other`` -- anything.
OUTPUT:
- ``False`` if ``other`` is a :class:`lattice polytope <LatticePolytopeClass>` equal to ``self``, ``True`` otherwise.
.. NOTE::
Two lattice polytopes are if they have the same vertices listed in the same order.
TESTS::
sage: p1 = LatticePolytope([(1,0), (0,1), (-1,-1)]) sage: p2 = LatticePolytope([(1,0), (0,1), (-1,-1)]) sage: p3 = LatticePolytope([(0,1), (1,0), (-1,-1)]) sage: p1 != p1 False sage: p1 != p2 False sage: p1 is p2 False sage: p1 != p3 True sage: p1 != 0 True """
def __reduce__(self): r""" Reduction function. Does not store data that can be relatively fast recomputed.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: o.vertices() == loads(o.dumps()).vertices() True """
def __setstate__(self, state): r""" Restores the state of pickled polytope.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: o.vertices() == loads(o.dumps()).vertices() True """
def _compute_embedding(self): r""" Compute embedding data for this polytope.
Useful only if the dimension of this polytope is not equal to its ambient dimension.
TESTS::
sage: p = LatticePolytope(([1], [2], [3])) sage: hasattr(p, "_sublattice_polytope") False sage: p._compute_embedding() sage: p._sublattice_polytope 1-d lattice polytope in 1-d lattice M """ return for point in points]) # In order to use facet normals obtained from subpolytopes, we # need the following (see Trac #9188). # Basis for the ambient space with spanned subspace in front # Let's represent it as columns of a matrix # Absolute value helps to keep normals "inner"
def _compute_facets(self): r""" Compute and cache equations of facets of ``self``.
TESTS::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p._compute_facets() sage: p._facet_normals N( 1, 1, 0), N( 1, -1, 0), N(-1, -1, 0), N(-1, 1, 0) in 3-d lattice N """ # Sort normals if facets are vertices and normals[0] * self.vertex(0) + constants[0] != 0): # vector(ZZ, constants) is slow self._facet_normals, compute_vertices=False)
def _compute_hodge_numbers(self): r""" Compute Hodge numbers for the current nef_partitions.
This function (currently) always raises an exception directing to use another way for computing Hodge numbers.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: o._compute_hodge_numbers() Traceback (most recent call last): ... NotImplementedError: use nef_partitions(hodge_numbers=True)! """
def _contains(self, point, region='whole polytope'): r""" Check if ``point`` is contained in ``self``.
This function is called by :meth:`__contains__` and :meth:`contains` to ensure the same call depth for warning messages.
INPUT:
- ``point`` -- an attempt will be made to convert it into a single element of the ambient space of ``self``; if it fails, ``False`` is returned
- ``region`` -- string; can be either ``'whole polytope'`` (default), ``'interior'``, or ``'relative interior'``
OUTPUT:
- ``True`` if ``point`` is contained in the specified ``region`` of ``self``, ``False`` otherwise
TESTS::
sage: p = lattice_polytope.cross_polytope(2) sage: p._contains((1,0)) True """ except TypeError as ex: if str(ex).endswith("have incompatible lattices"): warnings.warn("you have checked if a cone contains a point " "from an incompatible lattice, this is False", stacklevel=3) return False
raise ValueError("%s is an unknown region" % region) return False return False
def _embed(self, data): r""" Embed given point(s) into the ambient space of this polytope.
INPUT:
- ``data`` - point or matrix of points (as columns) in the affine subspace spanned by this polytope
OUTPUT: The same point(s) in the coordinates of the ambient space of this polytope.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: o._embed(o.vertices()) == o.vertices() True sage: m = matrix(ZZ, 3) sage: m[0, 0] = 1 sage: m[1, 1] = 1 sage: p = o.affine_transform(m) sage: p._embed((0,0)) M(-1, 0, 0) """ r = [M(self._embedding_matrix * point + self._shift_vector) for point in data] for point in r: point.set_immutable() return PointCollection(r, M) else: self._shift_vector)
def _latex_(self): r""" Return the latex representation of self.
OUTPUT:
- string
EXAMPLES:
Arbitrary lattice polytopes are printed as `\Delta^d`, where `d` is the (actual) dimension of the polytope::
sage: LatticePolytope([(1,1), (0,0)])._latex_() '\\Delta^{1}'
For 2- and 3-d reflexive polytopes the index in the internal database appears as a subscript::
sage: print(ReflexivePolytope(2, 3)._latex_()) \Delta^{2}_{3} """
def _palp(self, command, reduce_dimension=False): r""" Run ``command`` on vertices of this polytope.
Returns the output of ``command`` as a string.
.. note::
PALP cannot be called for polytopes that do not span the ambient space. If you specify ``reduce_dimension=True`` argument, PALP will be called for vertices of this polytope in some basis of the affine space it spans.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: o._palp("poly.x -f") 'M:7 6 N:27 8 Pic:17 Cor:0\n' sage: print(o._palp("nef.x -f -N -p")) # random time information M:27 8 N:7 6 codim=2 #part=5 H:[0] P:0 V:2 4 5 0sec 0cpu H:[0] P:2 V:3 4 5 0sec 0cpu H:[0] P:3 V:4 5 0sec 0cpu np=3 d:1 p:1 0sec 0cpu
sage: p = LatticePolytope([[1]]) sage: p._palp("poly.x -f") Traceback (most recent call last): ... ValueError: Cannot run "poly.x -f" for the zero-dimensional polytope! Polytope: 0-d lattice polytope in 1-d lattice M
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p._palp("poly.x -f") Traceback (most recent call last): ... ValueError: Cannot run PALP for a 2-dimensional polytope in a 3-dimensional space! sage: p._palp("poly.x -f", reduce_dimension=True) 'M:5 4 F:4\n' """ + "polytope!\nPolytope: %s") % (command, self)) "in a %d-dimensional space!") % (self.dim(), self.lattice_dim())) "!" in result or "failed." in result or "increase" in result or "Unable" in result): "Output:", result]
@cached_method def _PPL(self): r""" Return the Parma Polyhedra Library (PPL) representation of ``self``.
OUTPUT:
- :class:`~sage.libs.ppl.C_Polyhedron`
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o._PPL() A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 6 points sage: o._PPL().minimized_generators() Generator_System {point(-1/1, 0/1, 0/1), point(0/1, -1/1, 0/1), point(0/1, 0/1, -1/1), point(0/1, 0/1, 1/1), point(0/1, 1/1, 0/1), point(1/1, 0/1, 0/1)} sage: o._PPL().minimized_constraints() Constraint_System {x0-x1-x2+1>=0, x0+x1-x2+1>=0, x0+x1+x2+1>=0, x0-x1+x2+1>=0, -x0-x1+x2+1>=0, -x0-x1-x2+1>=0, -x0+x1-x2+1>=0, -x0+x1+x2+1>=0} """ [PPL_point(Linear_Expression(v, 0)) for v in self.vertices()]))
def _pullback(self, data): r""" Pull back given point(s) to the affine subspace spanned by this polytope.
INPUT:
- ``data`` -- rational point or matrix of points (as columns) in the ambient space
OUTPUT: The same point(s) in the coordinates of the affine subspace space spanned by this polytope.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: o._pullback(o.vertices().column_matrix()) == o.vertices().column_matrix() True sage: m = matrix(ZZ, 3) sage: m[0, 0] = 1 sage: m[1, 1] = 1 sage: p = o.affine_transform(m) sage: p._pullback((0, 0, 0)) [1, 0] """ r = [self._pullback(point) for point in data] for point in r: point.set_immutable() return PointCollection(r, self._sublattice) r = matrix([self._pullback(col) for col in data.columns(copy=False)]).transpose() return r
def _read_equations(self, data): r""" Read equations of facets/vertices of polar polytope from string or file.
TESTS:
For a reflexive polytope construct the polar polytope::
sage: p = LatticePolytope([(1,0), (0,1), (-1,-1)]) sage: p.vertices() M( 1, 0), M( 0, 1), M(-1, -1) in 2-d lattice M sage: s = p.poly_x("e") sage: print(s) 3 2 Vertices of P-dual <-> Equations of P 2 -1 -1 2 -1 -1 sage: "_polar" in p.__dict__ False sage: p._read_equations(s) sage: p._polar._vertices N( 2, -1), N(-1, 2), N(-1, -1) in 2-d lattice N
For a non-reflexive polytope cache facet equations::
sage: p = LatticePolytope([(1,0), (0,2), (-1,-3 )]) sage: p.vertices() M( 1, 0), M( 0, 2), M(-1, -3) in 2-d lattice M sage: "_facet_normals" in p.__dict__ False sage: "_facet_constants" in p.__dict__ False sage: s = p.poly_x("e") sage: print(s) 3 2 Equations of P 5 -1 2 -2 -1 2 -3 2 3 sage: p._read_equations(s) sage: p._facet_normals N( 5, -1), N(-2, -1), N(-3, 2) in 2-d lattice N sage: p._facet_constants (2, 2, 3) """ # If it is already known that this polytope is reflexive, its # polar (whose vertices are equations of facets of this one) # is already computed and there is no need to read equations # of facets of this polytope. Moreover, doing so can corrupt # data if this polytope was constructed as polar. Skip input. skip_palp_matrix(data) return read_palp_point_collection(data, N), compute_vertices=False) else:
def _read_nef_partitions(self, data): r""" Read nef-partitions of ``self`` from ``data``.
INPUT:
- ``data`` -- a string or a file.
OUTPUT:
- none.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: s = o.nef_x("-p -N -Lv") sage: print(s) # random time values M:27 8 N:7 6 codim=2 #part=5 3 6 Vertices in N-lattice: 1 0 0 -1 0 0 0 1 0 0 -1 0 0 0 1 0 0 -1 ------------------------------ 1 0 0 1 0 0 d=2 codim=2 0 1 0 0 1 0 d=2 codim=2 0 0 1 0 0 1 d=2 codim=2 P:0 V:2 4 5 (0 2) (1 1) (2 0) 0sec 0cpu P:2 V:3 4 5 (1 1) (1 1) (1 1) 0sec 0cpu P:3 V:4 5 (0 2) (1 1) (1 1) 0sec 0cpu np=3 d:1 p:1 0sec 0cpu
We make a copy of the octahedron since direct use of this function may destroy cache integrity and lead so strange effects in other doctests::
sage: o_copy = LatticePolytope(o.vertices()) sage: "_nef_partitions" in o_copy.__dict__ False sage: o_copy._read_nef_partitions(s) sage: o_copy._nef_partitions [ Nef-partition {0, 1, 3} U {2, 4, 5}, Nef-partition {0, 1, 2} U {3, 4, 5}, Nef-partition {0, 1, 2, 3} U {4, 5} ] """ raise RuntimeError("nef.x changed the order of vertices!") raise ValueError("more data expected!") # Set the stuff for h in line[start:end].split()) raise ValueError("""Wrong data format, cannot find "np="!""")
def _repr_(self): r""" Return a string representation of ``self``.
OUTPUT:
- a string.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: o._repr_() '3-d reflexive polytope in 3-d lattice M' """ except ValueError: pass else: parts.append("%d-d lattice" % self.lattice_dim()) else:
def _sort_faces(self, faces): r""" Return sorted (if necessary) ``faces`` as a tuple.
This function ensures that zero-dimensional faces are listed in agreement with the order of corresponding vertices and facets with facet normals.
INPUT:
- ``faces`` -- iterable of :class:`lattice polytopes <LatticePolytopeClass>`.
OUTPUT:
- :class:`tuple` of :class:`lattice polytopes <LatticePolytopeClass>`.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: # indirect doctest sage: for i, face in enumerate(o.faces(0)): ....: if face.vertex(0) != o.vertex(i): ....: print("Wrong order!") """ key=lambda f: f._ambient_vertex_indices)) hasattr(self, "_facet_normals"): # If we already have facet normals, sort according to them
@cached_method def adjacent(self): r""" Return faces adjacent to ``self`` in the ambient face lattice.
Two *distinct* faces `F_1` and `F_2` of the same face lattice are **adjacent** if all of the following conditions hold:
* `F_1` and `F_2` have the same dimension `d`;
* `F_1` and `F_2` share a facet of dimension `d-1`;
* `F_1` and `F_2` are facets of some face of dimension `d+1`, unless `d` is the dimension of the ambient structure.
OUTPUT:
- :class:`tuple` of :class:`lattice polytopes <LatticePolytopeClass>`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.adjacent() () sage: face = o.faces(1)[0] sage: face.adjacent() (1-d face of 3-d reflexive polytope in 3-d lattice M, 1-d face of 3-d reflexive polytope in 3-d lattice M, 1-d face of 3-d reflexive polytope in 3-d lattice M, 1-d face of 3-d reflexive polytope in 3-d lattice M) """
def affine_transform(self, a=1, b=0): r""" Return a*P+b, where P is this lattice polytope.
.. note::
#. While ``a`` and ``b`` may be rational, the final result must be a lattice polytope, i.e. all vertices must be integral.
#. If the transform (restricted to this polytope) is bijective, facial structure will be preserved, e.g. the first facet of the image will be spanned by the images of vertices which span the first facet of the original polytope.
INPUT:
- ``a`` - (default: 1) rational scalar or matrix
- ``b`` - (default: 0) rational scalar or vector, scalars are interpreted as vectors with the same components
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(2) sage: o.vertices() M( 1, 0), M( 0, 1), M(-1, 0), M( 0, -1) in 2-d lattice M sage: o.affine_transform(2).vertices() M( 2, 0), M( 0, 2), M(-2, 0), M( 0, -2) in 2-d lattice M sage: o.affine_transform(1,1).vertices() M(2, 1), M(1, 2), M(0, 1), M(1, 0) in 2-d lattice M sage: o.affine_transform(b=1).vertices() M(2, 1), M(1, 2), M(0, 1), M(1, 0) in 2-d lattice M sage: o.affine_transform(b=(1, 0)).vertices() M(2, 0), M(1, 1), M(0, 0), M(1, -1) in 2-d lattice M sage: a = matrix(QQ, 2, [1/2, 0, 0, 3/2]) sage: o.polar().vertices() N( 1, 1), N( 1, -1), N(-1, -1), N(-1, 1) in 2-d lattice N sage: o.polar().affine_transform(a, (1/2, -1/2)).vertices() M(1, 1), M(1, -2), M(0, -2), M(0, 1) in 2-d lattice M
While you can use rational transformation, the result must be integer::
sage: o.affine_transform(a) Traceback (most recent call last): ... ValueError: points [(1/2, 0), (0, 3/2), (-1/2, 0), (0, -3/2)] are not in 2-d lattice M! """ else: # Prevent long chains of "original-transform" r._original = self._original else:
def ambient(self): r""" Return the ambient structure of ``self``.
OUTPUT:
- lattice polytope containing ``self`` as a face.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.ambient() 3-d reflexive polytope in 3-d lattice M sage: o.ambient() is o True sage: face = o.faces(1)[0] sage: face 1-d face of 3-d reflexive polytope in 3-d lattice M sage: face.ambient() 3-d reflexive polytope in 3-d lattice M sage: face.ambient() is o True """
def ambient_facet_indices(self): r""" Return indices of facets of the ambient polytope containing ``self``.
OUTPUT:
- increasing :class:`tuple` of integers.
EXAMPLES:
The polytope itself is not contained in any of its facets::
sage: o = lattice_polytope.cross_polytope(3) sage: o.ambient_facet_indices() ()
But each of its other faces is contained in one or more facets::
sage: face = o.faces(1)[0] sage: face.ambient_facet_indices() (4, 5) sage: face.vertices() M(1, 0, 0), M(0, 1, 0) in 3-d lattice M sage: o.facets()[face.ambient_facet_indices()[0]].vertices() M(1, 0, 0), M(0, 1, 0), M(0, 0, -1) in 3-d lattice M """
@cached_method def ambient_point_indices(self): r""" Return indices of points of the ambient polytope contained in this one.
OUTPUT:
- :class:`tuple` of integers, the order corresponds to the order of points of this polytope.
EXAMPLES::
sage: cube = lattice_polytope.cross_polytope(3).polar() sage: face = cube.facets()[0] sage: face.ambient_point_indices() (4, 5, 6, 7, 8, 9, 10, 11, 12) sage: cube.points(face.ambient_point_indices()) == face.points() True """ return tuple(range(self.npoints()))
@cached_method def ambient_ordered_point_indices(self): r""" Return indices of points of the ambient polytope contained in this one.
OUTPUT:
- :class:`tuple` of integers such that ambient points in this order are geometrically ordered, e.g. for an edge points will appear from one end point to the other.
EXAMPLES::
sage: cube = lattice_polytope.cross_polytope(3).polar() sage: face = cube.facets()[0] sage: face.ambient_ordered_point_indices() (5, 8, 4, 9, 10, 11, 6, 12, 7) sage: cube.points(face.ambient_ordered_point_indices()) N(-1, -1, -1), N(-1, -1, 0), N(-1, -1, 1), N(-1, 0, -1), N(-1, 0, 0), N(-1, 0, 1), N(-1, 1, -1), N(-1, 1, 0), N(-1, 1, 1) in 3-d lattice N """ return tuple(range(self.npoints()))
def ambient_vertex_indices(self): r""" Return indices of vertices of the ambient structure generating ``self``.
OUTPUT:
- increasing :class:`tuple` of integers.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.ambient_vertex_indices() (0, 1, 2, 3, 4, 5) sage: face = o.faces(1)[0] sage: face.ambient_vertex_indices() (0, 1) """
@cached_method def boundary_point_indices(self): r""" Return indices of (relative) boundary lattice points of this polytope.
OUTPUT:
- increasing :class:`tuple` of integers.
EXAMPLES:
All points but the origin are on the boundary of this square::
sage: square = lattice_polytope.cross_polytope(2).polar() sage: square.points() N( 1, 1), N( 1, -1), N(-1, -1), N(-1, 1), N(-1, 0), N( 0, -1), N( 0, 0), N( 0, 1), N( 1, 0) in 2-d lattice N sage: square.boundary_point_indices() (0, 1, 2, 3, 4, 5, 7, 8)
For an edge the boundary is formed by the end points::
sage: face = square.edges()[0] sage: face.points() N(-1, -1), N(-1, 1), N(-1, 0) in 2-d lattice N sage: face.boundary_point_indices() (0, 1) """ for i, c in enumerate(self.distances().columns(copy=False)) if len(c.nonzero_positions()) < self.nfacets())
def boundary_points(self): r""" Return (relative) boundary lattice points of this polytope.
OUTPUT:
- a :class:`point collection <PointCollection>`.
EXAMPLES:
All points but the origin are on the boundary of this square::
sage: square = lattice_polytope.cross_polytope(2).polar() sage: square.boundary_points() N( 1, 1), N( 1, -1), N(-1, -1), N(-1, 1), N(-1, 0), N( 0, -1), N( 0, 1), N( 1, 0) in 2-d lattice N
For an edge the boundary is formed by the end points::
sage: face = square.edges()[0] sage: face.boundary_points() N(-1, -1), N(-1, 1) in 2-d lattice N """
def contains(self, *args): r""" Check if a given point is contained in ``self``.
INPUT:
- an attempt will be made to convert all arguments into a single element of the ambient space of ``self``; if it fails, ``False`` will be returned
OUTPUT:
- ``True`` if the given point is contained in ``self``, ``False`` otherwise
EXAMPLES::
sage: p = lattice_polytope.cross_polytope(2) sage: p.contains(p.lattice()(1,0)) True sage: p.contains((1,0)) True sage: p.contains(1,0) True sage: p.contains((2,0)) False """
@cached_method def dim(self): r""" Return the dimension of this polytope.
EXAMPLES:
We create a 3-dimensional octahedron and check its dimension::
sage: o = lattice_polytope.cross_polytope(3) sage: o.dim() 3
Now we create a 2-dimensional diamond in a 3-dimensional space::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p.dim() 2 sage: p.lattice_dim() 3 """
def distances(self, point=None): r""" Return the matrix of distances for this polytope or distances for the given point.
The matrix of distances m gives distances m[i,j] between the i-th facet (which is also the i-th vertex of the polar polytope in the reflexive case) and j-th point of this polytope.
If point is specified, integral distances from the point to all facets of this polytope will be computed.
EXAMPLES: The matrix of distances for a 3-dimensional octahedron::
sage: o = lattice_polytope.cross_polytope(3) sage: o.distances() [2 0 0 0 2 2 1] [2 2 0 0 0 2 1] [2 2 2 0 0 0 1] [2 0 2 0 2 0 1] [0 0 2 2 2 0 1] [0 0 0 2 2 2 1] [0 2 0 2 0 2 1] [0 2 2 2 0 0 1]
Distances from facets to the point (1,2,3)::
sage: o.distances([1,2,3]) (-3, 1, 7, 3, 1, -5, -1, 5)
It is OK to use RATIONAL coordinates::
sage: o.distances([1,2,3/2]) (-3/2, 5/2, 11/2, 3/2, -1/2, -7/2, 1/2, 7/2) sage: o.distances([1,2,sqrt(2)]) Traceback (most recent call last): ... TypeError: unable to convert sqrt(2) to an element of Rational Field
Now we create a non-spanning polytope::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p.distances() [2 2 0 0 1] [2 0 0 2 1] [0 0 2 2 1] [0 2 2 0 1] sage: p.distances((1/2, 3, 0)) (9/2, -3/2, -5/2, 7/2)
This point is not even in the affine subspace of the polytope::
sage: p.distances((1, 1, 1)) (3, 1, -1, 1) """ self.facet_constants()) for F, c in zip(self.facet_normals(), self.facet_constants())])
@cached_method def dual(self): r""" Return the dual face under face duality of polar reflexive polytopes.
This duality extends the correspondence between vertices and facets.
OUTPUT:
- a :class:`lattice polytope <LatticePolytopeClass>`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(4) sage: e = o.edges()[0]; e 1-d face of 4-d reflexive polytope in 4-d lattice M sage: ed = e.dual(); ed 2-d face of 4-d reflexive polytope in 4-d lattice N sage: ed.ambient() is e.ambient().polar() True sage: e.ambient_vertex_indices() == ed.ambient_facet_indices() True sage: e.ambient_facet_indices() == ed.ambient_vertex_indices() True """
@cached_method def dual_lattice(self): r""" Return the dual of the ambient lattice of ``self``.
OUTPUT:
- a lattice. If possible (that is, if :meth:`lattice` has a ``dual()`` method), the dual lattice is returned. Otherwise, `\ZZ^n` is returned, where `n` is the dimension of ``self``.
EXAMPLES::
sage: LatticePolytope([(1,0)]).dual_lattice() 2-d lattice N sage: LatticePolytope([], lattice=ZZ^3).dual_lattice() Ambient free module of rank 3 over the principal ideal domain Integer Ring """
def edges(self): r""" Return edges (faces of dimension 1) of ``self``.
OUTPUT:
- :class:`tuple` of :class:`lattice polytopes <LatticePolytopeClass>`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.edges() (1-d face of 3-d reflexive polytope in 3-d lattice M, ... 1-d face of 3-d reflexive polytope in 3-d lattice M) sage: len(o.edges()) 12 """
edges_lp = deprecated_function_alias(22122, edges)
@cached_method def face_lattice(self): r""" Return the face lattice of ``self``.
This lattice will have the empty polytope as the bottom and this polytope itself as the top.
OUTPUT:
- :class:`finite poset <sage.combinat.posets.posets.FinitePoset>` of :class:`lattice polytopes <LatticePolytopeClass>`.
EXAMPLES:
Let's take a look at the face lattice of a square::
sage: square = LatticePolytope([(0,0), (1,0), (1,1), (0,1)]) sage: L = square.face_lattice() sage: L Finite poset containing 10 elements with distinguished linear extension
To see all faces arranged by dimension, you can do this::
sage: for level in L.level_sets(): print(level) [-1-d face of 2-d lattice polytope in 2-d lattice M] [0-d face of 2-d lattice polytope in 2-d lattice M, 0-d face of 2-d lattice polytope in 2-d lattice M, 0-d face of 2-d lattice polytope in 2-d lattice M, 0-d face of 2-d lattice polytope in 2-d lattice M] [1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M] [2-d lattice polytope in 2-d lattice M]
For a particular face you can look at its actual vertices... ::
sage: face = L.level_sets()[1][0] sage: face.vertices() M(0, 0) in 2-d lattice M
... or you can see the index of the vertex of the original polytope that corresponds to the above one::
sage: face.ambient_vertex_indices() (0,) sage: square.vertex(0) M(0, 0)
An alternative to extracting faces from the face lattice is to use :meth:`faces` method::
sage: face is square.faces(dim=0)[0] True
The advantage of working with the face lattice directly is that you can (relatively easily) get faces that are related to the given one::
sage: face = L.level_sets()[1][0] sage: D = L.hasse_diagram() sage: D.neighbors(face) [1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M, -1-d face of 2-d lattice polytope in 2-d lattice M]
However, you can achieve some of this functionality using :meth:`facets`, :meth:`facet_of`, and :meth:`adjacent` methods::
sage: face = square.faces(0)[0] sage: face 0-d face of 2-d lattice polytope in 2-d lattice M sage: face.vertices() M(0, 0) in 2-d lattice M sage: face.facets() (-1-d face of 2-d lattice polytope in 2-d lattice M,) sage: face.facet_of() (1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M) sage: face.adjacent() (0-d face of 2-d lattice polytope in 2-d lattice M, 0-d face of 2-d lattice polytope in 2-d lattice M) sage: face.adjacent()[0].vertices() M(1, 0) in 2-d lattice M
Note that if ``p`` is a face of ``superp``, then the face lattice of ``p`` consists of (appropriate) faces of ``superp``::
sage: superp = LatticePolytope([(1,2,3,4), (5,6,7,8), ....: (1,2,4,8), (1,3,9,7)]) sage: superp.face_lattice() Finite poset containing 16 elements with distinguished linear extension sage: superp.face_lattice().top() 3-d lattice polytope in 4-d lattice M sage: p = superp.facets()[0] sage: p 2-d face of 3-d lattice polytope in 4-d lattice M sage: p.face_lattice() Finite poset containing 8 elements with distinguished linear extension sage: p.face_lattice().bottom() -1-d face of 3-d lattice polytope in 4-d lattice M sage: p.face_lattice().top() 2-d face of 3-d lattice polytope in 4-d lattice M sage: p.face_lattice().top() is p True """ # We need to compute face lattice on our own. if normal * vertex + self.facet_constant(j) == 0]
ambient_vertex_indices=vertices, ambient_facet_indices=facets)
vertex_to_facets, facet_to_vertices, LPFace, key = id(self)) else: # Get face lattice as a sublattice of the ambient one new_face._ambient_vertex_indices): else: # Last level is very easy to build, so we do it separately # even though the above cycle could do it too.
def faces(self, dim=None, codim=None): r""" Return faces of ``self`` of specified (co)dimension.
INPUT:
- ``dim`` -- integer, dimension of the requested faces;
- ``codim`` -- integer, codimension of the requested faces.
.. NOTE::
You can specify at most one parameter. If you don't give any, then all faces will be returned.
OUTPUT:
- if either ``dim`` or ``codim`` is given, the output will be a :class:`tuple` of :class:`lattice polytopes <LatticePolytopeClass>`;
- if neither ``dim`` nor ``codim`` is given, the output will be the :class:`tuple` of tuples as above, giving faces of all existing dimensions. If you care about inclusion relations between faces, consider using :meth:`face_lattice` or :meth:`adjacent`, :meth:`facet_of`, and :meth:`facets`.
EXAMPLES:
Let's take a look at the faces of a square::
sage: square = LatticePolytope([(0,0), (1,0), (1,1), (0,1)]) sage: square.faces() ((-1-d face of 2-d lattice polytope in 2-d lattice M,), (0-d face of 2-d lattice polytope in 2-d lattice M, 0-d face of 2-d lattice polytope in 2-d lattice M, 0-d face of 2-d lattice polytope in 2-d lattice M, 0-d face of 2-d lattice polytope in 2-d lattice M), (1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M), (2-d lattice polytope in 2-d lattice M,))
Its faces of dimension one (i.e., edges)::
sage: square.faces(dim=1) (1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M, 1-d face of 2-d lattice polytope in 2-d lattice M)
Its faces of codimension one are the same (also edges)::
sage: square.faces(codim=1) is square.faces(dim=1) True
Let's pick a particular face::
sage: face = square.faces(dim=1)[0]
Now you can look at the actual vertices of this face... ::
sage: face.vertices() M(0, 0), M(0, 1) in 2-d lattice M
... or you can see indices of the vertices of the original polytope that correspond to the above ones::
sage: face.ambient_vertex_indices() (0, 3) sage: square.vertices(face.ambient_vertex_indices()) M(0, 0), M(0, 1) in 2-d lattice M """ raise ValueError( "dimension and codimension cannot be specified together!") self.face_lattice().level_sets())) else:
faces_lp = deprecated_function_alias(22122, faces)
def facet_constant(self, i): r""" Return the constant in the ``i``-th facet inequality of this polytope.
This is equivalent to ``facet_constants()[i]``.
INPUT:
- ``i`` -- integer; the index of the facet
OUTPUT:
- integer -- the constant in the ``i``-th facet inequality.
.. SEEALSO::
:meth:`facet_constants`, :meth:`facet_normal`, :meth:`facet_normals`, :meth:`facets`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.facet_constant(0) 1 sage: o.facet_constant(0) == o.facet_constants()[0] True """
def facet_constants(self): r""" Return facet constants of ``self``.
Facet inequalities have form `n \cdot x + c \geq 0` where `n` is the inner normal and `c` is a constant.
OUTPUT:
- an integer vector
.. SEEALSO::
:meth:`facet_constant`, :meth:`facet_normal`, :meth:`facet_normals`, :meth:`facets`.
EXAMPLES:
For reflexive polytopes all constants are 1::
sage: o = lattice_polytope.cross_polytope(3) sage: o.vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M sage: o.facet_constants() (1, 1, 1, 1, 1, 1, 1, 1)
Here is an example of a 3-dimensional polytope in a 4-dimensional space with 3 facets containing the origin::
sage: p = LatticePolytope([(0,0,0,0), (1,1,1,3), ....: (1,-1,1,3), (-1,-1,1,3)]) sage: p.vertices() M( 0, 0, 0, 0), M( 1, 1, 1, 3), M( 1, -1, 1, 3), M(-1, -1, 1, 3) in 4-d lattice M sage: p.facet_constants() (0, 0, 3, 0) """
def facet_normal(self, i): r""" Return the inner normal to the ``i``-th facet of this polytope.
This is equivalent to ``facet_normals()[i]``.
INPUT:
- ``i`` -- integer; the index of the facet
OUTPUT:
- a vector
.. SEEALSO::
:meth:`facet_constant`, :meth:`facet_constants`, :meth:`facet_normals`, :meth:`facets`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.facet_normal(0) N(1, -1, -1) sage: o.facet_normal(0) is o.facet_normals()[0] True """
def facet_normals(self): r""" Return inner normals to the facets of ``self``.
If this polytope is not full-dimensional, facet normals will define this polytope in the affine subspace spanned by it.
OUTPUT:
- a :class:`point collection <PointCollection>` in the :meth:`dual_lattice` of ``self``.
.. SEEALSO::
:meth:`facet_constant`, :meth:`facet_constants`, :meth:`facet_normal`, :meth:`facets`.
EXAMPLES:
Normals to facets of an octahedron are vertices of a cube::
sage: o = lattice_polytope.cross_polytope(3) sage: o.vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M sage: o.facet_normals() N( 1, -1, -1), N( 1, 1, -1), N( 1, 1, 1), N( 1, -1, 1), N(-1, -1, 1), N(-1, -1, -1), N(-1, 1, -1), N(-1, 1, 1) in 3-d lattice N
Here is an example of a 3-dimensional polytope in a 4-dimensional space::
sage: p = LatticePolytope([(0,0,0,0), (1,1,1,3), ....: (1,-1,1,3), (-1,-1,1,3)]) sage: p.vertices() M( 0, 0, 0, 0), M( 1, 1, 1, 3), M( 1, -1, 1, 3), M(-1, -1, 1, 3) in 4-d lattice M sage: p.facet_normals() N( 0, 3, 0, 1), N( 1, -1, 0, 0), N( 0, 0, 0, -1), N(-3, 0, 0, 1) in 4-d lattice N sage: p.facet_constants() (0, 0, 3, 0)
Now we manually compute the distance matrix of this polytope. Since it is a simplex, each line (corresponding to a facet) should consist of zeros (indicating generating vertices of the corresponding facet) and a single positive number (since our normals are inner)::
sage: matrix([[n * v + c for v in p.vertices()] ....: for n, c in zip(p.facet_normals(), p.facet_constants())]) [0 6 0 0] [0 0 2 0] [3 0 0 0] [0 0 0 6] """
@cached_method def facet_of(self): r""" Return elements of the ambient face lattice having ``self`` as a facet.
OUTPUT:
- :class:`tuple` of :class:`lattice polytopes <LatticePolytopeClass>`.
EXAMPLES::
sage: square = LatticePolytope([(0,0), (1,0), (1,1), (0,1)]) sage: square.facet_of() () sage: face = square.faces(0)[0] sage: len(face.facet_of()) 2 sage: face.facet_of()[1] 1-d face of 2-d lattice polytope in 2-d lattice M """
def facets(self): r""" Return facets (faces of codimension 1) of ``self``.
OUTPUT:
- :class:`tuple` of :class:`lattice polytopes <LatticePolytopeClass>`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.facets() (2-d face of 3-d reflexive polytope in 3-d lattice M, ... 2-d face of 3-d reflexive polytope in 3-d lattice M) sage: len(o.facets()) 8 """
facets_lp = deprecated_function_alias(22122, facets)
# Dictionaries of normal forms _rp_dict = [None]*4
@cached_method def index(self): r""" Return the index of this polytope in the internal database of 2- or 3-dimensional reflexive polytopes. Databases are stored in the directory of the package.
.. note::
The first call to this function for each dimension can take a few seconds while the dictionary of all polytopes is constructed, but after that it is cached and fast.
:rtype: integer
EXAMPLES: We check what is the index of the "diamond" in the database::
sage: d = lattice_polytope.cross_polytope(2) sage: d.index() 3
Note that polytopes with the same index are not necessarily the same::
sage: d.vertices() M( 1, 0), M( 0, 1), M(-1, 0), M( 0, -1) in 2-d lattice M sage: lattice_polytope.ReflexivePolytope(2,3).vertices() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M
But they are in the same `GL(Z^n)` orbit and have the same normal form::
sage: d.normal_form() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M sage: lattice_polytope.ReflexivePolytope(2,3).normal_form() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M """ raise NotImplementedError("only reflexive polytopes can be indexed!") raise NotImplementedError("only 2- and 3-dimensional polytopes can be indexed!")
@cached_method def interior_point_indices(self): r""" Return indices of (relative) interior lattice points of this polytope.
OUTPUT:
- increasing :class:`tuple` of integers.
EXAMPLES:
The origin is the only interior point of this square::
sage: square = lattice_polytope.cross_polytope(2).polar() sage: square.points() N( 1, 1), N( 1, -1), N(-1, -1), N(-1, 1), N(-1, 0), N( 0, -1), N( 0, 0), N( 0, 1), N( 1, 0) in 2-d lattice N sage: square.interior_point_indices() (6,)
Its edges also have a single interior point each::
sage: face = square.edges()[0] sage: face.points() N(-1, -1), N(-1, 1), N(-1, 0) in 2-d lattice N sage: face.interior_point_indices() (2,) """ for i, c in enumerate(self.distances().columns(copy=False)) if len(c.nonzero_positions()) == self.nfacets())
def interior_points(self): r""" Return (relative) boundary lattice points of this polytope.
OUTPUT:
- a :class:`point collection <PointCollection>`.
EXAMPLES:
The origin is the only interior point of this square::
sage: square = lattice_polytope.cross_polytope(2).polar() sage: square.interior_points() N(0, 0) in 2-d lattice N
Its edges also have a single interior point each::
sage: face = square.edges()[0] sage: face.interior_points() N(-1, 0) in 2-d lattice N """
@cached_method def is_reflexive(self): r""" Return True if this polytope is reflexive.
EXAMPLES: The 3-dimensional octahedron is reflexive (and 4319 other 3-polytopes)::
sage: o = lattice_polytope.cross_polytope(3) sage: o.is_reflexive() True
But not all polytopes are reflexive::
sage: p = LatticePolytope([(1,0,0), (0,1,17), (-1,0,0), (0,-1,0)]) sage: p.is_reflexive() False
Only full-dimensional polytopes can be reflexive (otherwise the polar set is not a polytope at all, since it is unbounded)::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p.is_reflexive() False """ all(c == 1 for c in self.facet_constants())
def lattice(self): r""" Return the ambient lattice of ``self``.
OUTPUT:
- a lattice.
EXAMPLES::
sage: lattice_polytope.cross_polytope(3).lattice() 3-d lattice M """
def lattice_dim(self): r""" Return the dimension of the ambient lattice of ``self``.
OUTPUT:
- integer.
EXAMPLES::
sage: p = LatticePolytope([(1,0)]) sage: p.lattice_dim() 2 sage: p.dim() 0 """
def linearly_independent_vertices(self): r""" Return a maximal set of linearly independent vertices.
OUTPUT:
A tuple of vertex indices.
EXAMPLES::
sage: L = LatticePolytope([[0, 0], [-1, 1], [-1, -1]]) sage: L.linearly_independent_vertices() (1, 2) sage: L = LatticePolytope([[0, 0, 0]]) sage: L.linearly_independent_vertices() () sage: L = LatticePolytope([[0, 1, 0]]) sage: L.linearly_independent_vertices() (0,) """
def nef_partitions(self, keep_symmetric=False, keep_products=True, keep_projections=True, hodge_numbers=False): r""" Return 2-part nef-partitions of ``self``.
INPUT:
- ``keep_symmetric`` -- (default: ``False``) if ``True``, "-s" option will be passed to ``nef.x`` in order to keep symmetric partitions, i.e. partitions related by lattice automorphisms preserving ``self``;
- ``keep_products`` -- (default: ``True``) if ``True``, "-D" option will be passed to ``nef.x`` in order to keep product partitions, with corresponding complete intersections being direct products;
- ``keep_projections`` -- (default: ``True``) if ``True``, "-P" option will be passed to ``nef.x`` in order to keep projection partitions, i.e. partitions with one of the parts consisting of a single vertex;
- ``hodge_numbers`` -- (default: ``False``) if ``False``, "-p" option will be passed to ``nef.x`` in order to skip Hodge numbers computation, which takes a lot of time.
OUTPUT:
- a sequence of :class:`nef-partitions <NefPartition>`.
Type ``NefPartition?`` for definitions and notation.
EXAMPLES:
Nef-partitions of the 4-dimensional cross-polytope::
sage: p = lattice_polytope.cross_polytope(4) sage: p.nef_partitions() [ Nef-partition {0, 1, 4, 5} U {2, 3, 6, 7} (direct product), Nef-partition {0, 1, 2, 4} U {3, 5, 6, 7}, Nef-partition {0, 1, 2, 4, 5} U {3, 6, 7}, Nef-partition {0, 1, 2, 4, 5, 6} U {3, 7} (direct product), Nef-partition {0, 1, 2, 3} U {4, 5, 6, 7}, Nef-partition {0, 1, 2, 3, 4} U {5, 6, 7}, Nef-partition {0, 1, 2, 3, 4, 5} U {6, 7}, Nef-partition {0, 1, 2, 3, 4, 5, 6} U {7} (projection) ]
Now we omit projections::
sage: p.nef_partitions(keep_projections=False) [ Nef-partition {0, 1, 4, 5} U {2, 3, 6, 7} (direct product), Nef-partition {0, 1, 2, 4} U {3, 5, 6, 7}, Nef-partition {0, 1, 2, 4, 5} U {3, 6, 7}, Nef-partition {0, 1, 2, 4, 5, 6} U {3, 7} (direct product), Nef-partition {0, 1, 2, 3} U {4, 5, 6, 7}, Nef-partition {0, 1, 2, 3, 4} U {5, 6, 7}, Nef-partition {0, 1, 2, 3, 4, 5} U {6, 7} ]
Currently Hodge numbers cannot be computed for a given nef-partition::
sage: p.nef_partitions()[1].hodge_numbers() Traceback (most recent call last): ... NotImplementedError: use nef_partitions(hodge_numbers=True)!
But they can be obtained from ``nef.x`` for all nef-partitions at once. Partitions will be exactly the same::
sage: p.nef_partitions(hodge_numbers=True) # long time (2s on sage.math, 2011) [ Nef-partition {0, 1, 4, 5} U {2, 3, 6, 7} (direct product), Nef-partition {0, 1, 2, 4} U {3, 5, 6, 7}, Nef-partition {0, 1, 2, 4, 5} U {3, 6, 7}, Nef-partition {0, 1, 2, 4, 5, 6} U {3, 7} (direct product), Nef-partition {0, 1, 2, 3} U {4, 5, 6, 7}, Nef-partition {0, 1, 2, 3, 4} U {5, 6, 7}, Nef-partition {0, 1, 2, 3, 4, 5} U {6, 7}, Nef-partition {0, 1, 2, 3, 4, 5, 6} U {7} (projection) ]
Now it is possible to get Hodge numbers::
sage: p.nef_partitions(hodge_numbers=True)[1].hodge_numbers() (20,)
Since nef-partitions are cached, their Hodge numbers are accessible after the first request, even if you do not specify ``hodge_numbers=True`` anymore::
sage: p.nef_partitions()[1].hodge_numbers() (20,)
We illustrate removal of symmetric partitions on a diamond::
sage: p = lattice_polytope.cross_polytope(2) sage: p.nef_partitions() [ Nef-partition {0, 2} U {1, 3} (direct product), Nef-partition {0, 1} U {2, 3}, Nef-partition {0, 1, 2} U {3} (projection) ] sage: p.nef_partitions(keep_symmetric=True) [ Nef-partition {0, 1, 3} U {2} (projection), Nef-partition {0, 2, 3} U {1} (projection), Nef-partition {0, 3} U {1, 2}, Nef-partition {1, 2, 3} U {0} (projection), Nef-partition {1, 3} U {0, 2} (direct product), Nef-partition {2, 3} U {0, 1}, Nef-partition {0, 1, 2} U {3} (projection) ]
Nef-partitions can be computed only for reflexive polytopes::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (0,0,2), ....: (-1,0,0), (0,-1,0), (0,0,-1)]) sage: p.nef_partitions() Traceback (most recent call last): ... ValueError: The given polytope is not reflexive! Polytope: 3-d lattice polytope in 3-d lattice M """ + "Polytope: %s") % self) or keep_symmetric and oldkeys.find("-s") == -1 or not keep_symmetric and oldkeys.find("-s") != -1 or keep_projections and oldkeys.find("-P") == -1 or keep_products and oldkeys.find("-D") == -1): # Select only necessary partitions if (keep_projections or not p._is_projection) and (keep_products or not p._is_product)], cr=True, check=False)
def nef_x(self, keys): r""" Run nef.x with given ``keys`` on vertices of this polytope.
INPUT:
- ``keys`` - a string of options passed to nef.x. The key "-f" is added automatically.
OUTPUT: the output of nef.x as a string.
EXAMPLES: This call is used internally for computing nef-partitions::
sage: o = lattice_polytope.cross_polytope(3) sage: s = o.nef_x("-N -V -p") sage: s # output contains random time M:27 8 N:7 6 codim=2 #part=5 3 6 Vertices of P: 1 0 0 -1 0 0 0 1 0 0 -1 0 0 0 1 0 0 -1 P:0 V:2 4 5 0sec 0cpu P:2 V:3 4 5 0sec 0cpu P:3 V:4 5 0sec 0cpu np=3 d:1 p:1 0sec 0cpu """
def nfacets(self): r""" Return the number of facets of this polytope.
EXAMPLES: The number of facets of the 3-dimensional octahedron::
sage: o = lattice_polytope.cross_polytope(3) sage: o.nfacets() 8
The number of facets of an interval is 2::
sage: LatticePolytope(([1],[2])).nfacets() 2
Now consider a 2-dimensional diamond in a 3-dimensional space::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p.nfacets() 4 """
@cached_method def normal_form(self, algorithm="palp", permutation=False): r""" Return the normal form of vertices of ``self``.
Two full-dimensional lattice polytopes are in the same ``GL(\mathbb{Z})``-orbit if and only if their normal forms are the same. Normal form is not defined and thus cannot be used for polytopes whose dimension is smaller than the dimension of the ambient space.
The original algorithm was presented in [KS1998]_ and implemented in PALP. A modified version of the PALP algorithm is discussed in [GK2013]_ and available here as "palp_modified".
INPUT:
- ``algorithm`` -- (default: "palp") The algorithm which is used to compute the normal form. Options are:
* "palp" -- Run external PALP code, usually the fastest option.
* "palp_native" -- The original PALP algorithm implemented in sage. Currently considerably slower than PALP.
* "palp_modified" -- A modified version of the PALP algorithm which determines the maximal vertex-facet pairing matrix first and then computes its automorphisms, while the PALP algorithm does both things concurrently.
- ``permutation`` -- (default: ``False``) If ``True`` the permutation applied to vertices to obtain the normal form is returned as well. Note that the different algorithms may return different results that nevertheless lead to the same normal form.
OUTPUT:
- a :class:`point collection <PointCollection>` in the :meth:`lattice` of ``self`` or a tuple of it and a permutation.
EXAMPLES:
We compute the normal form of the "diamond"::
sage: d = LatticePolytope([(1,0), (0,1), (-1,0), (0,-1)]) sage: d.vertices() M( 1, 0), M( 0, 1), M(-1, 0), M( 0, -1) in 2-d lattice M sage: d.normal_form() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M
The diamond is the 3rd polytope in the internal database::
sage: d.index() 3 sage: d 2-d reflexive polytope #3 in 2-d lattice M
You can get it in its normal form (in the default lattice) as ::
sage: lattice_polytope.ReflexivePolytope(2, 3).vertices() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M
It is not possible to compute normal forms for polytopes which do not span the space::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p.normal_form() Traceback (most recent call last): ... ValueError: normal form is not defined for 2-d lattice polytope in 3-d lattice M
We can perform the same examples using other algorithms::
sage: o = lattice_polytope.cross_polytope(2) sage: o.normal_form(algorithm="palp_native") M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M
sage: o = lattice_polytope.cross_polytope(2) sage: o.normal_form(algorithm="palp_modified") M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M """ StringIO(self.poly_x("N")), M, permutation=permutation) else: raise ValueError('Algorithm must be palp, ' + 'palp_native, or palp_modified.') vertices, perm = result else:
def _palp_modified_normal_form(self, permutation=False): r""" Return the normal form of ``self`` using the modified PALP algorithm.
This is a helper function for :meth:`normal_form` and should not be called directly. The modified PALP algorithm can be faster than the native algorithm in case the automorphism group of the vertex-facet pairing matrix is large.
INPUT:
- ``permutation`` -- a Boolean, whether to return the permutation of the order of the vertices that was applied to obtain this matrix.
OUTPUT:
A matrix or a tuple of a matrix and a permutation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(2) sage: o.vertices() M( 1, 0), M( 0, 1), M(-1, 0), M( 0, -1) in 2-d lattice M sage: o._palp_modified_normal_form() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M sage: o._palp_modified_normal_form(permutation=True) (M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M, (3,4)) """ for k, p in enumerate(permutations)} else:
def _palp_native_normal_form(self, permutation=False): r""" Return the normal form of ``self`` using the native PALP algorithm implemented in Sage.
This is a helper function for :meth:`normal_form` and should not be called directly.
INPUT:
- ``permutation`` -- a Boolean, whether to return the permutation of the order of the vertices that was applied to obtain this matrix.
OUTPUT:
A matrix or a tuple of a matrix and a permutation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(2) sage: o.vertices() M( 1, 0), M( 0, 1), M(-1, 0), M( 0, -1) in 2-d lattice M sage: o._palp_native_normal_form() M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M sage: o._palp_native_normal_form(permutation=True) (M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M, (1,3,2,4)) """ else:
def _palp_PM_max(self, check=False): r""" Compute the permutation normal form of the vertex facet pairing matrix .
The permutation normal form of a matrix is defined as the lexicographic maximum under all permutations of its rows and columns. For more more detail, see also :meth:`~sage.matrix.matrix2.Matrix.permutation_normal_form`.
Instead of using the generic method for computing the permutation normal form, this method uses the PALP algorithm to compute the permutation normal form and its automorphisms concurrently.
INPUT:
- ``check`` -- Boolean (default: ``False``), whether to return the permutations leaving the maximal vertex-facet pairing matrix invariant.
OUTPUT:
A matrix or a tuple of a matrix and a dict whose values are the permutation group elements corresponding to the permutations that permute :meth:`vertices` such that the vertex-facet pairing matrix is maximal.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(2) sage: PM = o.vertex_facet_pairing_matrix() sage: PM_max = PM.permutation_normal_form() sage: PM_max == o._palp_PM_max() True sage: P2 = ReflexivePolytope(2, 0) sage: PM_max, permutations = P2._palp_PM_max(check=True) sage: PM_max [3 0 0] [0 3 0] [0 0 3] sage: list(permutations.values()) [[(1,2,3), (1,2,3)], [(1,3,2), (1,3,2)], [(1,3), (1,3)], [(1,2), (1,2)], [(), ()], [(2,3), (2,3)]] sage: PM_max.automorphisms_of_rows_and_columns() [((), ()), ((2,3), (2,3)), ((1,2), (1,2)), ((1,3,2), (1,3,2)), ((1,2,3), (1,2,3)), ((1,3), (1,3))] sage: PMs = [i._palp_PM_max(check=True) ....: for i in ReflexivePolytopes(2)] # long time sage: all(len(i) == len(j.automorphisms_of_rows_and_columns()) ....: for j, i in PMs) # long time True """
# and find all the ways of making the first row of PM_max # returns the index of max of any iterable
PGE(range(1, n_v + 1))]} [(PM.with_permuted_columns(permutations[0][1]))[0][i] for i in range(j, n_v)])
# Arrange other rows one by one and compare with first row # Error for k == 1 already! - permutations[0][1](first_row)[0]) # The largest elt of this row is smaller than largest elt # in 1st row, so nothing to do # otherwise: [PM.with_permuted_columns(permutations[n_s][1])[k][j] for j in range(i, n_v)]) * permutations[n_s][1] -permutations[0][1](first_row)[i]) # This row is smaller than 1st row, so nothing to do # This row is the same, so we have a symmetry! else: # This row is larger, so it becomes the first row and # the symmetries reset.
# Work out the restrictions the current permutations # place on other permutations as a automorphisms # of the first row # The array is such that: # S = [i, 1, ..., 1 (ith), j, i+1, ..., i+1 (jth), k ... ] # describes the "symmetry blocks" else:
# We determine the other rows of PM_max in turn by use of perms and # aut on previous rows. # Search for possible local permutations based off previous # global permutations. # number of local permutations associated with current global # We look for the line with the maximal entry in the first # subsymmetry block, i.e. we are allowed to swap elements # between 0 and S(0) *permutations_bar[n_p])[s] *permutations_bar[n_p])[s][0] else: *permutations_bar[n_p])[s][0] # We move to the next line # Maximal values agree, so possible symmetry else: # We found a greater maximal value for first entry. # It becomes our new reference: # Forget previous work done # Check if the permutations found just now work # with other elements # Now let us find out where the end of the # next symmetry block is: # Check through this block for each possible permutation # Find the largest value in this symmetry block *permutations_bar[s])[l] # Set reference and carry on to next permutation *permutations_bar[s])[l][c] else: *permutations_bar[s])[l][c] # The current case leads to a smaller matrix, # hence this case becomes our new reference # Update permutations # If the automorphisms are not already completely restricted, # update them # Take the old automorphisms and update by # the restrictions the last worked out # row imposes. S[c] = S[c - 1] S[S[c] - 1] += 1 else: # Now we have the perms, we construct PM_max using one of them else:
def npoints(self): r""" Return the number of lattice points of this polytope.
EXAMPLES: The number of lattice points of the 3-dimensional octahedron and its polar cube::
sage: o = lattice_polytope.cross_polytope(3) sage: o.npoints() 7 sage: cube = o.polar() sage: cube.npoints() 27 """
def nvertices(self): r""" Return the number of vertices of this polytope.
EXAMPLES: The number of vertices of the 3-dimensional octahedron and its polar cube::
sage: o = lattice_polytope.cross_polytope(3) sage: o.nvertices() 6 sage: cube = o.polar() sage: cube.nvertices() 8 """
@cached_method def origin(self): r""" Return the index of the origin in the list of points of self.
OUTPUT:
- integer if the origin belongs to this polytope, ``None`` otherwise.
EXAMPLES::
sage: p = lattice_polytope.cross_polytope(2) sage: p.origin() 4 sage: p.point(p.origin()) M(0, 0)
sage: p = LatticePolytope(([1],[2])) sage: p.points() M(1), M(2) in 1-d lattice M sage: print(p.origin()) None
Now we make sure that the origin of non-full-dimensional polytopes can be identified correctly (:trac:`10661`)::
sage: LatticePolytope([(1,0,0), (-1,0,0)]).origin() 2 """
def parent(self): """ Return the set of all lattice polytopes.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.parent() Set of all Lattice Polytopes """
def plot3d(self, show_facets=True, facet_opacity=0.5, facet_color=(0,1,0), facet_colors=None, show_edges=True, edge_thickness=3, edge_color=(0.5,0.5,0.5), show_vertices=True, vertex_size=10, vertex_color=(1,0,0), show_points=True, point_size=10, point_color=(0,0,1), show_vindices=None, vindex_color=(0,0,0), vlabels=None, show_pindices=None, pindex_color=(0,0,0), index_shift=1.1): r""" Return a 3d-plot of this polytope.
Polytopes with ambient dimension 1 and 2 will be plotted along x-axis or in xy-plane respectively. Polytopes of dimension 3 and less with ambient dimension 4 and greater will be plotted in some basis of the spanned space.
By default, everything is shown with more or less pretty combination of size and color parameters.
INPUT: Most of the parameters are self-explanatory:
- ``show_facets`` - (default:True)
- ``facet_opacity`` - (default:0.5)
- ``facet_color`` - (default:(0,1,0))
- ``facet_colors`` - (default:None) if specified, must be a list of colors for each facet separately, used instead of ``facet_color``
- ``show_edges`` - (default:True) whether to draw edges as lines
- ``edge_thickness`` - (default:3)
- ``edge_color`` - (default:(0.5,0.5,0.5))
- ``show_vertices`` - (default:True) whether to draw vertices as balls
- ``vertex_size`` - (default:10)
- ``vertex_color`` - (default:(1,0,0))
- ``show_points`` - (default:True) whether to draw other poits as balls
- ``point_size`` - (default:10)
- ``point_color`` - (default:(0,0,1))
- ``show_vindices`` - (default:same as show_vertices) whether to show indices of vertices
- ``vindex_color`` - (default:(0,0,0)) color for vertex labels
- ``vlabels`` - (default:None) if specified, must be a list of labels for each vertex, default labels are vertex indices
- ``show_pindices`` - (default:same as show_points) whether to show indices of other points
- ``pindex_color`` - (default:(0,0,0)) color for point labels
- ``index_shift`` - (default:1.1)) if 1, labels are placed exactly at the corresponding points. Otherwise the label position is computed as a multiple of the point position vector.
EXAMPLES: The default plot of a cube::
sage: c = lattice_polytope.cross_polytope(3).polar() sage: c.plot3d() Graphics3d Object
Plot without facets and points, shown without the frame::
sage: c.plot3d(show_facets=false,show_points=false).show(frame=False)
Plot with facets of different colors::
sage: c.plot3d(facet_colors=rainbow(c.nfacets(), 'rgbtuple')) Graphics3d Object
It is also possible to plot lower dimensional polytops in 3D (let's also change labels of vertices)::
sage: lattice_polytope.cross_polytope(2).plot3d(vlabels=["A", "B", "C", "D"]) Graphics3d Object
TESTS::
sage: p = LatticePolytope([[0,0,0],[0,1,1],[1,0,1],[1,1,0]]) sage: p.plot3d() Graphics3d Object """ raise ValueError("%d-dimensional polytopes can not be plotted in 3D!" % self.dim()) return self._sublattice_polytope.plot3d( show_facets, facet_opacity, facet_color, facet_colors, show_edges, edge_thickness, edge_color, show_vertices, vertex_size, vertex_color, show_points, point_size, point_color, show_vindices, vindex_color, vlabels, show_pindices, pindex_color, index_shift) else: for i in range(self.nvertices())] for i in range(self.nvertices(), self.npoints())] vertices, opacity=facet_opacity, rgbcolor=facet_color) vertices, opacity=facet_opacity, rgbcolor=c) pplot += line3d(vertices, thickness=edge_thickness, rgbcolor=edge_color) else: thickness=edge_thickness, rgbcolor=edge_color) # Compute the barycenter and shift text of labels away from it
def polyhedron(self): r""" Return the Polyhedron object determined by this polytope's vertices.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(2) sage: o.polyhedron() A 2-dimensional polyhedron in ZZ^2 defined as the convex hull of 4 vertices """
def show3d(self): """ Show a 3d picture of the polytope with default settings and without axes or frame.
See self.plot3d? for more details.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o.show3d() """
def point(self, i): r""" Return the i-th point of this polytope, i.e. the i-th column of the matrix returned by points().
EXAMPLES: First few points are actually vertices::
sage: o = lattice_polytope.cross_polytope(3) sage: o.vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M sage: o.point(1) M(0, 1, 0)
The only other point in the octahedron is the origin::
sage: o.point(6) M(0, 0, 0) sage: o.points() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1), M( 0, 0, 0) in 3-d lattice M """
def points(self, *args, **kwds): r""" Return all lattice points of ``self``.
INPUT:
- any arguments given will be passed on to the returned object.
OUTPUT:
- a :class:`point collection <PointCollection>`.
EXAMPLES:
Lattice points of the octahedron and its polar cube::
sage: o = lattice_polytope.cross_polytope(3) sage: o.points() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1), M( 0, 0, 0) in 3-d lattice M sage: cube = o.polar() sage: cube.points() N( 1, -1, -1), N( 1, 1, -1), N( 1, 1, 1), N( 1, -1, 1), N(-1, -1, 1), N(-1, -1, -1), N(-1, 1, -1), N(-1, 1, 1), N(-1, -1, 0), N(-1, 0, -1), N(-1, 0, 0), N(-1, 0, 1), N(-1, 1, 0), N( 0, -1, -1), N( 0, -1, 0), N( 0, -1, 1), N( 0, 0, -1), N( 0, 0, 0), N( 0, 0, 1), N( 0, 1, -1), N( 0, 1, 0), N( 0, 1, 1), N( 1, -1, 0), N( 1, 0, -1), N( 1, 0, 0), N( 1, 0, 1), N( 1, 1, 0) in 3-d lattice N
Lattice points of a 2-dimensional diamond in a 3-dimensional space::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p.points() M( 1, 0, 0), M( 0, 1, 0), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, 0) in 3-d lattice M
Only two of the above points:
sage: p.points(1, 3) M(0, 1, 0), M(0, -1, 0) in 3-d lattice M
We check that points of a zero-dimensional polytope can be computed::
sage: p = LatticePolytope([[1]]) sage: p.points() M(1) in 1-d lattice M """ else: M, [m[i, j] for i in range(M.rank())]) else:
def polar(self): r""" Return the polar polytope, if this polytope is reflexive.
EXAMPLES: The polar polytope to the 3-dimensional octahedron::
sage: o = lattice_polytope.cross_polytope(3) sage: cube = o.polar() sage: cube 3-d reflexive polytope in 3-d lattice N
The polar polytope "remembers" the original one::
sage: cube.polar() 3-d reflexive polytope in 3-d lattice M sage: cube.polar().polar() is cube True
Only reflexive polytopes have polars::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (0,0,2), ....: (-1,0,0), (0,-1,0), (0,0,-1)]) sage: p.polar() Traceback (most recent call last): ... ValueError: The given polytope is not reflexive! Polytope: 3-d lattice polytope in 3-d lattice M """ else: + "Polytope: %s") % self)
def poly_x(self, keys, reduce_dimension=False): r""" Run poly.x with given ``keys`` on vertices of this polytope.
INPUT:
- ``keys`` - a string of options passed to poly.x. The key "f" is added automatically.
- ``reduce_dimension`` - (default: False) if ``True`` and this polytope is not full-dimensional, poly.x will be called for the vertices of this polytope in some basis of the spanned affine space.
OUTPUT: the output of poly.x as a string.
EXAMPLES: This call is used for determining if a polytope is reflexive or not::
sage: o = lattice_polytope.cross_polytope(3) sage: print(o.poly_x("e")) 8 3 Vertices of P-dual <-> Equations of P -1 -1 1 1 -1 1 -1 1 1 1 1 1 -1 -1 -1 1 -1 -1 -1 1 -1 1 1 -1
Since PALP has limits on different parameters determined during compilation, the following code is likely to fail, unless you change default settings of PALP::
sage: BIG = lattice_polytope.cross_polytope(7) sage: BIG 7-d reflexive polytope in 7-d lattice M sage: BIG.poly_x("e") # possibly different output depending on your system Traceback (most recent call last): ... ValueError: Error executing 'poly.x -fe' for the given polytope! Output: Please increase POLY_Dmax to at least 7
You cannot call poly.x for polytopes that don't span the space (if you could, it would crush anyway)::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: p.poly_x("e") Traceback (most recent call last): ... ValueError: Cannot run PALP for a 2-dimensional polytope in a 3-dimensional space!
But if you know what you are doing, you can call it for the polytope in some basis of the spanned space::
sage: print(p.poly_x("e", reduce_dimension=True)) 4 2 Equations of P -1 1 0 1 1 2 -1 -1 0 1 -1 2 """
@cached_method def skeleton(self): r""" Return the graph of the one-skeleton of this polytope.
EXAMPLES::
sage: d = lattice_polytope.cross_polytope(2) sage: g = d.skeleton() sage: g Graph on 4 vertices sage: g.edges() [(0, 1, None), (0, 3, None), (1, 2, None), (2, 3, None)] """
def skeleton_points(self, k=1): r""" Return the increasing list of indices of lattice points in k-skeleton of the polytope (k is 1 by default).
EXAMPLES: We compute all skeleton points for the cube::
sage: o = lattice_polytope.cross_polytope(3) sage: c = o.polar() sage: c.skeleton_points() [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 11, 12, 13, 15, 19, 21, 22, 23, 25, 26]
The default was 1-skeleton::
sage: c.skeleton_points(k=1) [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 11, 12, 13, 15, 19, 21, 22, 23, 25, 26]
0-skeleton just lists all vertices::
sage: c.skeleton_points(k=0) [0, 1, 2, 3, 4, 5, 6, 7]
2-skeleton lists all points except for the origin (point #17)::
sage: c.skeleton_points(k=2) [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 18, 19, 20, 21, 22, 23, 24, 25, 26]
3-skeleton includes all points::
sage: c.skeleton_points(k=3) [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26]
It is OK to compute higher dimensional skeletons - you will get the list of all points::
sage: c.skeleton_points(k=100) [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26] """
def skeleton_show(self, normal=None): r"""Show the graph of one-skeleton of this polytope. Works only for polytopes in a 3-dimensional space.
INPUT:
- ``normal`` - a 3-dimensional vector (can be given as a list), which should be perpendicular to the screen. If not given, will be selected randomly (new each time and it may be far from "nice").
EXAMPLES: Show a pretty picture of the octahedron::
sage: o = lattice_polytope.cross_polytope(3) sage: o.skeleton_show([1,2,4])
Does not work for a diamond at the moment::
sage: d = lattice_polytope.cross_polytope(2) sage: d.skeleton_show() Traceback (most recent call last): ... NotImplementedError: skeleton view is implemented only in 3-d space """ normal = [ZZ.random_element(20),ZZ.random_element(20),ZZ.random_element(20)]
def traverse_boundary(self): r""" Return a list of indices of vertices of a 2-dimensional polytope in their boundary order.
Needed for plot3d function of polytopes.
EXAMPLES:
sage: p = lattice_polytope.cross_polytope(2).polar() sage: p.traverse_boundary() [3, 0, 1, 2] """ raise ValueError("Boundary can be traversed only for 2-polytopes!")
def vertex(self, i): r""" Return the i-th vertex of this polytope, i.e. the i-th column of the matrix returned by vertices().
EXAMPLES: Note that numeration starts with zero::
sage: o = lattice_polytope.cross_polytope(3) sage: o.vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M sage: o.vertex(3) M(-1, 0, 0) """
def vertex_facet_pairing_matrix(self): r""" Return the vertex facet pairing matrix `PM`.
Return a matrix whose the `i, j^\text{th}` entry is the height of the `j^\text{th}` vertex over the `i^\text{th}` facet. The ordering of the vertices and facets is as in :meth:`vertices` and :meth:`facets`.
EXAMPLES::
sage: L = lattice_polytope.cross_polytope(3) sage: L.vertex_facet_pairing_matrix() [2 0 0 0 2 2] [2 2 0 0 0 2] [2 2 2 0 0 0] [2 0 2 0 2 0] [0 0 2 2 2 0] [0 0 0 2 2 2] [0 2 0 2 0 2] [0 2 2 2 0 0] """ for n, c in zip(self.facet_normals(), self.facet_constants())])
def vertices(self, *args, **kwds): r""" Return vertices of ``self``.
INPUT:
- any arguments given will be passed on to the returned object.
OUTPUT:
- a :class:`point collection <PointCollection>`.
EXAMPLES:
Vertices of the octahedron and its polar cube are in dual lattices::
sage: o = lattice_polytope.cross_polytope(3) sage: o.vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M sage: cube = o.polar() sage: cube.vertices() N( 1, -1, -1), N( 1, 1, -1), N( 1, 1, 1), N( 1, -1, 1), N(-1, -1, 1), N(-1, -1, -1), N(-1, 1, -1), N(-1, 1, 1) in 3-d lattice N """ else:
def is_NefPartition(x): r""" Check if ``x`` is a nef-partition.
INPUT:
- ``x`` -- anything.
OUTPUT:
- ``True`` if ``x`` is a :class:`nef-partition <NefPartition>` and ``False`` otherwise.
EXAMPLES::
sage: from sage.geometry.lattice_polytope import is_NefPartition sage: is_NefPartition(1) False sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: is_NefPartition(np) True """
class NefPartition(SageObject, collections.Hashable): r""" Create a nef-partition.
INPUT:
- ``data`` -- a list of integers, the $i$-th element of this list must be the part of the $i$-th vertex of ``Delta_polar`` in this nef-partition;
- ``Delta_polar`` -- a :class:`lattice polytope <sage.geometry.lattice_polytope.LatticePolytopeClass>`;
- ``check`` -- by default the input will be checked for correctness, i.e. that ``data`` indeed specify a nef-partition. If you are sure that the input is correct, you can speed up construction via ``check=False`` option.
OUTPUT:
- a nef-partition of ``Delta_polar``.
Let $M$ and $N$ be dual lattices. Let $\Delta \subset M_\RR$ be a reflexive polytope with polar $\Delta^\circ \subset N_\RR$. Let $X_\Delta$ be the toric variety associated to the normal fan of $\Delta$. A **nef-partition** is a decomposition of the vertex set $V$ of $\Delta^\circ$ into a disjoint union $V = V_0 \sqcup V_1 \sqcup \dots \sqcup V_{k-1}$ such that divisors $E_i = \sum_{v\in V_i} D_v$ are Cartier (here $D_v$ are prime torus-invariant Weil divisors corresponding to vertices of $\Delta^\circ$). Equivalently, let $\nabla_i \subset N_\RR$ be the convex hull of vertices from $V_i$ and the origin. These polytopes form a nef-partition if their Minkowski sum $\nabla \subset N_\RR$ is a reflexive polytope.
The **dual nef-partition** is formed by polytopes $\Delta_i \subset M_\RR$ of $E_i$, which give a decomposition of the vertex set of $\nabla^\circ \subset M_\RR$ and their Minkowski sum is $\Delta$, i.e. the polar duality of reflexive polytopes switches convex hull and Minkowski sum for dual nef-partitions:
.. MATH::
\Delta^\circ &= \mathrm{Conv} \left(\nabla_0, \nabla_1, \dots, \nabla_{k-1}\right), \\ \nabla^{\phantom{\circ}} &= \nabla_0 + \nabla_1 + \dots + \nabla_{k-1}, \\ & \\ \Delta^{\phantom{\circ}} &= \Delta_0 + \Delta_1 + \dots + \Delta_{k-1}, \\ \nabla^\circ &= \mathrm{Conv} \left(\Delta_0, \Delta_1, \dots, \Delta_{k-1}\right).
One can also interpret the duality of nef-partitions as the duality of the associated cones. Below $\overline{M} = M \times \ZZ^k$ and $\overline{N} = N \times \ZZ^k$ are dual lattices.
The **Cayley polytope** $P \subset \overline{M}_\RR$ of a nef-partition is given by $P = \mathrm{Conv}(\Delta_0 \times e_0, \Delta_1 \times e_1, \ldots, \Delta_{k-1} \times e_{k-1})$, where $\{e_i\}_{i=0}^{k-1}$ is the standard basis of $\ZZ^k$. The **dual Cayley polytope** $P^* \subset \overline{N}_\RR$ is the Cayley polytope of the dual nef-partition.
The **Cayley cone** $C \subset \overline{M}_\RR$ of a nef-partition is the cone spanned by its Cayley polytope. The **dual Cayley cone** $C^\vee \subset \overline{M}_\RR$ is the usual dual cone of $C$. It turns out, that $C^\vee$ is spanned by $P^*$.
It is also possible to go back from the Cayley cone to the Cayley polytope, since $C$ is a reflexive Gorenstein cone supported by $P$: primitive integral ray generators of $C$ are contained in an affine hyperplane and coincide with vertices of $P$.
See Section 4.3.1 in [CK1999]_ and references therein for further details, or [BN2008]_ for a purely combinatorial approach.
EXAMPLES:
It is very easy to create a nef-partition for the octahedron, since for this polytope any decomposition of vertices is a nef-partition. We create a 3-part nef-partition with the 0-th and 1-st vertices belonging to the 0-th part (recall that numeration in Sage starts with 0), the 2-nd and 5-th vertices belonging to the 1-st part, and 3-rd and 4-th vertices belonging to the 2-nd part::
sage: o = lattice_polytope.cross_polytope(3) sage: np = NefPartition([0,0,1,2,2,1], o) sage: np Nef-partition {0, 1} U {2, 5} U {3, 4}
The octahedron plays the role of `\Delta^\circ` in the above description::
sage: np.Delta_polar() is o True
The dual nef-partition (corresponding to the "mirror complete intersection") gives decomposition of the vertex set of `\nabla^\circ`::
sage: np.dual() Nef-partition {0, 1, 2} U {3, 4} U {5, 6, 7} sage: np.nabla_polar().vertices() N(-1, -1, 0), N(-1, 0, 0), N( 0, -1, 0), N( 0, 0, -1), N( 0, 0, 1), N( 1, 0, 0), N( 0, 1, 0), N( 1, 1, 0) in 3-d lattice N
Of course, `\nabla^\circ` is `\Delta^\circ` from the point of view of the dual nef-partition::
sage: np.dual().Delta_polar() is np.nabla_polar() True sage: np.Delta(1).vertices() N(0, 0, -1), N(0, 0, 1) in 3-d lattice N sage: np.dual().nabla(1).vertices() N(0, 0, -1), N(0, 0, 1) in 3-d lattice N
Instead of constructing nef-partitions directly, you can request all 2-part nef-partitions of a given reflexive polytope (they will be computed using ``nef.x`` program from PALP)::
sage: o.nef_partitions() [ Nef-partition {0, 1, 3} U {2, 4, 5}, Nef-partition {0, 1, 3, 4} U {2, 5} (direct product), Nef-partition {0, 1, 2} U {3, 4, 5}, Nef-partition {0, 1, 2, 3} U {4, 5}, Nef-partition {0, 1, 2, 3, 4} U {5} (projection) ] """
def __init__(self, data, Delta_polar, check=True): r""" See :class:`NefPartition` for documentation.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: TestSuite(np).run() """ raise ValueError("nef-partitions can be constructed for reflexive " "polytopes ony!") raise ValueError("%s do not form a nef-partition!" % str(data))
def __eq__(self, other): r""" Compare ``self`` with ``other``.
INPUT:
- ``other`` -- anything.
OUTPUT:
- ``True`` if ``other`` is a :class:`nef-partition <NefPartition>` equal to ``self``, ``False`` otherwise.
.. NOTE::
Two nef-partitions are equal if they correspond to equal polytopes and their parts are the same, including their order.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np == np True sage: np == o.nef_partitions()[1] False sage: np2 = NefPartition(np._vertex_to_part, o) sage: np2 is np False sage: np2 == np True sage: np == 0 False """ and self._Delta_polar == other._Delta_polar and self._vertex_to_part == other._vertex_to_part)
def __hash__(self): r""" Return the hash of ``self``.
OUTPUT:
- an integer.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: hash(np) == hash(np) True """
def __ne__(self, other): r""" Compare ``self`` with ``other``.
INPUT:
- ``other`` -- anything.
OUTPUT:
- ``False`` if ``other`` is a :class:`nef-partition <NefPartition>` equal to ``self``, ``True`` otherwise.
.. NOTE::
Two nef-partitions are equal if they correspond to equal polytopes and their parts are the same, including their order.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np != np False sage: np != o.nef_partitions()[1] True sage: np2 = NefPartition(np._vertex_to_part, o) sage: np2 is np False sage: np2 != np False sage: np != 0 True """
def _latex_(self): r""" Return a LaTeX representation of ``self``.
OUTPUT:
- a string.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: latex(np) # indirect doctest \text{Nef-partition } \{0, 1, 3\} \sqcup \{2, 4, 5\} """ # We may or may not know the type of the partition result += r" \text{ (direct product)}" result += r" \text{ (projection)}" except AttributeError: pass
def _repr_(self): r""" Return a string representation of ``self``.
OUTPUT:
- a string.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: repr(np) # indirect doctest 'Nef-partition {0, 1, 3} U {2, 4, 5}' """ # We may or may not know the type of the partition
def Delta(self, i=None): r""" Return the polytope $\Delta$ or $\Delta_i$ corresponding to ``self``.
INPUT:
- ``i`` -- an integer. If not given, $\Delta$ will be returned.
OUTPUT:
- a :class:`lattice polytope <LatticePolytopeClass>`.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.Delta().polar() is o True sage: np.Delta().vertices() N( 1, -1, -1), N( 1, 1, -1), N( 1, 1, 1), N( 1, -1, 1), N(-1, -1, 1), N(-1, -1, -1), N(-1, 1, -1), N(-1, 1, 1) in 3-d lattice N sage: np.Delta(0).vertices() N(-1, -1, 0), N(-1, 0, 0), N( 1, 0, 0), N( 1, -1, 0) in 3-d lattice N """ else:
def Delta_polar(self): r""" Return the polytope $\Delta^\circ$ corresponding to ``self``.
OUTPUT:
- a :class:`lattice polytope <LatticePolytopeClass>`.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.Delta_polar() is o True """
def Deltas(self): r""" Return the polytopes $\Delta_i$ corresponding to ``self``.
OUTPUT:
- a tuple of :class:`lattice polytopes <LatticePolytopeClass>`.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.Delta().vertices() N( 1, -1, -1), N( 1, 1, -1), N( 1, 1, 1), N( 1, -1, 1), N(-1, -1, 1), N(-1, -1, -1), N(-1, 1, -1), N(-1, 1, 1) in 3-d lattice N sage: [Delta_i.vertices() for Delta_i in np.Deltas()] [N(-1, -1, 0), N(-1, 0, 0), N( 1, 0, 0), N( 1, -1, 0) in 3-d lattice N, N(0, 0, -1), N(0, 1, 1), N(0, 0, 1), N(0, 1, -1) in 3-d lattice N] sage: np.nabla_polar().vertices() N(-1, -1, 0), N( 1, -1, 0), N( 1, 0, 0), N(-1, 0, 0), N( 0, 1, -1), N( 0, 1, 1), N( 0, 0, 1), N( 0, 0, -1) in 3-d lattice N """
@cached_method def dual(self): r""" Return the dual nef-partition.
OUTPUT:
- a :class:`nef-partition <NefPartition>`.
See the class documentation for the definition.
ALGORITHM:
See Proposition 3.19 in [BN2008]_.
.. NOTE::
Automatically constructed dual nef-partitions will be ordered, i.e. vertex partition of `\nabla` will look like `\{0, 1, 2\} \sqcup \{3, 4, 5, 6\} \sqcup \{7, 8\}`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.dual() Nef-partition {0, 1, 2, 3} U {4, 5, 6, 7} sage: np.dual().Delta() is np.nabla() True sage: np.dual().nabla(0) is np.Delta(0) True """ # Delta and nabla are interchanged compared to [BN2008]_. # The order of vertices of this nabla_polar will be adjusted. reduce(minkowski_sum, (nabla.vertices() for nabla in self.nablas())), lattice=self._Delta_polar.lattice()).polar() # Make dual look "ordered", like {0,1,2} U {3,4,5,6} U {7,8}. compute_vertices=False) # If self is a valid nef-partition, the dual is as well.
def hodge_numbers(self): r""" Return Hodge numbers corresponding to ``self``.
OUTPUT:
- a tuple of integers (produced by ``nef.x`` program from PALP).
EXAMPLES:
Currently, you need to request Hodge numbers when you compute nef-partitions::
sage: p = lattice_polytope.cross_polytope(5) sage: np = p.nef_partitions()[0] # long time (4s on sage.math, 2011) sage: np.hodge_numbers() # long time Traceback (most recent call last): ... NotImplementedError: use nef_partitions(hodge_numbers=True)! sage: np = p.nef_partitions(hodge_numbers=True)[0] # long time (13s on sage.math, 2011) sage: np.hodge_numbers() # long time (19, 19) """ return self._hodge_numbers
def nabla(self, i=None): r""" Return the polytope $\nabla$ or $\nabla_i$ corresponding to ``self``.
INPUT:
- ``i`` -- an integer. If not given, $\nabla$ will be returned.
OUTPUT:
- a :class:`lattice polytope <LatticePolytopeClass>`.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.Delta_polar().vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M sage: np.nabla(0).vertices() M(-1, 0, 0), M( 1, 0, 0), M( 0, 1, 0) in 3-d lattice M sage: np.nabla().vertices() M(-1, 0, 1), M(-1, 0, -1), M( 1, 0, 1), M( 1, 0, -1), M( 0, 1, 1), M( 0, 1, -1), M( 1, -1, 0), M(-1, -1, 0) in 3-d lattice M """ else:
def nabla_polar(self): r""" Return the polytope $\nabla^\circ$ corresponding to ``self``.
OUTPUT:
- a :class:`lattice polytope <LatticePolytopeClass>`.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.nabla_polar().vertices() N(-1, -1, 0), N( 1, -1, 0), N( 1, 0, 0), N(-1, 0, 0), N( 0, 1, -1), N( 0, 1, 1), N( 0, 0, 1), N( 0, 0, -1) in 3-d lattice N sage: np.nabla_polar() is np.dual().Delta_polar() True """
def nablas(self): r""" Return the polytopes $\nabla_i$ corresponding to ``self``.
OUTPUT:
- a tuple of :class:`lattice polytopes <LatticePolytopeClass>`.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.Delta_polar().vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M sage: [nabla_i.vertices() for nabla_i in np.nablas()] [M(-1, 0, 0), M( 1, 0, 0), M( 0, 1, 0) in 3-d lattice M, M(0, -1, 0), M(0, 0, -1), M(0, 0, 1) in 3-d lattice M] """ [Delta_polar.vertex(j) for j in part] + origin, lattice=Delta_polar.lattice()) for part in self.parts())
def nparts(self): r""" Return the number of parts in ``self``.
OUTPUT:
- an integer.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.nparts() 2 """
def part(self, i, all_points=False): r""" Return the ``i``-th part of ``self``.
INPUT:
- ``i`` -- an integer
- ``all_points`` -- (default: False) whether to list all lattice points or just vertices
OUTPUT:
- a tuple of integers, indices of vertices (or all lattice points) of $\Delta^\circ$ belonging to $V_i$.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.part(0) (0, 1, 3) sage: np.part(0, all_points=True) (0, 1, 3) sage: np.dual().part(0) (0, 1, 2, 3) sage: np.dual().part(0, all_points=True) (0, 1, 2, 3, 8) """
@cached_method def parts(self, all_points=False): r""" Return all parts of ``self``.
INPUT:
- ``all_points`` -- (default: False) whether to list all lattice points or just vertices
OUTPUT:
- a tuple of tuples of integers. The $i$-th tuple contains indices of vertices (or all lattice points) of $\Delta^\circ$ belonging to $V_i$
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.parts() ((0, 1, 3), (2, 4, 5)) sage: np.parts(all_points=True) ((0, 1, 3), (2, 4, 5)) sage: np.dual().parts() ((0, 1, 2, 3), (4, 5, 6, 7)) sage: np.dual().parts(all_points=True) ((0, 1, 2, 3, 8), (4, 5, 6, 7, 10)) """ else:
def part_of(self, i): r""" Return the index of the part containing the ``i``-th vertex.
INPUT:
- ``i`` -- an integer.
OUTPUT:
- an integer $j$ such that the ``i``-th vertex of $\Delta^\circ$ belongs to $V_j$.
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: np = o.nef_partitions()[0] sage: np Nef-partition {0, 1, 3} U {2, 4, 5} sage: np.part_of(3) 0 sage: np.part_of(2) 1 """
@cached_method def part_of_point(self, i): r""" Return the index of the part containing the ``i``-th point.
INPUT:
- ``i`` -- an integer.
OUTPUT:
- an integer `j` such that the ``i``-th point of `\Delta^\circ` belongs to `\nabla_j`.
.. NOTE::
Since a nef-partition induces a partition on the set of boundary lattice points of `\Delta^\circ`, the value of `j` is well-defined for all `i` but the one that corresponds to the origin, in which case this method will raise a ``ValueError`` exception. (The origin always belongs to all `\nabla_j`.)
See :class:`nef-partition <NefPartition>` class documentation for definitions and notation.
EXAMPLES:
We consider a relatively complicated reflexive polytope #2252 (easily accessible in Sage as ``ReflexivePolytope(3, 2252)``, we create it here explicitly to avoid loading the whole database)::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (0,0,1), (0,1,-1), ....: (0,-1,1), (-1,1,0), (0,-1,-1), (-1,-1,0), (-1,-1,2)]) sage: np = p.nef_partitions()[0] sage: np Nef-partition {1, 2, 5, 7, 8} U {0, 3, 4, 6} sage: p.nvertices() 9 sage: p.npoints() 15
We see that the polytope has 6 more points in addition to vertices. One of them is the origin::
sage: p.origin() 14 sage: np.part_of_point(14) Traceback (most recent call last): ... ValueError: the origin belongs to all parts!
But the remaining 5 are partitioned by ``np``::
sage: [n for n in range(p.npoints()) ....: if p.origin() != n and np.part_of_point(n) == 0] [1, 2, 5, 7, 8, 9, 11, 13] sage: [n for n in range(p.npoints()) ....: if p.origin() != n and np.part_of_point(n) == 1] [0, 3, 4, 6, 10, 12] """
_palp_dimension = None
def _palp(command, polytopes, reduce_dimension=False): r""" Run ``command`` on vertices of given ``polytopes``.
Returns the name of the file containing the output of ``command``. You should delete it after using.
.. note::
PALP cannot be called for polytopes that do not span the ambient space. If you specify ``reduce_dimension=True`` argument, PALP will be called for vertices of this polytope in some basis of the affine space it spans.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: result_name = lattice_polytope._palp("poly.x -f", [o]) sage: f = open(result_name) sage: f.readlines() ['M:7 6 N:27 8 Pic:17 Cor:0\n'] sage: f.close() sage: os.remove(result_name)
sage: p = LatticePolytope([(1,0,0), (0,1,0), (-1,0,0), (0,-1,0)]) sage: lattice_polytope._palp("poly.x -f", [p]) Traceback (most recent call last): ValueError: Cannot run PALP for a 2-dimensional polytope in a 3-dimensional space!
sage: result_name = lattice_polytope._palp("poly.x -f", [p], reduce_dimension=True) sage: f = open(result_name) sage: f.readlines() ['M:5 4 F:4\n'] sage: f.close() sage: os.remove(result_name) """ raise ValueError(("Cannot run \"%s\" for the zero-dimensional " + "polytope!\nPolytope: %s") % (command, p)) "in a %d-dimensional space!") % (p.dim(), p.lattice_dim())) else: stdin=subprocess.PIPE, stdout=subprocess.PIPE, stderr=subprocess.PIPE, close_fds=True) raise RuntimeError(("Error executing \"%s\" for a polytope sequence!" + "\nOutput:\n%s") % (command, err)) except OSError: pass
def _palp_canonical_order(V, PM_max, permutations): r""" Compute the PALP normal form of the vertices V using auxiliary data computed elsewhere.
This is a helper function for :meth:`~sage.geometry.lattice_polytope.LatticePolytopeClass.normal_form` and should not be called directly.
Given a matrix of vertices, the maximal vertex-facet pairing matrix and the permutations realizing this matrix, apply the last part of the PALP algorithm and return the normal form.
INPUT:
- ``V`` -- :class:`point collection <PointCollection>`. The vertices.
- ``PM_max`` -- the maximal vertex-facet pairing matrix
- ``permutation`` -- the permutations of the vertices yielding ``PM_max``.
OUTPUT:
The PALP normal form as a :class:`point collection <PointCollection>`.
TESTS::
sage: L = lattice_polytope.cross_polytope(2) sage: V = L.vertices() sage: PM_max, permutations = L._palp_PM_max(check=True) sage: from sage.geometry.lattice_polytope import _palp_canonical_order sage: _palp_canonical_order(V, PM_max, permutations) (M( 1, 0), M( 0, 1), M( 0, -1), M(-1, 0) in 2-d lattice M, (1,3,2,4)) """ (M_max[j] == M_max[k] and S_max[j] < S_max[k]): k = j M_max[i], M_max[k] = M_max[k], M_max[i] S_max[i], S_max[k] = S_max[k], S_max[i] p_c = PermutationGroupElement((1 + i, 1 + k))*p_c # Create array of possible NFs. for k in permutations]
def _palp_convert_permutation(permutation): r""" Convert a permutation from PALPs notation to a PermutationGroupElement.
PALP specifies a permutation group element by its domain. Furthermore, it only supports permutations of up to 62 objects and labels these by `0 \dots 9`, `a \dots z`, and `A \dots Z`.
INPUT:
- ``permutation`` -- A string specifying a PALP style permutation.
OUTPUT:
A :class:`permutation group element <sage.groups.perm_gps.permgroup_element.PermmutationGroupElement>`.
EXAMPLES::
sage: from sage.geometry.lattice_polytope import _palp_convert_permutation sage: _palp_convert_permutation('1023') (1,2) sage: _palp_convert_permutation('0123456789bac') (11,12) """ else: elif o in range(65, 91): i = o - 28 else: raise ValueError('Cannot convert PALP index ' + i + ' to number.')
def _read_nef_x_partitions(data): r""" Read all nef-partitions for one polytope from a string or an open file.
``data`` should be an output of nef.x.
Returns the sequence of nef-partitions. Each nef-partition is given as a sequence of integers.
If there are no nef-partitions, returns the empty sequence. If the string is empty or EOF is reached, raises ValueError.
TESTS::
sage: o = lattice_polytope.cross_polytope(3) sage: s = o.nef_x("-N -p") sage: print(s) # random M:27 8 N:7 6 codim=2 #part=5 P:0 V:2 4 5 0sec 0cpu P:2 V:3 4 5 0sec 0cpu P:3 V:4 5 0sec 0cpu np=3 d:1 p:1 0sec 0cpu sage: lattice_polytope._read_nef_x_partitions(s) [[2, 4, 5], [3, 4, 5], [4, 5]] """ raise ValueError("Empty file!") # Compare the number of found partitions with np in data. raise ValueError("Found %d partitions, expected %d!" % (len(partitions), np)) else: raise ValueError("Wrong data format, cannot find \"np=\"!")
def _read_poly_x_incidences(data, dim): r""" Convert incidence data from binary numbers to sequences.
INPUT:
- ``data`` - an opened file with incidence information. The first line will be skipped, each consecutive line contains incidence information for all faces of one dimension, the first word of each line is a comment and is dropped.
- ``dim`` - dimension of the polytope.
OUTPUT: a sequence F, such that F[d][i] is a sequence of vertices or facets corresponding to the i-th d-dimensional face.
TESTS::
sage: p = lattice_polytope.cross_polytope(2) sage: result_name = lattice_polytope._palp("poly.x -fi", [p]) sage: with open(result_name) as f: ....: print(f.read()) Incidences as binary numbers [F-vector=(4 4)]: v[d][i]: sum_j Incidence(i'th dim-d-face, j-th vertex) x 2^j v[0]: 1000 0001 0100 0010 v[1]: 1001 1100 0011 0110 f[d][i]: sum_j Incidence(i'th dim-d-face, j-th facet) x 2^j f[0]: 0011 0101 1010 1100 f[1]: 0001 0010 0100 1000 sage: f = open(result_name) sage: l = f.readline() sage: lattice_polytope._read_poly_x_incidences(f, 2) [[[3], [0], [2], [1]], [[0, 3], [2, 3], [0, 1], [1, 2]]] sage: f.close() sage: os.remove(result_name) """ raise ValueError("Not enough data!")
def all_cached_data(polytopes): r""" Compute all cached data for all given ``polytopes`` and their polars.
This functions does it MUCH faster than member functions of ``LatticePolytope`` during the first run. So it is recommended to use this functions if you work with big sets of data. None of the polytopes in the given sequence should be constructed as the polar polytope to another one.
INPUT: a sequence of lattice polytopes.
EXAMPLES: This function has no output, it is just a fast way to work with long sequences of polytopes. Of course, you can use short sequences as well::
sage: o = lattice_polytope.cross_polytope(3) sage: lattice_polytope.all_cached_data([o]) """
def all_nef_partitions(polytopes, keep_symmetric=False): r""" Compute nef-partitions for all given ``polytopes``.
This functions does it MUCH faster than member functions of ``LatticePolytope`` during the first run. So it is recommended to use this functions if you work with big sets of data.
Note: member function ``is_reflexive`` will be called separately for each polytope. It is strictly recommended to call ``all_polars`` on the sequence of ``polytopes`` before using this function.
INPUT: a sequence of lattice polytopes.
EXAMPLES: This function has no output, it is just a fast way to work with long sequences of polytopes. Of course, you can use short sequences as well::
sage: o = lattice_polytope.cross_polytope(3) sage: lattice_polytope.all_nef_partitions([o]) sage: o.nef_partitions() [ Nef-partition {0, 1, 3} U {2, 4, 5}, Nef-partition {0, 1, 3, 4} U {2, 5} (direct product), Nef-partition {0, 1, 2} U {3, 4, 5}, Nef-partition {0, 1, 2, 3} U {4, 5}, Nef-partition {0, 1, 2, 3, 4} U {5} (projection) ]
You cannot use this function for non-reflexive polytopes::
sage: p = LatticePolytope([(1,0,0), (0,1,0), (0,0,2), ....: (-1,0,0), (0,-1,0), (0,0,-1)]) sage: lattice_polytope.all_nef_partitions([o, p]) Traceback (most recent call last): ... ValueError: nef-partitions can be computed for reflexive polytopes only """ keys += " -s" "polytopes only")
def all_points(polytopes): r""" Compute lattice points for all given ``polytopes``.
This functions does it MUCH faster than member functions of ``LatticePolytope`` during the first run. So it is recommended to use this functions if you work with big sets of data.
INPUT: a sequence of lattice polytopes.
EXAMPLES: This function has no output, it is just a fast way to work with long sequences of polytopes. Of course, you can use short sequences as well::
sage: o = lattice_polytope.cross_polytope(3) sage: lattice_polytope.all_points([o]) sage: o.points() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1), M( 0, 0, 0) in 3-d lattice M """ else: m = p._embed(read_palp_matrix(result)) if m.ncols() == nv: p._points = p.vertices() else: points = list(p.vertices()) for j in range(nv, m.ncols()): current = M.element_class( M, [m[i, j] for i in range(M.rank())]) current.set_immutable() points.append(current) p._points = PointCollection(points, M)
def all_polars(polytopes): r""" Compute polar polytopes for all reflexive and equations of facets for all non-reflexive ``polytopes``.
``all_facet_equations`` and ``all_polars`` are synonyms.
This functions does it MUCH faster than member functions of ``LatticePolytope`` during the first run. So it is recommended to use this functions if you work with big sets of data.
INPUT: a sequence of lattice polytopes.
EXAMPLES: This function has no output, it is just a fast way to work with long sequences of polytopes. Of course, you can use short sequences as well::
sage: o = lattice_polytope.cross_polytope(3) sage: lattice_polytope.all_polars([o]) sage: o.polar() 3-d reflexive polytope in 3-d lattice N """
# Synonym for the above function all_facet_equations = all_polars
def convex_hull(points): r""" Compute the convex hull of the given points.
.. note::
``points`` might not span the space. Also, it fails for large numbers of vertices in dimensions 4 or greater
INPUT:
- ``points`` - a list that can be converted into vectors of the same dimension over ZZ.
OUTPUT: list of vertices of the convex hull of the given points (as vectors).
EXAMPLES: Let's compute the convex hull of several points on a line in the plane::
sage: lattice_polytope.convex_hull([[1,2],[3,4],[5,6],[7,8]]) [(1, 2), (7, 8)] """ return [] return [p0] else:
def cross_polytope(dim): r""" Return a cross-polytope of the given dimension.
INPUT:
- ``dim`` -- an integer.
OUTPUT:
- a :class:`lattice polytope <LatticePolytopeClass>`.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: o 3-d reflexive polytope in 3-d lattice M sage: o.vertices() M( 1, 0, 0), M( 0, 1, 0), M( 0, 0, 1), M(-1, 0, 0), M( 0, -1, 0), M( 0, 0, -1) in 3-d lattice M """
def minkowski_sum(points1, points2): r""" Compute the Minkowski sum of two convex polytopes.
.. note::
Polytopes might not be of maximal dimension.
INPUT:
- ``points1, points2`` - lists of objects that can be converted into vectors of the same dimension, treated as vertices of two polytopes.
OUTPUT: list of vertices of the Minkowski sum, given as vectors.
EXAMPLES: Let's compute the Minkowski sum of two line segments::
sage: lattice_polytope.minkowski_sum([[1,0],[-1,0]],[[0,1],[0,-1]]) [(1, 1), (1, -1), (-1, 1), (-1, -1)] """
def positive_integer_relations(points): r""" Return relations between given points.
INPUT:
- ``points`` - lattice points given as columns of a matrix
OUTPUT: matrix of relations between given points with non-negative integer coefficients
EXAMPLES: This is a 3-dimensional reflexive polytope::
sage: p = LatticePolytope([(1,0,0), (0,1,0), ....: (-1,-1,0), (0,0,1), (-1,0,-1)]) sage: p.points() M( 1, 0, 0), M( 0, 1, 0), M(-1, -1, 0), M( 0, 0, 1), M(-1, 0, -1), M( 0, 0, 0) in 3-d lattice M
We can compute linear relations between its points in the following way::
sage: p.points().matrix().kernel().echelonized_basis_matrix() [ 1 0 0 1 1 0] [ 0 1 1 -1 -1 0] [ 0 0 0 0 0 1]
However, the above relations may contain negative and rational numbers. This function transforms them in such a way, that all coefficients are non-negative integers::
sage: lattice_polytope.positive_integer_relations(p.points().column_matrix()) [1 0 0 1 1 0] [1 1 1 0 0 0] [0 0 0 0 0 1]
sage: cm = ReflexivePolytope(2,1).vertices().column_matrix() sage: lattice_polytope.positive_integer_relations(cm) [2 1 1] """ # Find a non-negative linear combination of relations, # such that all components are non-negative and the i-th one is 1
# Use the new relation to remove negative entries in non-pivot columns # Get a new basis # Switch to integers
def read_all_polytopes(file_name): r""" Read all polytopes from the given file.
INPUT:
- ``file_name`` -- a string with the name of a file with VERTICES of polytopes.
OUTPUT:
- a sequence of polytopes.
EXAMPLES:
We use poly.x to compute two polar polytopes and read them::
sage: d = lattice_polytope.cross_polytope(2) sage: o = lattice_polytope.cross_polytope(3) sage: result_name = lattice_polytope._palp("poly.x -fe", [d, o]) sage: with open(result_name) as f: ....: print(f.read()) 4 2 Vertices of P-dual <-> Equations of P -1 1 1 1 -1 -1 1 -1 8 3 Vertices of P-dual <-> Equations of P -1 -1 1 1 -1 1 -1 1 1 1 1 1 -1 -1 -1 1 -1 -1 -1 1 -1 1 1 -1 sage: lattice_polytope.read_all_polytopes(result_name) [2-d reflexive polytope #14 in 2-d lattice M, 3-d reflexive polytope in 3-d lattice M] sage: os.remove(result_name) """
def read_palp_matrix(data, permutation=False): r""" Read and return an integer matrix from a string or an opened file.
First input line must start with two integers m and n, the number of rows and columns of the matrix. The rest of the first line is ignored. The next m lines must contain n numbers each.
If m>n, returns the transposed matrix. If the string is empty or EOF is reached, returns the empty matrix, constructed by ``matrix()``.
INPUT:
- ``data`` -- Either a string containing the filename or the file itself containing the output by PALP.
- ``permutation`` -- (default: ``False``) If ``True``, try to retrieve the permutation output by PALP. This parameter makes sense only when PALP computed the normal form of a lattice polytope.
OUTPUT:
A matrix or a tuple of a matrix and a permutation.
EXAMPLES::
sage: lattice_polytope.read_palp_matrix("2 3 comment \n 1 2 3 \n 4 5 6") [1 2 3] [4 5 6] sage: lattice_polytope.read_palp_matrix("3 2 Will be transposed \n 1 2 \n 3 4 \n 5 6") [1 3 5] [2 4 6] """ # If data is not a string, try to treat it as a file. return matrix() # In some cases there may be additional information to extract last_piece = first_line[-1] last_piece = last_piece.split('=') if last_piece[0] != 'perm': raise ValueError('PALP did not return a permutation.') p = _palp_convert_permutation(last_piece[1]) return (mat, p) else:
def set_palp_dimension(d): r""" Set the dimension for PALP calls to ``d``.
INPUT:
- ``d`` -- an integer from the list [4,5,6,11] or ``None``.
OUTPUT:
- none.
PALP has many hard-coded limits, which must be specified before compilation, one of them is dimension. Sage includes several versions with different dimension settings (which may also affect other limits and enable certain features of PALP). You can change the version which will be used by calling this function. Such a change is not done automatically for each polytope based on its dimension, since depending on what you are doing it may be necessary to use dimensions higher than that of the input polytope.
EXAMPLES:
Let's try to work with a 7-dimensional polytope::
sage: p = lattice_polytope.cross_polytope(7) sage: p._palp("poly.x -fv") Traceback (most recent call last): ... ValueError: Error executing 'poly.x -fv' for the given polytope! Output: Please increase POLY_Dmax to at least 7
However, we can work with this polytope by changing PALP dimension to 11::
sage: lattice_polytope.set_palp_dimension(11) sage: p._palp("poly.x -fv") '7 14 Vertices of P...'
Let's go back to default settings::
sage: lattice_polytope.set_palp_dimension(None) """ global _palp_dimension
def skip_palp_matrix(data, n=1): r""" Skip matrix data in a file.
INPUT:
- ``data`` - opened file with blocks of matrix data in the following format: A block consisting of m+1 lines has the number m as the first element of its first line.
- ``n`` - (default: 1) integer, specifies how many blocks should be skipped
If EOF is reached during the process, raises ValueError exception.
EXAMPLES: We create a file with vertices of the square and the cube, but read only the second set::
sage: d = lattice_polytope.cross_polytope(2) sage: o = lattice_polytope.cross_polytope(3) sage: result_name = lattice_polytope._palp("poly.x -fe", [d, o]) sage: with open(result_name) as f: ....: print(f.read()) 4 2 Vertices of P-dual <-> Equations of P -1 1 1 1 -1 -1 1 -1 8 3 Vertices of P-dual <-> Equations of P -1 -1 1 1 -1 1 -1 1 1 1 1 1 -1 -1 -1 1 -1 -1 -1 1 -1 1 1 -1 sage: f = open(result_name) sage: lattice_polytope.skip_palp_matrix(f) sage: lattice_polytope.read_palp_matrix(f) [-1 1 -1 1 -1 1 -1 1] [-1 -1 1 1 -1 -1 1 1] [ 1 1 1 1 -1 -1 -1 -1] sage: f.close() sage: os.remove(result_name) """ raise ValueError("There are not enough data to skip!") raise ValueError("There are not enough data to skip!")
def write_palp_matrix(m, ofile=None, comment="", format=None): r""" Write ``m`` into ``ofile`` in PALP format.
INPUT:
- ``m`` -- a matrix over integers or a :class:`point collection <PointCollection>`.
- ``ofile`` -- a file opened for writing (default: stdout)
- ``comment`` -- a string (default: empty) see output description
- ``format`` -- a format string used to print matrix entries.
OUTPUT:
- nothing is returned, output written to ``ofile`` has the format
* First line: number_of_rows number_of_columns comment * Next number_of_rows lines: rows of the matrix.
EXAMPLES::
sage: o = lattice_polytope.cross_polytope(3) sage: lattice_polytope.write_palp_matrix(o.vertices(), comment="3D Octahedron") 3 6 3D Octahedron 1 0 0 -1 0 0 0 1 0 0 -1 0 0 0 1 0 0 -1 sage: lattice_polytope.write_palp_matrix(o.vertices(), format="%4d") 3 6 1 0 0 -1 0 0 0 1 0 0 -1 0 0 0 1 0 0 -1 """ for i in range(m.nrows()) for j in range(m.ncols())) else: else: |